Added allele modes

protected/biopet-gatk-pipelines/src/main/scala/nl/lumc/sasc/biopet/pipelines/gatk/ShivaVariantcallingGatk.scala

@@ -12,7 +12,13 @@ class ShivaVariantcallingGatk(val root: Configurable) extends QScript with Shiva
   qscript =>
   def this() = this(null)
 
-  override def callers = new UnifiedGenotyper :: new HaplotypeCaller :: super.callers
+  override def callers = {
+    new HaplotypeCallerAllele ::
+      new UnifiedGenotyperAllele ::
+      new UnifiedGenotyper ::
+      new HaplotypeCaller ::
+      super.callers
+  }
 
   class HaplotypeCaller extends Variantcaller {
     val name = "haplotypecaller"
@@ -37,12 +43,46 @@ class ShivaVariantcallingGatk(val root: Configurable) extends QScript with Shiva
     def outputFile = new File(outputDir, namePrefix + "unifiedgenotyper.vcf.gz")
 
     def addJobs() {
-      val hc = new nl.lumc.sasc.biopet.extensions.gatk.broad.UnifiedGenotyper(qscript)
+      val ug = new nl.lumc.sasc.biopet.extensions.gatk.broad.UnifiedGenotyper(qscript)
+      ug.input_file = inputBams
+      ug.out = outputFile
+      add(ug)
+    }
+  }
+
+  class HaplotypeCallerAllele extends Variantcaller {
+    val name = "haplotypecaller_allele"
+    protected val defaultPrio = 5
+    protected val defaultUse = false
+
+    def outputFile = new File(outputDir, namePrefix + "haplotypecaller_allele.vcf.gz")
+
+    def addJobs() {
+      val hc = new nl.lumc.sasc.biopet.extensions.gatk.broad.HaplotypeCaller(qscript)
       hc.input_file = inputBams
       hc.out = outputFile
+      hc.alleles = config("input_alleles")
+      hc.genotyping_mode = org.broadinstitute.gatk.tools.walkers.genotyper.GenotypingOutputMode.GENOTYPE_GIVEN_ALLELES
       add(hc)
     }
   }
+
+  class UnifiedGenotyperAllele extends Variantcaller {
+    val name = "unifiedgenotyper_allele"
+    protected val defaultPrio = 6
+    protected val defaultUse = false
+
+    def outputFile = new File(outputDir, namePrefix + "unifiedgenotyper_allele.vcf.gz")
+
+    def addJobs() {
+      val ug = new nl.lumc.sasc.biopet.extensions.gatk.broad.UnifiedGenotyper(qscript)
+      ug.input_file = inputBams
+      ug.out = outputFile
+      ug.alleles = config("input_alleles")
+      ug.genotyping_mode = org.broadinstitute.gatk.tools.walkers.genotyper.GenotypingOutputMode.GENOTYPE_GIVEN_ALLELES
+      add(ug)
+    }
+  }
 }
 
 object ShivaVariantcallingGatk extends PipelineCommand
\ No newline at end of file

public/biopet-framework/src/main/scala/nl/lumc/sasc/biopet/extensions/gatk/CombineVariants.scala

@@ -18,7 +18,7 @@ class CombineVariants(val root: Configurable) extends Gatk {
   var outputFile: File = null
 
   var setKey: String = null
-  var rodPriorityList: List[String] = Nil
+  var rodPriorityList: String = null
   var minimumN: Int = config("minimumN", default = 1)
   var genotypeMergeOptions: Option[String] = config("genotypeMergeOptions")
 
@@ -45,6 +45,6 @@ class CombineVariants(val root: Configurable) extends Gatk {
     }).mkString +
     required("-o", outputFile) +
     optional("--setKey", setKey) +
-    (if (rodPriorityList.isEmpty) "" else optional("--rod_priority_list", rodPriorityList.mkString(","))) +
+    optional("--rod_priority_list", rodPriorityList) +
     optional("-genotypeMergeOptions", genotypeMergeOptions)
 }

public/shiva/src/main/scala/nl/lumc/sasc/biopet/pipelines/shiva/ShivaTrait.scala

@@ -7,6 +7,7 @@ import nl.lumc.sasc.biopet.core.summary.SummaryQScript
 import nl.lumc.sasc.biopet.core.MultiSampleQScript
 import nl.lumc.sasc.biopet.extensions.Ln
 import nl.lumc.sasc.biopet.extensions.picard.{ AddOrReplaceReadGroups, SamToFastq, MarkDuplicates }
+import nl.lumc.sasc.biopet.pipelines.bammetrics.BamMetrics
 import nl.lumc.sasc.biopet.pipelines.mapping.Mapping
 import scala.collection.JavaConversions._
 
@@ -182,15 +183,27 @@ trait ShivaTrait extends MultiSampleQScript with SummaryQScript {
     def addJobs(): Unit = {
       addPerLibJobs()
 
-      if (config("single_sample_variantcalling", default = false).asBoolean && preProcessBam.isDefined) {
-        val vc = makeVariantcalling(multisample = false)
-        vc.sampleId = Some(sampleId)
-        vc.outputDir = new File(sampleDir, "variantcalling")
-        vc.inputBams = preProcessBam.get :: Nil
-        vc.init
-        vc.biopetScript
-        addAll(vc.functions)
-        addSummaryQScript(vc)
+      if (preProcessBam.isDefined) {
+        val bamMetrics = new BamMetrics(root)
+        bamMetrics.sampleId = Some(sampleId)
+        bamMetrics.inputBam = preProcessBam.get
+        bamMetrics.outputDir = outputDir
+        bamMetrics.init
+        bamMetrics.biopetScript
+        addAll(bamMetrics.functions)
+        addSummaryQScript(bamMetrics)
+
+
+        if (config("single_sample_variantcalling", default = false).asBoolean) {
+          val vc = makeVariantcalling(multisample = false)
+          vc.sampleId = Some(sampleId)
+          vc.outputDir = new File(sampleDir, "variantcalling")
+          vc.inputBams = preProcessBam.get :: Nil
+          vc.init
+          vc.biopetScript
+          addAll(vc.functions)
+          addSummaryQScript(vc)
+        }
       }
     }
   }

public/shiva/src/main/scala/nl/lumc/sasc/biopet/pipelines/shiva/ShivaVariantcallingTrait.scala

@@ -5,7 +5,7 @@ import java.io.File
 import nl.lumc.sasc.biopet.core.{ PipelineCommand, SampleLibraryTag }
 import nl.lumc.sasc.biopet.core.summary.SummaryQScript
 import nl.lumc.sasc.biopet.extensions.gatk.CombineVariants
-import nl.lumc.sasc.biopet.tools.{ VcfFilter, MpileupToVcf }
+import nl.lumc.sasc.biopet.tools.{VcfStats, VcfFilter, MpileupToVcf}
 import nl.lumc.sasc.biopet.utils.ConfigUtils
 import org.broadinstitute.gatk.utils.commandline.Input
 
@@ -34,15 +34,22 @@ trait ShivaVariantcallingTrait extends SummaryQScript with SampleLibraryTag {
   def finalFile = new File(outputDir, namePrefix + "final.vcf.gz")
 
   def biopetScript: Unit = {
+    val callers = usedCallers
+
     val cv = new CombineVariants(qscript)
     cv.outputFile = finalFile
     cv.setKey = "VariantCaller"
-    for (caller <- usedCallers) {
+    for (caller <- callers.sortBy(_.prio)) {
       caller.addJobs()
       cv.addInput(caller.outputFile, caller.name)
     }
     add(cv)
 
+    val vcfStats = new VcfStats(qscript)
+    vcfStats.input = finalFile
+    vcfStats.setOutputDir(new File(outputDir, "vcfstats"))
+    add(vcfStats)
+
     addSummaryJobs
   }