/* (Programmed by Phil Chapman; revised 9-7-95; send comments or

   corrections to: pchapman@lamar.colostate.edu)

   This program calculates power for the F-tests in a two-factor

   experiment in a completely randomized design.  The user inputs

   the matrix of true treatment means into the matrix mu.  The

   program calculates power for the overall F-test for no trt effects,

   as well as the F-tests for main effects and interactions.

   The program can be modified if the design is blocked by changing

   the degrees of freedom for error (dfe) to (r*c-1)*(n-1).

*/

proc iml;

*Begin user inputs;

n1=2;             *Smallest number of observations per trt considered;

n2=20;            *Largest number of observations per trt considered;

n3=1;             *Increment in number of observations per trt;

alpha=0.05;       *Type I error rate;

sigma=3;           *sigma is the error standard deviation;

                   *The true means are input as the matrix mu.

                    In a 2 by 4 factorial mu will have 2 rows and

                    4 columns.  Note that a comma is used to

                    separate rows;

mu={1 2 3 4,

   5 6 7 8};

*End user inputs;

create power var{r c row col rcint trt n1 n2 n3 alpha sigma};

r=nrow(mu);

c=ncol(mu);

ra=j(c,1,1);

ca=j(1,r,1);

rtot=mu*ra;

ctot=ca*mu;

cor=sum(mu)**2/(r*c);

row=ssq(rtot)/c-cor;

col=ssq(ctot)/r-cor;

rcint=ssq(mu)-row-col-cor;

trt=ssq(mu)-cor;

print trt row col rcint;

append;

quit;

data power2;

  set power;

sig2=sigma*sigma;

dfrow=r-1;

dfcol=c-1;

dfint=(r-1)*(c-1);

dftrt=r*c-1;

do n=n1 to n2 by n3;

    dfe=r*c*(n-1);

critrow=finv(1-alpha,dfrow,dfe);

critcol=finv(1-alpha,dfcol,dfe);

critint=finv(1-alpha,dfint,dfe);

crittrt=finv(1-alpha,dftrt,dfe);

powrow=1-probf(critrow,dfrow,dfe,n*row/sig2);

powcol=1-probf(critcol,dfcol,dfe,n*col/sig2);

powint=1-probf(critint,dfint,dfe,n*rcint/sig2);

powtrt=1-probf(crittrt,dftrt,dfe,n*trt/sig2);

output;

end;

proc print;

  var r c n powrow powcol powint powtrt;

run;