Search in sources :

Example 1 with DoubleEquality

use of gdsc.core.utils.DoubleEquality in project GDSC-SMLM by aherbert.

the class EJMLLinearSolverTest method canSolveLinearEquationWithZerosInA.

@Test
public void canSolveLinearEquationWithZerosInA() {
    EJMLLinearSolver solver = new EJMLLinearSolver();
    DoubleEquality eq = new DoubleEquality(5e-3, 1e-16);
    solver.setEqual(eq);
    // Solves (one) linear equation, a x = b, for x[n]
    double[][] a = new double[][] { new double[] { 2, 0, -1, 0, 0, 0 }, new double[] { 0, 0, 0, 0, 0, 0 }, new double[] { -1, 0, 2, 0, 0, -1 }, new double[] { 0, 0, 0, 0, 0, 0 }, new double[] { 0, 0, 0, 0, 0, 0 }, new double[] { 0, 0, -1, 0, 0, 2 } };
    double[] b = new double[] { 3, 0, 3, 0, 0, 4 };
    // Expected solution
    double[] x = new double[] { 4.75, 0, 6.5, 0, 0, 5.25 };
    double[][] a_inv = new double[][] { new double[] { 0.75, 0, 0.5, 0, 0, 0.25 }, new double[] { 0, 0, 0, 0, 0, 0 }, new double[] { 0.5, 0, 1, 0, 0, 0.5 }, new double[] { 0, 0, 0, 0, 0, 0 }, new double[] { 0, 0, 0, 0, 0, 0 }, new double[] { 0.25, 0, 0.5, 0, 0, 0.75 } };
    boolean result = solver.solve(a, b);
    solver.invertLastA(a);
    Assert.assertTrue("Failed to solve", result);
    Assert.assertArrayEquals("Bad solution", x, b, 1e-4f);
    log("x = %s\n", Arrays.toString(b));
    for (int i = 0; i < b.length; i++) {
        log("a[%d] = %s\n", i, Arrays.toString(a[i]));
        Assert.assertArrayEquals("Bad inversion", a_inv[i], a[i], 1e-4f);
    }
}
Also used : DoubleEquality(gdsc.core.utils.DoubleEquality) Test(org.junit.Test)

Example 2 with DoubleEquality

use of gdsc.core.utils.DoubleEquality in project GDSC-SMLM by aherbert.

the class NonLinearFit method init.

private void init(StoppingCriteria sc, double maxRelativeError, double maxAbsoluteError) {
    setStoppingCriteria(sc);
    solver.setEqual(new DoubleEquality(maxRelativeError, maxAbsoluteError));
}
Also used : DoubleEquality(gdsc.core.utils.DoubleEquality)

Example 3 with DoubleEquality

use of gdsc.core.utils.DoubleEquality in project GDSC-SMLM by aherbert.

the class GradientCalculatorSpeedTest method gradientCalculatorComputesGradient.

private void gradientCalculatorComputesGradient(GradientCalculator calc) {
    int nparams = calc.nparams;
    Gaussian2DFunction func = new SingleEllipticalGaussian2DFunction(blockWidth, blockWidth);
    // Check the function is the correct size
    Assert.assertEquals(nparams, func.gradientIndices().length);
    int iter = 100;
    rdg = new RandomDataGenerator(new Well19937c(30051977));
    double[] beta = new double[nparams];
    double[] beta2 = new double[nparams];
    ArrayList<double[]> paramsList = new ArrayList<double[]>(iter);
    ArrayList<double[]> yList = new ArrayList<double[]>(iter);
    int[] x = createData(1, iter, paramsList, yList, true);
    double delta = 1e-3;
    DoubleEquality eq = new DoubleEquality(1e-3, 1e-3);
    for (int i = 0; i < paramsList.size(); i++) {
        double[] y = yList.get(i);
        double[] a = paramsList.get(i);
        double[] a2 = a.clone();
        //double s = 
        calc.evaluate(x, y, a, beta, func);
        for (int j = 0; j < nparams; j++) {
            double d = Precision.representableDelta(a[j], (a[j] == 0) ? 1e-3 : a[j] * delta);
            a2[j] = a[j] + d;
            double s1 = calc.evaluate(x, y, a2, beta2, func);
            a2[j] = a[j] - d;
            double s2 = calc.evaluate(x, y, a2, beta2, func);
            a2[j] = a[j];
            double gradient = (s1 - s2) / (2 * d);
            //System.out.printf("[%d,%d] %f  (%s %f+/-%f)  %f  ?=  %f\n", i, j, s, func.getName(j), a[j], d, beta[j],
            //		gradient);
            Assert.assertTrue("Not same gradient @ " + j, eq.almostEqualRelativeOrAbsolute(beta[j], gradient));
        }
    }
}
Also used : SingleEllipticalGaussian2DFunction(gdsc.smlm.function.gaussian.SingleEllipticalGaussian2DFunction) EllipticalGaussian2DFunction(gdsc.smlm.function.gaussian.EllipticalGaussian2DFunction) Gaussian2DFunction(gdsc.smlm.function.gaussian.Gaussian2DFunction) SingleFreeCircularGaussian2DFunction(gdsc.smlm.function.gaussian.SingleFreeCircularGaussian2DFunction) SingleFixedGaussian2DFunction(gdsc.smlm.function.gaussian.SingleFixedGaussian2DFunction) SingleNBFixedGaussian2DFunction(gdsc.smlm.function.gaussian.SingleNBFixedGaussian2DFunction) SingleCircularGaussian2DFunction(gdsc.smlm.function.gaussian.SingleCircularGaussian2DFunction) RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) SingleEllipticalGaussian2DFunction(gdsc.smlm.function.gaussian.SingleEllipticalGaussian2DFunction) ArrayList(java.util.ArrayList) DoubleEquality(gdsc.core.utils.DoubleEquality) Well19937c(org.apache.commons.math3.random.Well19937c)

Example 4 with DoubleEquality

use of gdsc.core.utils.DoubleEquality in project GDSC-SMLM by aherbert.

the class ErfTest method analyticErfGradientCorrectForErfApproximation.

@Test
public void analyticErfGradientCorrectForErfApproximation() {
    BaseErf erf = new Erf();
    int range = 7;
    int steps = 10000;
    double step = (double) range / steps;
    double delta = 1e-3;
    DoubleEquality eq = new DoubleEquality(5e-4, 1e-6);
    for (int i = 0; i < steps; i++) {
        double x = i * step;
        double x1 = x + Precision.representableDelta(x, delta);
        double x2 = x - Precision.representableDelta(x, delta);
        double o1 = erf.erf(x1);
        double o2 = erf.erf(x2);
        double delta2 = x1 - x2;
        double g = (o1 - o2) / delta2;
        double e = gdsc.smlm.function.Erf.dErf_dx(x);
        if (!eq.almostEqualRelativeOrAbsolute(e, g))
            Assert.assertTrue(x + " : " + e + " != " + g, false);
    }
}
Also used : DoubleEquality(gdsc.core.utils.DoubleEquality) Test(org.junit.Test)

Example 5 with DoubleEquality

use of gdsc.core.utils.DoubleEquality in project GDSC-SMLM by aherbert.

the class PoissonCalculatorTest method instanceAndFastMethodIsApproximatelyEqualToStaticMethod.

@Test
public void instanceAndFastMethodIsApproximatelyEqualToStaticMethod() {
    DoubleEquality eq = new DoubleEquality(3e-4, 0);
    RandomGenerator rg = new Well19937c(30051977);
    // Test for different x. The calculator approximation begins
    int n = 100;
    double[] u = new double[n];
    double[] x = new double[n];
    double e, o;
    for (double testx : new double[] { Math.nextAfter(PoissonCalculator.APPROXIMATION_X, -1), PoissonCalculator.APPROXIMATION_X, Math.nextUp(PoissonCalculator.APPROXIMATION_X), PoissonCalculator.APPROXIMATION_X * 1.1, PoissonCalculator.APPROXIMATION_X * 2, PoissonCalculator.APPROXIMATION_X * 10 }) {
        String X = Double.toString(testx);
        Arrays.fill(x, testx);
        PoissonCalculator pc = new PoissonCalculator(x);
        e = PoissonCalculator.maximumLogLikelihood(x);
        o = pc.getMaximumLogLikelihood();
        System.out.printf("[%s] Instance MaxLL = %g vs %g (error = %g)\n", X, e, o, DoubleEquality.relativeError(e, o));
        Assert.assertTrue("Instance Max LL not equal", eq.almostEqualRelativeOrAbsolute(e, o));
        o = PoissonCalculator.fastMaximumLogLikelihood(x);
        System.out.printf("[%s] Fast MaxLL = %g vs %g (error = %g)\n", X, e, o, DoubleEquality.relativeError(e, o));
        Assert.assertTrue("Fast Max LL not equal", eq.almostEqualRelativeOrAbsolute(e, o));
        // Generate data around the value
        for (int i = 0; i < n; i++) u[i] = x[i] + rg.nextDouble() - 0.5;
        e = PoissonCalculator.logLikelihood(u, x);
        o = pc.logLikelihood(u);
        System.out.printf("[%s] Instance LL = %g vs %g (error = %g)\n", X, e, o, DoubleEquality.relativeError(e, o));
        Assert.assertTrue("Instance LL not equal", eq.almostEqualRelativeOrAbsolute(e, o));
        o = PoissonCalculator.fastLogLikelihood(u, x);
        System.out.printf("[%s] Fast LL = %g vs %g (error = %g)\n", X, e, o, DoubleEquality.relativeError(e, o));
        Assert.assertTrue("Fast LL not equal", eq.almostEqualRelativeOrAbsolute(e, o));
        e = PoissonCalculator.logLikelihoodRatio(u, x);
        o = pc.getLogLikelihoodRatio(o);
        System.out.printf("[%s] Instance LLR = %g vs %g (error = %g)\n", X, e, o, DoubleEquality.relativeError(e, o));
        Assert.assertTrue("Instance LLR not equal", eq.almostEqualRelativeOrAbsolute(e, o));
    }
}
Also used : DoubleEquality(gdsc.core.utils.DoubleEquality) Well19937c(org.apache.commons.math3.random.Well19937c) RandomGenerator(org.apache.commons.math3.random.RandomGenerator) Test(org.junit.Test)

Aggregations

DoubleEquality (gdsc.core.utils.DoubleEquality)15 Well19937c (org.apache.commons.math3.random.Well19937c)7 Test (org.junit.Test)7 ArrayList (java.util.ArrayList)6 RandomDataGenerator (org.apache.commons.math3.random.RandomDataGenerator)6 Gaussian2DFunction (gdsc.smlm.function.gaussian.Gaussian2DFunction)3 EllipticalGaussian2DFunction (gdsc.smlm.function.gaussian.EllipticalGaussian2DFunction)2 ValueProcedure (gdsc.smlm.function.ValueProcedure)1 SingleCircularGaussian2DFunction (gdsc.smlm.function.gaussian.SingleCircularGaussian2DFunction)1 SingleEllipticalGaussian2DFunction (gdsc.smlm.function.gaussian.SingleEllipticalGaussian2DFunction)1 SingleFixedGaussian2DFunction (gdsc.smlm.function.gaussian.SingleFixedGaussian2DFunction)1 SingleFreeCircularGaussian2DFunction (gdsc.smlm.function.gaussian.SingleFreeCircularGaussian2DFunction)1 SingleNBFixedGaussian2DFunction (gdsc.smlm.function.gaussian.SingleNBFixedGaussian2DFunction)1 ErfGaussian2DFunction (gdsc.smlm.function.gaussian.erf.ErfGaussian2DFunction)1 SingleFreeCircularErfGaussian2DFunction (gdsc.smlm.function.gaussian.erf.SingleFreeCircularErfGaussian2DFunction)1 RandomGenerator (org.apache.commons.math3.random.RandomGenerator)1 DenseMatrix64F (org.ejml.data.DenseMatrix64F)1