[gnumeric] R.QHYPER: Fix problem at left edge of support.
- From: Morten Welinder <mortenw src gnome org>
- To: commits-list gnome org
- Cc:
- Subject: [gnumeric] R.QHYPER: Fix problem at left edge of support.
- Date: Thu, 3 Dec 2015 15:39:18 +0000 (UTC)
commit 0746e83fb186bb09ef14ac467a5d1a637becc25b
Author: Morten Welinder <terra gnome org>
Date: Thu Dec 3 10:38:14 2015 -0500
R.QHYPER: Fix problem at left edge of support.
ChangeLog | 5 +++
NEWS | 1 +
src/sf-dpq.c | 24 +++++++++++--
src/sstest.c | 103 ++++++++++++++++++++++++++++++++++++++--------------------
4 files changed, 94 insertions(+), 39 deletions(-)
---
diff --git a/ChangeLog b/ChangeLog
index 507d828..3b5ed07 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,8 @@
+2015-12-03 Morten Welinder <terra gnome org>
+
+ * src/sf-dpq.c (discpfuncinverter): Fix problem at left edge of
+ support. Fixes R.QHYPER(0.1,3,99,13)
+
2015-10-19 Morten Welinder <terra gnome org>
* src/func-builtin.c (gnumeric_table): Make sure to invalidate
diff --git a/NEWS b/NEWS
index 63b3e5f..c738b86 100644
--- a/NEWS
+++ b/NEWS
@@ -6,6 +6,7 @@ Morten:
* Add R.DGUMBEL, R.PGUMBEL, R.QGUMBEL.
* Fix conditional format problem. [#750271]
* Test also RANDWEIBULL and RANDLOGNORM.
+ * Fix problem with R.QHYPER
--------------------------------------------------------------------------
Gnumeric 1.12.24
diff --git a/src/sf-dpq.c b/src/sf-dpq.c
index 6bf9206..adf3f89 100644
--- a/src/sf-dpq.c
+++ b/src/sf-dpq.c
@@ -221,6 +221,7 @@ discpfuncinverter (gnm_float p, const gnm_float shape[],
{
gboolean have_xlow = gnm_finite (xlow);
gboolean have_xhigh = gnm_finite (xhigh);
+ gboolean check_left = TRUE;
gnm_float step;
int i;
@@ -253,16 +254,31 @@ discpfuncinverter (gnm_float p, const gnm_float shape[],
g_printerr ("x=%.20g e=%.20g\n", x0, ex0);
#endif
if (!lower_tail) ex0 = -ex0;
- if (ex0 <= 0)
- xlow = x0, have_xlow = TRUE;
- if (ex0 >= 0)
+
+ if (ex0 == 0)
+ return x0;
+ else if (ex0 < 0)
+ xlow = x0, have_xlow = TRUE, check_left = FALSE;
+ else if (ex0 > 0)
xhigh = x0, have_xhigh = TRUE, step = -gnm_abs (step);
if (i > 1 && have_xlow && have_xhigh) {
gnm_float xmid = gnm_floor ((xlow + xhigh) / 2);
if (xmid - xlow < 0.5 ||
- xmid - xlow < gnm_abs (xlow) * GNM_EPSILON)
+ xmid - xlow < gnm_abs (xlow) * GNM_EPSILON) {
+ if (check_left) {
+ /*
+ * The lower edge of the support might
+ * have a probability higher than what
+ * we are looking for.
+ */
+ gnm_float e = pfunc (xlow, shape, lower_tail, log_p) - p;
+ if (!lower_tail) e = -e;
+ if (e >= 0)
+ return xhigh = xlow;
+ }
return xhigh;
+ }
x0 = xmid;
} else {
gnm_float x1 = x0 + step;
diff --git a/src/sstest.c b/src/sstest.c
index ddd9ed5..fd8b198 100644
--- a/src/sstest.c
+++ b/src/sstest.c
@@ -563,6 +563,22 @@ rand_fractile_test (gnm_float const *vals, int N, int nf,
g_printerr (".\n");
}
+ for (i = 1; i < nf - 1; i++) {
+ if (!(fractiles[i] <= fractiles[i + 1])) {
+ g_printerr ("Severe fractile ordering problem.\n");
+ return FALSE;
+ }
+
+ if (probs && !(probs[i] <= probs[i + 1])) {
+ g_printerr ("Severe cumulative probabilities ordering problem.\n");
+ return FALSE;
+ }
+ }
+ if (probs && (probs[1] < 0 || probs[nf - 1] > 1)) {
+ g_printerr ("Severe cumulative probabilities range problem.\n");
+ return FALSE;
+ }
+
for (i = 0; i <= nf; i++)
fractilecount[i] = 0;
@@ -1002,6 +1018,8 @@ test_random_randnorm (int N)
char *expr;
gnm_float T;
int i;
+ gnm_float fractiles[10];
+ const int nf = G_N_ELEMENTS (fractiles);
expr = g_strdup_printf ("=RANDNORM(%.10" GNM_FORMAT_g ",%.10" GNM_FORMAT_g ")",
mean_target, var_target);
@@ -1029,6 +1047,12 @@ test_random_randnorm (int N)
ok = FALSE;
}
+ /* Fractile test */
+ for (i = 1; i < nf; i++)
+ fractiles[i] = qnorm (i / (double)nf, mean_target, gnm_sqrt (var_target), TRUE, FALSE);
+ if (!rand_fractile_test (vals, N, nf, fractiles, NULL))
+ ok = FALSE;
+
if (adtest < 0.05) {
g_printerr ("Anderson Darling Test rejected [%.10" GNM_FORMAT_g "]\n", adtest);
ok = FALSE;
@@ -1816,6 +1840,8 @@ test_random_randhyperg (int N)
gnm_float param_nb = gnm_floor (1 / (0.01 + gnm_pow (random_01 (), 4)));
gnm_float s = param_nr + param_nb;
gnm_float param_n = gnm_floor (random_01 () * (s + 1));
+ param_nr = 3, param_nb = 99, param_n = 13;
+
gnm_float mean_target = param_n * param_nr / s;
gnm_float var_target = s > 1
? mean_target * (param_nb / s) * (s - param_n) / (s - 1)
@@ -2361,47 +2387,54 @@ test_random (void)
const char *test_name = "test_random";
const int N = 20000;
const int High_N = 200000;
+ const char *single = g_getenv ("SSTEST_RANDOM");
mark_test_start (test_name);
- test_random_rand (N);
- test_random_randuniform (N);
- test_random_randnorm (High_N);
- test_random_randsnorm (High_N);
- test_random_randexp (N);
- test_random_randgamma (N);
- test_random_randbeta (N);
- test_random_randtdist (N);
- test_random_randfdist (N);
- test_random_randchisq (N);
- test_random_randcauchy (N);
-
- test_random_randbernoulli (N);
- test_random_randdiscrete (N);
- test_random_randbinom (N);
- test_random_randnegbinom (High_N);
- test_random_randhyperg (High_N);
- test_random_randbetween (N);
- test_random_randpoisson (High_N);
- test_random_randgeom (High_N);
- test_random_randlog (N);
- test_random_randweibull (N);
- test_random_randlognorm (N);
-
+#define CHECK1(NAME,C) \
+ do { if (!single || strcmp(single,#NAME) == 0) test_random_ ## NAME (C); } while (0)
+
+ /* Continuous */
+ CHECK1 (rand, N);
+ CHECK1 (randuniform, N);
+ CHECK1 (randbeta, N);
+ CHECK1 (randcauchy, N);
+ CHECK1 (randchisq, N);
+ CHECK1 (randexp, N);
+ CHECK1 (randfdist, N);
+ CHECK1 (randgamma, N);
+ CHECK1 (randlog, N);
+ CHECK1 (randlognorm, N);
+ CHECK1 (randnorm, High_N);
+ CHECK1 (randsnorm, High_N);
+ CHECK1 (randtdist, N);
+ CHECK1 (randweibull, N);
#if 0
- test_random_randexppow (N);
- test_random_randnormtail (N);
- test_random_randgumbel (N);
- test_random_randlandau (N);
- test_random_randlaplace (N);
- test_random_randlevy (N);
- test_random_randlogistic (N);
- test_random_randpareto (N);
- test_random_randrayleigh (N);
- test_random_randrayleightail (N);
- test_random_randstdist (N);
+ CHECK1 (randexppow, N);
+ CHECK1 (randgumbel, N);
+ CHECK1 (randlandau, N);
+ CHECK1 (randlaplace, N);
+ CHECK1 (randlevy, N);
+ CHECK1 (randlogistic, N);
+ CHECK1 (randnormtail, N);
+ CHECK1 (randpareto, N);
+ CHECK1 (randrayleigh, N);
+ CHECK1 (randrayleightail, N);
+ CHECK1 (randstdist, N);
#endif
+ /* Discrete */
+ CHECK1 (randbernoulli, N);
+ CHECK1 (randbetween, N);
+ CHECK1 (randbinom, N);
+ CHECK1 (randdiscrete, N);
+ CHECK1 (randgeom, High_N);
+ CHECK1 (randhyperg, High_N);
+ CHECK1 (randnegbinom, High_N);
+ CHECK1 (randpoisson, High_N);
+
+#undef CHECK1
+
mark_test_end (test_name);
}
[
Date Prev][
Date Next] [
Thread Prev][
Thread Next]
[
Thread Index]
[
Date Index]
[
Author Index]