Search in sources :

Example 6 with GradientCalculator

use of uk.ac.sussex.gdsc.smlm.fitting.nonlinear.gradient.GradientCalculator in project GDSC-SMLM by aherbert.

the class FisherInformationMatrixTest method createFisherInformationMatrix.

private static FisherInformationMatrix createFisherInformationMatrix(UniformRandomProvider rg, int columns, int zeroColumns) {
    final int maxx = 10;
    final int size = maxx * maxx;
    // Use a real Gaussian function here to compute the Fisher information.
    // The matrix may be sensitive to the type of equation used.
    int npeazeroColumnss = 1;
    Gaussian2DFunction fun = createFunction(maxx, npeazeroColumnss);
    while (fun.getNumberOfGradients() < columns) {
        npeazeroColumnss++;
        fun = createFunction(maxx, npeazeroColumnss);
    }
    final double[] a = new double[1 + npeazeroColumnss * Gaussian2DFunction.PARAMETERS_PER_PEAK];
    a[Gaussian2DFunction.BACKGROUND] = nextUniform(rg, 1, 5);
    for (int i = 0, j = 0; i < npeazeroColumnss; i++, j += Gaussian2DFunction.PARAMETERS_PER_PEAK) {
        a[j + Gaussian2DFunction.SIGNAL] = nextUniform(rg, 100, 300);
        // Non-overlapping peazeroColumnss otherwise the Crlb are poor
        a[j + Gaussian2DFunction.X_POSITION] = nextUniform(rg, 2 + i * 2, 4 + i * 2);
        a[j + Gaussian2DFunction.Y_POSITION] = nextUniform(rg, 2 + i * 2, 4 + i * 2);
        a[j + Gaussian2DFunction.X_SD] = nextUniform(rg, 1.5, 2);
        a[j + Gaussian2DFunction.Y_SD] = nextUniform(rg, 1.5, 2);
    }
    fun.initialise(a);
    final GradientCalculator calc = GradientCalculatorUtils.newCalculator(fun.getNumberOfGradients());
    double[][] matrixI = calc.fisherInformationMatrix(size, a, fun);
    // Reduce to the desired size
    matrixI = Arrays.copyOf(matrixI, columns);
    for (int i = 0; i < columns; i++) {
        matrixI[i] = Arrays.copyOf(matrixI[i], columns);
    }
    // Zero selected columns
    if (zeroColumns > 0) {
        final int[] zero = RandomUtils.sample(zeroColumns, columns, rg);
        for (final int i : zero) {
            for (int j = 0; j < columns; j++) {
                matrixI[i][j] = matrixI[j][i] = 0;
            }
        }
    }
    // Create matrix
    return new FisherInformationMatrix(matrixI, 1e-3);
}
Also used : Gaussian2DFunction(uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction) GradientCalculator(uk.ac.sussex.gdsc.smlm.fitting.nonlinear.gradient.GradientCalculator)

Aggregations

GradientCalculator (uk.ac.sussex.gdsc.smlm.fitting.nonlinear.gradient.GradientCalculator)6 Gaussian2DFunction (uk.ac.sussex.gdsc.smlm.function.gaussian.Gaussian2DFunction)3 UniformRandomProvider (org.apache.commons.rng.UniformRandomProvider)2 DenseMatrix64F (org.ejml.data.DenseMatrix64F)2 LocalList (uk.ac.sussex.gdsc.core.utils.LocalList)2 ValueProcedure (uk.ac.sussex.gdsc.smlm.function.ValueProcedure)2 TimingService (uk.ac.sussex.gdsc.test.utils.TimingService)2 FisherInformationMatrix (uk.ac.sussex.gdsc.smlm.fitting.FisherInformationMatrix)1 ExtendedNonLinearFunction (uk.ac.sussex.gdsc.smlm.function.ExtendedNonLinearFunction)1 NonLinearFunction (uk.ac.sussex.gdsc.smlm.function.NonLinearFunction)1 SingleFreeCircularGaussian2DFunction (uk.ac.sussex.gdsc.smlm.function.gaussian.SingleFreeCircularGaussian2DFunction)1