Example usage for org.apache.commons.math3.optim.nonlinear.scalar.noderiv NelderMeadSimplex NelderMeadSimplex

List of usage examples for org.apache.commons.math3.optim.nonlinear.scalar.noderiv NelderMeadSimplex NelderMeadSimplex

Introduction

In this page you can find the example usage for org.apache.commons.math3.optim.nonlinear.scalar.noderiv NelderMeadSimplex NelderMeadSimplex.

Prototype

public NelderMeadSimplex(final double[][] referenceSimplex) 

Source Link

Document

Build a Nelder-Mead simplex with default coefficients.

Usage

From source file:eu.crisis_economics.abm.markets.clearing.heterogeneous.NelderMeadClearingAlgorithm.java

@Override
public double applyToNetwork(final MixedClearingNetwork network) {
    Preconditions.checkNotNull(network);
    final SimplexOptimizer optimizer = new SimplexOptimizer(relErrorTarget, absErrorTarget);

    final ResidualCostFunction aggregateCostFunction = super.getResidualScalarCostFunction(network);
    final RealVector start = new ArrayRealVector(network.getNumberOfEdges());
    for (int i = 0; i < network.getNumberOfEdges(); ++i)
        start.setEntry(i, network.getEdges().get(i).getMaximumRateAdmissibleByBothParties());
    start.set(1.);/*  w ww .  j ava2  s. com*/

    final PointValuePair result = optimizer.optimize(new MaxEval(maximumEvaluations),
            new ObjectiveFunction(aggregateCostFunction), GoalType.MINIMIZE, new InitialGuess(start.toArray()),
            new NelderMeadSimplex(network.getNumberOfEdges()));

    final double residualCost = result.getValue();
    System.out.println("Network cleared: residual cost: " + residualCost + ".");

    return residualCost;
}

From source file:de.bund.bfr.math.MultivariateOptimization.java

@Override
public Result optimize(int nParameterSpace, int nOptimizations, boolean stopWhenSuccessful,
        Map<String, Double> minStartValues, Map<String, Double> maxStartValues, int maxIterations,
        DoubleConsumer progessListener, ExecutionContext exec) throws CanceledExecutionException {
    if (exec != null) {
        exec.checkCanceled();/*from   w  w w . ja v  a  2s  . c  o  m*/
    }

    progessListener.accept(0.0);

    List<ParamRange> ranges = MathUtils.getParamRanges(parameters, minStartValues, maxStartValues,
            nParameterSpace);

    ranges.set(parameters.indexOf(sdParam), new ParamRange(1.0, 1, 1.0));

    List<StartValues> startValuesList = MathUtils.createStartValuesList(ranges, nOptimizations,
            values -> optimizerFunction.value(Doubles.toArray(values)),
            progress -> progessListener.accept(0.5 * progress), exec);
    Result result = new Result();
    AtomicInteger currentIteration = new AtomicInteger();
    SimplexOptimizer optimizer = new SimplexOptimizer(new SimpleValueChecker(1e-10, 1e-10) {

        @Override
        public boolean converged(int iteration, PointValuePair previous, PointValuePair current) {
            if (super.converged(iteration, previous, current)) {
                return true;
            }

            return currentIteration.incrementAndGet() >= maxIterations;
        }
    });
    int count = 0;

    for (StartValues startValues : startValuesList) {
        if (exec != null) {
            exec.checkCanceled();
        }

        progessListener.accept(0.5 * count++ / startValuesList.size() + 0.5);

        try {
            PointValuePair optimizerResults = optimizer.optimize(new MaxEval(Integer.MAX_VALUE),
                    new MaxIter(maxIterations), new InitialGuess(Doubles.toArray(startValues.getValues())),
                    new ObjectiveFunction(optimizerFunction), GoalType.MAXIMIZE,
                    new NelderMeadSimplex(parameters.size()));
            double logLikelihood = optimizerResults.getValue() != null ? optimizerResults.getValue()
                    : Double.NaN;

            if (result.logLikelihood == null || logLikelihood > result.logLikelihood) {
                result = getResults(optimizerResults);

                if (result.logLikelihood == 0.0 || stopWhenSuccessful) {
                    break;
                }
            }
        } catch (TooManyEvaluationsException | TooManyIterationsException | ConvergenceException e) {
        }
    }

    return result;
}

From source file:hulo.localization.models.obs.GaussianProcessLDPLMean.java

static PointValuePair minimize(MultivariateFunction func, double[] pointInit, double[] diffPointInit) {
    SimplexOptimizer optimizer = new SimplexOptimizer(1e-5, 1e-10);
    NelderMeadSimplex simplex = new NelderMeadSimplex(diffPointInit);
    ObjectiveFunction objFunc = new ObjectiveFunction(func);
    InitialGuess initGuess = new InitialGuess(pointInit);
    PointValuePair pair = optimizer.optimize(maxEval, objFunc, GoalType.MINIMIZE, initGuess, simplex);
    return pair;//from w ww. j a  va2 s  .  c  o m
}

From source file:edu.cmu.tetrad.sem.GeneralizedSemEstimator.java

private double[] optimize(MultivariateFunction function, double[] values, int optimizer) {
    PointValuePair pair = null;/*from w  w  w.ja  va2s. co m*/

    if (optimizer == 1) {
        //            0.01, 0.000001
        //2.0D * FastMath.ulp(1.0D), 1e-8
        MultivariateOptimizer search = new PowellOptimizer(1e-7, 1e-7);
        pair = search.optimize(new InitialGuess(values), new ObjectiveFunction(function), GoalType.MINIMIZE,
                new MaxEval(100000));
    } else if (optimizer == 2) {
        MultivariateOptimizer search = new SimplexOptimizer(1e-7, 1e-7);
        pair = search.optimize(new InitialGuess(values), new ObjectiveFunction(function), GoalType.MINIMIZE,
                new MaxEval(100000), new NelderMeadSimplex(values.length));
    } else if (optimizer == 3) {
        int dim = values.length;
        int additionalInterpolationPoints = 0;
        final int numIterpolationPoints = 2 * dim + 1 + additionalInterpolationPoints;

        BOBYQAOptimizer search = new BOBYQAOptimizer(numIterpolationPoints);
        pair = search.optimize(new MaxEval(100000), new ObjectiveFunction(function), GoalType.MINIMIZE,
                new InitialGuess(values), SimpleBounds.unbounded(dim));
    } else if (optimizer == 4) {
        MultivariateOptimizer search = new CMAESOptimizer(3000000, .05, false, 0, 0, new MersenneTwister(),
                false, new SimplePointChecker<PointValuePair>(0.5, 0.5));
        pair = search.optimize(new MaxEval(30000), new ObjectiveFunction(function), GoalType.MINIMIZE,
                new InitialGuess(values), new CMAESOptimizer.Sigma(new double[values.length]),
                new CMAESOptimizer.PopulationSize(1000));
    } else if (optimizer == 5) {
        //            0.01, 0.000001
        //2.0D * FastMath.ulp(1.0D), 1e-8
        MultivariateOptimizer search = new PowellOptimizer(.05, .05);
        pair = search.optimize(new InitialGuess(values), new ObjectiveFunction(function), GoalType.MINIMIZE,
                new MaxEval(100000));
    } else if (optimizer == 6) {
        MultivariateOptimizer search = new PowellOptimizer(1e-7, 1e-7);
        pair = search.optimize(new InitialGuess(values), new ObjectiveFunction(function), GoalType.MAXIMIZE,
                new MaxEval(10000));
    }

    return pair.getPoint();
}

From source file:gdsc.smlm.ij.plugins.pcpalm.PCPALMMolecules.java

private double[] optimiseSimplex(float[] x, float[] y, double[] initialSolution) {
    // Simplex optimisation
    SkewNormalMultivariateFunction sn2 = new SkewNormalMultivariateFunction(initialSolution);
    sn2.addData(x, y);// w w w  .j a  va2  s .  c o  m
    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;
}

From source file:nl.systemsgenetics.cellTypeSpecificAlleleSpecificExpression.BetaBinomialTest.java

public BetaBinomialTest(ArrayList<IndividualSnpData> all_individuals) throws Exception {
    boolean debug = true;
    //basic information, get the zero instance, was checked as the same in ReadAsLinesIntoIndividualSNPdata
    snpName = all_individuals.get(0).getSnpName();
    chromosome = all_individuals.get(0).getChromosome();
    position = all_individuals.get(0).getPosition();
    genotype = all_individuals.get(0).genotype;

    ArrayList<IndividualSnpData> het_individuals;
    het_individuals = UtilityMethods.isolateHeterozygotesFromIndividualSnpData(all_individuals);

    numberOfHets = het_individuals.size();

    hetSampleNames = new ArrayList<String>();
    asRef = new ArrayList<Integer>();
    asAlt = new ArrayList<Integer>();
    asNo = new ArrayList<Integer>();

    dispersion = new ArrayList<Double>();

    int total_overlap = 0;

    for (IndividualSnpData temp_het : het_individuals) {
        //Do nothing if there is no data in het_individuals

        hetSampleNames.add(temp_het.getSampleName());

        asRef.add(temp_het.getRefNum());
        asAlt.add(temp_het.getAltNum());
        asNo.add(temp_het.getNoNum());//from   w  w  w  .  ja  va  2 s . c  om

        dispersion.add(temp_het.getDispersion());

        //this is used to check if we will continue with calculations.
        total_overlap += temp_het.getRefNum() + temp_het.getAltNum();

    }

    if ((total_overlap >= GlobalVariables.minReads) && (numberOfHets >= GlobalVariables.minHets)) {
        // There is data to perform the binomial test, perform it.       
        if (GlobalVariables.verbosity >= 100) {
            System.out.println();
            System.out.println("---- Starting beta binomial LRT test estimate ----");
            System.out.println("SNP name: " + snpName);
            System.out.println("at: chr" + chromosome + ":" + position);
            System.out.println("--------------------------------------------------");

            System.out.println("debug:");
            System.out.println("Num of hets: " + Integer.toString(numberOfHets));
            System.out.println(het_individuals.get(0).getSnpName());
            System.out.println(total_overlap);

            System.out.println("asRef:       " + asRef.toString());
            System.out.println("asAlt:       " + asAlt.toString());
            System.out.println("dispersion:  " + dispersion.toString());
            System.out.println("Starting Null estimation");
        }

        Integer[] asRefArray = asRef.toArray(new Integer[asRef.size()]);
        Integer[] asAltArray = asAlt.toArray(new Integer[asAlt.size()]);
        Double[] dispArray = dispersion.toArray(new Double[dispersion.size()]);

        BetaBinomNullLikelihood betaBinomNull;
        betaBinomNull = new BetaBinomNullLikelihood(asRefArray, asAltArray, dispArray);

        nullLogLik = betaBinomNull.value(new double[] { 0.5 });

        if (GlobalVariables.verbosity >= 100) {
            System.out.println("Starting Alt estimation");
        }

        BetaBinomAltLikelihood betaBinomAlt;
        betaBinomAlt = new BetaBinomAltLikelihood(asRefArray, asAltArray, dispArray);
        NelderMeadSimplex simplex;
        simplex = new NelderMeadSimplex(2);
        SimplexOptimizer optimizer = new SimplexOptimizer(GlobalVariables.simplexThreshold,
                GlobalVariables.simplexThreshold); //numbers are to which precision you want it to be done.
        PointValuePair solutionAlt = optimizer.optimize(new ObjectiveFunction(betaBinomAlt),
                new MaxEval(GlobalVariables.maximumIterations), simplex, GoalType.MINIMIZE,
                new InitialGuess(new double[] { 0.5, 0.5 }), new SearchInterval(-1000.0, 1000.0));

        double[] valueAlt = solutionAlt.getPoint();

        altLogLik = betaBinomAlt.value(valueAlt);
        iterations = optimizer.getIterations();
        alphaParam = valueAlt[0];
        betaParam = valueAlt[1];
        binomRatio = valueAlt[0] / (valueAlt[0] + valueAlt[1]);

        //chi squared statistic is determined based on both null and alt loglikelihoods.
        chiSq = 2.0 * (nullLogLik - altLogLik);

        //determine P value based on distribution
        ChiSquaredDistribution distribution = new ChiSquaredDistribution(1);
        pVal = 1 - distribution.cumulativeProbability(chiSq);

        if (GlobalVariables.verbosity >= 10) {
            System.out.println("\n--- Beta Binomial Test Statistics: ---");
            System.out.println("LogLik of Alt converged to a threshold of "
                    + Double.toString(GlobalVariables.simplexThreshold));
            System.out.println("\tAlpha parameter:      " + Double.toString(valueAlt[0]));
            System.out.println("\tBeta parameter:       " + Double.toString(valueAlt[1]));
            System.out.println("\tIterations to converge:           " + Integer.toString(iterations) + "\n");
            System.out.println("\tNull log likelihood:  " + Double.toString(nullLogLik));
            System.out.println("\tAlt log likelihood:   " + Double.toString(altLogLik) + "\n");
            System.out.println("\tChisq statistic:      " + Double.toString(chiSq));
            System.out.println("\tP value:              " + Double.toString(pVal));
            System.out.println("------------------------------------------");
            //TODO, I want to format this properly, but not necessary

        }

        testPerformed = true;

        //This below is to compare everything to python created by WASP
        //Will put the backup back into python, just to make sure.

        if (GlobalVariables.verbosity >= 100) {
            System.out.println("loglikRef = " + Double.toString(altLogLik));
            System.out.println("loglik = 0");
            for (int i = 0; i < asRef.size(); i++) {

                System.out.println("loglik += AS_betabinom_loglike(" + "[math.log("
                        + Double.toString(alphaParam) + ") - math.log(" + Double.toString(alphaParam) + " + "
                        + Double.toString(betaParam) + "), " + "math.log(" + Double.toString(betaParam)
                        + ") - math.log(" + Double.toString(alphaParam) + " + " + Double.toString(betaParam)
                        + ")]," + Double.toString(dispersion.get(i)) + "," + Double.toString(asRef.get(i)) + ","
                        + Double.toString(asAlt.get(i)) + "," + "0.980198, 0.005)");
            }
            System.out.println("if allclose(-1.0 * loglik, loglikRef): correctTest += 1");

            //I find exactly the correct values after running it through WASP implementation in python, after one check.
            //Using this as a foundation stone to write some tests:
            for (int i = 0; i < asRef.size(); i++) {
                String command = "Assert.assertEquals( " + "likelihoodFunctions.BetaBinomLogLik("
                        + String.format("%.15g", dispArray[i]) + "," + String.format("%.15g", alphaParam) + ","
                        + String.format("%.15g", betaParam) + "," + Integer.toString(asRefArray[i]) + ","
                        + Integer.toString(asAltArray[i]) + ")" + ", "
                        + Double.toString(LikelihoodFunctions.BetaBinomLogLik(dispArray[i], alphaParam,
                                betaParam, asRefArray[i], asAltArray[i]))
                        + ");";
                System.out.println(command);
            }
        }

    }
}

From source file:nl.systemsgenetics.cellTypeSpecificAlleleSpecificExpression.BetaBinomOverdispInSample.java

@SuppressWarnings({ "empty-statement", "empty-statement" })
public BetaBinomOverdispInSample(String filename) throws FileNotFoundException, IOException {

    sampleName = filename;/*from w  w w  .j a v a2 s  .c  om*/

    if (GlobalVariables.verbosity >= 10) {
        System.out.println("\n--- Starting beta binomial dispersion estimate ---");
        System.out.println("AS File: " + sampleName);
        System.out.println("--------------------------------------------------");
    }

    ArrayList<IndividualSnpData> allSnps = new ArrayList<IndividualSnpData>();

    File file = new File(filename);
    FileReader fr = new FileReader(file);
    BufferedReader br = new BufferedReader(fr);
    String line;
    while ((line = br.readLine()) != null) {
        allSnps.add(new IndividualSnpData(file.getAbsolutePath(), line));
    }
    br.close();
    fr.close();

    ArrayList<IndividualSnpData> hets;
    hets = UtilityMethods.isolateHeterozygotesFromIndividualSnpData(allSnps);

    int numOfHets = hets.size();

    int[] asRef = new int[numOfHets];
    int[] asAlt = new int[numOfHets];

    int totalRef = 0;
    int totalAlt = 0;

    for (int i = 0; i < numOfHets; i++) {
        asRef[i] = hets.get(i).getRefNum();
        asAlt[i] = hets.get(i).getAltNum();

        totalRef += asRef[i];
        totalAlt += asAlt[i];

    }
    if (GlobalVariables.verbosity >= 10) {
        System.out.println("sample loaded");

        System.out.println("\ttotal_ref = " + totalRef);
        System.out.println("\ttotal_alt = " + totalAlt);
    }

    BetaBinomLikelihoodForOverdispersion betaBinom = new BetaBinomLikelihoodForOverdispersion(asRef, asAlt);
    NelderMeadSimplex simplex;
    simplex = new NelderMeadSimplex(1);
    SimplexOptimizer optimizer = new SimplexOptimizer(GlobalVariables.simplexThreshold,
            GlobalVariables.simplexThreshold); //numbers are to which precision you want it to be done.
    PointValuePair solution = optimizer.optimize(new ObjectiveFunction(betaBinom),
            new MaxEval(GlobalVariables.maximumIterations), simplex, GoalType.MINIMIZE,
            new InitialGuess(new double[] { 0.5 }), new SearchInterval(0.0, 1.0));

    overdispersion = solution.getFirst();
    if (GlobalVariables.verbosity >= 10) {
        System.out.println("Log likelihood converged to a threshold of "
                + Double.toString(GlobalVariables.simplexThreshold));
        System.out.println("\tDispersion sigma: " + Double.toString(overdispersion[0]));
        System.out.println("\tLog likelihood:   " + Double.toString(betaBinom.value(overdispersion)));
    }

}

From source file:nl.systemsgenetics.cellTypeSpecificAlleleSpecificExpression.CTSBetaBinomialTest.java

public CTSBetaBinomialTest(ArrayList<IndividualSnpData> all_individuals) throws Exception {

    //basic information, get the zero instance.
    snpName = all_individuals.get(0).getSnpName();
    chromosome = all_individuals.get(0).getChromosome();
    position = all_individuals.get(0).getPosition();

    //isolate heterozygote individuals.
    ArrayList<IndividualSnpData> het_individuals = UtilityMethods
            .isolateHeterozygotesFromIndividualSnpData(all_individuals);

    numberOfHets = het_individuals.size();

    hetSampleNames = new ArrayList<String>();
    asRef = new ArrayList<Integer>();
    asAlt = new ArrayList<Integer>();
    asNo = new ArrayList<Integer>();

    dispersion = new ArrayList<Double>();
    cellProp = new ArrayList<Double>();
    int total_overlap = 0;

    //Get the basic data without doing any tests.
    for (IndividualSnpData temp_het : het_individuals) {
        //Do nothing if there is no data in het_individuals

        hetSampleNames.add(temp_het.getSampleName());

        asRef.add(temp_het.getRefNum());
        asAlt.add(temp_het.getAltNum());
        asNo.add(temp_het.getNoNum());/*  w  w w.j a  va 2  s .  c  o  m*/

        dispersion.add(temp_het.getDispersion());
        cellProp.add(temp_het.getCellTypeProp());

        //this is used to check if we will continue with calculations.
        //BASED on the minHets and MinReads
        total_overlap += temp_het.getRefNum() + temp_het.getAltNum();
    }

    //Check if we do a test.
    if ((total_overlap >= GlobalVariables.minReads) && (numberOfHets >= GlobalVariables.minHets)) {
        // There is data to perform the binomial test, perform it.       

        if (GlobalVariables.verbosity >= 100) {
            System.out.println();

            System.out.println("SNP name: " + snpName);
            System.out.println("at: chr" + chromosome + ":" + position);
            System.out.println("------------------------------------------------------");
            System.out.println("Num of hets: " + Integer.toString(numberOfHets));
            System.out.println(het_individuals.get(0).getSnpName());
            System.out.println(total_overlap);
            System.out.println("asRef:       " + asRef.toString());
            System.out.println("asAlt:       " + asAlt.toString());
            System.out.println("dispersion:  " + dispersion.toString());
            System.out.println("cellProp:    " + cellProp.toString());
            System.out.println("Starting non-CTS Beta binomial Null estimation");

        }

        /*
        Going to do this in a two step procedure:
                
        1. Determine a   2 parameter null the same way as what is done 
                in the beta binomial alt, with starting values:
                This will be the actual null model in the CTS beta binomial output
        2. Determine the cell type specific proportion, but now we use a
                cell type specific model as the alternative.
                        
                
        */

        Integer[] asRefArray = asRef.toArray(new Integer[asRef.size()]);
        Integer[] asAltArray = asAlt.toArray(new Integer[asAlt.size()]);
        Double[] dispArray = dispersion.toArray(new Double[dispersion.size()]);
        Double[] cellPropArray = cellProp.toArray(new Double[cellProp.size()]);

        try {
            BetaBinomAltLikelihood betaBinomNull;
            betaBinomNull = new BetaBinomAltLikelihood(asRefArray, asAltArray, dispArray);
            NelderMeadSimplex simplex;
            simplex = new NelderMeadSimplex(2);
            SimplexOptimizer optimizer = new SimplexOptimizer(GlobalVariables.simplexThreshold,
                    GlobalVariables.simplexThreshold); //numbers are to which precision you want it to be done.
            PointValuePair solutionNull = optimizer.optimize(new ObjectiveFunction(betaBinomNull),
                    new MaxEval(GlobalVariables.maximumIterations), simplex, GoalType.MINIMIZE,
                    new InitialGuess(new double[] { 0.5, 0.5 }), new SearchInterval(-1000.0, 1000.0));

            double[] valueNull = solutionNull.getPoint();

            nullLogLik = betaBinomNull.value(valueNull);
            nulliterations = optimizer.getIterations();
            NullAlphaParam = valueNull[0];
            NullBetaParam = valueNull[1];
            NullbinomRatio = valueNull[0] / (valueNull[0] + valueNull[1]);

            if (GlobalVariables.verbosity >= 100) {
                System.out.println("LogLik of Null converged to a threshold of "
                        + Double.toString(GlobalVariables.simplexThreshold));
                System.out.println("\tNull Alpha parameter:      " + Double.toString(valueNull[0]));
                System.out.println("\tNull Beta parameter:       " + Double.toString(valueNull[1]));
                System.out.println("\tIterations to converge:    " + Integer.toString(nulliterations) + "\n");
            }

            //CHECK WHAT THE version does in terms of loglik.
            CTSbetaBinomialAltLikelihoodVersion2 CTSbetaBinomAlt;
            CTSbetaBinomAlt = new CTSbetaBinomialAltLikelihoodVersion2(asRefArray, asAltArray, dispArray,
                    cellPropArray);

            // Going to make some initialGuesses, because sometimes 
            // this falls into a local minimum with the MLE, when only the 
            // null parameters are used as a starting point.

            double nonCTSprop;
            nonCTSprop = valueNull[0] / (valueNull[0] + valueNull[1]);

            ArrayList<InitialGuess> GuessList = new ArrayList<InitialGuess>();

            GuessList.add(new InitialGuess(new double[] { 0.0, nonCTSprop }));
            GuessList.add(new InitialGuess(new double[] { 0.25, 0.5 }));
            GuessList.add(new InitialGuess(new double[] { 0.75, 0.5 }));

            simplex = new NelderMeadSimplex(2);

            ArrayList<double[]> OptimizerResults = new ArrayList<double[]>();
            Double[] OptimizedLogLik = new Double[8];

            int i = 0;
            int lowestIndices = -1;
            double lowestLogLik = 0;

            for (InitialGuess IGuess : GuessList) {

                PointValuePair solutionAlt = optimizer.optimize(new ObjectiveFunction(CTSbetaBinomAlt),
                        new MaxEval(500), simplex, GoalType.MINIMIZE, IGuess, new SearchInterval(0, 1));

                double[] valueAlt = solutionAlt.getPoint();

                OptimizerResults.add(valueAlt);
                OptimizedLogLik[i] = CTSbetaBinomAlt.value(valueAlt);

                //determine lowest loglik.
                if (i == 0) {
                    lowestIndices = 0;
                    lowestLogLik = OptimizedLogLik[i];
                } else if (OptimizedLogLik[i] < lowestLogLik) {
                    lowestIndices = i;
                    lowestLogLik = OptimizedLogLik[i];
                }

                if (GlobalVariables.verbosity >= 100) {
                    System.out.println("\nAlt Loglik convergence of starting coordinates"
                            + Arrays.toString(IGuess.getInitialGuess()));
                    System.out.println("\tFinal parameters:     " + Arrays.toString(valueAlt));
                    System.out.println("\tLogLik:               " + Double.toString(OptimizedLogLik[i]));
                }

                i++;
            }

            //Now select the lowest LogLik.

            double[] bestParams = OptimizerResults.get(lowestIndices);
            altLogLik = OptimizedLogLik[lowestIndices];

            binomRatioCellType = bestParams[0];
            binomRatioResidual = bestParams[1];

            //chi squared statistic is determined based on both null and alt loglikelihoods.
            chiSq = LikelihoodFunctions.ChiSqFromLogLik(nullLogLik, altLogLik);

            //determine P value based on distribution

            pVal = LikelihoodFunctions.determinePvalFrom1DFchiSq(chiSq);

            if (GlobalVariables.verbosity >= 10) {
                System.out.println("\n--- Starting cell type specific beta binomial LRT test estimate ---");
                System.out.println("LogLik of converged to a threshold of "
                        + Double.toString(GlobalVariables.simplexThreshold) + "\n");
                System.out.println(
                        "\tCellType Binomial ratio:       " + Double.toString(binomRatioCellType) + "\n");
                System.out.println(
                        "\tResidual Binomial ratio:       " + Double.toString(binomRatioResidual) + "\n");
                System.out.println("\tNull log likelihood:           " + Double.toString(nullLogLik));
                System.out.println("\tAlt log likelihood:            " + Double.toString(altLogLik) + "\n");
                System.out.println("\tChisq statistic:               " + Double.toString(chiSq));
                System.out.println("\tP value:                       " + Double.toString(pVal));
                System.out.println("-----------------------------------------------------------------------");

            }

            testConverged = true;

        } catch (TooManyEvaluationsException e) {

            if (GlobalVariables.verbosity >= 1) {
                System.out.println("WARNING: Did not converge to a solution for SNP " + snpName
                        + "in cell type specific beta binomial");
                System.out.println("         After " + Integer.toString(GlobalVariables.maximumIterations)
                        + " iterations.");
                System.out.println("         Continue-ing with the next.");
            }

        }

        //Finally test was done, so we say the test was performed.
        testPerformed = true;

    }

}

From source file:org.dawnsci.plotting.tools.powdercheck.PowderCheckJob.java

private List<PowderCheckResult> fitPeaksToTrace(final Dataset xIn, final Dataset yIn, Dataset baselineIn) {

    resultList.clear();//from w w w . j a v a 2s  .  co m

    List<HKL> spacings = CalibrationFactory.getCalibrationStandards().getCalibrant().getHKLs();
    final double[] qVals = new double[spacings.size()];

    for (int i = 0; i < spacings.size(); i++) {
        if (xAxis == XAxis.ANGLE)
            qVals[i] = 2 * Math.toDegrees(Math.asin((metadata.getDiffractionCrystalEnvironment().getWavelength()
                    / (2 * spacings.get(i).getDNano() * 10))));
        else
            qVals[i] = (Math.PI * 2) / (spacings.get(i).getDNano() * 10);
    }

    double qMax = xIn.max().doubleValue();
    double qMin = xIn.min().doubleValue();

    List<Double> qList = new ArrayList<Double>();

    int count = 0;

    for (double q : qVals) {
        if (q > qMax || q < qMin)
            continue;
        count++;
        qList.add(q);
    }

    double minPeak = Collections.min(qList);
    double maxPeak = Collections.max(qList);

    int minXidx = ROISliceUtils.findPositionOfClosestValueInAxis(xIn, minPeak) - EDGE_PIXEL_NUMBER;
    int maxXidx = ROISliceUtils.findPositionOfClosestValueInAxis(xIn, maxPeak) + EDGE_PIXEL_NUMBER;

    int maxSize = xIn.getSize();

    minXidx = minXidx < 0 ? 0 : minXidx;
    maxXidx = maxXidx > maxSize - 1 ? maxSize - 1 : maxXidx;

    final Dataset x = xIn.getSlice(new int[] { minXidx }, new int[] { maxXidx }, null);
    final Dataset y = yIn.getSlice(new int[] { minXidx }, new int[] { maxXidx }, null);
    y.setName("Fit");
    Dataset baseline = baselineIn.getSlice(new int[] { minXidx }, new int[] { maxXidx }, null);

    List<APeak> peaks = Generic1DFitter.fitPeaks(x, y, Gaussian.class, count + 10);

    List<PowderCheckResult> initResults = new ArrayList<PowderCheckResult>();

    CompositeFunction cf = new CompositeFunction();

    for (APeak peak : peaks)
        cf.addFunction(peak);

    double limit = findMatchLimit(qList, cf);

    while (cf.getNoOfFunctions() != 0 && !qList.isEmpty())
        findMatches(initResults, qList, cf, limit);

    final CompositeFunction cfFinal = compositeFunctionFromResults(initResults);

    double[] initParam = new double[cfFinal.getFunctions().length * 3];

    {
        int i = 0;
        for (IFunction func : cfFinal.getFunctions()) {
            initParam[i++] = func.getParameter(0).getValue();
            initParam[i++] = func.getParameter(1).getValue();
            initParam[i++] = func.getParameter(2).getValue();
        }
    }

    final Dataset yfit = DatasetFactory.zeros(x, Dataset.FLOAT64);

    MultivariateOptimizer opt = new SimplexOptimizer(REL_TOL, ABS_TOL);

    MultivariateFunction fun = new MultivariateFunction() {

        @Override
        public double value(double[] arg0) {

            int j = 0;
            for (IFunction func : cfFinal.getFunctions()) {

                double[] p = func.getParameterValues();
                p[0] = arg0[j++];
                p[1] = arg0[j++];
                p[2] = arg0[j++];
                func.setParameterValues(p);
            }

            for (int i = 0; i < yfit.getSize(); i++) {
                yfit.set(cfFinal.val(x.getDouble(i)), i);
            }

            return y.residual(yfit);
        }
    };

    opt.optimize(new InitialGuess(initParam), GoalType.MINIMIZE, new ObjectiveFunction(fun),
            new MaxEval(MAX_EVAL), new NelderMeadSimplex(initParam.length));

    Dataset fit = Maths.add(yfit, baseline);
    fit.setName("Fit");
    Dataset residual = Maths.subtract(y, yfit);
    residual.setName("Residual");

    system.updatePlot1D(x, Arrays.asList(new IDataset[] { fit, residual }), null);
    setPlottingSystemAxes();
    for (int i = 0; i < cfFinal.getNoOfFunctions(); i++) {
        resultList.add(new PowderCheckResult(cfFinal.getFunction(i), initResults.get(i).getCalibrantQValue()));
    }

    return resultList;

}

From source file:org.hawkular.datamining.forecast.models.AbstractModelOptimizer.java

protected void optimize(double[] initialGuess, MultivariateFunctionMappingAdapter costFunction) {

    // Nelder-Mead Simplex
    SimplexOptimizer optimizer = new SimplexOptimizer(0.0001, 0.0001);
    PointValuePair unBoundedResult = optimizer.optimize(GoalType.MINIMIZE, new MaxIter(MAX_ITER),
            new MaxEval(MAX_EVAL), new InitialGuess(initialGuess), new ObjectiveFunction(costFunction),
            new NelderMeadSimplex(initialGuess.length));

    result = costFunction.unboundedToBounded(unBoundedResult.getPoint());
}