diff --git a/java/src/org/broadinstitute/sting/oneoffprojects/walkers/association/modules/MismatchRate.java b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/association/modules/MismatchRate.java new file mode 100755 index 000000000..465c2ffe9 --- /dev/null +++ b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/association/modules/MismatchRate.java @@ -0,0 +1,91 @@ +package org.broadinstitute.sting.oneoffprojects.walkers.association.modules; + +import org.broadinstitute.sting.gatk.datasources.sample.Sample; +import org.broadinstitute.sting.oneoffprojects.walkers.association.MapExtender; +import org.broadinstitute.sting.oneoffprojects.walkers.association.RegionalAssociationWalker; +import org.broadinstitute.sting.oneoffprojects.walkers.association.statistics.casecontrol.UStatistic; +import org.broadinstitute.sting.utils.BaseUtils; +import org.broadinstitute.sting.utils.MathUtils; +import org.broadinstitute.sting.utils.collections.Pair; +import org.broadinstitute.sting.utils.pileup.ReadBackedPileup; + +import java.util.*; + +/** + * Created by IntelliJ IDEA. + * User: chartl + * Date: 3/8/11 + * Time: 2:32 PM + * To change this template use File | Settings | File Templates. + */ +public class MismatchRate extends UStatistic { + + private Map sampleStats = new HashMap(); + private int currentRefBase = 0; + + public void init(RegionalAssociationWalker walker) { + Set samples = walker.getSamples(); + for ( Sample s : samples ) { + if ( s.hasProperty("mismatch_rate.mean") && s.hasProperty("mismatch_rate.std") ) { + double mn = (Double) s.getProperty("mismatch_rate.mean"); + double std = (Double) s.getProperty("mismatch_rate.std"); + sampleStats.put(s,new Pair(mn,std)); + } else { + sampleStats.put(s,new MathUtils.RunningAverage()); + } + } + } + + @Override + public Map mapLocus(MapExtender extender) { + currentRefBase = BaseUtils.simpleBaseToBaseIndex(extender.getReferenceContext().getBase()); + Map pileups = extender.getReadFilteredPileup(); + Map maps = new HashMap(pileups.size()); + for ( Map.Entry samPileup : pileups.entrySet() ) { + maps.put(samPileup.getKey(),map(samPileup.getKey(),samPileup.getValue())); + } + + return maps; + } + + public Collection map(Sample sample, ReadBackedPileup pileup) { + Object stats = sampleStats.get(sample); + double mn; + double std; + if ( stats instanceof Pair ) { + mn = ((Pair)stats).first; + std = ((Pair)stats).second; + } else { + MathUtils.RunningAverage ra = (MathUtils.RunningAverage) stats; + mn = ra.mean(); + std = ra.stddev(); + if ( std <= 0.0 ) { + std = 1.0; + } + ra.add(pileup.size()); + } + + return Arrays.asList((Number) ((calcMMR(pileup) - mn) / std)); + } + + public double calcMMR(ReadBackedPileup rbp) { + int[] counts = rbp.getBaseCounts(); + int total = 0; + int nonref = 0; + for ( int base : new int[]{0,1,2,3} ) { + total += counts[base]; + if ( base != currentRefBase ) { + nonref += counts[base]; + } + } + + return ((double)nonref)/total; + } + + // note: this is to satisfy the interface, and is never called due to override + public Collection map(ReadBackedPileup pileup) { return null; } + + public int getWindowSize() { return 100; } + public int slideByValue() { return 10; } + public boolean usePreviouslySeenReads() { return true; } +} diff --git a/java/src/org/broadinstitute/sting/oneoffprojects/walkers/association/modules/ReferenceMismatches.java b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/association/modules/ReferenceMismatches.java new file mode 100755 index 000000000..424f6bb8c --- /dev/null +++ b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/association/modules/ReferenceMismatches.java @@ -0,0 +1,48 @@ +package org.broadinstitute.sting.oneoffprojects.walkers.association.modules; + +import org.broadinstitute.sting.gatk.datasources.sample.Sample; +import org.broadinstitute.sting.oneoffprojects.walkers.association.MapExtender; +import org.broadinstitute.sting.oneoffprojects.walkers.association.statistics.casecontrol.ZStatistic; +import org.broadinstitute.sting.utils.BaseUtils; +import org.broadinstitute.sting.utils.collections.Pair; +import org.broadinstitute.sting.utils.pileup.ReadBackedPileup; + +import java.util.Map; + +/** + * Created by IntelliJ IDEA. + * User: chartl + * Date: 3/8/11 + * Time: 1:17 PM + * To change this template use File | Settings | File Templates. + */ +public class ReferenceMismatches extends ZStatistic { + + final static int[] BASE_INDEX = {0,1,2,3}; + + int currentRefBase = 0; + + @Override + public Map mapLocus(MapExtender extender) { + currentRefBase = BaseUtils.simpleBaseToBaseIndex(extender.getReferenceContext().getBase()); + return super.mapLocus(extender); + } + + public Pair map(ReadBackedPileup rbp) { + int[] counts = rbp.getBaseCounts(); + int total = 0; + int nonref = 0; + for ( int base : BASE_INDEX ) { + total += counts[base]; + if ( base != currentRefBase ) { + nonref += counts[base]; + } + } + + return new Pair(nonref,total); + } + + public int getWindowSize() { return 100; } + public int slideByValue() { return 25; } + public boolean usePreviouslySeenReads() { return true; } +}