Search in sources :

Example 6 with DoubleDoubleBiPredicate

use of uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate in project GDSC-SMLM by aherbert.

the class GaussianKernelTest method canComputeDownscaleGaussianKernelIncRangeIncScale.

@Test
void canComputeDownscaleGaussianKernelIncRangeIncScale() {
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    for (int i = 0; i < 5; i++) {
        final double s = 0.33 * (1 << i);
        final GaussianKernel k = new GaussianKernel(s);
        for (int range = 3; range < 5; range++) {
            for (int scale = 1; scale <= 5; scale++) {
                for (final boolean edgeCorrection : new boolean[] { true, false }) {
                    final double[] e = GaussianKernel.makeGaussianKernel(s / scale, range, edgeCorrection);
                    final double[] o = k.getDownscaleGaussianKernel(scale, range, edgeCorrection);
                    if (MathUtils.isPow2(scale)) {
                        Assertions.assertArrayEquals(e, o);
                    } else {
                        TestAssertions.assertArrayTest(e, o, predicate);
                    }
                }
            }
        }
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) Test(org.junit.jupiter.api.Test)

Example 7 with DoubleDoubleBiPredicate

use of uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate in project GDSC-SMLM by aherbert.

the class GaussianKernelTest method canComputeDownscaleGaussianKernelIncRangeDecScale.

@Test
void canComputeDownscaleGaussianKernelIncRangeDecScale() {
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    for (int i = 0; i < 5; i++) {
        final double s = 0.33 * (1 << i);
        final GaussianKernel k = new GaussianKernel(s);
        for (int range = 3; range < 5; range++) {
            for (int scale = 5; scale >= 1; scale--) {
                for (final boolean edgeCorrection : new boolean[] { true, false }) {
                    final double[] e = GaussianKernel.makeGaussianKernel(s / scale, range, edgeCorrection);
                    final double[] o = k.getDownscaleGaussianKernel(scale, range, edgeCorrection);
                    if (MathUtils.isPow2(scale)) {
                        Assertions.assertArrayEquals(e, o);
                    } else {
                        TestAssertions.assertArrayTest(e, o, predicate);
                    }
                }
            }
        }
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) Test(org.junit.jupiter.api.Test)

Example 8 with DoubleDoubleBiPredicate

use of uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate in project GDSC-SMLM by aherbert.

the class GaussianKernelTest method canComputeDownscaleGaussianKernelDecScaleIncRange.

@Test
void canComputeDownscaleGaussianKernelDecScaleIncRange() {
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    for (int i = 0; i < 5; i++) {
        final double s = 0.33 * (1 << i);
        final GaussianKernel k = new GaussianKernel(s);
        for (int scale = 5; scale >= 1; scale--) {
            for (int range = 3; range < 5; range++) {
                for (final boolean edgeCorrection : new boolean[] { true, false }) {
                    final double[] e = GaussianKernel.makeGaussianKernel(s / scale, range, edgeCorrection);
                    final double[] o = k.getDownscaleGaussianKernel(scale, range, edgeCorrection);
                    if (MathUtils.isPow2(scale)) {
                        Assertions.assertArrayEquals(e, o);
                    } else {
                        TestAssertions.assertArrayTest(e, o, predicate);
                    }
                }
            }
        }
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) Test(org.junit.jupiter.api.Test)

Example 9 with DoubleDoubleBiPredicate

use of uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate in project GDSC-SMLM by aherbert.

the class GradientCalculatorSpeedTest method gradientCalculatorComputesSameOutputWithBias.

@SeededTest
void gradientCalculatorComputesSameOutputWithBias(RandomSeed seed) {
    final Gaussian2DFunction func = new SingleEllipticalGaussian2DFunction(blockWidth, blockWidth);
    final int nparams = func.getNumberOfGradients();
    final GradientCalculator calc = new GradientCalculator(nparams);
    final int n = func.size();
    final int iter = 50;
    final ArrayList<double[]> paramsList = new ArrayList<>(iter);
    final ArrayList<double[]> yList = new ArrayList<>(iter);
    final ArrayList<double[][]> alphaList = new ArrayList<>(iter);
    final ArrayList<double[]> betaList = new ArrayList<>(iter);
    final ArrayList<double[]> xList = new ArrayList<>(iter);
    // Manipulate the background
    final double defaultBackground = background;
    final boolean report = logger.isLoggable(Level.INFO);
    try {
        background = 1e-2;
        createData(RngUtils.create(seed.getSeed()), 1, iter, paramsList, yList, true);
        final EjmlLinearSolver solver = new EjmlLinearSolver(1e-5, 1e-6);
        for (int i = 0; i < paramsList.size(); i++) {
            final double[] y = yList.get(i);
            final double[] a = paramsList.get(i);
            final double[][] alpha = new double[nparams][nparams];
            final double[] beta = new double[nparams];
            calc.findLinearised(n, y, a, alpha, beta, func);
            alphaList.add(alpha);
            betaList.add(beta.clone());
            for (int j = 0; j < nparams; j++) {
                if (Math.abs(beta[j]) < 1e-6) {
                    logger.info(FunctionUtils.getSupplier("[%d] Tiny beta %s %g", i, func.getGradientParameterName(j), beta[j]));
                }
            }
            // Solve
            if (!solver.solve(alpha, beta)) {
                throw new AssertionError();
            }
            xList.add(beta);
        // System.out.println(Arrays.toString(beta));
        }
        final double[][] alpha = new double[nparams][nparams];
        final double[] beta = new double[nparams];
        final Statistics[] rel = new Statistics[nparams];
        final Statistics[] abs = new Statistics[nparams];
        for (int i = 0; i < nparams; i++) {
            rel[i] = new Statistics();
            abs[i] = new Statistics();
        }
        final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
        // for (double b : new double[] { -500, -100, -10, -1, -0.1, 0.1, 1, 10, 100, 500 })
        for (final double b : new double[] { -10, -1, -0.1, 0.1, 1, 10 }) {
            if (report) {
                for (int i = 0; i < nparams; i++) {
                    rel[i].reset();
                    abs[i].reset();
                }
            }
            for (int i = 0; i < paramsList.size(); i++) {
                final double[] y = add(yList.get(i), b);
                final double[] a = paramsList.get(i).clone();
                a[0] += b;
                calc.findLinearised(n, y, a, alpha, beta, func);
                final double[][] alpha2 = alphaList.get(i);
                final double[] beta2 = betaList.get(i);
                final double[] x2 = xList.get(i);
                TestAssertions.assertArrayTest(beta2, beta, predicate, "Beta");
                TestAssertions.assertArrayTest(alpha2, alpha, predicate, "Alpha");
                // Solve
                solver.solve(alpha, beta);
                Assertions.assertArrayEquals(x2, beta, 1e-10, "X");
                if (report) {
                    for (int j = 0; j < nparams; j++) {
                        rel[j].add(DoubleEquality.relativeError(x2[j], beta[j]));
                        abs[j].add(Math.abs(x2[j] - beta[j]));
                    }
                }
            }
            if (report) {
                for (int i = 0; i < nparams; i++) {
                    logger.info(FunctionUtils.getSupplier("Bias = %.2f : %s : Rel %g +/- %g: Abs %g +/- %g", b, func.getGradientParameterName(i), rel[i].getMean(), rel[i].getStandardDeviation(), abs[i].getMean(), abs[i].getStandardDeviation()));
                }
            }
        }
    } finally {
        background = defaultBackground;
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) SingleEllipticalGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.SingleEllipticalGaussian2DFunction) ArrayList(java.util.ArrayList) Statistics(uk.ac.sussex.gdsc.core.utils.Statistics) SingleFreeCircularGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.SingleFreeCircularGaussian2DFunction) SingleEllipticalGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.SingleEllipticalGaussian2DFunction) SingleNbFixedGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.SingleNbFixedGaussian2DFunction) EllipticalGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.EllipticalGaussian2DFunction) SingleFixedGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.SingleFixedGaussian2DFunction) Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction) SingleCircularGaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.SingleCircularGaussian2DFunction) EjmlLinearSolver(uk.ac.sussex.gdsc.smlm.fitting.linear.EjmlLinearSolver) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Example 10 with DoubleDoubleBiPredicate

use of uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate in project GDSC-SMLM by aherbert.

the class PoissonGradientProcedureTest method gradientProcedureComputesSameAsGradientCalculator.

private void gradientProcedureComputesSameAsGradientCalculator(RandomSeed seed, int nparams) {
    final int iter = 10;
    final ArrayList<double[]> paramsList = new ArrayList<>(iter);
    createFakeParams(RngUtils.create(seed.getSeed()), nparams, iter, paramsList);
    final int n = blockWidth * blockWidth;
    final FakeGradientFunction func = new FakeGradientFunction(blockWidth, nparams);
    final GradientCalculator calc = GradientCalculatorUtils.newCalculator(nparams, false);
    // Create messages
    final IntArrayFormatSupplier msgOal = getMessage(nparams, "[%d] Observations: Not same alpha linear @ %d");
    final IntArrayFormatSupplier msgOam = getMessage(nparams, "[%d] Observations: Not same alpha matrix @ %d");
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    for (int i = 0; i < paramsList.size(); i++) {
        final PoissonGradientProcedure p = PoissonGradientProcedureUtils.create(func);
        p.computeFisherInformation(paramsList.get(i));
        final double[][] m = calc.fisherInformationMatrix(n, paramsList.get(i), func);
        // Not exactly the same ...
        final double[] al = p.getLinear();
        TestAssertions.assertArrayTest(al, new DenseMatrix64F(m).data, predicate, msgOal.set(1, i));
        final double[][] am = p.getMatrix();
        TestAssertions.assertArrayTest(am, m, predicate, msgOam.set(1, i));
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) ArrayList(java.util.ArrayList) IntArrayFormatSupplier(uk.ac.sussex.gdsc.test.utils.functions.IntArrayFormatSupplier) FakeGradientFunction(uk.ac.sussex.gdsc.smlm.function.FakeGradientFunction) DenseMatrix64F(org.ejml.data.DenseMatrix64F)

Aggregations

DoubleDoubleBiPredicate (uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate)63 SeededTest (uk.ac.sussex.gdsc.test.junit5.SeededTest)27 Test (org.junit.jupiter.api.Test)22 UniformRandomProvider (org.apache.commons.rng.UniformRandomProvider)12 PoissonDistribution (org.apache.commons.math3.distribution.PoissonDistribution)6 TDoubleArrayList (gnu.trove.list.array.TDoubleArrayList)4 ArrayList (java.util.ArrayList)4 MixtureMultivariateGaussianDistribution (uk.ac.sussex.gdsc.smlm.math3.distribution.fitting.MultivariateGaussianMixtureExpectationMaximization.MixtureMultivariateGaussianDistribution)4 MultivariateGaussianDistribution (uk.ac.sussex.gdsc.smlm.math3.distribution.fitting.MultivariateGaussianMixtureExpectationMaximization.MixtureMultivariateGaussianDistribution.MultivariateGaussianDistribution)4 UnivariateFunction (org.apache.commons.math3.analysis.UnivariateFunction)3 SimpsonIntegrator (org.apache.commons.math3.analysis.integration.SimpsonIntegrator)3 ArrayRealVector (org.apache.commons.math3.linear.ArrayRealVector)3 RealMatrix (org.apache.commons.math3.linear.RealMatrix)3 RealVector (org.apache.commons.math3.linear.RealVector)3 ParameterizedTest (org.junit.jupiter.params.ParameterizedTest)3 BigDecimal (java.math.BigDecimal)2 MathContext (java.math.MathContext)2 MixtureMultivariateNormalDistribution (org.apache.commons.math3.distribution.MixtureMultivariateNormalDistribution)2 MultivariateNormalDistribution (org.apache.commons.math3.distribution.MultivariateNormalDistribution)2 MultivariateJacobianFunction (org.apache.commons.math3.fitting.leastsquares.MultivariateJacobianFunction)2