Search in sources :

Example 6 with CopyNumberTriState

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());
    }
}
Also used : IndexRange(org.broadinstitute.hellbender.utils.IndexRange) Target(org.broadinstitute.hellbender.tools.exome.Target) CopyNumberTriState(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState)

Example 7 with CopyNumberTriState

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());
    }
}
Also used : IndexRange(org.broadinstitute.hellbender.utils.IndexRange) Target(org.broadinstitute.hellbender.tools.exome.Target) CopyNumberTriState(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState)

Example 8 with CopyNumberTriState

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);
}
Also used : Allele(htsjdk.variant.variantcontext.Allele) CopyNumberTriState(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState) ArrayList(java.util.ArrayList) Test(org.testng.annotations.Test)

Example 9 with CopyNumberTriState

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);
        }
    }
}
Also used : Genotype(htsjdk.variant.variantcontext.Genotype) IOException(java.io.IOException) Target(org.broadinstitute.hellbender.tools.exome.Target) CopyNumberTriState(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState) HiddenStateSegment(org.broadinstitute.hellbender.utils.hmm.segmentation.HiddenStateSegment) HiddenStateSegmentRecord(org.broadinstitute.hellbender.utils.hmm.segmentation.HiddenStateSegmentRecord) SimpleInterval(org.broadinstitute.hellbender.utils.SimpleInterval) UserException(org.broadinstitute.hellbender.exceptions.UserException)

Example 10 with CopyNumberTriState

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;
}
Also used : Target(org.broadinstitute.hellbender.tools.exome.Target) CopyNumberTriState(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState) HiddenStateSegmentRecord(org.broadinstitute.hellbender.utils.hmm.segmentation.HiddenStateSegmentRecord) ArrayList(java.util.ArrayList) List(java.util.List) File(java.io.File) Test(org.testng.annotations.Test)

Aggregations

CopyNumberTriState (org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriState)24 Target (org.broadinstitute.hellbender.tools.exome.Target)20 HiddenStateSegmentRecord (org.broadinstitute.hellbender.utils.hmm.segmentation.HiddenStateSegmentRecord)10 Test (org.testng.annotations.Test)10 Genotype (htsjdk.variant.variantcontext.Genotype)8 IndexRange (org.broadinstitute.hellbender.utils.IndexRange)8 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)8 Allele (htsjdk.variant.variantcontext.Allele)6 File (java.io.File)6 IOException (java.io.IOException)6 ArrayList (java.util.ArrayList)6 VariantContext (htsjdk.variant.variantcontext.VariantContext)4 htsjdk.variant.vcf (htsjdk.variant.vcf)4 java.util (java.util)4 Collectors (java.util.stream.Collectors)4 IntStream (java.util.stream.IntStream)4 StandardArgumentDefinitions (org.broadinstitute.hellbender.cmdline.StandardArgumentDefinitions)4 TargetArgumentCollection (org.broadinstitute.hellbender.tools.exome.TargetArgumentCollection)4 CopyNumberTriStateAllele (org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriStateAllele)4 CopyNumberTriStateTransitionProbabilityCache (org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriStateTransitionProbabilityCache)4