List of usage examples for org.apache.commons.math3.optim SimplePointChecker SimplePointChecker
public SimplePointChecker(final double relativeThreshold, final double absoluteThreshold)
From source file:net.sf.tweety.math.opt.solver.ApacheCommonsNonLinearConjugateGradientOptimizer.java
@Override public Map<Variable, Term> solve(ConstraintSatisfactionProblem problem) throws GeneralMathException { // only optimization problems if (!(problem instanceof OptimizationProblem)) throw new IllegalArgumentException("Only optimization problems allowed for this solver."); OptimizationProblem p = (OptimizationProblem) problem; // no constraints allowed if (!p.isEmpty()) throw new IllegalArgumentException( "Only optimization problems without constraints allowed for this solver."); final Term target = p.getTargetFunction(); final List<Variable> vars = new ArrayList<Variable>(target.getVariables()); MultivariateFunction acTarget = new MultivariateFunction() { @Override/*from w ww . j a va 2 s.c o m*/ public double value(double[] arg0) { return target.replaceAllTerms(arg0, vars).doubleValue(); } }; final Term[] targetGradient = new Term[vars.size()]; for (int i = 0; i < vars.size(); i++) targetGradient[i] = target.derive(vars.get(i)); MultivariateVectorFunction acTargetGradient = new MultivariateVectorFunction() { @Override public double[] value(double[] arg0) throws IllegalArgumentException { double[] result = new double[arg0.length]; for (int i = 0; i < arg0.length; i++) result[i] = targetGradient[i].replaceAllTerms(arg0, vars).doubleValue(); return result; } }; // create solver NonLinearConjugateGradientOptimizer optimizer = new NonLinearConjugateGradientOptimizer( NonLinearConjugateGradientOptimizer.Formula.FLETCHER_REEVES, new SimplePointChecker<PointValuePair>(this.precision, this.precision)); double[] s = new double[vars.size()]; for (int i = 0; i < vars.size(); i++) s[i] = 0.5; PointValuePair val = optimizer.optimize(new ObjectiveFunction(acTarget), new ObjectiveFunctionGradient(acTargetGradient), new InitialGuess(s), p.getType() == OptimizationProblem.MAXIMIZE ? GoalType.MAXIMIZE : GoalType.MINIMIZE, new MaxEval(this.maxEval)); Map<Variable, Term> result = new HashMap<Variable, Term>(); for (int i = 0; i < vars.size(); i++) result.put(vars.get(i), new FloatConstant(val.getPoint()[i])); return result; }
From source file:net.sf.tweety.math.opt.solver.ApacheCommonsCMAESOptimizer.java
@Override public Map<Variable, Term> solve(ConstraintSatisfactionProblem problem) throws GeneralMathException { // only optimization problems if (!(problem instanceof OptimizationProblem)) throw new IllegalArgumentException("Only optimization problems allowed for this solver."); OptimizationProblem p = (OptimizationProblem) problem; // only box constraints allowed (so far) if (!p.isEmpty()) throw new IllegalArgumentException( "Only optimization problems with box constraints on variables allowed for this solver (no other constraints."); final List<Variable> vars = new ArrayList<Variable>(p.getTargetFunction().getVariables()); double[] lb = new double[vars.size()]; double[] ub = new double[vars.size()]; double[] s = new double[vars.size()]; double[] sigma = new double[vars.size()]; for (int i = 0; i < vars.size(); i++) { lb[i] = vars.get(i).getLowerBound(); ub[i] = vars.get(i).getUpperBound(); s[i] = (lb[i] + ub[i]) / 2;/*from w ww.j a v a 2s.c om*/ sigma[i] = ub[i] - lb[i]; } final Term targetFunction = p.getTargetFunction(); MultivariateFunction target = new MultivariateFunction() { @Override public double value(double[] arg0) { return targetFunction.replaceAllTerms(arg0, vars).doubleValue(); } }; // construct solver CMAESOptimizer optimizer = new CMAESOptimizer(this.maxIterations, this.stopFitness, this.isActiveCMA, this.diagonalOnly, this.checkFeasableCount, new JDKRandomGenerator(), false, new SimplePointChecker<PointValuePair>(this.precision, this.precision)); PointValuePair val = optimizer.optimize(new CMAESOptimizer.Sigma(sigma), new ObjectiveFunction(target), new InitialGuess(s), p.getType() == OptimizationProblem.MAXIMIZE ? GoalType.MAXIMIZE : GoalType.MINIMIZE, new MaxEval(this.maxIterations), new SimpleBounds(lb, ub), new CMAESOptimizer.PopulationSize(this.populationSize)); Map<Variable, Term> result = new HashMap<Variable, Term>(); for (int i = 0; i < vars.size(); i++) result.put(vars.get(i), new FloatConstant(val.getPoint()[i])); return result; }
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 v a 2s . 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.utils.Cell_Outliner.java
/** * Find an ellipse that optimises the fit to the polygon detected edges. * //from ww w . jav a 2 s. co m * @param roi * @param params * @param weightMap * @param angle * @return */ private double[] fitPolygonalCell(PolygonRoi roi, double[] params, FloatProcessor weightMap, FloatProcessor angle) { // Get an estimate of the starting parameters using the current polygon double[] startPoint = params; startPoint = estimateStartPoint(roi, weightMap.getWidth(), weightMap.getHeight()); int maxEval = 2000; final DifferentiableEllipticalFitFunction func = new DifferentiableEllipticalFitFunction(roi, weightMap); double relativeThreshold = 100 * Precision.EPSILON; double absoluteThreshold = 100 * Precision.SAFE_MIN; ConvergenceChecker<PointVectorValuePair> checker = new SimplePointChecker<PointVectorValuePair>( relativeThreshold, absoluteThreshold); double initialStepBoundFactor = 10; double costRelativeTolerance = 1e-10; double parRelativeTolerance = 1e-10; double orthoTolerance = 1e-10; double threshold = Precision.SAFE_MIN; LevenbergMarquardtOptimizer optimiser = new LevenbergMarquardtOptimizer(initialStepBoundFactor, checker, costRelativeTolerance, parRelativeTolerance, orthoTolerance, threshold); try { PointVectorValuePair solution = optimiser.optimize(new MaxIter(maxEval), new MaxEval(Integer.MAX_VALUE), new ModelFunctionJacobian(new MultivariateMatrixFunction() { public double[][] value(double[] point) throws IllegalArgumentException { return func.jacobian(point); } }), new ModelFunction(func), new Target(func.calculateTarget()), new Weight(func.calculateWeights()), new InitialGuess(startPoint)); if (debug) IJ.log(String.format("Eval = %d (Iter = %d), RMS = %f", optimiser.getEvaluations(), optimiser.getIterations(), optimiser.getRMS())); return solution.getPointRef(); } catch (Exception e) { IJ.log("Failed to find an elliptical solution, defaulting to polygon"); e.printStackTrace(); } return null; }
From source file:org.compevol.ssgd.MaximumLikelihood.java
public MaximumLikelihood(final MultivariateFunction likelihood, final Parameter variables) { this.likelihood = likelihood; this.variables = variables; initial = new double[variables.getDimension()]; Arrays.fill(initial, 1);//ww w . ja v a 2s . c om optimizer = new CMAESOptimizer(Integer.MAX_VALUE, 0.0, true, 0, 8096, new RandomGenerator() { @Override public void setSeed(int i) { throw new UnsupportedOperationException(); } @Override public void setSeed(int[] ints) { throw new UnsupportedOperationException(); } @Override public void setSeed(long l) { throw new UnsupportedOperationException(); } @Override public void nextBytes(byte[] bytes) { MathUtils.nextBytes(bytes); } @Override public int nextInt() { return MathUtils.nextInt(); } @Override public int nextInt(int i) { return MathUtils.nextInt(i); } @Override public long nextLong() { return MathUtils.nextLong(); } @Override public boolean nextBoolean() { return MathUtils.nextBoolean(); } @Override public float nextFloat() { return MathUtils.nextFloat(); } @Override public double nextDouble() { return MathUtils.nextDouble(); } @Override public double nextGaussian() { return MathUtils.nextGaussian(); } }, true, new SimplePointChecker<PointValuePair>(MachineAccuracy.SQRT_EPSILON, MachineAccuracy.EPSILON)); }
From source file:uk.ac.diamond.scisoft.analysis.diffraction.FittingUtils.java
/** * @param opt// ww w.j av a2s . co m * @param n * @return optimizer */ public static MultivariateOptimizer createOptimizer(Optimizer opt, int n) { switch (opt) { case BOBYQA: return new BOBYQAOptimizer(n + 2); case CMAES: return new CMAESOptimizer(MAX_ITER, 0., true, 0, 10, seed == null ? new Well19937c() : new Well19937c(seed), false, new SimplePointChecker<PointValuePair>(REL_TOL, ABS_TOL)); case Simplex: default: return new SimplexOptimizer(new SimplePointChecker<PointValuePair>(REL_TOL * 1e3, ABS_TOL * 1e3)); } }
From source file:uk.ac.diamond.scisoft.analysis.optimize.ApacheOptimizer.java
private MultivariateOptimizer createOptimizer() { SimplePointChecker<PointValuePair> checker = new SimplePointChecker<PointValuePair>(REL_TOL, ABS_TOL); switch (optimizer) { case CONJUGATE_GRADIENT: return new NonLinearConjugateGradientOptimizer(Formula.POLAK_RIBIERE, checker); case BOBYQA://from w w w . j av a2s . co m return new BOBYQAOptimizer(n + 2); case CMAES: return new CMAESOptimizer(MAX_ITER, 0., true, 0, 10, seed == null ? new Well19937c() : new Well19937c(seed), false, new SimplePointChecker<PointValuePair>(REL_TOL, ABS_TOL)); case POWELL: return new PowellOptimizer(REL_TOL, ABS_TOL, checker); case SIMPLEX_MD: case SIMPLEX_NM: return new SimplexOptimizer(checker); default: throw new IllegalStateException("Should not be called"); } }
From source file:uk.ac.diamond.scisoft.ncd.calibration.rcp.views.SaxsQAxisCalibration.java
@Override protected void runJavaCommand() { try {//ww w . j a v a 2s . c o m String saxsDet = getDetectorName(); if (saxsDet == null) { throw new IllegalArgumentException(NcdMessages.NO_SAXS_DETECTOR); } Amount<Energy> energy = getEnergy(); if (energy == null) { throw new IllegalArgumentException(NLS.bind(NcdMessages.NO_ENERGY_DATA, "calibration settings")); } Amount<Length> pxSize = getPixel(); if (pxSize == null) { throw new IllegalArgumentException(NcdMessages.NO_SAXS_PIXEL); } final boolean runRefinement = beamRefineButton.getSelection(); IPlottingSystem plotSystem = PlottingFactory.getPlottingSystem(GUI_PLOT_NAME); Collection<IRegion> sectorRegions = plotSystem.getRegions(RegionType.SECTOR); if (sectorRegions == null || sectorRegions.isEmpty()) { throw new IllegalArgumentException(NcdMessages.NO_SEC_DATA); } if (sectorRegions.size() > 1) { throw new IllegalArgumentException(NcdMessages.NO_SEC_SUPPORT); } final SectorROI sroi = (SectorROI) sectorRegions.iterator().next().getROI(); if (runRefinement) { IImageTrace trace = (IImageTrace) plotSystem.getTraces().iterator().next(); final Dataset dataset = (Dataset) trace.getData(); final Dataset mask = (Dataset) trace.getMask(); final MultivariateFunctionWithMonitor beamOffset = new MultivariateFunctionWithMonitor(dataset, mask, sroi); beamOffset.addSourceProviders(service); int cmaesLambda = uk.ac.diamond.scisoft.ncd.core.rcp.Activator.getDefault().getPreferenceStore() .getInt(NcdPreferences.CMAESlambda); double cmaesInputSigmaPref = uk.ac.diamond.scisoft.ncd.core.rcp.Activator.getDefault() .getPreferenceStore().getInt(NcdPreferences.CMAESsigma); double[] cmaesInputSigma = new double[] { cmaesInputSigmaPref, cmaesInputSigmaPref }; int cmaesMaxIterations = uk.ac.diamond.scisoft.ncd.core.rcp.Activator.getDefault() .getPreferenceStore().getInt(NcdPreferences.CMAESmaxiteration); int cmaesCheckerPref = uk.ac.diamond.scisoft.ncd.core.rcp.Activator.getDefault() .getPreferenceStore().getInt(NcdPreferences.CMAESchecker); SimplePointChecker<PointValuePair> cmaesChecker = new SimplePointChecker<PointValuePair>(1e-4, 1.0 / cmaesCheckerPref); beamOffset.configureOptimizer(cmaesLambda == 0 ? null : cmaesLambda, cmaesInputSigmaPref == 0 ? null : cmaesInputSigma, cmaesMaxIterations == 0 ? null : cmaesMaxIterations, null, cmaesCheckerPref == 0 ? null : cmaesChecker); beamOffset.setInitPeaks(peaks); beamOffset.optimize(sroi.getPoint()); } else { peaksSourceProvider.putPeaks(peaks); } } catch (Exception e) { Status status = new Status(IStatus.ERROR, Activator.PLUGIN_ID, e.getMessage()); StatusManager.getManager().handle(status, StatusManager.SHOW); return; } }
From source file:uk.ac.diamond.scisoft.ncd.core.data.plots.GuinierPlotData.java
public Object[] getGuinierPlotParameters(IDataset data, IDataset axis) { Dataset guinierData = getSaxsPlotDataset(data, axis); Dataset guinierAxis = getSaxsPlotAxis(axis); /* int sliceSize = data.getSize() / 10; int[] step = new int[] {sliceSize}; IndexIterator dataSliceIter = guinierData.getSliceIterator(null, null, null); IndexIterator axisSliceIter = guinierAxis.getSliceIterator(null, null, null); // Iterate over data slices//from w ww. ja va 2 s.com int[] minPosition = new int[] {-1}; double minError = Double.MAX_VALUE; double slope = Double.NaN; double intercept = Double.NaN; Map<DoublePoint, Double> clusterInputMap = new HashMap<DoublePoint, Double>(); while (dataSliceIter.hasNext() && axisSliceIter.hasNext()) { SimpleRegression regression = new SimpleRegression(); int[] pos = dataSliceIter.getPos(); // Iterate over data values for linear regression IndexIterator dataIter = new SliceIterator( guinierData.getShape(), guinierData.getSize(), pos, step); pos = axisSliceIter.getPos(); IndexIterator axisIter = new SliceIterator( guinierAxis.getShape(), guinierAxis.getSize(), pos, step); int points = 0; while (dataIter.hasNext() && axisIter.hasNext()) { double dataVal = guinierData.getDouble(dataIter.getPos()); double axisVal = guinierAxis.getDouble(axisIter.getPos()); regression.addData(axisVal, dataVal); points++; } if (points == sliceSize) { regression.regress(); double err = regression.getMeanSquareError(); if (err < minError) { minError = err; minPosition = Arrays.copyOf(pos, pos.length); slope = regression.getSlope(); intercept = regression.getIntercept(); double I0 = Math.exp(intercept); double Rg = Math.sqrt(-3.0*slope); System.out.println(" Min Pos : " + Arrays.toString(minPosition)); System.out.println(" Min Error : " + Double.toString(minError)); System.out.println(" Slope : " + Double.toString(slope) + " Intercept : " + Double.toString(intercept)); System.out.println(" I(0) : " + Double.toString(I0) + " Rg : " + Double.toString(Rg)); clusterInputMap.put(new DoublePoint(minPosition), minError); } } else { break; } } DBSCANClusterer<DoublePoint> clusterer = new DBSCANClusterer<DoublePoint>(5, 5); List<Cluster<DoublePoint>> clusterResults = clusterer.cluster(clusterInputMap.keySet()); // output the clusters for (int i = 0; i < clusterResults.size(); i++) { System.out.println("Cluster " + i); double[] minPoint = null; double minVal = Double.MAX_VALUE; for (DoublePoint point : clusterResults.get(i).getPoints()) { System.out.println(Arrays.toString(point.getPoint())); Double val = clusterInputMap.get(point); if (val < minVal) { minVal = val; minPoint = Arrays.copyOf(point.getPoint(), point.getPoint().length); } minVal = (val < minVal ? val : minVal); } System.out.println("Min cluster point : " + Arrays.toString(minPoint)); System.out.println("Min cluster value : " + Double.toString(minVal)); System.out.println(); } */ ConvergenceChecker<PointValuePair> cmaesChecker = new SimplePointChecker<PointValuePair>(1e-6, 1e-8); RandomDataGenerator rnd = new RandomDataGenerator(); CMAESOptimizer optimizer = new CMAESOptimizer(cmaesMaxIterations, 0.0, true, 0, cmaesCheckFeasableCount, rnd.getRandomGenerator(), false, cmaesChecker); GuinierLineFitFunction function = new GuinierLineFitFunction(guinierData, guinierAxis); Amount<Dimensionless> I0 = Amount.valueOf(Double.NaN, Double.NaN, Dimensionless.UNIT); Amount<Dimensionless> Rg = Amount.valueOf(Double.NaN, Double.NaN, Dimensionless.UNIT); double[] qvals = new double[] { Double.NaN, Double.NaN }; double q0 = guinierAxis.getDouble(0); double qMin = guinierAxis.getDouble(1); double qMax = guinierAxis.getDouble(guinierAxis.getSize() - 1); double[] startPosition = new double[] { guinierAxis.getDouble(0), guinierAxis.getDouble(GuinierLineFitFunction.MIN_POINTS) }; double[] cmaesInputSigma = new double[] { (qMin - q0) * 0.1, qMax * 0.1 }; try { final PointValuePair res = optimizer.optimize(new MaxEval(cmaesMaxIterations), new ObjectiveFunction(function), GoalType.MAXIMIZE, new CMAESOptimizer.PopulationSize(cmaesLambda), new CMAESOptimizer.Sigma(cmaesInputSigma), new SimpleBounds(new double[] { q0, q0 }, new double[] { qMin, qMax }), new InitialGuess(startPosition)); qvals = res.getPoint(); function.value(qvals); I0 = getI0(function.regression); Rg = getRg(function.regression); System.out.println("Final Result"); String msg = StringUtils.join(new String[] { " I(0) ", I0.toString(), " Rg ", Rg.toString() }, " : "); System.out.println(msg); msg = StringUtils.join(new String[] { "Slice", ArrayUtils.toString(res.getPoint()), "R", Double.toString(function.regression.getR()) }, " : "); System.out.println(msg); /* // Run Monte-Carlo simulation to generate error estimates Rg values //double finalR = function.regression.getR(); int maxSample = 10000; int minSample = 10; int totalSample = 100000; int counter = 0; int totalCounter = 0; GuinierLineFitFunction mcFunction = new GuinierLineFitFunction(guinierData, guinierAxis); DescriptiveStatistics statsR = new DescriptiveStatistics(); List<Pair<Double, Amount<Dimensionless>>> listI0 = new ArrayList<Pair<Double,Amount<Dimensionless>>>(); List<Pair<Double, Amount<Dimensionless>>> listRg = new ArrayList<Pair<Double,Amount<Dimensionless>>>(); while ((counter < maxSample && totalCounter < totalSample) || (counter < minSample && totalCounter >= totalSample)) { double q1 = rnd.nextUniform(q0, qMin); double q2 = rnd.nextUniform(q0, qMax); if (!(q2 > q1)) { continue; } totalCounter++; mcFunction.value(new double[] {q1, q2}); double tmpR = Math.abs(mcFunction.regression.getR()); //boolean equalsR = Precision.equalsWithRelativeTolerance(tmpR, finalR, 0.1); if (!(Double.isNaN(tmpR) || Double.isInfinite(tmpR))) { statsR.addValue(tmpR); Amount<Dimensionless> tmpI0 = getI0(mcFunction.regression); Amount<Dimensionless> tmpRg = getRg(mcFunction.regression); if (Double.isNaN(tmpI0.getEstimatedValue()) || Double.isInfinite(tmpI0.getEstimatedValue()) || Double.isNaN(tmpRg.getEstimatedValue()) || Double.isInfinite(tmpRg.getEstimatedValue())) { continue; } listI0.add(new Pair<Double, Amount<Dimensionless>>(tmpR, tmpI0)); listRg.add(new Pair<Double, Amount<Dimensionless>>(tmpR, tmpRg)); counter++; } } double threshold = statsR.getPercentile(90); //double threshold = 0.95*statsR.getMax(); SummaryStatistics statsI0 = new SummaryStatistics(); SummaryStatistics statsRg = new SummaryStatistics(); for (Pair<Double, Amount<Dimensionless>> tmpVal : listRg) { if (tmpVal.getFirst() > threshold) { statsRg.addValue(tmpVal.getSecond().getEstimatedValue()); } } for (Pair<Double, Amount<Dimensionless>> tmpVal : listI0) { if (tmpVal.getFirst() > threshold) { statsI0.addValue(tmpVal.getSecond().getEstimatedValue()); } } double meanI0 = statsI0.getMean(); double stdI0 = statsI0.getStandardDeviation(); I0 = Amount.valueOf(meanI0, stdI0, Dimensionless.UNIT); double meanRg = statsRg.getMean(); double stdRg = statsRg.getStandardDeviation(); Rg = Amount.valueOf(meanRg, stdRg, Dimensionless.UNIT); String msg = StringUtils.join(new String[] { "Monte-Carlo Rg", Rg.toString() }, " : "); System.out.println(msg); */ } catch (MaxCountExceededException e) { System.out.println("Maximum counts exceeded"); return null; } return new Object[] { I0, Rg, qvals[0], qvals[1] }; }