# This file is posted as power.ssc n _ 2:50 power _ 1 - pf(qf(.95,2,3*n-3),2,3*n-3,n/2) # Plot the power curve par(mfrow=c(1,1),fin=c(5,5)) plot(c(0,50), c(0,1), type="n", xlab="Sample Size",ylab="Power", main="Power versus Sample Size\n F-test for equal means") lines(n, power, type="l",lty=1) # Create a loop to serch for the # smallest sample size to achieve # the specified power level. Here # we create the model matrix for # the one way anova needed to # evaluate the noncentrality # parameter nmax=50 alpha= .05 power= .90 # Enter the C matrix for the null # hypothesis of CB=0 C <- matrix(c(0, 1, 0, -1, 0, 0, 1, -1),2,4,byrow=T) # Enter deviations from the null hypothesis as # multiples of the stadard deviation for the # random errors. d <- matrix(c(0.5, 1), 2,1,byrow=T) # Enter the portion of the model matrix # that will be replicated as the sample # size is increased X <- cbind(matrix(1,3,1),diag(rep(1,3))) # Set the initial sample size n<-2 X<- rbind(X,X) # Create a function for computing power powerfun<-function(n,X,C,d,alpha){ nc<-t(d)%*%solve(C%*%ginverse(t(X)%*%X)%*%t(C))%*%d df1 <- nrow(C) df2 <- nrow(X)-qr(X)\$rank 1-pf(qf(.95,df1,df2),df1,df2,nc) } # Enter the loop to evaluate the # required sample size npower<-0 while (nnmax you have # not found the desired sample size and you # must increase nmax