Examples with Gradient1Function uk.ac.sussex.gdsc.smlm.function.Gradient1Function used on opensource projects

Example 1 with Gradient1Function

use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.

the class WPoissonGradientProcedureTest method gradientProcedureUnrolledComputesSameAsGradientProcedure.

private void gradientProcedureUnrolledComputesSameAsGradientProcedure(RandomSeed seed, int nparams, boolean precomputed) {
    final int iter = 10;
    final ArrayList<double[]> paramsList = new ArrayList<>(iter);
    final UniformRandomProvider rng = RngUtils.create(seed.getSeed());
    createFakeParams(rng, nparams, iter, paramsList);
    final Gradient1Function func = new FakeGradientFunction(blockWidth, nparams);
    final double[] v = (precomputed) ? dataCache.computeIfAbsent(seed, WPoissonGradientProcedureTest::createData) : null;
    final IntArrayFormatSupplier msg = getMessage(nparams, "[%d] Observations: Not same linear @ %d");
    for (int i = 0; i < paramsList.size(); i++) {
        final double[] y = createFakeData(rng);
        final WPoissonGradientProcedure p1 = new WPoissonGradientProcedure(y, v, func);
        p1.computeFisherInformation(paramsList.get(i));
        final WPoissonGradientProcedure p2 = WPoissonGradientProcedureUtils.create(y, v, func);
        p2.computeFisherInformation(paramsList.get(i));
        // Exactly the same ...
        Assertions.assertArrayEquals(p1.getLinear(), p2.getLinear(), msg.set(1, i));
    }
}

Example 2 with Gradient1Function

use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.

the class UnivariateLikelihoodFisherInformationCalculatorTest method canComputePerPixelPoissonGaussianApproximationFisherInformation.

private static void canComputePerPixelPoissonGaussianApproximationFisherInformation(UniformRandomProvider rng) {
    // Create function
    final Gaussian2DFunction func = GaussianFunctionFactory.create2D(1, 10, 10, GaussianFunctionFactory.FIT_ERF_CIRCLE, null);
    final double[] params = new double[1 + Gaussian2DFunction.PARAMETERS_PER_PEAK];
    params[Gaussian2DFunction.BACKGROUND] = nextUniform(rng, 0.1, 0.3);
    params[Gaussian2DFunction.SIGNAL] = nextUniform(rng, 100, 300);
    params[Gaussian2DFunction.X_POSITION] = nextUniform(rng, 4, 6);
    params[Gaussian2DFunction.Y_POSITION] = nextUniform(rng, 4, 6);
    params[Gaussian2DFunction.X_SD] = nextUniform(rng, 1, 1.3);
    Gradient1Function f1 = func;
    FisherInformation[] fi;
    // Get a per-pixel variance
    final double[] var = new double[func.size()];
    fi = new FisherInformation[var.length];
    for (int i = var.length; i-- > 0; ) {
        var[i] = 0.9 + 0.2 * rng.nextDouble();
        fi[i] = new PoissonGaussianApproximationFisherInformation(Math.sqrt(var[i]));
    }
    f1 = (Gradient1Function) OffsetFunctionFactory.wrapFunction(func, var);
    // This introduces a dependency on a different package, and relies on that
    // computing the correct answer. However that code predates this and so the
    // test ensures that the FisherInformationCalculator functions correctly.
    final PoissonGradientProcedure p1 = PoissonGradientProcedureUtils.create(f1);
    p1.computeFisherInformation(params);
    final double[] e = p1.getLinear();
    final FisherInformationCalculator calc = new UnivariateLikelihoodFisherInformationCalculator(func, fi);
    final FisherInformationMatrix I = calc.compute(params);
    final double[] o = I.getMatrix().data;
    TestAssertions.assertArrayTest(e, o, TestHelper.doublesAreClose(1e-6, 0));
}

Example 3 with Gradient1Function

use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.

the class LseLvmSteppingFunctionSolver method computeDeviationsAndValues.

@Override
protected void computeDeviationsAndValues(double[] parametersVariance, double[] fx) {
    Gradient1Function f1 = (Gradient1Function) this.function;
    // Capture the y-values if necessary
    if (fx != null && fx.length == f1.size()) {
        f1 = new Gradient2FunctionValueStore(f1, fx);
    }
    final LsqVarianceGradientProcedure p = createVarianceProcedure(f1);
    if (p.variance(null) == LsqVarianceGradientProcedure.STATUS_OK) {
        setDeviations(parametersVariance, p.variance);
    }
}

Example 4 with Gradient1Function

use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.

the class SteppingFunctionSolver method computeValue.

@Override
protected boolean computeValue(double[] y, double[] fx, double[] a) {
    // If the fx array is not null then wrap the gradient function.
    // Compute the value and the wrapper will store the values appropriately.
    // Then reset the gradient function.
    // Note: If a sub class wraps the function with weights
    // then the weights will not be stored in the function value.
    // Only the value produced by the original function is stored:
    // Wrapped (+weights) < FunctionStore < Function
    // However if the base function is already wrapped then this will occur:
    // Wrapped (+weights) < FunctionStore < Wrapped (+precomputed) < Function
    gradientIndices = function.gradientIndices();
    if (fx != null && fx.length == ((Gradient1Function) function).size()) {
        final GradientFunction tmp = function;
        function = new Gradient1FunctionStore((Gradient1Function) function, fx, null);
        lastY = prepareFunctionValue(y, a);
        value = computeFunctionValue(a);
        function = tmp;
    } else {
        lastY = prepareFunctionValue(y, a);
        value = computeFunctionValue(a);
    }
    return true;
}

Example 5 with Gradient1Function

use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.

the class LsqVarianceGradientProcedureTest method gradientProcedureIsFasterUnrolledThanGradientProcedure.

private void gradientProcedureIsFasterUnrolledThanGradientProcedure(RandomSeed seed, final int nparams, final boolean precomputed) {
    Assumptions.assumeTrue(TestSettings.allow(TestComplexity.MEDIUM));
    final int iter = 100;
    final ArrayList<double[]> paramsList = new ArrayList<>(iter);
    createFakeParams(RngUtils.create(seed.getSeed()), nparams, iter, paramsList);
    // Remove the timing of the function call by creating a dummy function
    final FakeGradientFunction f = new FakeGradientFunction(blockWidth, nparams);
    final Gradient1Function func = (precomputed) ? OffsetGradient1Function.wrapGradient1Function(f, SimpleArrayUtils.newArray(f.size(), 0.1, 1.3)) : f;
    final IndexSupplier msg = new IndexSupplier(1, "M ", null);
    for (int i = 0; i < paramsList.size(); i++) {
        final LsqVarianceGradientProcedure p1 = new LsqVarianceGradientProcedure(func);
        p1.variance(paramsList.get(i));
        p1.variance(paramsList.get(i));
        final LsqVarianceGradientProcedure p2 = LsqVarianceGradientProcedureUtils.create(func);
        p2.variance(paramsList.get(i));
        p2.variance(paramsList.get(i));
        // Check they are the same
        Assertions.assertArrayEquals(p1.variance, p2.variance, msg.set(0, i));
    }
    // Realistic loops for an optimisation
    final int loops = 15;
    // Run till stable timing
    final Timer t1 = new Timer() {

        @Override
        void run() {
            for (int i = 0, k = 0; i < paramsList.size(); i++) {
                final LsqVarianceGradientProcedure p1 = new LsqVarianceGradientProcedure(func);
                for (int j = loops; j-- > 0; ) {
                    p1.variance(paramsList.get(k++ % iter));
                }
            }
        }
    };
    final long time1 = t1.getTime();
    final Timer t2 = new Timer(t1.loops) {

        @Override
        void run() {
            for (int i = 0, k = 0; i < paramsList.size(); i++) {
                final LsqVarianceGradientProcedure p2 = LsqVarianceGradientProcedureUtils.create(func);
                for (int j = loops; j-- > 0; ) {
                    p2.variance(paramsList.get(k++ % iter));
                }
            }
        }
    };
    final long time2 = t2.getTime();
    logger.log(TestLogUtils.getTimingRecord("precomputed=" + precomputed + " Standard " + nparams, time1, "Unrolled", time2));
}