Examples with UnivariateObjectiveFunction org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction used on opensource projects

Example 1 with UnivariateObjectiveFunction

use of org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction in project GDSC-SMLM by aherbert.

the class FIRE method findMin.

private UnivariatePointValuePair findMin(UnivariatePointValuePair current, UnivariateOptimizer o, UnivariateFunction f, double qValue, double factor) {
    try {
        BracketFinder bracket = new BracketFinder();
        bracket.search(f, GoalType.MINIMIZE, qValue * factor, qValue / factor);
        UnivariatePointValuePair next = o.optimize(GoalType.MINIMIZE, new MaxEval(3000), new SearchInterval(bracket.getLo(), bracket.getHi(), bracket.getMid()), new UnivariateObjectiveFunction(f));
        if (next == null)
            return current;
        //System.out.printf("LineMin [%.1f]  %f = %f\n", factor, next.getPoint(), next.getValue());
        if (current != null)
            return (next.getValue() < current.getValue()) ? next : current;
        return next;
    } catch (Exception e) {
        return current;
    }
}

Example 2 with UnivariateObjectiveFunction

use of org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction in project gatk-protected by broadinstitute.

the class PosteriorSummaryUtils method calculatePosteriorMode.

/**
     * Given a list of posterior samples, returns an estimate of the posterior mode (using
     * mllib kernel density estimation in {@link KernelDensity} and {@link BrentOptimizer}).
     * Note that estimate may be poor if number of samples is small (resulting in poor kernel density estimation),
     * or if posterior is not unimodal (or is sufficiently pathological otherwise). If the samples contain
     * {@link Double#NaN}, {@link Double#NaN} will be returned.
     * @param samples   posterior samples, cannot be {@code null} and number of samples must be greater than 0
     * @param ctx       {@link JavaSparkContext} used by {@link KernelDensity} for mllib kernel density estimation
     */
public static double calculatePosteriorMode(final List<Double> samples, final JavaSparkContext ctx) {
    Utils.nonNull(samples);
    Utils.validateArg(samples.size() > 0, "Number of samples must be greater than zero.");
    //calculate sample min, max, mean, and standard deviation
    final double sampleMin = Collections.min(samples);
    final double sampleMax = Collections.max(samples);
    final double sampleMean = new Mean().evaluate(Doubles.toArray(samples));
    final double sampleStandardDeviation = new StandardDeviation().evaluate(Doubles.toArray(samples));
    //if samples are all the same or contain NaN, can simply return mean
    if (sampleStandardDeviation == 0. || Double.isNaN(sampleMean)) {
        return sampleMean;
    }
    //use Silverman's rule to set bandwidth for kernel density estimation from sample standard deviation
    //see https://en.wikipedia.org/wiki/Kernel_density_estimation#Practical_estimation_of_the_bandwidth
    final double bandwidth = SILVERMANS_RULE_CONSTANT * sampleStandardDeviation * Math.pow(samples.size(), SILVERMANS_RULE_EXPONENT);
    //use kernel density estimation to approximate posterior from samples
    final KernelDensity pdf = new KernelDensity().setSample(ctx.parallelize(samples, 1)).setBandwidth(bandwidth);
    //use Brent optimization to find mode (i.e., maximum) of kernel-density-estimated posterior
    final BrentOptimizer optimizer = new BrentOptimizer(RELATIVE_TOLERANCE, RELATIVE_TOLERANCE * (sampleMax - sampleMin));
    final UnivariateObjectiveFunction objective = new UnivariateObjectiveFunction(f -> pdf.estimate(new double[] { f })[0]);
    //search for mode within sample range, start near sample mean
    final SearchInterval searchInterval = new SearchInterval(sampleMin, sampleMax, sampleMean);
    return optimizer.optimize(objective, GoalType.MAXIMIZE, searchInterval, BRENT_MAX_EVAL).getPoint();
}

Example 3 with UnivariateObjectiveFunction

use of org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction in project gatk-protected by broadinstitute.

the class CNLOHCaller method optimizeIt.

private double optimizeIt(final Function<Double, Double> objectiveFxn, final SearchInterval searchInterval) {
    final MaxEval BRENT_MAX_EVAL = new MaxEval(1000);
    final double RELATIVE_TOLERANCE = 0.001;
    final double ABSOLUTE_TOLERANCE = 0.001;
    final BrentOptimizer OPTIMIZER = new BrentOptimizer(RELATIVE_TOLERANCE, ABSOLUTE_TOLERANCE);
    final UnivariateObjectiveFunction objective = new UnivariateObjectiveFunction(x -> objectiveFxn.apply(x));
    return OPTIMIZER.optimize(objective, GoalType.MAXIMIZE, searchInterval, BRENT_MAX_EVAL).getPoint();
}

Example 4 with UnivariateObjectiveFunction

use of org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction in project imagingbook-common by imagingbook.

the class FourierDescriptor method getStartPointPhase.

/**
 * Calculates the 'canonical' start point. This version uses
 * (a) a coarse search for a global maximum of fp() and subsequently
 * (b) a numerical optimization using Brent's method
 * (implemented with Apache Commons Math).
 *
 * @param Mp number of Fourier coefficient pairs
 * @return start point phase
 */
public double getStartPointPhase(int Mp) {
    Mp = Math.min(Mp, (G.length - 1) / 2);
    UnivariateFunction fp = new TargetFunction(Mp);
    // search for the global maximum in coarse steps
    double cmax = Double.NEGATIVE_INFINITY;
    int kmax = -1;
    // number of steps over 180 degrees
    int K = 25;
    for (int k = 0; k < K; k++) {
        // phase to evaluate
        final double phi = Math.PI * k / K;
        final double c = fp.value(phi);
        if (c > cmax) {
            cmax = c;
            kmax = k;
        }
    }
    // optimize using previous and next point as the bracket.
    double minPhi = Math.PI * (kmax - 1) / K;
    double maxPhi = Math.PI * (kmax + 1) / K;
    UnivariateOptimizer optimizer = new BrentOptimizer(1E-4, 1E-6);
    int maxIter = 20;
    UnivariatePointValuePair result = optimizer.optimize(new MaxEval(maxIter), new UnivariateObjectiveFunction(fp), GoalType.MAXIMIZE, new SearchInterval(minPhi, maxPhi));
    double phi0 = result.getPoint();
    // the canonical start point phase
    return phi0;
}

Example 5 with UnivariateObjectiveFunction

use of org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction in project GDSC-SMLM by aherbert.

the class PoissonGammaGaussianFisherInformation method findMaximum.

/**
 * Find the maximum of the integrated function A^2/P. P is the Poisson-Gamma convolution, A is the
 * partial gradient.
 *
 * <p>When both A and P are convolved with a Gaussian kernel, the integral of this function - 1 is
 * the Fisher information.
 *
 * <p>This method is used to determine the maximum of the function. The integral of A should equal
 * 1. The range can be determined using a fraction of the maximum, or integrating until the sum is
 * 1.
 *
 * @param theta the Poisson mean
 * @param rel Relative threshold
 * @return [max,max value,evaluations]
 */
public double[] findMaximum(final double theta, double rel) {
    if (theta < MIN_MEAN) {
        throw new IllegalArgumentException("Poisson mean must be positive");
    }
    final double dirac = PoissonGammaFunction.dirac(theta);
    final UnivariateFunction f = new UnivariateFunction() {

        double[] gradient = new double[1];

        @Override
        public double value(double x) {
            // Return F without the dirac. This is so that the maximum of the function
            // is known for relative height analysis.
            double value;
            if (x == 0) {
                gradient[0] = dirac / gain;
                value = gradient[0] * theta;
            } else {
                value = PoissonGammaFunction.unscaledPoissonGammaPartial(x, theta, gain, gradient);
            }
            return getF(value, gradient[0]);
        }
    };
    // Determine the extent of the function: A^2/P
    // Without convolution this will have a maximum that is above, and
    // converges to [Poisson mean * amplification].
    final double mean = theta * gain;
    final int iterations = 10;
    final QuickConvergenceChecker checker = new QuickConvergenceChecker(iterations);
    final BrentOptimizer opt = new BrentOptimizer(rel, Double.MIN_VALUE, checker);
    UnivariatePointValuePair pair;
    try {
        final double start = (theta < 1) ? 0 : mean;
        pair = opt.optimize(new SearchInterval(0, mean * 1.5, start), GoalType.MAXIMIZE, new MaxEval(50000), new UnivariateObjectiveFunction(f));
    } catch (final TooManyEvaluationsException ex) {
        // This should not occur with the quick checker fixed iterations
        // - just get the current best
        pair = checker.getBest();
    }
    final double[] max = new double[3];
    max[0] = pair.getPoint();
    max[1] = pair.getValue();
    max[2] = opt.getEvaluations();
    return max;
}