一般化線形モデルを用いた予測

nag_glm_predict (g02gpc) Example Program Data

Training Data
Nag_Reci
Nag_MeanInclude
Nag_FALSE
Nag_FALSE
5 1 0.0 0           : slink,smean,st_ioffset,st_weight,t_n,m,scale,print_iter
1.0 25.0            : t_x,y
2.0 10.0
3.0  6.0
4.0  4.0
5.0  3.0
1                   : sx

Prediction Data
2 Nag_TRUE
Nag_FALSE Nag_FALSE : n,svfobs,soffset,sweight
32.0                : x
18.0

nag_glm_predict (g02gpc) Example Program Results

Residual sum of squares =       0.3872, Degrees of freedom = 3.000000

       Estimate     Standard error

        -0.0239         0.0028
         0.0638         0.0026

   I      ETA         SE(ETA)      Predicted  SE(Predicted)

   1)    2.01807       0.08168       0.49552       0.35981
   2)    1.12472       0.04476       0.88911       0.36098

/* nag_glm_predict (g02gpc) Example Program.
 *
 * CLL6I261D/CLL6I261DL Version.
 *
 * Copyright 2017 Numerical Algorithms Group.
 *
 * Mark 26.1, 2017.
 */
/* Pre-processor includes */
#include <stdio.h>
#include <math.h>
#include <ctype.h>
#include <nag.h>
#include <nag_stdlib.h>
#include <nagg02.h>

#define T_X(I, J) t_x[(I) *t_tdx + J]
#define X(I, J)   x[(I) *tdx + J]
int main(void)
{
  /* Integer scalar and array declarations */
  Integer i, ip, j, m, n, t_n, tdx, t_tdx, print_iter;
  Integer exit_status = 0, tdv, rank, lx, lt_x, lv;
  Integer *sx = 0;
  /* nAG structures */
  Nag_Link link;
  Nag_IncludeMean mean;
  Nag_Boolean vfobs, weight, t_weight, ioffset, t_ioffset;
  Nag_Distributions errfn;
  NagError fail;
  /* Character scalar and array declarations */
  char sioffset[40], st_ioffset[40], sweight[40], st_weight[40];
  char slink[40], smean[40], svfobs[40];
  /* Double scalar and array declarations */
  double rss, scale, ex_power, df;
  double *b = 0, *cov = 0, *eta = 0, *offset = 0, *t_offset = 0;
  double *pred = 0, *se = 0, *seeta = 0, *sepred = 0, *binom_t = 0;
  double *v = 0, *wt = 0, *x = 0, *y = 0, *t_x = 0, *t_wt = 0;
  /* Set control parameters */
  double eps = 0.000001;
  double tol = 0.00005;
  Integer max_iter = 10;

  /* Initialize the error structure */
  INIT_FAIL(fail);

  printf("nag_glm_predict (g02gpc) Example Program Results\n");

  /* Skip headings in data file */
  scanf("%*[^\n] ");
  scanf("%*[^\n] ");
  /* Read in training data for model that will be used for prediction */
  scanf("%39s %39s %39s %39s %ld %ld %lf %ld"
        "%*[^\n] ", slink, smean, st_ioffset, st_weight, &t_n, &m, &scale,
        &print_iter);
  /*
   * nag_enum_name_to_value (x04nac).
   * Converts nAG enum member name to value
   */
  link = (Nag_Link) nag_enum_name_to_value(slink);
  mean = (Nag_IncludeMean) nag_enum_name_to_value(smean);
  t_ioffset = (Nag_Boolean) nag_enum_name_to_value(st_ioffset);
  t_weight = (Nag_Boolean) nag_enum_name_to_value(st_weight);

  t_tdx = m;
  lt_x = t_tdx * t_n;

  /* Allocate memory */
  if (t_weight) {
    if (!(t_wt = nAG_ALLOC(t_n, double)))
    {
      printf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }
  }
  if (t_ioffset) {
    if (!(t_offset = nAG_ALLOC(t_n, double)))
    {
      printf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }
  }
  if (!(t_x = nAG_ALLOC(lt_x, double)) ||
      !(y = nAG_ALLOC(t_n, double)) || !(sx = nAG_ALLOC(m, Integer)))
  {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }

  /* Read in the data */
  for (i = 0; i < t_n; i++) {
    for (j = 0; j < m; j++)
      scanf("%lf", &T_X(i, j));
    scanf("%lf", &y[i]);
    if (t_ioffset)
      scanf("%lf", &t_offset[i]);
    if (t_weight)
      scanf("%lf", &t_wt[i]);
    scanf("%*[^\n] ");
  }

  for (j = 0; j < m; j++)
    scanf("%ld%*[^\n] ", &sx[j]);

  if (link == Nag_Expo)
    scanf("%lf%*[^\n] ", &ex_power);
  else
    ex_power = 0.0;

  /* Calculate ip */
  ip = 0;
  for (j = 0; j < m; j++)
    if (sx[j] > 0)
      ip++;
  if (mean == Nag_MeanInclude)
    ip++;

  tdv = ip + 6;
  lv = tdv * t_n;

  if (!(b = nAG_ALLOC(ip, double)) ||
      !(v = nAG_ALLOC(lv, double)) ||
      !(se = nAG_ALLOC(ip, double)) ||
      !(cov = nAG_ALLOC(ip * (ip + 1) / 2, double)))
  {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }

  /* Call nag_glm_normal (g02gac) to fit model to training data */
  nag_glm_normal(link, mean, t_n, t_x, t_tdx, m, sx, ip, y, t_wt, t_offset,
                 &scale, ex_power, &rss, &df, b, &rank, se, cov, v, tdv,
                 tol, max_iter, print_iter, "", eps, &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_glm_normal (g02gac).\n%s\n", fail.message);
    exit_status = 1;
    goto END;
  }

  /* Display parameter estimates for training data */
  printf("\nResidual sum of squares = %12.4g, Degrees of freedom = %2f\n\n",
         rss, df);
  printf("       Estimate     Standard error\n\n");
  for (i = 0; i < ip; i++)
    printf(" %14.4f %14.4f\n", b[i], se[i]);
  printf("\n");

  /* Skip second lot of headings in data file */
  scanf("%*[^\n] ");

  /* Read in data to predict from and check array sizes */
  scanf("%ld %39s %39s %39s%*[^\n] ", &n, svfobs, sioffset,
        sweight);
  /*
   * nag_enum_name_to_value (x04nac).
   * Converts nAG enum member name to value
   */
  vfobs = (Nag_Boolean) nag_enum_name_to_value(svfobs);
  ioffset = (Nag_Boolean) nag_enum_name_to_value(sioffset);
  weight = (Nag_Boolean) nag_enum_name_to_value(sweight);

  if (weight) {
    if (!(wt = nAG_ALLOC(n, double)))
    {
      printf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }
  }
  if (ioffset) {
    if (!(offset = nAG_ALLOC(n, double)))
    {
      printf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }
  }

  tdx = m;
  lx = tdx * n;

  if (!(x = nAG_ALLOC(lx, double)) ||
      !(eta = nAG_ALLOC(n, double)) ||
      !(seeta = nAG_ALLOC(n, double)) ||
      !(pred = nAG_ALLOC(n, double)) || !(sepred = nAG_ALLOC(n, double)))
  {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }
  for (i = 0; i < n; i++) {
    for (j = 0; j < m; j++)
      scanf("%lf", &X(i, j));
    if (offset)
      scanf("%lf", &offset[i]);
    if (weight)
      scanf("%lf", &wt[i]);
    scanf("%*[^\n] ");
  }

  /* Using nag_glm_normal (g02gac) to fit training model, so error structure
     is normal */
  errfn = Nag_Normal;

  /* Call nag_glm_predict (g02gpc) to calculate predictions */
  nag_glm_predict(errfn, link, mean, n, x, tdx, m, sx, ip, binom_t, offset,
                  wt, scale, ex_power, b, cov, vfobs, eta, seeta, pred,
                  sepred, &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_glm_predict (g02gpc).\n%s\n", fail.message);
    exit_status = 1;
    goto END;
  }

  /* Display predicted values */
  printf("   I      ETA         SE(ETA)      Predicted  SE(Predicted)\n");
  printf("\n");
  for (i = 0; i < n; i++) {
    printf(" %3ld) %10.5f    %10.5f    %10.5f    %10.5f\n", i + 1,
           eta[i], seeta[i], pred[i], sepred[i]);
  }

END:
  nAG_FREE(t_wt);
  nAG_FREE(t_x);
  nAG_FREE(y);
  nAG_FREE(sx);
  nAG_FREE(b);
  nAG_FREE(v);
  nAG_FREE(se);
  nAG_FREE(cov);
  nAG_FREE(wt);
  nAG_FREE(x);
  nAG_FREE(offset);
  nAG_FREE(eta);
  nAG_FREE(seeta);
  nAG_FREE(pred);
  nAG_FREE(sepred);

  return exit_status;
}