# PlotPoweR.txt:  Plotting power functions with R
corr = 0.1; n = seq(from=0,to=3000,by=10)
fmaxinf = -log(1-corr^2); ncp = n*fmaxinf; ncp[1:20]
critval = qchisq(0.95,5)
Power = 1 - pchisq(critval,5,ncp) 
plot(n,Power,type="l",ylim=c(0,1))
titlestring = expression(paste(  "Power of the Chi-squared test for fit with ",
              rho, " = Corr(", e[1], ",", e[2], ")"  ))
title(titlestring)
# Add lines for other values of corr
corr = 0.15
fmaxinf = -log(1-corr^2); ncp = n*fmaxinf
Power = 1 - pchisq(critval,5,ncp) 
lines(n,Power,lty=2)
corr = 0.20
fmaxinf = -log(1-corr^2); ncp = n*fmaxinf
Power = 1 - pchisq(critval,5,ncp) 
lines(n,Power,lty=3)
# Make line labels
x1 = c(1500,2250) ; y1 = c(0.40,0.40) ; lines(x1,y1,lty=1)
x2 = c(1500,2250) ; y2 = c(0.30,0.30) ; lines(x2,y2,lty=2)
x3 = c(1500,2250) ; y3 = c(0.20,0.20) ; lines(x3,y3,lty=3)
rhoeq = expression(paste(rho,' = 0.10')); text(2500,0.4,rhoeq)
rhoeq = expression(paste(rho,' = 0.15')); text(2500,0.3,rhoeq)
rhoeq = expression(paste(rho,' = 0.20')); text(2500,0.2,rhoeq)