Skip to content

Commit

Permalink
[ADAM-1569] Support log space probabilities in PhredUtils for extreme…
Browse files Browse the repository at this point in the history 
…ly high Phred scores.

Resolves #1569. Can support phred scores >1,000.
  • Loading branch information
@fnothaft
fnothaft committed Jun 22, 2017
1 parent 706c428 commit 412d092
Show file tree
Hide file tree
Showing 7 changed files with 88 additions and 31 deletions.
16 changes: 8 additions & 8 deletions adam-core/src/main/scala/org/bdgenomics/adam/converters/VariantContextConverter.scala
View file
Original file line number Diff line number Diff line change
Expand Up @@ -303,7 +303,7 @@ class VariantContextConverter(
// is the last allele the non-ref allele?
val alleles = vc.getAlternateAlleles.toSeq
val referenceModelIndex = if (alleles.nonEmpty && alleles.last == NON_REF_ALLELE) {
Some(alleles.length - 1)
Some(alleles.length)
} else {
None
}
Expand Down Expand Up @@ -818,7 +818,7 @@ class VariantContextConverter(
val pl = g.getPL
try {
val likelihoods = gIndices.map(idx => {
jFloat(PhredUtils.phredToLogProbability(pl(idx)))
jFloat(PhredUtils.phredToLogProbability(pl(idx)).toFloat)
}).toList
gb.setGenotypeLikelihoods(likelihoods)
} catch {
Expand All @@ -838,7 +838,7 @@ class VariantContextConverter(
if (g.hasPL) {
val pl = g.getPL
gb.setNonReferenceLikelihoods(gIndices.map(idx => {
jFloat(PhredUtils.phredToLogProbability(pl(idx)))
jFloat(PhredUtils.phredToLogProbability(pl(idx)).toFloat)
}).toList)
} else {
gb
Expand Down Expand Up @@ -979,12 +979,12 @@ class VariantContextConverter(
val array = Array.fill(elements) { Int.MinValue }

(0 to copyNumber).foreach(idx => {
array(idx) = PhredUtils.logProbabilityToPhred(gls.get(idx))
array(idx) = PhredUtils.logProbabilityToPhred(gls.get(idx).toDouble)
})

var cnIdx = copyNumber + 1
(1 to copyNumber).foreach(idx => {
array(cnIdx + 1) = PhredUtils.logProbabilityToPhred(nls.get(idx))
array(cnIdx + 1) = PhredUtils.logProbabilityToPhred(nls.get(idx).toDouble)
cnIdx += (copyNumber - idx)
})

Expand All @@ -1004,11 +1004,11 @@ class VariantContextConverter(
}
gb.noPL
} else {
gb.PL(nls.map(l => PhredUtils.logProbabilityToPhred(l))
gb.PL(nls.map(l => PhredUtils.logProbabilityToPhred(l.toDouble))
.toArray)
}
} else if (nls.isEmpty) {
gb.PL(gls.map(l => PhredUtils.logProbabilityToPhred(l))
gb.PL(gls.map(l => PhredUtils.logProbabilityToPhred(l.toDouble))
.toArray)
} else {
gb.PL(nonRefPls(gls, nls))
Expand Down Expand Up @@ -1663,7 +1663,7 @@ class VariantContextConverter(
val coreWithOptNonRefs = nonRefIndex.fold(convertedCore)(nonRefAllele => {

// non-ref pl indices
val nrIndices = GenotypeLikelihoods.getPLIndecesOfAlleles(alleleIdx, nonRefAllele)
val nrIndices = GenotypeLikelihoods.getPLIndecesOfAlleles(0, nonRefAllele)

formatNonRefGenotypeLikelihoods(g, convertedCore, nrIndices)
})
Expand Down
29 changes: 23 additions & 6 deletions adam-core/src/main/scala/org/bdgenomics/adam/util/PhredUtils.scala
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,15 @@
*/
package org.bdgenomics.adam.util

import scala.math.{ exp, pow, log, log10 }
import scala.math.{
exp,
expm1,
pow,
log,
log1p,
log10,
round
}

/**
* Helper singleton for converting Phred scores to/from probabilities.
Expand All @@ -36,11 +44,20 @@ object PhredUtils extends Serializable {
phredToErrorProbabilityCache.map { p => 1.0 - p }
}

private val MLOG10_DIV10 = -log(10.0) / 10.0
private val M10_DIV_LOG10 = -10.0 / log(10.0)

/**
* @param phred The phred score to convert.
* @return The input phred score as a log success probability.
*/
def phredToLogProbability(phred: Int): Float = log(phredToSuccessProbability(phred)).toFloat
def phredToLogProbability(phred: Int): Double = {
if (phred < 255) {
log(phredToSuccessProbability(phred)).toFloat
} else {
log1p(-exp(phred * MLOG10_DIV10))
}
}

/**
* @param phred The phred score to convert.
Expand Down Expand Up @@ -82,11 +99,11 @@ object PhredUtils extends Serializable {
/**
* @param p A log-scaled success probability to conver to phred.
* @return Returns this probability as a Phred score. If the log value is 0.0,
* we clip the phred score to 256.
* we clip the phred score to Int.MaxValue.
*/
def logProbabilityToPhred(p: Float): Int = if (p == 0.0) {
256
def logProbabilityToPhred(p: Double): Int = if (p == 0.0) {
Int.MaxValue
} else {
successProbabilityToPhred(exp(p))
round(M10_DIV_LOG10 * log(-expm1(p))).toInt
}
}
8 changes: 4 additions & 4 deletions adam-core/src/test/resources/gvcf_multiallelic/multiallelic.vcf
View file
Original file line number Diff line number Diff line change
Expand Up @@ -134,8 +134,8 @@
##reference=file:///dlmp/misc-data/reference/ggpsReferences/genomes/withMT/hs37d5.chr.fa
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA12878i
chr22 16157521 . C <NON_REF> . . END=16157602 GT:DP:GQ:MIN_DP:PL 0/0:4:9:2:0,0,45
chr22 16157603 . G C,<NON_REF> . . DP=1;MLEAC=1;MLEAF=0.5;MQ=59.0;MQ0=0 GT:AD:DP:GQ:PL 1/1:0,1:1:3:41,3,0,-2147483648,3,41
chr22 16157603 . G C,<NON_REF> . . DP=1;MLEAC=1;MLEAF=0.5;MQ=59.0;MQ0=0 GT:AD:DP:GQ:PL 1/1:0,1:1:3:41,3,0,-2147483648,41,41
chr22 16157604 . G <NON_REF> . . END=16157639 GT:DP:GQ:MIN_DP:PL 0/0:2:0:1:0,0,0
chr22 18030096 . TAAA TA,<NON_REF> . . BaseQRankSum=-0.133;ClippingRankSum=-1.438;DP=114;MLEAC=1;MLEAF=0.5;MQ=69.72;MQ0=0;MQRankSum=-0.686;ReadPosRankSum=-0.013 GT:AD:DP:GQ:PL 1/.:13,17:50:86:256,86,256,-2147483648,0,256
chr22 18030096 . TAAA T,<NON_REF> . . BaseQRankSum=-0.133;ClippingRankSum=-1.438;DP=114;MLEAC=0;MLEAF=0.0;MQ=69.72;MQ0=0;MQRankSum=-0.686;ReadPosRankSum=-0.013 GT:AD:DP:GQ:PL ./.:13,3:50:86:256,256,256,-2147483648,256,256
chr22 18030096 . TAAA TAA,<NON_REF> . . BaseQRankSum=-0.133;ClippingRankSum=-1.438;DP=114;MLEAC=1;MLEAF=0.5;MQ=69.72;MQ0=0;MQRankSum=-0.686;ReadPosRankSum=-0.013 GT:AD:DP:GQ:PL ./1:13,17:50:86:256,137,256,-2147483648,256,256
chr22 18030096 . TAAA T,<NON_REF> . . BaseQRankSum=-0.133;ClippingRankSum=-1.438;DP=114;MLEAC=0;MLEAF=0.0;MQ=69.72;MQ0=0;MQRankSum=-0.686;ReadPosRankSum=-0.013 GT:AD:DP:GQ:PL ./.:13,3:50:86:2147483647,2147483647,2147483647,-2147483648,2147483647,2147483647
chr22 18030096 . TAAA TA,<NON_REF> . . BaseQRankSum=-0.133;ClippingRankSum=-1.438;DP=114;MLEAC=1;MLEAF=0.5;MQ=69.72;MQ0=0;MQRankSum=-0.686;ReadPosRankSum=-0.013 GT:AD:DP:GQ:PL 1/.:13,17:50:86:2147483647,86,2147483647,-2147483648,2147483647,2147483647
chr22 18030096 . TAAA TAA,<NON_REF> . . BaseQRankSum=-0.133;ClippingRankSum=-1.438;DP=114;MLEAC=1;MLEAF=0.5;MQ=69.72;MQ0=0;MQRankSum=-0.686;ReadPosRankSum=-0.013 GT:AD:DP:GQ:PL ./1:13,17:50:86:2147483647,137,281,-2147483648,2147483647,2147483647
12 changes: 6 additions & 6 deletions adam-core/src/test/resources/sorted.lex.vcf
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,9 +66,9 @@
##contig=<ID=2,length=249250621>
##contig=<ID=13,length=249250621>
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA12878 NA12891 NA12892
1 14397 . CTGT C . IndelQD AC=2;AF=0.333;AN=6;BaseQRankSum=1.8;ClippingRankSum=0.138;DP=69;FS=7.786;MLEAC=2;MLEAF=0.333;MQ=26.84;MQ0=0;MQRankSum=-1.906;QD=1.55;ReadPosRankSum=0.384 GT:AD:DP:FT:GQ:PL 0/1:16,4:20:rd:99:120,0,256 0/1:8,2:10:dp;rd:60:60,0,256 0/0:39,0:39:PASS:99:0,116,256
1 14522 . G A . VQSRTrancheSNP99.95to100.00 AC=2;AF=0.333;AN=6;BaseQRankSum=2.044;ClippingRankSum=-2.196;DP=48;FS=13.179;MLEAC=2;MLEAF=0.333;MQ=25.89;MQ0=0;MQRankSum=-0.063;QD=8.87;ReadPosRankSum=0.952;VQSLOD=-3.333;culprit=MQ GT:AD:DP:FT:GQ:PL 0/1:10,5:15:dp:99:99,0,256 0/1:2,5:7:dp;rd:34:128,0,34 0/0:26,0:26:PASS:78:0,78,256
1 63735 rs201888535 CCTA C . PASS AC=1;AF=0.167;AN=6;BaseQRankSum=1.138;ClippingRankSum=0.448;DB;DP=176;FS=13.597;MLEAC=1;MLEAF=0.167;MQ=31.06;MQ0=0;MQRankSum=0.636;QD=9.98;ReadPosRankSum=-1.18 GT:AD:DP:FT:GQ:PL 0/0:27,0:27:PASS:79:0,79,256 0/0:40,0:40:PASS:99:0,117,256 0/1:23,74:97:rd:99:256,0,256
13 752721 rs3131972 A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL 1/1:0,27:27:PASS:81:256,81,0 1/1:0,19:19:dp:57:256,57,0 1/1:0,22:22:PASS:66:256,66,0
13 752791 . A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL:SB 1/1:0,27:27:PASS:81:256,81,0:0,1,2,3 1/1:0,19:19:dp:57:256,57,0:4,5,6,7 1/1:0,22:22:PASS:66:256,66,0:2,3,4,5
2 19190 . GC G . PASS AC=3;AF=0.5;AN=6;BaseQRankSum=4.157;ClippingRankSum=3.666;DP=74;FS=37.037;MLEAC=3;MLEAF=0.5;MQ=22.26;MQ0=0;MQRankSum=0.195;QD=16.04;ReadPosRankSum=-4.072 GT:AD:DP:FT:GQ:PL 0/1:8,14:22:PASS:99:256,0,256 0/1:18,13:31:PASS:99:256,0,256 0/1:5,15:20:rd:99:256,0,107
1 14397 . CTGT C . IndelQD AC=2;AF=0.333;AN=6;BaseQRankSum=1.8;ClippingRankSum=0.138;DP=69;FS=7.786;MLEAC=2;MLEAF=0.333;MQ=26.84;MQ0=0;MQRankSum=-1.906;QD=1.55;ReadPosRankSum=0.384 GT:AD:DP:FT:GQ:PL 0/1:16,4:20:rd:99:120,0,2147483647 0/1:8,2:10:dp;rd:60:60,0,414 0/0:39,0:39:PASS:99:0,116,2147483647
1 14522 . G A . VQSRTrancheSNP99.95to100.00 AC=2;AF=0.333;AN=6;BaseQRankSum=2.044;ClippingRankSum=-2.196;DP=48;FS=13.179;MLEAC=2;MLEAF=0.333;MQ=25.89;MQ0=0;MQRankSum=-0.063;QD=8.87;ReadPosRankSum=0.952;VQSLOD=-3.333;culprit=MQ GT:AD:DP:FT:GQ:PL 0/1:10,5:15:dp:99:99,0,2147483647 0/1:2,5:7:dp;rd:34:128,0,34 0/0:26,0:26:PASS:78:0,78,2147483647
1 63735 rs201888535 CCTA C . PASS AC=1;AF=0.167;AN=6;BaseQRankSum=1.138;ClippingRankSum=0.448;DB;DP=176;FS=13.597;MLEAC=1;MLEAF=0.167;MQ=31.06;MQ0=0;MQRankSum=0.636;QD=9.98;ReadPosRankSum=-1.18 GT:AD:DP:FT:GQ:PL 0/0:27,0:27:PASS:79:0,79,2147483647 0/0:40,0:40:PASS:99:0,117,2147483647 0/1:23,74:97:rd:99:2147483647,0,2147483647
13 752721 rs3131972 A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL 1/1:0,27:27:PASS:81:2147483647,81,0 1/1:0,19:19:dp:57:2147483647,57,0 1/1:0,22:22:PASS:66:2147483647,66,0
13 752791 . A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL:SB 1/1:0,27:27:PASS:81:2147483647,81,0:0,1,2,3 1/1:0,19:19:dp:57:2147483647,57,0:4,5,6,7 1/1:0,22:22:PASS:66:2147483647,66,0:2,3,4,5
2 19190 . GC G . PASS AC=3;AF=0.5;AN=6;BaseQRankSum=4.157;ClippingRankSum=3.666;DP=74;FS=37.037;MLEAC=3;MLEAF=0.5;MQ=22.26;MQ0=0;MQRankSum=0.195;QD=16.04;ReadPosRankSum=-4.072 GT:AD:DP:FT:GQ:PL 0/1:8,14:22:PASS:99:416,0,2147483647 0/1:18,13:31:PASS:99:353,0,2147483647 0/1:5,15:20:rd:99:2147483647,0,107
12 changes: 6 additions & 6 deletions adam-core/src/test/resources/sorted.vcf
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,9 +66,9 @@
##contig=<ID=2,length=249250621>
##contig=<ID=13,length=249250621>
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA12878 NA12891 NA12892
1 14397 . CTGT C . IndelQD AC=2;AF=0.333;AN=6;BaseQRankSum=1.8;ClippingRankSum=0.138;DP=69;FS=7.786;MLEAC=2;MLEAF=0.333;MQ=26.84;MQ0=0;MQRankSum=-1.906;QD=1.55;ReadPosRankSum=0.384 GT:AD:DP:FT:GQ:PL 0/1:16,4:20:rd:99:120,0,256 0/1:8,2:10:dp;rd:60:60,0,256 0/0:39,0:39:PASS:99:0,116,256
1 14522 . G A . VQSRTrancheSNP99.95to100.00 AC=2;AF=0.333;AN=6;BaseQRankSum=2.044;ClippingRankSum=-2.196;DP=48;FS=13.179;MLEAC=2;MLEAF=0.333;MQ=25.89;MQ0=0;MQRankSum=-0.063;QD=8.87;ReadPosRankSum=0.952;VQSLOD=-3.333;culprit=MQ GT:AD:DP:FT:GQ:PL 0/1:10,5:15:dp:99:99,0,256 0/1:2,5:7:dp;rd:34:128,0,34 0/0:26,0:26:PASS:78:0,78,256
1 63735 rs201888535 CCTA C . PASS AC=1;AF=0.167;AN=6;BaseQRankSum=1.138;ClippingRankSum=0.448;DB;DP=176;FS=13.597;MLEAC=1;MLEAF=0.167;MQ=31.06;MQ0=0;MQRankSum=0.636;QD=9.98;ReadPosRankSum=-1.18 GT:AD:DP:FT:GQ:PL 0/0:27,0:27:PASS:79:0,79,256 0/0:40,0:40:PASS:99:0,117,256 0/1:23,74:97:rd:99:256,0,256
2 19190 . GC G . PASS AC=3;AF=0.5;AN=6;BaseQRankSum=4.157;ClippingRankSum=3.666;DP=74;FS=37.037;MLEAC=3;MLEAF=0.5;MQ=22.26;MQ0=0;MQRankSum=0.195;QD=16.04;ReadPosRankSum=-4.072 GT:AD:DP:FT:GQ:PL 0/1:8,14:22:PASS:99:256,0,256 0/1:18,13:31:PASS:99:256,0,256 0/1:5,15:20:rd:99:256,0,107
13 752721 rs3131972 A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL 1/1:0,27:27:PASS:81:256,81,0 1/1:0,19:19:dp:57:256,57,0 1/1:0,22:22:PASS:66:256,66,0
13 752791 . A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL:SB 1/1:0,27:27:PASS:81:256,81,0:0,1,2,3 1/1:0,19:19:dp:57:256,57,0:4,5,6,7 1/1:0,22:22:PASS:66:256,66,0:2,3,4,5
1 14397 . CTGT C . IndelQD AC=2;AF=0.333;AN=6;BaseQRankSum=1.8;ClippingRankSum=0.138;DP=69;FS=7.786;MLEAC=2;MLEAF=0.333;MQ=26.84;MQ0=0;MQRankSum=-1.906;QD=1.55;ReadPosRankSum=0.384 GT:AD:DP:FT:GQ:PL 0/1:16,4:20:rd:99:120,0,2147483647 0/1:8,2:10:dp;rd:60:60,0,414 0/0:39,0:39:PASS:99:0,116,2147483647
1 14522 . G A . VQSRTrancheSNP99.95to100.00 AC=2;AF=0.333;AN=6;BaseQRankSum=2.044;ClippingRankSum=-2.196;DP=48;FS=13.179;MLEAC=2;MLEAF=0.333;MQ=25.89;MQ0=0;MQRankSum=-0.063;QD=8.87;ReadPosRankSum=0.952;VQSLOD=-3.333;culprit=MQ GT:AD:DP:FT:GQ:PL 0/1:10,5:15:dp:99:99,0,2147483647 0/1:2,5:7:dp;rd:34:128,0,34 0/0:26,0:26:PASS:78:0,78,2147483647
1 63735 rs201888535 CCTA C . PASS AC=1;AF=0.167;AN=6;BaseQRankSum=1.138;ClippingRankSum=0.448;DB;DP=176;FS=13.597;MLEAC=1;MLEAF=0.167;MQ=31.06;MQ0=0;MQRankSum=0.636;QD=9.98;ReadPosRankSum=-1.18 GT:AD:DP:FT:GQ:PL 0/0:27,0:27:PASS:79:0,79,2147483647 0/0:40,0:40:PASS:99:0,117,2147483647 0/1:23,74:97:rd:99:2147483647,0,2147483647
2 19190 . GC G . PASS AC=3;AF=0.5;AN=6;BaseQRankSum=4.157;ClippingRankSum=3.666;DP=74;FS=37.037;MLEAC=3;MLEAF=0.5;MQ=22.26;MQ0=0;MQRankSum=0.195;QD=16.04;ReadPosRankSum=-4.072 GT:AD:DP:FT:GQ:PL 0/1:8,14:22:PASS:99:416,0,2147483647 0/1:18,13:31:PASS:99:353,0,2147483647 0/1:5,15:20:rd:99:2147483647,0,107
13 752721 rs3131972 A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL 1/1:0,27:27:PASS:81:2147483647,81,0 1/1:0,19:19:dp:57:2147483647,57,0 1/1:0,22:22:PASS:66:2147483647,66,0
13 752791 . A G . PASS AC=6;AF=1.0;AN=6;DB;DP=69;FS=0.0;MLEAC=6;MLEAF=1.0;MQ=60.0;MQ0=0;POSITIVE_TRAIN_SITE;QD=31.67;VQSLOD=18.94;culprit=QD GT:AD:DP:FT:GQ:PL:SB 1/1:0,27:27:PASS:81:2147483647,81,0:0,1,2,3 1/1:0,19:19:dp:57:2147483647,57,0:4,5,6,7 1/1:0,22:22:PASS:66:2147483647,66,0:2,3,4,5
5 changes: 4 additions & 1 deletion adam-core/src/test/scala/org/bdgenomics/adam/converters/VariantContextConverterSuite.scala
View file
Original file line number Diff line number Diff line change
Expand Up @@ -344,13 +344,15 @@ class VariantContextConverterSuite extends ADAMFunSuite {
assert(adamGT1.getAlternateReadDepth === 2)
assert(adamGT1.getGenotypeLikelihoods
.map(f => f: scala.Float)
.map(_.toDouble)
.map(PhredUtils.logProbabilityToPhred)
.sameElements(List(59, 0, 256)))
.sameElements(List(59, 0, Int.MaxValue)))

assert(adamGT2.getAlleles.sameElements(List(GenotypeAllele.OTHER_ALT, GenotypeAllele.ALT)))
assert(adamGT2.getAlternateReadDepth === 3)
assert(adamGT2.getGenotypeLikelihoods
.map(f => f: scala.Float)
.map(_.toDouble)
.map(PhredUtils.logProbabilityToPhred)
.sameElements(List(59, 1, 102)))

Expand Down Expand Up @@ -391,6 +393,7 @@ class VariantContextConverterSuite extends ADAMFunSuite {
assert(adamGT.getGenotypeLikelihoods.isEmpty)
assert(adamGT.getNonReferenceLikelihoods
.map(f => f: scala.Float)
.map(_.toDouble)
.map(PhredUtils.logProbabilityToPhred)
.sameElements(List(0, 1, 2)))
}
Expand Down
37 changes: 37 additions & 0 deletions adam-core/src/test/scala/org/bdgenomics/adam/util/PhredUtilsSuite.scala
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,37 @@
/**
* Licensed to Big Data Genomics (BDG) under one
* or more contributor license agreements. See the NOTICE file
* distributed with this work for additional information
* regarding copyright ownership. The BDG licenses this file
* to you under the Apache License, Version 2.0 (the
* "License"); you may not use this file except in compliance
* with the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.bdgenomics.adam.util

import org.scalatest.FunSuite

class PhredUtilsSuite extends FunSuite {

test("convert low phred score to log and back") {
val logP = PhredUtils.phredToLogProbability(10)
val phred = PhredUtils.logProbabilityToPhred(logP)
assert(phred === 10)
}

test("convert high phred score to log and back") {
val logP = PhredUtils.phredToLogProbability(1000)
println(logP)
val phred = PhredUtils.logProbabilityToPhred(logP)
println(phred)
assert(phred === 1000)
}
}

0 comments on commit 412d092

Please sign in to comment.