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 pspp_linreg_cache * c, int n_vals)
  35 {
  36   int i;
  37   int j;
  38   const struct pspp_linreg_coeff **found;
  39   const struct pspp_linreg_coeff *coe;
  40   double result;
  41   double tmp;
  42
  43   if (predictors == NULL || vals == NULL || c == NULL)
  44     {
  45       return GSL_NAN;
  46     }
  47   if (c->coeff == NULL)
  48     {
  49       /* The stupid model: just guess the mean. */
  50       return c->depvar_mean;
  51     }
  52   found = xnmalloc (c->n_coeffs, sizeof (*found));
  53   *found = c->coeff;
  54   result = c->coeff->estimate;  /* Intercept. */
  55
  56   /*
  57     The loops guard against the possibility that the caller passed us
  58     inadequate information, such as too few or too many values, or
  59     a redundant list of variable names.
  60    */
  61   for (j = 0; j < n_vals; j++)
  62     {
  63       coe = pspp_linreg_get_coeff (c, predictors[j], vals[j]);
  64       i = 1;
  65       while (found[i] == coe && i < c->n_coeffs)
  66         {
  67           i++;
  68         }
  69       if (i < c->n_coeffs)
  70         {
  71           found[i] = coe;
  72           tmp = pspp_linreg_coeff_get_est (coe);
  73           if (predictors[j]->type == NUMERIC)
  74             {
  75               tmp *= vals[j]->f;
  76             }
  77           result += tmp;
  78         }
  79     }
  80   free (found);
  81
  82   return result;
  83 }
  84
  85 double
  86 pspp_linreg_residual (const struct variable **predictors,
  87                       const union value **vals,
  88                       const union value *obs,
  89                       const pspp_linreg_cache * c, int n_vals)
  90 {
  91   double pred;
  92   double result;
  93
  94   if (predictors == NULL || vals == NULL || c == NULL || obs == NULL)
  95     {
  96       return GSL_NAN;
  97     }
  98   pred = pspp_linreg_predict (predictors, vals, c, n_vals);
  99
 100   result = gsl_isnan (pred) ? GSL_NAN : (obs->f - pred);
 101   return result;
 102 }