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

Example 41 with Allele

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

the class VariantDataManager method writeOutRecalibrationTable.

public void writeOutRecalibrationTable(final VariantContextWriter recalWriter, final SAMSequenceDictionary seqDictionary) {
    // we need to sort in coordinate order in order to produce a valid VCF
    Collections.sort(data, VariantDatum.getComparator(seqDictionary));
    // create dummy alleles to be used
    List<Allele> alleles = Arrays.asList(Allele.create("N", true), Allele.create("<VQSR>", false));
    for (final VariantDatum datum : data) {
        if (VRAC.useASannotations)
            //use the alleles to distinguish between multiallelics in AS mode
            alleles = Arrays.asList(datum.referenceAllele, datum.alternateAllele);
        VariantContextBuilder builder = new VariantContextBuilder("VQSR", datum.loc.getContig(), datum.loc.getStart(), datum.loc.getEnd(), alleles);
        builder.attribute(VCFConstants.END_KEY, datum.loc.getEnd());
        builder.attribute(GATKVCFConstants.VQS_LOD_KEY, String.format("%.4f", datum.lod));
        builder.attribute(GATKVCFConstants.CULPRIT_KEY, (datum.worstAnnotation != -1 ? annotationKeys.get(datum.worstAnnotation) : "NULL"));
        if (datum.atTrainingSite)
            builder.attribute(GATKVCFConstants.POSITIVE_LABEL_KEY, true);
        if (datum.atAntiTrainingSite)
            builder.attribute(GATKVCFConstants.NEGATIVE_LABEL_KEY, true);
        recalWriter.add(builder.make());
    }
}

Example 42 with Allele

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

the class GenotypeUtils method computeDiploidGenotypeCounts.

/**
     * Returns a triple of ref/het/hom genotype "counts".
     *
     * The exact meaning of the count is dependent on the rounding behavior.
     * if {@code roundContributionFromEachGenotype}: the counts are discrete integer counts of the most probable genotype for each {@link Genotype}
     * else: they are the sum of the normalized likelihoods of each genotype and will not be integers
     *
     * Skips non-diploid genotypes.
     *
     *
     * @param vc the VariantContext that the {@link Genotype}s originated from, non-null
     * @param genotypes a GenotypesContext containing genotypes to count, these must be a subset of {@code vc.getGenotypes()}, non-null
     * @param roundContributionFromEachGenotype if this is true, the normalized likelihood from each genotype will be rounded before
     *                                          adding to the total count
     */
public static GenotypeCounts computeDiploidGenotypeCounts(final VariantContext vc, final GenotypesContext genotypes, final boolean roundContributionFromEachGenotype) {
    Utils.nonNull(vc, "vc");
    Utils.nonNull(genotypes, "genotypes");
    final boolean doMultiallelicMapping = !vc.isBiallelic();
    int idxAA = 0, idxAB = 1, idxBB = 2;
    double refCount = 0;
    double hetCount = 0;
    double homCount = 0;
    for (final Genotype g : genotypes) {
        if (!isDiploidWithLikelihoods(g)) {
            continue;
        }
        // Genotype::getLikelihoods returns a new array, so modification in-place is safe
        final double[] normalizedLikelihoods = MathUtils.normalizeFromLog10ToLinearSpace(g.getLikelihoods().getAsVector());
        if (doMultiallelicMapping) {
            if (g.isHetNonRef()) {
                //all likelihoods go to homCount
                homCount++;
                continue;
            }
            //get alternate allele for each sample
            final Allele a1 = g.getAllele(0);
            final Allele a2 = g.getAllele(1);
            if (a2.isNonReference()) {
                final int[] idxVector = vc.getGLIndecesOfAlternateAllele(a2);
                idxAA = idxVector[0];
                idxAB = idxVector[1];
                idxBB = idxVector[2];
            } else //I expect hets to be reference first, but there are no guarantees (e.g. phasing)
            if (a1.isNonReference()) {
                final int[] idxVector = vc.getGLIndecesOfAlternateAllele(a1);
                idxAA = idxVector[0];
                idxAB = idxVector[1];
                idxBB = idxVector[2];
            }
        }
        if (roundContributionFromEachGenotype) {
            refCount += MathUtils.fastRound(normalizedLikelihoods[idxAA]);
            hetCount += MathUtils.fastRound(normalizedLikelihoods[idxAB]);
            homCount += MathUtils.fastRound(normalizedLikelihoods[idxBB]);
        } else {
            refCount += normalizedLikelihoods[idxAA];
            hetCount += normalizedLikelihoods[idxAB];
            homCount += normalizedLikelihoods[idxBB];
        }
    }
    return new GenotypeCounts(refCount, hetCount, homCount);
}

Example 43 with Allele

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

the class CoverageUnitTest method testLikelihoods.

@Test
public void testLikelihoods() {
    final Allele REF = Allele.create("T", true);
    final Allele ALT = Allele.create("A");
    final int refDepth = 5;
    final int altDepth = 3;
    final List<GATKRead> refReads = IntStream.range(0, refDepth).mapToObj(n -> makeRead()).collect(Collectors.toList());
    final List<GATKRead> altReads = IntStream.range(0, altDepth).mapToObj(n -> makeRead()).collect(Collectors.toList());
    final ReadLikelihoods<Allele> likelihoods = AnnotationArtificialData.makeLikelihoods("sample1", refReads, altReads, -100.0, -100.0, REF, ALT);
    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext = null;
    final Map<String, Object> annotate = new Coverage().annotate(referenceContext, vc, likelihoods);
    Assert.assertEquals(annotate.size(), 1, "size");
    Assert.assertEquals(annotate.keySet(), Collections.singleton(VCFConstants.DEPTH_KEY), "annots");
    final int m = altDepth + refDepth;
    Assert.assertEquals(annotate.get(VCFConstants.DEPTH_KEY), String.valueOf(m));
}

Example 44 with Allele

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

the class CoverageUnitTest method makeVC.

private VariantContext makeVC() {
    final GenotypesContext testGC = GenotypesContext.create(2);
    final Allele refAllele = Allele.create("A", true);
    final Allele altAllele = Allele.create("T");
    return (new VariantContextBuilder()).alleles(Arrays.asList(refAllele, altAllele)).chr("1").start(15L).stop(15L).genotypes(testGC).make();
}

Example 45 with Allele

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

the class MappingQualityZeroUnitTest method testLikelihoods.

@Test
public void testLikelihoods() {
    final Allele REF = Allele.create("T", true);
    final Allele ALT = Allele.create("A");
    final int refDepth = 5;
    final int altDepth = 3;
    final int refMQ = 10;
    final int altMQ = 0;
    final List<GATKRead> refReads = IntStream.range(0, refDepth).mapToObj(i -> makeRead(refMQ)).collect(Collectors.toList());
    final List<GATKRead> altReads = IntStream.range(0, altDepth).mapToObj(i -> makeRead(altMQ)).collect(Collectors.toList());
    final ReadLikelihoods<Allele> likelihoods = AnnotationArtificialData.makeLikelihoods("sample1", refReads, altReads, -1.0, -1.0, REF, ALT);
    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext = null;
    final Map<String, Object> annotate = new MappingQualityZero().annotate(referenceContext, vc, likelihoods);
    Assert.assertEquals(annotate.size(), 1, "size");
    Assert.assertEquals(annotate.keySet(), Collections.singleton(VCFConstants.MAPPING_QUALITY_ZERO_KEY), "annots");
    Assert.assertEquals(annotate.get(VCFConstants.MAPPING_QUALITY_ZERO_KEY), String.valueOf(altDepth));
}