zzh
/
gatk-3.8
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Adding Convey HC-1 HMM acceleration

This commit is contained in:
sathibault 2013-05-08 11:01:20 -05:00
parent 803f666fd5
commit d79b5f0931
3 changed files with 103 additions and 10 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
protected/java/src/org/broadinstitute/sting/gatk/walkers/haplotypecaller/LikelihoodCalculationEngine.java
View File

@ -90,7 +90,10 @@ public class LikelihoodCalculationEngine {
pairHMM = new Log10PairHMM(false);
break;
case LOGLESS_CACHING:
pairHMM = new LoglessPairHMM();
if (CnyPairHMM.isAvailable())
pairHMM = new CnyPairHMM();
else
pairHMM = new LoglessPairHMM();
break;
default:
throw new UserException.BadArgumentValue("pairHMM", "Specified pairHMM implementation is unrecognized or incompatible with the HaplotypeCaller. Acceptable options are ORIGINAL, EXACT, CACHING, and LOGLESS_CACHING.");
@ -151,6 +154,8 @@ public class LikelihoodCalculationEngine {
private PerReadAlleleLikelihoodMap computeReadLikelihoods( final List<Haplotype> haplotypes, final List<GATKSAMRecord> reads) {
// first, a little set up to get copies of the Haplotypes that are Alleles (more efficient than creating them each time)
final BatchPairHMM batchPairHMM = (pairHMM instanceof BatchPairHMM) ? (BatchPairHMM)pairHMM : null;
final Vector<GATKSAMRecord> batchedReads = new Vector<GATKSAMRecord>(reads.size());
final int numHaplotypes = haplotypes.size();
final Map<Haplotype, Allele> alleleVersions = new HashMap<Haplotype, Allele>(numHaplotypes);
for ( final Haplotype haplotype : haplotypes ) {
@ -177,16 +182,29 @@ public class LikelihoodCalculationEngine {
readQuals[kkk] = ( readQuals[kkk] < (byte) 18 ? QualityUtils.MIN_USABLE_Q_SCORE : readQuals[kkk] );
}
for( int jjj = 0; jjj < numHaplotypes; jjj++ ) {
final Haplotype haplotype = haplotypes.get(jjj);
final boolean isFirstHaplotype = jjj == 0;
final double log10l = pairHMM.computeReadLikelihoodGivenHaplotypeLog10(haplotype.getBases(),
read.getReadBases(), readQuals, readInsQuals, readDelQuals, overallGCP, isFirstHaplotype);
perReadAlleleLikelihoodMap.add(read, alleleVersions.get(haplotype), log10l);
}
if ( batchPairHMM != null ) {
batchPairHMM.batchAdd(haplotypes, read.getReadBases(), readQuals, readInsQuals, readDelQuals, overallGCP);
batchedReads.add(read);
} else {
for( int jjj = 0; jjj < numHaplotypes; jjj++ ) {
final Haplotype haplotype = haplotypes.get(jjj);
final boolean isFirstHaplotype = jjj == 0;
final double log10l = pairHMM.computeReadLikelihoodGivenHaplotypeLog10(haplotype.getBases(),
read.getReadBases(), readQuals, readInsQuals, readDelQuals, overallGCP, isFirstHaplotype);
perReadAlleleLikelihoodMap.add(read, alleleVersions.get(haplotype), log10l);
}
}
}
if ( batchPairHMM != null ) {
for( final GATKSAMRecord read : batchedReads ) {
final double[] likelihoods = batchPairHMM.batchResult();
for( int jjj = 0; jjj < numHaplotypes; jjj++ ) {
final Haplotype haplotype = haplotypes.get(jjj);
perReadAlleleLikelihoodMap.add(read, alleleVersions.get(haplotype), likelihoods[jjj]);
}
}
}
return perReadAlleleLikelihoodMap;
}

59
protected/java/src/org/broadinstitute/sting/utils/pairhmm/CnyPairHMM.java 100644
View File

@ -0,0 +1,59 @@
package org.broadinstitute.sting.utils.pairhmm;
import java.util.*;
import org.broadinstitute.sting.utils.haplotype.Haplotype;
public final class CnyPairHMM extends PairHMM implements BatchPairHMM {
private static class HmmInput {
public List<Haplotype> haplotypes;
public byte[] readBases;
public byte[] readQuals;
public byte[] insertionGOP;
public byte[] deletionGOP;
public byte[] overallGCP;
};
private static boolean loaded = false;
private List<HmmInput> pending = new LinkedList<HmmInput>();
static public boolean isAvailable() {
return true;
}
public void batchAdd(final List<Haplotype> haplotypes,
final byte[] readBases,
final byte[] readQuals,
final byte[] insertionGOP,
final byte[] deletionGOP,
final byte[] overallGCP) {
HmmInput test=new HmmInput();
test.haplotypes=haplotypes;
test.readBases=readBases;
test.readQuals=readQuals;
test.insertionGOP=insertionGOP;
test.deletionGOP=deletionGOP;
test.overallGCP=overallGCP;
pending.add(test);
}
public double[] batchResult() {
HmmInput test=pending.remove(0);
double[] results=new double[test.haplotypes.size()];
for (int i=0; i<results.length; i++) {
results[i]=0.0;
}
return results;
}
protected double subComputeReadLikelihoodGivenHaplotypeLog10( final byte[] haplotypeBases,
final byte[] readBases,
final byte[] readQuals,
final byte[] insertionGOP,
final byte[] deletionGOP,
final byte[] overallGCP,
final int hapStartIndex,
final boolean recacheReadValues ) {
return 0.0;
}
}

16
public/java/src/org/broadinstitute/sting/utils/pairhmm/BatchPairHMM.java 100644
View File

@ -0,0 +1,16 @@
package org.broadinstitute.sting.utils.pairhmm;
import java.util.List;
import org.broadinstitute.sting.utils.haplotype.Haplotype;
public interface BatchPairHMM {
public void batchAdd(final List<Haplotype> haplotypes,
final byte[] readBases,
final byte[] readQuals,
final byte[] insertionGOP,
final byte[] deletionGOP,
final byte[] overallGCP);
public double[] batchResult();
}