Skip to content

B
biopet.biopet
Commit 8db9bc07 authored by pjvan_thof's avatar pjvan_thof
Browse files
Options

Adding unit tests to MuTect2

parent d22fc2f7
Hide whitespace changes
Inline Side-by-side
Showing with 80 additions and 42 deletions
shiva/src/main/scala/nl/lumc/sasc/biopet/pipelines/shiva/variantcallers/somatic/MuTect2.scala
View file @ 8db9bc07
......@@ -3,7 +3,7 @@ package nl.lumc.sasc.biopet.pipelines.shiva.variantcallers.somatic
import nl.lumc.sasc.biopet.extensions.bcftools.BcftoolsReheader
import nl.lumc.sasc.biopet.extensions.gatk.{BqsrGather, CombineVariants}
import nl.lumc.sasc.biopet.extensions._
import nl.lumc.sasc.biopet.utils.IoUtils
import nl.lumc.sasc.biopet.utils.{IoUtils, Logging}
import nl.lumc.sasc.biopet.utils.config.Configurable
import org.broadinstitute.gatk.queue.extensions.gatk.TaggedFile
......@@ -19,6 +19,7 @@ class MuTect2(val parent: Configurable) extends SomaticVariantCaller {
def biopetScript(): Unit = {
if (tnPairs.isEmpty) Logging.addError("No tumor-normal found in conig")
val outputFiles = for (pair <- tnPairs) yield {
val bqsrFile =
......@@ -38,23 +39,22 @@ class MuTect2(val parent: Configurable) extends SomaticVariantCaller {
new File(outputDir, s"${pair.tumorSample}-${pair.normalSample}.$name.vcf.gz")
val muTect2 = new gatk.MuTect2(this)
muTect2.input_file :+= TaggedFile(inputBams(pair.tumorSample), "tumor")
muTect2.input_file :+= TaggedFile(inputBams(pair.normalSample), "normal")
inputBams.get(pair.tumorSample).foreach(muTect2.input_file :+= TaggedFile(_, "tumor"))
inputBams.get(pair.normalSample).foreach(muTect2.input_file :+= TaggedFile(_, "normal"))
muTect2.BQSR = bqsrFile
if (runConEst) {
val namePrefix = outputFile.getAbsolutePath.stripSuffix(".vcf.gz")
val contEstOutput: File = new File(s"$namePrefix.contamination.txt")
val contEst = gatk.ContEst(this,
inputBams(pair.tumorSample),
inputBams(pair.normalSample),
contEstOutput)
val contEst = new gatk.ContEst(this)
inputBams.get(pair.tumorSample).foreach(contEst.input_file :+= _)
inputBams.get(pair.normalSample).foreach(contEst.input_file :+= _)
contEst.output = new File(s"$namePrefix.contamination.txt")
contEst.BQSR = bqsrFile
add(contEst)
val contaminationPerSample: File = new File(s"$namePrefix.contamination.short.txt")
val awk: Awk = Awk(this, "BEGIN{OFS=\"\\t\"}{if($1 != \"name\") print $1,$4;}")
awk.input = contEstOutput
awk.input = contEst.output
add(awk > contaminationPerSample)
muTect2.contaminationFile = Some(contaminationPerSample)
......@@ -71,6 +71,7 @@ class MuTect2(val parent: Configurable) extends SomaticVariantCaller {
muTect2.BQSR = Some(gather.outputBqsrFile)
}
outputDir.mkdirs()
val renameFile = new File(outputDir, s".rename.${pair.tumorSample}-${pair.normalSample}.txt")
IoUtils.writeLinesToFile(renameFile,
List(
......
shiva/src/test/scala/nl/lumc/sasc/biopet/pipelines/shiva/ShivaSvCallingTest.scala
View file @ 8db9bc07
......@@ -55,7 +55,7 @@ class ShivaSvCallingTest extends TestNGSuite with Matchers {
private var dirs: List[File] = Nil
@DataProvider(name = "shivaSvCallingOptions")
def shivaSvCallingOptions = {
def shivaSvCallingOptions: Array[Array[AnyVal]] = {
val bool = Array(true, false)
(for (bams <- 0 to 3;
delly <- bool;
......@@ -69,7 +69,7 @@ class ShivaSvCallingTest extends TestNGSuite with Matchers {
delly: Boolean,
clever: Boolean,
breakdancer: Boolean,
pindel: Boolean) = {
pindel: Boolean): Unit = {
val outputDir = ShivaSvCallingTest.outputDir
dirs :+= outputDir
val callers: ListBuffer[String] = ListBuffer()
......@@ -95,7 +95,7 @@ class ShivaSvCallingTest extends TestNGSuite with Matchers {
pipeline.script()
val summaryCallers =
pipeline.summarySettings.get("sv_callers").get.asInstanceOf[List[String]]
pipeline.summarySettings("sv_callers").asInstanceOf[List[String]]
if (delly) assert(summaryCallers.contains("delly"))
else assert(!summaryCallers.contains("delly"))
if (clever) assert(summaryCallers.contains("clever"))
......@@ -123,7 +123,7 @@ class ShivaSvCallingTest extends TestNGSuite with Matchers {
}
@DataProvider(name = "dellyOptions")
def dellyOptions = {
def dellyOptions: Array[Array[AnyVal]] = {
val bool = Array(true, false)
for (del <- bool;
dup <- bool;
......@@ -188,21 +188,21 @@ class ShivaSvCallingTest extends TestNGSuite with Matchers {
pipeline.script()
val summaryCallers: List[String] =
pipeline.summarySettings.get("sv_callers").get.asInstanceOf[List[String]]
pipeline.summarySettings("sv_callers").asInstanceOf[List[String]]
assert(summaryCallers.contains("delly"))
assert(summaryCallers.contains("clever"))
assert(summaryCallers.contains("breakdancer"))
}
// remove temporary run directory all tests in the class have been run
@AfterClass def removeTempOutputDir() = {
@AfterClass def removeTempOutputDir(): Unit = {
dirs.foreach(FileUtils.deleteDirectory)
}
}
object ShivaSvCallingTest {
def outputDir = Files.createTempDir()
val inputDir = Files.createTempDir()
def outputDir: File = Files.createTempDir()
val inputDir: File = Files.createTempDir()
private def inputTouch(name: String): File = {
val file = new File(outputDir, name).getAbsoluteFile
......
shiva/src/test/scala/nl/lumc/sasc/biopet/pipelines/shiva/ShivaTest.scala
View file @ 8db9bc07
......@@ -49,7 +49,7 @@ trait ShivaTestTrait extends TestNGSuite with Matchers {
}
@DataProvider(name = "shivaOptions")
def shivaOptions = {
def shivaOptions: Array[Array[Any]] = {
for (s1 <- sample1; s2 <- sample2;
realign <- realignProvider; baseRecalibration <- baseRecalibrationProvider)
yield Array("", s1, s2, realign, baseRecalibration)
......@@ -121,7 +121,7 @@ trait ShivaTestTrait extends TestNGSuite with Matchers {
numberSamples
else 0)
pipeline.functions.count(_.isInstanceOf[BaseRecalibrator]) shouldBe (if (dbsnp && baseRecalibration)
(numberLibs * 2)
numberLibs * 2
else 0)
pipeline.functions.count(_.isInstanceOf[PrintReads]) shouldBe (if (dbsnp && baseRecalibration && usePrintReads)
numberLibs
......@@ -212,9 +212,9 @@ class ShivaWithAnnotationTest extends ShivaTestTrait {
}
object ShivaTest {
def outputDir = Files.createTempDir()
def outputDir: File = Files.createTempDir()
val inputDir = Files.createTempDir()
val inputDir: File = Files.createTempDir()
def inputTouch(name: String): String = {
val file = new File(inputDir, name)
......
shiva/src/test/scala/nl/lumc/sasc/biopet/pipelines/shiva/ShivaVariantcallingTest.scala
View file @ 8db9bc07
......@@ -25,16 +25,12 @@ import com.google.common.io.Files
import nl.lumc.sasc.biopet.core.BiopetPipe
import nl.lumc.sasc.biopet.extensions.Freebayes
import nl.lumc.sasc.biopet.extensions.bcftools.{BcftoolsCall, BcftoolsMerge}
import nl.lumc.sasc.biopet.extensions.gatk.{
CombineVariants,
GenotypeConcordance,
HaplotypeCaller,
UnifiedGenotyper
}
import nl.lumc.sasc.biopet.extensions.gatk._
import nl.lumc.sasc.biopet.utils.config.Config
import nl.lumc.sasc.biopet.extensions.tools.{MpileupToVcf, VcfFilter, VcfStats}
import nl.lumc.sasc.biopet.extensions.vt.{VtDecompose, VtNormalize}
import nl.lumc.sasc.biopet.utils.ConfigUtils
import nl.lumc.sasc.biopet.pipelines.shiva.variantcallers.somatic.TumorNormalPair
import nl.lumc.sasc.biopet.utils.{ConfigUtils, Logging}
import org.apache.commons.io.FileUtils
import org.broadinstitute.gatk.queue.QSettings
import org.scalatest.Matchers
......@@ -60,6 +56,7 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
}
def raw: Boolean = false
def mutect2: Boolean = false
def bcftools: Boolean = false
def bcftools_singlesample: Boolean = false
def haplotypeCallerGvcf: Boolean = false
......@@ -72,17 +69,23 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
def referenceVcf: Option[File] = None
def roiBedFiles: List[File] = Nil
def ampliconBedFile: Option[File] = None
def runContest: Option[Boolean] = None
def normalize = false
def decompose = false
def bamRange: List[Int] = (0 to 2).toList
def tumorNormalPairs: List[TumorNormalPair] = Nil
@DataProvider(name = "shivaVariantcallingOptions")
def shivaVariantcallingOptions = {
(for (bams <- 0 to 2)
def shivaVariantcallingOptions: Array[Array[Any]] = {
(for (bams <- bamRange)
yield
Array[Any](
bams,
raw,
mutect2,
bcftools,
bcftools_singlesample,
unifiedGenotyper,
......@@ -91,7 +94,8 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
haplotypeCallerAllele,
unifiedGenotyperAllele,
freebayes,
varscanCnsSinglesample
varscanCnsSinglesample,
tumorNormalPairs
)).toArray
}
......@@ -100,6 +104,7 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
@Test(dataProvider = "shivaVariantcallingOptions")
def testShivaVariantcalling(bams: Int,
raw: Boolean,
mutect2: Boolean,
bcftools: Boolean,
bcftoolsSinglesample: Boolean,
unifiedGenotyper: Boolean,
......@@ -108,10 +113,13 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
haplotypeCallerAllele: Boolean,
unifiedGenotyperAllele: Boolean,
freebayes: Boolean,
varscanCnsSinglesample: Boolean) = {
varscanCnsSinglesample: Boolean,
tumorNormalPairs: List[TumorNormalPair]): Unit = {
val outputDir = ShivaVariantcallingTest.outputDir
dirs :+= outputDir
val callers: ListBuffer[String] = ListBuffer()
if (raw) callers.append("raw")
if (mutect2) callers.append("mutect2")
if (bcftools) callers.append("bcftools")
if (bcftoolsSinglesample) callers.append("bcftools_singlesample")
if (unifiedGenotyper) callers.append("unifiedgenotyper")
......@@ -121,25 +129,34 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
if (haplotypeCaller) callers.append("haplotypecaller")
if (freebayes) callers.append("freebayes")
if (varscanCnsSinglesample) callers.append("varscan_cns_singlesample")
val map = Map(
val sampleTags: Map[String, Any] = tumorNormalPairs.foldLeft(Map[String, Any]()) {
case (m, pair) =>
val tag = Map(
"samples" -> Map(pair.tumorSample -> Map(
"tags" -> Map("type" -> "tumor", "normal" -> pair.normalSample))))
ConfigUtils.mergeMaps(m, tag)
}
val map = sampleTags ++ Map(
"variantcallers" -> callers.toList,
"execute_vt_normalize" -> normalize,
"execute_vt_decompose" -> decompose,
"regions_of_interest" -> roiBedFiles.map(_.getAbsolutePath)
) ++ referenceVcf.map("reference_vcf" -> _) ++ ampliconBedFile.map(
"amplicon_bed" -> _.getAbsolutePath)
"amplicon_bed" -> _.getAbsolutePath) ++ runContest.map("run_contest" -> _)
val pipeline = initPipeline(map, outputDir)
pipeline.inputBams = (for (n <- 1 to bams)
yield n.toString -> ShivaVariantcallingTest.inputTouch("bam_" + n + ".bam")).toMap
yield
s"sample_${n.toString}" -> ShivaVariantcallingTest.inputTouch("bam_" + n + ".bam")).toMap
val illegalArgumentException = pipeline.inputBams.isEmpty || callers.isEmpty
val illegalStateException = mutect2 && bams == 1
if (illegalArgumentException) intercept[IllegalArgumentException] {
pipeline.script()
}
if (!illegalArgumentException) {
} else if (illegalStateException) intercept[IllegalStateException] {
pipeline.script()
} else {
pipeline.script()
val pipesJobs = pipeline.functions
......@@ -173,6 +190,11 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
pipeline.functions.count(_.isInstanceOf[GenotypeConcordance]) shouldBe (if (referenceVcf.isDefined)
1 + callers.size
else 0)
pipesJobs.count(_.isInstanceOf[MuTect2]) shouldBe (if (mutect2) tumorNormalPairs.size else 0)
pipeline.functions.count(_.isInstanceOf[ContEst]) shouldBe (if (mutect2 && runContest
.getOrElse(false))
tumorNormalPairs.size
else 0)
pipeline.summarySettings
.get("variantcallers")
......@@ -182,10 +204,11 @@ trait ShivaVariantcallingTestTrait extends TestNGSuite with Matchers {
pipeline.summarySettings.get("regions_of_interest") shouldBe Some(
roiBedFiles.map(_.getAbsolutePath))
}
Logging.errors.clear()
}
// remove temporary run directory all tests in the class have been run
@AfterClass def removeTempOutputDir() = {
@AfterClass def removeTempOutputDir(): Unit = {
dirs.foreach(FileUtils.deleteDirectory)
}
}
......@@ -206,6 +229,19 @@ class ShivaVariantcallingAllTest extends ShivaVariantcallingTestTrait {
class ShivaVariantcallingRawTest extends ShivaVariantcallingTestTrait {
override def raw: Boolean = true
}
class ShivaVariantcallingMuTect2Test extends ShivaVariantcallingTestTrait {
override def mutect2: Boolean = true
override def haplotypeCaller: Boolean = true
override def tumorNormalPairs: List[TumorNormalPair] =
TumorNormalPair("sample_1", "sample_2") :: Nil
}
class ShivaVariantcallingMuTect2ContestTest extends ShivaVariantcallingTestTrait {
override def mutect2: Boolean = true
override def haplotypeCaller: Boolean = true
override def runContest = Some(true)
override def tumorNormalPairs: List[TumorNormalPair] =
TumorNormalPair("sample_1", "sample_2") :: Nil
}
class ShivaVariantcallingBcftoolsTest extends ShivaVariantcallingTestTrait {
override def bcftools: Boolean = true
}
......@@ -256,8 +292,8 @@ class ShivaVariantcallingAmpliconTest extends ShivaVariantcallingTestTrait {
}
object ShivaVariantcallingTest {
def outputDir = Files.createTempDir()
val inputDir = Files.createTempDir()
def outputDir: File = Files.createTempDir()
val inputDir: File = Files.createTempDir()
def inputTouch(name: String): File = {
val file = new File(inputDir, name).getAbsoluteFile
......@@ -292,6 +328,7 @@ object ShivaVariantcallingTest {
"tabix" -> Map("exe" -> "test"),
"input_alleles" -> "test.vcf.gz",
"varscan_jar" -> "test",
"vt" -> Map("exe" -> "test")
"vt" -> Map("exe" -> "test"),
"popfile" -> "test"
)
}
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment