Yamsvp.scala 5.52 KB
Newer Older
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.
 */
Wai Yi Leung's avatar
Wai Yi Leung committed
16
17
18
19
20
21
/*
 * Structural variation calling
 */

package nl.lumc.sasc.biopet.pipelines.yamsvp

22
23
import java.io.File

Wai Yi Leung's avatar
Wai Yi Leung committed
24
import nl.lumc.sasc.biopet.core.config.Configurable
Peter van 't Hof's avatar
Peter van 't Hof committed
25
import nl.lumc.sasc.biopet.core.{ BiopetQScript, MultiSampleQScript, PipelineCommand }
Wai Yi Leung's avatar
Wai Yi Leung committed
26

27
import nl.lumc.sasc.biopet.extensions.Ln
28
29
import nl.lumc.sasc.biopet.extensions.breakdancer.Breakdancer
import nl.lumc.sasc.biopet.extensions.clever.CleverCaller
30
import nl.lumc.sasc.biopet.extensions.igvtools.IGVToolsCount
31
import nl.lumc.sasc.biopet.extensions.sambamba.{ SambambaMerge, SambambaMarkdup }
32
33
//import nl.lumc.sasc.biopet.extensions.pindel.Pindel
import nl.lumc.sasc.biopet.extensions.delly.Delly
Wai Yi Leung's avatar
Wai Yi Leung committed
34
import nl.lumc.sasc.biopet.pipelines.bammetrics.BamMetrics
35

Wai Yi Leung's avatar
Wai Yi Leung committed
36
37
38
39
import nl.lumc.sasc.biopet.pipelines.mapping.Mapping

import org.broadinstitute.gatk.queue.QScript
import org.broadinstitute.gatk.queue.function._
Wai Yi Leung's avatar
Wai Yi Leung committed
40
41
import org.broadinstitute.gatk.queue.engine.JobRunInfo

42
43
class Yamsvp(val root: Configurable) extends QScript with MultiSampleQScript {
  qscript =>
Wai Yi Leung's avatar
Wai Yi Leung committed
44
45
  def this() = this(null)

Peter van 't Hof's avatar
Peter van 't Hof committed
46
  var reference: File = config("reference", required = true)
Wai Yi Leung's avatar
Wai Yi Leung committed
47

48
49
  def makeSample(id: String) = new Sample(id)
  class Sample(sampleId: String) extends AbstractSample(sampleId) {
Wai Yi Leung's avatar
Wai Yi Leung committed
50

51
52
    val alignmentDir: String = sampleDir + "alignment/"
    val svcallingDir: String = sampleDir + "svcalls/"
Wai Yi Leung's avatar
Wai Yi Leung committed
53

54
55
    def makeLibrary(id: String) = new Library(id)
    class Library(libraryId: String) extends AbstractLibrary(libraryId) {
Wai Yi Leung's avatar
Wai Yi Leung committed
56

57
      //      val runDir: String = alignmentDir + "run_" + libraryId + "/"
58

59
60
61
      val mapping = new Mapping(qscript)
      mapping.libraryId = libraryId
      mapping.sampleId = sampleId
62

63
64
65
66
      protected def addJobs(): Unit = {
        mapping.input_R1 = config("R1", required = true)
        mapping.input_R2 = config("R2", required = true)
        mapping.outputDir = libDir
67

68
69
70
71
        mapping.init
        mapping.biopetScript
        qscript.addAll(mapping.functions)
      }
Wai Yi Leung's avatar
Wai Yi Leung committed
72
    }
73
74
75
76
77
78
79
80
81
    protected def addJobs(): Unit = {
      addLibsJobs()
      val libraryBamfiles = libraries.map(_._2.mapping.finalBamFile).toList

      val bamFile: File = if (libraryBamfiles.size == 1) {
        val alignmentlink = Ln(qscript, libraryBamfiles.head,
          alignmentDir + sampleId + ".merged.bam", true)
        alignmentlink.isIntermediate = true
        add(alignmentlink)
82
        alignmentlink.out
83
      } else if (libraryBamfiles.size > 1) {
84
        val mergeSamFiles = new SambambaMerge(qscript)
85
        mergeSamFiles.input = libraryBamfiles
86
87
88
        mergeSamFiles.output = sampleDir + sampleId + ".merged.bam"
        mergeSamFiles.isIntermediate = true
        add(mergeSamFiles)
89
90
        mergeSamFiles.output
      } else null
Peter van 't Hof's avatar
Peter van 't Hof committed
91

92
93
      val bamMarkDup = SambambaMarkdup(qscript, bamFile)
      add(bamMarkDup)
94

95
      addAll(BamMetrics(qscript, bamMarkDup.output, alignmentDir + "metrics" + File.separator).functions)
Wai Yi Leung's avatar
Wai Yi Leung committed
96

97
98
99
      // create an IGV TDF file
      val tdfCount = IGVToolsCount(qscript, bamMarkDup.output, config("genome_name", default = "hg19"))
      add(tdfCount)
100

101
102
      /// bamfile will be used as input for the SV callers. First run Clever
      //    val cleverVCF : File = sampleDir + "/" + sampleID + ".clever.vcf"
103

104
105
106
      val cleverDir = svcallingDir + sampleId + ".clever/"
      val clever = CleverCaller(qscript, bamMarkDup.output, qscript.reference, svcallingDir, cleverDir)
      add(clever)
107

108
109
      val clever_vcf = Ln(qscript, clever.outputvcf, svcallingDir + sampleId + ".clever.vcf", relative = true)
      add(clever_vcf)
110

111
112
113
      val breakdancerDir = svcallingDir + sampleId + ".breakdancer/"
      val breakdancer = Breakdancer(qscript, bamMarkDup.output, qscript.reference, breakdancerDir)
      addAll(breakdancer.functions)
114

115
116
      val bd_vcf = Ln(qscript, breakdancer.outputvcf, svcallingDir + sampleId + ".breakdancer.vcf", relative = true)
      add(bd_vcf)
117

118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
      val dellyDir = svcallingDir + sampleId + ".delly/"
      val delly = Delly(qscript, bamMarkDup.output, dellyDir)
      addAll(delly.functions)

      val delly_vcf = Ln(qscript, delly.outputvcf, svcallingDir + sampleId + ".delly.vcf", relative = true)
      add(delly_vcf)

      // for pindel we should use per library config collected into one config file
      //    val pindelDir = svcallingDir + sampleID + ".pindel/"
      //    val pindel = Pindel(qscript, analysisBam, this.reference, pindelDir)
      //    sampleOutput.vcf += ("pindel" -> List(pindel.outputvcf))
      //    addAll(pindel.functions)
      //
      //    val pindel_vcf = Ln(qscript, pindel.outputvcf, svcallingDir + sampleID + ".pindel.vcf", relative = true)
      //    add(pindel_vcf)
      //
    }
  }
136

137
138
  def init() {
  }
139

140
141
142
  def biopetScript() {
    logger.info("Starting YAM SV Pipeline")
    addSamplesJobs
Wai Yi Leung's avatar
Wai Yi Leung committed
143
  }
144

145
146
  override def onExecutionDone(jobs: Map[QFunction, JobRunInfo], success: Boolean) {
    logger.info("YAM SV Pipeline has run .......................")
Wai Yi Leung's avatar
Wai Yi Leung committed
147
148
149
  }
}

150
object Yamsvp extends PipelineCommand