src/language/stats/regression-export.h

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 2005 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 /*
  18   Exported C code for a regression model. The EXPORT subcommand causes PSPP
  19   to save a model as a small C program. This file contains some of the code
  20   of that saved program.
  21  */
  22 #ifndef REG_EXPORT_COMMENTS_H
  23 #define REG_EXPORT_COMMENTS_H
  24 const char reg_header[] = "#ifndef REG_EXPORT_COMMENTS_H\n#define REG_EXPORT_COMMENTS_H\n"
  25 "double pspp_reg_estimate (const double *, const char *[]);\n\n"
  26 "double pspp_reg_variance (const double *var_vals, const char *[]);\n\n"
  27 "double pspp_reg_confidence_interval_U "
  28 "(const double *var_vals, const char *var_names[], double p);\n\n"
  29 "double pspp_reg_confidence_interval_L "
  30 "(const double *var_vals, const char *var_names[], double p);\n\n"
  31 "double pspp_reg_prediction_interval_U "
  32 "(const double *var_vals, const char *var_names[], double p);\n\n"
  33 "double pspp_reg_prediction_interval_L "
  34 "(const double *var_vals, const char *var_names[], double p);\n"
  35 "#endif\n";
  36
  37 const char reg_preamble[] =  "/*\n   This program contains functions which return estimates\n"
  38 "   and confidence intervals for a linear model. The EXPORT subcommand\n"
  39 "   of the REGRESSION procedure of GNU PSPP generated this program.\n*/\n\n"
  40 "#include <string.h>\n#include <math.h>\n#define PSPP_REG_MAXLEN 1024\n\n";
  41
  42 const char reg_mean_cmt[] =  "/*\n   Estimate the mean of Y, the dependent variable for\n"
  43 "   the linear model of the form \n\n"
  44 "      Y = b0 + b1 * X1 + b2 * X2 + ... + bk * Xk + error\n\n"
  45 "   where X1, ..., Xk are the independent variables\n"
  46 "   whose values are stored in var_vals and whose names, \n"
  47 "   as known by PSPP, are stored in var_names. The estimated \n"
  48 "   regression coefficients (i.e., the estimates of b0,...,bk) \n"
  49 "   are stored in model_coeffs.\n*/\n";
  50
  51 const char reg_getvar[] = "{\n\t\tj = pspp_reg_getvar (var_names[i]);\n"
  52 "\t\testimate += var_vals[j] * model_coeffs[j];\n"
  53 "\t}\n\t\n\treturn estimate;\n}\n\n"
  54 "/*\n    Variance of an estimated mean of this form:\n\t"
  55 "Y = b0 + b1 * X1 + ... + bk * Xk\n    where X1,...Xk are the dependent variables,"
  56 " stored in\n    var_vals and b0,...,bk are the estimated regression coefficients.\n*/\n"
  57 "double\npspp_reg_variance (const double *var_vals, "
  58 "const char *var_names[])\n{\n\t";
  59
  60 const char reg_export_t_quantiles_1[] = "/*\n    Quantiles for the T distribution.\n*/\n"
  61 "static double\npspp_reg_t_quantile "
  62 "(double prob)\n{\n\n\tint i;\n\tdouble quantiles[] = {\n\t\t";
  63
  64 const char reg_export_t_quantiles_2[] = "i = (int) 100.0 * prob;\n\treturn quantiles[i];\n}\n";
  65
  66 const char reg_variance[] = "double result = 0.0;\n\n\tfor(i = 0; i < n_vars; i++)\n\t"
  67 "{\n\t\tj = pspp_reg_getvar (var_names[i]);\n\t\t"
  68 "unshuffled_vals[j] = var_vals[i];\n\t}\n\t"
  69 "for (i = 0; i < n_vars; i++)\n\t"
  70 "{\n\t\tresult += cov[i][i] * unshuffled_vals[i] * unshuffled_vals[i];\n\t\t"
  71 "for (j = i + 1; j < n_vars; j++)\n\t\t{\n\t\t\t"
  72 "result += 2.0 * cov[i][j] * unshuffled_vals[i] * unshuffled_vals[j];"
  73 "\n\t\t}\n\t}\n\treturn result;\n}\n";
  74
  75 const char reg_export_confidence_interval[] = "/*\n    Upper confidence limit for an "
  76 "estimated mean b0 + b1 * X1 + ... + bk * Xk.\n    The confidence interval is a "
  77 "100 * p percent confidence interval.\n*/\n"
  78 "double pspp_reg_confidence_interval_U "
  79 "(const double *var_vals, const char *var_names[], double p)\n{\n\t"
  80 "double result;\n\t"
  81 "result = sqrt (pspp_reg_variance (var_vals, var_names));\n\t"
  82 "result *= pspp_reg_t_quantile ((1.0 + p) / 2.0);\n\t"
  83 "result += pspp_reg_estimate (var_vals, var_names);\n\treturn result;\n}\n"
  84 "/*\n    Lower confidence limit for an "
  85 "estimated mean b0 + b1 * X1 + ... + bk * Xk.\n    The confidence interval is a "
  86 "100 * p percent confidence interval.\n*/\n"
  87 "double pspp_reg_confidence_interval_L "
  88 "(const double *var_vals, const char *var_names[], double p)\n{\n\t"
  89 "double result;\n\t"
  90 "result = -sqrt (pspp_reg_variance (var_vals, var_names));\n\t"
  91 "result *= pspp_reg_t_quantile ((1.0 + p) / 2.0);\n\t"
  92 "result += pspp_reg_estimate (var_vals, var_names);\n\treturn result;\n}\n";
  93
  94 const char reg_export_prediction_interval_1[] = "/*\n    Upper prediction limit for a "
  95 "predicted value b0 + b1 * X1 + ... + bk * Xk.\n    The prediction interval is a "
  96 "100 * p percent prediction interval.\n*/\n"
  97 "double pspp_reg_prediction_interval_U "
  98 "(const double *var_vals, const char *var_names[], double p)\n{\n\t"
  99 "double result;\n\tresult = sqrt (";
 100
 101 const char reg_export_prediction_interval_2[] = " + pspp_reg_variance (var_vals, var_names));\n"
 102 "\tresult *= pspp_reg_t_quantile ((1.0 + p) / 2.0);\n\t"
 103 "result += pspp_reg_estimate (var_vals, var_names);\n\treturn result;\n}\n"
 104 "/*\n    Lower prediction limit for a "
 105 "predicted value b0 + b1 * X1 + ... + bk * Xk.\n    The prediction interval is a "
 106 "100 * p percent prediction interval.\n*/\n"
 107 "double pspp_reg_prediction_interval_L "
 108 "(const double *var_vals, const char *var_names[], double p)\n{\n\t"
 109 "double result;\n\t"
 110 "result = -sqrt (";
 111
 112 const char reg_export_prediction_interval_3[] = " + pspp_reg_variance (var_vals, var_names));"
 113 "\n\tresult *= pspp_reg_t_quantile ((1.0 + p) / 2.0);\n\t"
 114 "result += pspp_reg_estimate (var_vals, var_names);\n\treturn result;\n}\n";
 115
 116 /*
 117   Change categorical values to binary vectors. The routine will use
 118   an encoding in which a categorical variable with n values is mapped
 119   to a vector with n-1 entries. Value 0 is mapped to the zero vector,
 120   value 1 is mapped to a vector whose first entry is 1 and all others are
 121   0, etc. For example, if a variable can have 'a', 'b' or 'c' as values,
 122   then the value 'a' will be encoded as (0,0), 'b' as (1,0) and 'c' as
 123   (0,1). If the design matrix used to create the model used a different
 124   encoding, then the function pspp_reg_categorical_encode () will return
 125   a vector which does not match its categorical value in the model.
 126  */
 127 const char reg_export_categorical_encode_1[] = "struct pspp_reg_categorical_variable\n"
 128 "{\n\tchar * name;\n\tsize_t n_vals;\n\tchar *values[1024];\n};\n\n"
 129 "/*\n   This function returns the binary vector which corresponds to the value\n"
 130 "   of the categorical variable stored in 'value'. The name of the variable is\n"
 131 "   stored in the 'var' argument. Notice the values stored in the\n"
 132 "   pspp_categorical_variable structures all end with a space character.\n"
 133 "   That means the values of the categorical variables you pass to any function\n"
 134 "   in this program should also end with a space character.\n*/\n"
 135 "static\ndouble * pspp_reg_get_value_vector (char *var, char *value)\n{\n\tdouble *result;\n\t"
 136 "int i;\n\t";
 137
 138 const char reg_export_categorical_encode_2[] = "int v_index = 0;\n\t"
 139 "while (v_index < n_vars && strncmp (var, varlist[i]->name, PSPP_REG_MAXLEN) != 0)\n\t{\n\t\t"
 140 "v_index++;\n\t}\n\tresult = (double *) malloc (varlist[v_index]->n_vals * sizeof (*result));\n\t"
 141 "for (i = 0; i < varlist[v_index]->n_vals; i++)\n\t{\n\t\t"
 142 "if (strncmp ( (varlist[v_index]->values)[i], value, PSPP_REG_MAXLEN) == 0)\n\t\t{\n\t\t\t"
 143 "result[i] = 1.0;\n\t\t}\n\t\telse result[i] = 0.0;\n\t}\n\n\t"
 144 "return result;\n}\n\n";
 145 #endif