use of org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState in project gatk-protected by broadinstitute.
the class XHMMSegmentCallerIntegrationTest method assertSampleSegmentsCoordinates.
private void assertSampleSegmentsCoordinates(List<HiddenStateSegmentRecord<CopyNumberTriState, Target>> sampleRecords, TargetCollection<Target> targets) {
for (final HiddenStateSegmentRecord<CopyNumberTriState, Target> record : sampleRecords) {
final IndexRange range = targets.indexRange(record.getSegment());
Assert.assertTrue(range.size() > 0);
Assert.assertEquals(record.getSegment().getContig(), targets.location(range.from).getContig());
Assert.assertEquals(record.getSegment().getStart(), targets.location(range.from).getStart());
Assert.assertEquals(record.getSegment().getEnd(), targets.location(range.to - 1).getEnd());
}
}
use of org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState in project gatk-protected by broadinstitute.
the class XHMMSegmentCallerIntegrationTest method assertOutputHasConsistentNumberOfTargets.
private void assertOutputHasConsistentNumberOfTargets(final List<HiddenStateSegmentRecord<CopyNumberTriState, Target>> outputRecords, final TargetCollection<Target> targets) {
for (final HiddenStateSegmentRecord<CopyNumberTriState, Target> nextRecord : outputRecords) {
final IndexRange indexRange = targets.indexRange(nextRecord.getSegment());
Assert.assertEquals(indexRange.to - indexRange.from, nextRecord.getSegment().getTargetCount());
}
}
use of org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState in project gatk by broadinstitute.
the class CopyNumberTriStateUnitTest method testAlleleListsLogic.
@Test
public void testAlleleListsLogic() {
Assert.assertEquals(CopyNumberTriState.ALL_ALLELES.size(), CopyNumberTriState.values().length);
Assert.assertEquals(CopyNumberTriState.ALL_ALLELES.subList(1, CopyNumberTriState.ALL_ALLELES.size()), CopyNumberTriState.ALTERNATIVE_ALLELES);
Assert.assertFalse(CopyNumberTriState.ALTERNATIVE_ALLELES.stream().anyMatch(Allele::isReference));
Assert.assertTrue(CopyNumberTriState.ALL_ALLELES.get(0).isReference());
Assert.assertEquals(CopyNumberTriState.ALL_ALLELES.get(0), CopyNumberTriState.NEUTRAL.allele);
final List<Allele> expectedOrder = new ArrayList<>();
for (final CopyNumberTriState state : CopyNumberTriState.values()) {
if (state == CopyNumberTriState.NEUTRAL) {
continue;
}
expectedOrder.add(state.allele);
}
Assert.assertEquals(CopyNumberTriState.ALTERNATIVE_ALLELES, expectedOrder);
}
use of org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState in project gatk by broadinstitute.
the class ConvertGSVariantsToSegments method apply.
@Override
public void apply(final VariantContext variant, final ReadsContext readsContext, final ReferenceContext referenceContext, final FeatureContext featureContext) {
final SimpleInterval interval = new SimpleInterval(variant);
final int targetCount = targets.indexRange(interval).size();
final int[] callCounts = new int[CopyNumberTriState.values().length];
for (final Genotype genotype : variant.getGenotypes().iterateInSampleNameOrder()) {
final String sample = genotype.getSampleName();
final double mean = doubleFrom(genotype.getExtendedAttribute(GS_COPY_NUMBER_FRACTION));
final int copyNumber = intFrom(genotype.getExtendedAttribute(GS_COPY_NUMBER_FORMAT));
final CopyNumberTriState call = copyNumber == neutralCopyNumber ? CopyNumberTriState.NEUTRAL : (copyNumber < neutralCopyNumber) ? CopyNumberTriState.DELETION : CopyNumberTriState.DUPLICATION;
callCounts[call.ordinal()]++;
final double[] probs = doubleArrayFrom(genotype.getExtendedAttribute(GS_COPY_NUMBER_POSTERIOR));
final double log10PostQualCall = calculateLog10CallQuality(probs, call);
final double log10PostQualNonRef = calculateLog10CallQualityNonRef(probs);
final double phredProbCall = -10.0 * log10PostQualCall;
final double phredProbNonRef = -10.0 * log10PostQualNonRef;
final HiddenStateSegment<CopyNumberTriState, Target> segment = new HiddenStateSegment<>(interval, targetCount, mean, // GS VCF does not contain any stddev or var estimate for coverage fraction.
0.0, call, // GS does not provide an EQ, we approximate it to be the 1 - sum of all call compatible CN corresponding posterior probs
phredProbCall, // GS does not provide a SQ, we leave is a NaN.
Double.NaN, // GS does not provide a START Q.
Double.NaN, // GS does not provide a END Q.
Double.NaN, phredProbNonRef);
final HiddenStateSegmentRecord<CopyNumberTriState, Target> record = new HiddenStateSegmentRecord<>(sample, segment);
try {
outputWriter.writeRecord(record);
} catch (final IOException ex) {
throw new UserException.CouldNotCreateOutputFile(outputFile, ex);
}
}
}
use of org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState in project gatk by broadinstitute.
the class XHMMSegmentCallerIntegrationTest method testRunCommandLine.
//TODO: this test used to contain a test of concordance with XHMM. It no longer does that because our model has
//TODO: diverged from XHMM's. Eventually the right thing to do is use the simulateChain() method to generate
//TODO: simulated data for some artificial set of CNV segments and to test concordance with those segments.
//TODO: however we still use XHMM's emission model, which is both not generative and quite silly. Once we
//TODO: have a generative model of coverage we can modify simulateChain() accordingly and then write a concordance
//TODO: test here. Until then, we do not have an integration test but we do have our ongoing evaluations, which
//TODO: show the superiority of our modifications versus the original XHMM model.
@Test(dataProvider = "simulatedChainData")
public File testRunCommandLine(final XHMMData chain) {
final File inputFile = writeChainInTempFile(chain);
final File outputFile = createTempFile("output", ".tab");
runCommandLine(chain, inputFile, outputFile);
Assert.assertTrue(outputFile.exists());
final TargetCollection<Target> targets = TargetArgumentCollection.readTargetCollection(REALISTIC_TARGETS_FILE);
final List<HiddenStateSegmentRecord<CopyNumberTriState, Target>> outputRecords = readOutputRecords(outputFile);
assertOutputIsInOrder(outputRecords, targets);
assertOutputHasConsistentNumberOfTargets(outputRecords, targets);
final Map<String, List<HiddenStateSegmentRecord<CopyNumberTriState, Target>>> outputRecordsBySample = splitOutputRecordBySample(outputRecords);
assertSampleNames(outputRecordsBySample.keySet(), chain);
for (final List<HiddenStateSegmentRecord<CopyNumberTriState, Target>> sampleRecords : outputRecordsBySample.values()) {
assertSampleSegmentsCoverAllTargets(sampleRecords, targets);
assertSampleSegmentsCoordinates(sampleRecords, targets);
}
return outputFile;
}
Aggregations