pspp_coeff_var_to_coeff: Guard against a null pointer in coefs[i]->v_info.

[pspp-builds.git] / src / language / stats / glm.q

diff --git a/src/language/stats/glm.q b/src/language/stats/glm.q
index c015afa7ec5784ab5a46c8006b560ce28040f750..ddc2a78ba00388e8e507a2f32ce6e753aba25b77 100644 (file)
--- a/src/language/stats/glm.q
+++ b/src/language/stats/glm.q
@@ -47,7 +47,7 @@
 #include "xalloc.h"
 #include "gettext.h"
 
-#define GLM_LARGE_DATA 1000
+#define GLM_LARGE_DATA 10000
 
 /* (headers) */
 
@@ -151,13 +151,18 @@ glm_custom_dependent (struct lexer *lexer, struct dataset *ds,
   return 1;
 }
 
+/*
+  COV is the covariance matrix for variables included in the
+  model. That means the dependent variable is in there, too.
+ */
 static void
-coeff_init (pspp_linreg_cache * c, struct design_matrix *dm)
+coeff_init (pspp_linreg_cache * c, struct design_matrix *cov)
 {
-  c->coeff = xnmalloc (dm->m->size2 + 1, sizeof (*c->coeff));
+  c->coeff = xnmalloc (cov->m->size2, sizeof (*c->coeff));
+  c->n_coeffs = cov->m->size2 - 1;
   c->coeff[0] = xmalloc (sizeof (*(c->coeff[0])));     /* The first coefficient is the intercept. */
   c->coeff[0]->v_info = NULL;  /* Intercept has no associated variable. */
-  pspp_coeff_init (c->coeff + 1, dm);
+  pspp_coeff_init (c->coeff + 1, cov);
 }
 
 /*
@@ -257,7 +262,7 @@ run_glm (struct casereader *input,
   int n_indep = 0;
   struct ccase c;
   const struct variable **indep_vars;
-  struct variable **all_vars;
+  const struct variable **all_vars;
   struct design_matrix *X;
   struct moments_var **mom;
   struct casereader *reader;
@@ -307,9 +312,9 @@ run_glm (struct casereader *input,
 
   reader = casereader_clone (input);
   reader = casereader_create_filter_missing (reader, indep_vars, n_indep,
-                                            MV_ANY, NULL);
+                                            MV_ANY, NULL, NULL);
   reader = casereader_create_filter_missing (reader, v_dependent, 1,
-                                            MV_ANY, NULL);
+                                            MV_ANY, NULL, NULL);
   n_data = data_pass_one (casereader_clone (reader),
                          (const struct variable **) all_vars, n_all_vars,
                          mom);
@@ -339,6 +344,16 @@ run_glm (struct casereader *input,
                }
            }
        }
+      model = pspp_linreg_cache_alloc (v_dependent, indep_vars, n_data, n_indep);
+      /*
+       For large data sets, use QR decomposition.
+      */
+      if (n_data > sqrt (n_indep) && n_data > GLM_LARGE_DATA)
+       {
+         model->method = PSPP_LINREG_QR;
+       }
+      coeff_init (model, X);
+      pspp_linreg_with_cov (X, model);
       casereader_destroy (reader);
       for (i = 0; i < n_all_vars; i++)
        {
@@ -357,6 +372,7 @@ run_glm (struct casereader *input,
     }
   free (indep_vars);
   free (lopts.get_indep_mean_std);
+  pspp_linreg_cache_free (model);
   casereader_destroy (input);
 
   return true;