Examples with VariantContext htsjdk.variant.variantcontext.VariantContext used on opensource projects

Example 56 with VariantContext

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

the class GermlineProbabilityCalculator method calculateAnnotations.

public static Map<String, Object> calculateAnnotations(List<VariantContext> germlineResourceVariants, final List<Allele> altAlleles, final double[] tumorLog10Odds, final Optional<double[]> normalLog10Odds, final double afOfAllelesNotInGermlineResource, final double log10PriorProbOfSomaticEvent) {
    final double[] normalLog10OddsOrFlat = normalLog10Odds.orElseGet(() -> MathUtils.applyToArray(tumorLog10Odds, x -> 0));
    final Optional<VariantContext> germlineVC = germlineResourceVariants.isEmpty() ? Optional.empty() : // assume only one VC per site
    Optional.of(germlineResourceVariants.get(0));
    final double[] populationAlleleFrequencies = getGermlineAltAlleleFrequencies(altAlleles, germlineVC, afOfAllelesNotInGermlineResource);
    // note the minus sign required because Mutect has the convention that this is log odds of allele *NOT* being in the normal
    final double[] germlineLog10Posteriors = new IndexRange(0, altAlleles.size()).mapToDouble(n -> log10PosteriorProbabilityOfGermlineVariant(-normalLog10OddsOrFlat[n], tumorLog10Odds[n], populationAlleleFrequencies[n], log10PriorProbOfSomaticEvent));
    return ImmutableMap.of(POPULATION_AF_VCF_ATTRIBUTE, populationAlleleFrequencies, GERMLINE_POSTERIORS_VCF_ATTRIBUTE, germlineLog10Posteriors);
}

Example 57 with VariantContext

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

the class CalculateMixingFractions method apply.

@Override
public void apply(final VariantContext vc, final ReadsContext readsContext, final ReferenceContext refContext, final FeatureContext fc) {
    if (!isBiallelicSingletonHetSnp(vc)) {
        return;
    }
    final Optional<String> variantSample = StreamSupport.stream(vc.getGenotypes().spliterator(), false).filter(genotype -> genotype.isHet()).map(genotype -> genotype.getSampleName()).findFirst();
    if (!variantSample.isPresent()) {
        return;
    }
    final List<GATKRead> reads = new ArrayList<>();
    final List<Integer> offsets = new ArrayList<>();
    for (final GATKRead read : readsContext) {
        if (read.failsVendorQualityCheck()) {
            continue;
        }
        final AlignmentStateMachine asm = new AlignmentStateMachine(read);
        while (asm.stepForwardOnGenome() != null && asm.getGenomePosition() < vc.getStart()) {
        }
        if (asm.getGenomePosition() == vc.getStart()) {
            reads.add(read);
            offsets.add(asm.getReadOffset());
        }
    }
    final ReadPileup pileup = new ReadPileup(vc, reads, offsets);
    final byte altBase = vc.getAlternateAllele(0).getBases()[0];
    final long altCount = StreamSupport.stream(pileup.spliterator(), false).filter(pe -> pe.getBase() == altBase).count();
    final long totalCount = pileup.size();
    sampleCounts.get(variantSample.get()).addCounts(altCount, totalCount);
}

Example 58 with VariantContext

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

the class ReadThreadingAssembler method composeGivenHaplotypes.

/**
     * Create the list of artificial GGA-mode haplotypes by injecting each of the provided alternate alleles into the reference haplotype
     *
     * @param refHaplotype the reference haplotype
     * @param givenHaplotypes the list of alternate alleles in VariantContexts
     * @param activeRegionWindow the window containing the reference haplotype
     *
     * @return a non-null list of haplotypes
     */
private static List<Haplotype> composeGivenHaplotypes(final Haplotype refHaplotype, final List<VariantContext> givenHaplotypes, final SimpleInterval activeRegionWindow) {
    Utils.nonNull(refHaplotype, "the reference haplotype cannot be null");
    Utils.nonNull(givenHaplotypes, "given haplotypes cannot be null");
    Utils.nonNull(activeRegionWindow, "active region window cannot be null");
    Utils.validateArg(activeRegionWindow.size() == refHaplotype.length(), "inconsistent reference haplotype and active region window");
    final Set<Haplotype> returnHaplotypes = new LinkedHashSet<>();
    final int activeRegionStart = refHaplotype.getAlignmentStartHapwrtRef();
    for (final VariantContext compVC : givenHaplotypes) {
        Utils.validateArg(GATKVariantContextUtils.overlapsRegion(compVC, activeRegionWindow), " some variant provided does not overlap with active region window");
        for (final Allele compAltAllele : compVC.getAlternateAlleles()) {
            final Haplotype insertedRefHaplotype = refHaplotype.insertAllele(compVC.getReference(), compAltAllele, activeRegionStart + compVC.getStart() - activeRegionWindow.getStart(), compVC.getStart());
            if (insertedRefHaplotype != null) {
                // can be null if the requested allele can't be inserted into the haplotype
                returnHaplotypes.add(insertedRefHaplotype);
            }
        }
    }
    return new ArrayList<>(returnHaplotypes);
}

Example 59 with VariantContext

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

the class ASEReadCounter method apply.

@Override
public void apply(AlignmentContext alignmentContext, ReferenceContext referenceContext, FeatureContext featureContext) {
    final String contig = alignmentContext.getContig();
    final long position = alignmentContext.getPosition();
    final char refAllele = (char) referenceContext.getBase();
    final List<VariantContext> VCs = featureContext.getValues(variants);
    if (VCs != null && VCs.size() > 1) {
        throw new UserException("More then one variant context at position: " + contig + ":" + position);
    }
    if (VCs == null || VCs.isEmpty()) {
        return;
    }
    final VariantContext vc = VCs.get(0);
    if (!vc.isBiallelic()) {
        logger.warn("Ignoring site: cannot run ASE on non-biallelic sites: " + vc.toString());
        return;
    }
    if (vc.getHetCount() < 1) {
        logger.warn("Ignoring site: variant is not het at postion: " + contig + ":" + position);
        return;
    }
    if (vc.getNAlleles() == 1 || vc.getAlternateAllele(0).getBases().length == 0) {
        throw new UserException("The file of variant sites must contain heterozygous sites and cannot be a GVCF file containing <NON_REF> alleles.");
    }
    final char altAllele = (char) vc.getAlternateAllele(0).getBases()[0];
    final String siteID = vc.getID();
    final ReadPileup pileup = filterPileup(alignmentContext.getBasePileup(), countType);
    // count up the depths of all and QC+ bases
    final String line = calculateLineForSite(pileup, siteID, refAllele, altAllele);
    if (line != null) {
        outputStream.println(line);
    }
}

Example 60 with VariantContext

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

the class FeatureDataSourceUnitTest method testGetSequenceDictionary.

@Test
public void testGetSequenceDictionary() {
    try (FeatureDataSource<VariantContext> featureSource = new FeatureDataSource<>(QUERY_TEST_VCF, "CustomName")) {
        final SAMSequenceDictionary dict = featureSource.getSequenceDictionary();
        Assert.assertEquals(dict.size(), 4);
        Assert.assertEquals(dict.getSequences().stream().map(s -> s.getSequenceName()).collect(Collectors.toList()), Arrays.asList("1", "2", "3", "4"));
    }
}