Fixes for SLOD: 1) didn't work properly for multi-allelics (randomly chose an allele, possibly one that wasn't genotyped in the full context); 2) in cases when there were more alt alleles than the max allowed and the user is calculating SB, we would recompute the best alt alleles(s); 3) for some reason, we were recomputing the LOD for the full context when we'd already done that. Given that this passes integration tests on my end, this should be the last commit before the release.
This commit is contained in:
parent
1ee46e5c06
commit
b4749757f8
|
|
@ -38,6 +38,7 @@ import org.broadinstitute.sting.utils.pileup.ReadBackedPileup;
|
||||||
import org.broadinstitute.sting.utils.variantcontext.Allele;
|
import org.broadinstitute.sting.utils.variantcontext.Allele;
|
||||||
import org.broadinstitute.sting.utils.variantcontext.VariantContext;
|
import org.broadinstitute.sting.utils.variantcontext.VariantContext;
|
||||||
|
|
||||||
|
import java.util.List;
|
||||||
import java.util.Map;
|
import java.util.Map;
|
||||||
|
|
||||||
|
|
||||||
|
|
@ -81,19 +82,19 @@ public abstract class GenotypeLikelihoodsCalculationModel implements Cloneable {
|
||||||
* @param contexts stratified alignment contexts
|
* @param contexts stratified alignment contexts
|
||||||
* @param contextType stratified context type
|
* @param contextType stratified context type
|
||||||
* @param priors priors to use for GLs
|
* @param priors priors to use for GLs
|
||||||
* @param alternateAlleleToUse the alternate allele to use, null if not set
|
* @param alternateAllelesToUse the alternate allele to use, null if not set
|
||||||
* @param useBAQedPileup should we use the BAQed pileup or the raw one?
|
* @param useBAQedPileup should we use the BAQed pileup or the raw one?
|
||||||
* @param locParser Genome Loc Parser
|
* @param locParser Genome Loc Parser
|
||||||
* @return variant context where genotypes are no-called but with GLs
|
* @return variant context where genotypes are no-called but with GLs
|
||||||
*/
|
*/
|
||||||
public abstract VariantContext getLikelihoods(RefMetaDataTracker tracker,
|
public abstract VariantContext getLikelihoods(final RefMetaDataTracker tracker,
|
||||||
ReferenceContext ref,
|
final ReferenceContext ref,
|
||||||
Map<String, AlignmentContext> contexts,
|
final Map<String, AlignmentContext> contexts,
|
||||||
AlignmentContextUtils.ReadOrientation contextType,
|
final AlignmentContextUtils.ReadOrientation contextType,
|
||||||
GenotypePriors priors,
|
final GenotypePriors priors,
|
||||||
Allele alternateAlleleToUse,
|
final List<Allele> alternateAllelesToUse,
|
||||||
boolean useBAQedPileup,
|
final boolean useBAQedPileup,
|
||||||
GenomeLocParser locParser);
|
final GenomeLocParser locParser);
|
||||||
|
|
||||||
|
|
||||||
protected int getFilteredDepth(ReadBackedPileup pileup) {
|
protected int getFilteredDepth(ReadBackedPileup pileup) {
|
||||||
|
|
|
||||||
|
|
@ -284,13 +284,14 @@ public class IndelGenotypeLikelihoodsCalculationModel extends GenotypeLikelihood
|
||||||
|
|
||||||
private final static EnumSet<VariantContext.Type> allowableTypes = EnumSet.of(VariantContext.Type.INDEL, VariantContext.Type.MIXED);
|
private final static EnumSet<VariantContext.Type> allowableTypes = EnumSet.of(VariantContext.Type.INDEL, VariantContext.Type.MIXED);
|
||||||
|
|
||||||
public VariantContext getLikelihoods(RefMetaDataTracker tracker,
|
public VariantContext getLikelihoods(final RefMetaDataTracker tracker,
|
||||||
ReferenceContext ref,
|
final ReferenceContext ref,
|
||||||
Map<String, AlignmentContext> contexts,
|
final Map<String, AlignmentContext> contexts,
|
||||||
AlignmentContextUtils.ReadOrientation contextType,
|
final AlignmentContextUtils.ReadOrientation contextType,
|
||||||
GenotypePriors priors,
|
final GenotypePriors priors,
|
||||||
Allele alternateAlleleToUse,
|
final List<Allele> alternateAllelesToUse,
|
||||||
boolean useBAQedPileup, GenomeLocParser locParser) {
|
final boolean useBAQedPileup,
|
||||||
|
final GenomeLocParser locParser) {
|
||||||
|
|
||||||
if (tracker == null)
|
if (tracker == null)
|
||||||
return null;
|
return null;
|
||||||
|
|
|
||||||
|
|
@ -64,7 +64,7 @@ public class SNPGenotypeLikelihoodsCalculationModel extends GenotypeLikelihoodsC
|
||||||
final Map<String, AlignmentContext> contexts,
|
final Map<String, AlignmentContext> contexts,
|
||||||
final AlignmentContextUtils.ReadOrientation contextType,
|
final AlignmentContextUtils.ReadOrientation contextType,
|
||||||
final GenotypePriors priors,
|
final GenotypePriors priors,
|
||||||
final Allele alternateAlleleToUse,
|
final List<Allele> alternateAllelesToUse,
|
||||||
final boolean useBAQedPileup,
|
final boolean useBAQedPileup,
|
||||||
final GenomeLocParser locParser) {
|
final GenomeLocParser locParser) {
|
||||||
|
|
||||||
|
|
@ -95,8 +95,8 @@ public class SNPGenotypeLikelihoodsCalculationModel extends GenotypeLikelihoodsC
|
||||||
}
|
}
|
||||||
|
|
||||||
// find the alternate allele(s) that we should be using
|
// find the alternate allele(s) that we should be using
|
||||||
if ( alternateAlleleToUse != null ) {
|
if ( alternateAllelesToUse != null ) {
|
||||||
alleles.add(alternateAlleleToUse);
|
alleles.addAll(alternateAllelesToUse);
|
||||||
} else if ( useAlleleFromVCF ) {
|
} else if ( useAlleleFromVCF ) {
|
||||||
final VariantContext vc = UnifiedGenotyperEngine.getVCFromAllelesRod(tracker, ref, ref.getLocus(), true, logger, UAC.alleles);
|
final VariantContext vc = UnifiedGenotyperEngine.getVCFromAllelesRod(tracker, ref, ref.getLocus(), true, logger, UAC.alleles);
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -237,14 +237,14 @@ public class UnifiedGenotyperEngine {
|
||||||
// ---------------------------------------------------------------------------------------------------------
|
// ---------------------------------------------------------------------------------------------------------
|
||||||
|
|
||||||
// private method called by both UnifiedGenotyper and UGCalcLikelihoods entry points into the engine
|
// private method called by both UnifiedGenotyper and UGCalcLikelihoods entry points into the engine
|
||||||
private VariantContext calculateLikelihoods(RefMetaDataTracker tracker, ReferenceContext refContext, Map<String, AlignmentContext> stratifiedContexts, AlignmentContextUtils.ReadOrientation type, Allele alternateAlleleToUse, boolean useBAQedPileup, final GenotypeLikelihoodsCalculationModel.Model model) {
|
private VariantContext calculateLikelihoods(RefMetaDataTracker tracker, ReferenceContext refContext, Map<String, AlignmentContext> stratifiedContexts, AlignmentContextUtils.ReadOrientation type, List<Allele> alternateAllelesToUse, boolean useBAQedPileup, final GenotypeLikelihoodsCalculationModel.Model model) {
|
||||||
|
|
||||||
// initialize the data for this thread if that hasn't been done yet
|
// initialize the data for this thread if that hasn't been done yet
|
||||||
if ( glcm.get() == null ) {
|
if ( glcm.get() == null ) {
|
||||||
glcm.set(getGenotypeLikelihoodsCalculationObject(logger, UAC));
|
glcm.set(getGenotypeLikelihoodsCalculationObject(logger, UAC));
|
||||||
}
|
}
|
||||||
|
|
||||||
return glcm.get().get(model).getLikelihoods(tracker, refContext, stratifiedContexts, type, getGenotypePriors(model), alternateAlleleToUse, useBAQedPileup && BAQEnabledOnCMDLine, genomeLocParser);
|
return glcm.get().get(model).getLikelihoods(tracker, refContext, stratifiedContexts, type, getGenotypePriors(model), alternateAllelesToUse, useBAQedPileup && BAQEnabledOnCMDLine, genomeLocParser);
|
||||||
}
|
}
|
||||||
|
|
||||||
private VariantCallContext generateEmptyContext(RefMetaDataTracker tracker, ReferenceContext ref, Map<String, AlignmentContext> stratifiedContexts, AlignmentContext rawContext) {
|
private VariantCallContext generateEmptyContext(RefMetaDataTracker tracker, ReferenceContext ref, Map<String, AlignmentContext> stratifiedContexts, AlignmentContext rawContext) {
|
||||||
|
|
@ -398,16 +398,14 @@ public class UnifiedGenotyperEngine {
|
||||||
//final boolean DEBUG_SLOD = false;
|
//final boolean DEBUG_SLOD = false;
|
||||||
|
|
||||||
// the overall lod
|
// the overall lod
|
||||||
VariantContext vcOverall = calculateLikelihoods(tracker, refContext, stratifiedContexts, AlignmentContextUtils.ReadOrientation.COMPLETE, vc.getAlternateAllele(0), false, model);
|
|
||||||
clearAFarray(AFresult.log10AlleleFrequencyLikelihoods);
|
|
||||||
clearAFarray(AFresult.log10AlleleFrequencyPosteriors);
|
|
||||||
afcm.get().getLog10PNonRef(vcOverall, getAlleleFrequencyPriors(model), AFresult);
|
|
||||||
//double overallLog10PofNull = AFresult.log10AlleleFrequencyPosteriors[0];
|
//double overallLog10PofNull = AFresult.log10AlleleFrequencyPosteriors[0];
|
||||||
double overallLog10PofF = MathUtils.log10sumLog10(AFresult.log10AlleleFrequencyPosteriors[0], 0);
|
double overallLog10PofF = MathUtils.log10sumLog10(AFresult.log10AlleleFrequencyPosteriors[0], 0);
|
||||||
//if ( DEBUG_SLOD ) System.out.println("overallLog10PofF=" + overallLog10PofF);
|
//if ( DEBUG_SLOD ) System.out.println("overallLog10PofF=" + overallLog10PofF);
|
||||||
|
|
||||||
|
List<Allele> alternateAllelesToUse = builder.make().getAlternateAlleles();
|
||||||
|
|
||||||
// the forward lod
|
// the forward lod
|
||||||
VariantContext vcForward = calculateLikelihoods(tracker, refContext, stratifiedContexts, AlignmentContextUtils.ReadOrientation.FORWARD, vc.getAlternateAllele(0), false, model);
|
VariantContext vcForward = calculateLikelihoods(tracker, refContext, stratifiedContexts, AlignmentContextUtils.ReadOrientation.FORWARD, alternateAllelesToUse, false, model);
|
||||||
clearAFarray(AFresult.log10AlleleFrequencyLikelihoods);
|
clearAFarray(AFresult.log10AlleleFrequencyLikelihoods);
|
||||||
clearAFarray(AFresult.log10AlleleFrequencyPosteriors);
|
clearAFarray(AFresult.log10AlleleFrequencyPosteriors);
|
||||||
afcm.get().getLog10PNonRef(vcForward, getAlleleFrequencyPriors(model), AFresult);
|
afcm.get().getLog10PNonRef(vcForward, getAlleleFrequencyPriors(model), AFresult);
|
||||||
|
|
@ -417,7 +415,7 @@ public class UnifiedGenotyperEngine {
|
||||||
//if ( DEBUG_SLOD ) System.out.println("forwardLog10PofNull=" + forwardLog10PofNull + ", forwardLog10PofF=" + forwardLog10PofF);
|
//if ( DEBUG_SLOD ) System.out.println("forwardLog10PofNull=" + forwardLog10PofNull + ", forwardLog10PofF=" + forwardLog10PofF);
|
||||||
|
|
||||||
// the reverse lod
|
// the reverse lod
|
||||||
VariantContext vcReverse = calculateLikelihoods(tracker, refContext, stratifiedContexts, AlignmentContextUtils.ReadOrientation.REVERSE, vc.getAlternateAllele(0), false, model);
|
VariantContext vcReverse = calculateLikelihoods(tracker, refContext, stratifiedContexts, AlignmentContextUtils.ReadOrientation.REVERSE, alternateAllelesToUse, false, model);
|
||||||
clearAFarray(AFresult.log10AlleleFrequencyLikelihoods);
|
clearAFarray(AFresult.log10AlleleFrequencyLikelihoods);
|
||||||
clearAFarray(AFresult.log10AlleleFrequencyPosteriors);
|
clearAFarray(AFresult.log10AlleleFrequencyPosteriors);
|
||||||
afcm.get().getLog10PNonRef(vcReverse, getAlleleFrequencyPriors(model), AFresult);
|
afcm.get().getLog10PNonRef(vcReverse, getAlleleFrequencyPriors(model), AFresult);
|
||||||
|
|
|
||||||
Loading…
Reference in New Issue