Skip to content

Commit e1c633c

Browse files
AllinCottrellAllinCottrell
committed
try an experiment with aggregate() using matrix data
1 parent 133dcfa commit e1c633c

File tree

4 files changed

+178
-58
lines changed

4 files changed

+178
-58
lines changed

‎lib/src/geneval.c

+7-1
Original file line numberDiff line numberDiff line change
@@ -13604,6 +13604,7 @@ static NODE *eval_3args_func (NODE *l, NODE *m, NODE *r,
1360413604
int *ylist = NULL;
1360513605
gretl_matrix *xm = NULL;
1360613606
gretl_matrix *ym = NULL;
13607+
int matrix_version = 0;
1360713608
int data_args = 1;
1360813609
int mat_args = 0;
1360913610

@@ -13641,6 +13642,9 @@ static NODE *eval_3args_func (NODE *l, NODE *m, NODE *r,
1364113642
if (mat_args == 1 && data_args == 2) {
1364213643
/* matrix arguments: can't be mixed with series/list */
1364313644
p->err = E_TYPES;
13645+
} else if (mat_args == 2 && ym->cols == 1) {
13646+
/* experimental */
13647+
matrix_version = 1;
1364413648
} else if (mat_args > 0) {
1364513649
/* convert to dataset (with borrowed Z) */
1364613650
dset = mdset = matrix_dset_plus_lists(xm, ym, &xlist,
@@ -13650,7 +13654,9 @@ static NODE *eval_3args_func (NODE *l, NODE *m, NODE *r,
1365013654
if (!p->err && ylist != NULL) {
1365113655
p->err = aggregate_discrete_check(ylist, dset);
1365213656
}
13653-
if (!p->err) {
13657+
if (!p->err && matrix_version) {
13658+
A = matrix_aggregate(xm, ym, fnname, &p->err);
13659+
} else if (!p->err) {
1365413660
A = aggregate_by(x, y, xlist, ylist, fnname, dset, &p->err);
1365513661
}
1365613662
if (mdset != NULL) {

‎lib/src/genfuncs.c

+165-56
Original file line numberDiff line numberDiff line change
@@ -7139,7 +7139,8 @@ static gretl_matrix *real_aggregate_by (const double *x,
71397139
gretl_matrix **listvals;
71407140
gretl_matrix *yvals;
71417141
int n = sample_size(dset);
7142-
int match, skipnull = 0;
7142+
int match;
7143+
int skipnull = 0;
71437144
int countcol = 1;
71447145
int maxcases;
71457146
int ny, nx, mcols;
@@ -7382,6 +7383,7 @@ static void aggr_add_colnames (gretl_matrix *m,
73827383
*/
73837384

73847385
static fncall *get_user_aggrby_call (const char *s,
7386+
int nparam,
73857387
double *tmp,
73867388
int *err)
73877389
{
@@ -7393,9 +7395,8 @@ static fncall *get_user_aggrby_call (const char *s,
73937395
*err = E_INVARG;
73947396
} else {
73957397
GretlType rt = user_func_get_return_type(uf);
7396-
int np = fn_n_params(uf);
73977398

7398-
if (rt != GRETL_TYPE_DOUBLE || np != 1) {
7399+
if (rt != GRETL_TYPE_DOUBLE || fn_n_params(uf) != nparam) {
73997400
*err = E_INVARG;
74007401
}
74017402
}
@@ -7404,16 +7405,77 @@ static fncall *get_user_aggrby_call (const char *s,
74047405
fc = user_func_get_fncall(uf);
74057406
if (fc == NULL) {
74067407
*err = E_ALLOC;
7407-
} else {
7408+
} else if (tmp != NULL) {
74087409
*err = push_function_arg(fc, NULL, NULL,
74097410
GRETL_TYPE_SERIES,
74107411
tmp);
7412+
} else {
7413+
; /* not ready yet! */
74117414
}
74127415
}
74137416

74147417
return fc;
74157418
}
74167419

7420+
static int get_aggregator (const char *fname,
7421+
double (**dbuiltin) (int, int, const double *),
7422+
int (**ibuiltin) (int, int, const double *))
7423+
{
7424+
int f = function_lookup(fname);
7425+
7426+
switch (f) {
7427+
case F_SUM:
7428+
*dbuiltin = gretl_sum;
7429+
break;
7430+
case F_SUMALL:
7431+
*dbuiltin = series_sum_all;
7432+
break;
7433+
case F_MEAN:
7434+
*dbuiltin = gretl_mean;
7435+
break;
7436+
case F_SD:
7437+
*dbuiltin = gretl_stddev;
7438+
break;
7439+
case F_VCE:
7440+
*dbuiltin = gretl_variance;
7441+
break;
7442+
case F_SST:
7443+
*dbuiltin = gretl_sst;
7444+
break;
7445+
case F_SKEWNESS:
7446+
*dbuiltin = gretl_skewness;
7447+
break;
7448+
case F_KURTOSIS:
7449+
*dbuiltin = gretl_kurtosis;
7450+
break;
7451+
case F_MIN:
7452+
*dbuiltin = gretl_min;
7453+
break;
7454+
case F_MAX:
7455+
*dbuiltin = gretl_max;
7456+
break;
7457+
case F_MEDIAN:
7458+
*dbuiltin = gretl_median;
7459+
break;
7460+
case F_GINI:
7461+
*dbuiltin = gretl_gini;
7462+
break;
7463+
case F_NOBS:
7464+
*dbuiltin = series_get_nobs;
7465+
break;
7466+
case F_ISCONST:
7467+
*ibuiltin = gretl_isconst;
7468+
break;
7469+
case F_ISDUMMY:
7470+
*ibuiltin = gretl_isdummy;
7471+
break;
7472+
default:
7473+
break;
7474+
}
7475+
7476+
return *dbuiltin != NULL || *ibuiltin != NULL;
7477+
}
7478+
74177479
/**
74187480
* aggregate_by:
74197481
* @x: data array.
@@ -7462,57 +7524,7 @@ gretl_matrix *aggregate_by (const double *x,
74627524
}
74637525

74647526
if (!just_count) {
7465-
int f = function_lookup(fname);
7466-
7467-
switch (f) {
7468-
case F_SUM:
7469-
dbuiltin = gretl_sum;
7470-
break;
7471-
case F_SUMALL:
7472-
dbuiltin = series_sum_all;
7473-
break;
7474-
case F_MEAN:
7475-
dbuiltin = gretl_mean;
7476-
break;
7477-
case F_SD:
7478-
dbuiltin = gretl_stddev;
7479-
break;
7480-
case F_VCE:
7481-
dbuiltin = gretl_variance;
7482-
break;
7483-
case F_SST:
7484-
dbuiltin = gretl_sst;
7485-
break;
7486-
case F_SKEWNESS:
7487-
dbuiltin = gretl_skewness;
7488-
break;
7489-
case F_KURTOSIS:
7490-
dbuiltin = gretl_kurtosis;
7491-
break;
7492-
case F_MIN:
7493-
dbuiltin = gretl_min;
7494-
break;
7495-
case F_MAX:
7496-
dbuiltin = gretl_max;
7497-
break;
7498-
case F_MEDIAN:
7499-
dbuiltin = gretl_median;
7500-
break;
7501-
case F_GINI:
7502-
dbuiltin = gretl_gini;
7503-
break;
7504-
case F_NOBS:
7505-
dbuiltin = series_get_nobs;
7506-
break;
7507-
case F_ISCONST:
7508-
ibuiltin = gretl_isconst;
7509-
break;
7510-
case F_ISDUMMY:
7511-
ibuiltin = gretl_isdummy;
7512-
break;
7513-
default:
7514-
break;
7515-
}
7527+
get_aggregator(fname, &dbuiltin, &ibuiltin);
75167528
}
75177529

75187530
n = sample_size(dset);
@@ -7524,7 +7536,7 @@ gretl_matrix *aggregate_by (const double *x,
75247536
if (tmp == NULL) {
75257537
*err = E_ALLOC;
75267538
} else if (dbuiltin == NULL && ibuiltin == NULL) {
7527-
fc = get_user_aggrby_call(fname, tmp, err);
7539+
fc = get_user_aggrby_call(fname, 1, tmp, err);
75287540
}
75297541
}
75307542

@@ -7553,6 +7565,103 @@ gretl_matrix *aggregate_by (const double *x,
75537565
return m;
75547566
}
75557567

7568+
gretl_matrix *matrix_aggregate (const gretl_matrix *X,
7569+
const gretl_matrix *y,
7570+
const char *func,
7571+
int *err)
7572+
{
7573+
double (*dbuiltin) (int, int, const double *) = NULL;
7574+
int (*ibuiltin) (int, int, const double *) = NULL;
7575+
gretl_matrix *Tmp;
7576+
gretl_matrix *M;
7577+
gretl_matrix *ret;
7578+
double *z, *zsave;
7579+
double by = 0;
7580+
double rkj;
7581+
int nx = X->cols;
7582+
int nby = 0;
7583+
int t1, t2, count;
7584+
int i, j, k;
7585+
7586+
if (X == NULL || y == NULL || X->rows != y->rows) {
7587+
*err = E_INVARG;
7588+
} else if (!get_aggregator(func, &dbuiltin, &ibuiltin)) {
7589+
// fc = get_user_aggrby_call(func, 3, NULL, err);
7590+
*err = E_INVARG;
7591+
}
7592+
7593+
if (*err) {
7594+
return NULL;
7595+
}
7596+
7597+
/* allocate space for y ~ X */
7598+
Tmp = gretl_matrix_alloc(X->rows, 1 + X->cols);
7599+
if (Tmp == NULL) {
7600+
*err = E_ALLOC;
7601+
return NULL;
7602+
}
7603+
7604+
/* make M = (y ~ X), sorted by y */
7605+
memcpy(Tmp->val, y->val, X->rows * sizeof(double));
7606+
memcpy(Tmp->val + X->rows, X->val, X->rows * nx * sizeof(double));
7607+
M = gretl_matrix_sort_by_column(Tmp, 0, err);
7608+
gretl_matrix_free(Tmp);
7609+
if (*err) {
7610+
return NULL;
7611+
}
7612+
7613+
/* determine the number of 'by' values */
7614+
nby = 0;
7615+
for (i=0; i<M->rows; i++) {
7616+
by = M->val[i];
7617+
if (i == 0 || by != M->val[i-1]) {
7618+
nby++;
7619+
}
7620+
}
7621+
7622+
/* allocate the return matrix */
7623+
ret = gretl_zero_matrix_new(nby, 2 + nx);
7624+
7625+
z = M->val + M->rows;
7626+
t2 = t1 = 0;
7627+
k = 0;
7628+
7629+
for (i=0; i<M->rows; i++) {
7630+
by = M->val[i];
7631+
if (i > 0) {
7632+
if (by > M->val[i-1] || i == M->rows - 1) {
7633+
if (i == M->rows - 1) {
7634+
t2++;
7635+
}
7636+
/* calculate */
7637+
gretl_matrix_set(ret, k, 0, M->val[t1]);
7638+
count = t2 - t1 + 1;
7639+
gretl_matrix_set(ret, k, 1, count);
7640+
zsave = z;
7641+
for (j=0; j<nx; j++) {
7642+
if (dbuiltin) {
7643+
rkj = (*dbuiltin)(t1, t2, z);
7644+
} else {
7645+
rkj = (double) (*ibuiltin)(t1, t2, z);
7646+
}
7647+
gretl_matrix_set(ret, k, j+2, rkj);
7648+
z += M->rows;
7649+
}
7650+
z = zsave + count;
7651+
k++;
7652+
t1 = t2 = 0;
7653+
} else {
7654+
/* increment sample range */
7655+
t2++;
7656+
}
7657+
}
7658+
}
7659+
7660+
gretl_matrix_free(M);
7661+
7662+
return ret;
7663+
}
7664+
75567665
static int monthly_or_quarterly_dates (const DATASET *dset,
75577666
int pd, double sd0, int T,
75587667
double *x)

‎lib/src/genfuncs.h

+5
Original file line numberDiff line numberDiff line change
@@ -237,6 +237,11 @@ gretl_matrix *aggregate_by (const double *x,
237237
DATASET *dset,
238238
int *err);
239239

240+
gretl_matrix *matrix_aggregate (const gretl_matrix *X,
241+
const gretl_matrix *y,
242+
const char *func,
243+
int *err);
244+
240245
int fill_dataset_dates_series (const DATASET *dset, double *x);
241246

242247
int fill_day_of_week_array (double *dow,

‎lib/src/gretl_matrix.c

+1-1
Original file line numberDiff line numberDiff line change
@@ -15370,7 +15370,7 @@ gretl_matrix *gretl_matrix_sort_by_column (const gretl_matrix *m,
1537015370
return NULL;
1537115371
}
1537215372

15373-
a = gretl_matrix_copy(m);
15373+
a = gretl_matrix_alloc(m->rows, m->cols);
1537415374
if (a == NULL) {
1537515375
free(rs);
1537615376
*err = E_ALLOC;

0 commit comments

Comments
 (0)