Search in sources :

Example 1 with ReadLikelihoods

use of org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods in project gatk by broadinstitute.

the class RMSMappingQualityUnitTest method testLikelihoods.

@Test
public void testLikelihoods() {
    final Allele REF = Allele.create("A", true);
    final Allele ALT = Allele.create("C", true);
    final int[] MQs = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, QualityUtils.MAPPING_QUALITY_UNAVAILABLE };
    final List<Integer> MQsList = Arrays.asList(ArrayUtils.toObject(MQs));
    //MQ 255 are excluded from the calculations, we test it here.
    final List<Integer> MQsListOK = new ArrayList<>(MQsList);
    //NOTE: if we just call remove(i), Java thinks i is an index.
    //A workaround for this overloading bogosity to to call removeAll and pass a collection
    //(casting i to (Object) would work too but it's more error prone)
    MQsListOK.removeAll(Collections.singleton(QualityUtils.MAPPING_QUALITY_UNAVAILABLE));
    final List<GATKRead> reads = Arrays.stream(MQs).mapToObj(mq -> AnnotationArtificialData.makeRead(20, mq)).collect(Collectors.toList());
    final ReadLikelihoods<Allele> likelihoods = AnnotationArtificialData.makeLikelihoods("sample1", reads, -1.0, REF, ALT);
    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext = null;
    final Map<String, Object> annotate = new RMSMappingQuality().annotate(referenceContext, vc, likelihoods);
    Assert.assertEquals(annotate.size(), 1, "size");
    Assert.assertEquals(annotate.keySet(), Collections.singleton(GATKVCFConstants.RAW_RMS_MAPPING_QUALITY_KEY), "annots");
    //only those are MQ0
    final double rms = MathUtils.sumOfSquares(MQsListOK);
    Assert.assertNull(annotate.get(VCFConstants.RMS_MAPPING_QUALITY_KEY));
    Assert.assertEquals(annotate.get(GATKVCFConstants.RAW_RMS_MAPPING_QUALITY_KEY), String.format("%.2f", rms));
}
Also used : java.util(java.util) AS_RMSMappingQuality(org.broadinstitute.hellbender.tools.walkers.annotator.allelespecific.AS_RMSMappingQuality) ImmutableMap(com.google.common.collect.ImmutableMap) GATKVCFConstants(org.broadinstitute.hellbender.utils.variant.GATKVCFConstants) ArrayUtils(org.apache.commons.lang3.ArrayUtils) QualityUtils(org.broadinstitute.hellbender.utils.QualityUtils) Test(org.testng.annotations.Test) htsjdk.variant.variantcontext(htsjdk.variant.variantcontext) GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) Collectors(java.util.stream.Collectors) Sets(com.google.cloud.dataflow.sdk.repackaged.com.google.common.collect.Sets) ReadLikelihoods(org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods) MathUtils(org.broadinstitute.hellbender.utils.MathUtils) UserException(org.broadinstitute.hellbender.exceptions.UserException) Assert(org.testng.Assert) AlleleList(org.broadinstitute.hellbender.utils.genotyper.AlleleList) IndexedAlleleList(org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) VCFConstants(htsjdk.variant.vcf.VCFConstants) GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) AS_RMSMappingQuality(org.broadinstitute.hellbender.tools.walkers.annotator.allelespecific.AS_RMSMappingQuality) Test(org.testng.annotations.Test)

Example 2 with ReadLikelihoods

use of org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods in project gatk by broadinstitute.

the class RMSMappingQualityUnitTest method testLikelihoodsEmpty_AS.

@Test
public void testLikelihoodsEmpty_AS() throws Exception {
    final List<GATKRead> reads = Collections.emptyList();
    final Map<String, List<GATKRead>> readsBySample = ImmutableMap.of("sample1", reads);
    final org.broadinstitute.hellbender.utils.genotyper.SampleList sampleList = new IndexedSampleList(Arrays.asList("sample1"));
    final AlleleList<Allele> alleleList = new IndexedAlleleList<>(Arrays.asList(Allele.create("A")));
    final ReadLikelihoods<Allele> likelihoods = new ReadLikelihoods<>(sampleList, alleleList, readsBySample);
    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext = null;
    final Map<String, Object> annotate = new AS_RMSMappingQuality().annotate(referenceContext, vc, likelihoods);
    final String[] split = ((String) annotate.get(GATKVCFConstants.AS_RAW_RMS_MAPPING_QUALITY_KEY)).split(AS_RMSMappingQuality.SPLIT_DELIM);
    Assert.assertEquals(split.length, 2);
    Assert.assertEquals(split[0], String.format("%.2f", 0.0));
    Assert.assertEquals(split[1], String.format("%.2f", 0.0));
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) IndexedAlleleList(org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList) ReadLikelihoods(org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods) AS_RMSMappingQuality(org.broadinstitute.hellbender.tools.walkers.annotator.allelespecific.AS_RMSMappingQuality) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) AlleleList(org.broadinstitute.hellbender.utils.genotyper.AlleleList) IndexedAlleleList(org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) Test(org.testng.annotations.Test)

Example 3 with ReadLikelihoods

use of org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods in project gatk by broadinstitute.

the class RMSMappingQualityUnitTest method testLikelihoods_AS.

@Test
public void testLikelihoods_AS() {
    final Allele REF = Allele.create("A", true);
    final Allele ALT = Allele.create("C");
    final int[] MQs = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, QualityUtils.MAPPING_QUALITY_UNAVAILABLE };
    final List<Integer> MQsList = Arrays.asList(ArrayUtils.toObject(MQs));
    //MQ 255 are excluded from the calculations, we test it here.
    final List<Integer> MQsListOK = new ArrayList<>(MQsList);
    //NOTE: if we just call remove(i), Java thinks i is an index.
    //A workaround for this overloading bogosity to to call removeAll and pass a collection
    //(casting i to (Object) would work too but it's more error prone)
    MQsListOK.removeAll(Collections.singleton(QualityUtils.MAPPING_QUALITY_UNAVAILABLE));
    final List<GATKRead> reads = Arrays.stream(MQs).mapToObj(mq -> AnnotationArtificialData.makeRead(30, mq)).collect(Collectors.toList());
    final ReadLikelihoods<Allele> likelihoods = AnnotationArtificialData.makeLikelihoods("sample1", reads, -10.0, REF, ALT);
    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext = null;
    final Map<String, Object> annotate = new AS_RMSMappingQuality().annotateRawData(referenceContext, vc, likelihoods);
    Assert.assertEquals(annotate.size(), 1, "size");
    Assert.assertEquals(annotate.keySet(), Collections.singleton(GATKVCFConstants.AS_RAW_RMS_MAPPING_QUALITY_KEY), "annots");
    //only those are MQ0
    final double rms = MathUtils.sumOfSquares(MQsListOK);
    final String[] split = ((String) annotate.get(GATKVCFConstants.AS_RAW_RMS_MAPPING_QUALITY_KEY)).split(AS_RMSMappingQuality.SPLIT_DELIM);
    Assert.assertEquals(split.length, 2);
    Assert.assertEquals(split[0], String.format("%.2f", rms));
    Assert.assertEquals(split[1], String.format("%.2f", 0.0));
}
Also used : java.util(java.util) AS_RMSMappingQuality(org.broadinstitute.hellbender.tools.walkers.annotator.allelespecific.AS_RMSMappingQuality) ImmutableMap(com.google.common.collect.ImmutableMap) GATKVCFConstants(org.broadinstitute.hellbender.utils.variant.GATKVCFConstants) ArrayUtils(org.apache.commons.lang3.ArrayUtils) QualityUtils(org.broadinstitute.hellbender.utils.QualityUtils) Test(org.testng.annotations.Test) htsjdk.variant.variantcontext(htsjdk.variant.variantcontext) GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) Collectors(java.util.stream.Collectors) Sets(com.google.cloud.dataflow.sdk.repackaged.com.google.common.collect.Sets) ReadLikelihoods(org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods) MathUtils(org.broadinstitute.hellbender.utils.MathUtils) UserException(org.broadinstitute.hellbender.exceptions.UserException) Assert(org.testng.Assert) AlleleList(org.broadinstitute.hellbender.utils.genotyper.AlleleList) IndexedAlleleList(org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) VCFConstants(htsjdk.variant.vcf.VCFConstants) GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) AS_RMSMappingQuality(org.broadinstitute.hellbender.tools.walkers.annotator.allelespecific.AS_RMSMappingQuality) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) Test(org.testng.annotations.Test)

Example 4 with ReadLikelihoods

use of org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods in project gatk by broadinstitute.

the class RMSMappingQualityUnitTest method testEmptyLikelihoods.

@Test
public void testEmptyLikelihoods() throws Exception {
    final List<GATKRead> reads = Collections.emptyList();
    final Map<String, List<GATKRead>> readsBySample = ImmutableMap.of("sample1", reads);
    final org.broadinstitute.hellbender.utils.genotyper.SampleList sampleList = new IndexedSampleList(Arrays.asList("sample1"));
    final AlleleList<Allele> alleleList = new IndexedAlleleList<>(Arrays.asList(Allele.create("A")));
    final ReadLikelihoods<Allele> likelihoods = new ReadLikelihoods<>(sampleList, alleleList, readsBySample);
    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext = null;
    final Map<String, Object> annotate = new RMSMappingQuality().annotate(referenceContext, vc, likelihoods);
    Assert.assertTrue(annotate.isEmpty());
}
Also used : GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) IndexedAlleleList(org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList) ReadLikelihoods(org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) AlleleList(org.broadinstitute.hellbender.utils.genotyper.AlleleList) IndexedAlleleList(org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList) IndexedSampleList(org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList) AS_RMSMappingQuality(org.broadinstitute.hellbender.tools.walkers.annotator.allelespecific.AS_RMSMappingQuality) Test(org.testng.annotations.Test)

Example 5 with ReadLikelihoods

use of org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods in project gatk by broadinstitute.

the class AS_StrandBiasTest method getStrandCountsFromLikelihoodMap.

/**
     Allocate and fill a 2x2 strand contingency table.  In the end, it'll look something like this:
     *             fw      rc
     *   allele1   #       #
     *   allele2   #       #
     * @return a 2x2 contingency table
     */
public void getStrandCountsFromLikelihoodMap(final VariantContext vc, final ReadLikelihoods<Allele> likelihoods, final ReducibleAnnotationData<List<Integer>> perAlleleValues, final int minCount) {
    if (likelihoods == null || vc == null) {
        return;
    }
    final Allele ref = vc.getReference();
    final List<Allele> allAlts = vc.getAlternateAlleles();
    for (final String sample : likelihoods.samples()) {
        final ReducibleAnnotationData<List<Integer>> sampleTable = new AlleleSpecificAnnotationData<>(vc.getAlleles(), null);
        likelihoods.bestAlleles(sample).stream().filter(ba -> ba.isInformative()).forEach(ba -> updateTable(ba.allele, ba.read, ref, allAlts, sampleTable));
        if (passesMinimumThreshold(sampleTable, minCount)) {
            combineAttributeMap(sampleTable, perAlleleValues);
        }
    }
}
Also used : Allele(htsjdk.variant.variantcontext.Allele) GenotypesContext(htsjdk.variant.variantcontext.GenotypesContext) Logger(org.apache.log4j.Logger) ReadLikelihoods(org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods) java.util(java.util) StrandBiasTest(org.broadinstitute.hellbender.tools.walkers.annotator.StrandBiasTest) GATKVCFConstants(org.broadinstitute.hellbender.utils.variant.GATKVCFConstants) VCFInfoHeaderLine(htsjdk.variant.vcf.VCFInfoHeaderLine) GATKVCFHeaderLines(org.broadinstitute.hellbender.utils.variant.GATKVCFHeaderLines) VariantContext(htsjdk.variant.variantcontext.VariantContext) ReferenceContext(org.broadinstitute.hellbender.engine.ReferenceContext) GATKRead(org.broadinstitute.hellbender.utils.read.GATKRead) Allele(htsjdk.variant.variantcontext.Allele)

Aggregations

ReadLikelihoods (org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods)22 ReferenceContext (org.broadinstitute.hellbender.engine.ReferenceContext)17 GATKRead (org.broadinstitute.hellbender.utils.read.GATKRead)17 Allele (htsjdk.variant.variantcontext.Allele)11 Test (org.testng.annotations.Test)11 GATKVCFConstants (org.broadinstitute.hellbender.utils.variant.GATKVCFConstants)10 VariantContext (htsjdk.variant.variantcontext.VariantContext)9 java.util (java.util)9 Collectors (java.util.stream.Collectors)9 AlleleList (org.broadinstitute.hellbender.utils.genotyper.AlleleList)9 IndexedAlleleList (org.broadinstitute.hellbender.utils.genotyper.IndexedAlleleList)9 IndexedSampleList (org.broadinstitute.hellbender.utils.genotyper.IndexedSampleList)9 Assert (org.testng.Assert)7 htsjdk.variant.variantcontext (htsjdk.variant.variantcontext)6 VCFConstants (htsjdk.variant.vcf.VCFConstants)6 List (java.util.List)6 GenotypesContext (htsjdk.variant.variantcontext.GenotypesContext)5 ImmutableMap (com.google.common.collect.ImmutableMap)4 TextCigarCodec (htsjdk.samtools.TextCigarCodec)4 Genotype (htsjdk.variant.variantcontext.Genotype)4