bin.pwr <- function(n,p0=1/2,RR=NULL,alpha=.05) { # plot power vs p for H0: p=p0, default p0=1/2 # RR if present is the chosen rejection region # alpha is the target significance level if(is.null(RR)){ X = 0:n db = dbinom(X,n,p0) sdb = sort(db) cb = cumsum(sdb) m = max((1:(n+1))[cb <= alpha]) R = X[db <= sdb[m]] } else{ R=RR } Alpha <- sum(dbinom(R,n,p0)) p <- seq(.01,.99,.01) P <- p for(i in 1:length(p)) { P[i] <- sum(dbinom(R,n,p[i])) } plot(p,P,ylab="Power",type="l") title(paste("Power curve for n=",n," alpha =",round(Alpha,4))) }