yydy323b.sas

高老师的经典多元统计分析教材源代码非常经典使用

/*  yydy323b.sas */

 proc iml;
   n1=20;n2=20;n3=20;
   n=n1+n2+n3;k=3; p=4;
   x1={
      260  75 40 18 , 200  72 34 17 , 240  87 45 18 , 170  65 39 17 ,
      270 110 39 24 , 205 130 34 23 , 190  69 27 15 , 200  46 45 15 ,
      250 117 21 20 , 200 107 28 20 , 225 130 36 11 , 210 125 26 17 ,
      170  64 31 14 , 270  76 33 13 , 190  60 34 16 , 280  81 20 18 ,
      310 119 25 15 , 270  57 31  8 , 250  67 31 14 , 260 135 39 29 } ;
   x2={
       310 122 30 21 ,310  60 35 18 ,190  40 27 15 ,225  65 34 16 ,
       170  65 37 16 ,210  82 31 17 ,280  67 37 18 ,210  38 36 17 ,
       280  65 30 23 ,200  76 40 17 ,200  76 39 20 ,280  94 26 11 ,
       190  60 33 17 ,295  55 30 16 ,270 125 24 21 ,280 120 32 18 ,
       240  62 32 20 ,280  69 29 20 ,370  70 30 20 ,280  40 37 17 };
   x3={
       320  64 39 17 ,260  59 37 11 ,360  88 28 26 ,295 100 36 12 ,
       270  65 32 21 ,380 114 36 21 ,240  55 42 10 ,260  55 34 20 ,
       260 110 29 20 ,295  73 33 21 ,240 114 38 18 ,310 103 32 18 ,
       330 112 21 11 ,345 127 24 20 ,250  62 22 16 ,260  59 21 19 ,
       225 100 34 30 ,345 120 36 18 ,360 107 25 23 ,250 117 36 16 };
   xx=x1//x2//x3;
   ln={[20] 1} ; lnn={[60] 1};
   x10=(ln*x1)/n1; print x10;
   x20=(ln*x2)/n1; print x20;
   x30=(ln*x3)/n1; print x30;
   xx0=(lnn*xx)/n; print xx0;
   mm1=i(n1)-j(n1,n1,1)/n1;
   mm=i(n)-j(n,n,1)/n;
   a1=x1`*mm1*x1; print a1;
   a2=x2`*mm1*x2; print a2;
   a3=x3`*mm1*x3; print a3;
   tt=xx`*mm*xx;  print tt;
   a=a1+a2+a3;   print a;
   da=det(a);
   dt=det(tt); a0=da/dt;
   print da dt a0;
   b=sqrt(a0);
   f=(n-k-p+1)*(1-b)/(b*p);
   df1=2*p;df2=2*(n-k-p+1);
   p0=1-probf(f,df1,df2);
   print f p0;
   f1=(tt[1,1]-a[1,1])*(n-k)/((k-1)*a[1,1]);
   p1=1-probf(f1,k-1,n-k);
   fa1=finv(0.99,k-1,n-k); fa5=finv(0.95,k-1,n-k);
   print fa1 fa5 f1 p1;
 run;