MappingTest.scala 7.62 KB
Newer Older
bow's avatar
bow committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
/**
 * Biopet is built on top of GATK Queue for building bioinformatic
 * pipelines. It is mainly intended to support LUMC SHARK cluster which is running
 * SGE. But other types of HPC that are supported by GATK Queue (such as PBS)
 * should also be able to execute Biopet tools and pipelines.
 *
 * Copyright 2014 Sequencing Analysis Support Core - Leiden University Medical Center
 *
 * Contact us at: sasc@lumc.nl
 *
 * A dual licensing mode is applied. The source code within this project that are
 * not part of GATK Queue is freely available for non-commercial use under an AGPL
 * license; For commercial users or users who do not want to follow the AGPL
 * license, please contact us to obtain a separate license.
 */
16
17
package nl.lumc.sasc.biopet.pipelines.mapping

Peter van 't Hof's avatar
Peter van 't Hof committed
18
import java.io.{ File, FileOutputStream }
19

20
import com.google.common.io.Files
Peter van 't Hof's avatar
Peter van 't Hof committed
21
import nl.lumc.sasc.biopet.utils.config.Config
Peter van 't Hof's avatar
Peter van 't Hof committed
22
import nl.lumc.sasc.biopet.extensions._
Peter van 't Hof's avatar
Peter van 't Hof committed
23
24
import nl.lumc.sasc.biopet.extensions.bwa.{ BwaAln, BwaMem, BwaSampe, BwaSamse }
import nl.lumc.sasc.biopet.extensions.picard.{ AddOrReplaceReadGroups, MarkDuplicates, MergeSamFiles, SortSam }
Peter van 't Hof's avatar
Peter van 't Hof committed
25
import nl.lumc.sasc.biopet.pipelines.flexiprep.Fastqc
26
import nl.lumc.sasc.biopet.extensions.tools.{ FastqSync, SeqStat }
Peter van 't Hof's avatar
Peter van 't Hof committed
27
import nl.lumc.sasc.biopet.utils.ConfigUtils
28
29
30
31
32
33
import org.apache.commons.io.FileUtils
import org.broadinstitute.gatk.queue.QSettings
import org.scalatest.Matchers
import org.scalatest.testng.TestNGSuite
import org.testng.annotations.{ AfterClass, DataProvider, Test }

34
/**
Peter van 't Hof's avatar
Peter van 't Hof committed
35
36
 * Test class for [[Mapping]]
 *
37
38
39
40
41
42
43
44
 * Created by pjvan_thof on 2/12/15.
 */
class MappingTest extends TestNGSuite with Matchers {
  def initPipeline(map: Map[String, Any]): Mapping = {
    new Mapping {
      override def configName = "mapping"
      override def globalConfig = new Config(map)
      qSettings = new QSettings
Peter van 't Hof's avatar
Peter van 't Hof committed
45
      qSettings.runName = "test"
46
47
48
    }
  }

Peter van 't Hof's avatar
Peter van 't Hof committed
49
  @DataProvider(name = "mappingOptions")
50
  def mappingOptions = {
Peter van 't Hof's avatar
Peter van 't Hof committed
51
    val aligners = Array("bwa-mem", "bwa-aln", "star", "star-2pass", "bowtie", "stampy", "gsnap", "tophat")
52
53
54
55
    val paired = Array(true, false)
    val chunks = Array(1, 5, 10, 100)
    val skipMarkDuplicates = Array(true, false)
    val skipFlexipreps = Array(true, false)
56
    val zipped = Array(true, false)
57
58
59
60
61
62

    for (
      aligner <- aligners;
      pair <- paired;
      chunk <- chunks;
      skipMarkDuplicate <- skipMarkDuplicates;
63
64
65
      skipFlexiprep <- skipFlexipreps;
      zipped <- zipped
    ) yield Array(aligner, pair, chunk, skipMarkDuplicate, skipFlexiprep, zipped)
66
67
68
  }

  @Test(dataProvider = "mappingOptions")
69
70
71
72
  def testMapping(aligner: String, paired: Boolean, chunks: Int,
                  skipMarkDuplicate: Boolean,
                  skipFlexiprep: Boolean,
                  zipped: Boolean) = {
73
74
75
76
77
    val map = ConfigUtils.mergeMaps(Map("output_dir" -> MappingTest.outputDir,
      "aligner" -> aligner,
      "number_chunks" -> chunks,
      "skip_markduplicates" -> skipMarkDuplicate,
      "skip_flexiprep" -> skipFlexiprep
78
    ), Map(MappingTest.executables.toSeq: _*))
79
80
    val mapping: Mapping = initPipeline(map)

81
    if (zipped) {
Peter van 't Hof's avatar
Peter van 't Hof committed
82
83
      mapping.input_R1 = MappingTest.r1Zipped
      if (paired) mapping.input_R2 = Some(MappingTest.r2Zipped)
84
    } else {
Peter van 't Hof's avatar
Peter van 't Hof committed
85
86
      mapping.input_R1 = MappingTest.r1
      if (paired) mapping.input_R2 = Some(MappingTest.r2)
87
    }
88
89
    mapping.sampleId = Some("1")
    mapping.libId = Some("1")
90
91
92
93
    mapping.script()

    //Flexiprep
    mapping.functions.count(_.isInstanceOf[Fastqc]) shouldBe (if (skipFlexiprep) 0 else if (paired) 4 else 2)
Peter van 't Hof's avatar
Peter van 't Hof committed
94
95
96
97
98
99
100
    //mapping.functions.count(_.isInstanceOf[Zcat]) shouldBe (if (!zipped || (chunks > 1 && skipFlexiprep)) 0 else if (paired) 2 else 1)
    //mapping.functions.count(_.isInstanceOf[SeqStat]) shouldBe ((if (skipFlexiprep) 0 else if (paired) 2 else 1) * (chunks + 1))
    //mapping.functions.count(_.isInstanceOf[SeqtkSeq]) shouldBe ((if (skipFlexiprep) 0 else if (paired) 2 else 1) * chunks)
    //mapping.functions.count(_.isInstanceOf[Cutadapt]) shouldBe ((if (skipFlexiprep) 0 else if (paired) 2 else 1) * chunks)
    //mapping.functions.count(_.isInstanceOf[FastqSync]) shouldBe ((if (skipFlexiprep) 0 else if (paired && !skipFlexiprep) 1 else 0) * chunks)
    //mapping.functions.count(_.isInstanceOf[Sickle]) shouldBe ((if (skipFlexiprep) 0 else 1) * chunks)
    //mapping.functions.count(_.isInstanceOf[Gzip]) shouldBe (if (skipFlexiprep) 0 else if (paired) 2 else 1)
101
102

    //aligners
Peter van 't Hof's avatar
Peter van 't Hof committed
103
    //mapping.functions.count(_.isInstanceOf[BwaMem]) shouldBe ((if (aligner == "bwa-mem") 1 else 0) * chunks)
Peter van 't Hof's avatar
Peter van 't Hof committed
104
105
106
    mapping.functions.count(_.isInstanceOf[BwaAln]) shouldBe ((if (aligner == "bwa-aln") if (paired) 2 else 1 else 0) * chunks)
    mapping.functions.count(_.isInstanceOf[BwaSampe]) shouldBe ((if (aligner == "bwa-aln") if (paired) 1 else 0 else 0) * chunks)
    mapping.functions.count(_.isInstanceOf[BwaSamse]) shouldBe ((if (aligner == "bwa-aln") if (paired) 0 else 1 else 0) * chunks)
107
108
109
110
111
112
    mapping.functions.count(_.isInstanceOf[Star]) shouldBe ((if (aligner == "star") 1 else if (aligner == "star-2pass") 3 else 0) * chunks)
    mapping.functions.count(_.isInstanceOf[Bowtie]) shouldBe ((if (aligner == "bowtie") 1 else 0) * chunks)
    mapping.functions.count(_.isInstanceOf[Stampy]) shouldBe ((if (aligner == "stampy") 1 else 0) * chunks)

    // Sort sam or replace readgroup
    val sort = aligner match {
Peter van 't Hof's avatar
Peter van 't Hof committed
113
      case "bwa-mem" | "bwa-aln" | "stampy" => "sortsam"
114
      case "star" | "star-2pass" | "bowtie" | "gsnap" | "tophat" => "replacereadgroups"
Peter van 't Hof's avatar
Peter van 't Hof committed
115
      case _ => throw new IllegalArgumentException("aligner: " + aligner + " does not exist")
116
117
    }

bow's avatar
bow committed
118
    if (aligner != "tophat") { // FIXME
Peter van 't Hof's avatar
Peter van 't Hof committed
119
      //mapping.functions.count(_.isInstanceOf[SortSam]) shouldBe ((if (sort == "sortsam") 1 else 0) * chunks)
bow's avatar
bow committed
120
121
122
123
      mapping.functions.count(_.isInstanceOf[AddOrReplaceReadGroups]) shouldBe ((if (sort == "replacereadgroups") 1 else 0) * chunks)
      mapping.functions.count(_.isInstanceOf[MergeSamFiles]) shouldBe (if (skipMarkDuplicate && chunks > 1) 1 else 0)
      mapping.functions.count(_.isInstanceOf[MarkDuplicates]) shouldBe (if (skipMarkDuplicate) 0 else 1)
    }
124
  }
125
126
127
128
129

  // remove temporary run directory all tests in the class have been run
  @AfterClass def removeTempOutputDir() = {
    FileUtils.deleteDirectory(MappingTest.outputDir)
  }
130
131
132
}

object MappingTest {
133
  val outputDir = Files.createTempDir()
Peter van 't Hof's avatar
Peter van 't Hof committed
134
135
136
137
138
139
140
141
142
143
  new File(outputDir, "input").mkdirs()

  val r1 = new File(outputDir, "input" + File.separator + "R1.fq")
  Files.touch(r1)
  val r2 = new File(outputDir, "input" + File.separator + "R2.fq")
  Files.touch(r2)
  val r1Zipped = new File(outputDir, "input" + File.separator + "R1.fq.gz")
  Files.touch(r1Zipped)
  val r2Zipped = new File(outputDir, "input" + File.separator + "R2.fq.gz")
  Files.touch(r2Zipped)
144

145
146
147
148
149
150
151
152
153
154
  private def copyFile(name: String): Unit = {
    val is = getClass.getResourceAsStream("/" + name)
    val os = new FileOutputStream(new File(outputDir, name))
    org.apache.commons.io.IOUtils.copy(is, os)
    os.close()
  }

  copyFile("ref.fa")
  copyFile("ref.dict")
  copyFile("ref.fa.fai")
Peter van 't Hof's avatar
Peter van 't Hof committed
155
156
  copyFile("ref.1.bt2")
  copyFile("ref.1.ebwt")
157

158
  val executables = Map(
Peter van 't Hof's avatar
Peter van 't Hof committed
159
    "reference_fasta" -> (outputDir + File.separator + "ref.fa"),
160
    "db" -> "test",
Peter van 't Hof's avatar
Peter van 't Hof committed
161
    "bowtie_index" -> (outputDir + File.separator + "ref"),
162
163
    "fastqc" -> Map("exe" -> "test"),
    "seqtk" -> Map("exe" -> "test"),
164
165
    "gsnap" -> Map("exe" -> "test"),
    "tophat" -> Map("exe" -> "test"),
166
    "sickle" -> Map("exe" -> "test"),
167
    "cutadapt" -> Map("exe" -> "test"),
168
169
170
    "bwa" -> Map("exe" -> "test"),
    "star" -> Map("exe" -> "test"),
    "bowtie" -> Map("exe" -> "test"),
171
    "stampy" -> Map("exe" -> "test", "genome" -> "test", "hash" -> "test"),
172
173
    "samtools" -> Map("exe" -> "test"),
    "md5sum" -> Map("exe" -> "test")
174
175
  )
}