use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.
the class LsqVarianceGradientProcedureTest method gradientProcedureUnrolledComputesSameAsGradientProcedure.
private void gradientProcedureUnrolledComputesSameAsGradientProcedure(RandomSeed seed, int nparams, boolean precomputed) {
final int iter = 10;
final ArrayList<double[]> paramsList = new ArrayList<>(iter);
createFakeParams(RngUtils.create(seed.getSeed()), nparams, iter, paramsList);
Gradient1Function func = new FakeGradientFunction(blockWidth, nparams);
if (precomputed) {
func = OffsetGradient1Function.wrapGradient1Function(func, SimpleArrayUtils.newArray(func.size(), 0.1, 1.3));
}
final IntArrayFormatSupplier msg = new IntArrayFormatSupplier("[%d] Observations: Not same variance @ %d", 2);
msg.set(0, nparams);
for (int i = 0; i < paramsList.size(); i++) {
final LsqVarianceGradientProcedure p1 = new LsqVarianceGradientProcedure(func);
p1.variance(paramsList.get(i));
final LsqVarianceGradientProcedure p2 = LsqVarianceGradientProcedureUtils.create(func);
p2.variance(paramsList.get(i));
// Exactly the same ...
Assertions.assertArrayEquals(p1.variance, p2.variance, msg.set(1, i));
}
}
use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.
the class PoissonGradientProcedureTest method gradientProcedureUnrolledComputesSameAsGradientProcedure.
private void gradientProcedureUnrolledComputesSameAsGradientProcedure(RandomSeed seed, int nparams, boolean precomputed) {
final int iter = 10;
final ArrayList<double[]> paramsList = new ArrayList<>(iter);
createFakeParams(RngUtils.create(seed.getSeed()), nparams, iter, paramsList);
Gradient1Function func = new FakeGradientFunction(blockWidth, nparams);
if (precomputed) {
func = OffsetGradient1Function.wrapGradient1Function(func, SimpleArrayUtils.newArray(func.size(), 0.1, 1.3));
}
// 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");
for (int i = 0; i < paramsList.size(); i++) {
final PoissonGradientProcedure p1 = new PoissonGradientProcedure(func);
p1.computeFisherInformation(paramsList.get(i));
final PoissonGradientProcedure p2 = PoissonGradientProcedureUtils.create(func);
p2.computeFisherInformation(paramsList.get(i));
// Exactly the same ...
Assertions.assertArrayEquals(p1.getLinear(), p2.getLinear(), msgOal.set(1, i));
final double[][] am1 = p1.getMatrix();
final double[][] am2 = p2.getMatrix();
Assertions.assertArrayEquals(am1, am2, msgOam.set(1, i));
}
}
use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.
the class FastMleSteppingFunctionSolver method computeFunctionFisherInformationMatrix.
@Override
protected FisherInformationMatrix computeFunctionFisherInformationMatrix(double[] y, double[] a) {
// The fisher information is that for a Poisson process.
// We must wrap the gradient function if weights are present.
Gradient1Function f1 = (Gradient1Function) function;
if (obsVariances != null) {
f1 = OffsetGradient1Function.wrapGradient1Function(f1, obsVariances);
}
final PoissonGradientProcedure p = PoissonGradientProcedureUtils.create(f1);
p.computeFisherInformation(a);
if (p.isNaNGradients()) {
throw new FunctionSolverException(FitStatus.INVALID_GRADIENTS);
}
return new FisherInformationMatrix(p.getLinear(), p.numberOfGradients);
}
use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.
the class LsqLvmGradientProcedureTest method gradientProcedureLinearIsFasterThanGradientProcedureMatrix.
private void gradientProcedureLinearIsFasterThanGradientProcedureMatrix(RandomSeed seed, final int nparams, final BaseLsqLvmGradientProcedureFactory factory1, final BaseLsqLvmGradientProcedureFactory factory2, boolean doAssert) {
Assumptions.assumeTrue(TestSettings.allow(TestComplexity.MEDIUM));
final int iter = 100;
final ArrayList<double[]> paramsList = new ArrayList<>(iter);
final ArrayList<double[]> yList = new ArrayList<>(iter);
createData(RngUtils.create(seed.getSeed()), 1, iter, paramsList, yList);
// Remove the timing of the function call by creating a dummy function
final Gradient1Function func = new FakeGradientFunction(blockWidth, nparams);
// Create messages
final IndexSupplier msgA = new IndexSupplier(1, "A ", null);
final IndexSupplier msgB = new IndexSupplier(1, "B ", null);
for (int i = 0; i < paramsList.size(); i++) {
final BaseLsqLvmGradientProcedure p1 = factory1.createProcedure(yList.get(i), func);
p1.gradient(paramsList.get(i));
p1.gradient(paramsList.get(i));
final BaseLsqLvmGradientProcedure p2 = factory2.createProcedure(yList.get(i), func);
p2.gradient(paramsList.get(i));
p2.gradient(paramsList.get(i));
// Check they are the same
Assertions.assertArrayEquals(p1.getAlphaLinear(), p2.getAlphaLinear(), msgA.set(0, i));
Assertions.assertArrayEquals(p1.beta, p2.beta, msgB.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 BaseLsqLvmGradientProcedure p = factory1.createProcedure(yList.get(i), func);
for (int j = loops; j-- > 0; ) {
p.gradient(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 BaseLsqLvmGradientProcedure p2 = factory2.createProcedure(yList.get(i), func);
for (int j = loops; j-- > 0; ) {
p2.gradient(paramsList.get(k++ % iter));
}
}
}
};
final long time2 = t2.getTime();
final LogRecord record = TestLogUtils.getTimingRecord(factory1.getClass().getSimpleName() + nparams, time1, factory2.getClass().getSimpleName(), time2);
if (!doAssert && record.getLevel() == TestLevel.TEST_FAILURE) {
record.setLevel(TestLevel.TEST_WARNING);
}
logger.log(record);
}
use of uk.ac.sussex.gdsc.smlm.function.Gradient1Function in project GDSC-SMLM by aherbert.
the class LvmGradientProcedureTest method gradientProcedureIsFasterUnrolledThanGradientProcedure.
private void gradientProcedureIsFasterUnrolledThanGradientProcedure(RandomSeed seed, final int nparams, final Type type, final boolean precomputed) {
Assumptions.assumeTrue(TestSettings.allow(TestComplexity.MEDIUM));
final int iter = 100;
final ArrayList<double[]> paramsList = new ArrayList<>(iter);
final ArrayList<double[]> yList = new ArrayList<>(iter);
createData(RngUtils.create(seed.getSeed()), 1, iter, paramsList, yList);
// Remove the timing of the function call by creating a dummy function
final FakeGradientFunction fgf = new FakeGradientFunction(blockWidth, nparams);
final Gradient1Function func;
if (precomputed) {
final double[] b = SimpleArrayUtils.newArray(fgf.size(), 0.1, 1.3);
func = OffsetGradient1Function.wrapGradient1Function(fgf, b);
} else {
func = fgf;
}
final FastLog fastLog = type == Type.FAST_LOG_MLE ? getFastLog() : null;
final IntArrayFormatSupplier msgA = new IntArrayFormatSupplier("A [%d]", 1);
final IntArrayFormatSupplier msgB = new IntArrayFormatSupplier("B [%d]", 1);
for (int i = 0; i < paramsList.size(); i++) {
final LvmGradientProcedure p1 = createProcedure(type, yList.get(i), func, fastLog);
p1.gradient(paramsList.get(i));
p1.gradient(paramsList.get(i));
final LvmGradientProcedure p2 = LvmGradientProcedureUtils.create(yList.get(i), func, type, fastLog);
p2.gradient(paramsList.get(i));
p2.gradient(paramsList.get(i));
// Check they are the same
Assertions.assertArrayEquals(p1.getAlphaLinear(), p2.getAlphaLinear(), msgA.set(0, i));
Assertions.assertArrayEquals(p1.beta, p2.beta, msgB.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 LvmGradientProcedure p1 = createProcedure(type, yList.get(i), func, fastLog);
for (int j = loops; j-- > 0; ) {
p1.gradient(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 LvmGradientProcedure p2 = LvmGradientProcedureUtils.create(yList.get(i), func, type, fastLog);
for (int j = loops; j-- > 0; ) {
p2.gradient(paramsList.get(k++ % iter));
}
}
}
};
final long time2 = t2.getTime();
logger.log(TestLogUtils.getTimingRecord(new TimingResult(() -> String.format("%s, Precomputed=%b : Standard", type, precomputed), time1), new TimingResult(() -> String.format("Unrolled %d", nparams), time2)));
}
Aggregations