use of org.apache.commons.math3.fitting.leastsquares.MultivariateJacobianFunction in project GDSC-SMLM by aherbert.
the class TrackPopulationAnalysisTest method canComputeFbmDiffusionFunction.
@Test
void canComputeFbmDiffusionFunction() {
final int size = 10;
final double delta = 1e-6;
final DoubleDoubleBiPredicate test = TestHelper.doublesAreClose(1e-5);
for (final double t : new double[] { 0.8, 1, 1.2 }) {
final MultivariateJacobianFunction f = new FbmDiffusionFunction(size, t);
for (final double d : new double[] { 0.8, 0.9, 1, 1.1, 1.2 }) {
for (final double s : new double[] { 2, 20 }) {
for (final double alpha : new double[] { 1e-2, 0.8, 0.9, 1, 1.1, 1.2, 2.0 }) {
// Check the value and Jacobian
final RealVector point = new ArrayRealVector(new double[] { d, s, alpha }, false);
final Pair<RealVector, RealMatrix> p = f.value(point);
final double[] value = p.getFirst().toArray();
Assertions.assertEquals(size, value.length);
final double a1 = alpha + 1;
final double a2 = alpha + 2;
final double ta = Math.pow(t, alpha);
for (int n = 1; n <= size; n++) {
// MSD = [4Dt^a / (a+2)(a+1)] * [(n+1)^(a+2) + (n-1)^(a+2) - 2n^(a+2)]
// - [8Dt^a / (a+2)(a+1)] + 4s^2
// assume t=1.
final double msd = (4 * d * ta / (a2 * a1)) * (Math.pow(n + 1, a2) + Math.pow(n - 1, a2) - 2 * Math.pow(n, a2)) - 8 * d * ta / (a2 * a1) + 4 * s * s;
TestAssertions.assertTest(msd, value[n - 1], test, "value");
}
// Columns of the Jacobian
final double[] dfda1 = p.getSecond().getColumn(0);
final double[] dfdb1 = p.getSecond().getColumn(1);
final double[] dfdc1 = p.getSecond().getColumn(2);
point.setEntry(0, d - delta);
RealVector v1 = f.value(point).getFirst();
point.setEntry(0, d + delta);
RealVector v2 = f.value(point).getFirst();
final double[] dfda = v2.subtract(v1).mapDivide(2 * delta).toArray();
point.setEntry(0, d);
point.setEntry(1, s - delta);
v1 = f.value(point).getFirst();
point.setEntry(1, s + delta);
v2 = f.value(point).getFirst();
final double[] dfdb = v2.subtract(v1).mapDivide(2 * delta).toArray();
point.setEntry(1, s);
point.setEntry(2, alpha - delta);
v1 = f.value(point).getFirst();
point.setEntry(2, alpha + delta);
v2 = f.value(point).getFirst();
final double[] dfdc = v2.subtract(v1).mapDivide(2 * delta).toArray();
// Element-by-element relative error
TestAssertions.assertArrayTest(dfda, dfda1, test, "jacobian dfda");
TestAssertions.assertArrayTest(dfdb, dfdb1, test, "jacobian dfdb");
TestAssertions.assertArrayTest(dfdc, dfdc1, test, "jacobian dfdc");
}
}
}
}
}
use of org.apache.commons.math3.fitting.leastsquares.MultivariateJacobianFunction in project GDSC-SMLM by aherbert.
the class TrackPopulationAnalysisTest method canComputeBrownianDiffusionFunction1.
@Test
void canComputeBrownianDiffusionFunction1() {
final int size = 10;
final double delta = 1e-6;
final DoubleDoubleBiPredicate test = TestHelper.doublesAreClose(1e-5);
for (final double t : new double[] { 0.8, 1, 1.2 }) {
final MultivariateJacobianFunction f = new BrownianDiffusionFunction(size, t);
for (final double d : new double[] { 0.8, 0.9, 1, 1.1, 1.2 }) {
for (final double s : new double[] { 2, 20 }) {
// Check the value and Jacobian
final RealVector point = new ArrayRealVector(new double[] { d, s }, false);
final Pair<RealVector, RealMatrix> p = f.value(point);
final double[] value = p.getFirst().toArray();
Assertions.assertEquals(size, value.length);
for (int n = 1; n <= size; n++) {
// MSD = 4Dt * (n - 1/3) + 4s^2
final double msd = 4 * d * t * (n - 1.0 / 3) + 4 * s * s;
TestAssertions.assertTest(msd, value[n - 1], test, "value");
}
// Columns of the Jacobian
final double[] dfda1 = p.getSecond().getColumn(0);
final double[] dfdb1 = p.getSecond().getColumn(1);
point.setEntry(0, d - delta);
RealVector v1 = f.value(point).getFirst();
point.setEntry(0, d + delta);
RealVector v2 = f.value(point).getFirst();
final double[] dfda = v2.subtract(v1).mapDivide(2 * delta).toArray();
point.setEntry(0, d);
point.setEntry(1, s - delta);
v1 = f.value(point).getFirst();
point.setEntry(1, s + delta);
v2 = f.value(point).getFirst();
final double[] dfdb = v2.subtract(v1).mapDivide(2 * delta).toArray();
// Element-by-element relative error
TestAssertions.assertArrayTest(dfda, dfda1, test, "jacobian dfda");
TestAssertions.assertArrayTest(dfdb, dfdb1, test, "jacobian dfdb");
}
}
}
}
Aggregations