pintos-os.org Git - pspp/blob - src/language/stats/glm.c

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 2010, 2011 Free Software Foundation, Inc.
   3
   4    This program is free software: you can redistribute it and/or modify
   5    it under the terms of the GNU General Public License as published by
   6    the Free Software Foundation, either version 3 of the License, or
   7    (at your option) any later version.
   8
   9    This program is distributed in the hope that it will be useful,
  10    but WITHOUT ANY WARRANTY; without even the implied warranty of
  11    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  12    GNU General Public License for more details.
  13
  14    You should have received a copy of the GNU General Public License
  15    along with this program.  If not, see <http://www.gnu.org/licenses/>. */
  16
  17 #include <config.h>
  18
  19 #include <gsl/gsl_cdf.h>
  20 #include <gsl/gsl_matrix.h>
  21 #include <math.h>
  22
  23 #include "data/case.h"
  24 #include "data/casegrouper.h"
  25 #include "data/casereader.h"
  26 #include "data/dataset.h"
  27 #include "data/dictionary.h"
  28 #include "data/format.h"
  29 #include "data/value.h"
  30 #include "language/command.h"
  31 #include "language/dictionary/split-file.h"
  32 #include "language/lexer/lexer.h"
  33 #include "language/lexer/value-parser.h"
  34 #include "language/lexer/variable-parser.h"
  35 #include "libpspp/ll.h"
  36 #include "libpspp/message.h"
  37 #include "libpspp/misc.h"
  38 #include "libpspp/taint.h"
  39 #include "linreg/sweep.h"
  40 #include "math/categoricals.h"
  41 #include "math/covariance.h"
  42 #include "math/interaction.h"
  43 #include "math/moments.h"
  44 #include "output/tab.h"
  45
  46 #include "gettext.h"
  47 #define _(msgid) gettext (msgid)
  48
  49 struct glm_spec
  50 {
  51   size_t n_dep_vars;
  52   const struct variable **dep_vars;
  53
  54   size_t n_factor_vars;
  55   const struct variable **factor_vars;
  56
  57   /* In the current implementation, design_vars will
  58      normally be the same as factor_vars.
  59      This will change once interactions, nested variables
  60      and repeated measures become involved.
  61   */
  62   size_t n_design_vars;
  63   const struct variable **design_vars;
  64
  65   size_t n_interactions;
  66   struct interaction **interactions;
  67
  68   enum mv_class exclude;
  69
  70   /* The weight variable */
  71   const struct variable *wv;
  72
  73   const struct dictionary *dict;
  74
  75   bool intercept;
  76
  77   double alpha;
  78 };
  79
  80 struct glm_workspace
  81 {
  82   double total_ssq;
  83   struct moments *totals;
  84
  85   struct categoricals *cats;
  86
  87   /*
  88      Sums of squares due to different variables. Element 0 is the SSE
  89      for the entire model. For i > 0, element i is the SS due to
  90      variable i.
  91    */
  92   gsl_vector *ssq;
  93 };
  94
  95 static void output_glm (const struct glm_spec *,
  96                         const struct glm_workspace *ws);
  97 static void run_glm (struct glm_spec *cmd, struct casereader *input,
  98                      const struct dataset *ds);
  99
 100
 101 static bool parse_design_spec (struct lexer *lexer, struct glm_spec *glm);
 102
 103
 104 int
 105 cmd_glm (struct lexer *lexer, struct dataset *ds)
 106 {
 107   int i;
 108   struct const_var_set *factors = NULL;
 109   struct glm_spec glm;
 110   bool design = false;
 111   glm.dict = dataset_dict (ds);
 112   glm.n_dep_vars = 0;
 113   glm.n_factor_vars = 0;
 114   glm.n_design_vars = 0;
 115   glm.n_interactions = 0;
 116   glm.interactions = NULL;
 117   glm.dep_vars = NULL;
 118   glm.factor_vars = NULL;
 119   glm.design_vars = NULL;
 120   glm.exclude = MV_ANY;
 121   glm.intercept = true;
 122   glm.wv = dict_get_weight (glm.dict);
 123   glm.alpha = 0.05;
 124
 125   if (!parse_variables_const (lexer, glm.dict,
 126                               &glm.dep_vars, &glm.n_dep_vars,
 127                               PV_NO_DUPLICATE | PV_NUMERIC))
 128     goto error;
 129
 130   lex_force_match (lexer, T_BY);
 131
 132   if (!parse_variables_const (lexer, glm.dict,
 133                               &glm.factor_vars, &glm.n_factor_vars,
 134                               PV_NO_DUPLICATE | PV_NUMERIC))
 135     goto error;
 136
 137   if (glm.n_dep_vars > 1)
 138     {
 139       msg (ME, _("Multivariate analysis is not yet implemented"));
 140       return CMD_FAILURE;
 141     }
 142
 143   factors =
 144     const_var_set_create_from_array (glm.factor_vars, glm.n_factor_vars);
 145
 146   while (lex_token (lexer) != T_ENDCMD)
 147     {
 148       lex_match (lexer, T_SLASH);
 149
 150       if (lex_match_id (lexer, "MISSING"))
 151         {
 152           lex_match (lexer, T_EQUALS);
 153           while (lex_token (lexer) != T_ENDCMD
 154                  && lex_token (lexer) != T_SLASH)
 155             {
 156               if (lex_match_id (lexer, "INCLUDE"))
 157                 {
 158                   glm.exclude = MV_SYSTEM;
 159                 }
 160               else if (lex_match_id (lexer, "EXCLUDE"))
 161                 {
 162                   glm.exclude = MV_ANY;
 163                 }
 164               else
 165                 {
 166                   lex_error (lexer, NULL);
 167                   goto error;
 168                 }
 169             }
 170         }
 171       else if (lex_match_id (lexer, "INTERCEPT"))
 172         {
 173           lex_match (lexer, T_EQUALS);
 174           while (lex_token (lexer) != T_ENDCMD
 175                  && lex_token (lexer) != T_SLASH)
 176             {
 177               if (lex_match_id (lexer, "INCLUDE"))
 178                 {
 179                   glm.intercept = true;
 180                 }
 181               else if (lex_match_id (lexer, "EXCLUDE"))
 182                 {
 183                   glm.intercept = false;
 184                 }
 185               else
 186                 {
 187                   lex_error (lexer, NULL);
 188                   goto error;
 189                 }
 190             }
 191         }
 192       else if (lex_match_id (lexer, "CRITERIA"))
 193         {
 194           lex_match (lexer, T_EQUALS);
 195           if (lex_match_id (lexer, "ALPHA"))
 196             {
 197               if (lex_force_match (lexer, T_LPAREN))
 198                 {
 199                   if (! lex_force_num (lexer))
 200                     {
 201                       lex_error (lexer, NULL);
 202                       goto error;
 203                     }
 204
 205                   glm.alpha = lex_number (lexer);
 206                   lex_get (lexer);
 207                   if ( ! lex_force_match (lexer, T_RPAREN))
 208                     {
 209                       lex_error (lexer, NULL);
 210                       goto error;
 211                     }
 212                 }
 213             }
 214           else
 215             {
 216               lex_error (lexer, NULL);
 217               goto error;
 218             }
 219         }
 220       else if (lex_match_id (lexer, "METHOD"))
 221         {
 222           lex_match (lexer, T_EQUALS);
 223           if ( !lex_force_match_id (lexer, "SSTYPE"))
 224             {
 225               lex_error (lexer, NULL);
 226               goto error;
 227             }
 228
 229           if ( ! lex_force_match (lexer, T_LPAREN))
 230             {
 231               lex_error (lexer, NULL);
 232               goto error;
 233             }
 234
 235           if ( ! lex_force_int (lexer))
 236             {
 237               lex_error (lexer, NULL);
 238               goto error;
 239             }
 240
 241           if (3 != lex_integer (lexer))
 242             {
 243               msg (ME, _("Only type 3 sum of squares are currently implemented"));
 244               goto error;
 245             }
 246
 247           lex_get (lexer);
 248
 249           if ( ! lex_force_match (lexer, T_RPAREN))
 250             {
 251               lex_error (lexer, NULL);
 252               goto error;
 253             }
 254         }
 255       else if (lex_match_id (lexer, "DESIGN"))
 256         {
 257           lex_match (lexer, T_EQUALS);
 258
 259           if (! parse_design_spec (lexer, &glm))
 260             goto error;
 261
 262           if ( glm.n_interactions == 0)
 263             {
 264               msg (ME, _("One or more design  variables must be given"));
 265               goto error;
 266             }
 267
 268           design = true;
 269         }
 270       else
 271         {
 272           lex_error (lexer, NULL);
 273           goto error;
 274         }
 275     }
 276
 277   if ( ! design )
 278     {
 279       lex_error (lexer, _("/DESIGN is mandatory in GLM"));
 280       goto error;
 281     }
 282
 283   {
 284     struct casegrouper *grouper;
 285     struct casereader *group;
 286     bool ok;
 287
 288     grouper = casegrouper_create_splits (proc_open (ds), glm.dict);
 289     while (casegrouper_get_next_group (grouper, &group))
 290       run_glm (&glm, group, ds);
 291     ok = casegrouper_destroy (grouper);
 292     ok = proc_commit (ds) && ok;
 293   }
 294
 295   const_var_set_destroy (factors);
 296   free (glm.factor_vars);
 297   for (i = 0 ; i < glm.n_interactions; ++i)
 298     interaction_destroy (glm.interactions[i]);
 299   free (glm.interactions);
 300   free (glm.dep_vars);
 301   free (glm.design_vars);
 302
 303
 304   return CMD_SUCCESS;
 305
 306 error:
 307
 308   const_var_set_destroy (factors);
 309   free (glm.factor_vars);
 310   for (i = 0 ; i < glm.n_interactions; ++i)
 311     interaction_destroy (glm.interactions[i]);
 312
 313   free (glm.interactions);
 314   free (glm.dep_vars);
 315   free (glm.design_vars);
 316
 317   return CMD_FAILURE;
 318 }
 319
 320 static void get_ssq (struct covariance *, gsl_vector *,
 321                      const struct glm_spec *);
 322
 323 static bool
 324 not_dropped (size_t j, size_t * dropped, size_t n_dropped)
 325 {
 326   size_t i;
 327
 328   for (i = 0; i < n_dropped; i++)
 329     {
 330       if (j == dropped[i])
 331         return false;
 332     }
 333   return true;
 334 }
 335
 336 static void
 337 get_ssq (struct covariance *cov, gsl_vector * ssq, const struct glm_spec *cmd)
 338 {
 339   const struct variable **vars;
 340   gsl_matrix *small_cov = NULL;
 341   gsl_matrix *cm = covariance_calculate_unnormalized (cov);
 342   size_t i;
 343   size_t j;
 344   size_t k;
 345   size_t n;
 346   size_t m;
 347   size_t *dropped;
 348   size_t n_dropped;
 349
 350   dropped = xcalloc (covariance_dim (cov), sizeof (*dropped));
 351   vars = xcalloc (covariance_dim (cov), sizeof (*vars));
 352   covariance_get_var_indices (cov, vars);
 353
 354   for (k = 0; k < cmd->n_design_vars; k++)
 355     {
 356       n_dropped = 0;
 357       for (i = 1; i < covariance_dim (cov); i++)
 358         {
 359           if (vars[i] == cmd->design_vars[k])
 360             {
 361               dropped[n_dropped++] = i;
 362             }
 363         }
 364       small_cov =
 365         gsl_matrix_alloc (cm->size1 - n_dropped, cm->size2 - n_dropped);
 366       gsl_matrix_set (small_cov, 0, 0, gsl_matrix_get (cm, 0, 0));
 367       n = 0;
 368       m = 0;
 369       for (i = 0; i < cm->size1; i++)
 370         {
 371           if (not_dropped (i, dropped, n_dropped))
 372             {
 373               m = 0;
 374               for (j = 0; j < cm->size2; j++)
 375                 {
 376                   if (not_dropped (j, dropped, n_dropped))
 377                     {
 378                       gsl_matrix_set (small_cov, n, m,
 379                                       gsl_matrix_get (cm, i, j));
 380                       m++;
 381                     }
 382                 }
 383               n++;
 384             }
 385         }
 386       reg_sweep (small_cov, 0);
 387       gsl_vector_set (ssq, k + 1,
 388                       gsl_matrix_get (small_cov, 0, 0)
 389                       - gsl_vector_get (ssq, 0));
 390       gsl_matrix_free (small_cov);
 391     }
 392
 393   free (dropped);
 394   free (vars);
 395   gsl_matrix_free (cm);
 396 }
 397
 398 //static  void dump_matrix (const gsl_matrix *m);
 399
 400 static void
 401 run_glm (struct glm_spec *cmd, struct casereader *input,
 402          const struct dataset *ds)
 403 {
 404   bool warn_bad_weight = true;
 405   int v;
 406   struct taint *taint;
 407   struct dictionary *dict = dataset_dict (ds);
 408   struct casereader *reader;
 409   struct ccase *c;
 410
 411   struct glm_workspace ws;
 412   struct covariance *cov;
 413
 414   ws.cats = categoricals_create (cmd->interactions, cmd->n_interactions,
 415                                  cmd->wv, cmd->exclude,
 416                                  NULL, NULL, NULL, NULL);
 417
 418   cov = covariance_2pass_create (cmd->n_dep_vars, cmd->dep_vars,
 419                                  ws.cats, cmd->wv, cmd->exclude);
 420
 421
 422   c = casereader_peek (input, 0);
 423   if (c == NULL)
 424     {
 425       casereader_destroy (input);
 426       return;
 427     }
 428   output_split_file_values (ds, c);
 429   case_unref (c);
 430
 431   taint = taint_clone (casereader_get_taint (input));
 432
 433   ws.totals = moments_create (MOMENT_VARIANCE);
 434
 435   for (reader = casereader_clone (input);
 436        (c = casereader_read (reader)) != NULL; case_unref (c))
 437     {
 438       double weight = dict_get_case_weight (dict, c, &warn_bad_weight);
 439
 440       for (v = 0; v < cmd->n_dep_vars; ++v)
 441         moments_pass_one (ws.totals, case_data (c, cmd->dep_vars[v])->f,
 442                           weight);
 443
 444       covariance_accumulate_pass1 (cov, c);
 445     }
 446   casereader_destroy (reader);
 447
 448   categoricals_done (ws.cats);
 449
 450   for (reader = input;
 451        (c = casereader_read (reader)) != NULL; case_unref (c))
 452     {
 453       double weight = dict_get_case_weight (dict, c, &warn_bad_weight);
 454
 455       for (v = 0; v < cmd->n_dep_vars; ++v)
 456         moments_pass_two (ws.totals, case_data (c, cmd->dep_vars[v])->f,
 457                           weight);
 458
 459       covariance_accumulate_pass2 (cov, c);
 460     }
 461   casereader_destroy (reader);
 462
 463   {
 464     gsl_matrix *cm = covariance_calculate_unnormalized (cov);
 465
 466     //    dump_matrix (cm);
 467
 468     ws.total_ssq = gsl_matrix_get (cm, 0, 0);
 469
 470     reg_sweep (cm, 0);
 471
 472     /*
 473       Store the overall SSE.
 474     */
 475     ws.ssq = gsl_vector_alloc (cm->size1);
 476     gsl_vector_set (ws.ssq, 0, gsl_matrix_get (cm, 0, 0));
 477     get_ssq (cov, ws.ssq, cmd);
 478     //    dump_matrix (cm);
 479
 480     gsl_matrix_free (cm);
 481   }
 482
 483   if (!taint_has_tainted_successor (taint))
 484     output_glm (cmd, &ws);
 485
 486   gsl_vector_free (ws.ssq);
 487
 488   covariance_destroy (cov);
 489   moments_destroy (ws.totals);
 490
 491   taint_destroy (taint);
 492 }
 493
 494 static void
 495 output_glm (const struct glm_spec *cmd, const struct glm_workspace *ws)
 496 {
 497   const struct fmt_spec *wfmt =
 498     cmd->wv ? var_get_print_format (cmd->wv) : &F_8_0;
 499
 500   double n_total, mean;
 501   double df_corr = 0.0;
 502   double mse = 0;
 503
 504   int f;
 505   int r;
 506   const int heading_columns = 1;
 507   const int heading_rows = 1;
 508   struct tab_table *t;
 509
 510   const int nc = 6;
 511   int nr = heading_rows + 4 + cmd->n_interactions;
 512   if (cmd->intercept)
 513     nr++;
 514
 515   t = tab_create (nc, nr);
 516   tab_title (t, _("Tests of Between-Subjects Effects"));
 517
 518   tab_headers (t, heading_columns, 0, heading_rows, 0);
 519
 520   tab_box (t, TAL_2, TAL_2, -1, TAL_1, 0, 0, nc - 1, nr - 1);
 521
 522   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
 523   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
 524
 525   tab_text (t, 0, 0, TAB_CENTER | TAT_TITLE, _("Source"));
 526
 527   /* TRANSLATORS: The parameter is a roman numeral */
 528   tab_text_format (t, 1, 0, TAB_CENTER | TAT_TITLE,
 529                    _("Type %s Sum of Squares"), "III");
 530   tab_text (t, 2, 0, TAB_CENTER | TAT_TITLE, _("df"));
 531   tab_text (t, 3, 0, TAB_CENTER | TAT_TITLE, _("Mean Square"));
 532   tab_text (t, 4, 0, TAB_CENTER | TAT_TITLE, _("F"));
 533   tab_text (t, 5, 0, TAB_CENTER | TAT_TITLE, _("Sig."));
 534
 535   moments_calculate (ws->totals, &n_total, &mean, NULL, NULL, NULL);
 536
 537   if (cmd->intercept)
 538     df_corr += 1.0;
 539
 540   for (f = 0; f < cmd->n_interactions; ++f)
 541     df_corr += categoricals_n_count (ws->cats, f) - 1.0;
 542
 543   mse = gsl_vector_get (ws->ssq, 0) / (n_total - df_corr);
 544
 545   r = heading_rows;
 546   tab_text (t, 0, r, TAB_LEFT | TAT_TITLE, _("Corrected Model"));
 547
 548   r++;
 549
 550   if (cmd->intercept)
 551     {
 552       const double intercept = pow2 (mean * n_total) / n_total;
 553       const double df = 1.0;
 554       const double F = intercept / df / mse;
 555       tab_text (t, 0, r, TAB_LEFT | TAT_TITLE, _("Intercept"));
 556       tab_double (t, 1, r, 0, intercept, NULL);
 557       tab_double (t, 2, r, 0, 1.00, wfmt);
 558       tab_double (t, 3, r, 0, intercept / df, NULL);
 559       tab_double (t, 4, r, 0, F, NULL);
 560       tab_double (t, 5, r, 0, gsl_cdf_fdist_Q (F, df, n_total - df_corr),
 561                   NULL);
 562       r++;
 563     }
 564
 565   for (f = 0; f < cmd->n_interactions; ++f)
 566     {
 567       struct string str = DS_EMPTY_INITIALIZER;
 568       const double df = categoricals_n_count (ws->cats, f) - 1.0;
 569       const double ssq = gsl_vector_get (ws->ssq, f + 1);
 570       const double F = ssq / df / mse;
 571       interaction_to_string (cmd->interactions[f], &str);
 572       tab_text (t, 0, r, TAB_LEFT | TAT_TITLE, ds_cstr (&str));
 573       ds_destroy (&str);
 574
 575       tab_double (t, 1, r, 0, ssq, NULL);
 576       tab_double (t, 2, r, 0, df, wfmt);
 577       tab_double (t, 3, r, 0, ssq / df, NULL);
 578       tab_double (t, 4, r, 0, F, NULL);
 579
 580       tab_double (t, 5, r, 0, gsl_cdf_fdist_Q (F, df, n_total - df_corr),
 581                   NULL);
 582       r++;
 583     }
 584
 585   {
 586     /* Corrected Model */
 587     const double df = df_corr - 1.0;
 588     const double ssq = ws->total_ssq - gsl_vector_get (ws->ssq, 0);
 589     const double F = ssq / df / mse;
 590     tab_double (t, 1, heading_rows, 0, ssq, NULL);
 591     tab_double (t, 2, heading_rows, 0, df, wfmt);
 592     tab_double (t, 3, heading_rows, 0, ssq / df, NULL);
 593     tab_double (t, 4, heading_rows, 0, F, NULL);
 594
 595     tab_double (t, 5, heading_rows, 0,
 596                 gsl_cdf_fdist_Q (F, df, n_total - df_corr), NULL);
 597   }
 598
 599   {
 600     const double df = n_total - df_corr;
 601     const double ssq = gsl_vector_get (ws->ssq, 0);
 602     const double mse = ssq / df;
 603     tab_text (t, 0, r, TAB_LEFT | TAT_TITLE, _("Error"));
 604     tab_double (t, 1, r, 0, ssq, NULL);
 605     tab_double (t, 2, r, 0, df, wfmt);
 606     tab_double (t, 3, r++, 0, mse, NULL);
 607   }
 608
 609   if (cmd->intercept)
 610     {
 611       const double intercept = pow2 (mean * n_total) / n_total;
 612       const double ssq = intercept + ws->total_ssq;
 613
 614       tab_text (t, 0, r, TAB_LEFT | TAT_TITLE, _("Total"));
 615       tab_double (t, 1, r, 0, ssq, NULL);
 616       tab_double (t, 2, r, 0, n_total, wfmt);
 617
 618       r++;
 619     }
 620
 621   tab_text (t, 0, r, TAB_LEFT | TAT_TITLE, _("Corrected Total"));
 622
 623
 624   tab_double (t, 1, r, 0, ws->total_ssq, NULL);
 625   tab_double (t, 2, r, 0, n_total - 1.0, wfmt);
 626
 627   tab_submit (t);
 628 }
 629
 630 #if 0
 631 static void
 632 dump_matrix (const gsl_matrix * m)
 633 {
 634   size_t i, j;
 635   for (i = 0; i < m->size1; ++i)
 636     {
 637       for (j = 0; j < m->size2; ++j)
 638         {
 639           double x = gsl_matrix_get (m, i, j);
 640           printf ("%.3f ", x);
 641         }
 642       printf ("\n");
 643     }
 644   printf ("\n");
 645 }
 646 #endif
 647
 648
 649 \f
 650
 651 /* Match a variable.
 652    If the match succeeds, the variable will be placed in VAR.
 653    Returns true if successful */
 654 static bool
 655 lex_match_variable (struct lexer *lexer, const struct glm_spec *glm, const struct variable **var)
 656 {
 657   if (lex_token (lexer) !=  T_ID)
 658     return false;
 659
 660   *var = parse_variable_const  (lexer, glm->dict);
 661
 662   if ( *var == NULL)
 663     return false;
 664   return true;
 665 }
 666
 667 /* An interaction is a variable followed by {*, BY} followed by an interaction */
 668 static bool
 669 parse_design_interaction (struct lexer *lexer, struct glm_spec *glm, struct interaction **iact)
 670 {
 671   const struct variable *v = NULL;
 672   assert (iact);
 673
 674   switch  (lex_next_token (lexer, 1))
 675     {
 676     case T_ENDCMD:
 677     case T_SLASH:
 678     case T_COMMA:
 679     case T_ID:
 680     case T_BY:
 681     case T_ASTERISK:
 682       break;
 683     default:
 684       return false;
 685       break;
 686     }
 687
 688   if (! lex_match_variable (lexer, glm, &v))
 689     {
 690       interaction_destroy (*iact);
 691       *iact = NULL;
 692       return false;
 693     }
 694
 695   assert (v);
 696
 697   if ( *iact == NULL)
 698     *iact = interaction_create (v);
 699   else
 700     interaction_add_variable (*iact, v);
 701
 702   if ( lex_match (lexer, T_ASTERISK) || lex_match (lexer, T_BY))
 703     {
 704       lex_error (lexer, "Interactions are not yet implemented"); return false;
 705       return parse_design_interaction (lexer, glm, iact);
 706     }
 707
 708   glm->n_design_vars++;
 709   glm->design_vars = xrealloc (glm->design_vars, sizeof (*glm->design_vars) * glm->n_design_vars);
 710   glm->design_vars[glm->n_design_vars - 1] = v;
 711
 712   return true;
 713 }
 714
 715 static bool
 716 parse_nested_variable (struct lexer *lexer, struct glm_spec *glm)
 717 {
 718   const struct variable *v = NULL;
 719   if ( ! lex_match_variable (lexer, glm, &v))
 720     return false;
 721
 722   if (lex_match (lexer, T_LPAREN))
 723     {
 724       if ( ! parse_nested_variable (lexer, glm))
 725         return false;
 726
 727       if ( ! lex_force_match (lexer, T_RPAREN))
 728         return false;
 729     }
 730
 731   lex_error (lexer, "Nested variables are not yet implemented"); return false;
 732   return true;
 733 }
 734
 735 /* A design term is an interaction OR a nested variable */
 736 static bool
 737 parse_design_term (struct lexer *lexer, struct glm_spec *glm)
 738 {
 739   struct interaction *iact = NULL;
 740   if (parse_design_interaction (lexer, glm, &iact))
 741     {
 742       /* Interaction parsing successful.  Add to list of interactions */
 743       glm->interactions = xrealloc (glm->interactions, sizeof *glm->interactions * ++glm->n_interactions);
 744       glm->interactions[glm->n_interactions - 1] = iact;
 745       return true;
 746     }
 747
 748   if ( parse_nested_variable (lexer, glm))
 749     return true;
 750
 751   return false;
 752 }
 753
 754
 755
 756 /* Parse a complete DESIGN specification.
 757    A design spec is a design term, optionally followed by a comma,
 758    and another design spec.
 759 */
 760 static bool
 761 parse_design_spec (struct lexer *lexer, struct glm_spec *glm)
 762 {
 763   if  (lex_token (lexer) == T_ENDCMD || lex_token (lexer) == T_SLASH)
 764     return true;
 765
 766   if ( ! parse_design_term (lexer, glm))
 767     return false;
 768
 769   lex_match (lexer, T_COMMA);
 770
 771   return parse_design_spec (lexer, glm);
 772 }
 773