Java GenotypeBuilder.noPL方法代码示例

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
/**
 * Assign genotypes (GTs) to the samples in the Variant Context greedily based on the PLs
 *
 * @param newLikelihoods the PL array
 * @param allelesToUse   the list of alleles to choose from (corresponding to the PLs)
 * @param numChromosomes Number of chromosomes per pool
 */
private void assignGenotype(final GenotypeBuilder gb,
                            final double[] newLikelihoods,
                            final List<Allele> allelesToUse,
                            final int numChromosomes) {
    final int numNewAltAlleles = allelesToUse.size() - 1;

    // find the genotype with maximum likelihoods
    final int PLindex = numNewAltAlleles == 0 ? 0 : MathUtils.maxElementIndex(newLikelihoods);
    final GenotypeLikelihoodCalculator calculator = GenotypeLikelihoodCalculators.getInstance(numChromosomes, allelesToUse.size());
    final GenotypeAlleleCounts alleleCounts = calculator.genotypeAlleleCountsAt(PLindex);

    gb.alleles(alleleCounts.asAlleleList(allelesToUse));

    // remove PLs if necessary
    if (newLikelihoods.length > MAX_LENGTH_FOR_POOL_PL_LOGGING)
        gb.noPL();

    // TODO - deprecated so what is the appropriate method to call?
    if (numNewAltAlleles > 0)
        gb.log10PError(GenotypeLikelihoods.getGQLog10FromLikelihoods(PLindex, newLikelihoods));
}
 

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
/**
 * Create the new GenotypesContext with the subsetted PLs and ADs
 *
 * @param originalGs               the original GenotypesContext
 * @param vc                       the original VariantContext
 * @param allelesToUse             the actual alleles to use with the new Genotypes
 * @param likelihoodIndexesToUse   the indexes in the PL to use given the allelesToUse (@see #determineLikelihoodIndexesToUse())
 * @param assignGenotypes          assignment strategy for the (subsetted) PLs
 * @return a new non-null GenotypesContext
 */
private static GenotypesContext createGenotypesWithSubsettedLikelihoods(final GenotypesContext originalGs,
                                                                        final VariantContext vc,
                                                                        final List<Allele> allelesToUse,
                                                                        final List<Integer> likelihoodIndexesToUse,
                                                                        final GenotypeAssignmentMethod assignGenotypes) {
    // the new genotypes to create
    final GenotypesContext newGTs = GenotypesContext.create(originalGs.size());

    // make sure we are seeing the expected number of likelihoods per sample
    final int expectedNumLikelihoods = GenotypeLikelihoods.numLikelihoods(vc.getNAlleles(), 2);

    // the samples
    final List<String> sampleIndices = originalGs.getSampleNamesOrderedByName();

    // create the new genotypes
    for ( int k = 0; k < originalGs.size(); k++ ) {
        final Genotype g = originalGs.get(sampleIndices.get(k));
        final GenotypeBuilder gb = new GenotypeBuilder(g);

        // create the new likelihoods array from the alleles we are allowed to use
        double[] newLikelihoods;
        if ( !g.hasLikelihoods() ) {
            // we don't have any likelihoods, so we null out PLs and make G ./.
            newLikelihoods = null;
            gb.noPL();
        } else {
            final double[] originalLikelihoods = g.getLikelihoods().getAsVector();
            if ( likelihoodIndexesToUse == null ) {
                newLikelihoods = originalLikelihoods;
            } else if ( originalLikelihoods.length != expectedNumLikelihoods ) {
                newLikelihoods = null;
            } else {
                newLikelihoods = new double[likelihoodIndexesToUse.size()];
                int newIndex = 0;
                for ( final int oldIndex : likelihoodIndexesToUse )
                    newLikelihoods[newIndex++] = originalLikelihoods[oldIndex];

                // might need to re-normalize
                newLikelihoods = MathUtils.normalizeFromLog10(newLikelihoods, false, true);
            }

            if ( newLikelihoods == null || likelihoodsAreUninformative(newLikelihoods) )
                gb.noPL();
            else
                gb.PL(newLikelihoods);
        }

        updateGenotypeAfterSubsetting(g.getAlleles(), gb, assignGenotypes, newLikelihoods, allelesToUse);
        newGTs.add(gb.make());
    }

    return fixADFromSubsettedAlleles(newGTs, vc, allelesToUse);
}
 

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
/**
 * Cleans up genotype-level annotations that need to be updated.
 * 1. move MIN_DP to DP if present
 * 2. propagate DP to AD if not present
 * 3. remove SB if present
 * 4. change the PGT value from "0|1" to "1|1" for homozygous variant genotypes
 *
 * @param VC            the VariantContext with the Genotypes to fix
 * @param createRefGTs  if true we will also create proper hom ref genotypes since we assume the site is monomorphic
 * @return a new set of Genotypes
 */
private List<Genotype> cleanupGenotypeAnnotations(final VariantContext VC, final boolean createRefGTs) {
    final GenotypesContext oldGTs = VC.getGenotypes();
    final List<Genotype> recoveredGs = new ArrayList<>(oldGTs.size());
    for ( final Genotype oldGT : oldGTs ) {
        final Map<String, Object> attrs = new HashMap<>(oldGT.getExtendedAttributes());

        final GenotypeBuilder builder = new GenotypeBuilder(oldGT);
        int depth = oldGT.hasDP() ? oldGT.getDP() : 0;

        // move the MIN_DP to DP
        if ( oldGT.hasExtendedAttribute("MIN_DP") ) {
            depth = Integer.parseInt((String)oldGT.getAnyAttribute("MIN_DP"));
            builder.DP(depth);
            attrs.remove("MIN_DP");
        }

        // remove SB
        attrs.remove("SB");

        // update PGT for hom vars
        if ( oldGT.isHomVar() && oldGT.hasExtendedAttribute(HaplotypeCaller.HAPLOTYPE_CALLER_PHASING_GT_KEY) ) {
            attrs.put(HaplotypeCaller.HAPLOTYPE_CALLER_PHASING_GT_KEY, "1|1");
        }

        // create AD if it's not there
        if ( !oldGT.hasAD() && VC.isVariant() ) {
            final int[] AD = new int[VC.getNAlleles()];
            AD[0] = depth;
            builder.AD(AD);
        }

        if ( createRefGTs ) {
            final int ploidy = oldGT.getPloidy();
            final List<Allele> refAlleles = Collections.nCopies(ploidy,VC.getReference());

            //keep 0 depth samples as no-call
            if (depth > 0) {
                builder.alleles(refAlleles);
            }

            // also, the PLs are technically no longer usable
            builder.noPL();
        }

        recoveredGs.add(builder.noAttributes().attributes(attrs).make());
    }
    return recoveredGs;
}
 

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
/**
 * From a given variant context, extract a given subset of alleles, and update genotype context accordingly,
 * including updating the PL's, and assign genotypes accordingly
 *
 * @param vc              variant context with alleles and genotype likelihoods
 * @param defaultPloidy   ploidy to assume in case that {@code vc} does not contain that information
 *                        for a sample.
 * @param allelesToUse    alleles to subset
 * @param assignGenotypes true: assign hard genotypes, false: leave as no-call
 * @return GenotypesContext with new PLs
 */
public GenotypesContext subsetAlleles(final VariantContext vc, final int defaultPloidy,
                                      final List<Allele> allelesToUse,
                                      final boolean assignGenotypes) {
    // the genotypes with PLs
    final GenotypesContext oldGTs = vc.getGenotypes();

    // samples
    final List<String> sampleIndices = oldGTs.getSampleNamesOrderedByName();

    // the new genotypes to create
    final GenotypesContext newGTs = GenotypesContext.create();

    // we need to determine which of the alternate alleles (and hence the likelihoods) to use and carry forward
    final int numOriginalAltAlleles = vc.getAlternateAlleles().size();
    final int numNewAltAlleles = allelesToUse.size() - 1;


    // create the new genotypes
    for (int k = 0; k < oldGTs.size(); k++) {
        final Genotype g = oldGTs.get(sampleIndices.get(k));
        final int declaredPloidy = g.getPloidy();
        final int ploidy = declaredPloidy <= 0 ? defaultPloidy : declaredPloidy;
        if (!g.hasLikelihoods()) {
            newGTs.add(GenotypeBuilder.create(g.getSampleName(), GATKVariantContextUtils.noCallAlleles(ploidy)));
            continue;
        }

        // create the new likelihoods array from the alleles we are allowed to use
        final double[] originalLikelihoods = g.getLikelihoods().getAsVector();
        double[] newLikelihoods;

        // Optimization: if # of new alt alleles = 0 (pure ref call), keep original likelihoods so we skip normalization
        // and subsetting
        if (numOriginalAltAlleles == numNewAltAlleles || numNewAltAlleles == 0) {
            newLikelihoods = originalLikelihoods;
        } else {
            newLikelihoods = GeneralPloidyGenotypeLikelihoods.subsetToAlleles(originalLikelihoods, ploidy, vc.getAlleles(), allelesToUse);

            // might need to re-normalize
            newLikelihoods = MathUtils.normalizeFromLog10(newLikelihoods, false, true);
        }

        // if there is no mass on the (new) likelihoods, then just no-call the sample
        if (MathUtils.sum(newLikelihoods) > GATKVariantContextUtils.SUM_GL_THRESH_NOCALL) {
            newGTs.add(GenotypeBuilder.create(g.getSampleName(), GATKVariantContextUtils.noCallAlleles(ploidy)));
        } else {
            final GenotypeBuilder gb = new GenotypeBuilder(g);

            if (numNewAltAlleles == 0)
                gb.noPL();
            else
                gb.PL(newLikelihoods);

            // if we weren't asked to assign a genotype, then just no-call the sample
            if (!assignGenotypes || MathUtils.sum(newLikelihoods) > GATKVariantContextUtils.SUM_GL_THRESH_NOCALL)
                gb.alleles(GATKVariantContextUtils.noCallAlleles(ploidy));
            else
                assignGenotype(gb, newLikelihoods, allelesToUse, ploidy);
            newGTs.add(gb.make());
        }
    }

    return newGTs;

}
 

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
private List<Genotype> cleanupGenotypeAnnotations(final VariantContext VC, final boolean createRefGTs) {
	final GenotypesContext oldGTs = VC.getGenotypes();
	final List<Genotype> recoveredGs = new ArrayList<>(oldGTs.size());
	for (final Genotype oldGT : oldGTs) {
		final Map<String, Object> attrs = new HashMap<>(oldGT.getExtendedAttributes());

		final GenotypeBuilder builder = new GenotypeBuilder(oldGT);
		int depth = oldGT.hasDP() ? oldGT.getDP() : 0;

		// move the MIN_DP to DP
		if (oldGT.hasExtendedAttribute(GaeaVCFConstants.MIN_DP_FORMAT_KEY)) {
			depth = Integer.parseInt((String) oldGT.getAnyAttribute(GaeaVCFConstants.MIN_DP_FORMAT_KEY));
			builder.DP(depth);
			attrs.remove(GaeaVCFConstants.MIN_DP_FORMAT_KEY);
		}

		// move the GQ to RGQ
		if (createRefGTs && oldGT.hasGQ()) {
			builder.noGQ();
			attrs.put(GaeaVCFConstants.REFERENCE_GENOTYPE_QUALITY, oldGT.getGQ());
		}

		// remove SB
		attrs.remove(GaeaVCFConstants.STRAND_BIAS_BY_SAMPLE_KEY);

		// update PGT for hom vars
		if (oldGT.isHomVar() && oldGT.hasExtendedAttribute(GaeaVCFConstants.HAPLOTYPE_CALLER_PHASING_GT_KEY)) {
			attrs.put(GaeaVCFConstants.HAPLOTYPE_CALLER_PHASING_GT_KEY, "1|1");
		}

		// create AD if it's not there
		if (!oldGT.hasAD() && VC.isVariant()) {
			final int[] AD = new int[VC.getNAlleles()];
			AD[0] = depth;
			builder.AD(AD);
		}

		if (createRefGTs) {
			final int ploidy = oldGT.getPloidy();
			final List<Allele> refAlleles = Collections.nCopies(ploidy, VC.getReference());

			// keep 0 depth samples and 0 GQ samples as no-call
			if (depth > 0 && oldGT.hasGQ() && oldGT.getGQ() > 0) {
				builder.alleles(refAlleles);
			}

			// also, the PLs are technically no longer usable
			builder.noPL();
		}

		recoveredGs.add(builder.noAttributes().attributes(attrs).make());
	}
	return recoveredGs;
}
 

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
@Override
@Requires("vc != null && allelesToUse != null")
public GenotypesContext subsetAlleles(VariantContext vc, int defaultPloidy, List<Allele> allelesToUse, boolean assignGenotypes) {
    // the genotypes with PLs
    final GenotypesContext oldGTs = vc.getGenotypes();

    // samples
    final List<String> sampleIndices = oldGTs.getSampleNamesOrderedByName();

    // the new genotypes to create
    final GenotypesContext newGTs = GenotypesContext.create();

    // we need to determine which of the alternate alleles (and hence the likelihoods) to use and carry forward
    final int numOriginalAltAlleles = vc.getAlternateAlleles().size();
    final int numNewAltAlleles = allelesToUse.size() - 1;


    // create the new genotypes
    for ( int k = 0; k < oldGTs.size(); k++ ) {
        final Genotype g = oldGTs.get(sampleIndices.get(k));
        final int declaredPloidy = g.getPloidy();
        final int ploidy = declaredPloidy <= 0 ? defaultPloidy : declaredPloidy;
        if ( !g.hasLikelihoods() ) {
            newGTs.add(GenotypeBuilder.create(g.getSampleName(),GaeaGvcfVariantContextUtils.noCallAlleles(ploidy)));
            continue;
        }

        // create the new likelihoods array from the alleles we are allowed to use
        final double[] originalLikelihoods = g.getLikelihoods().getAsVector();
        double[] newLikelihoods;

        // Optimization: if # of new alt alleles = 0 (pure ref call), keep original likelihoods so we skip normalization
        // and subsetting
        if ( numOriginalAltAlleles == numNewAltAlleles || numNewAltAlleles == 0) {
            newLikelihoods = originalLikelihoods;
        } else {
            newLikelihoods = GeneralPloidyGenotypeLikelihoods.subsetToAlleles(originalLikelihoods, ploidy, vc.getAlleles(), allelesToUse);

            // might need to re-normalize
            newLikelihoods = GvcfMathUtils.normalizeFromLog10(newLikelihoods, false, true);
        }

        // if there is no mass on the (new) likelihoods, then just no-call the sample
        if ( GvcfMathUtils.sum(newLikelihoods) > GaeaGvcfVariantContextUtils.SUM_GL_THRESH_NOCALL ) {
            newGTs.add(GenotypeBuilder.create(g.getSampleName(), GaeaGvcfVariantContextUtils.noCallAlleles(ploidy)));
        } else {
            final GenotypeBuilder gb = new GenotypeBuilder(g);
            final String sampleName = g.getSampleName();

            if ( numNewAltAlleles == 0 )
                gb.noPL();
            else
                gb.PL(newLikelihoods);

            // if we weren't asked to assign a genotype, then just no-call the sample
            if ( !assignGenotypes || GvcfMathUtils.sum(newLikelihoods) > GaeaGvcfVariantContextUtils.SUM_GL_THRESH_NOCALL )
                gb.alleles(GaeaGvcfVariantContextUtils.noCallAlleles(ploidy));
            else
                assignGenotype(gb, vc, sampleName, newLikelihoods, allelesToUse, ploidy);
            newGTs.add(gb.make());
        }
    }

    return GaeaGvcfVariantContextUtils.fixADFromSubsettedAlleles(newGTs, vc, allelesToUse);
}
 

import htsjdk.variant.variantcontext.GenotypeBuilder; //导入方法依赖的package包/类
/**
 * Strips PLs from the specified GenotypeBuilder if their number exceeds the
 * maximum allowed. Corresponding counters are updated.
 * 
 * @param gb
 *            the GenotypeBuilder to modify
 * @param vc
 *            the VariantContext
 * @param sampleName
 *            the sample name
 * @param newLikelihoods
 *            the PL array
 */
protected void removePLsIfMaxNumPLValuesExceeded(final GenotypeBuilder gb, final VariantContext vc,
		final String sampleName, final double[] newLikelihoods) {
	final int numPLValuesFound = newLikelihoods.length;
	if (numPLValuesFound > maxNumPLValues) {
		numTimesMaxNumPLValuesExceeded++;
		gb.noPL();
		if (numPLValuesFound > maxNumPLValuesObserved) {
			maxNumPLValuesObserved = numPLValuesFound;
		}
	}
}