Examples with DoubleDoubleBiPredicate uk.ac.sussex.gdsc.test.api.function.DoubleDoubleBiPredicate used on opensource projects

Example 11 with DoubleDoubleBiPredicate

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

the class LogFactorialTest method testLogFactorialDouble.

/**
 * Test the factorial of a fractional number against Commons Math logGamma(1+n).
 */
@SeededTest
void testLogFactorialDouble(RandomSeed seed) {
    final UniformRandomProvider rng = RngUtils.create(seed.getSeed());
    final DoubleDoubleBiPredicate tol = TestHelper.doublesAreClose(1e-14);
    for (int i = 0; i < 200; i++) {
        final double n = rng.nextDouble() * 200;
        final double expected = n <= 1.5 ? Gamma.logGamma1p(n) : Gamma.logGamma(1 + n);
        TestAssertions.assertTest(expected, LogFactorial.value(n), tol, () -> Double.toString(n));
    }
}

Example 12 with DoubleDoubleBiPredicate

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

the class PoissonCalculatorTest method canComputeLogLikelihoodRatio.

private static void canComputeLogLikelihoodRatio(RandomSeed seed, BaseNonLinearFunction nlf) {
    logger.log(TestLogUtils.getRecord(LOG_LEVEL, nlf.name));
    final int n = maxx * maxx;
    final double[] a = new double[] { 1 };
    // Simulate Poisson process
    nlf.initialise(a);
    final UniformRandomProvider rng = RngUtils.create(seed.getSeed());
    final double[] x = new double[n];
    final double[] u = new double[n];
    for (int i = 0; i < n; i++) {
        u[i] = nlf.eval(i);
        if (u[i] > 0) {
            x[i] = GdscSmlmTestUtils.createPoissonSampler(rng, u[i]).sample();
        }
    }
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-10, 0);
    double ll = PoissonCalculator.logLikelihood(u, x);
    final double mll = PoissonCalculator.maximumLogLikelihood(x);
    double llr = -2 * (ll - mll);
    double llr2 = PoissonCalculator.logLikelihoodRatio(u, x);
    logger.log(TestLogUtils.getRecord(LOG_LEVEL, "llr=%f, llr2=%f", llr, llr2));
    TestAssertions.assertTest(llr, llr2, predicate, "Log-likelihood ratio");
    final double[] op = new double[x.length];
    for (int i = 0; i < n; i++) {
        op[i] = PoissonCalculator.maximumLikelihood(x[i]);
    }
    // TestSettings.setLogLevel(uk.ac.sussex.gdsc.smlm.TestLog.Level.FINE);
    final int df = n - 1;
    final ChiSquaredDistributionTable table = ChiSquaredDistributionTable.createUpperTailed(0.05, df);
    final ChiSquaredDistributionTable table2 = ChiSquaredDistributionTable.createUpperTailed(0.001, df);
    if (logger.isLoggable(LOG_LEVEL)) {
        logger.log(TestLogUtils.getRecord(LOG_LEVEL, "Chi2 = %f (q=%.3f), %f (q=%.3f)  %f %b  %f", table.getCrititalValue(df), table.getSignificanceValue(), table2.getCrititalValue(df), table2.getSignificanceValue(), ChiSquaredDistributionTable.computeQValue(24, 2), ChiSquaredDistributionTable.createUpperTailed(0.05, 2).reject(24, 2), ChiSquaredDistributionTable.getChiSquared(1e-6, 2)));
    }
    final TDoubleArrayList list = new TDoubleArrayList();
    final int imin = 5;
    final int imax = 15;
    for (int i = imin; i <= imax; i++) {
        a[0] = (double) i / 10;
        nlf.initialise(a);
        for (int j = 0; j < n; j++) {
            u[j] = nlf.eval(j);
        }
        ll = PoissonCalculator.logLikelihood(u, x);
        list.add(ll);
        llr = PoissonCalculator.logLikelihoodRatio(u, x);
        BigDecimal product = new BigDecimal(1);
        double ll2 = 0;
        for (int j = 0; j < n; j++) {
            final double p1 = PoissonCalculator.likelihood(u[j], x[j]);
            ll2 += Math.log(p1);
            final double ratio = p1 / op[j];
            product = product.multiply(new BigDecimal(ratio));
        }
        llr2 = -2 * Math.log(product.doubleValue());
        final double p = ChiSquaredDistributionTable.computePValue(llr, df);
        final double q = ChiSquaredDistributionTable.computeQValue(llr, df);
        logger.log(TestLogUtils.getRecord(LOG_LEVEL, "a=%f, ll=%f, ll2=%f, llr=%f, llr2=%f, product=%s, p=%f, q=%f " + "(reject=%b @ %.3f, reject=%b @ %.3f)", a[0], ll, ll2, llr, llr2, product.round(new MathContext(4)).toString(), p, q, table.reject(llr, df), table.getSignificanceValue(), table2.reject(llr, df), table2.getSignificanceValue()));
        // too small to store in a double.
        if (product.doubleValue() > 0) {
            TestAssertions.assertTest(ll, ll2, predicate, "Log-likelihood");
            TestAssertions.assertTest(llr, llr2, predicate, "Log-likelihood ratio");
        }
    }
    // Find max using quadratic fit
    final double[] data = list.toArray();
    int index = SimpleArrayUtils.findMaxIndex(data);
    final double maxa = (double) (imin + index) / 10;
    double fita = maxa;
    try {
        if (index == 0) {
            index++;
        }
        if (index == data.length - 1) {
            index--;
        }
        final int i1 = index - 1;
        final int i2 = index;
        final int i3 = index + 1;
        fita = QuadraticUtils.findMinMax((double) (imin + i1) / 10, data[i1], (double) (imin + i2) / 10, data[i2], (double) (imin + i3) / 10, data[i3]);
    } catch (final DataException ex) {
    // Ignore
    }
    // Allow a tolerance as the random data may alter the p-value computation.
    // Should allow it to be less than 2 increment either side of the answer.
    logger.log(TestLogUtils.getRecord(LOG_LEVEL, "max fit = %g => %g", maxa, fita));
    Assertions.assertEquals(1, fita, 0.199, "max");
}

Example 13 with DoubleDoubleBiPredicate

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

the class PoissonCalculatorTest method canComputeFastLog_FastLikelihoodForIntegerData.

@Test
void canComputeFastLog_FastLikelihoodForIntegerData() {
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-4, 0);
    final FastLog fastLog = FastLogFactory.getFastLog();
    for (final double u : photons) {
        final PoissonDistribution pd = new PoissonDistribution(u);
        for (int x = 0; x < 100; x++) {
            double expected = pd.probability(x);
            double observed = PoissonCalculator.fastLikelihood(u, x, fastLog);
            if (expected > 1e-100) {
                TestAssertions.assertTest(expected, observed, predicate);
            }
            expected = pd.logProbability(x);
            observed = PoissonCalculator.fastLogLikelihood(u, x, fastLog);
            TestAssertions.assertTest(expected, observed, predicate);
        }
    }
}

Example 14 with DoubleDoubleBiPredicate

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

the class PoissonCalculatorTest method canComputeFastLikelihoodForIntegerData.

@Test
void canComputeFastLikelihoodForIntegerData() {
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-4, 0);
    for (final double u : photons) {
        final PoissonDistribution pd = new PoissonDistribution(u);
        for (int x = 0; x < 100; x++) {
            double expected = pd.probability(x);
            double observed = PoissonCalculator.fastLikelihood(u, x);
            if (expected > 1e-100) {
                TestAssertions.assertTest(expected, observed, predicate);
            }
            expected = pd.logProbability(x);
            observed = PoissonCalculator.fastLogLikelihood(u, x);
            TestAssertions.assertTest(expected, observed, predicate);
        }
    }
}

Example 15 with DoubleDoubleBiPredicate

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

the class PoissonGammaFunctionTest method probabilityMatchesLogProbability.

private static void probabilityMatchesLogProbability(final double gain, double mu) {
    final PoissonGammaFunction f = PoissonGammaFunction.createWithAlpha(1.0 / gain);
    // Evaluate an initial range.
    // Gaussian should have >99% within +/- s
    // Poisson will have mean mu with a variance mu.
    // At large mu it is approximately normal so use 3 sqrt(mu) for the range added to the mean
    // Note: The input s parameter is after-gain so adjust.
    final int[] range = PoissonGaussianFunctionTest.getRange(gain, mu, 0);
    final int min = range[0];
    final int max = range[1];
    // Note: The input mu parameter is pre-gain.
    final double e = mu;
    final Supplier<String> msg = () -> String.format("g=%f, mu=%f", gain, mu);
    final DoubleDoubleBiPredicate predicate = TestHelper.doublesAreClose(1e-6, 0);
    for (int x = min; x <= max; x++) {
        final double p = f.likelihood(x, e);
        if (p == 0) {
            continue;
        }
        final double logP = f.logLikelihood(x, e);
        TestAssertions.assertTest(Math.log(p), logP, predicate, msg);
    }
}