src/math/linreg/predict.c

   1 /*
   2    lib/linreg/predict.c
   3
   4    Copyright (C) 2005 Free Software Foundation, Inc. Written by Jason H. Stover.
   5
   6    This program is free software; you can redistribute it and/or modify it under
   7    the terms of the GNU General Public License as published by the Free
   8    Software Foundation; either version 2 of the License, or (at your option)
   9    any later version.
  10
  11    This program is distributed in the hope that it will be useful, but WITHOUT
  12    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  13    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License for
  14    more details.
  15
  16    You should have received a copy of the GNU General Public License along with
  17    this program; if not, write to the Free Software Foundation, Inc., 51
  18    Franklin Street, Fifth Floor, Boston, MA 02111-1307, USA.
  19  */
  20
  21 #include <math/linreg/linreg.h>
  22 #include <math/linreg/coefficient.h>
  23 #include <gl/xalloc.h>
  24
  25 /*
  26   Predict the value of the dependent variable with the
  27   new set of predictors. PREDICTORS must point to a list
  28   of variables, each of whose values are stored in VALS,
  29   in the same order.
  30  */
  31 double
  32 pspp_linreg_predict (const struct variable **predictors,
  33                      const union value **vals,
  34                      const void *c_, int n_vals)
  35 {
  36   const pspp_linreg_cache *c = c_;
  37   int i;
  38   int j;
  39   const struct pspp_linreg_coeff **found;
  40   const struct pspp_linreg_coeff *coe;
  41   double result;
  42   double tmp;
  43
  44   if (predictors == NULL || vals == NULL || c == NULL)
  45     {
  46       return GSL_NAN;
  47     }
  48   if (c->coeff == NULL)
  49     {
  50       /* The stupid model: just guess the mean. */
  51       return c->depvar_mean;
  52     }
  53   found = xnmalloc (c->n_coeffs, sizeof (*found));
  54   *found = c->coeff;
  55   result = c->coeff->estimate;  /* Intercept. */
  56
  57   /*
  58     The loops guard against the possibility that the caller passed us
  59     inadequate information, such as too few or too many values, or
  60     a redundant list of variable names.
  61    */
  62   for (j = 0; j < n_vals; j++)
  63     {
  64       coe = pspp_linreg_get_coeff (c, predictors[j], vals[j]);
  65       i = 1;
  66       while (found[i] == coe && i < c->n_coeffs)
  67         {
  68           i++;
  69         }
  70       if (i < c->n_coeffs)
  71         {
  72           found[i] = coe;
  73           tmp = pspp_linreg_coeff_get_est (coe);
  74           if (predictors[j]->type == NUMERIC)
  75             {
  76               tmp *= vals[j]->f;
  77             }
  78           result += tmp;
  79         }
  80     }
  81   free (found);
  82
  83   return result;
  84 }
  85
  86 double
  87 pspp_linreg_residual (const struct variable **predictors,
  88                       const union value **vals,
  89                       const union value *obs,
  90                       const void *c, int n_vals)
  91 {
  92   double pred;
  93   double result;
  94
  95   if (predictors == NULL || vals == NULL || c == NULL || obs == NULL)
  96     {
  97       return GSL_NAN;
  98     }
  99   pred = pspp_linreg_predict (predictors, vals, c, n_vals);
 100
 101   result = gsl_isnan (pred) ? GSL_NAN : (obs->f - pred);
 102   return result;
 103 }