From 162f02edbb6e7b8eb382f049d598c46a92afb381 Mon Sep 17 00:00:00 2001 From: Peter van 't Hof <p.j.van_t_hof@lumc.nl> Date: Mon, 1 Feb 2016 14:33:53 +0100 Subject: [PATCH] Added test on main method --- .../lumc/sasc/biopet/tools/BaseCounter.scala | 42 +++++++++--------- .../src/test/resources/chrQ.refflat | 1 + .../biopet-tools/src/test/resources/empty.bai | Bin 0 -> 104 bytes .../biopet-tools/src/test/resources/empty.bam | Bin 0 -> 189 bytes .../sasc/biopet/tools/BaseCounterTest.scala | 20 +++++++++ 5 files changed, 43 insertions(+), 20 deletions(-) create mode 100644 public/biopet-tools/src/test/resources/chrQ.refflat create mode 100644 public/biopet-tools/src/test/resources/empty.bai create mode 100644 public/biopet-tools/src/test/resources/empty.bam diff --git a/public/biopet-tools/src/main/scala/nl/lumc/sasc/biopet/tools/BaseCounter.scala b/public/biopet-tools/src/main/scala/nl/lumc/sasc/biopet/tools/BaseCounter.scala index 756e187a6..7a2b0907f 100644 --- a/public/biopet-tools/src/main/scala/nl/lumc/sasc/biopet/tools/BaseCounter.scala +++ b/public/biopet-tools/src/main/scala/nl/lumc/sasc/biopet/tools/BaseCounter.scala @@ -56,7 +56,7 @@ object BaseCounter extends ToolCommand { val overlapGenes = groupGenesOnOverlap(geneReader.getAll) logger.info("Start reading bamFile") - val counts = (for (genes <- overlapGenes.values.flatten.par) yield bamToGeneCount(cmdArgs.bamFile, genes)).toList + val counts = (for (genes <- overlapGenes.values.flatten.par) yield runThread(cmdArgs.bamFile, genes)).toList logger.info("Done reading bamFile") writeGeneCounts(counts.flatMap(_.geneCounts), cmdArgs.outputDir, cmdArgs.prefix) @@ -281,7 +281,7 @@ object BaseCounter extends ToolCommand { nonStrandedMetaExonCounts: List[(String, RegionCount)], strandedMetaExonCounts: List[(String, RegionCount)]) - def bamToGeneCount(bamFile: File, genes: List[Gene]): ThreadOutput = { + private[tools] def runThread(bamFile: File, genes: List[Gene]): ThreadOutput = { val counts = genes.map(gene => gene -> new GeneCount(gene)).toMap val bamReader = SamReaderFactory.makeDefault().open(bamFile) @@ -304,25 +304,27 @@ object BaseCounter extends ToolCommand { } def createMetaExonCounts(genes: List[Gene]): List[(String, RegionCount)] = { - val regions = genes.map(gene => gene.getName -> generateMergedExonRegions(gene).sorted) - val chr = genes.head.getContig - val begin = regions.map(_._2.allRecords.head.start).min - val end = regions.map(_._2.allRecords.last.end).max - - val posibleEnds = (regions.flatMap(_._2.allRecords.map(_.end)) ++ regions.flatMap(_._2.allRecords.map(_.start))).distinct.sorted - - def mergeRegions(newBegin: Int, output: List[(String, RegionCount)] = Nil): List[(String, RegionCount)] = { - val newEnds = posibleEnds.filter(_ > newBegin) - if (newBegin > end || newEnds.isEmpty) output - else { - val newEnd = newEnds.min - val record = BedRecord(chr, newBegin, newEnd) - val names = regions.filter(_._2.overlapWith(record).nonEmpty).map(_._1) - if (names.nonEmpty) mergeRegions(newEnd, (names.mkString(","), new RegionCount(record.start + 1, record.end)) :: output) - else mergeRegions(newEnd, output) + if (genes.nonEmpty) { + val regions = genes.map(gene => gene.getName -> generateMergedExonRegions(gene).sorted) + val chr = genes.head.getContig + val begin = regions.map(_._2.allRecords.head.start).min + val end = regions.map(_._2.allRecords.last.end).max + + val posibleEnds = (regions.flatMap(_._2.allRecords.map(_.end)) ++ regions.flatMap(_._2.allRecords.map(_.start))).distinct.sorted + + def mergeRegions(newBegin: Int, output: List[(String, RegionCount)] = Nil): List[(String, RegionCount)] = { + val newEnds = posibleEnds.filter(_ > newBegin) + if (newBegin > end || newEnds.isEmpty) output + else { + val newEnd = newEnds.min + val record = BedRecord(chr, newBegin, newEnd) + val names = regions.filter(_._2.overlapWith(record).nonEmpty).map(_._1) + if (names.nonEmpty) mergeRegions(newEnd, (names.mkString(","), new RegionCount(record.start + 1, record.end)) :: output) + else mergeRegions(newEnd, output) + } } - } - mergeRegions(begin) + mergeRegions(begin) + } else Nil } def bamRecordBasesOverlap(samRecord: SAMRecord, start: Int, end: Int): Int = { diff --git a/public/biopet-tools/src/test/resources/chrQ.refflat b/public/biopet-tools/src/test/resources/chrQ.refflat new file mode 100644 index 000000000..dfcc51100 --- /dev/null +++ b/public/biopet-tools/src/test/resources/chrQ.refflat @@ -0,0 +1 @@ +geneA t1 chrQ + 200 500 225 475 3 200,300,400 250,350,500 diff --git a/public/biopet-tools/src/test/resources/empty.bai b/public/biopet-tools/src/test/resources/empty.bai new file mode 100644 index 0000000000000000000000000000000000000000..e4e22396b4877369154096123eaf25e8bb47503a GIT binary patch literal 104 ucmZ>A^kigWU|;}YPay^dMj%}S#2~N_LNRzvg2<zbgVZ7eWWC7ZAT|JCaRmMV literal 0 HcmV?d00001 diff --git a/public/biopet-tools/src/test/resources/empty.bam b/public/biopet-tools/src/test/resources/empty.bam new file mode 100644 index 0000000000000000000000000000000000000000..bd4073d80cbf73c689e97d2814f9b1fd00d63895 GIT binary patch literal 189 zcmb2|=3rp}f&Xj_PR>jW3mA%vzNGRbCnOYnD0s;8d9%?K<EK2|wtJmBukUs8jJJ;G z(;^POCp_Oad!5tZdh&%)%YcDN<)8u(mHYXNYcn%DpG$~eJ?Y9LAw9JPC)YUeTso_; zfbTg=qh%S7K!IRKWR6X%M5Bqg!2fop2{98J1rIYR++D%cq_t{QR@SRmQc{lZ0v)O) X7(`S2e(eLgOCHUE(hSUC_ksujDyu*2 literal 0 HcmV?d00001 diff --git a/public/biopet-tools/src/test/scala/nl/lumc/sasc/biopet/tools/BaseCounterTest.scala b/public/biopet-tools/src/test/scala/nl/lumc/sasc/biopet/tools/BaseCounterTest.scala index 2ec0b6080..09a7fde18 100644 --- a/public/biopet-tools/src/test/scala/nl/lumc/sasc/biopet/tools/BaseCounterTest.scala +++ b/public/biopet-tools/src/test/scala/nl/lumc/sasc/biopet/tools/BaseCounterTest.scala @@ -1,5 +1,9 @@ package nl.lumc.sasc.biopet.tools +import java.io.File +import java.nio.file.Paths + +import com.google.common.io.Files import htsjdk.samtools.{SAMReadGroupRecord, SAMSequenceRecord, SAMLineParser, SAMFileHeader} import org.scalatest.Matchers import org.scalatest.testng.TestNGSuite @@ -164,6 +168,22 @@ class BaseCounterTest extends TestNGSuite with Matchers { assert(ab.exists(x => x._1 == "geneB" && x._2.start == 201 && x._2.end == 210)) assert(ab.exists(x => x._1 == "geneB" && x._2.start == 221 && x._2.end == 250)) } + + private def resourcePath(p: String): String = { + Paths.get(getClass.getResource(p).toURI).toString + } + + @Test + def testMain: Unit = { + val outputDir = Files.createTempDir() + outputDir.deleteOnExit() + val prefix = "test" + val bamFile = new File(resourcePath("/empty.bam")) + val refflat = new File(resourcePath("/chrQ.refflat")) + main(Array("-o", outputDir.getAbsolutePath, "-p", prefix, + "-b", bamFile.getAbsolutePath, "-r", refflat.getAbsolutePath)) + outputDir.list().size shouldBe 34 + } } object BaseCounterTest { -- GitLab