Search in sources :

Example 11 with InitialGuess

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

the class PCPALMMolecules method fitGaussian.

private double[] fitGaussian(float[] x, float[] y) {
    WeightedObservedPoints obs = new WeightedObservedPoints();
    for (int i = 0; i < x.length; i++) obs.add(x[i], y[i]);
    Collection<WeightedObservedPoint> observations = obs.toList();
    GaussianCurveFitter fitter = GaussianCurveFitter.create().withMaxIterations(2000);
    GaussianCurveFitter.ParameterGuesser guess = new GaussianCurveFitter.ParameterGuesser(observations);
    double[] initialGuess = null;
    try {
        initialGuess = guess.guess();
        return fitter.withStartPoint(initialGuess).fit(observations);
    } catch (TooManyEvaluationsException e) {
        // Use the initial estimate
        return initialGuess;
    } catch (Exception e) {
        // Just in case there is another exception type, or the initial estimate failed
        return null;
    }
}
Also used : WeightedObservedPoints(org.apache.commons.math3.fitting.WeightedObservedPoints) TooManyEvaluationsException(org.apache.commons.math3.exception.TooManyEvaluationsException) WeightedObservedPoint(org.apache.commons.math3.fitting.WeightedObservedPoint) WeightedObservedPoint(org.apache.commons.math3.fitting.WeightedObservedPoint) ClusterPoint(gdsc.core.clustering.ClusterPoint) GaussianCurveFitter(org.apache.commons.math3.fitting.GaussianCurveFitter) TooManyEvaluationsException(org.apache.commons.math3.exception.TooManyEvaluationsException)

Example 12 with InitialGuess

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

the class PCPALMMolecules method optimiseSimplex.

private double[] optimiseSimplex(float[] x, float[] y, double[] initialSolution) {
    // Simplex optimisation
    SkewNormalMultivariateFunction sn2 = new SkewNormalMultivariateFunction(initialSolution);
    sn2.addData(x, y);
    NelderMeadSimplex simplex = new NelderMeadSimplex(4);
    SimplexOptimizer opt = new SimplexOptimizer(1e-6, 1e-10);
    PointValuePair solution = opt.optimize(new MaxEval(1000), new InitialGuess(initialSolution), simplex, new ObjectiveFunction(sn2), GoalType.MINIMIZE);
    double[] skewParameters2 = solution.getPointRef();
    return skewParameters2;
}
Also used : MaxEval(org.apache.commons.math3.optim.MaxEval) InitialGuess(org.apache.commons.math3.optim.InitialGuess) NelderMeadSimplex(org.apache.commons.math3.optim.nonlinear.scalar.noderiv.NelderMeadSimplex) SimplexOptimizer(org.apache.commons.math3.optim.nonlinear.scalar.noderiv.SimplexOptimizer) ObjectiveFunction(org.apache.commons.math3.optim.nonlinear.scalar.ObjectiveFunction) PointValuePair(org.apache.commons.math3.optim.PointValuePair)

Example 13 with InitialGuess

use of org.apache.commons.math3.optim.InitialGuess in project vcell by virtualcell.

the class FitBleachSpotOp method fitToGaussian.

static GaussianFitResults fitToGaussian(double init_center_i, double init_center_j, double init_radius2, FloatImage image) {
    // 
    // do some optimization on the image (fitting to a Gaussian)
    // set initial guesses from ROI operation.
    // 
    ISize imageSize = image.getISize();
    final int num_i = imageSize.getX();
    final int num_j = imageSize.getY();
    final float[] floatPixels = image.getFloatPixels();
    // 
    // initial guess based on previous fit of ROI
    // do gaussian fit in index space for center and standard deviation (later to translate it back to world coordinates)
    // 
    final int window_size = (int) Math.sqrt(init_radius2) * 4;
    // final int window_min_i = 0;       // (int) Math.max(0, Math.floor(init_center_i - window_size/2));
    // final int window_max_i = num_i-1; // (int) Math.min(num_i-1, Math.ceil(init_center_i + window_size/2));
    // final int window_min_j = 0;       // (int) Math.max(0, Math.floor(init_center_j - window_size/2));
    // final int window_max_j = num_j-1; // (int) Math.min(num_j-1, Math.ceil(init_center_j + window_size/2));
    final int window_min_i = (int) Math.max(0, Math.floor(init_center_i - window_size / 2));
    final int window_max_i = (int) Math.min(num_i - 1, Math.ceil(init_center_i + window_size / 2));
    final int window_min_j = (int) Math.max(0, Math.floor(init_center_j - window_size / 2));
    final int window_max_j = (int) Math.min(num_j - 1, Math.ceil(init_center_j + window_size / 2));
    final int PARAM_INDEX_CENTER_I = 0;
    final int PARAM_INDEX_CENTER_J = 1;
    final int PARAM_INDEX_K = 2;
    final int PARAM_INDEX_HIGH = 3;
    final int PARAM_INDEX_RADIUS_SQUARED = 4;
    final int NUM_PARAMETERS = 5;
    double[] initParameters = new double[NUM_PARAMETERS];
    initParameters[PARAM_INDEX_CENTER_I] = init_center_i;
    initParameters[PARAM_INDEX_CENTER_J] = init_center_j;
    initParameters[PARAM_INDEX_HIGH] = 1.0;
    initParameters[PARAM_INDEX_K] = 10;
    initParameters[PARAM_INDEX_RADIUS_SQUARED] = init_radius2;
    PowellOptimizer optimizer = new PowellOptimizer(1e-4, 1e-1);
    MultivariateFunction func = new MultivariateFunction() {

        @Override
        public double value(double[] point) {
            double center_i = point[PARAM_INDEX_CENTER_I];
            double center_j = point[PARAM_INDEX_CENTER_J];
            double high = point[PARAM_INDEX_HIGH];
            double K = point[PARAM_INDEX_K];
            double radius2 = point[PARAM_INDEX_RADIUS_SQUARED];
            double error2 = 0;
            for (int j = window_min_j; j <= window_max_j; j++) {
                // double y = j - center_j;
                double y = j;
                for (int i = window_min_i; i <= window_max_i; i++) {
                    // double x = i - center_i;
                    double x = i;
                    double modelValue = high - FastMath.exp(-K * FastMath.exp(-2 * (x * x + y * y) / radius2));
                    double imageValue = floatPixels[j * num_i + i];
                    double error = modelValue - imageValue;
                    error2 += error * error;
                }
            }
            System.out.println(new GaussianFitResults(center_i, center_j, radius2, K, high, error2));
            return error2;
        }
    };
    PointValuePair pvp = optimizer.optimize(new ObjectiveFunction(func), new InitialGuess(initParameters), new MaxEval(100000), GoalType.MINIMIZE);
    double[] fittedParamValues = pvp.getPoint();
    double fitted_center_i = fittedParamValues[PARAM_INDEX_CENTER_I];
    double fitted_center_j = fittedParamValues[PARAM_INDEX_CENTER_J];
    double fitted_radius2 = fittedParamValues[PARAM_INDEX_RADIUS_SQUARED];
    double fitted_K = fittedParamValues[PARAM_INDEX_K];
    double fitted_high = fittedParamValues[PARAM_INDEX_HIGH];
    double objectiveFunctionValue = pvp.getValue();
    return new GaussianFitResults(fitted_center_i, fitted_center_j, fitted_radius2, fitted_K, fitted_high, objectiveFunctionValue);
}
Also used : MultivariateFunction(org.apache.commons.math3.analysis.MultivariateFunction) InitialGuess(org.apache.commons.math3.optim.InitialGuess) MaxEval(org.apache.commons.math3.optim.MaxEval) ISize(org.vcell.util.ISize) ObjectiveFunction(org.apache.commons.math3.optim.nonlinear.scalar.ObjectiveFunction) PowellOptimizer(org.apache.commons.math3.optim.nonlinear.scalar.noderiv.PowellOptimizer) PointValuePair(org.apache.commons.math3.optim.PointValuePair)

Aggregations

InitialGuess (org.apache.commons.math3.optim.InitialGuess)11 MaxEval (org.apache.commons.math3.optim.MaxEval)11 PointValuePair (org.apache.commons.math3.optim.PointValuePair)11 ObjectiveFunction (org.apache.commons.math3.optim.nonlinear.scalar.ObjectiveFunction)11 TooManyEvaluationsException (org.apache.commons.math3.exception.TooManyEvaluationsException)10 TooManyIterationsException (org.apache.commons.math3.exception.TooManyIterationsException)5 OptimizationData (org.apache.commons.math3.optim.OptimizationData)5 SimpleBounds (org.apache.commons.math3.optim.SimpleBounds)5 CMAESOptimizer (org.apache.commons.math3.optim.nonlinear.scalar.noderiv.CMAESOptimizer)5 ConvergenceException (org.apache.commons.math3.exception.ConvergenceException)4 CustomPowellOptimizer (org.apache.commons.math3.optim.nonlinear.scalar.noderiv.CustomPowellOptimizer)4 MultivariateFunction (org.apache.commons.math3.analysis.MultivariateFunction)3 SimpleValueChecker (org.apache.commons.math3.optim.SimpleValueChecker)3 NelderMeadSimplex (org.apache.commons.math3.optim.nonlinear.scalar.noderiv.NelderMeadSimplex)3 RandomGenerator (org.apache.commons.math3.random.RandomGenerator)3 MultivariateMatrixFunction (org.apache.commons.math3.analysis.MultivariateMatrixFunction)2 GaussianCurveFitter (org.apache.commons.math3.fitting.GaussianCurveFitter)2 WeightedObservedPoint (org.apache.commons.math3.fitting.WeightedObservedPoint)2 WeightedObservedPoints (org.apache.commons.math3.fitting.WeightedObservedPoints)2 LeastSquaresBuilder (org.apache.commons.math3.fitting.leastsquares.LeastSquaresBuilder)2