Program Version: 1.1
Date: June 18, 1996}

randomaccess;

{Assign Starting Values}
     maxv:=20;                 {max. no of independent variables}
     redim var(maxv);
     step:=20;               {no. of step}

{Get parameter Function}
     Function GetBracketsForSearch(k)
     begin
      dim k(2);
      k(1):=0;
      k(2):=0.02;
      iret:=0;
      if DisplayNumericInputBox ('Range of the shrinkage parameter k',
       'Lower limit:|Upper limit:', k)=0 then begin
        iret:=1;
        goto thatsit;
      end;
      if k(1)<0.0 then begin
       DisplayMessageBox (MB_ICONSTOP, 'Invalid lower limit',
        'Lower limit value must >= 0.0');
       iret:=1;
      end;
      if k(2)<k(1) then begin
       DisplayMessageBox (MB_ICONSTOP, 'Invalid upper lambda',
        'Upper limit value must > Lower limit value');
       iret:=1;
      end;
      if (abs(k(1))>1) or (abs(k(2))>1) then begin
       DisplayMessageBox (MB_ICONSTOP, 'Invalid limit values',
        'User-specified values too large.');
       iret:=1;
      end;
thatsit:;
      GetBracketsForSearch:=iret;
     end;

{Ridge Regression Subroutines}
     Sub ridge (x,y,b,ixtx,hat,shrink,nv,nc);
     begin
      dim b(nv);
      dim x(nc,nv);
      dim xt(nv,nc);
      dim xtx(nv,nv);
      dim ixtx(nv,nv);
      dim hat(nc,nc);
      dim r1(nv,nc);
      dim r2(nc,nv);
      dim y(nc);

      MatrixTranspose (x, xt);
      MatrixMultiply (xt,x,xtx);
      for i:=1 to nv do xtx(i,i):=xtx(i,i)+shrink;
      MatrixInverse (xtx, ixtx);
      MatrixMultiply (ixtx,xt,r1);
      MatrixMultiply (r1,y,b);
      MatrixMultiply (x,ixtx,r2);
      MatrixMultiply (r2,xt,hat);
     end;

{Convert the scaling estimates to un-scaling data estimates Subroutine}
     Sub convert (tm,si,b,ymean,n,bt);
     begin
      dim tm(n);
      dim si(n);
      dim b(n);
      dim bt(n+1);

      temp:=0;
      for i:=1 to n do bt(i+1):=b(i)/si(i);
      for i:=1 to n do temp:=temp+b(i)*tm(i)/si(i);
      bt(1):=ymean-temp;
     end;

{Choose Variables}
     if SelectVariables2 ('Select Variables',
           1, 10, dv, i,'Dependent variable:',
           2, maxv-1-1, var, n, 'Independent variables:')=0 then stop;

{Allocate Memory}
     redim k(2);
     redim b(n);
     redim x(ncases,n);
     redim y(ncases);
     redim mat(n,n);
     redim hat(ncases,ncases);
     redim predicted(ncases);
     redim residuals(ncases);
     redim     result(step+1,n+5);
     redim tm(n);
     redim si(n);
     redim bt(n+1);


{ get values ("brackets") for k}
     if GetBracketsForSearch(k) then stop;
     diff:=(k(2)-k(1))/step;

{Copy Data to Array}
     for i:=1 to n do MatrixCopy (data,1, var(i), ncases,1,x,1,i);
     MatrixCopy (data,1,dv,ncases,1,y,1,1);

{Center and Scale the regressor}
     for i:=1 to n do begin
          for j:=1 to ncases do tm(i):=tm(i)+x(j,i)/ncases;
          for j:=1 to ncases do si(i):=si(i)+(x(j,i)-tm(i))**2;
          si(i):=si(i)**0.5;
          for j:=1 to ncases do
               x(j,i):=(x(j,i)-tm(i))/si(i);
     end;

{Compute the y mean}
     ValMean (y,1,ncases,ymean);

{Compute MSE & Trace Statistics from OLS Regression}
     ridge (x,y,b,mat,hat,0.0,n,ncases);
     convert(tm,si,b,ymean,n,bt);
     for i:=1 to ncases do begin;
      predicted(i):=ymean;
      for j:=1 to n do
        predicted(i):=predicted(i)+b(j)*x(i,j);
      residuals(i):=y(i)-predicted(i);
      sse:=sse+residuals(i)^2;
     end;
     dfe:=ncases-n-1;
     mse:=sse/dfe;
     MatrixTrace (hat, trh);
     Ck:=(sse/mse)-ncases+2+2*trh;
     GCV:=sse/((ncases-1-trh)**2);

{Put the result into the result matrix}
     result(1,1):=0;
     result(1,2):=Ck;
     result(1,3):=GCV;
     result(1,4):=trh;
     for i:=1 to n+1 do result(1,i+4):=bt(i);

{Iterative Loop}
     for count:=1 to step do begin
     sse:=0;
     shrink:=k(1)+count*diff;
     ridge (x,y,b,mat,hat,shrink,n,ncases);
     convert(tm,si,b,ymean,n,bt);

{Compute SSE & Trace Statistics}
     for i:=1 to ncases do begin;
      predicted(i):=ymean;
      for j:=1 to n do
        predicted(i):=predicted(i)+b(j)*x(i,j);
      residuals(i):=y(i)-predicted(i);
      sse:=sse+residuals(i)^2;
     end;
     MatrixTrace (hat, trh);
     Ck:=(sse/mse)-ncases+2+2*trh;
     GCV:=sse/((ncases-1-trh)**2);

{Put the result into the result matrix}
     result(count+1,1):=shrink;
     result(count+1,2):=Ck;
     result(count+1,3):=GCV;
     result(count+1,4):=trh;
     for i:=1 to n+1 do result(count+1,i+4):=bt(i);

     end;

{Display Ridge Regression Trace Statistics}
     line01$:='Trace Statistics and The Coefficients for different k |';
     line02$:='Note: The min. Ck row is in red colour';
      kname01$:='k | Ck | GCV | df-trace | Intercept';
     for i:=1 to n do kname01$:=kname01$ + ' | ' + varname(var(i));
     shandle:=NewScrollsheet (step+1, n+5, result, line01$+line02$, '', kname01$);

{Search and Highlight the Min. Ck row}
     redim vk(step+1);
     redim kk(step+1);
     MatrixGetColumn (result,2,vk);
     MatrixGetColumn (result,1,kk);
     ValMin (vk,1,step+1,min_ck);
     for i:=1 to step+1 do begin
          if vk(i)=min_ck then begin
               for j:=1 to n+3 do
               ScrollsheetSetHilite (shandle,i,j,1);
          end;
     end;

{Plot of Ck against k}
       Graph1:=Newgraph (LINEPLOT, 'Plot of Ck against k',
          'The value of Ck','The shrinkage parameter k', step+1, kk, vk);
       GraphSetPlotPointStyle (Graph1, 1, ON, P_filled_square, ?Size, ?Color);

{Plot of GCV against k}
        MatrixGetColumn (result,3,vk);
       Graph2:=Newgraph (LINEPLOT, 'Plot of GCV against k',
          'The value of GCV','The shrinkage parameter k', step+1, kk, vk);
       GraphSetPlotPointStyle (Graph2, 1, ON, P_filled_square, ?Size, ?Color);

{Plot of DF-Trace}
        MatrixGetColumn (result,4,vk);
       Graph2:=Newgraph (LINEPLOT, 'DF-Trace Plot',
          'DF value','The shrinkage parameter k', step+1, kk, vk);
       GraphSetPlotPointStyle (Graph2, 1, ON, P_filled_square, ?Size, ?Color);

2300 East 14th Street, Tulsa, OK 74104
Phone: (918) 749-1119; Fax: (918) 749-2217
e-mail: info@statsoft.com

©Copyright StatSoft, Inc., 1984-2004.
StatSoft, StatSoft logo, STATISTICA, SEWSS, SEDAS, Data Miner, SEPATH and GTrees are trademarks of StatSoft, Inc.