Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Walk by reads instead of by reference #2

Open
wants to merge 5 commits into
base: master
Choose a base branch
from
Open

Walk by reads instead of by reference #2

Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
a3b9806
walk by reads
Jun 27, 2013
980f7cc
fixed sorting error (Comparison method violates its general contract)
Jun 27, 2013
2cf2153
limit whitespace diffs
Jun 27, 2013
1c40251
added back walk by reads (removed accidentally when trying to limit
Jun 27, 2013
8705402
fixed pileup sort order (put highest quality 1st)
Jun 27, 2013
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
72 changes: 36 additions & 36 deletions src/org/broadinstitute/cga/tools/gatk/walkers/cancer/mutect/MuTect.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@

@PartitionBy(PartitionType.LOCUS)
@Reference(window=@Window(start=-1* MuTect.REFERENCE_HALF_WINDOW_LENGTH,stop= MuTect.REFERENCE_HALF_WINDOW_LENGTH))
@By(DataSource.REFERENCE)
@By(DataSource.READS)
public class MuTect extends LocusWalker<Integer, Integer> implements TreeReducible<Integer> {
public static final int REFERENCE_HALF_WINDOW_LENGTH = 150;
public static final String BAM_TAG_TUMOR = "tumor";
Expand Down Expand Up @@ -139,22 +139,15 @@ public class MuTect extends LocusWalker<Integer, Integer> implements TreeReducib
private NormalPowerCalculator normalDbSNPSitePowerCalculator;
private TumorPowerCalculator strandArtifactPowerCalculator;

private static class PileupComparatorByAltRefQual implements Comparator<PileupElement> {
private byte alt;
private static class PileupComparatorByQual implements Comparator<PileupElement> {


private PileupComparatorByAltRefQual(byte alt) {
this.alt = alt;
}

public int compare(PileupElement o1, PileupElement o2) {
// if the bases are the same, the higher quality score comes first
if (o1.getBase() == o2.getBase()) {
if (o1.getQual() == o2.getQual()) { return 0; }
return (o1.getQual() > o2.getQual())?-1:1;
} else {
return (o1.getBase() == alt)?-1:1;
}




public int compare(PileupElement o1, PileupElement o2) {
return o2.getQual() - o1.getQual();
}
}

Expand Down Expand Up @@ -481,27 +474,34 @@ public Integer map(final RefMetaDataTracker tracker, final ReferenceContext ref,
VariableAllelicRatioGenotypeLikelihoods contaminantLikelihoods
= new VariableAllelicRatioGenotypeLikelihoods(upRef, contaminantF);

List<PileupElement> peList = new ArrayList<PileupElement>(tumorReadPile.finalPileup.depthOfCoverage());
for(PileupElement pe : tumorReadPile.finalPileup) {
peList.add(pe);
}

Collections.sort(peList, new PileupComparatorByAltRefQual((byte)altAllele));
int readsToKeep = (int) (peList.size() * contaminantAlternateFraction);

for(PileupElement pe : peList) {
byte base = pe.getBase();
if (pe.getBase() == altAllele) {
// if we've retained all we need, then turn the remainder of alts to ref
if (readsToKeep == 0) {
base = (byte) upRef;
} else {
readsToKeep--;
}
}

contaminantLikelihoods.add(base, pe.getQual());
}
List<PileupElement> refList = new ArrayList<PileupElement>();
List<PileupElement> altList = new ArrayList<PileupElement>();

for (PileupElement element : tumorReadPile.finalPileup) {
if (element.getBase() == altAllele) {
altList.add(element);
} else {
refList.add(element);
}
}
Collections.sort(altList, new PileupComparatorByQual());

int readsToKeep = (int) ((altList.size() + refList.size()) * contaminantAlternateFraction);

for (PileupElement pe : altList) {
byte base = pe.getBase();
// if we've retained all we need, then turn the
// remainder of alts to ref
if (readsToKeep == 0) {
base = (byte) upRef;
} else {
readsToKeep--;
}
contaminantLikelihoods.add(base, pe.getQual());
}
for (PileupElement pe : refList) {
contaminantLikelihoods.add((byte) upRef, pe.getQual());
}
double[] refHetHom = LocusReadPile.extractRefHetHom(contaminantLikelihoods, upRef, altAllele);
double contaminantLod = refHetHom[1] - refHetHom[0];
candidate.setContaminantLod(contaminantLod);
Expand Down