Search in sources :

Example 6 with GenomeLoc

use of org.broadinstitute.hellbender.utils.GenomeLoc in project gatk-protected by broadinstitute.

the class ReferenceConfidenceModelUnitTest method checkOverlapping.

private void checkOverlapping(final int pos, Collection<VariantContext> calls, final VariantContext expected) {
    final GenomeLoc loc = parser.createGenomeLoc(parser.getSequenceDictionary().getSequences().get(0).getSequenceName(), pos, pos);
    final VariantContext actual = model.getOverlappingVariantContext(loc, calls);
    Assert.assertEquals(actual, expected);
}
Also used : VariantContext(htsjdk.variant.variantcontext.VariantContext) GenomeLoc(org.broadinstitute.hellbender.utils.GenomeLoc)

Example 7 with GenomeLoc

use of org.broadinstitute.hellbender.utils.GenomeLoc in project gatk by broadinstitute.

the class ReferenceConfidenceModelUnitTest method checkOverlapping.

private void checkOverlapping(final int pos, Collection<VariantContext> calls, final VariantContext expected) {
    final GenomeLoc loc = parser.createGenomeLoc(parser.getSequenceDictionary().getSequences().get(0).getSequenceName(), pos, pos);
    final VariantContext actual = model.getOverlappingVariantContext(loc, calls);
    Assert.assertEquals(actual, expected);
}
Also used : VariantContext(htsjdk.variant.variantcontext.VariantContext) GenomeLoc(org.broadinstitute.hellbender.utils.GenomeLoc)

Example 8 with GenomeLoc

use of org.broadinstitute.hellbender.utils.GenomeLoc in project gatk by broadinstitute.

the class ReferenceConfidenceModelUnitTest method checkReferenceModelResult.

private void checkReferenceModelResult(final RefConfData data, final List<VariantContext> contexts, final List<Integer> expectedDPs, final List<VariantContext> calls) {
    Assert.assertNotNull(contexts);
    final SimpleInterval loc = data.getActiveRegion().getExtendedSpan();
    final List<Boolean> seenBP = new ArrayList<>(Collections.nCopies(data.getActiveRegion().getSpan().size(), false));
    for (int i = 0; i < loc.size(); i++) {
        final GenomeLoc curPos = parser.createGenomeLoc(loc.getContig(), loc.getStart() + i);
        final VariantContext call = model.getOverlappingVariantContext(curPos, calls);
        final VariantContext refModel = model.getOverlappingVariantContext(curPos, contexts);
        if (!data.getActiveRegion().getSpan().contains(curPos)) {
            // part of the extended interval, but not the full interval
            Assert.assertNull(refModel);
            continue;
        }
        if (call != null) {
            if (call.isVariant() && refModel.getType() == VariantContext.Type.SYMBOLIC) {
                //Assert.assertEquals(refModel, call, "Should have found call " + call + " but found " + refModel + " instead");
                // must be a deletion.
                Assert.assertTrue(call.getReference().length() > 1);
                // the deletion must not start at the same position
                Assert.assertTrue(call.getStart() < refModel.getStart());
                Assert.assertEquals(call.getReference().getBaseString().substring(refModel.getStart() - call.getStart(), refModel.getStart() - call.getStart() + 1), refModel.getReference().getBaseString(), // the reference must be the same.
                "" + data.getRefHap());
                // No confidence in the reference hom-ref call across the deletion
                Assert.assertTrue(refModel.getGenotype(0).getGQ() <= 0);
                // the reference and the lonelly <NON_REF>
                Assert.assertEquals(refModel.getAlleles().size(), 2);
                Assert.assertEquals(refModel.getAlleles().get(1), GATKVCFConstants.NON_REF_SYMBOLIC_ALLELE);
            } else {
                Assert.assertEquals(refModel, call, "Should have found call " + call + " but found " + refModel + " instead");
            }
        } else {
            final int expectedDP = expectedDPs.get(curPos.getStart() - data.getActiveRegion().getSpan().getStart());
            Assert.assertEquals(refModel.getStart(), loc.getStart() + i);
            Assert.assertEquals(refModel.getEnd(), loc.getStart() + i);
            Assert.assertFalse(refModel.hasLog10PError());
            Assert.assertEquals(refModel.getAlternateAlleles().size(), 1);
            Assert.assertEquals(refModel.getAlternateAllele(0), GATKVCFConstants.NON_REF_SYMBOLIC_ALLELE);
            Assert.assertTrue(refModel.hasGenotype(sample));
            final Genotype g = refModel.getGenotype(sample);
            Assert.assertTrue(g.hasAD());
            Assert.assertTrue(g.hasDP());
            Assert.assertEquals(g.getDP(), expectedDP);
            Assert.assertTrue(g.hasGQ());
            Assert.assertTrue(g.hasPL());
        }
        final VariantContext vc = call == null ? refModel : call;
        if (curPos.getStart() == vc.getStart()) {
            for (int pos = vc.getStart(); pos <= vc.getEnd(); pos++) {
                final int j = pos - data.getActiveRegion().getSpan().getStart();
                Assert.assertFalse(seenBP.get(j));
                seenBP.set(j, true);
            }
        }
    }
    for (int i = 0; i < seenBP.size(); i++) {
        Assert.assertEquals((boolean) seenBP.get(i), true);
    }
}
Also used : VariantContext(htsjdk.variant.variantcontext.VariantContext) Genotype(htsjdk.variant.variantcontext.Genotype) SimpleInterval(org.broadinstitute.hellbender.utils.SimpleInterval) GenomeLoc(org.broadinstitute.hellbender.utils.GenomeLoc)

Example 9 with GenomeLoc

use of org.broadinstitute.hellbender.utils.GenomeLoc in project gatk by broadinstitute.

the class BandPassActivityProfileUnitTest method testBandPass.

@Test(dataProvider = "BandPassBasicTest")
public void testBandPass(final int start, final boolean precedingIsActive, final int nPrecedingSites, final int bandPassSize, final double sigma) {
    final BandPassActivityProfile profile = new BandPassActivityProfile(null, MAX_PROB_PROPAGATION_DISTANCE, ACTIVE_PROB_THRESHOLD, bandPassSize, sigma, false, header);
    final int expectedBandSize = bandPassSize * 2 + 1;
    Assert.assertEquals(profile.getFilteredSize(), bandPassSize, "Wrong filter size");
    Assert.assertEquals(profile.getSigma(), sigma, "Wrong sigma");
    Assert.assertEquals(profile.getBandSize(), expectedBandSize, "Wrong expected band size");
    final String contig = genomeLocParser.getSequenceDictionary().getSequences().get(0).getSequenceName();
    final double precedingProb = precedingIsActive ? 1.0 : 0.0;
    for (int i = 0; i < nPrecedingSites; i++) {
        final GenomeLoc loc = genomeLocParser.createGenomeLoc(contig, i + start);
        final ActivityProfileState state = new ActivityProfileState(new SimpleInterval(loc), precedingProb);
        profile.add(state);
    }
    final GenomeLoc nextLoc = genomeLocParser.createGenomeLoc(contig, nPrecedingSites + start);
    profile.add(new ActivityProfileState(new SimpleInterval(nextLoc), 1.0));
    if (!precedingIsActive && nPrecedingSites >= bandPassSize && bandPassSize < start) {
        // we have enough space that all probs fall on the genome
        final double[] probs = profile.getProbabilitiesAsArray();
        Assert.assertEquals(MathUtils.sum(probs), 1.0 * (nPrecedingSites * precedingProb + 1), 1e-3, "Activity profile doesn't sum to number of non-zero prob states");
    }
}
Also used : SimpleInterval(org.broadinstitute.hellbender.utils.SimpleInterval) GenomeLoc(org.broadinstitute.hellbender.utils.GenomeLoc) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 10 with GenomeLoc

use of org.broadinstitute.hellbender.utils.GenomeLoc in project gatk by broadinstitute.

the class BandPassActivityProfileUnitTest method testBandPassComposition.

@Test(dataProvider = "BandPassComposition")
public void testBandPassComposition(final int bandPassSize, final int integrationLength) {
    final int start = 1;
    final BandPassActivityProfile profile = new BandPassActivityProfile(null, MAX_PROB_PROPAGATION_DISTANCE, ACTIVE_PROB_THRESHOLD, bandPassSize, BandPassActivityProfile.DEFAULT_SIGMA, header);
    final double[] rawActiveProbs = new double[integrationLength + bandPassSize * 2];
    // add a buffer so that we can get all of the band pass values
    final String contig = genomeLocParser.getSequenceDictionary().getSequences().get(0).getSequenceName();
    int pos = start;
    int rawProbsOffset = 0;
    for (int i = 0; i < bandPassSize; i++) {
        final GenomeLoc loc = genomeLocParser.createGenomeLoc(contig, pos++);
        final ActivityProfileState state = new ActivityProfileState(new SimpleInterval(loc), 0.0);
        profile.add(state);
        rawActiveProbs[rawProbsOffset++] = 0.0;
        rawActiveProbs[rawActiveProbs.length - rawProbsOffset] = 0.0;
    }
    for (int i = 0; i < integrationLength; i++) {
        final GenomeLoc nextLoc = genomeLocParser.createGenomeLoc(contig, pos++);
        profile.add(new ActivityProfileState(new SimpleInterval(nextLoc), 1.0));
        rawActiveProbs[rawProbsOffset++] = 1.0;
        for (int j = 0; j < profile.size(); j++) {
            Assert.assertTrue(profile.getStateList().get(j).isActiveProb() >= 0.0, "State probability < 0 at " + j);
            Assert.assertTrue(profile.getStateList().get(j).isActiveProb() <= 1.0 + 1e-3, "State probability > 1 at " + j);
        }
    }
    final double[] expectedProbs = bandPassInOnePass(profile, rawActiveProbs);
    for (int j = 0; j < profile.size(); j++) {
        Assert.assertEquals(profile.getStateList().get(j).isActiveProb(), expectedProbs[j], "State probability not expected at " + j);
    }
}
Also used : SimpleInterval(org.broadinstitute.hellbender.utils.SimpleInterval) GenomeLoc(org.broadinstitute.hellbender.utils.GenomeLoc) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Aggregations

GenomeLoc (org.broadinstitute.hellbender.utils.GenomeLoc)13 VariantContext (htsjdk.variant.variantcontext.VariantContext)6 BaseTest (org.broadinstitute.hellbender.utils.test.BaseTest)6 Test (org.testng.annotations.Test)6 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)5 UnvalidatingGenomeLoc (org.broadinstitute.hellbender.utils.UnvalidatingGenomeLoc)4 Genotype (htsjdk.variant.variantcontext.Genotype)2 DataProvider (org.testng.annotations.DataProvider)2 ArrayList (java.util.ArrayList)1