Skip to content

Commit fc248df

Browse files
authored
Prioritize het calls when merging clustered SVs (#9058)
1 parent cf88d08 commit fc248df

File tree

3 files changed

+174
-48
lines changed

3 files changed

+174
-48
lines changed

src/main/java/org/broadinstitute/hellbender/tools/sv/cluster/CanonicalSVCollapser.java

Lines changed: 28 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -122,30 +122,35 @@ public enum FlagFieldLogic {
122122
/**
123123
* Comparators used for picking the representative genotype for a given sample
124124
*/
125+
// Priotize non-ref over ref
125126
final Comparator<Genotype> genotypeIsNonRefComparator = (o1, o2) -> {
126127
final long count1 = Math.min(1, o1.getAlleles().stream().filter(Allele::isNonReference).filter(Allele::isCalled).count());
127128
final long count2 = Math.min(1, o2.getAlleles().stream().filter(Allele::isNonReference).filter(Allele::isCalled).count());
128129
return Long.compare(count1, count2);
129130
};
130131

132+
// Priotize fewer ALT alleles over more. When applied after non-ref comparator, hom-ref genotypes will not be encountered.
131133
final Comparator<Genotype> genotypeNonRefCountComparator = (o1, o2) -> {
132134
final long count1 = o1.getAlleles().stream().filter(Allele::isNonReference).filter(Allele::isCalled).count();
133135
final long count2 = o2.getAlleles().stream().filter(Allele::isNonReference).filter(Allele::isCalled).count();
134-
return Long.compare(count1, count2);
136+
return Long.compare(count2, count1);
135137
};
136138

139+
// Priotize called genotypes
137140
final Comparator<Genotype> genotypeCalledComparator = (o1, o2) -> {
138141
final long count1 = o1.getAlleles().stream().filter(Allele::isCalled).count();
139142
final long count2 = o2.getAlleles().stream().filter(Allele::isCalled).count();
140143
return Long.compare(count1, count2);
141144
};
142145

146+
// Priotize higher quality
143147
final Comparator<Genotype> genotypeQualityComparator = (o1, o2) -> {
144148
final int quality1 = VariantContextGetters.getAttributeAsInt(o1, VCFConstants.GENOTYPE_QUALITY_KEY, 0);
145149
final int quality2 = VariantContextGetters.getAttributeAsInt(o2, VCFConstants.GENOTYPE_QUALITY_KEY, 0);
146150
return Integer.compare(quality1, quality2);
147151
};
148152

153+
// Priotize higher depth genotyping quality
149154
final Comparator<Genotype> genotypeCopyNumberQualityComparator = new Comparator<Genotype>() {
150155
@Override
151156
public int compare(Genotype o1, Genotype o2) {
@@ -155,14 +160,33 @@ public int compare(Genotype o1, Genotype o2) {
155160
}
156161
};
157162

163+
// Priotize depth genotypes closer to reference
158164
final Comparator<Genotype> genotypeCopyNumberComparator = new Comparator<Genotype>() {
159165
@Override
160166
public int compare(Genotype o1, Genotype o2) {
161167
final int expectedQualityNumber1 = VariantContextGetters.getAttributeAsInt(o1, GATKSVVCFConstants.EXPECTED_COPY_NUMBER_FORMAT, 0);
162168
final int copyNumber1 = VariantContextGetters.getAttributeAsInt(o1, GATKSVVCFConstants.COPY_NUMBER_FORMAT, 0);
163169
final int expectedQualityNumber2 = VariantContextGetters.getAttributeAsInt(o2, GATKSVVCFConstants.EXPECTED_COPY_NUMBER_FORMAT, 0);
164170
final int copyNumber2 = VariantContextGetters.getAttributeAsInt(o2, GATKSVVCFConstants.COPY_NUMBER_FORMAT, 0);
165-
return Double.compare(Math.abs(expectedQualityNumber1 - copyNumber1), Math.abs(expectedQualityNumber2 - copyNumber2));
171+
return Double.compare(Math.abs(expectedQualityNumber2 - copyNumber2), Math.abs(expectedQualityNumber1 - copyNumber1));
172+
}
173+
};
174+
175+
// Priotize DEL over DUP as final tiebreaker
176+
final Comparator<Genotype> genotypeDelOverDupComparator = new Comparator<Genotype>() {
177+
@Override
178+
public int compare(Genotype o1, Genotype o2) {
179+
final int expectedCN1 = VariantContextGetters.getAttributeAsInt(o1, GATKSVVCFConstants.EXPECTED_COPY_NUMBER_FORMAT, 0);
180+
final boolean isDel1 = VariantContextGetters.getAttributeAsInt(o1, GATKSVVCFConstants.COPY_NUMBER_FORMAT, expectedCN1) < expectedCN1;
181+
final int expectedCN2 = VariantContextGetters.getAttributeAsInt(o2, GATKSVVCFConstants.EXPECTED_COPY_NUMBER_FORMAT, 0);
182+
final boolean isDel2 = VariantContextGetters.getAttributeAsInt(o2, GATKSVVCFConstants.COPY_NUMBER_FORMAT, expectedCN2) < expectedCN2;
183+
if (isDel1 && !isDel2) {
184+
return 1;
185+
} else if (isDel2 && !isDel1) {
186+
return -1;
187+
} else {
188+
return 0;
189+
}
166190
}
167191
};
168192

@@ -461,7 +485,8 @@ protected Genotype getRepresentativeGenotype(final Collection<Genotype> genotype
461485
.thenComparing(genotypeQualityComparator)
462486
.thenComparing(genotypeNonRefCountComparator)
463487
.thenComparing(genotypeCopyNumberQualityComparator)
464-
.thenComparing(genotypeCopyNumberComparator)).get();
488+
.thenComparing(genotypeCopyNumberComparator)
489+
.thenComparing(genotypeDelOverDupComparator)).get();
465490
}
466491

467492

0 commit comments

Comments
 (0)