/*Example 13.4.1  Intraclass Correlation Test */
proc iml;
x1={72 60 56 41 32 30 39 42 37 33 32 63 54 47
    91 56 79 81 78 46 39 32 60 35 39 50 43 48};
x2={66 53 57 29 32 35 39 43 40 29 30 45 46 51
    79 68 65 80 55 38 35 30 50 37 36 34 37 54};
x3={76 66 64 36 35 34 31 31 31 27 34 74 60 52
   100 47 70 68 67 37 34 30 67 48 39 37 39 57};
x4={77 63 58 38 36 26 27 25 25 36 28 63 52 45
    75 50 61 58 60 38 37 32 54 39 31 40 50 43};
x=x1`||x2`||x3`||x4`;
n=28;p=4;
xb=(x`*J(n,n,1/n))`;
s=(x-xb) `*(x-xb)/(n-1);
sum=0;
do i=1 to p;
do j=i+1 to p;
sum=sum+s[i,j];
end;
end;
R=2/(p*(p-1))*sum*p/trace(s);
c=J(p,1,1);
lambda=det(p*s)/(c`*s*c*exp((p-1)*log((p*trace(s)-c`*s*c)/(p-1))));
Q=-(n-1-p*(p+1)**2*(2*p-3)/(6*(p-1)*(P**2+p-4)))*log(lambda);
g=p*(p+1)/2-2;
chi=cinv(1-0.05,g);
print  R lambda, Q chi;
if Q>chi then print "Q > CHI: reject";
         else print "Q < CHI: accept";
quit;
                  R          LAMBDA
          0.8584975       0.4634634
                  Q             CHI
          19.428637       15.507313
                 Q > CHI: reject