@@ -121,7 +121,7 @@ object DownsampleVcf extends LazyLogging {
121
121
* @param epsilon the error rate for genotyping
122
122
* @return a new `Likelihood` that has the likelihoods of AA, AB, and BB
123
123
*/
124
- def biallelic (alleleDepthA : Int , alleleDepthB : Int , epsilon : Double = 0.01 ): Likelihoods = {
124
+ def biallelic (alleleDepthA : Int , alleleDepthB : Int , epsilon : Double = 0.01 ): IndexedSeq [ Double ] = {
125
125
val aGivenAA = log10(1 - epsilon)
126
126
val aGivenBB = log10(epsilon)
127
127
val aGivenAB = log10((1 - epsilon) / 2 )
@@ -130,7 +130,7 @@ object DownsampleVcf extends LazyLogging {
130
130
val rawGlBB = ((alleleDepthA * aGivenBB) + (alleleDepthB * aGivenAA))
131
131
val rawGlAB = ((alleleDepthA + alleleDepthB) * aGivenAB)
132
132
133
- Likelihoods ( 2 , IndexedSeq (rawGlAA, rawGlAB, rawGlBB) )
133
+ IndexedSeq (rawGlAA, rawGlAB, rawGlBB)
134
134
}
135
135
136
136
/** Computes the likelihoods for each possible genotype given a sequence of read depths for any
@@ -140,30 +140,28 @@ object DownsampleVcf extends LazyLogging {
140
140
* @return a new `Likelihood` that has the log likelihoods of all possible genotypes in the
141
141
* order specified in VFC spec for the GL/PL tags.
142
142
*/
143
- def generalized (alleleDepths : IndexedSeq [Int ], epsilon : Double = 0.01 ): Likelihoods = {
143
+ def generalized (alleleDepths : IndexedSeq [Int ], epsilon : Double = 0.01 ): IndexedSeq [ Double ] = {
144
144
val numAlleles = alleleDepths.length
145
145
// probabilities associated with each possible genotype for a pair of alleles
146
146
val logProbs : Array [Double ] = Array (
147
147
math.log10(epsilon),
148
148
math.log10((1 - epsilon) / 2 ),
149
149
math.log10(1 - epsilon)
150
150
)
151
- // raw genotype log-likelihoods
152
- Likelihoods (
153
- numAlleles,
154
- (0 until numAlleles).flatMap(b =>
155
- (0 to b).map(a =>
156
- (0 until numAlleles).map(allele =>
157
- logProbs(Array (a, b).count(_ == allele)) * alleleDepths(allele)
158
- ).sum
159
- )
151
+ // compute genotype log-likelihoods
152
+ (0 until numAlleles).flatMap(b =>
153
+ (0 to b).map(a =>
154
+ (0 until numAlleles).map(allele =>
155
+ logProbs(Array (a, b).count(_ == allele)) * alleleDepths(allele)
156
+ ).sum
160
157
)
161
158
)
162
159
}
163
160
164
161
def apply (alleleDepths : IndexedSeq [Int ], epsilon : Double = 0.01 ): Likelihoods = {
165
- require(alleleDepths.length >= 2 , " at least two alleles are required to calculate genotype likelihoods" )
166
- generalized(alleleDepths, epsilon)
162
+ val numAlleles = alleleDepths.length
163
+ require(numAlleles >= 2 , " at least two alleles are required to calculate genotype likelihoods" )
164
+ Likelihoods (numAlleles, generalized(alleleDepths, epsilon))
167
165
}
168
166
}
169
167
0 commit comments