Search in sources :

Example 41 with Genotype

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

the class OrientationBiasFiltererUnitTest method testAnnotateVariantContextWithPreprocessingValues.

@Test
public void testAnnotateVariantContextWithPreprocessingValues() {
    final FeatureDataSource<VariantContext> featureDataSource = new FeatureDataSource<>(new File(smallM2Vcf));
    SortedSet<Transition> relevantTransitions = new TreeSet<>();
    relevantTransitions.add(Transition.transitionOf('G', 'T'));
    final Map<Transition, Double> preAdapterQFakeScoreMap = new HashMap<>();
    final double amGTPreAdapterQ = 20.0;
    // preAdapterQ suppression will do nothing.
    preAdapterQFakeScoreMap.put(Transition.transitionOf('G', 'T'), amGTPreAdapterQ);
    for (final VariantContext vc : featureDataSource) {
        final VariantContext updatedVariantContext = OrientationBiasFilterer.annotateVariantContextWithPreprocessingValues(vc, relevantTransitions, preAdapterQFakeScoreMap);
        final Genotype genotypeTumor = updatedVariantContext.getGenotype("TUMOR");
        final Genotype genotypeNormal = updatedVariantContext.getGenotype("NORMAL");
        Assert.assertNotEquals(genotypeTumor.getExtendedAttribute(OrientationBiasFilterConstants.FOB, VCFConstants.EMPTY_ALLELE), VCFConstants.EMPTY_ALLELE);
        Assert.assertNotEquals(genotypeTumor.getExtendedAttribute(OrientationBiasFilterConstants.P_ARTIFACT_FIELD_NAME, VCFConstants.EMPTY_ALLELE), VCFConstants.EMPTY_ALLELE);
        assertArtifact(amGTPreAdapterQ, genotypeTumor, Transition.transitionOf('G', 'T'));
        // The NORMAL is always ref/ref in the example file.
        assertNormal(genotypeNormal);
    }
}
Also used : Transition(org.broadinstitute.hellbender.tools.picard.analysis.artifacts.Transition) VariantContext(htsjdk.variant.variantcontext.VariantContext) Genotype(htsjdk.variant.variantcontext.Genotype) File(java.io.File) FeatureDataSource(org.broadinstitute.hellbender.engine.FeatureDataSource) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 42 with Genotype

use of htsjdk.variant.variantcontext.Genotype 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 43 with Genotype

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

the class FilterByOrientationBiasIntegrationTest method testRun.

@Test
public void testRun() throws IOException {
    final File outputFile = File.createTempFile("ob_", ".vcf");
    final List<String> arguments = new ArrayList<>();
    arguments.add("-" + FilterByOrientationBias.PRE_ADAPTER_METRICS_DETAIL_FILE_SHORT_NAME);
    arguments.add(preAdapterQFile);
    arguments.add("-" + StandardArgumentDefinitions.VARIANT_SHORT_NAME);
    arguments.add(smallM2VcfMore);
    arguments.add("-" + StandardArgumentDefinitions.OUTPUT_SHORT_NAME);
    arguments.add(outputFile.getAbsolutePath());
    runCommandLine(arguments);
    Assert.assertTrue(outputFile.exists());
    final File summaryFile = new File(outputFile.getAbsolutePath() + FilterByOrientationBias.SUMMARY_FILE_SUFFIX);
    Assert.assertTrue(summaryFile.exists());
    final List<VariantContext> variantContexts = new ArrayList<>();
    final FeatureDataSource<VariantContext> featureDataSource = new FeatureDataSource<>(outputFile);
    for (final VariantContext vc : featureDataSource) {
        variantContexts.add(vc);
    }
    Assert.assertEquals(variantContexts.size(), 11);
    Assert.assertTrue(FileUtils.sizeOf(outputFile) > 0);
    Assert.assertTrue(FileUtils.sizeOf(summaryFile) > 0);
    boolean is_variant_context_tested = false;
    //  Also, make sure that the variant context has the filter as well.  Not just the genotypes.
    for (final VariantContext vc : variantContexts) {
        final Genotype tumorGenotype = vc.getGenotype("TUMOR");
        Assert.assertTrue((tumorGenotype.getFilters() == null) || (tumorGenotype.getFilters().contains(OrientationBiasFilterConstants.IS_ORIENTATION_BIAS_CUT)) || !OrientationBiasUtils.isGenotypeInTransitionWithComplement(tumorGenotype, Transition.transitionOf('G', 'T')));
        // If we see a filtered genotype, make sure the variant context was filtered as well.
        if ((tumorGenotype.getFilters() != null) && (tumorGenotype.getFilters().contains(OrientationBiasFilterConstants.IS_ORIENTATION_BIAS_CUT))) {
            Assert.assertTrue(vc.getFilters().contains(OrientationBiasFilterConstants.IS_ORIENTATION_BIAS_CUT));
            is_variant_context_tested = true;
        }
        final Genotype normalGenotype = vc.getGenotype("NORMAL");
        Assert.assertTrue((normalGenotype.getFilters() == null) || normalGenotype.getFilters().equals(VCFConstants.UNFILTERED) || normalGenotype.getFilters().equals(VCFConstants.PASSES_FILTERS_v4));
    }
    Assert.assertTrue(is_variant_context_tested, "Unit test may be broken.  Should have tested that variant context contained filter as well as genotype fields.");
    final List<OrientationSampleTransitionSummary> summaries = OrientationBiasUtils.readOrientationBiasSummaryTable(summaryFile);
    Assert.assertEquals(summaries.size(), 2);
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("NORMAL")).count(), 1);
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("TUMOR")).count(), 1);
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("NORMAL")).map(s -> s.getArtifactMode()).filter(am -> am.equals(Transition.GtoT)).count(), 1);
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("NORMAL")).map(s -> s.getArtifactModeComplement()).filter(am -> am.equals(Transition.CtoA)).count(), 1);
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("TUMOR")).map(s -> s.getArtifactMode()).filter(am -> am.equals(Transition.GtoT)).count(), 1);
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("TUMOR")).map(s -> s.getArtifactModeComplement()).filter(am -> am.equals(Transition.CtoA)).count(), 1);
    Assert.assertEquals(summaries.stream().filter(s -> s.getArtifactModeComplement().equals(s.getArtifactMode().complement())).count(), summaries.size());
    Assert.assertEquals(summaries.stream().filter(s -> s.getSample().equals("TUMOR")).map(s -> s.getArtifactModeComplement()).filter(am -> am.equals(Transition.CtoA)).count(), 1);
}
Also used : Genotype(htsjdk.variant.variantcontext.Genotype) FileUtils(org.apache.commons.io.FileUtils) StandardArgumentDefinitions(org.broadinstitute.hellbender.cmdline.StandardArgumentDefinitions) Test(org.testng.annotations.Test) IOException(java.io.IOException) CommandLineProgramTest(org.broadinstitute.hellbender.CommandLineProgramTest) OrientationBiasFilterConstants(org.broadinstitute.hellbender.tools.exome.orientationbiasvariantfilter.OrientationBiasFilterConstants) File(java.io.File) OrientationSampleTransitionSummary(org.broadinstitute.hellbender.tools.exome.orientationbiasvariantfilter.OrientationSampleTransitionSummary) ArrayList(java.util.ArrayList) List(java.util.List) Assert(org.testng.Assert) FeatureDataSource(org.broadinstitute.hellbender.engine.FeatureDataSource) VariantContext(htsjdk.variant.variantcontext.VariantContext) Transition(org.broadinstitute.hellbender.tools.picard.analysis.artifacts.Transition) OrientationBiasUtils(org.broadinstitute.hellbender.tools.exome.orientationbiasvariantfilter.OrientationBiasUtils) VCFConstants(htsjdk.variant.vcf.VCFConstants) ArrayList(java.util.ArrayList) VariantContext(htsjdk.variant.variantcontext.VariantContext) Genotype(htsjdk.variant.variantcontext.Genotype) OrientationSampleTransitionSummary(org.broadinstitute.hellbender.tools.exome.orientationbiasvariantfilter.OrientationSampleTransitionSummary) File(java.io.File) FeatureDataSource(org.broadinstitute.hellbender.engine.FeatureDataSource) Test(org.testng.annotations.Test) CommandLineProgramTest(org.broadinstitute.hellbender.CommandLineProgramTest)

Example 44 with Genotype

use of htsjdk.variant.variantcontext.Genotype in project hmftools by hartwigmedical.

the class StrelkaPostProcess method simplifyVariant.

@NotNull
static VariantContext simplifyVariant(@NotNull final VariantContext variant, @NotNull final String sampleName) {
    // MIVO: force GT to 0/1 even for variants with multiple alts
    final List<Allele> outputVariantAlleles = variant.getAlleles().subList(0, 2);
    final Genotype genotype = new GenotypeBuilder(sampleName, outputVariantAlleles).DP(getDP(variant)).AD(getAD(variant)).make();
    return new VariantContextBuilder(variant).genotypes(genotype).make();
}
Also used : Allele(htsjdk.variant.variantcontext.Allele) VariantContextBuilder(htsjdk.variant.variantcontext.VariantContextBuilder) Genotype(htsjdk.variant.variantcontext.Genotype) GenotypeBuilder(htsjdk.variant.variantcontext.GenotypeBuilder) NotNull(org.jetbrains.annotations.NotNull)

Example 45 with Genotype

use of htsjdk.variant.variantcontext.Genotype in project hmftools by hartwigmedical.

the class StrelkaPostProcess method readAf.

private static double readAf(@NotNull final VariantContext variant, final int tierIndex, @NotNull Function<Allele, String> alleleKey) {
    final Genotype tumorGenotype = variant.getGenotype(TUMOR_GENOTYPE);
    final String altField = tumorGenotype.getExtendedAttribute(alleleKey.apply(variant.getAlternateAllele(0))).toString();
    final String refField = tumorGenotype.getExtendedAttribute(alleleKey.apply(variant.getReference())).toString();
    final String[] altFieldParts = altField.split(SEPARATOR);
    final String[] refFieldParts = refField.split(SEPARATOR);
    final double altReads = Double.parseDouble(altFieldParts[tierIndex]);
    final double refReads = Double.parseDouble(refFieldParts[tierIndex]);
    final double total = altReads + refReads;
    if (total == 0) {
        return 0;
    }
    return altReads / total;
}
Also used : Genotype(htsjdk.variant.variantcontext.Genotype)

Aggregations

Genotype (htsjdk.variant.variantcontext.Genotype)150 VariantContext (htsjdk.variant.variantcontext.VariantContext)97 Allele (htsjdk.variant.variantcontext.Allele)82 ArrayList (java.util.ArrayList)54 VariantContextBuilder (htsjdk.variant.variantcontext.VariantContextBuilder)52 GenotypeBuilder (htsjdk.variant.variantcontext.GenotypeBuilder)51 File (java.io.File)48 VCFHeader (htsjdk.variant.vcf.VCFHeader)46 IOException (java.io.IOException)45 Collectors (java.util.stream.Collectors)42 Test (org.testng.annotations.Test)37 HashSet (java.util.HashSet)35 VariantContextWriter (htsjdk.variant.variantcontext.writer.VariantContextWriter)29 List (java.util.List)29 VCFHeaderLine (htsjdk.variant.vcf.VCFHeaderLine)27 VCFFormatHeaderLine (htsjdk.variant.vcf.VCFFormatHeaderLine)25 SAMSequenceDictionaryProgress (com.github.lindenb.jvarkit.util.picard.SAMSequenceDictionaryProgress)23 HashMap (java.util.HashMap)23 java.util (java.util)22 Set (java.util.Set)22