33 |
#include "gettext.h" |
#include "gettext.h" |
34 |
#include "lexer.h" |
#include "lexer.h" |
35 |
#include <linreg/pspp_linreg.h> |
#include <linreg/pspp_linreg.h> |
36 |
|
#include "missing-values.h" |
37 |
#include "tab.h" |
#include "tab.h" |
38 |
#include "var.h" |
#include "var.h" |
39 |
#include "vfm.h" |
#include "vfm.h" |
500 |
size_t k; |
size_t k; |
501 |
size_t n_data = 0; |
size_t n_data = 0; |
502 |
size_t row; |
size_t row; |
503 |
|
size_t case_num; |
504 |
int n_indep; |
int n_indep; |
505 |
int j = 0; |
int j = 0; |
506 |
|
/* |
507 |
|
Keep track of the missing cases. |
508 |
|
*/ |
509 |
|
int *is_missing_case; |
510 |
const union value *val; |
const union value *val; |
511 |
struct casereader *r; |
struct casereader *r; |
512 |
struct casereader *r2; |
struct casereader *r2; |
520 |
pspp_linreg_opts lopts; |
pspp_linreg_opts lopts; |
521 |
|
|
522 |
n_data = casefile_get_case_cnt (cf); |
n_data = casefile_get_case_cnt (cf); |
523 |
|
|
524 |
|
is_missing_case = xnmalloc (n_data, sizeof (*is_missing_case)); |
525 |
|
for (i = 0; i < n_data; i++) |
526 |
|
is_missing_case[i] = 0; |
527 |
|
|
528 |
n_indep = cmd.n_variables - cmd.n_dependent; |
n_indep = cmd.n_variables - cmd.n_dependent; |
529 |
indep_vars = xnmalloc (n_indep, sizeof *indep_vars); |
indep_vars = xnmalloc (n_indep, sizeof *indep_vars); |
530 |
|
|
|
Y = gsl_vector_alloc (n_data); |
|
531 |
lopts.get_depvar_mean_std = 1; |
lopts.get_depvar_mean_std = 1; |
532 |
lopts.get_indep_mean_std = xnmalloc (n_indep, sizeof (int)); |
lopts.get_indep_mean_std = xnmalloc (n_indep, sizeof (int)); |
533 |
|
|
|
lcache = pspp_linreg_cache_alloc (n_data, n_indep); |
|
|
lcache->indep_means = gsl_vector_alloc (n_indep); |
|
|
lcache->indep_std = gsl_vector_alloc (n_indep); |
|
534 |
|
|
535 |
/* |
/* |
536 |
Read from the active file. The first pass encodes categorical |
Read from the active file. The first pass encodes categorical |
537 |
variables. |
variables and drops cases with missing values. |
538 |
*/ |
*/ |
539 |
ca = cr_recoded_cat_ar_create (cmd.n_variables, cmd.v_variables); |
ca = cr_recoded_cat_ar_create (cmd.n_variables, cmd.v_variables); |
540 |
for (r = casefile_get_reader (cf); |
for (r = casefile_get_reader (cf); |
541 |
casereader_read (r, &c); case_destroy (&c)) |
casereader_read (r, &c); case_destroy (&c)) |
542 |
{ |
{ |
543 |
|
row = casereader_cnum (r) - 1; |
544 |
for (i = 0; i < ca->n_vars; i++) |
for (i = 0; i < ca->n_vars; i++) |
545 |
{ |
{ |
546 |
v = (*(ca->a + i))->v; |
v = (*(ca->a + i))->v; |
547 |
val = case_data (&c, v->fv); |
val = case_data (&c, v->fv); |
548 |
cr_value_update (*(ca->a + i), val); |
cr_value_update (*(ca->a + i), val); |
549 |
} |
} |
550 |
|
for (i = 0; i < cmd.n_variables; i++) |
551 |
|
{ |
552 |
|
v = cmd.v_variables[i]; |
553 |
|
val = case_data (&c, v->fv); |
554 |
|
if (mv_is_value_missing (&v->miss, val)) |
555 |
|
{ |
556 |
|
n_data--; |
557 |
|
is_missing_case[row] = 1; |
558 |
|
} |
559 |
|
} |
560 |
} |
} |
561 |
|
Y = gsl_vector_alloc (n_data); |
562 |
cr_create_value_matrices (ca); |
cr_create_value_matrices (ca); |
563 |
X = |
X = |
564 |
design_matrix_create (n_indep, (const struct variable **) cmd.v_variables, |
design_matrix_create (n_indep, (const struct variable **) cmd.v_variables, |
565 |
ca, n_data); |
ca, n_data); |
566 |
|
lcache = pspp_linreg_cache_alloc (n_data, n_indep); |
567 |
|
lcache->indep_means = gsl_vector_alloc (n_indep); |
568 |
|
lcache->indep_std = gsl_vector_alloc (n_indep); |
569 |
|
|
570 |
/* |
/* |
571 |
The second pass creates the design matrix. |
The second pass creates the design matrix. |
572 |
*/ |
*/ |
573 |
|
row = 0; |
574 |
for (r2 = casefile_get_reader (cf); casereader_read (r2, &c); |
for (r2 = casefile_get_reader (cf); casereader_read (r2, &c); |
575 |
case_destroy (&c)) |
case_destroy (&c)) |
576 |
/* Iterate over the cases. */ |
/* Iterate over the cases. */ |
577 |
{ |
{ |
578 |
k = 0; |
k = 0; |
579 |
row = casereader_cnum (r2) - 1; |
case_num = casereader_cnum (r2) - 1; |
580 |
for (i = 0; i < cmd.n_variables; ++i) /* Iterate over the variables |
if (!is_missing_case[case_num]) |
|
for the current case. |
|
|
*/ |
|
581 |
{ |
{ |
582 |
v = cmd.v_variables[i]; |
for (i = 0; i < cmd.n_variables; ++i) /* Iterate over the variables |
583 |
val = case_data (&c, v->fv); |
for the current case. |
584 |
/* |
*/ |
|
Independent/dependent variable separation. The |
|
|
'variables' subcommand specifies a varlist which contains |
|
|
both dependent and independent variables. The dependent |
|
|
variables are specified with the 'dependent' |
|
|
subcommand. We need to separate the two. |
|
|
*/ |
|
|
if (is_depvar (i)) |
|
|
{ |
|
|
if (v->type != NUMERIC) |
|
|
{ |
|
|
msg (SE, gettext ("Dependent variable must be numeric.")); |
|
|
pspp_reg_rc = CMD_FAILURE; |
|
|
return; |
|
|
} |
|
|
lcache->depvar = (const struct var *) v; |
|
|
gsl_vector_set (Y, row, val->f); |
|
|
} |
|
|
else |
|
585 |
{ |
{ |
586 |
if (v->type == ALPHA) |
v = cmd.v_variables[i]; |
587 |
|
val = case_data (&c, v->fv); |
588 |
|
/* |
589 |
|
Independent/dependent variable separation. The |
590 |
|
'variables' subcommand specifies a varlist which contains |
591 |
|
both dependent and independent variables. The dependent |
592 |
|
variables are specified with the 'dependent' |
593 |
|
subcommand. We need to separate the two. |
594 |
|
*/ |
595 |
|
if (is_depvar (i)) |
596 |
{ |
{ |
597 |
rc = cr_var_to_recoded_categorical (v, ca); |
if (v->type != NUMERIC) |
598 |
design_matrix_set_categorical (X, row, v, val, rc); |
{ |
599 |
|
msg (SE, gettext ("Dependent variable must be numeric.")); |
600 |
|
pspp_reg_rc = CMD_FAILURE; |
601 |
|
return; |
602 |
|
} |
603 |
|
lcache->depvar = (const struct var *) v; |
604 |
|
gsl_vector_set (Y, row, val->f); |
605 |
} |
} |
606 |
else if (v->type == NUMERIC) |
else |
607 |
{ |
{ |
608 |
design_matrix_set_numeric (X, row, v, val); |
if (v->type == ALPHA) |
609 |
|
{ |
610 |
|
rc = cr_var_to_recoded_categorical (v, ca); |
611 |
|
design_matrix_set_categorical (X, row, v, val, rc); |
612 |
|
} |
613 |
|
else if (v->type == NUMERIC) |
614 |
|
{ |
615 |
|
design_matrix_set_numeric (X, row, v, val); |
616 |
|
} |
617 |
|
|
618 |
|
indep_vars[k] = i; |
619 |
|
k++; |
620 |
|
lopts.get_indep_mean_std[i] = 1; |
621 |
} |
} |
|
|
|
|
indep_vars[k] = i; |
|
|
k++; |
|
|
lopts.get_indep_mean_std[i] = 1; |
|
622 |
} |
} |
623 |
|
row++; |
624 |
} |
} |
625 |
} |
} |
626 |
/* |
/* |
627 |
Now that we know the number of coefficients, allocate space |
Now that we know the number of coefficients, allocate space |
628 |
and store pointers to the variables that correspond to the |
and store pointers to the variables that correspond to the |
629 |
coefficients. |
coefficients. |
630 |
*/ |
*/ |
631 |
lcache->coeff = xnmalloc (X->m->size2 + 1, sizeof (*lcache->coeff)); |
lcache->coeff = xnmalloc (X->m->size2 + 1, sizeof (*lcache->coeff)); |
632 |
for (i = 0; i < X->m->size2; i++) |
for (i = 0; i < X->m->size2; i++) |
646 |
pspp_linreg_cache_free (lcache); |
pspp_linreg_cache_free (lcache); |
647 |
free (lopts.get_indep_mean_std); |
free (lopts.get_indep_mean_std); |
648 |
free (indep_vars); |
free (indep_vars); |
649 |
|
free (is_missing_case); |
650 |
casereader_destroy (r); |
casereader_destroy (r); |
651 |
return; |
return; |
652 |
} |
} |