ShivaReport.scala 8.77 KB
Newer Older
1
2
package nl.lumc.sasc.biopet.pipelines.shiva

3
import java.io.{ PrintWriter, File }
4

5
import nl.lumc.sasc.biopet.core.config.Configurable
Peter van 't Hof's avatar
Peter van 't Hof committed
6
import nl.lumc.sasc.biopet.core.report._
7
import nl.lumc.sasc.biopet.core.summary.{ SummaryValue, Summary }
8
import nl.lumc.sasc.biopet.extensions.rscript.StackedBarPlot
9
import nl.lumc.sasc.biopet.pipelines.bammetrics.BammetricsReport
10
11
12
13
14
import nl.lumc.sasc.biopet.pipelines.flexiprep.FlexiprepReport

/**
 * Created by pjvan_thof on 3/30/15.
 */
15
16
17
18
class ShivaReport(val root: Configurable) extends ReportBuilderExtension {
  val builder = ShivaReport
}

19
/** Object for report generation for Shiva pipeline */
20
object ShivaReport extends MultisampleReportBuilder {
21

22
  /** Root page for the shiva report */
23
  def indexPage = {
24
    val regions = regionsPage
25
    ReportPage(
26
      List("Samples" -> generateSamplesPage(pageArgs)) ++
27
28
        (if (regions.isDefined) Map(regions.get) else Map()) ++
        Map("Files" -> filesPage,
29
          "Versions" -> ReportPage(List(), List((
30
31
32
            "Executables" -> ReportSection("/nl/lumc/sasc/biopet/core/report/executables.ssp"
            ))), Map())
        ),
33
34
35
      List(
        "Report" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/shiva/shivaFront.ssp"),
        "Variantcalling" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/shiva/sampleVariants.ssp",
Peter van 't Hof's avatar
Peter van 't Hof committed
36
          Map("showPlot" -> true, "showTable" -> false)),
37
        "Alignment" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/alignmentSummary.ssp",
Peter van 't Hof's avatar
Peter van 't Hof committed
38
          Map("sampleLevel" -> true, "showPlot" -> true, "showTable" -> false)
39
        ),
40
        "Insert Size" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/insertSize.ssp",
Peter van 't Hof's avatar
Peter van 't Hof committed
41
          Map("sampleLevel" -> true, "showPlot" -> true, "showTable" -> false)),
Peter van 't Hof's avatar
Peter van 't Hof committed
42
43
        "Whole genome coverage" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/wgsHistogram.ssp",
          Map("sampleLevel" -> true, "showPlot" -> true, "showTable" -> false)),
44
        "QC reads" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepReadSummary.ssp",
Peter van 't Hof's avatar
Peter van 't Hof committed
45
          Map("showPlot" -> true, "showTable" -> false)),
46
        "QC bases" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepBaseSummary.ssp",
Peter van 't Hof's avatar
Peter van 't Hof committed
47
          Map("showPlot" -> true, "showTable" -> false))
48
49
50
51
      ),
      pageArgs
    )
  }
52

53
54
55
  //TODO: Add variants per target
  /** Generate a page with all target coverage stats */
  def regionsPage: Option[(String, ReportPage)] = {
56
57
58
59
60
61
    val roi = summary.getValue("shiva", "settings", "regions_of_interest")
    val amplicon = summary.getValue("shiva", "settings", "amplicon_bed")

    var regionPages: Map[String, ReportPage] = Map()

    def createPage(name: String, amplicon: Boolean = false): ReportPage = {
Peter van 't Hof's avatar
Peter van 't Hof committed
62
      ReportPage(
63
        List(),
Peter van 't Hof's avatar
Peter van 't Hof committed
64
65
66
        List("Coverage" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/covstatsMultiTable.ssp")),
        Map("target" -> name)
      )
67
68
69
70
71
72
73
74
75
76
77
78
79
    }

    amplicon match {
      case Some(x: String) => regionPages += (x + " (Amplicon)") -> createPage(x, true)
      case _               =>
    }

    roi match {
      case Some(x: String)  => regionPages += x -> createPage(x, false)
      case Some(x: List[_]) => x.foreach(x => regionPages += x.toString -> createPage(x.toString, false))
      case _                =>
    }

Peter van 't Hof's avatar
Peter van 't Hof committed
80
    if (regionPages.nonEmpty) Some("Regions" -> ReportPage(
81
      List(),
Peter van 't Hof's avatar
Peter van 't Hof committed
82
83
84
85
86
      regionPages.map(p => (p._1 -> ReportSection(
        "/nl/lumc/sasc/biopet/pipelines/bammetrics/covstatsMultiTable.ssp",
        Map("target" -> p._1.stripSuffix(" (Amplicon)"))
      ))).toList.sortBy(_._1),
      Map())
87
88
    )
    else None
89
90
  }

91
92
  /** Files page, can be used general or at sample level */
  def filesPage: ReportPage = ReportPage(List(), List(
93
94
    "Input fastq files" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepInputfiles.ssp"),
    "After QC fastq files" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepOutputfiles.ssp"),
95
96
97
    "Bam files per lib" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/mapping/outputBamfiles.ssp", Map("sampleLevel" -> false)),
    "Preprocessed bam files" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/mapping/outputBamfiles.ssp",
      Map("pipelineName" -> "shiva", "fileTag" -> "preProcessBam")),
98
99
100
    "VCF files" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/shiva/outputVcfFiles.ssp", Map("sampleId" -> None))
  ), Map())

101
102
  /** Single sample page */
  def samplePage(sampleId: String, args: Map[String, Any]): ReportPage = {
103
    ReportPage(List(
104
      "Libraries" -> generateLibraryPage(args),
105
      "Alignment" -> BammetricsReport.bamMetricsPage(summary, Some(sampleId), None),
106
      "Files" -> filesPage
107
    ), List(
108
109
      "Alignment" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/alignmentSummary.ssp",
        if (summary.libraries(sampleId).size > 1) Map("showPlot" -> true) else Map()),
110
      "Preprocessing" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/alignmentSummary.ssp", Map("sampleLevel" -> true)),
Peter van 't Hof's avatar
Peter van 't Hof committed
111
      "Variantcalling" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/shiva/sampleVariants.ssp"),
112
113
      "QC reads" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepReadSummary.ssp"),
      "QC bases" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepBaseSummary.ssp")
114
115
116
    ), args)
  }

117
118
  /** Library page */
  def libraryPage(sampleId: String, libId: String, args: Map[String, Any]): ReportPage = {
119
    ReportPage(List(
120
      "Alignment" -> BammetricsReport.bamMetricsPage(summary, Some(sampleId), Some(libId)),
121
122
123
      "QC" -> FlexiprepReport.flexiprepPage
    ), List(
      "Alignment" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/bammetrics/alignmentSummary.ssp"),
124
125
126
      "QC reads" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepReadSummary.ssp"),
      "QC bases" -> ReportSection("/nl/lumc/sasc/biopet/pipelines/flexiprep/flexiprepBaseSummary.ssp")
    ), args)
127
128
  }

129
  /** Name of the report */
130
  def reportName = "Shiva Report"
131

132
133
134
135
136
137
138
139
  /**
   * Generate a stackbar plot for found variants
   * @param outputDir OutputDir for the tsv and png file
   * @param prefix Prefix of the tsv and png file
   * @param summary Summary class
   * @param libraryLevel Default false, when set true plot will be based on library stats instead of sample stats
   * @param sampleId Default it selects all sampples, when sample is giving it limits to selected sample
   */
140
  def variantSummaryPlot(outputDir: File,
141
142
143
144
                         prefix: String,
                         summary: Summary,
                         libraryLevel: Boolean = false,
                         sampleId: Option[String] = None): Unit = {
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
    val tsvFile = new File(outputDir, prefix + ".tsv")
    val pngFile = new File(outputDir, prefix + ".png")
    val tsvWriter = new PrintWriter(tsvFile)
    if (libraryLevel) tsvWriter.print("Library") else tsvWriter.print("Sample")
    tsvWriter.println("\tHomVar\tHet\tHomRef\tNoCall")

    def getLine(summary: Summary, sample: String, lib: Option[String] = None): String = {
      val homVar = new SummaryValue(List("shivavariantcalling", "stats", "multisample-vcfstats-final", "genotype", "HomVar"),
        summary, Some(sample), lib).value.getOrElse(0).toString.toLong
      val homRef = new SummaryValue(List("shivavariantcalling", "stats", "multisample-vcfstats-final", "genotype", "HomRef"),
        summary, Some(sample), lib).value.getOrElse(0).toString.toLong
      val noCall = new SummaryValue(List("shivavariantcalling", "stats", "multisample-vcfstats-final", "genotype", "NoCall"),
        summary, Some(sample), lib).value.getOrElse(0).toString.toLong
      val het = new SummaryValue(List("shivavariantcalling", "stats", "multisample-vcfstats-final", "genotype", "Het"),
        summary, Some(sample), lib).value.getOrElse(0).toString.toLong
      val sb = new StringBuffer()
      if (lib.isDefined) sb.append(sample + "-" + lib.get + "\t") else sb.append(sample + "\t")
      sb.append(homVar + "\t")
      sb.append(het + "\t")
      sb.append(homRef + "\t")
      sb.append(noCall)
      sb.toString
    }

    if (libraryLevel) {
170
171
172
173
      for (
        sample <- summary.samples if (sampleId.isEmpty || sample == sampleId.get);
        lib <- summary.libraries(sample)
      ) {
174
175
176
177
178
179
180
181
182
183
184
185
186
187
        tsvWriter.println(getLine(summary, sample, Some(lib)))
      }
    } else {
      for (sample <- summary.samples if (sampleId.isEmpty || sample == sampleId.get)) {
        tsvWriter.println(getLine(summary, sample))
      }
    }

    tsvWriter.close()

    val plot = new StackedBarPlot(null)
    plot.input = tsvFile
    plot.output = pngFile
    plot.ylabel = Some("VCF records")
188
189
190
    if (libraryLevel) {
      plot.width = Some(200 + (summary.libraries.filter(s => sampleId.getOrElse(s._1) == s._1).foldLeft(0)(_ + _._2.size) * 10))
    } else plot.width = Some(200 + (summary.samples.filter(s => sampleId.getOrElse(s) == s).size * 10))
191
192
    plot.runLocal()
  }
193
}