Search in sources :

Example 31 with ReadPileup

use of org.broadinstitute.hellbender.utils.pileup.ReadPileup in project gatk by broadinstitute.

the class ReferenceConfidenceModelUnitTest method testCalcNIndelInformativeReads.

@Test(dataProvider = "CalcNIndelInformativeReadsData")
public void testCalcNIndelInformativeReads(final String readBases, final String ref, final int maxIndelSize, final List<Integer> expected) {
    final byte qual = (byte) 30;
    final byte[] quals = Utils.dupBytes(qual, readBases.length());
    for (int i = 0; i < readBases.getBytes().length; i++) {
        final GATKRead read = ArtificialReadUtils.createArtificialRead(readBases.getBytes(), quals, readBases.length() + "M");
        final SimpleInterval loc = new SimpleInterval("20", i + 1, i + 1);
        final ReadPileup pileup = new ReadPileup(loc, Collections.singletonList(read), i);
        final int actual = model.calcNIndelInformativeReads(pileup, i, ref.getBytes(), maxIndelSize);
        Assert.assertEquals(actual, (int) expected.get(i), "failed at position " + i);
    }
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) ReadPileup(org.broadinstitute.hellbender.utils.pileup.ReadPileup) SimpleInterval(org.broadinstitute.hellbender.utils.SimpleInterval) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 32 with ReadPileup

use of org.broadinstitute.hellbender.utils.pileup.ReadPileup in project gatk by broadinstitute.

the class AlignmentContextUnitTest method test1Sample2Readgroups.

@Test
public void test1Sample2Readgroups() throws Exception {
    final SAMReadGroupRecord readGroupOne = new SAMReadGroupRecord("rg1");
    readGroupOne.setSample("sample1");
    final SAMReadGroupRecord readGroupTwo = new SAMReadGroupRecord("rg2");
    readGroupTwo.setSample("sample1");
    final SAMFileHeader header = ArtificialReadUtils.createArtificialSamHeader(1, 1, 1000);
    header.addReadGroup(readGroupOne);
    header.addReadGroup(readGroupTwo);
    final Locatable loc = new SimpleInterval("chr1", 1, 1);
    final GATKRead read1 = ArtificialReadUtils.createArtificialRead(header, "read1", 0, 1, 10);
    read1.setReadGroup(readGroupOne.getId());
    final GATKRead read2 = ArtificialReadUtils.createArtificialRead(header, "read2", 0, 1, 10);
    read2.setReadGroup(readGroupTwo.getId());
    final GATKRead read3 = ArtificialReadUtils.createArtificialRead(header, "read3", 0, 1, 10);
    read3.setReadGroup(readGroupOne.getId());
    final GATKRead read4 = ArtificialReadUtils.createArtificialRead(header, "read4", 0, 1, 10);
    read4.setReadGroup(readGroupTwo.getId());
    final GATKRead read5 = ArtificialReadUtils.createArtificialRead(header, "read5", 0, 1, 10);
    read5.setReadGroup(readGroupTwo.getId());
    final GATKRead read6 = ArtificialReadUtils.createArtificialRead(header, "read6", 0, 1, 10);
    read6.setReadGroup(readGroupOne.getId());
    final GATKRead read7 = ArtificialReadUtils.createArtificialRead(header, "read7", 0, 1, 10);
    read7.setReadGroup(readGroupOne.getId());
    final ReadPileup pileup = new ReadPileup(loc, Arrays.asList(read1, read2, read3, read4, read5, read6, read7), 1);
    final AlignmentContext ac = new AlignmentContext(loc, pileup);
    Assert.assertSame(ac.getBasePileup(), pileup);
    Assert.assertEquals(ac.getContig(), loc.getContig());
    Assert.assertEquals(ac.getEnd(), loc.getEnd());
    Assert.assertEquals(ac.getLocation(), loc);
    Assert.assertEquals(ac.getPosition(), loc.getStart());
    Assert.assertEquals(ac.getStart(), loc.getStart());
    Assert.assertEquals(ac.hasPileupBeenDownsampled(), false);
    Assert.assertEquals(ac.size(), pileup.size());
    Assert.assertSame(ac.stratify(AlignmentContext.ReadOrientation.COMPLETE), ac, "Complete");
    final AlignmentContext acFWD = ac.stratify(AlignmentContext.ReadOrientation.FORWARD);
    Assert.assertEquals(acFWD.getLocation(), loc, "Forward Loc");
    Assert.assertEquals((Iterable<?>) acFWD.getBasePileup(), (Iterable<?>) pileup, "Forward Pileup");
    final AlignmentContext acREV = ac.stratify(AlignmentContext.ReadOrientation.REVERSE);
    AlignmentContext emptyAC = new AlignmentContext(loc, new ReadPileup(loc));
    Assert.assertEquals(acREV.getLocation(), loc, "Reverse Loc");
    Assert.assertEquals((Iterable<?>) acREV.getBasePileup(), (Iterable<?>) emptyAC.getBasePileup(), "Reverse pileup");
    Assert.assertNotNull(ac.toString());
    final Map<String, AlignmentContext> bySample = ac.splitContextBySampleName(header);
    Assert.assertEquals(bySample.size(), 1);
    Assert.assertEquals(bySample.keySet(), ReadUtils.getSamplesFromHeader(header));
    final AlignmentContext firstAC = bySample.values().iterator().next();
    Assert.assertEquals(firstAC.getLocation(), ac.getLocation());
    Assert.assertEquals(firstAC.getBasePileup(), ac.getBasePileup());
    final Map<String, AlignmentContext> bySampleAssume1 = ac.splitContextBySampleName("sample1", header);
    Assert.assertEquals(bySampleAssume1.size(), 1);
    Assert.assertEquals(bySampleAssume1.keySet(), ReadUtils.getSamplesFromHeader(header));
    final AlignmentContext firstACAssume1 = bySampleAssume1.values().iterator().next();
    Assert.assertEquals(firstACAssume1.getLocation(), ac.getLocation());
    Assert.assertEquals(firstACAssume1.getBasePileup(), ac.getBasePileup());
    final Map<String, AlignmentContext> stringAlignmentContextMap = AlignmentContext.splitContextBySampleName(pileup, header);
    Assert.assertEquals(stringAlignmentContextMap.keySet(), Collections.singleton("sample1"));
    Assert.assertEquals(stringAlignmentContextMap.get("sample1").getLocation(), loc);
    Assert.assertEquals(stringAlignmentContextMap.get("sample1").getBasePileup(), pileup);
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) ReadPileup(org.broadinstitute.hellbender.utils.pileup.ReadPileup) SAMReadGroupRecord(htsjdk.samtools.SAMReadGroupRecord) SimpleInterval(org.broadinstitute.hellbender.utils.SimpleInterval) SAMFileHeader(htsjdk.samtools.SAMFileHeader) Locatable(htsjdk.samtools.util.Locatable) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 33 with ReadPileup

use of org.broadinstitute.hellbender.utils.pileup.ReadPileup in project gatk by broadinstitute.

the class LocusIteratorByStateUnitTest method testLIBS_ComplexPileupTests.

@Test(enabled = true, dataProvider = "LIBS_ComplexPileupTests")
public void testLIBS_ComplexPileupTests(final int nReadsPerLocus, final int nLoci, final int nSamples, final boolean keepReads, final boolean grabReadsAfterEachCycle, final int downsampleTo) {
    final int readLength = 10;
    final boolean downsample = downsampleTo != -1;
    final DownsamplingMethod downsampler = downsample ? new DownsamplingMethod(DownsampleType.BY_SAMPLE, downsampleTo, null) : new DownsamplingMethod(DownsampleType.NONE, null, null);
    final ArtificialBAMBuilder bamBuilder = new ArtificialBAMBuilder(header.getSequenceDictionary(), nReadsPerLocus, nLoci);
    bamBuilder.createAndSetHeader(nSamples).setReadLength(readLength).setAlignmentStart(1);
    final List<GATKRead> reads = bamBuilder.makeReads();
    final LocusIteratorByState li;
    li = new LocusIteratorByState(new FakeCloseableIterator<>(reads.iterator()), downsampler, keepReads, bamBuilder.getSamples(), bamBuilder.getHeader(), true);
    final Set<GATKRead> seenSoFar = new LinkedHashSet<>();
    final Set<GATKRead> keptReads = new LinkedHashSet<>();
    int bpVisited = 0;
    while (li.hasNext()) {
        bpVisited++;
        final AlignmentContext alignmentContext = li.next();
        final ReadPileup p = alignmentContext.getBasePileup();
        AssertWellOrderedPileup(p);
        if (downsample) {
        // just not a safe test
        //Assert.assertTrue(p.getNumberOfElements() <= maxDownsampledCoverage * nSamples, "Too many reads at locus after downsampling");
        } else {
            final int minPileupSize = nReadsPerLocus * nSamples;
            Assert.assertTrue(p.size() >= minPileupSize);
        }
        // the number of reads starting here
        int nReadsStartingHere = 0;
        for (final GATKRead read : p.getReads()) if (read.getStart() == alignmentContext.getPosition())
            nReadsStartingHere++;
        // we can have no more than maxDownsampledCoverage per sample
        final int maxCoveragePerLocus = downsample ? downsampleTo : nReadsPerLocus;
        Assert.assertTrue(nReadsStartingHere <= maxCoveragePerLocus * nSamples);
        seenSoFar.addAll(p.getReads());
        if (keepReads && grabReadsAfterEachCycle) {
            final List<GATKRead> locusReads = li.transferReadsFromAllPreviousPileups();
            if (downsample) {
                // with downsampling we might have some reads here that were downsampled away
                // in the pileup.  We want to ensure that no more than the max coverage per sample is added
                Assert.assertTrue(locusReads.size() >= nReadsStartingHere);
                Assert.assertTrue(locusReads.size() <= maxCoveragePerLocus * nSamples);
            } else {
                Assert.assertEquals(locusReads.size(), nReadsStartingHere);
            }
            keptReads.addAll(locusReads);
            // check that all reads we've seen so far are in our keptReads
            for (final GATKRead read : seenSoFar) {
                Assert.assertTrue(keptReads.contains(read), "A read that appeared in a pileup wasn't found in the kept reads: " + read);
            }
        }
        if (!keepReads)
            Assert.assertTrue(li.getReadsFromAllPreviousPileups().isEmpty(), "Not keeping reads but the underlying list of reads isn't empty");
    }
    if (keepReads && !grabReadsAfterEachCycle)
        keptReads.addAll(li.transferReadsFromAllPreviousPileups());
    if (!downsample) {
        // downsampling may drop loci
        final int expectedBpToVisit = nLoci + readLength - 1;
        Assert.assertEquals(bpVisited, expectedBpToVisit, "Didn't visit the expected number of bp");
    }
    if (keepReads) {
        // check we have the right number of reads
        final int totalReads = nLoci * nReadsPerLocus * nSamples;
        if (!downsample) {
            // downsampling may drop reads
            Assert.assertEquals(keptReads.size(), totalReads, "LIBS didn't keep the right number of reads during the traversal");
            // check that the order of reads is the same as in our read list
            for (int i = 0; i < reads.size(); i++) {
                final GATKRead inputRead = reads.get(i);
                final GATKRead keptRead = reads.get(i);
                Assert.assertSame(keptRead, inputRead, "Input reads and kept reads differ at position " + i);
            }
        } else {
            Assert.assertTrue(keptReads.size() <= totalReads, "LIBS didn't keep the right number of reads during the traversal");
        }
        // check uniqueness
        final Set<String> readNames = new LinkedHashSet<>();
        for (final GATKRead read : keptReads) {
            Assert.assertFalse(readNames.contains(read.getName()), "Found duplicate reads in the kept reads");
            readNames.add(read.getName());
        }
        // check that all reads we've seen are in our keptReads
        for (final GATKRead read : seenSoFar) {
            Assert.assertTrue(keptReads.contains(read), "A read that appeared in a pileup wasn't found in the kept reads: " + read);
        }
        if (!downsample) {
            // check that every read in the list of keep reads occurred at least once in one of the pileups
            for (final GATKRead keptRead : keptReads) {
                Assert.assertTrue(seenSoFar.contains(keptRead), "There's a read " + keptRead + " in our keptReads list that never appeared in any pileup");
            }
        }
    }
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) ReadPileup(org.broadinstitute.hellbender.utils.pileup.ReadPileup) ArtificialBAMBuilder(org.broadinstitute.hellbender.utils.read.ArtificialBAMBuilder) AlignmentContext(org.broadinstitute.hellbender.engine.AlignmentContext) DownsamplingMethod(org.broadinstitute.hellbender.utils.downsampling.DownsamplingMethod) Test(org.testng.annotations.Test)

Example 34 with ReadPileup

use of org.broadinstitute.hellbender.utils.pileup.ReadPileup in project gatk by broadinstitute.

the class LocusIteratorByStateUnitTest method testUnmappedAndAllIReadsPassThrough.

@Test(enabled = false)
public void testUnmappedAndAllIReadsPassThrough() {
    final int readLength = 10;
    final GATKRead mapped1 = ArtificialReadUtils.createArtificialRead(header, "mapped1", 0, 1, readLength);
    final GATKRead mapped2 = ArtificialReadUtils.createArtificialRead(header, "mapped2", 0, 1, readLength);
    final GATKRead unmapped = ArtificialReadUtils.createArtificialRead(header, "unmapped", 0, 1, readLength);
    final GATKRead allI = ArtificialReadUtils.createArtificialRead(header, "allI", 0, 1, readLength);
    unmapped.setIsUnmapped();
    unmapped.setCigar("*");
    allI.setCigar(readLength + "I");
    final List<GATKRead> reads = Arrays.asList(mapped1, unmapped, allI, mapped2);
    // create the iterator by state with the fake reads and fake records
    final LocusIteratorByState li;
    li = makeLIBS(reads, DownsamplingMethod.NONE, true, header);
    Assert.assertTrue(li.hasNext());
    AlignmentContext context = li.next();
    ReadPileup pileup = context.getBasePileup();
    Assert.assertEquals(pileup.size(), 2, "Should see only 2 reads in pileup, even with unmapped and all I reads");
    final List<GATKRead> rawReads = li.transferReadsFromAllPreviousPileups();
    Assert.assertEquals(rawReads, reads, "Input and transferred read lists should be the same, and include the unmapped and all I reads");
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) AlignmentContext(org.broadinstitute.hellbender.engine.AlignmentContext) ReadPileup(org.broadinstitute.hellbender.utils.pileup.ReadPileup) Test(org.testng.annotations.Test)

Example 35 with ReadPileup

use of org.broadinstitute.hellbender.utils.pileup.ReadPileup in project gatk by broadinstitute.

the class LocusIteratorByStateUnitTest method testAllIReadsPassThrough.

@Test
public void testAllIReadsPassThrough() {
    final int readLength = 10;
    final GATKRead mapped1 = ArtificialReadUtils.createArtificialRead(header, "mapped1", 0, 1, readLength);
    final GATKRead mapped2 = ArtificialReadUtils.createArtificialRead(header, "mapped2", 0, 1, readLength);
    final GATKRead allI = ArtificialReadUtils.createArtificialRead(header, "allI", 0, 1, readLength);
    allI.setCigar(readLength + "I");
    final List<GATKRead> reads = Arrays.asList(mapped1, allI, mapped2);
    // create the iterator by state with the fake reads and fake records
    final LocusIteratorByState li;
    li = makeLIBS(reads, DownsamplingMethod.NONE, true, header);
    Assert.assertTrue(li.hasNext());
    final AlignmentContext context = li.next();
    final ReadPileup pileup = context.getBasePileup();
    Assert.assertEquals(pileup.size(), 2, "Should see only 2 reads in pileup, even with all I reads");
    final List<GATKRead> rawReads = li.transferReadsFromAllPreviousPileups();
    Assert.assertEquals(rawReads, reads, "Input and transferred read lists should be the same, and include and all I reads");
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) AlignmentContext(org.broadinstitute.hellbender.engine.AlignmentContext) ReadPileup(org.broadinstitute.hellbender.utils.pileup.ReadPileup) Test(org.testng.annotations.Test)

Aggregations

ReadPileup (org.broadinstitute.hellbender.utils.pileup.ReadPileup)36 GATKRead (org.broadinstitute.hellbender.utils.read.GATKRead)24 Test (org.testng.annotations.Test)19 AlignmentContext (org.broadinstitute.hellbender.engine.AlignmentContext)14 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)13 SAMFileHeader (htsjdk.samtools.SAMFileHeader)11 Locatable (htsjdk.samtools.util.Locatable)11 BaseTest (org.broadinstitute.hellbender.utils.test.BaseTest)9 PileupElement (org.broadinstitute.hellbender.utils.pileup.PileupElement)7 java.util (java.util)5 VariantContext (htsjdk.variant.variantcontext.VariantContext)4 Collectors (java.util.stream.Collectors)4 SAMReadGroupRecord (htsjdk.samtools.SAMReadGroupRecord)3 FeatureContext (org.broadinstitute.hellbender.engine.FeatureContext)3 ReferenceContext (org.broadinstitute.hellbender.engine.ReferenceContext)3 VisibleForTesting (com.google.common.annotations.VisibleForTesting)2 CigarOperator (htsjdk.samtools.CigarOperator)2 htsjdk.variant.variantcontext (htsjdk.variant.variantcontext)2 htsjdk.variant.vcf (htsjdk.variant.vcf)2 VCFConstants (htsjdk.variant.vcf.VCFConstants)2