Search in sources :

Example 16 with Allele

use of htsjdk.variant.variantcontext.Allele 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 17 with Allele

use of htsjdk.variant.variantcontext.Allele in project gatk-protected by broadinstitute.

the class GermlineProbabilityCalculatorUnitTest method testGetGermlineAltAlleleFrequencies.

@Test
public void testGetGermlineAltAlleleFrequencies() {
    final double defaultAF = 0.001;
    final double nonDefaultAF1 = 0.1;
    final double nonDefaultAF2 = 0.01;
    final Allele Aref = Allele.create("A", true);
    final Allele C = Allele.create("C");
    final Allele G = Allele.create("G");
    final Allele T = Allele.create("T");
    final String source = "SOURCE";
    final int start = 1;
    final int stop = 1;
    //biallelic, vc has the same alt allele
    final List<Allele> altAlleles1 = Arrays.asList(C);
    final VariantContext vc1 = new VariantContextBuilder(source, "1", start, stop, Arrays.asList(Aref, C)).attribute(VCFConstants.ALLELE_FREQUENCY_KEY, new double[] { nonDefaultAF1 }).make();
    final double[] af1 = GermlineProbabilityCalculator.getGermlineAltAlleleFrequencies(altAlleles1, Optional.of(vc1), defaultAF);
    Assert.assertEquals(af1.length, altAlleles1.size());
    Assert.assertEquals(af1[0], nonDefaultAF1, 0.00001);
    //biallelic, vc has different alt allele
    final List<Allele> altAlleles2 = Arrays.asList(C);
    final VariantContext vc2 = new VariantContextBuilder(source, "1", start, stop, Arrays.asList(Aref, G)).attribute(VCFConstants.ALLELE_FREQUENCY_KEY, new double[] { nonDefaultAF1 }).make();
    final double[] af2 = GermlineProbabilityCalculator.getGermlineAltAlleleFrequencies(altAlleles2, Optional.of(vc2), defaultAF);
    Assert.assertEquals(af2.length, altAlleles2.size());
    Assert.assertEquals(af2[0], defaultAF, 0.00001);
    //triallelic, same alt alleles
    final List<Allele> altAlleles3 = Arrays.asList(C, G);
    final VariantContext vc3 = new VariantContextBuilder(source, "1", start, stop, Arrays.asList(Aref, C, G)).attribute(VCFConstants.ALLELE_FREQUENCY_KEY, new double[] { nonDefaultAF1, nonDefaultAF2 }).make();
    final double[] af3 = GermlineProbabilityCalculator.getGermlineAltAlleleFrequencies(altAlleles3, Optional.of(vc3), defaultAF);
    Assert.assertEquals(af3.length, altAlleles3.size());
    Assert.assertEquals(af3[0], nonDefaultAF1, 0.00001);
    Assert.assertEquals(af3[1], nonDefaultAF2, 0.00001);
    //triallelic, same alt alleles in different order
    final List<Allele> altAlleles4 = Arrays.asList(C, G);
    final VariantContext vc4 = new VariantContextBuilder(source, "1", start, stop, Arrays.asList(Aref, G, C)).attribute(VCFConstants.ALLELE_FREQUENCY_KEY, new double[] { nonDefaultAF1, nonDefaultAF2 }).make();
    final double[] af4 = GermlineProbabilityCalculator.getGermlineAltAlleleFrequencies(altAlleles4, Optional.of(vc4), defaultAF);
    Assert.assertEquals(af4.length, altAlleles4.size());
    Assert.assertEquals(af4[0], nonDefaultAF2, 0.00001);
    Assert.assertEquals(af4[1], nonDefaultAF1, 0.00001);
    //triallelic, only one allele in common
    final List<Allele> altAlleles5 = Arrays.asList(C, G);
    final VariantContext vc5 = new VariantContextBuilder(source, "1", start, stop, Arrays.asList(Aref, C, T)).attribute(VCFConstants.ALLELE_FREQUENCY_KEY, new double[] { nonDefaultAF1, nonDefaultAF2 }).make();
    final double[] af5 = GermlineProbabilityCalculator.getGermlineAltAlleleFrequencies(altAlleles5, Optional.of(vc5), defaultAF);
    Assert.assertEquals(af5.length, altAlleles5.size());
    Assert.assertEquals(af5[0], nonDefaultAF1, 0.00001);
    Assert.assertEquals(af5[1], defaultAF, 0.00001);
}
Also used : Allele(htsjdk.variant.variantcontext.Allele) VariantContextBuilder(htsjdk.variant.variantcontext.VariantContextBuilder) VariantContext(htsjdk.variant.variantcontext.VariantContext) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 18 with Allele

use of htsjdk.variant.variantcontext.Allele in project gatk by broadinstitute.

the class EvaluateCopyNumberTriStateCalls method buildAndAnnotateTruthOverlappingGenotype.

private Genotype buildAndAnnotateTruthOverlappingGenotype(final String sample, final TargetCollection<Target> targets, final Genotype truthGenotype, final int truthCopyNumber, final CopyNumberTriStateAllele truthAllele, final List<Pair<VariantContext, Genotype>> calls) {
    final Set<CopyNumberTriStateAllele> calledAlleles = calls.stream().map(pair -> CopyNumberTriStateAllele.valueOf(pair.getRight().getAllele(0))).collect(Collectors.toSet());
    final Allele calledAllele = calledAlleles.size() == 1 ? calledAlleles.iterator().next() : Allele.NO_CALL;
    final GenotypeBuilder builder = new GenotypeBuilder(sample);
    // Set the call allele.
    builder.alleles(Collections.singletonList(calledAllele));
    // Set the truth allele.
    builder.attribute(VariantEvaluationContext.TRUTH_GENOTYPE_KEY, CopyNumberTriStateAllele.ALL_ALLELES.indexOf(truthAllele));
    // Annotate the genotype with the number of calls.
    builder.attribute(VariantEvaluationContext.CALLED_SEGMENTS_COUNT_KEY, calls.size());
    // When there is more than one qualified type of event we indicate how many.
    builder.attribute(VariantEvaluationContext.CALLED_ALLELE_COUNTS_KEY, CopyNumberTriStateAllele.ALL_ALLELES.stream().mapToInt(allele -> (int) calls.stream().filter(pair -> pair.getRight().getAllele(0).equals(allele, true)).count()).toArray());
    // Calculate the length in targets of the call as the sum across all calls.
    builder.attribute(VariantEvaluationContext.CALLED_TARGET_COUNT_KEY, calls.stream().mapToInt(pair -> getTargetCount(targets, pair.getLeft(), pair.getRight())).sum());
    // Calculate call quality-- if there is more than one overlapping call we take the maximum qual one.
    builder.attribute(VariantEvaluationContext.CALL_QUALITY_KEY, calls.stream().mapToDouble(pair -> GATKProtectedVariantContextUtils.calculateGenotypeQualityFromPLs(pair.getRight())).max().orElse(0.0));
    // Calculate the truth copy fraction.
    builder.attribute(VariantEvaluationContext.TRUTH_COPY_FRACTION_KEY, truthGenotype.getExtendedAttribute(GS_COPY_NUMBER_FRACTION_KEY));
    // Calculate the truth call quality.
    final double truthQuality = calculateTruthQuality(truthGenotype, truthCopyNumber);
    builder.attribute(VariantEvaluationContext.TRUTH_QUALITY_KEY, truthQuality);
    // Set genotype filters:
    final boolean truthPassQualityMinimum = truthQuality >= filterArguments.minimumTruthSegmentQuality;
    builder.filter(truthPassQualityMinimum ? EvaluationFilter.PASS : EvaluationFilter.LowQuality.acronym);
    // Calculate the evaluation class (TP, FN, etc.). Only if there is actually either a truth or a call that is not ref.
    if (calledAlleles.contains(CopyNumberTriStateAllele.DEL) || calledAlleles.contains(CopyNumberTriStateAllele.DUP) || truthAllele != CopyNumberTriStateAllele.REF) {
        final EvaluationClass evaluationClass;
        if (calledAlleles.isEmpty() || (calledAlleles.size() == 1 && calledAlleles.contains(CopyNumberTriStateAllele.REF))) {
            evaluationClass = EvaluationClass.FALSE_NEGATIVE;
        } else if (calledAlleles.size() == 1) {
            evaluationClass = calledAlleles.contains(truthAllele) ? EvaluationClass.TRUE_POSITIVE : truthAllele == CopyNumberTriStateAllele.REF ? EvaluationClass.FALSE_POSITIVE : /* else */
            EvaluationClass.DISCORDANT_POSITIVE;
        } else {
            evaluationClass = truthAllele == CopyNumberTriStateAllele.REF ? EvaluationClass.FALSE_POSITIVE : EvaluationClass.MIXED_POSITIVE;
        }
        builder.attribute(VariantEvaluationContext.EVALUATION_CLASS_KEY, evaluationClass.acronym);
    }
    return builder.make();
}
Also used : Genotype(htsjdk.variant.variantcontext.Genotype) DocumentedFeature(org.broadinstitute.barclay.help.DocumentedFeature) Allele(htsjdk.variant.variantcontext.Allele) htsjdk.variant.vcf(htsjdk.variant.vcf) CommandLineProgramProperties(org.broadinstitute.barclay.argparser.CommandLineProgramProperties) java.util(java.util) CopyNumberProgramGroup(org.broadinstitute.hellbender.cmdline.programgroups.CopyNumberProgramGroup) Argument(org.broadinstitute.barclay.argparser.Argument) VariantContextWriterBuilder(htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder) StandardArgumentDefinitions(org.broadinstitute.hellbender.cmdline.StandardArgumentDefinitions) ArgumentCollection(org.broadinstitute.barclay.argparser.ArgumentCollection) TargetArgumentCollection(org.broadinstitute.hellbender.tools.exome.TargetArgumentCollection) Function(java.util.function.Function) Pair(org.apache.commons.lang3.tuple.Pair) StreamSupport(java.util.stream.StreamSupport) TargetCollection(org.broadinstitute.hellbender.tools.exome.TargetCollection) org.broadinstitute.hellbender.utils(org.broadinstitute.hellbender.utils) Locatable(htsjdk.samtools.util.Locatable) CommandLineProgram(org.broadinstitute.hellbender.cmdline.CommandLineProgram) GenotypeBuilder(htsjdk.variant.variantcontext.GenotypeBuilder) IOException(java.io.IOException) CopyNumberTriStateAllele(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriStateAllele) Collectors(java.util.stream.Collectors) ImmutablePair(org.apache.commons.lang3.tuple.ImmutablePair) File(java.io.File) HMMPostProcessor(org.broadinstitute.hellbender.utils.hmm.segmentation.HMMPostProcessor) XHMMSegmentGenotyper(org.broadinstitute.hellbender.tools.exome.germlinehmm.xhmm.XHMMSegmentGenotyper) Stream(java.util.stream.Stream) UserException(org.broadinstitute.hellbender.exceptions.UserException) VariantContextWriter(htsjdk.variant.variantcontext.writer.VariantContextWriter) Target(org.broadinstitute.hellbender.tools.exome.Target) XHMMSegmentCaller(org.broadinstitute.hellbender.tools.exome.germlinehmm.xhmm.XHMMSegmentCaller) VariantContext(htsjdk.variant.variantcontext.VariantContext) BufferedReader(java.io.BufferedReader) FileReader(java.io.FileReader) CopyNumberTriStateAllele(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriStateAllele) Allele(htsjdk.variant.variantcontext.Allele) CopyNumberTriStateAllele(org.broadinstitute.hellbender.tools.exome.germlinehmm.CopyNumberTriStateAllele) GenotypeBuilder(htsjdk.variant.variantcontext.GenotypeBuilder)

Example 19 with Allele

use of htsjdk.variant.variantcontext.Allele in project gatk by broadinstitute.

the class AlleleSpecificAnnotationDataUnitTest method testTwoRef.

@Test(expectedExceptions = IllegalArgumentException.class)
public void testTwoRef() throws Exception {
    final Allele Aref = Allele.create("A", true);
    final Allele Tref = Allele.create("T", true);
    final List<Allele> alleles = Arrays.asList(Aref, Tref);
    String rawData = "1|2";
    new AlleleSpecificAnnotationData<>(alleles, rawData);
}
Also used : Allele(htsjdk.variant.variantcontext.Allele) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 20 with Allele

use of htsjdk.variant.variantcontext.Allele in project gatk by broadinstitute.

the class AlleleSpecificAnnotationDataUnitTest method testCreate.

@Test
public void testCreate() throws Exception {
    final Allele Aref = Allele.create("A", true);
    final Allele T = Allele.create("T", false);
    final List<Allele> alleles = Arrays.asList(Aref, T);
    String rawData = "1|2";
    final AlleleSpecificAnnotationData<Integer> asad = new AlleleSpecificAnnotationData<>(alleles, rawData);
    Assert.assertEquals(asad.getAlleles(), alleles);
    Assert.assertEquals(asad.getRefAllele(), Aref);
    Assert.assertNull(asad.getAttribute(Aref));
    Assert.assertNull(asad.getAttribute(T));
    Assert.assertEquals(asad.getRawData(), rawData);
    final Map<Allele, Integer> map = new HashMap<>();
    map.put(Aref, 10);
    map.put(T, 11);
    asad.setAttributeMap(map);
    Assert.assertEquals(asad.getAttribute(Aref), (Integer) 10);
    Assert.assertEquals(asad.getAttribute(T), (Integer) 11);
}
Also used : Allele(htsjdk.variant.variantcontext.Allele) HashMap(java.util.HashMap) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Aggregations

Allele (htsjdk.variant.variantcontext.Allele)91 Test (org.testng.annotations.Test)48 VariantContext (htsjdk.variant.variantcontext.VariantContext)44 VariantContextBuilder (htsjdk.variant.variantcontext.VariantContextBuilder)27 Genotype (htsjdk.variant.variantcontext.Genotype)26 BaseTest (org.broadinstitute.hellbender.utils.test.BaseTest)22 java.util (java.util)19 Collectors (java.util.stream.Collectors)16 IntStream (java.util.stream.IntStream)14 ReferenceContext (org.broadinstitute.hellbender.engine.ReferenceContext)13 CommandLineProgramTest (org.broadinstitute.hellbender.CommandLineProgramTest)12 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)12 File (java.io.File)11 StandardArgumentDefinitions (org.broadinstitute.hellbender.cmdline.StandardArgumentDefinitions)11 ReadLikelihoods (org.broadinstitute.hellbender.utils.genotyper.ReadLikelihoods)11 GATKRead (org.broadinstitute.hellbender.utils.read.GATKRead)11 GenotypeBuilder (htsjdk.variant.variantcontext.GenotypeBuilder)10 VCFConstants (htsjdk.variant.vcf.VCFConstants)10 IOException (java.io.IOException)10 Target (org.broadinstitute.hellbender.tools.exome.Target)10