BedtoolsCoverageToCounts.scala 2.78 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.
 */
16
package nl.lumc.sasc.biopet.tools
17

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

20
import nl.lumc.sasc.biopet.core.config.Configurable
Peter van 't Hof's avatar
Peter van 't Hof committed
21
import nl.lumc.sasc.biopet.core.{ ToolCommand, ToolCommandFuntion }
22
import org.broadinstitute.gatk.utils.commandline.{ Input, Output }
Peter van 't Hof's avatar
Peter van 't Hof committed
23

Peter van 't Hof's avatar
Peter van 't Hof committed
24
import scala.collection.{ mutable, SortedMap }
25 26
import scala.io.Source

27
class BedtoolsCoverageToCounts(val root: Configurable) extends ToolCommandFuntion {
28 29 30 31
  javaMainClass = getClass.getName

  @Input(doc = "Input fasta", shortName = "input", required = true)
  var input: File = _
Peter van 't Hof's avatar
Peter van 't Hof committed
32

33 34
  @Output(doc = "Output tag library", shortName = "output", required = true)
  var output: File = _
Peter van 't Hof's avatar
Peter van 't Hof committed
35

Peter van 't Hof's avatar
Peter van 't Hof committed
36
  override val defaultCoreMemory = 3.0
Peter van 't Hof's avatar
Peter van 't Hof committed
37 38 39

  override def commandLine = super.commandLine +
    required("-I", input) +
40 41 42
    required("-o", output)
}

43
object BedtoolsCoverageToCounts extends ToolCommand {
Peter van 't Hof's avatar
Peter van 't Hof committed
44
  case class Args(input: File = null, output: File = null) extends AbstractArgs
45 46

  class OptParser extends AbstractOptParser {
Peter van 't Hof's avatar
Peter van 't Hof committed
47
    opt[File]('I', "input") required () valueName "<file>" action { (x, c) =>
Peter van 't Hof's avatar
Peter van 't Hof committed
48 49
      c.copy(input = x)
    }
Peter van 't Hof's avatar
Peter van 't Hof committed
50
    opt[File]('o', "output") required () unbounded () valueName "<file>" action { (x, c) =>
Peter van 't Hof's avatar
Peter van 't Hof committed
51 52
      c.copy(output = x)
    }
53
  }
Peter van 't Hof's avatar
Peter van 't Hof committed
54

55 56 57 58
  /**
   * @param args the command line arguments
   */
  def main(args: Array[String]): Unit = {
59 60
    val argsParser = new OptParser
    val commandArgs: Args = argsParser.parse(args, Args()) getOrElse sys.exit(1)
Peter van 't Hof's avatar
Peter van 't Hof committed
61

62
    if (!commandArgs.input.exists) throw new IllegalStateException("Input file not found, file: " + commandArgs.input)
Peter van 't Hof's avatar
Peter van 't Hof committed
63

Peter van 't Hof's avatar
Peter van 't Hof committed
64 65
    val counts: mutable.Map[String, Long] = mutable.Map()
    for (line <- Source.fromFile(commandArgs.input).getLines()) {
66 67 68 69 70 71
      val values = line.split("\t")
      val gene = values(3)
      val count = values(6).toLong
      if (counts.contains(gene)) counts(gene) += count
      else counts += gene -> count
    }
Peter van 't Hof's avatar
Peter van 't Hof committed
72 73 74

    val sortedCounts: SortedMap[String, Long] = SortedMap(counts.toArray: _*)

75
    val writer = new PrintWriter(commandArgs.output)
Peter van 't Hof's avatar
Peter van 't Hof committed
76
    for ((seq, count) <- sortedCounts) {
77 78
      if (count > 0) writer.println(seq + "\t" + count)
    }
Peter van 't Hof's avatar
Peter van 't Hof committed
79
    writer.close()
80 81
  }
}