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

Simple bug fix to allow multiple records at same site when genotyping given alleles. Takes only the first record (respecting filters, SNP type, etc), and issues a warning if there is more than one valid record at a site 

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@5789 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
depristo 2011-05-11 14:17:14 +00:00
parent 16db86e6cb
commit 9a1d0d7076
1 changed files with 18 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

26
java/src/org/broadinstitute/sting/gatk/walkers/genotyper/SNPGenotypeLikelihoodsCalculationModel.java
View File

@ -73,18 +73,28 @@ public class SNPGenotypeLikelihoodsCalculationModel extends GenotypeLikelihoodsC
if ( alternateAlleleToUse != null ) {
bestAlternateAllele = alternateAlleleToUse.getBases()[0];
} else if ( useAlleleFromVCF ) {
final VariantContext vcInput = tracker.getVariantContext(ref, "alleles", null, ref.getLocus(), true);
if ( vcInput == null || vcInput.isFiltered() || !vcInput.isSNP() )
return null;
if ( !vcInput.isSNP() ) {
logger.info("Record at position " + ref.getLocus() + " is not a SNP; skipping...");
return null;
VariantContext vc = null;
// search for usable record
for( final VariantContext vc_input : tracker.getVariantContexts(ref, "alleles", null, ref.getLocus(), true, false) ) {
if ( vc_input != null && ! vc_input.isFiltered() && vc_input.isSNP() ) {
if ( vc == null ) {
vc = vc_input;
} else {
logger.warn("Multiple valid VCF records detected at site " + ref.getLocus() + ", only considering alleles from first record only");
}
}
}
if ( !vcInput.isBiallelic() ) {
// ignore places where we don't have a variant
if ( vc == null )
return null;
if ( !vc.isBiallelic() ) {
// for multi-allelic sites go back to the reads and find the most likely alternate allele
initializeBestAlternateAllele(refBase, contexts);
} else {
bestAlternateAllele = vcInput.getAlternateAllele(0).getBases()[0];
bestAlternateAllele = vc.getAlternateAllele(0).getBases()[0];
}
} else {
initializeBestAlternateAllele(refBase, contexts);