Search in sources :

Example 26 with Haplotype

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

the class PairHMMLikelihoodCalculationEngineUnitTest method testComputeLikelihoods.

@Test
public void testComputeLikelihoods() {
    final LikelihoodEngineArgumentCollection LEAC = new LikelihoodEngineArgumentCollection();
    PairHMMLikelihoodCalculationEngine.writeLikelihoodsToFile = true;
    final ReadLikelihoodCalculationEngine lce = new PairHMMLikelihoodCalculationEngine((byte) SAMUtils.MAX_PHRED_SCORE, new PairHMMNativeArguments(), PairHMM.Implementation.LOGLESS_CACHING, MathUtils.logToLog10(QualityUtils.qualToErrorProbLog10(LEAC.phredScaledGlobalReadMismappingRate)), PairHMMLikelihoodCalculationEngine.PCRErrorModel.CONSERVATIVE);
    final Map<String, List<GATKRead>> perSampleReadList = new HashMap<>();
    final int n = 10;
    final GATKRead read1 = ArtificialReadUtils.createArtificialRead(TextCigarCodec.decode(n + "M"));
    read1.setMappingQuality(60);
    final String sample1 = "sample1";
    perSampleReadList.put(sample1, Arrays.asList(read1));
    final SampleList samples = new IndexedSampleList(sample1);
    final AssemblyResultSet assemblyResultSet = new AssemblyResultSet();
    final byte[] bases = Strings.repeat("A", n + 1).getBytes();
    final Haplotype hap1 = new Haplotype(bases, true);
    hap1.setGenomeLocation(read1);
    assemblyResultSet.add(hap1);
    final byte[] basesModified = bases;
    //different bases
    basesModified[5] = 'C';
    final Haplotype hap2 = new Haplotype(basesModified, false);
    //use same loc
    hap2.setGenomeLocation(read1);
    assemblyResultSet.add(hap2);
    final ReadLikelihoods<Haplotype> likes = lce.computeReadLikelihoods(assemblyResultSet, samples, perSampleReadList);
    final LikelihoodMatrix<Haplotype> mtx = likes.sampleMatrix(0);
    Assert.assertEquals(mtx.numberOfAlleles(), 2);
    Assert.assertEquals(mtx.numberOfReads(), 1);
    final double v1 = mtx.get(0, 0);
    final double v2 = mtx.get(1, 0);
    Assert.assertTrue(v1 > v2, "matching haplotype should have a higher likelihood");
    lce.close();
    new File(PairHMMLikelihoodCalculationEngine.LIKELIHOODS_FILENAME).delete();
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) PairHMMNativeArguments(org.broadinstitute.gatk.nativebindings.pairhmm.PairHMMNativeArguments) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) SampleList(org.broadinstitute.hellbender.utils.genotyper.SampleList) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) Haplotype(org.broadinstitute.hellbender.utils.haplotype.Haplotype) File(java.io.File) SampleList(org.broadinstitute.hellbender.utils.genotyper.SampleList) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 27 with Haplotype

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

the class ReferenceConfidenceModelUnitTest method testRefConfidenceBasic.

@Test(dataProvider = "RefConfidenceData")
public void testRefConfidenceBasic(final int nReads, final int extension) {
    final RefConfData data = new RefConfData("ACGTAACCGGTT", extension);
    final List<Haplotype> haplotypes = Arrays.asList(data.getRefHap());
    final List<VariantContext> calls = Collections.emptyList();
    for (int i = 0; i < nReads; i++) {
        data.getActiveRegion().add(data.makeRead(0, data.getRefLength()));
    }
    final ReadLikelihoods<Haplotype> likelihoods = createDummyStratifiedReadMap(data.getRefHap(), samples, data.getActiveRegion());
    final PloidyModel ploidyModel = new HomogeneousPloidyModel(samples, 2);
    final IndependentSampleGenotypesModel genotypingModel = new IndependentSampleGenotypesModel();
    final List<Integer> expectedDPs = Collections.nCopies(data.getActiveRegion().getSpan().size(), nReads);
    final List<VariantContext> contexts = model.calculateRefConfidence(data.getRefHap(), haplotypes, data.getPaddedRefLoc(), data.getActiveRegion(), likelihoods, ploidyModel, calls);
    checkReferenceModelResult(data, contexts, expectedDPs, calls);
}
Also used : IndependentSampleGenotypesModel(org.broadinstitute.hellbender.tools.walkers.genotyper.IndependentSampleGenotypesModel) HomogeneousPloidyModel(org.broadinstitute.hellbender.tools.walkers.genotyper.HomogeneousPloidyModel) VariantContext(htsjdk.variant.variantcontext.VariantContext) HomogeneousPloidyModel(org.broadinstitute.hellbender.tools.walkers.genotyper.HomogeneousPloidyModel) PloidyModel(org.broadinstitute.hellbender.tools.walkers.genotyper.PloidyModel) Haplotype(org.broadinstitute.hellbender.utils.haplotype.Haplotype) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 28 with Haplotype

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

the class ReadThreadingAssemblerUnitTest method testAssemblyWithVariant.

private void testAssemblyWithVariant(final ReadThreadingAssembler assembler, final byte[] refBases, final SimpleInterval loc, final int nReadsToUse, final VariantContext site) {
    final String preRef = new String(refBases).substring(0, site.getStart());
    final String postRef = new String(refBases).substring(site.getEnd() + 1, refBases.length);
    final byte[] altBases = (preRef + site.getAlternateAllele(0).getBaseString() + postRef).getBytes();
    //        logger.warn("ref " + new String(refBases));
    //        logger.warn("alt " + new String(altBases));
    final List<GATKRead> reads = new LinkedList<>();
    for (int i = 0; i < nReadsToUse; i++) {
        final byte[] bases = altBases.clone();
        final byte[] quals = Utils.dupBytes((byte) 30, altBases.length);
        final String cigar = altBases.length + "M";
        final GATKRead read = ArtificialReadUtils.createArtificialRead(header, loc.getContig(), loc.getContig(), loc.getStart(), bases, quals, cigar);
        reads.add(read);
    }
    final Haplotype refHaplotype = new Haplotype(refBases, true);
    final Haplotype altHaplotype = new Haplotype(altBases, false);
    final List<Haplotype> haplotypes = assemble(assembler, refBases, loc, reads);
    Assert.assertEquals(haplotypes, Arrays.asList(refHaplotype, altHaplotype));
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) Haplotype(org.broadinstitute.hellbender.utils.haplotype.Haplotype) KBestHaplotype(org.broadinstitute.hellbender.tools.walkers.haplotypecaller.graphs.KBestHaplotype)

Example 29 with Haplotype

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

the class ReadThreadingAssemblerUnitTest method testAssembleRef.

@Test(dataProvider = "AssembleIntervalsData")
public void testAssembleRef(final ReadThreadingAssembler assembler, final SimpleInterval loc, final int nReadsToUse) {
    final byte[] refBases = seq.getSubsequenceAt(loc.getContig(), loc.getStart(), loc.getEnd()).getBases();
    final List<GATKRead> reads = new LinkedList<>();
    for (int i = 0; i < nReadsToUse; i++) {
        final byte[] bases = refBases.clone();
        final byte[] quals = Utils.dupBytes((byte) 30, refBases.length);
        final String cigar = refBases.length + "M";
        final GATKRead read = ArtificialReadUtils.createArtificialRead(header, loc.getContig(), loc.getContig(), loc.getStart(), bases, quals, cigar);
        reads.add(read);
    }
    // TODO -- generalize to all assemblers
    final Haplotype refHaplotype = new Haplotype(refBases, true);
    final List<Haplotype> haplotypes = assemble(assembler, refBases, loc, reads);
    Assert.assertEquals(haplotypes, Collections.singletonList(refHaplotype));
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) Haplotype(org.broadinstitute.hellbender.utils.haplotype.Haplotype) KBestHaplotype(org.broadinstitute.hellbender.tools.walkers.haplotypecaller.graphs.KBestHaplotype) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 30 with Haplotype

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

the class ReadThreadingAssemblerUnitTest method assemble.

private List<Haplotype> assemble(final ReadThreadingAssembler assembler, final byte[] refBases, final SimpleInterval loc, final List<GATKRead> reads) {
    final Haplotype refHaplotype = new Haplotype(refBases, true);
    final Cigar c = new Cigar();
    c.add(new CigarElement(refHaplotype.getBases().length, CigarOperator.M));
    refHaplotype.setCigar(c);
    final AssemblyRegion activeRegion = new AssemblyRegion(loc, null, true, 0, header);
    activeRegion.addAll(reads);
    //        logger.warn("Assembling " + activeRegion + " with " + engine);
    final AssemblyResultSet assemblyResultSet = assembler.runLocalAssembly(activeRegion, refHaplotype, refBases, loc, Collections.<VariantContext>emptyList(), null, header);
    return assemblyResultSet.getHaplotypeList();
}
Also used : Cigar(htsjdk.samtools.Cigar) AssemblyRegion(org.broadinstitute.hellbender.engine.AssemblyRegion) AssemblyResultSet(org.broadinstitute.hellbender.tools.walkers.haplotypecaller.AssemblyResultSet) Haplotype(org.broadinstitute.hellbender.utils.haplotype.Haplotype) KBestHaplotype(org.broadinstitute.hellbender.tools.walkers.haplotypecaller.graphs.KBestHaplotype) CigarElement(htsjdk.samtools.CigarElement)

Aggregations

Haplotype (org.broadinstitute.hellbender.utils.haplotype.Haplotype)88 Test (org.testng.annotations.Test)31 GATKRead (org.broadinstitute.hellbender.utils.read.GATKRead)28 BaseTest (org.broadinstitute.hellbender.utils.test.BaseTest)28 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)12 DataProvider (org.testng.annotations.DataProvider)10 VariantContext (htsjdk.variant.variantcontext.VariantContext)9 KBestHaplotype (org.broadinstitute.hellbender.tools.walkers.haplotypecaller.graphs.KBestHaplotype)8 ArrayList (java.util.ArrayList)7 Cigar (htsjdk.samtools.Cigar)6 AssemblyRegion (org.broadinstitute.hellbender.engine.AssemblyRegion)6 HomogeneousPloidyModel (org.broadinstitute.hellbender.tools.walkers.genotyper.HomogeneousPloidyModel)6 IndependentSampleGenotypesModel (org.broadinstitute.hellbender.tools.walkers.genotyper.IndependentSampleGenotypesModel)6 PloidyModel (org.broadinstitute.hellbender.tools.walkers.genotyper.PloidyModel)6 ReadThreadingGraph (org.broadinstitute.hellbender.tools.walkers.haplotypecaller.readthreading.ReadThreadingGraph)6 SampleList (org.broadinstitute.hellbender.utils.genotyper.SampleList)6 EventMap (org.broadinstitute.hellbender.utils.haplotype.EventMap)6 CigarElement (htsjdk.samtools.CigarElement)4 Allele (htsjdk.variant.variantcontext.Allele)4 File (java.io.File)4