Examples with RandomGenerator org.apache.commons.math3.random.RandomGenerator used on opensource projects

Example 21 with RandomGenerator

use of org.apache.commons.math3.random.RandomGenerator in project gatk by broadinstitute.

the class HMM method generateHiddenStateChain.

default default List<S> generateHiddenStateChain(final List<T> positions) {
    final RandomGenerator rg = RandomGeneratorFactory.createRandomGenerator(new Random(RANDOM_SEED_FOR_CHAIN_GENERATION));
    final List<S> hiddenStates = hiddenStates();
    final List<S> result = new ArrayList<>(positions.size());
    final S initialState = GATKProtectedMathUtils.randomSelect(hiddenStates, s -> Math.exp(logPriorProbability(s, positions.get(0))), rg);
    result.add(initialState);
    IntStream.range(1, positions.size()).forEach(n -> result.add(GATKProtectedMathUtils.randomSelect(hiddenStates, s -> Math.exp(logTransitionProbability(result.get(n - 1), positions.get(n - 1), s, positions.get(n))), rg)));
    return result;
}

Example 22 with RandomGenerator

use of org.apache.commons.math3.random.RandomGenerator in project gatk by broadinstitute.

the class AdaptiveMetropolisSampler method sample.

public double sample(final RandomGenerator rng, final Function<Double, Double> logPDF) {
    Utils.nonNull(rng);
    Utils.nonNull(logPDF);
    final AbstractRealDistribution normal = new NormalDistribution(rng, 0, 1);
    final double proposal = xCurrent + stepSize * normal.sample();
    final double acceptanceProbability = (proposal < lowerBound || upperBound < proposal) ? 0 : Math.min(1, Math.exp(logPDF.apply(proposal) - logPDF.apply(xCurrent)));
    //adjust stepSize larger/smaller to decrease/increase the acceptance rate
    final double correctionFactor = (acceptanceProbability - optimalAcceptanceRate) * adjustmentRate * (timeScale / (timeScale + iteration));
    stepSize *= Math.exp(correctionFactor);
    iteration++;
    return rng.nextDouble() < acceptanceProbability ? proposal : xCurrent;
}

Example 23 with RandomGenerator

use of org.apache.commons.math3.random.RandomGenerator in project gatk-protected by broadinstitute.

the class CoverageModelParameters method generateRandomModel.

/**
     * Generates random coverage model parameters.
     *
     * @param targetList list of targets
     * @param numLatents number of latent variables
     * @param seed random seed
     * @param randomMeanLogBiasStandardDeviation std of mean log bias (mean is set to 0)
     * @param randomBiasCovariatesStandardDeviation std of bias covariates (mean is set to 0)
     * @param randomMaxUnexplainedVariance max value of unexplained variance (samples are taken from a uniform
     *                                     distribution [0, {@code randomMaxUnexplainedVariance}])
     * @param initialBiasCovariatesARDCoefficients initial row vector of ARD coefficients
     * @return an instance of {@link CoverageModelParameters}
     */
public static CoverageModelParameters generateRandomModel(final List<Target> targetList, final int numLatents, final long seed, final double randomMeanLogBiasStandardDeviation, final double randomBiasCovariatesStandardDeviation, final double randomMaxUnexplainedVariance, final INDArray initialBiasCovariatesARDCoefficients) {
    Utils.validateArg(numLatents >= 0, "Dimension of the bias space must be non-negative");
    Utils.validateArg(randomBiasCovariatesStandardDeviation >= 0, "Standard deviation of random bias covariates" + " must be non-negative");
    Utils.validateArg(randomMeanLogBiasStandardDeviation >= 0, "Standard deviation of random mean log bias" + " must be non-negative");
    Utils.validateArg(randomMaxUnexplainedVariance >= 0, "Max random unexplained variance must be non-negative");
    Utils.validateArg(initialBiasCovariatesARDCoefficients == null || numLatents > 0 && initialBiasCovariatesARDCoefficients.length() == numLatents, "If ARD is enabled, the dimension" + " of the bias latent space must be positive and match the length of ARD coeffecient vector");
    final boolean biasCovariatesEnabled = numLatents > 0;
    final int numTargets = targetList.size();
    final RandomGenerator rng = RandomGeneratorFactory.createRandomGenerator(new Random(seed));
    /* Gaussian random for mean log bias */
    final INDArray initialMeanLogBias = Nd4j.create(getNormalRandomNumbers(numTargets, 0, randomMeanLogBiasStandardDeviation, rng), new int[] { 1, numTargets });
    /* Uniform random for unexplained variance */
    final INDArray initialUnexplainedVariance = Nd4j.create(getUniformRandomNumbers(numTargets, 0, randomMaxUnexplainedVariance, rng), new int[] { 1, numTargets });
    final INDArray initialMeanBiasCovariates;
    if (biasCovariatesEnabled) {
        /* Gaussian random for bias covariates */
        initialMeanBiasCovariates = Nd4j.create(getNormalRandomNumbers(numTargets * numLatents, 0, randomBiasCovariatesStandardDeviation, rng), new int[] { numTargets, numLatents });
    } else {
        initialMeanBiasCovariates = null;
    }
    return new CoverageModelParameters(targetList, initialMeanLogBias, initialUnexplainedVariance, initialMeanBiasCovariates, initialBiasCovariatesARDCoefficients);
}

Example 24 with RandomGenerator

use of org.apache.commons.math3.random.RandomGenerator in project gatk-protected by broadinstitute.

the class AlleleFractionSegmenterUnitTest method generateCounts.

//visible for testing joint segmentation
protected static AllelicCountCollection generateCounts(final List<Double> minorAlleleFractionSequence, final List<SimpleInterval> positions, final RandomGenerator rng, final AlleleFractionGlobalParameters trueParams) {
    //translate to ApacheCommons' parametrization of the gamma distribution
    final GammaDistribution biasGenerator = getGammaDistribution(trueParams, rng);
    final double outlierProbability = trueParams.getOutlierProbability();
    final AllelicCountCollection counts = new AllelicCountCollection();
    for (int n = 0; n < minorAlleleFractionSequence.size(); n++) {
        counts.add(generateAllelicCount(minorAlleleFractionSequence.get(n), positions.get(n), rng, biasGenerator, outlierProbability));
    }
    return counts;
}

Example 25 with RandomGenerator

use of org.apache.commons.math3.random.RandomGenerator in project gatk-protected by broadinstitute.

the class AlleleFractionSegmenterUnitTest method testSegmentation.

@Test
public void testSegmentation() {
    final RandomGenerator rng = RandomGeneratorFactory.createRandomGenerator(new Random(563));
    final List<Double> trueWeights = Arrays.asList(0.2, 0.5, 0.3);
    final List<Double> trueMinorAlleleFractions = Arrays.asList(0.12, 0.32, 0.5);
    final double trueMemoryLength = 1e5;
    final AlleleFractionGlobalParameters trueParams = new AlleleFractionGlobalParameters(1.0, 0.01, 0.01);
    final AlleleFractionHMM trueModel = new AlleleFractionHMM(trueMinorAlleleFractions, trueWeights, trueMemoryLength, AllelicPanelOfNormals.EMPTY_PON, trueParams);
    // randomly set positions
    final int chainLength = 10000;
    final List<SimpleInterval> positions = CopyRatioSegmenterUnitTest.randomPositions("chr1", chainLength, rng, trueMemoryLength / 4);
    final List<Integer> trueStates = trueModel.generateHiddenStateChain(positions);
    final List<Double> truthMinorFractions = trueStates.stream().map(trueModel::getMinorAlleleFraction).collect(Collectors.toList());
    final AllelicCountCollection counts = generateCounts(truthMinorFractions, positions, rng, trueParams);
    final AlleleFractionSegmenter segmenter = new AlleleFractionSegmenter(10, counts, AllelicPanelOfNormals.EMPTY_PON);
    final List<ModeledSegment> segments = segmenter.getModeledSegments();
    final double[] segmentMinorFractions = segments.stream().flatMap(s -> Collections.nCopies((int) s.getTargetCount(), s.getSegmentMean()).stream()).mapToDouble(x -> x).toArray();
    final double averageMinorFractionError = IntStream.range(0, truthMinorFractions.size()).mapToDouble(n -> Math.abs(segmentMinorFractions[n] - truthMinorFractions.get(n))).average().getAsDouble();
    Assert.assertEquals(averageMinorFractionError, 0, 0.01);
}