I agree that this test doesn't is not calculated by PROC GLM.
Here is my approach. Since the maksimum of the log-likelihood is just -n*log(determinant(E))/2+constant, where E is the error-matrix, it is just to let proc glm output the error matrix into a dataset. Do this for both model, and then calculate the likelihood ratio value manually. The pvalue is calculated by using that the 2 x difference in log-likelihood is chi-square distributed. I think there is a slightly better approach where it is used that some transformation of the ratio is exact F-distributed, but I cant find this in my notes from univerversity.
here is an example. It is neccessary to take determinant, so I therefor use the matrix-package here (Matrix-Package). Alternatively it can be done with PROC IML.
Of some reason this site will not show my code correct, so I therefore also attach it in a separate file
data simulation;
do gp1=1 to 4;
do gp2=1 to 4;
do i=1 to 5;
y1=rannor(-1);
y2=rannor(-1);
y3=rannor(-1);
output;
end;
end;
end;
run;
ods output errorsscp=sscp0;
proc glm data=simulation;
class gp1 gp2;
model y1 y2 y3 = gp1 gp2;
manova h=gp1 gp2/printe ;
run;
quit;
ods output errorsscp=sscp1;
proc glm data=simulation;
model y1 y2 y3 = ;
manova /printe;
run;
quit;
data _NULL_; array sigma0{3,3} _temporary_ ; array sigma1{3,3} _temporary ; dsid0= open('sscp0'); dsid1= open('sscp1'); do i=1 to 3; rc= fetchobs(dsid0,i); rc= fetchobs(dsid1,i); do j=1 to 3; sigma0[i,j]=getvarn(dsid0,2+j); sigma1[i,j]=getvarn(dsid1,2+j); end; end; call show(sigma0); put; call show(sigma1); d0=determinant(sigma1); d1=determinant(sigma0); likelihoodratio=80*(log(d0)-log(d1));/*n=80 in this example*/
put d0= d1= likelihoodratio=; pvalue=sdf('chisquare',likelihoodratio,3*( ((4-1)+(4-1)+1)-1));/*third argument is the difference in parameters*/; put pvalue= pvalue6.4; run;
... View more
