/****************************************************** * Wechslerr Adult Intelligence Scale subtest scores * * From Morrison, Appendix B * * * * These results are a bit different from what Morrison* * reports. Either there's a typo in the data or in * * the reporting of results....it's probably in the * * data (I double checked with what's printed), * ******************************************************/ data IQ; input sub senile $ information simiarlities arithmetic picture; a1=1; a2=0; if senile='n' then a2=1; datalines; 1 n 7 5 9 8 2 n 8 8 5 6 3 n 16 18 11 9 4 n 8 3 7 9 5 n 6 3 13 9 6 n 11 8 10 10 7 n 12 7 9 8 8 n 8 11 9 3 9 n 14 12 11 4 10 n 13 13 13 6 11 n 13 9 9 9 12 n 13 10 15 7 13 n 14 11 12 8 14 n 15 11 11 10 15 n 13 10 16 9 16 n 10 5 8 6 17 n 10 3 7 7 18 n 17 13 13 7 19 n 10 6 10 7 20 n 10 10 15 8 21 n 14 7 11 5 22 n 16 11 12 11 23 n 10 7 14 6 24 n 10 10 9 6 25 n 10 7 10 10 26 n 7 6 5 9 27 n 15 12 10 6 28 n 17 15 15 8 29 n 16 13 16 9 30 n 13 10 17 8 31 n 13 10 17 10 32 n 19 12 16 10 33 n 19 12 17 11 34 n 13 10 7 8 35 n 15 11 12 8 36 n 16 9 11 11 37 n 14 13 14 9 38 y 9 5 10 8 39 y 10 0 6 2 40 y 8 9 11 1 41 y 13 7 14 9 42 y 4 0 4 0 43 y 4 0 6 0 44 y 11 9 9 8 45 y 5 3 3 6 46 y 9 7 8 6 47 y 7 2 6 4 48 y 12 10 14 3 49 y 13 12 11 10 proc sort; by senile; run; proc means; by senile; var information simiarlities arithmetic picture; output out=xbar mean=x1 x2 x3 x4; run; proc print data=xbar; run; proc glm data=iq; class senile; model information simiarlities arithmetic picture = senile ; contrast 'Parallel' senile 1 -1; manova h=senile m=(1 -1 0 0, 0 1 -1 0, 0 0 1 -1) /printe printh; lsmeans senile; title 'Question1: Are profiles parallel?'; run; proc glm data=iq;; class senile; model information simiarlities arithmetic picture = senile ; manova h=senile m=(1 1 1 1) / printe printh; title 'Question2: Are profiles coincident? MANOVA (Hotelling T^2)'; run; data showit; set iq; sum = information +simiarlities+ arithmetic +picture ; run; proc glm data=showit; class senile; model sum = senile; title1 'Question 2: Are profiles coincident?'; title2 'ANOVA on sums'; run; /***** We should NOT test whether level, but I will to demonstrate how to do it *******/ proc glm data=iq;; class senile; model information simiarlities arithmetic picture = ; manova h=intercept m=(1 -1 0 0, 0 1 -1 0, 0 0 1 -1) /printe printh; title 'Question 3: Are profiles Flat?'; run; /******************** IML **************************/ proc iml; /************ Handy module ***************/ start samplestats(X,Xbar,W,S,R,n); n = nrow(X); one = J(n,1); Xbar = X`*one/n; W = (X - one*Xbar`)` * (X - one*Xbar`); S = W/(n-1); Dsqrt = sqrt(diag(S)); R = inv(Dsqrt)*S*inv(Dsqrt); Finish samplestats; /**********************************************/ use iq; read all var{information simiarlities arithmetic picture} into X; read all var{information simiarlities arithmetic picture} where(senile='y') into Xy; read all var{information simiarlities arithmetic picture} where(senile='n') into Xn; run samplestats(Xy,Xbar1,W1,S1,R,n1); run samplestats(Xn,Xbar2,W2,S2,R,n2); print 'Profile Analysis: Senile or not'; C = { 1 -1 0 0, 0 1 -1 0, 0 0 1 -1}; print C; Spool = (W1+W2)/(n1+n2-2); diff = Xbar1-Xbar2; Cdiff = C*diff; p = 4; Tsq = ((n1*n2)/(n1+n2))*Cdiff`*inv(C*Spool*C`)*Cdiff; F = ((n1+n2-p)/((n1+n2-2)*(p-1))) * Tsq; dfnum = p-1; dfden = n1+n2-p; pvalue = 1 - cdf('F',F,dfnum,dfden); print p n1 n2 Spool; print 'Ho(1): parallel profiles C*mu_1 = C*mu_2'; print diff , Cdiff, Tsq dfnum dfden F pvalue; /* Profiles are parrallel, so now test coincidence */ one = J(p,1); tden = one`*(Xbar1-Xbar2); tnum = ((1/n1)+(1/n2))*one`*Spool*one; t = abs(tden/sqrt(tnum)); df = n1+n2-2; pvalue = 2*(1 - cdf('t',t,df)); F = t**2; pvalueF = 1 - cdf('F',F,1,df); print 'Ho(2): profiles are the same or 1`*(mu_1-mu_2)'; print tden tnum df t pvalue, F pvalueF ; /* Profiles are coincident, now test whether flat or level */ p=4; Print 'Ho(3): profiles are flat or C*(mu1 + mu2) = 0'; dfnum = p-1; dfden = n1+n2-p+1; run samplestats(X,Xbar,W,S,R,n); T2 = (n1+n2)*Xbar`*C`*inv(C*S*C`)*C*Xbar; F =( (n - p + 1)/((n-1)*(p-1)) ) * T2; pvalue = 1 - cdf('F',F,dfnum,dfden); Print 'Ho(3): profiles are flat or C*(mu1 + mu2) = 0'; print n Xbar T2 F dfnum dfden pvalue;