Example usage for org.apache.commons.math.optimization.direct NelderMead NelderMead

List of usage examples for org.apache.commons.math.optimization.direct NelderMead NelderMead

Introduction

In this page you can find the example usage for org.apache.commons.math.optimization.direct NelderMead NelderMead.

Prototype

public NelderMead() 

Source Link

Document

Build a Nelder-Mead optimizer with default coefficients.

Usage

From source file:elkfed.main.AutoTune.java

public static void main(String[] args) {
    try {/*www  . j a  va2 s  .  c o  m*/
        CorefExperimentDocument doc;
        doc = CorefExperimentDocument.Factory.parse(new FileInputStream(args[0]));
        Evaluator eval = new Evaluator(doc);
        double[][] samples = new double[eval.getNumParameters() + 1][];
        for (int i = 0; i < eval.getNumParameters() + 1; i++) {
            samples[i] = eval.sample();
        }
        DirectSearchOptimizer optimizer = new NelderMead();
        try {
            optimizer.setStartConfiguration(samples);
            optimizer.setConvergenceChecker(new StupidConvergenceChecker());
            optimizer.optimize(eval, GoalType.MINIMIZE, samples[0]);
        } catch (ConvergenceException ex) {
            ex.printStackTrace();
        }
        System.out.println("*** Best values:");
        eval.saveBest(System.out);
    } catch (Exception ex) {
        ex.printStackTrace();
        System.exit(1);
    }
}

From source file:ch.algotrader.option.SABR.java

/**
 * Perfors a SABR calibartion based on specified volatilities.
 *
 * @return SABRSmileVO The SABR smile//from  w  w  w.  j  a v a2 s . c om
 */
public static SABRSmileVO calibrate(final Double[] strikes, final Double[] volatilities, final double atmVol,
        final double forward, final double years) throws SABRException {

    MultivariateRealFunction estimateRhoAndVol = x -> {

        double r = x[0];
        double v = x[1];
        double alpha = findAlpha(forward, forward, atmVol, years, beta, x[0], x[1]);
        double sumErrors = 0;

        for (int i = 0; i < volatilities.length; i++) {

            double modelVol = vol(forward, strikes[i], years, alpha, beta, r, v);
            sumErrors += Math.pow(modelVol - volatilities[i], 2);
        }

        if (Math.abs(r) > 1) {
            sumErrors = 1e100;
        }

        return sumErrors;
    };

    NelderMead nelderMead = new NelderMead();
    RealPointValuePair result;
    try {
        result = nelderMead.optimize(estimateRhoAndVol, GoalType.MINIMIZE, new double[] { -0.5, 2.6 });
    } catch (MathException ex) {
        throw new SABRException(ex.getMessage(), ex);
    }

    double rho = result.getPoint()[0];
    double volVol = result.getPoint()[1];

    SABRSmileVO params = new SABRSmileVO();
    params.setYears(years);
    params.setRho(rho);
    params.setVolVol(volVol);
    params.setAlpha(findAlpha(forward, forward, atmVol, years, beta, rho, volVol));
    params.setAtmVol(atmVol);

    return params;
}

From source file:com.opengamma.analytics.math.minimization.NelderMeadDownhillSimplexMinimizer.java

/**
 * {@inheritDoc}//from  w  w  w.jav a2 s . c  o  m
 */
@Override
public DoubleMatrix1D minimize(final Function1D<DoubleMatrix1D, Double> function,
        final DoubleMatrix1D startPosition) {
    Validate.notNull(function, "function");
    Validate.notNull(startPosition, "start position");
    final MultivariateRealOptimizer optimizer = new NelderMead();
    final MultivariateRealFunction commonsFunction = CommonsMathWrapper.wrapMultivariate(function);
    try {
        return new DoubleMatrix1D(CommonsMathWrapper
                .unwrap(optimizer.optimize(commonsFunction, MINIMIZER, startPosition.getData())));
    } catch (final ConvergenceException e) {
        throw new MathException(e);
    } catch (final FunctionEvaluationException e) {
        throw new MathException(e);
    }
}

From source file:edu.valelab.GaussianFit.ZCalibrator.java

/**
 * Creates fitFunctionWx_ and fitFunctionWy_ based on data in data_
 * //from  w  ww .  ja va  2  s.  co m
 * 
 */
public void fitFunction() throws FunctionEvaluationException, OptimizationException {

    NelderMead nmx = new NelderMead();
    SimpleScalarValueChecker convergedChecker_ = new SimpleScalarValueChecker(1e-6, -1);

    double[][] wxData = getDataAsArray(0);
    MultiVariateZCalibrationFunction mvcx = new MultiVariateZCalibrationFunction(wxData);

    double[] params0_ = new double[5]; // initial estimates:
    params0_[0] = 37; // TODO: better estimate for c
    params0_[1] = 200; // Estimate for w0
    params0_[2] = 10; // TODO: better estimate for d
    params0_[3] = 1; // TODO: better estimate for A
    params0_[4] = 1; // TODO: better estimate for B

    nmx.setStartConfiguration(params0_);
    nmx.setConvergenceChecker(convergedChecker_);
    nmx.setMaxIterations(maxIterations_);

    double[] paramsOut = { 0.0 };

    RealPointValuePair result = nmx.optimize(mvcx, GoalType.MINIMIZE, params0_);
    paramsOut = result.getPoint();

    //for (int i = 0; i < paramsOut.length; i++) {
    //  System.out.println("Result " + i + " value: " + (int) paramsOut[i]);
    //}

    // write fit result to Results Table:
    ResultsTable res = new ResultsTable();
    res.incrementCounter();
    res.addValue("c", paramsOut[0]);
    res.addValue("w0", paramsOut[1]);
    res.addValue("d", paramsOut[2]);
    res.addValue("A", paramsOut[3]);
    res.addValue("B", paramsOut[4]);

    fitFunctionWx_ = paramsOut;

    double[][] yxData = getDataAsArray(1);
    MultiVariateZCalibrationFunction yvcx = new MultiVariateZCalibrationFunction(yxData);

    nmx.setStartConfiguration(params0_);

    result = nmx.optimize(yvcx, GoalType.MINIMIZE, params0_);
    paramsOut = result.getPoint();

    System.out.println("Y:");

    //for (int i = 0; i < paramsOut.length; i++) {
    //  System.out.println("Result " + i + " value: " + (int) paramsOut[i]);
    //}

    res.incrementCounter();
    res.addValue("c", paramsOut[0]);
    res.addValue("w0", paramsOut[1]);
    res.addValue("d", paramsOut[2]);
    res.addValue("A", paramsOut[3]);
    res.addValue("B", paramsOut[4]);

    res.show("Fit Parameters");

    fitFunctionWy_ = paramsOut;

    plotFitFunctions();

}

From source file:edu.valelab.gaussianfit.fitting.ZCalibrator.java

/**
 * Creates fitFunctionWx_ and fitFunctionWy_ based on data in data_
 * /*from   w w  w.  j a  va  2  s.  c  om*/
 * 
 * @throws org.apache.commons.math.FunctionEvaluationException
 * @throws org.apache.commons.math.optimization.OptimizationException
 */
public void fitFunction() throws FunctionEvaluationException, OptimizationException {

    NelderMead nmx = new NelderMead();
    SimpleScalarValueChecker convergedChecker_ = new SimpleScalarValueChecker(1e-6, -1);

    double[][] wxData = getDataAsArray(0);
    MultiVariateZCalibrationFunction mvcx = new MultiVariateZCalibrationFunction(wxData);

    double[] params0_ = new double[5]; // initial estimates:
    params0_[0] = 37; // TODO: better estimate for c
    params0_[1] = 200; // Estimate for w0
    params0_[2] = 10; // TODO: better estimate for d
    params0_[3] = 1; // TODO: better estimate for A
    params0_[4] = 1; // TODO: better estimate for B

    nmx.setStartConfiguration(params0_);
    nmx.setConvergenceChecker(convergedChecker_);
    nmx.setMaxIterations(maxIterations_);

    double[] paramsOut;

    RealPointValuePair result = nmx.optimize(mvcx, GoalType.MINIMIZE, params0_);
    paramsOut = result.getPoint();

    // write fit result to Results Table:
    ResultsTable res = new ResultsTable();
    res.incrementCounter();
    res.addValue("c", paramsOut[0]);
    res.addValue("w0", paramsOut[1]);
    res.addValue("d", paramsOut[2]);
    res.addValue("A", paramsOut[3]);
    res.addValue("B", paramsOut[4]);

    fitFunctionWx_ = paramsOut;

    double[][] yxData = getDataAsArray(1);
    MultiVariateZCalibrationFunction yvcx = new MultiVariateZCalibrationFunction(yxData);

    nmx.setStartConfiguration(params0_);

    result = nmx.optimize(yvcx, GoalType.MINIMIZE, params0_);
    paramsOut = result.getPoint();

    res.incrementCounter();
    res.addValue("c", paramsOut[0]);
    res.addValue("w0", paramsOut[1]);
    res.addValue("d", paramsOut[2]);
    res.addValue("A", paramsOut[3]);
    res.addValue("B", paramsOut[4]);

    res.show("Fit Parameters");

    fitFunctionWy_ = paramsOut;

    plotFitFunctions();

}

From source file:endrov.nucAutoJH.FitGaussian.java

private static double[] fitGaussian2D_(EvPixels p, double sigmaInit, final double midxInit,
        final double midyInit) {
    //sigma00, sigma01, sigma11, mu_x, mu_y, c 

    p = p.getReadOnly(EvPixelsType.DOUBLE);
    final double[] arrPixels = p.getArrayDouble();
    final int w = p.getWidth();
    final int h = p.getHeight();

    int extent = (int) Math.round(3 * sigmaInit);
    extent = Math.max(extent, 2);

    final int sx = Math.max(0, (int) (midxInit - extent));
    final int ex = Math.min(w, (int) (midxInit + extent + 1)); //+1 to the right?
    final int sy = Math.max(0, (int) (midyInit - extent));
    final int ey = Math.min(h, (int) (midyInit + extent + 1));

    double minIntensity = Double.MAX_VALUE;
    double maxIntensity = Double.MIN_VALUE;
    for (int y = sy; y < ey; y++) {
        int base = y * w;
        double dy2 = y - midyInit;
        dy2 = dy2 * dy2;/*www.  ja  v a2  s . c om*/
        for (int x = sx; x < ex; x++) {
            double dx2 = x - midxInit;
            dx2 = dx2 * dx2;
            double t = arrPixels[base + x];
            //if(dx2+dy2<=extent*extent)
            {
                if (t < minIntensity)
                    minIntensity = t;
                if (t > maxIntensity)
                    maxIntensity = t;
            }
        }
    }

    //double[] weights=new double[]{1};
    double[] startPoint = new double[] { sigmaInit, 0, sigmaInit, midxInit, midyInit, minIntensity,
            maxIntensity - minIntensity };
    //double[] output=new double[startPoint.length];

    try {
        MultivariateRealFunction func = new MultivariateRealFunction() {
            //      opt.optimize(

            public double value(double[] arg) throws FunctionEvaluationException, IllegalArgumentException {
                double sigma00 = arg[0];
                double sigma01 = arg[1];
                double sigma11 = arg[2];
                double mu0 = arg[3];
                double mu1 = arg[4];
                double C = arg[5];
                double D = arg[6];

                double sumError = 0;

                Matrix2d sigma = new Matrix2d(sigma00, sigma01, sigma01, sigma11);
                Matrix2d sigmaInv = new Matrix2d();
                sigma.invert(sigmaInv);
                double sigmaDet = sigma.determinant();
                double front = 1.0 / (2 * Math.PI * Math.sqrt(sigmaDet));
                //System.out.println("front: "+front);
                //System.out.println("sigma inv "+sigmaInv);

                if (mu0 < sx || mu0 > ex)
                    sumError += 1000000;
                if (mu1 < sy || mu1 > ey)
                    sumError += 1000000;
                if (sigma00 < 1)
                    sumError += 1000000;
                //if(sigma01<0)               sumError+=1000000;
                if (sigma11 < 1)
                    sumError += 1000000;
                if (D <= 0)
                    sumError += 1000000;

                for (int y = sy; y < ey; y++) {
                    int base = y * w;
                    double dy2 = y - midyInit;
                    dy2 = dy2 * dy2;
                    for (int x = sx; x < ex; x++) {
                        double dx2 = x - midxInit;
                        dx2 = dx2 * dx2;
                        double thisReal = arrPixels[base + x];
                        //                  if(dx2+dy2<=extent*extent)
                        {
                            //               DoubleMatrix2D sigma=new DenseDoubleMatrix2D(new double[][]{{sigma00,sigma01},{sigma01,sigma11}});
                            //double sigmaDet=sigma00*sigma11-sigma01*sigma01;

                            double dx0 = x - mu0;
                            double dx1 = y - mu1;

                            //http://en.wikipedia.org/wiki/Multivariate_normal_distribution

                            Vector2d vX = new Vector2d(dx0, dx1);
                            Vector2d op = new Vector2d(vX);
                            sigmaInv.transform(op);
                            double upper = -0.5 * op.dot(vX);
                            double exp = Math.exp(upper);

                            //System.out.println("front "+front+" "+exp+" C "+C+" thisreal"+thisReal+" upper "+upper);

                            if (upper > -0.4)
                                exp = 1;
                            else
                                exp = Math.max(0, 1 + upper + 0.4);

                            /*
                            if(exp<0.7)
                               exp=0;
                            else
                               exp=1;
                            */

                            double thisExpected = D * front * exp + C;
                            double diff = thisExpected - thisReal;
                            sumError += diff * diff;

                        }
                    }
                }

                //System.out.println(sigma00+"\t"+sigma01+"\t"+sigma11+"\tC"+C+"\tmu "+mu0+","+mu1+"\terr "+sumError);
                return sumError;
                //            return new double[]{sumError};
            }

        };

        NelderMead opt = new NelderMead();
        //LevenbergMarquardtOptimizer opt=new LevenbergMarquardtOptimizer();
        opt.setMaxIterations(10000);
        RealPointValuePair pair = opt.optimize(func, GoalType.MINIMIZE, startPoint);

        int numit = opt.getIterations();
        System.out.println("#it " + numit);
        System.out.println("err " + func.value(pair.getPointRef()));
        return pair.getPointRef();

        //         for(int i=0;i<startPoint.length;i++)
        //         System.out.println("i: "+i+"  "+output[i]);
        //, output, weights, startPoint);
    }
    /*
    catch (MaxIterationsExceededException e)
       {
       System.out.println("max it reached");
               
       }*/
    catch (Exception e) {
        e.printStackTrace();
    }

    //Maybe this is a bad point?
    System.out.println("max it reached");
    return startPoint;
    //      return output;
}

From source file:endrov.typeLineageAutoNucJH.FitGaussian.java

private static double[] fitGaussian2D_(EvPixels p, double sigmaInit, final double midxInit,
        final double midyInit) {
    //sigma00, sigma01, sigma11, mu_x, mu_y, c 

    p = p.getReadOnly(EvPixelsType.DOUBLE);
    final double[] arrPixels = p.getArrayDouble();
    final int w = p.getWidth();
    final int h = p.getHeight();

    int extent = (int) Math.round(3 * sigmaInit);
    extent = Math.max(extent, 2);

    final int sx = Math.max(0, (int) (midxInit - extent));
    final int ex = Math.min(w, (int) (midxInit + extent + 1)); //+1 to the right?
    final int sy = Math.max(0, (int) (midyInit - extent));
    final int ey = Math.min(h, (int) (midyInit + extent + 1));

    double minIntensity = Double.MAX_VALUE;
    double maxIntensity = -Double.MAX_VALUE;
    for (int y = sy; y < ey; y++) {
        int base = y * w;
        double dy2 = y - midyInit;
        dy2 = dy2 * dy2;/*from  w ww  .  j a  v  a  2  s .  co m*/
        for (int x = sx; x < ex; x++) {
            double dx2 = x - midxInit;
            dx2 = dx2 * dx2;
            double t = arrPixels[base + x];
            //if(dx2+dy2<=extent*extent)
            {
                if (t < minIntensity)
                    minIntensity = t;
                if (t > maxIntensity)
                    maxIntensity = t;
            }
        }
    }

    //double[] weights=new double[]{1};
    double[] startPoint = new double[] { sigmaInit, 0, sigmaInit, midxInit, midyInit, minIntensity,
            maxIntensity - minIntensity };
    //double[] output=new double[startPoint.length];

    try {
        MultivariateRealFunction func = new MultivariateRealFunction() {
            //      opt.optimize(

            public double value(double[] arg) throws FunctionEvaluationException, IllegalArgumentException {
                double sigma00 = arg[0];
                double sigma01 = arg[1];
                double sigma11 = arg[2];
                double mu0 = arg[3];
                double mu1 = arg[4];
                double C = arg[5];
                double D = arg[6];

                double sumError = 0;

                Matrix2d sigma = new Matrix2d(sigma00, sigma01, sigma01, sigma11);
                Matrix2d sigmaInv = new Matrix2d();
                sigma.invert(sigmaInv);
                double sigmaDet = sigma.determinant();
                double front = 1.0 / (2 * Math.PI * Math.sqrt(sigmaDet));
                //System.out.println("front: "+front);
                //System.out.println("sigma inv "+sigmaInv);

                if (mu0 < sx || mu0 > ex)
                    sumError += 1000000;
                if (mu1 < sy || mu1 > ey)
                    sumError += 1000000;
                if (sigma00 < 1)
                    sumError += 1000000;
                //if(sigma01<0)               sumError+=1000000;
                if (sigma11 < 1)
                    sumError += 1000000;
                if (D <= 0)
                    sumError += 1000000;

                for (int y = sy; y < ey; y++) {
                    int base = y * w;
                    double dy2 = y - midyInit;
                    dy2 = dy2 * dy2;
                    for (int x = sx; x < ex; x++) {
                        double dx2 = x - midxInit;
                        dx2 = dx2 * dx2;
                        double thisReal = arrPixels[base + x];
                        //                  if(dx2+dy2<=extent*extent)
                        {
                            //               DoubleMatrix2D sigma=new DenseDoubleMatrix2D(new double[][]{{sigma00,sigma01},{sigma01,sigma11}});
                            //double sigmaDet=sigma00*sigma11-sigma01*sigma01;

                            double dx0 = x - mu0;
                            double dx1 = y - mu1;

                            //http://en.wikipedia.org/wiki/Multivariate_normal_distribution

                            Vector2d vX = new Vector2d(dx0, dx1);
                            Vector2d op = new Vector2d(vX);
                            sigmaInv.transform(op);
                            double upper = -0.5 * op.dot(vX);
                            double exp = Math.exp(upper);

                            //System.out.println("front "+front+" "+exp+" C "+C+" thisreal"+thisReal+" upper "+upper);

                            if (upper > -0.4)
                                exp = 1;
                            else
                                exp = Math.max(0, 1 + upper + 0.4);

                            /*
                            if(exp<0.7)
                               exp=0;
                            else
                               exp=1;
                            */

                            double thisExpected = D * front * exp + C;
                            double diff = thisExpected - thisReal;
                            sumError += diff * diff;

                        }
                    }
                }

                //System.out.println(sigma00+"\t"+sigma01+"\t"+sigma11+"\tC"+C+"\tmu "+mu0+","+mu1+"\terr "+sumError);
                return sumError;
                //            return new double[]{sumError};
            }

        };

        NelderMead opt = new NelderMead();
        //LevenbergMarquardtOptimizer opt=new LevenbergMarquardtOptimizer();
        opt.setMaxIterations(10000);
        RealPointValuePair pair = opt.optimize(func, GoalType.MINIMIZE, startPoint);

        int numit = opt.getIterations();
        System.out.println("#it " + numit);
        System.out.println("err " + func.value(pair.getPointRef()));
        return pair.getPointRef();

        //         for(int i=0;i<startPoint.length;i++)
        //         System.out.println("i: "+i+"  "+output[i]);
        //, output, weights, startPoint);
    }
    /*
    catch (MaxIterationsExceededException e)
       {
       System.out.println("max it reached");
               
       }*/
    catch (Exception e) {
        e.printStackTrace();
    }

    //Maybe this is a bad point?
    System.out.println("max it reached");
    return startPoint;
    //      return output;
}

From source file:edu.valelab.gaussianfit.fitting.ZCalibrator.java

/**
 * Use the fitfunction to estimate the z position given width in x and y
 * //from  w  w w  .ja v a2  s  .c o m
 * minimize the distance D in sqrt wx and sqrt wy space
 * D = sqrt (  square (sqrt wx - sqrt wx, calib) + sqr(sqrt wy - sqrt w, calib) )
 * 
 * 
 * @param wx - width in x
 * @param wy - width in y
 * @return - calculated z position
 */

public double getZ(double wx, double wy) {
    if (!hasFitFunctions())
        return 0.0;

    NelderMead nmx = new NelderMead();
    SimpleScalarValueChecker convergedChecker_ = new SimpleScalarValueChecker(1e-6, -1);

    MultiVariateZFunction mz = new MultiVariateZFunction(fitFunctionWx_, fitFunctionWy_, wx, wy);

    double[] params0_ = new double[1]; // initial estimates:
    params0_[0] = 15; // TODO: Need the middle z value of the stack here!!!

    nmx.setStartConfiguration(params0_);
    nmx.setConvergenceChecker(convergedChecker_);
    nmx.setMaxIterations(maxIterations_);

    double[] paramsOut = { 0.0 };

    try {
        RealPointValuePair result = nmx.optimize(mz, GoalType.MINIMIZE, params0_);
        paramsOut = result.getPoint();
    } catch (java.lang.OutOfMemoryError e) {
        throw (e);
    } catch (FunctionEvaluationException e) {
        ij.IJ.log(" " + e.toString());
    } catch (OptimizationException e) {
        ij.IJ.log(" " + e.toString());
    } catch (IllegalArgumentException e) {
        ij.IJ.log(" " + e.toString());
    }

    return paramsOut[0];
}

From source file:edu.valelab.GaussianFit.ZCalibrator.java

/**
 * Use the fitfunction to estimate the z position given width in x and y
 * //  w  w  w. j a  v  a  2s .c  om
 * minimize the distance D in sqrt wx and sqrt wy space
 * D = sqrt (  square (sqrt wx - sqrt wx, calib) + sqr(sqrt wy - sqrt w, calib) )
 * 
 * 
 */

public double getZ(double wx, double wy) {
    if (!hasFitFunctions())
        return 0.0;

    NelderMead nmx = new NelderMead();
    SimpleScalarValueChecker convergedChecker_ = new SimpleScalarValueChecker(1e-6, -1);

    MultiVariateZFunction mz = new MultiVariateZFunction(fitFunctionWx_, fitFunctionWy_, wx, wy);

    double[] params0_ = new double[1]; // initial estimates:
    params0_[0] = 15; // TODO: Need the middle z value of the stack here!!!

    nmx.setStartConfiguration(params0_);
    nmx.setConvergenceChecker(convergedChecker_);
    nmx.setMaxIterations(maxIterations_);

    double[] paramsOut = { 0.0 };

    try {
        RealPointValuePair result = nmx.optimize(mz, GoalType.MINIMIZE, params0_);
        paramsOut = result.getPoint();
    } catch (java.lang.OutOfMemoryError e) {
        throw (e);
    } catch (Exception e) {
        ij.IJ.log(" " + e.toString());
    }

    return paramsOut[0];
}

From source file:org.renjin.primitives.optimize.Optimizations.java

/**
 * General-purpose optimization based on NelderMead, quasi-Newton and conjugate-gradient algorithms.
 * It includes an option for box-constrained
 * optimization and simulated annealing.
 *
 * @param par initial parameters//  ww  w.ja  v  a 2s  .co  m
 * @param fn
 * @param gradientFunction
 * @param method
 * @param controlParameters
 * @param lower
 * @param upper
 * @return
 */
@Primitive
public static ListVector optim(@Current Context context, @Current Environment rho, DoubleVector par,
        Function fn, SEXP gradientFunction, String method, ListVector controlParameters, DoubleVector lower,
        DoubleVector upper) {

    MultivariateRealClosure g = new MultivariateRealClosure(context, rho, fn);

    if (method.equals("Nelder-Mead")) {

        NelderMead optimizer = new NelderMead();
        try {
            RealPointValuePair res = optimizer.optimize(g, GoalType.MINIMIZE, par.toDoubleArray());
            ListVector.Builder result = new ListVector.Builder();
            result.add(new DoubleArrayVector(res.getPoint()));
            result.add(new DoubleArrayVector(res.getValue()));
            result.add(new IntArrayVector(IntVector.NA, IntVector.NA));
            result.add(new IntArrayVector(0));
            result.add(Null.INSTANCE);
            return result.build();

        } catch (FunctionEvaluationException e) {
            throw new EvalException(e);
        } catch (OptimizationException e) {
            throw new EvalException(e);
        }
    } else {
        throw new EvalException("method '" + method + "' not implemented.");
    }
}