ユーザ提供のRidge回帰パラメータを用いたRidge回帰

nag_regsn_ridge (g02kbc) Example Program Data
20 3 16 5 Nag_OrigPara : n, m, lh, lpec, wantb
1 1 1 : isx
Nag_LOOCV Nag_GCV Nag_EUV Nag_FPE Nag_BIC : pec
  19.5  43.1  29.1  11.9
  24.7  49.8  28.2  22.8
  30.7  51.9  37.0  18.7
  29.8  54.3  31.1  20.1
  19.1  42.2  30.9  12.9
  25.6  53.9  23.7  21.7
  31.4  58.5  27.6  27.1
  27.9  52.1  30.6  25.4
  22.1  49.9  23.2  21.3
  25.5  53.5  24.8  19.3
  31.1  56.6  30.0  25.4
  30.4  56.7  28.3  27.2
  18.7  46.5  23.0  11.7
  19.7  44.2  28.6  17.8
  14.6  42.7  21.3  12.8
  29.5  54.4  30.1  23.9
  27.7  55.3  25.7  22.6
  30.2  58.6  24.6  25.4
  22.7  48.2  27.1  14.8
  25.2  51.0  27.5  21.1 : End of observations
0.0 0.002 0.004 0.006 0.008 0.010 0.012 0.014 0.016 0.018 0.020 0.022
0.024 0.026 0.028 0.030 : Ridge co-efficients

nag_regsn_ridge (g02kbc) Example Program Results
Number of parameters used =          4
Effective number of parameters (NEP):
   Ridge   
   Coeff.  NEP
     0.0000     4.0000
     0.0020     3.2634
     0.0040     3.1475
     0.0060     3.0987
     0.0080     3.0709
     0.0100     3.0523
     0.0120     3.0386
     0.0140     3.0278
     0.0160     3.0189
     0.0180     3.0112
     0.0200     3.0045
     0.0220     2.9984
     0.0240     2.9928
     0.0260     2.9876
     0.0280     2.9828
     0.0300     2.9782

Parameter Estimates (Original scalings)
  Ridge  
   Coeff.    Intercept           1          2          3 
    0.0000  117.0847     4.3341    -2.8568    -2.1861 
    0.0020   22.2748     1.4644    -0.4012    -0.6738 
    0.0040    7.7209     1.0229    -0.0242    -0.4408 
    0.0060    1.8363     0.8437     0.1282    -0.3460 
    0.0080   -1.3396     0.7465     0.2105    -0.2944 
    0.0100   -3.3219     0.6853     0.2618    -0.2619 
    0.0120   -4.6734     0.6432     0.2968    -0.2393 
    0.0140   -5.6511     0.6125     0.3222    -0.2228 
    0.0160   -6.3891     0.5890     0.3413    -0.2100 
    0.0180   -6.9642     0.5704     0.3562    -0.1999 
    0.0200   -7.4236     0.5554     0.3681    -0.1916 
    0.0220   -7.7978     0.5429     0.3779    -0.1847 
    0.0240   -8.1075     0.5323     0.3859    -0.1788 
    0.0260   -8.3673     0.5233     0.3926    -0.1737 
    0.0280   -8.5874     0.5155     0.3984    -0.1693 
    0.0300   -8.7758     0.5086     0.4033    -0.1653 

Variance Inflation Factors
  Ridge  
  Coeff.           1          2          3 
    0.0000  708.8429   564.3434   104.6060 
    0.0020   50.5592    40.4483     8.2797 
    0.0040   16.9816    13.7247     3.3628 
    0.0060    8.5033     6.9764     2.1185 
    0.0080    5.1472     4.3046     1.6238 
    0.0100    3.4855     2.9813     1.3770 
    0.0120    2.5434     2.2306     1.2356 
    0.0140    1.9581     1.7640     1.1463 
    0.0160    1.5698     1.4541     1.0859 
    0.0180    1.2990     1.2377     1.0428 
    0.0200    1.1026     1.0805     1.0105 
    0.0220    0.9556     0.9627     0.9855 
    0.0240    0.8427     0.8721     0.9655 
    0.0260    0.7541     0.8007     0.9491 
    0.0280    0.6832     0.7435     0.9353 
    0.0300    0.6257     0.6969     0.9235 

Prediction error criterion
  Ridge  
  Coeff.           1          2          3          4          5

    0.0000    8.0368     7.6879     6.1503     7.3804     8.6052
    0.0020    7.5464     7.4238     6.2124     7.2261     8.2355
    0.0040    7.5575     7.4520     6.2793     7.2675     8.2515
    0.0060    7.5656     7.4668     6.3100     7.2876     8.2611
    0.0080    7.5701     7.4749     6.3272     7.2987     8.2661
    0.0100    7.5723     7.4796     6.3381     7.3053     8.2685
    0.0120    7.5732     7.4823     6.3455     7.3095     8.2695
    0.0140    7.5734     7.4838     6.3508     7.3122     8.2696
    0.0160    7.5731     7.4845     6.3548     7.3140     8.2691
    0.0180    7.5724     7.4848     6.3578     7.3151     8.2683
    0.0200    7.5715     7.4847     6.3603     7.3158     8.2671
    0.0220    7.5705     7.4843     6.3623     7.3161     8.2659
    0.0240    7.5694     7.4838     6.3639     7.3162     8.2645
    0.0260    7.5682     7.4832     6.3654     7.3162     8.2630
    0.0280    7.5669     7.4825     6.3666     7.3161     8.2615
    0.0300    7.5657     7.4818     6.3677     7.3159     8.2600

Key:
      1  Leave one out cross-validation
      2  Generalized cross-validation
      3  Unbiased estimate of variance
      4  Final prediction error
      5  Bayesian information criterion

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

int main(void)
{
  /*Integer scalar and array declarations */
  Integer exit_status = 0;
  Integer i, ip, ip1, j, lh, lpec, m, n, pl;
  Integer pdb, pdpe, pdvf, pdx;
  Integer *isx = 0;
  /*Double scalar and array declarations */
  double *b = 0, *h = 0, *nep = 0, *pe = 0, *vf = 0, *x = 0, *y = 0;
  /*Character scalar and array declarations */
  char spec[40], swantb[40];
  /*nAG Types */
  Nag_OrderType order;
  Nag_ParaOption wantb;
  Nag_VIFOption wantvf;
  Nag_PredictError *pec = 0;
  NagError fail;

  INIT_FAIL(fail);

  printf("%s\n", "nag_regsn_ridge (g02kbc) Example Program Results");
  /* Skip heading in data file */
  scanf("%*[^\n] ");
  /* Read in the problem size information */
  scanf("%ld%ld%ld%ld%39s%*[^\n] ",
        &n, &m, &lh, &lpec, swantb);
  wantb = (Nag_ParaOption) nag_enum_name_to_value(swantb);
  if (!(isx = nAG_ALLOC(m, Integer)))
  {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }
  /* Read in the ISX flags */
  for (i = 0; i < m; i++)
    scanf("%ld ", &isx[i]);
  scanf("%*[^\n] ");
  /* Total number of variables */
  ip = 0;
  for (j = 0; j < m; j++) {
    if (isx[j] == 1)
      ip = ip + 1;
  }
#ifdef nAG_COLUMN_MAJOR
  pdb = ip + 1;
#define B(I, J)  b[(J-1)*pdb + I-1]
  pdpe = lpec;
#define PE(I, J) pe[(J-1)*pdpe + I-1]
  pdvf = ip;
#define VF(I, J) vf[(J-1)*pdvf + I-1]
  pdx = n;
#define X(I, J)  x[(J-1)*pdx + I-1]
  order = Nag_ColMajor;
#else
  pdb = lh;
#define B(I, J)  b[(I-1)*pdb + J-1]
  pdpe = lh;
#define PE(I, J) pe[(I-1)*pdpe + J-1]
  pdvf = lh;
#define VF(I, J) vf[(I-1)*pdvf + J-1]
  pdx = m;
#define X(I, J)  x[(I-1)*pdx + J-1]
  order = Nag_RowMajor;
#endif
  if (!(b = nAG_ALLOC(pdb * (order == Nag_RowMajor ? (ip + 1) : lh), double))
      || !(h = nAG_ALLOC(lh, double)) || !(nep = nAG_ALLOC(lh, double))
      || !(pe = nAG_ALLOC(pdpe * (order == Nag_RowMajor ? lpec : lh), double))
      || !(vf = nAG_ALLOC(pdvf * (order == Nag_RowMajor ? ip : lh), double))
      || !(x = nAG_ALLOC(pdx * (order == Nag_RowMajor ? n : m), double))
      || !(y = nAG_ALLOC(n, double))
      || !(pec = nAG_ALLOC(lpec, Nag_PredictError)))
  {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }
  /* Read in the data */
  if (lpec > 0) {
    for (i = 0; i < lpec; i++) {
      scanf("%39s ", spec);
      pec[i] = (Nag_PredictError) nag_enum_name_to_value(spec);
    }
    scanf("%*[^\n] ");
  }
  for (i = 1; i <= n; i++) {
    for (j = 1; j <= m; j++)
      scanf("%lf ", &X(i, j));
    scanf("%lf ", &y[i - 1]);
  }
  scanf("%*[^\n] ");
  /*  Read in the ridge coefficients */
  for (i = 0; i < lh; i++)
    scanf("%lf ", &h[i]);
  scanf("%*[^\n] ");
  /* Output the variance inflation factors and parameter estimates */
  wantvf = Nag_WantVIF;
  /* Run the analysis */
  /*
   * nag_regsn_ridge (g02kbc)
   * Ridge regression
   */
  nag_regsn_ridge(order, n, m, x, pdx, isx, ip, y, lh, h, nep, wantb, b, pdb,
                  wantvf, vf, pdvf, lpec, pec, pe, pdpe, &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_regsn_ridge (g02kbc).\n%s\n", fail.message);
    exit_status = 1;
    goto END;
  }
  /* Output results */
  ip1 = ip - 1;
  /* Summaries */
  printf("%s%10ld\n", "Number of parameters used = ", ip + 1);
  printf("%s\n", "Effective number of parameters (NEP):");
  printf("%s\n", "   Ridge   ");
  printf("%s%s\n", "   Coeff.  ", "NEP");
  for (i = 1; i <= lh; i++)
    printf(" %10.4f %10.4f\n", h[i - 1], nep[i - 1]);
  /* Parameter estimates */
  if (wantb != Nag_NoPara) {
    printf("\n");
    if (wantb == Nag_OrigPara) {
      printf("%s\n", "Parameter Estimates (Original scalings)");
    }
    else {
      printf("%s\n", "Parameter Estimates (Standarised)");
    }
    pl = MIN(ip, 4);
    printf("%s\n", "  Ridge  ");
    printf("%s  ", "   Coeff.  ");
    printf("%s ", "Intercept ");
    for (i = 1; i <= pl; i++)
      printf("%10ld%s", i, i % 4 ? " " : "\n");
    printf("\n");
    if (pl < ip1) {
      for (i = pl + 1; i <= ip1; i++)
        printf("%10ld%s", i, i % 5 ? " " : "\n");
      printf("\n");
    }
    pl = MIN(ip + 1, 5);
    for (i = 1; i <= lh; i++) {
      printf("%10.4f", h[i - 1]);
      for (j = 1; j <= pl; j++)
        printf("%10.4f%s", B(j, i), j % 5 ? " " : "\n");
      printf("\n");
      if (pl < ip) {
        for (j = pl + 1; j <= ip; j++)
          printf("%10.4f%s", B(j, i), j % 5 ? " " : "\n");
        printf("\n");
      }
    }
  }
  /* Variance inflation factors */
  if (wantvf != Nag_NoVIF) {
    printf("\n");
    printf("%s\n", "Variance Inflation Factors");
    pl = MIN(ip, 5);
    printf("%s\n", "  Ridge  ");
    printf("%s", "  Coeff.  ");
    for (i = 1; i <= pl; i++)
      printf("%10ld%s", i, i % 5 ? " " : "\n");
    printf("\n");
    if (pl < ip) {
      for (i = pl + 1; i <= ip; i++)
        printf("%10ld%s", i, i % 5 ? " " : "\n");
      printf("\n");
    }
    for (i = 1; i <= lh; i++) {
      printf("%10.4f", h[i - 1]);
      for (j = 1; j <= pl; j++)
        printf("%10.4f%s", VF(j, i), j % 5 ? " " : "\n");
      printf("\n");
      if (pl < ip) {
        for (j = pl + 1; j <= ip; j++)
          printf("%10.4f%s", VF(j, i), j % 5 ? " " : "\n");
        printf("\n");
      }
    }
  }
  /* Prediction error criterion */
  if (lpec > 0) {
    printf("\n");
    printf("%s\n", "Prediction error criterion");
    pl = MIN(lpec, 5);
    printf("%s\n", "  Ridge  ");
    printf("%s", "  Coeff.  ");
    for (i = 1; i <= pl; i++)
      printf("%10ld%s", i, i % 5 ? " " : "\n");
    printf("\n");
    if (pl < lpec) {
      for (i = pl + 1; i <= lpec; i++)
        printf("%10ld%s", i, i % 5 ? " " : "\n");
      printf("\n");
    }
    for (i = 1; i <= lh; i++) {
      printf("%10.4f", h[i - 1]);
      for (j = 1; j <= pl; j++)
        printf("%10.4f%s", PE(j, i), j % 5 ? " " : "\n");
      if (pl < ip) {
        for (j = pl + 1; j <= ip; j++)
          printf("%10.4f%s", PE(j, i), j % 5 ? " " : "\n");
      }
    }
    printf("\n");
    printf("%s\n", "Key:");
    for (i = 1; i <= lpec; i++) {
      if (pec[i - 1] == Nag_LOOCV) {
        printf("  %5ld  Leave one out cross-validation\n", i);
      }
      else if (pec[i - 1] == Nag_GCV) {
        printf("  %5ld  Generalized cross-validation\n", i);
      }
      else if (pec[i - 1] == Nag_EUV) {
        printf("  %5ld  Unbiased estimate of variance\n", i);
      }
      else if (pec[i - 1] == Nag_FPE) {
        printf("  %5ld  Final prediction error\n", i);
      }
      else if (pec[i - 1] == Nag_BIC) {
        printf("  %5ld  Bayesian information criterion\n", i);
      }
    }
  }

END:
  nAG_FREE(b);
  nAG_FREE(h);
  nAG_FREE(nep);
  nAG_FREE(pe);
  nAG_FREE(vf);
  nAG_FREE(x);
  nAG_FREE(y);
  nAG_FREE(isx);
  nAG_FREE(pec);

  return exit_status;
}