Search in sources :

Example 91 with UniformRandomProvider

use of org.apache.commons.rng.UniformRandomProvider in project GDSC-SMLM by aherbert.

the class MultivariateGaussianMixtureExpectationMaximizationTest method testExpectationMaximizationSpeedWithDifferentNumberOfComponents.

/**
 * Test the speed of implementations of the expectation maximization algorithm with a mixture of n
 * 2D Gaussian distributions.
 *
 * @param seed the seed
 */
@SpeedTag
@SeededTest
void testExpectationMaximizationSpeedWithDifferentNumberOfComponents(RandomSeed seed) {
    Assumptions.assumeTrue(TestSettings.allow(TestComplexity.HIGH));
    // Create data
    final UniformRandomProvider rng = RngUtils.create(seed.getSeed());
    for (int n = 2; n <= 4; n++) {
        final double[][][] data = new double[10][][];
        for (int i = 0; i < data.length; i++) {
            final double[] sampleWeights = createWeights(n, rng);
            final double[][] sampleMeans = create(n, 2, rng, -5, 5);
            final double[][] sampleStdDevs = create(n, 2, rng, 1, 10);
            final double[] sampleCorrelations = create(n, rng, -0.9, 0.9);
            data[i] = createData2d(1000, rng, sampleWeights, sampleMeans, sampleStdDevs, sampleCorrelations);
        }
        final int numComponents = n;
        // Time initial estimation and fitting
        final TimingService ts = new TimingService();
        ts.execute(new FittingSpeedTask("Commons n=" + n + " 2D", data) {

            @Override
            Object run(double[][] data) {
                final MultivariateNormalMixtureExpectationMaximization fitter = new MultivariateNormalMixtureExpectationMaximization(data);
                fitter.fit(MultivariateNormalMixtureExpectationMaximization.estimate(data, numComponents));
                return fitter.getLogLikelihood();
            }
        });
        ts.execute(new FittingSpeedTask("GDSC n=" + n + " 2D", data) {

            @Override
            Object run(double[][] data) {
                final MultivariateGaussianMixtureExpectationMaximization fitter = new MultivariateGaussianMixtureExpectationMaximization(data);
                fitter.fit(MultivariateGaussianMixtureExpectationMaximization.estimate(data, numComponents));
                return fitter.getLogLikelihood();
            }
        });
        if (logger.isLoggable(Level.INFO)) {
            logger.info(ts.getReport());
        }
        // More than twice as fast
        Assertions.assertTrue(ts.get(-1).getMean() < ts.get(-2).getMean() / 2);
    }
}
Also used : UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) TimingService(uk.ac.sussex.gdsc.test.utils.TimingService) MultivariateNormalMixtureExpectationMaximization(org.apache.commons.math3.distribution.fitting.MultivariateNormalMixtureExpectationMaximization) SpeedTag(uk.ac.sussex.gdsc.test.junit5.SpeedTag) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Example 92 with UniformRandomProvider

use of org.apache.commons.rng.UniformRandomProvider in project GDSC-SMLM by aherbert.

the class PerPixelCameraModelTest method canCropAndConvertDataWithCropBounds.

@SeededTest
void canCropAndConvertDataWithCropBounds(RandomSeed seed) {
    final PerPixelCameraModelTestData data = (PerPixelCameraModelTestData) dataCache.computeIfAbsent(seed, PerPixelCameraModelTest::createData);
    final PerPixelCameraModel model = new PerPixelCameraModel(w, h, data.bias, data.gain, data.variance);
    final UniformRandomProvider rand = RngUtils.create(seed.getSeed());
    final ImageExtractor ie = ImageExtractor.wrap(data.bias, w, h);
    for (int j = 0; j < 10; j++) {
        final Rectangle bounds = getBounds(rand, ie);
        checkConversion(data, bounds, model.crop(bounds, false));
    }
}
Also used : Rectangle(java.awt.Rectangle) UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) ImageExtractor(uk.ac.sussex.gdsc.core.utils.ImageExtractor) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Example 93 with UniformRandomProvider

use of org.apache.commons.rng.UniformRandomProvider in project GDSC-SMLM by aherbert.

the class PerPixelCameraModelTest method canGetMeanVarianceData.

private static void canGetMeanVarianceData(RandomSeed seed, boolean initialise, boolean normalised) {
    final PerPixelCameraModelTestData data = (PerPixelCameraModelTestData) dataCache.computeIfAbsent(seed, PerPixelCameraModelTest::createData);
    final PerPixelCameraModel model = createModel(data, initialise);
    final UniformRandomProvider rand = RngUtils.create(seed.getSeed());
    final ImageExtractor ie = ImageExtractor.wrap(data.bias, w, h);
    for (int i = 0; i < 10; i++) {
        final Rectangle bounds = getBounds(rand, ie);
        final float[] v = (normalised) ? model.getNormalisedVariance(bounds) : model.getVariance(bounds);
        final double e = MathUtils.sum(v) / v.length;
        final double o = (normalised) ? model.getMeanNormalisedVariance(bounds) : model.getMeanVariance(bounds);
        Assertions.assertEquals(e, o);
    }
}
Also used : Rectangle(java.awt.Rectangle) UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) ImageExtractor(uk.ac.sussex.gdsc.core.utils.ImageExtractor)

Example 94 with UniformRandomProvider

use of org.apache.commons.rng.UniformRandomProvider in project GDSC-SMLM by aherbert.

the class Gaussian2DPeakResultHelperTest method canComputeMeanSignalUsingP.

@SeededTest
void canComputeMeanSignalUsingP(RandomSeed seed) {
    final UniformRandomProvider rg = RngUtils.create(seed.getSeed());
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    for (int i = 0; i < 10; i++) {
        final double intensity = rg.nextDouble() * 100;
        final double sx = rg.nextDouble() * 2;
        final double sy = rg.nextDouble() * 2;
        final double p = rg.nextDouble();
        double expected = intensity * p / (Math.PI * MathUtils.pow2(Gaussian2DPeakResultHelper.inverseCumulative2D(p)) * sx * sy);
        double observed = Gaussian2DPeakResultHelper.getMeanSignalUsingP(intensity, sx, sy, p);
        assertEquals(expected, observed, predicate);
        // Test fixed versions verse dynamic
        expected = Gaussian2DPeakResultHelper.getMeanSignalUsingP(intensity, sx, sy, 0.5);
        observed = Gaussian2DPeakResultHelper.getMeanSignalUsingP05(intensity, sx, sy);
        assertEquals(expected, observed, predicate);
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Example 95 with UniformRandomProvider

use of org.apache.commons.rng.UniformRandomProvider in project GDSC-SMLM by aherbert.

the class Gaussian2DPeakResultHelperTest method canComputeMeanSignalUsingR.

@SeededTest
void canComputeMeanSignalUsingR(RandomSeed seed) {
    final UniformRandomProvider rg = RngUtils.create(seed.getSeed());
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    for (int i = 0; i < 10; i++) {
        final double intensity = rg.nextDouble() * 100;
        final double sx = rg.nextDouble() * 2;
        final double sy = rg.nextDouble() * 2;
        final double r = rg.nextDouble() * 5;
        assertEquals(intensity * Gaussian2DPeakResultHelper.cumulative2D(r) / (Math.PI * r * r * sx * sy), Gaussian2DPeakResultHelper.getMeanSignalUsingR(intensity, sx, sy, r), predicate);
        // Test fixed versions verse dynamic
        assertEquals(Gaussian2DPeakResultHelper.getMeanSignalUsingR(intensity, sx, sy, 1), Gaussian2DPeakResultHelper.getMeanSignalUsingR1(intensity, sx, sy), predicate);
        assertEquals(Gaussian2DPeakResultHelper.getMeanSignalUsingR(intensity, sx, sy, 2), Gaussian2DPeakResultHelper.getMeanSignalUsingR2(intensity, sx, sy), predicate);
    }
}
Also used : DoubleDoubleBiPredicate(uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate) UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Aggregations

UniformRandomProvider (org.apache.commons.rng.UniformRandomProvider)213 SeededTest (uk.ac.sussex.gdsc.test.junit5.SeededTest)145 SharedStateContinuousSampler (org.apache.commons.rng.sampling.distribution.SharedStateContinuousSampler)17 TimingService (uk.ac.sussex.gdsc.test.utils.TimingService)14 Rectangle (java.awt.Rectangle)13 DoubleDoubleBiPredicate (uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate)13 SpeedTag (uk.ac.sussex.gdsc.test.junit5.SpeedTag)12 TDoubleArrayList (gnu.trove.list.array.TDoubleArrayList)10 ArrayList (java.util.ArrayList)10 NormalizedGaussianSampler (org.apache.commons.rng.sampling.distribution.NormalizedGaussianSampler)9 StoredDataStatistics (uk.ac.sussex.gdsc.core.utils.StoredDataStatistics)8 CalibrationWriter (uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter)8 ContinuousSampler (org.apache.commons.rng.sampling.distribution.ContinuousSampler)6 ImageExtractor (uk.ac.sussex.gdsc.core.utils.ImageExtractor)6 LocalList (uk.ac.sussex.gdsc.core.utils.LocalList)6 Gaussian2DFunction (uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction)6 FloatProcessor (ij.process.FloatProcessor)5 ErfGaussian2DFunction (uk.ac.sussex.gdsc.smlm.function.gaussian.erf.ErfGaussian2DFunction)5 TimingResult (uk.ac.sussex.gdsc.test.utils.TimingResult)5 File (java.io.File)4