From 0cdadffe14b1b635bb842547192add2c9de3c156 Mon Sep 17 00:00:00 2001 From: Ryan Poplin Date: Thu, 10 May 2012 13:09:03 -0400 Subject: [PATCH] Committing the best of the frantic pre-CSHL experiments: Better algorithm for partioning reads amongst the alleles they support. Require the read's original alignment to actually overlap the variant. QD uses the non-informative reads when calculating D. More HC-specific annotations for potential use in a statistical filtering strategy. Increasing the minimum kmer length in the assembly graphs. Misc minor bug fixes. --- .../sting/gatk/iterators/LocusIteratorByState.java | 8 ++++---- .../sting/gatk/walkers/annotator/QualByDepth.java | 3 +-- .../sting/utils/fragments/FragmentUtils.java | 3 +++ 3 files changed, 8 insertions(+), 6 deletions(-) diff --git a/public/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java b/public/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java index 8b9674353..c71373cd1 100755 --- a/public/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java +++ b/public/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java @@ -451,12 +451,12 @@ public class LocusIteratorByState extends LocusIterator { nextAlignmentContext = new AlignmentContext(loc, new ReadBackedExtendedEventPileupImpl(loc, fullExtendedEventPileup), hasBeenSampled); } else { // this is a regular event pileup (not extended) - GenomeLoc location = getLocation(); - Map fullPileup = new HashMap(); + final GenomeLoc location = getLocation(); + final Map fullPileup = new HashMap(); boolean hasBeenSampled = false; for (final String sample : samples) { - Iterator iterator = readStates.iterator(sample); - List pile = new ArrayList(readStates.size(sample)); + final Iterator iterator = readStates.iterator(sample); + final List pile = new ArrayList(readStates.size(sample)); hasBeenSampled |= location.getStart() <= readStates.getDownsamplingExtent(sample); size = 0; // number of elements in this sample's pileup diff --git a/public/java/src/org/broadinstitute/sting/gatk/walkers/annotator/QualByDepth.java b/public/java/src/org/broadinstitute/sting/gatk/walkers/annotator/QualByDepth.java index 51b834bd2..44a0498d4 100755 --- a/public/java/src/org/broadinstitute/sting/gatk/walkers/annotator/QualByDepth.java +++ b/public/java/src/org/broadinstitute/sting/gatk/walkers/annotator/QualByDepth.java @@ -86,8 +86,7 @@ public class QualByDepth extends InfoFieldAnnotation implements StandardAnnotati continue; for ( final Map.Entry> alleleBin : alleleBins.entrySet() ) { - if ( !alleleBin.getKey().equals(Allele.NO_CALL) ) - depth += alleleBin.getValue().size(); + depth += alleleBin.getValue().size(); } } diff --git a/public/java/src/org/broadinstitute/sting/utils/fragments/FragmentUtils.java b/public/java/src/org/broadinstitute/sting/utils/fragments/FragmentUtils.java index 858f7a2ae..fcab47618 100644 --- a/public/java/src/org/broadinstitute/sting/utils/fragments/FragmentUtils.java +++ b/public/java/src/org/broadinstitute/sting/utils/fragments/FragmentUtils.java @@ -166,6 +166,9 @@ public class FragmentUtils { if( firstReadQuals[iii] > MIN_QUAL_BAD_OVERLAP && secondReadQuals[iii-firstReadStop] > MIN_QUAL_BAD_OVERLAP && firstReadBases[iii] != secondReadBases[iii-firstReadStop] ) { return overlappingPair;// high qual bases don't match exactly, probably indel in only one of the fragments, so don't merge them } + if( firstReadQuals[iii] < MIN_QUAL_BAD_OVERLAP && secondReadQuals[iii-firstReadStop] < MIN_QUAL_BAD_OVERLAP ) { + return overlappingPair; // both reads have low qual bases in the overlap region so don't merge them because don't know what is going on + } bases[iii] = ( firstReadQuals[iii] > secondReadQuals[iii-firstReadStop] ? firstReadBases[iii] : secondReadBases[iii-firstReadStop] ); quals[iii] = ( firstReadQuals[iii] > secondReadQuals[iii-firstReadStop] ? firstReadQuals[iii] : secondReadQuals[iii-firstReadStop] ); }