(* Multivariate Analysis Package - Multiple Linear Regression Module
   Copyright 1985 Douglas L. Anderton.  This program may be freely
   circulated so long as it is not sold for profit and any charge does
   not exceed costs of reproduction *)

Program Regress(Input,Output);
Const
   N=30;
Type
   SUBS=1..N;
   RVEC = Array [SUBS] Of Real;
   NBYN = Array [SUBS] Of RVEC;
   IVEC = Array [SUBS] Of Integer;
   S8 = Array [SUBS] Of String[8];
Var
   dfile, ofile : Text;
   sel : IVEC;
   mean, stdev, b, beta : RVEC;
   varn : S8;
   cor : NBYN;
   i, j, nv, dv, ot : SUBS;
   det, rsq, multr, ftest, tss, ssr, sse, seb, ttest, nc, df1, df2 : Real;
   error : Boolean;

Procedure openfiles(Var dfile, ofile:Text; Var nv, dv, ot:SUBS);
Var
   dfl, ofl:String[12];
Begin
   ClrScr; Writeln(' *** REGRESS: MULTIPLE LINEAR REGRESSION ***');
   Writeln; Write('Name of the CORREL file? ');
   Readln(dfl); Assign(dfile,dfl); Reset(dfile);
   Write('Name of the output file? ');
   Readln(ofl); Assign(ofile,ofl); Rewrite(ofile);
   ot:=0;
   If (ofl='CON:') Or (ofl='con:') Then ot:=1;
   If (ofl='LST:') Or (ofl='lst:') Then ot:=2;
   If (ot=2) Then
      Begin
      Writeln(ofile,'Multivariate Analysis Package (1.6) - ',
         'Program: REGRESS, Datafile: ',dfl); Writeln(ofile);
      End;
   Write('How many variables in CORREL matrix? '); Readln(nv);
   Write('Number of variables to use in REGRESS? '); Readln(dv);
   End; (* Of openfiles *)

Procedure selectvar2(Var sel:IVEC; Var dv:SUBS);
Var
   i:SUBS;
Begin
   Writeln;
   Write('Column number for dependent variable? ');
   Readln(sel[dv]);
   For i:=1 To dv-1 Do
     Begin
     Write('Column number for independent variable ',i,'? '); Readln(sel[i]);
     End;
   End; (* Of selectvar *)

{$I GETCOR.LIB}
{$I REGGAUSS.LIB}

Begin (* main *)
   openfiles(dfile, ofile, nv, dv, ot);
   selectvar2(sel, dv);
   (* read correlation matrix matrix *)
   getcor(cor, mean, stdev, varn, nc, sel, nv, dv, dfile);
   (* solve R(xx)B(xy)=R(xy) with gaussian elimination *)
   Gausslice(dv-1,cor,cor[dv],b,error);
   If error Then Begin
     Writeln(ofile,'Error in Gausslice - ill conditioned correlation matrix');
     bdos(0);
     End;
   (* multiply for unstandardized coefficients *)
   For i:=1 To dv-1 Do beta[i] := b[i]*(stdev[dv]/stdev[i]);
   beta[dv] := mean[dv];
   For i:=1 To dv-1 Do beta[dv] := beta[dv]-(beta[i]*mean[i]);
   (* significance tests and anova *)
   rsq := 0.0; For i:=1 To dv-1 Do rsq := rsq + (b[i]*cor[dv,i]);
   multr := sqrt(rsq);
   df1 := dv-1; df2 := nc-dv; ftest := (rsq*df2)/((-rsq+1.0)*df1);
   tss := sqr(stdev[dv])*(nc-1.0); ssr := rsq*tss; sse := tss-ssr;
   (* output *)
   ClrScr;
   Writeln(ofile,'Summary of Regression on ',varn[dv]); Writeln(ofile);
   Writeln(ofile,'R-Square: ',rsq:8:6,'  Multiple R:  ',multr:8:6);
   Writeln(ofile,'F test:   ',ftest:8:4,'  df1 = ',df1:3:0,' df2 = ',df2:10:0);
   Writeln(ofile);
   Writeln(ofile,'          Sums of Squares       Mean Squares           df');
   Writeln(ofile,'Total      ',tss:16:4,' ':22,(nc-1):10:0);
   Writeln(ofile,'Regression ',ssr:16:4,'   ',ssr/df1:16:5,'   ',df1:10:0);
   Writeln(ofile,'Error      ',sse:16:4,'   ',sse/df2:16:5,'   ',df2:10:0);
   Writeln(ofile);
   If ot=1 Then
      Begin
      Write('- Press any key to continue -'); While Not KeyPressed Do;
      ClrScr;
      End;
   Writeln(ofile,'Coefficents of Regression on ',varn[dv]);
   Write(ofile,' ':9,'Beta':12,'Std Beta':12,'Std Error':13,'t test':10);
   Writeln(ofile,'df':12,'corr':8);Writeln(ofile);
   Writeln(ofile,'Intercept',beta[dv]:13:6);
   { put inverse diagonals in unused means vector }
   InverseDiagonals(dv-1,cor,mean);
   For i:=1 To dv-1 Do
      Begin
      seb := sse/df2; rsq := 1.0-(1.0/mean[i]);
      ssr := sqr(stdev[i])*(nc-1.0); seb := sqrt(seb/(ssr*(-rsq+1.0)));
      ttest := beta[i]/seb;
      Write(ofile,varn[i]:8,' ',beta[i]:13:6,' ',b[i]:11:6,'  ',seb:11:6);
      Writeln(ofile,' ',ttest:10:4,' ',df2:10:0,' ',cor[dv,i]:6:4);
      End;
   Writeln(ofile);
   If ot=1 Then
      Begin
      Write('- Press any key to exit -'); While Not KeyPressed Do;
      ClrScr;
      End;
   Close(ofile); Close(dfile);
   Assign(dfile,'MAPSTAT.COM'); Execute(dfile);
End.