diff --git a/src/main/scala/com/fulcrumgenomics/umi/DuplexConsensusCaller.scala b/src/main/scala/com/fulcrumgenomics/umi/DuplexConsensusCaller.scala index d402cfe40..00434b71c 100644 --- a/src/main/scala/com/fulcrumgenomics/umi/DuplexConsensusCaller.scala +++ b/src/main/scala/com/fulcrumgenomics/umi/DuplexConsensusCaller.scala @@ -340,11 +340,16 @@ class DuplexConsensusCaller(override val readNamePrefix: String, private[umi] def duplexConsensus(ab: Option[VanillaConsensusRead], ba: Option[VanillaConsensusRead], sourceReads: Seq[SourceRead]): Option[DuplexConsensusRead] = { - (ab, ba) match { + // Calculate the length that we'll create a duplex read as "the shorter of the available reads", + // and then filter each read to make sure it has some coverage in that section of the read + val len = min(ab.map(_.length).getOrElse(Int.MaxValue), ba.map(_.length).getOrElse(Int.MaxValue)) + val abX = ab.filter(_.depths.iterator.take(len).exists(_ > 0)) + val baX = ba.filter(_.depths.iterator.take(len).exists(_ > 0)) + + (abX, baX) match { case (Some(a), None) => Some(DuplexConsensusRead(id=a.id, a.bases, a.quals, a.errors, a, None)) case (None, Some(b)) => Some(DuplexConsensusRead(id=b.id, b.bases, b.quals, b.errors, b, None)) case (Some(a), Some(b)) => - val len = min(a.length, b.length) val id = a.id val bases = new Array[Byte](len) val quals = new Array[Byte](len) diff --git a/src/test/scala/com/fulcrumgenomics/umi/DuplexConsensusCallerTest.scala b/src/test/scala/com/fulcrumgenomics/umi/DuplexConsensusCallerTest.scala index 25047dcbc..e51fad2bd 100644 --- a/src/test/scala/com/fulcrumgenomics/umi/DuplexConsensusCallerTest.scala +++ b/src/test/scala/com/fulcrumgenomics/umi/DuplexConsensusCallerTest.scala @@ -28,8 +28,9 @@ import com.fulcrumgenomics.testing.SamBuilder.{Minus, Plus} import com.fulcrumgenomics.testing.{SamBuilder, UnitSpec} import com.fulcrumgenomics.umi.ConsensusTags.PerRead.{AbRawReadErrorRate, BaRawReadErrorRate, RawReadErrorRate} import com.fulcrumgenomics.util.NumericTypes.PhredScore +import org.scalatest.OptionValues -class DuplexConsensusCallerTest extends UnitSpec { +class DuplexConsensusCallerTest extends UnitSpec with OptionValues { // Function to create a caller without so much writing def caller(q: Int = 10, pre: Int = DuplexConsensusCaller.ErrorRatePreUmi, post: Int = DuplexConsensusCaller.ErrorRatePostUmi, minReads: Seq[Int] = Seq(1)) = new DuplexConsensusCaller(readNamePrefix="test", minInputBaseQuality=q.toByte, errorRatePreUmi=pre.toByte, errorRatePostUmi=post.toByte, @@ -468,4 +469,29 @@ class DuplexConsensusCallerTest extends UnitSpec { } } } + "DuplexConsensusCaller.duplexConsensus" should "swap AB and BA reads when AB reads have zero depth across the minimum length" in { + val cc = caller(minReads=Seq(1, 1, 0)) + val ab = VanillaConsensusRead( + id = "ab", + bases = "NNN".getBytes, + quals = "III".getBytes, + depths = Array(0, 0, 0), + errors = Array(0, 0, 0) + ) + val ba = ab.copy(depths=Array(1, 1, 1)) + + // ab has zero depth, so it should be swapped with BA + { + val read = cc.duplexConsensus(ab=Some(ab), ba=Some(ba), sourceReads=Seq.empty).value + read.abConsensus shouldBe ba + read.baConsensus.isEmpty shouldBe true + } + + // ba has zero depth, so it should be undefined. + { + val read = cc.duplexConsensus(ab=Some(ba), ba=Some(ab), sourceReads=Seq.empty).value + read.abConsensus shouldBe ba + read.baConsensus.isEmpty shouldBe true + } + } }