summaryrefslogtreecommitdiffstats
path: root/test/unit/regress.c
diff options
context:
space:
mode:
Diffstat (limited to '')
-rw-r--r--test/unit/regress.c119
1 files changed, 119 insertions, 0 deletions
diff --git a/test/unit/regress.c b/test/unit/regress.c
new file mode 100644
index 0000000..f47d1c4
--- /dev/null
+++ b/test/unit/regress.c
@@ -0,0 +1,119 @@
+/*
+ **********************************************************************
+ * Copyright (C) Miroslav Lichvar 2017
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of version 2 of the GNU General Public License as
+ * published by the Free Software Foundation.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License along
+ * with this program; if not, write to the Free Software Foundation, Inc.,
+ * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ *
+ **********************************************************************
+ */
+#include <regress.c>
+#include "test.h"
+
+#define POINTS 64
+
+void
+test_unit(void)
+{
+ double x[POINTS], x2[POINTS], y[POINTS], w[POINTS];
+ double b0, b1, b2, s2, sb0, sb1, slope, slope2, intercept, sd, median;
+ double xrange, yrange, wrange, x2range;
+ int i, j, n, m, c1, c2, c3, runs, best_start, dof;
+
+ for (n = 3; n <= POINTS; n++) {
+ for (i = 0; i < 200; i++) {
+ slope = TST_GetRandomDouble(-0.1, 0.1);
+ intercept = TST_GetRandomDouble(-1.0, 1.0);
+ sd = TST_GetRandomDouble(1e-6, 1e-4);
+ slope2 = (random() % 2 ? 1 : -1) * TST_GetRandomDouble(0.1, 0.5);
+
+ DEBUG_LOG("iteration %d n=%d intercept=%e slope=%e sd=%e",
+ i, n, intercept, slope, sd);
+
+ for (j = 0; j < n; j++) {
+ x[j] = -j;
+ y[j] = intercept + slope * x[j] + (j % 2 ? 1 : -1) * TST_GetRandomDouble(1e-6, sd);
+ w[j] = TST_GetRandomDouble(1.0, 2.0);
+ x2[j] = (y[j] - intercept - slope * x[j]) / slope2;
+ }
+
+ RGR_WeightedRegression(x, y, w, n, &b0, &b1, &s2, &sb0, &sb1);
+ DEBUG_LOG("WR b0=%e b1=%e s2=%e sb0=%e sb1=%e", b0, b1, s2, sb0, sb1);
+ TEST_CHECK(fabs(b0 - intercept) < sd + 1e-3);
+ TEST_CHECK(fabs(b1 - slope) < sd);
+
+ if (RGR_FindBestRegression(x, y, w, n, 0, 3, &b0, &b1, &s2, &sb0, &sb1,
+ &best_start, &runs, &dof)) {
+ DEBUG_LOG("BR b0=%e b1=%e s2=%e sb0=%e sb1=%e runs=%d bs=%d dof=%d",
+ b0, b1, s2, sb0, sb1, runs, best_start, dof);
+
+ TEST_CHECK(fabs(b0 - intercept) < sd + 1e-3);
+ TEST_CHECK(fabs(b1 - slope) < sd);
+ }
+
+ if (RGR_MultipleRegress(x, x2, y, n, &b2)) {
+ DEBUG_LOG("MR b2=%e", b2);
+ TEST_CHECK(fabs(b2 - slope2) < 1e-6);
+ }
+
+ for (j = 0; j < n / 7; j++)
+ y[random() % n] += 100 * sd;
+
+ if (RGR_FindBestRobustRegression(x, y, n, 1e-8, &b0, &b1, &runs, &best_start)) {
+ DEBUG_LOG("BRR b0=%e b1=%e runs=%d bs=%d", b0, b1, runs, best_start);
+
+ TEST_CHECK(fabs(b0 - intercept) < sd + 1e-2);
+ TEST_CHECK(fabs(b1 - slope) < 5.0 * sd);
+ }
+
+ for (j = 0; j < n; j++)
+ x[j] = random() % 4 * TST_GetRandomDouble(-1000, 1000);
+
+ median = RGR_FindMedian(x, n);
+
+ for (j = c1 = c2 = c3 = 0; j < n; j++) {
+ if (x[j] < median)
+ c1++;
+ if (x[j] > median)
+ c3++;
+ else
+ c2++;
+ }
+
+ TEST_CHECK(c1 + c2 >= c3 && c1 <= c2 + c3);
+
+ xrange = TST_GetRandomDouble(1e-6, pow(10.0, random() % 10));
+ yrange = random() % 3 * TST_GetRandomDouble(0.0, pow(10.0, random() % 10));
+ wrange = random() % 3 * TST_GetRandomDouble(0.0, pow(10.0, random() % 10));
+ x2range = random() % 3 * TST_GetRandomDouble(0.0, pow(10.0, random() % 10));
+ m = random() % n;
+
+ for (j = 0; j < n; j++) {
+ x[j] = (j ? x[j - 1] : 0.0) + TST_GetRandomDouble(1e-6, xrange);
+ y[j] = TST_GetRandomDouble(-yrange, yrange);
+ w[j] = 1.0 + TST_GetRandomDouble(0.0, wrange);
+ x2[j] = TST_GetRandomDouble(-x2range, x2range);
+ }
+
+ RGR_WeightedRegression(x, y, w, n, &b0, &b1, &s2, &sb0, &sb1);
+
+ if (RGR_FindBestRegression(x + m, y + m, w, n - m, m, 3, &b0, &b1, &s2, &sb0, &sb1,
+ &best_start, &runs, &dof))
+ ;
+ if (RGR_MultipleRegress(x, x2, y, n, &b2))
+ ;
+ if (RGR_FindBestRobustRegression(x, y, n, 1e-8, &b0, &b1, &runs, &best_start))
+ ;
+ }
+ }
+}