Search in sources :

Example 21 with Well19937c

use of org.apache.commons.math3.random.Well19937c in project GDSC-SMLM by aherbert.

the class FastMLEGradient2ProcedureTest method gradientProcedureComputesSameLogLikelihoodAsMLEGradientCalculator.

private void gradientProcedureComputesSameLogLikelihoodAsMLEGradientCalculator(int nparams) {
    int iter = 10;
    rdg = new RandomDataGenerator(new Well19937c(30051977));
    ArrayList<double[]> paramsList = new ArrayList<double[]>(iter);
    ArrayList<double[]> yList = new ArrayList<double[]>(iter);
    createFakeData(nparams, iter, paramsList, yList);
    FakeGradientFunction func = new FakeGradientFunction(blockWidth, nparams);
    MLEGradientCalculator calc = (MLEGradientCalculator) GradientCalculatorFactory.newCalculator(nparams, true);
    String name = String.format("[%d]", nparams);
    for (int i = 0; i < paramsList.size(); i++) {
        FastMLEGradient2Procedure p = FastMLEGradient2ProcedureFactory.createUnrolled(yList.get(i), func);
        double s = p.computeLogLikelihood(paramsList.get(i));
        double s2 = calc.logLikelihood(yList.get(i), paramsList.get(i), func);
        // Virtually the same ...
        Assert.assertEquals(name + " Result: Not same @ " + i, s, s2, Math.abs(s) * 1e-5);
    }
}
Also used : RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) ArrayList(java.util.ArrayList) Well19937c(org.apache.commons.math3.random.Well19937c) FakeGradientFunction(gdsc.smlm.function.FakeGradientFunction)

Example 22 with Well19937c

use of org.apache.commons.math3.random.Well19937c in project GDSC-SMLM by aherbert.

the class FastMLEGradient2ProcedureTest method gradientProcedureUnrolledComputesSameAsGradientProcedure.

private void gradientProcedureUnrolledComputesSameAsGradientProcedure(int nparams) {
    int iter = 10;
    rdg = new RandomDataGenerator(new Well19937c(30051977));
    ArrayList<double[]> paramsList = new ArrayList<double[]>(iter);
    ArrayList<double[]> yList = new ArrayList<double[]>(iter);
    createFakeData(nparams, iter, paramsList, yList);
    FakeGradientFunction func = new FakeGradientFunction(blockWidth, nparams);
    FastMLEGradient2Procedure p1, p2;
    String name = String.format("[%d]", nparams);
    for (int i = 0; i < paramsList.size(); i++) {
        p1 = new FastMLEGradient2Procedure(yList.get(i), func);
        p2 = FastMLEGradient2ProcedureFactory.createUnrolled(yList.get(i), func);
        double[] a = paramsList.get(i);
        double ll1 = p1.computeLogLikelihood(a);
        double ll2 = p2.computeLogLikelihood(a);
        Assert.assertEquals(name + " LL: Not same @ " + i, ll1, ll2, 0);
        p1 = new FastMLEGradient2Procedure(yList.get(i), func);
        p2 = FastMLEGradient2ProcedureFactory.createUnrolled(yList.get(i), func);
        p1.computeFirstDerivative(a);
        p2.computeFirstDerivative(a);
        Assert.assertArrayEquals(name + " first derivative value: Not same @ " + i, p1.u, p2.u, 0);
        Assert.assertArrayEquals(name + " first derivative: Not same @ " + i, p1.d1, p2.d1, 0);
        p1 = new FastMLEGradient2Procedure(yList.get(i), func);
        p2 = FastMLEGradient2ProcedureFactory.createUnrolled(yList.get(i), func);
        p1.computeSecondDerivative(a);
        p2.computeSecondDerivative(a);
        Assert.assertArrayEquals(name + " update value: Not same @ " + i, p1.u, p2.u, 0);
        Assert.assertArrayEquals(name + " update: Not same d1 @ " + i, p1.d1, p2.d1, 0);
        Assert.assertArrayEquals(name + " update: Not same d2 @ " + i, p1.d2, p2.d2, 0);
    }
}
Also used : RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) ArrayList(java.util.ArrayList) Well19937c(org.apache.commons.math3.random.Well19937c) FakeGradientFunction(gdsc.smlm.function.FakeGradientFunction)

Example 23 with Well19937c

use of org.apache.commons.math3.random.Well19937c in project GDSC-SMLM by aherbert.

the class LSQLVMGradientProcedureTest method gradientProcedureComputesSameAsGradientCalculator.

private void gradientProcedureComputesSameAsGradientCalculator(int nparams, BaseLSQLVMGradientProcedureFactory factory) {
    int iter = 10;
    rdg = new RandomDataGenerator(new Well19937c(30051977));
    double[][] alpha = new double[nparams][nparams];
    double[] beta = new double[nparams];
    ArrayList<double[]> paramsList = new ArrayList<double[]>(iter);
    ArrayList<double[]> yList = new ArrayList<double[]>(iter);
    int[] x = createFakeData(nparams, iter, paramsList, yList);
    int n = x.length;
    FakeGradientFunction func = new FakeGradientFunction(blockWidth, nparams);
    GradientCalculator calc = GradientCalculatorFactory.newCalculator(nparams, false);
    String name = factory.getClass().getSimpleName();
    for (int i = 0; i < paramsList.size(); i++) {
        BaseLSQLVMGradientProcedure p = factory.createProcedure(yList.get(i), func);
        p.gradient(paramsList.get(i));
        double s = p.value;
        double s2 = calc.findLinearised(n, yList.get(i), paramsList.get(i), alpha, beta, func);
        // Exactly the same ...
        Assert.assertEquals(name + " Result: Not same @ " + i, s, s2, 0);
        Assert.assertArrayEquals(name + " Observations: Not same beta @ " + i, p.beta, beta, 0);
        double[] al = p.getAlphaLinear();
        Assert.assertArrayEquals(name + " Observations: Not same alpha @ " + i, al, new DenseMatrix64F(alpha).data, 0);
        double[][] am = p.getAlphaMatrix();
        for (int j = 0; j < nparams; j++) Assert.assertArrayEquals(name + " Observations: Not same alpha @ " + i, am[j], alpha[j], 0);
    }
}
Also used : RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) ArrayList(java.util.ArrayList) Well19937c(org.apache.commons.math3.random.Well19937c) DenseMatrix64F(org.ejml.data.DenseMatrix64F) FakeGradientFunction(gdsc.smlm.function.FakeGradientFunction)

Example 24 with Well19937c

use of org.apache.commons.math3.random.Well19937c in project GDSC-SMLM by aherbert.

the class LSQLVMGradientProcedureTest method gradientProcedureComputesGradient.

private void gradientProcedureComputesGradient(ErfGaussian2DFunction func) {
    int nparams = func.getNumberOfGradients();
    int[] indices = func.gradientIndices();
    int iter = 100;
    rdg = new RandomDataGenerator(new Well19937c(30051977));
    ArrayList<double[]> paramsList = new ArrayList<double[]>(iter);
    ArrayList<double[]> yList = new ArrayList<double[]>(iter);
    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();
        BaseLSQLVMGradientProcedure p = LSQLVMGradientProcedureFactory.create(y, func);
        p.gradient(a);
        //double s = p.ssx;
        double[] beta = p.beta.clone();
        for (int j = 0; j < nparams; j++) {
            int k = indices[j];
            double d = Precision.representableDelta(a[k], (a[k] == 0) ? 1e-3 : a[k] * delta);
            a2[k] = a[k] + d;
            p.value(a2);
            double s1 = p.value;
            a2[k] = a[k] - d;
            p.value(a2);
            double s2 = p.value;
            a2[k] = a[k];
            // Apply a factor of -2 to compute the actual gradients:
            // See Numerical Recipes in C++, 2nd Ed. Equation 15.5.6 for Nonlinear Models 
            beta[j] *= -2;
            double gradient = (s1 - s2) / (2 * d);
            //System.out.printf("[%d,%d] %f  (%s %f+/-%f)  %f  ?=  %f\n", i, k, s, func.getName(k), a[k], d, beta[j],
            //		gradient);
            Assert.assertTrue("Not same gradient @ " + j, eq.almostEqualRelativeOrAbsolute(beta[j], gradient));
        }
    }
}
Also used : RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) ArrayList(java.util.ArrayList) DoubleEquality(gdsc.core.utils.DoubleEquality) Well19937c(org.apache.commons.math3.random.Well19937c)

Example 25 with Well19937c

use of org.apache.commons.math3.random.Well19937c in project GDSC-SMLM by aherbert.

the class LVMGradientProcedureTest method gradientProcedureComputesSameAsGradientCalculator.

private void gradientProcedureComputesSameAsGradientCalculator(int nparams, boolean mle) {
    int iter = 10;
    rdg = new RandomDataGenerator(new Well19937c(30051977));
    double[][] alpha = new double[nparams][nparams];
    double[] beta = new double[nparams];
    ArrayList<double[]> paramsList = new ArrayList<double[]>(iter);
    ArrayList<double[]> yList = new ArrayList<double[]>(iter);
    int[] x = createFakeData(nparams, iter, paramsList, yList);
    int n = x.length;
    FakeGradientFunction func = new FakeGradientFunction(blockWidth, nparams);
    GradientCalculator calc = GradientCalculatorFactory.newCalculator(nparams, mle);
    String name = String.format("[%d] %b", nparams, mle);
    for (int i = 0; i < paramsList.size(); i++) {
        LVMGradientProcedure p = LVMGradientProcedureFactory.create(yList.get(i), func, (mle) ? Type.MLE : Type.LSQ);
        p.gradient(paramsList.get(i));
        double s = p.value;
        double s2 = calc.findLinearised(n, yList.get(i), paramsList.get(i), alpha, beta, func);
        // Exactly the same ...
        Assert.assertEquals(name + " Result: Not same @ " + i, s, s2, 0);
        Assert.assertArrayEquals(name + " Observations: Not same beta @ " + i, p.beta, beta, 0);
        double[] al = p.getAlphaLinear();
        Assert.assertArrayEquals(name + " Observations: Not same alpha @ " + i, al, new DenseMatrix64F(alpha).data, 0);
        double[][] am = p.getAlphaMatrix();
        for (int j = 0; j < nparams; j++) Assert.assertArrayEquals(name + " Observations: Not same alpha @ " + i, am[j], alpha[j], 0);
    }
}
Also used : RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) ArrayList(java.util.ArrayList) Well19937c(org.apache.commons.math3.random.Well19937c) DenseMatrix64F(org.ejml.data.DenseMatrix64F) FakeGradientFunction(gdsc.smlm.function.FakeGradientFunction)

Aggregations

Well19937c (org.apache.commons.math3.random.Well19937c)66 RandomDataGenerator (org.apache.commons.math3.random.RandomDataGenerator)41 ArrayList (java.util.ArrayList)31 RandomGenerator (org.apache.commons.math3.random.RandomGenerator)26 FakeGradientFunction (gdsc.smlm.function.FakeGradientFunction)17 Test (org.junit.Test)17 DoubleEquality (gdsc.core.utils.DoubleEquality)7 Gaussian2DFunction (gdsc.smlm.function.gaussian.Gaussian2DFunction)7 DenseMatrix64F (org.ejml.data.DenseMatrix64F)6 StoredDataStatistics (gdsc.core.utils.StoredDataStatistics)5 PointValuePair (org.apache.commons.math3.optim.PointValuePair)5 TimingService (gdsc.core.test.TimingService)4 Statistics (gdsc.core.utils.Statistics)4 Gradient1Function (gdsc.smlm.function.Gradient1Function)4 ValueProcedure (gdsc.smlm.function.ValueProcedure)4 ErfGaussian2DFunction (gdsc.smlm.function.gaussian.erf.ErfGaussian2DFunction)4 TooManyEvaluationsException (org.apache.commons.math3.exception.TooManyEvaluationsException)4 SimpleValueChecker (org.apache.commons.math3.optim.SimpleValueChecker)4 PseudoRandomGenerator (gdsc.core.utils.PseudoRandomGenerator)3 TurboList (gdsc.core.utils.TurboList)3