Java Genotype.isNonInformative方法代码示例

import htsjdk.variant.variantcontext.Genotype; //导入方法依赖的package包/类
/**
 * Returns a new Genotype with the PLs of the multi-allelic original reduced to a bi-allelic case
 *
 * This is handled in the following way:
 *
 * Suppose we have for a A/B/C site the following GLs:
 *
 * AA AB BB AC BC CC
 *
 * and we want to get the bi-allelic GLs for X/B, where X is everything not B
 *
 * XX = AA + AC + CC (since X = A or C)
 * XB = AB + BC
 * BB = BB
 *
 * @param original the original multi-allelic genotype
 * @param altIndex the index of the alt allele we wish to keep in the bialleic case -- with ref == 0
 * @param nAlts the total number of alt alleles
 * @return a new biallelic genotype with appropriate PLs
 */
@Deprecated
protected Genotype combineGLs(final Genotype original, final int altIndex, final int nAlts ) {
    if ( original.isNonInformative() )
        return new GenotypeBuilder(original).PL(BIALLELIC_NON_INFORMATIVE_PLS).alleles(BIALLELIC_NOCALL).make();

    if ( altIndex < 1 || altIndex > nAlts ) throw new IllegalStateException("altIndex must be between 1 and nAlts " + nAlts);

    final double[] normalizedPr = MathUtils.normalizeFromLog10(GenotypeLikelihoods.fromPLs(original.getPL()).getAsVector());
    final double[] biAllelicPr = new double[3];

    for ( int index = 0; index < normalizedPr.length; index++ ) {
        final GenotypeLikelihoods.GenotypeLikelihoodsAllelePair pair = GenotypeLikelihoods.getAllelePair(index);

        if ( pair.alleleIndex1 == altIndex ) {
            if ( pair.alleleIndex2 == altIndex )
                // hom-alt case
                biAllelicPr[2] = normalizedPr[index];
            else
                // het-alt case
                biAllelicPr[1] += normalizedPr[index];
        } else {
            if ( pair.alleleIndex2 == altIndex )
                // het-alt case
                biAllelicPr[1] += normalizedPr[index];
            else
                // hom-non-alt
                biAllelicPr[0] += normalizedPr[index];
        }
    }

    final double[] GLs = new double[3];
    for ( int i = 0; i < GLs.length; i++ ) GLs[i] = Math.log10(biAllelicPr[i]);

    return new GenotypeBuilder(original).PL(GLs).alleles(BIALLELIC_NOCALL).make();
}
 

import htsjdk.variant.variantcontext.Genotype; //导入方法依赖的package包/类
/**
 * Returns a new Genotype with the PLs of the multi-allelic original reduced to a bi-allelic case.
 *
 * <p>Uses the log-sum-exp trick in order to work well with very low PLs</p>
 *
 * <p>This is handled in the following way:</p>
 *
 * <p>Suppose we have for a A/B/C site the following GLs:</p>
 *
 * <p>AA AB BB AC BC CC</p>
 *
 * <p>and we want to get the bi-allelic GLs for X/B, where X is everything not B</p>
 *
 * <p>XX = AA + AC + CC (since X = A or C)<br/>
 * XB = AB + BC                           <br/>
 * BB = BB                                     <br/>
 * </p>
 * <p>
 *     This implementation use the log sum trick in order to avoid numeric inestability.
 * </p>
 *
 * @param original the original multi-allelic genotype
 * @param altIndex the index of the alt allele we wish to keep in the bialleic case -- with ref == 0
 * @param nAlts the total number of alt alleles
 * @return a new biallelic genotype with appropriate PLs
 */
protected Genotype combineGLsPrecise(final Genotype original, final int altIndex, final int nAlts ) {

    if ( original.isNonInformative() )
        return new GenotypeBuilder(original).PL(BIALLELIC_NON_INFORMATIVE_PLS).alleles(BIALLELIC_NOCALL).make();

    if ( altIndex < 1 || altIndex > nAlts ) throw new IllegalStateException("altIndex must be between 1 and nAlts " + nAlts);

    final int[] pls = original.getPL();

    final int nAlleles = nAlts + 1;

    final int plCount = pls.length;

    double BB = 0;
    final double[] XBvalues = new double[nAlleles - 1];
    final double[] XXvalues = new double[plCount - nAlleles];

    int xbOffset = 0;
    int xxOffset = 0;
    for ( int index = 0; index < plCount; index++ ) {
        final GenotypeLikelihoods.GenotypeLikelihoodsAllelePair pair = GenotypeLikelihoods.getAllelePair(index);
        int i = pair.alleleIndex1;
        int j = pair.alleleIndex2;
        if (i == j) {
            if (i == altIndex) BB = PHRED_2_LOG10_COEFF * pls[index]; else XXvalues[xxOffset++] = PHRED_2_LOG10_COEFF * pls[index];
        } else if (i == altIndex || j == altIndex)
            XBvalues[xbOffset++] = PHRED_2_LOG10_COEFF * pls[index];
        else
            XXvalues[xxOffset++] = PHRED_2_LOG10_COEFF * pls[index];
    }

    final double XB = MathUtils.log10sumLog10(XBvalues);
    final double XX = MathUtils.log10sumLog10(XXvalues);

    final double[] GLs = new double[] { XX, XB, BB};
    return new GenotypeBuilder(original).PL(GLs).alleles(BIALLELIC_NOCALL).make();
}
 

import htsjdk.variant.variantcontext.Genotype; //导入方法依赖的package包/类
private VariantContext makeCombinedAltAllelesVariantContext(final VariantContext vc) {
    final int nAltAlleles = vc.getNAlleles() - 1;

    if ( nAltAlleles == 1 )
        return vc;
    else {
        final VariantContextBuilder vcb = new VariantContextBuilder(vc);
        final Allele reference = vcb.getAlleles().get(0);
        vcb.alleles(Arrays.asList(reference, GaeaVCFConstants.NON_REF_SYMBOLIC_ALLELE));
        final int genotypeCount = GenotypeLikelihoodCalculators.genotypeCount(2, vc.getNAlleles());
        final double[] hetLikelihoods = new double[vc.getNAlleles() - 1];
        final double[] homAltLikelihoods = new double[genotypeCount - hetLikelihoods.length - 1];
        final double[] newLikelihoods = new double[3];
        final List<Genotype> newGenotypes = new ArrayList<>(vc.getNSamples());
        for (final Genotype oldGenotype : vc.getGenotypes()) {
            final GenotypeBuilder gb = new GenotypeBuilder(oldGenotype);
            final List<Allele> oldAlleles = oldGenotype.getAlleles();
            if (oldAlleles != null) {
                final List<Allele> newAlleles = new ArrayList<>(oldAlleles.size());
                for (int i = 0; i < oldAlleles.size(); i++) {
                    final Allele oldAllele = oldAlleles.get(i);
                    if (oldAllele.isReference())
                        newAlleles.add(reference);
                    else if (oldAllele.isNoCall())
                        newAlleles.add(Allele.NO_CALL);
                    else
                        newAlleles.add(GaeaVCFConstants.NON_REF_SYMBOLIC_ALLELE);
                }
                gb.alleles(newAlleles);
            }
            if (oldGenotype.isNonInformative())
                gb.PL(BIALLELIC_NON_INFORMATIVE_PLS);
            else if (combineAltAlleleLikelihoods(oldGenotype, genotypeCount, newLikelihoods, hetLikelihoods, homAltLikelihoods))
                gb.PL(newLikelihoods);
            newGenotypes.add(gb.make());
        }
        return vcb.genotypes(newGenotypes).make();
    }
}
 

import htsjdk.variant.variantcontext.Genotype; //导入方法依赖的package包/类
/**
 * Returns a new Genotype with the PLs of the multi-allelic original reduced to a bi-allelic case.
 *
 * <p>Uses the log-sum-exp trick in order to work well with very low PLs</p>
 *
 * <p>This is handled in the following way:</p>
 *
 * <p>Suppose we have for a A/B/C site the following GLs:</p>
 *
 * <p>AA AB BB AC BC CC</p>
 *
 * <p>and we want to get the bi-allelic GLs for X/B, where X is everything not B</p>
 *
 * <p>XX = AA + AC + CC (since X = A or C)<br/>
 * XB = AB + BC                           <br/>
 * BB = BB                                     <br/>
 * </p>
 * <p>
 *     This implementation use the log sum trick in order to avoid numeric inestability.
 * </p>
 *
 * @param original the original multi-allelic genotype
 * @param altIndex the index of the alt allele we wish to keep in the bialleic case -- with ref == 0
 * @param nAlts the total number of alt alleles
 * @return a new biallelic genotype with appropriate PLs
 */
@Requires({"original.hasLikelihoods()"})
@Ensures({"result.hasLikelihoods()", "result.getPL().length == 3"})
protected Genotype combineGLsPrecise(final Genotype original, final int altIndex, final int nAlts ) {

    if ( original.isNonInformative() )
        return new GenotypeBuilder(original).PL(BIALLELIC_NON_INFORMATIVE_PLS).alleles(BIALLELIC_NOCALL).make();

    if ( altIndex < 1 || altIndex > nAlts ) throw new IllegalStateException("altIndex must be between 1 and nAlts " + nAlts);

    final int[] pls = original.getPL();

    final int nAlleles = nAlts + 1;

    final int plCount = pls.length;

    double BB = 0;
    final double[] XBvalues = new double[nAlleles - 1];
    final double[] XXvalues = new double[plCount - nAlleles];

    int xbOffset = 0;
    int xxOffset = 0;
    for ( int index = 0; index < plCount; index++ ) {
        final GenotypeLikelihoods.GenotypeLikelihoodsAllelePair pair = GenotypeLikelihoods.getAllelePair(index);
        int i = pair.alleleIndex1;
        int j = pair.alleleIndex2;
        if (i == j) {
          if (i == altIndex) BB = PHRED_2_LOG10_COEFF * pls[index]; else XXvalues[xxOffset++] = PHRED_2_LOG10_COEFF * pls[index];
        } else if (i == altIndex || j == altIndex)
          XBvalues[xbOffset++] = PHRED_2_LOG10_COEFF * pls[index];
        else
          XXvalues[xxOffset++] = PHRED_2_LOG10_COEFF * pls[index];
    }

    final double XB = GvcfMathUtils.log10sumLog10(XBvalues);
    final double XX = GvcfMathUtils.log10sumLog10(XXvalues);

    final double[] GLs = new double[] { XX, XB, BB};
    return new GenotypeBuilder(original).PL(GLs).alleles(BIALLELIC_NOCALL).make();
}
 

import htsjdk.variant.variantcontext.Genotype; //导入方法依赖的package包/类
/**
 * Returns a new Genotype with the PLs of the multi-allelic original reduced to a bi-allelic case.
 *
 * <p>Uses the log-sum-exp trick in order to work well with very low PLs</p>
 *
 * <p>This is handled in the following way:</p>
 *
 * <p>Suppose we have for a A/B/C site the following GLs:</p>
 *
 * <p>AA AB BB AC BC CC</p>
 *
 * <p>and we want to get the bi-allelic GLs for X/B, where X is everything not B</p>
 *
 * <p>XX = AA + AC + CC (since X = A or C)<br/>
 * XB = AB + BC                           <br/>
 * BB = BB                                <br/>
 * </p>
 * <p>
 *     This implementation uses the log-sum-exp trick in order to avoid numeric instability (underflow).
 * </p>
 *
 * @param original the original multi-allelic genotype
 * @param alleleIndex the index of the alt allele we wish to keep in the bialleic case -- with ref == 0
 * @param numberOfAlleles the total number of alleles (alternatives + the reference).
 * @return a new biallelic genotype with appropriate PLs
 */
@Requires({"original.hasLikelihoods() && alleleIndex >= 0"})
@Ensures({"result.hasLikelihoods()"})
private Genotype combineGLs(final Genotype original, final int defaultPloidy, final int alleleIndex, final int numberOfAlleles ) {

    final int declaredPloidy = original.getPloidy();
    final int ploidy = declaredPloidy <= 0 ? defaultPloidy : declaredPloidy;
    if ( original.isNonInformative() )
        return new GenotypeBuilder(original).PL(biallelicNonInformativePls(ploidy)).alleles(biallelicNoCall(ploidy)).make();

    final int[] pls = original.getPL();

    final GenotypeLikelihoodCalculator calculator = GenotypeLikelihoodCalculators.getInstance(ploidy, numberOfAlleles);
    final double[] newPLs = new double[ploidy + 1];
    Arrays.fill(newPLs, Double.NEGATIVE_INFINITY);
    for (int i = 0; i < pls.length; i++) {
        final GenotypeAlleleCounts alleleCounts = calculator.genotypeAlleleCountsAt(i);
        final int alleleCount = alleleCounts.alleleCountFor(alleleIndex);
        final int newPLIndex = alleleIndex == 0 ? ploidy - alleleCount : alleleCount;
        newPLs[newPLIndex] = GvcfMathUtils.approximateLog10SumLog10(newPLs[newPLIndex], -.1 * pls[i]);
    }

    return new GenotypeBuilder(original).PL(newPLs).alleles(biallelicNoCall(ploidy)).make();
}