Changeset: bfb1f607de02 for MonetDB
URL: http://dev.monetdb.org/hg/MonetDB?cmd=changeset;node=bfb1f607de02
Modified Files:
clients/Tests/exports.stable.out
gdk/gdk_aggr.c
gdk/gdk_calc.h
monetdb5/modules/kernel/aggr.c
monetdb5/modules/kernel/aggr.mal
monetdb5/modules/kernel/aggr.mal.sh
monetdb5/modules/kernel/algebra.mx
Branch: default
Log Message:
Implemented standard deviation (sample and population) in single scan.
This implementation is currently not yet used by SQL.
See bug 3178.
diffs (truncated from 883 to 300 lines):
diff --git a/clients/Tests/exports.stable.out b/clients/Tests/exports.stable.out
--- a/clients/Tests/exports.stable.out
+++ b/clients/Tests/exports.stable.out
@@ -90,6 +90,8 @@ BAT *BATcalcorcst(BAT *b, const ValRecor
BAT *BATcalcrsh(BAT *b1, BAT *b2, BAT *s, int abort_on_error);
BAT *BATcalcrshcst(BAT *b, const ValRecord *v, BAT *s, int abort_on_error);
BAT *BATcalcsign(BAT *b, BAT *s);
+dbl BATcalcstdev_population(dbl *avgp, BAT *b);
+dbl BATcalcstdev_sample(dbl *avgp, BAT *b);
BAT *BATcalcsub(BAT *b1, BAT *b2, BAT *s, int tp, int abort_on_error);
BAT *BATcalcsubcst(BAT *b, const ValRecord *v, BAT *s, int tp, int
abort_on_error);
BAT *BATcalcxor(BAT *b1, BAT *b2, BAT *s);
@@ -124,6 +126,8 @@ BAT *BATgroupmedian(BAT *b, BAT *g, BAT
BAT *BATgroupmin(BAT *b, BAT *g, BAT *e, BAT *s, int tp, int skip_nils, int
abort_on_error);
BAT *BATgroupprod(BAT *b, BAT *g, BAT *e, BAT *s, int tp, int skip_nils, int
abort_on_error);
BAT *BATgroupsize(BAT *b, BAT *g, BAT *e, BAT *s, int tp, int skip_nils, int
abort_on_error);
+BAT *BATgroupstdev_population(BAT *b, BAT *g, BAT *e, BAT *s, int tp, int
skip_nils, int abort_on_error);
+BAT *BATgroupstdev_sample(BAT *b, BAT *g, BAT *e, BAT *s, int tp, int
skip_nils, int abort_on_error);
BAT *BATgroupsum(BAT *b, BAT *g, BAT *e, BAT *s, int tp, int skip_nils, int
abort_on_error);
BUN BATgrows(BAT *b);
BAT *BAThash(BAT *b, BUN masksize);
@@ -621,6 +625,10 @@ str AGGRprod3_lng(bat *retval, bat *bid,
str AGGRprod3_sht(bat *retval, bat *bid, bat *gid, bat *eid);
str AGGRprod3_wrd(bat *retval, bat *bid, bat *gid, bat *eid);
str AGGRsize2(bat *retval, bat *bid, bat *eid);
+str AGGRstdev2_dbl(bat *retval, bat *bid, bat *eid);
+str AGGRstdev3_dbl(bat *retval, bat *bid, bat *gid, bat *eid);
+str AGGRstdevp2_dbl(bat *retval, bat *bid, bat *eid);
+str AGGRstdevp3_dbl(bat *retval, bat *bid, bat *gid, bat *eid);
str AGGRsubavg_dbl(bat *retval, bat *bid, bat *gid, bat *eid, int *skip_nils,
int *abort_on_error);
str AGGRsubavgcand_dbl(bat *retval, bat *bid, bat *gid, bat *eid, bat *sid,
int *skip_nils, int *abort_on_error);
str AGGRsubcount(bat *retval, bat *bid, bat *gid, bat *eid, int *skip_nils);
@@ -649,6 +657,10 @@ str AGGRsubprodcand_int(bat *retval, bat
str AGGRsubprodcand_lng(bat *retval, bat *bid, bat *gid, bat *eid, bat *sid,
int *skip_nils, int *abort_on_error);
str AGGRsubprodcand_sht(bat *retval, bat *bid, bat *gid, bat *eid, bat *sid,
int *skip_nils, int *abort_on_error);
str AGGRsubprodcand_wrd(bat *retval, bat *bid, bat *gid, bat *eid, bat *sid,
int *skip_nils, int *abort_on_error);
+str AGGRsubstdev_dbl(bat *retval, bat *bid, bat *gid, bat *eid, int
*skip_nils, int *abort_on_error);
+str AGGRsubstdevcand_dbl(bat *retval, bat *bid, bat *gid, bat *eid, bat *sid,
int *skip_nils, int *abort_on_error);
+str AGGRsubstdevp_dbl(bat *retval, bat *bid, bat *gid, bat *eid, int
*skip_nils, int *abort_on_error);
+str AGGRsubstdevpcand_dbl(bat *retval, bat *bid, bat *gid, bat *eid, bat *sid,
int *skip_nils, int *abort_on_error);
str AGGRsubsum_bte(bat *retval, bat *bid, bat *gid, bat *eid, int *skip_nils,
int *abort_on_error);
str AGGRsubsum_dbl(bat *retval, bat *bid, bat *gid, bat *eid, int *skip_nils,
int *abort_on_error);
str AGGRsubsum_flt(bat *retval, bat *bid, bat *gid, bat *eid, int *skip_nils,
int *abort_on_error);
@@ -792,13 +804,8 @@ str ALGslice_wrd(int *ret, bat *bid, wrd
str ALGsplit(int *result, int *bid);
str ALGssort(int *result, int *bid);
str ALGssort_rev(int *result, int *bid);
-str ALGstdev_bte(dbl *res, int *bid);
-str ALGstdev_dbl(dbl *res, int *bid);
-str ALGstdev_flt(dbl *res, int *bid);
-str ALGstdev_int(dbl *res, int *bid);
-str ALGstdev_lng(dbl *res, int *bid);
-str ALGstdev_sht(dbl *res, int *bid);
-str ALGstdev_wrd(dbl *res, int *bid);
+str ALGstdev(dbl *res, int *bid);
+str ALGstdevp(dbl *res, int *bid);
str ALGsubselect1(bat *result, bat *bid, const void *low, const void *high,
const bit *li, const bit *hi, const bit *anti);
str ALGsubselect2(bat *result, bat *bid, bat *sid, const void *low, const void
*high, const bit *li, const bit *hi, const bit *anti);
str ALGsubslice_wrd(int *ret, bat *bid, wrd *start, wrd *end);
diff --git a/gdk/gdk_aggr.c b/gdk/gdk_aggr.c
--- a/gdk/gdk_aggr.c
+++ b/gdk/gdk_aggr.c
@@ -2114,3 +2114,294 @@ BATgroupmedian(BAT *b, BAT *g, BAT *e, B
BBPunfix(bn->batCacheid);
return NULL;
}
+
+/* ---------------------------------------------------------------------- */
+/* standard deviation (both biases and non-biased) */
+
+#define AGGR_STDEV(TYPE) \
+ do { \
+ TYPE x; \
+ for (i = 0; i < cnt; i++) { \
+ x = ((const TYPE *) values)[i]; \
+ if (x == TYPE##_nil) \
+ continue; \
+ n++; \
+ delta = (dbl) x - mean; \
+ mean += delta / n; \
+ m2 += delta * ((dbl) x - mean); \
+ } \
+ } while (0)
+
+static dbl
+calcstdev(dbl *avgp, const void *values, BUN cnt, int tp, int issample)
+{
+ BUN n = 0, i;
+ dbl mean = 0;
+ dbl m2 = 0;
+ dbl delta;
+
+ assert(issample == 0 || issample == 1);
+
+ switch (ATOMstorage(tp)) {
+ case TYPE_bte:
+ AGGR_STDEV(bte);
+ break;
+ case TYPE_sht:
+ AGGR_STDEV(sht);
+ break;
+ case TYPE_int:
+ AGGR_STDEV(int);
+ break;
+ case TYPE_lng:
+ AGGR_STDEV(lng);
+ break;
+ case TYPE_flt:
+ AGGR_STDEV(flt);
+ break;
+ case TYPE_dbl:
+ AGGR_STDEV(dbl);
+ break;
+ default:
+ return dbl_nil;
+ }
+ if (n <= (BUN) issample) {
+ if (avgp)
+ *avgp = dbl_nil;
+ return dbl_nil;
+ }
+ if (avgp)
+ *avgp = mean;
+ return sqrt(m2 / (n - issample));
+}
+
+dbl
+BATcalcstdev_population(dbl *avgp, BAT *b)
+{
+ return calcstdev(avgp, (const void *) Tloc(b, BUNfirst(b)),
+ BATcount(b), b->ttype, 0);
+}
+
+dbl
+BATcalcstdev_sample(dbl *avgp, BAT *b)
+{
+ return calcstdev(avgp, (const void *) Tloc(b, BUNfirst(b)),
+ BATcount(b), b->ttype, 1);
+}
+
+#undef AGGR_STDEV
+#define AGGR_STDEV(TYPE) \
+ do { \
+ const TYPE *vals = (const TYPE *) Tloc(b, BUNfirst(b)); \
+ for (i = start; i < end; i++, vals++) { \
+ if (cand) { \
+ if (i < *cand - b->hseqbase) { \
+ if (gids) \
+ gids++; \
+ continue; \
+ } \
+ assert(i == *cand - b->hseqbase); \
+ if (++cand == candend) \
+ end = i + 1; \
+ } \
+ if (gids == NULL || \
+ (*gids >= min && *gids <= max)) { \
+ gid = gids ? *gids - min : (oid) i; \
+ if (*vals == TYPE##_nil) { \
+ if (!skip_nils) \
+ cnts[gid] = BUN_NONE; \
+ } else if (cnts[gid] != BUN_NONE) { \
+ cnts[gid]++; \
+ delta[gid] = (dbl) *vals - mean[gid]; \
+ mean[gid] += delta[gid] / cnts[gid]; \
+ m2[gid] += delta[gid] * ((dbl) *vals -
mean[gid]); \
+ } \
+ } \
+ if (gids) \
+ gids++; \
+ } \
+ for (i = 0; i < ngrp; i++) { \
+ if (cnts[i] == 0 || cnts[i] == BUN_NONE) { \
+ dbls[i] = dbl_nil; \
+ mean[i] = dbl_nil; \
+ nils++; \
+ } else if (cnts[i] == 1) { \
+ dbls[i] = 0; \
+ } else { \
+ dbls[i] = m2[i] / (cnts[i] - issample); \
+ } \
+ } \
+ } while (0)
+
+/* Calculate group standard deviation (population (i.e. biased) or
+ * sample (i.e. non-biased)) with optional candidates list.
+ *
+ * Note that this helper function is prepared to return two BATs: one
+ * (as return value) with the standard deviation per group, and one
+ * (as return argument) with the average per group. This isn't
+ * currently used since it doesn't fit into the mold of grouped
+ * aggregates. */
+static BAT *
+dogroupstdev(BAT **avgb, BAT *b, BAT *g, BAT *e, BAT *s, int tp,
+ int skip_nils, int issample, const char *func)
+{
+ const oid *gids;
+ oid gid;
+ oid min, max;
+ BUN i, ngrp;
+ BUN nils = 0;
+ BUN *cnts = NULL;
+ dbl *dbls, *mean, *delta, *m2;
+ BAT *bn = NULL;
+ BUN start, end, cnt;
+ const oid *cand = NULL, *candend = NULL;
+ const char *err;
+
+ assert(tp == TYPE_dbl);
+ (void) tp; /* compatibility (with other BATgroup*
+ * functions) argument */
+
+ if ((err = BATgroupaggrinit(b, g, e, s, &min, &max, &ngrp, &start, &end,
+ &cnt, &cand, &candend)) != NULL) {
+ GDKerror("%s: %s\n", func, err);
+ return NULL;
+ }
+ if (g == NULL) {
+ GDKerror("%s: b and g must be aligned\n", func);
+ return NULL;
+ }
+
+ if (BATcount(b) == 0 || ngrp == 0) {
+ /* trivial: no products, so return bat aligned with g
+ * with nil in the tail */
+ bn = BATconstant(TYPE_dbl, &dbl_nil, ngrp);
+ BATseqbase(bn, ngrp == 0 ? 0 : min);
+ return bn;
+ }
+
+ if ((e == NULL ||
+ (BATcount(e) == BATcount(b) && e->hseqbase == b->hseqbase)) &&
+ (BATtdense(g) || (g->tkey && g->T->nonil))) {
+ /* trivial: singleton groups, so all results are equal
+ * to zero (population) or nil (sample) */
+ dbl v = issample ? dbl_nil : 0;
+ bn = BATconstant(TYPE_dbl, &v, ngrp);
+ BATseqbase(bn, ngrp == 0 ? 0 : min);
+ return bn;
+ }
+
+ delta = GDKmalloc(ngrp * sizeof(dbl));
+ m2 = GDKmalloc(ngrp * sizeof(dbl));
+ cnts = GDKzalloc(ngrp * sizeof(BUN));
+ if (avgb) {
+ if ((*avgb = BATnew(TYPE_void, TYPE_dbl, ngrp)) == NULL) {
+ mean = NULL;
+ goto alloc_fail;
+ }
+ mean = (dbl *) Tloc(*avgb, BUNfirst(*avgb));
+ } else {
+ mean = GDKmalloc(ngrp * sizeof(dbl));
+ }
+ if (mean == NULL || delta == NULL || m2 == NULL || cnts == NULL)
+ goto alloc_fail;
+
+ bn = BATnew(TYPE_void, TYPE_dbl, ngrp);
+ if (bn == NULL)
+ goto alloc_fail;
+ dbls = (dbl *) Tloc(bn, BUNfirst(bn));
+
+ for (i = 0; i < ngrp; i++) {
+ mean[i] = 0;
+ delta[i] = 0;
+ m2[i] = 0;
+ }
+
+ if (BATtdense(g))
+ gids = NULL;
+ else
+ gids = (const oid *) Tloc(g, BUNfirst(g) + start);
+
+ switch (ATOMstorage(b->ttype)) {
+ case TYPE_bte:
+ AGGR_STDEV(bte);
+ break;
+ case TYPE_sht:
+ AGGR_STDEV(sht);
+ break;
+ case TYPE_int:
+ AGGR_STDEV(int);
+ break;
+ case TYPE_lng:
+ AGGR_STDEV(lng);
+ break;
+ case TYPE_flt:
+ AGGR_STDEV(flt);
+ break;
+ case TYPE_dbl:
+ AGGR_STDEV(dbl);
+ break;
+ default:
+ if (avgb)
+ BBPreclaim(*avgb);
+ else
+ GDKfree(mean);
+ GDKfree(delta);
+ GDKfree(m2);
+ GDKfree(cnts);
+ BBPunfix(bn->batCacheid);
_______________________________________________
checkin-list mailing list
[email protected]
http://mail.monetdb.org/mailman/listinfo/checkin-list
