##### R-code for making Figure 3.8.

postscript("fig38.ps",width=7.2,height=7.2,horizontal=F)

par(mfrow=c(2,2), mar=c(4,4,1,2))

alpha = 0.01
tpi = seq(0.01,0.99,0.01)

m = 25
talpha25 = numeric()

for(i in 1:length(tpi)){

U=m*tpi[i]+qnorm(0.995)*sqrt(m*tpi[i]*(1-tpi[i]))
L=m*tpi[i]-qnorm(0.995)*sqrt(m*tpi[i]*(1-tpi[i]))-0.0000001

talpha25[i] = pbinom(L,m,tpi[i])+(1-pbinom(U,m,tpi[i]))

}

plot(tpi,talpha25,type="l",lty=1,xlab=expression(pi),
     ylab=expression(tilde(alpha)),mgp=c(2,1,0),xlim=c(0,1), 
     ylim=c(0,0.04), cex=0.8)
lines(tpi,rep(alpha,length(tpi)),lty=2,cex=0.8)
title(xlab="(a)",cex=0.9)

m = 50
talpha50 = numeric()

for(i in 1:length(tpi)){

U=m*tpi[i]+qnorm(0.995)*sqrt(m*tpi[i]*(1-tpi[i]))
L=m*tpi[i]-qnorm(0.995)*sqrt(m*tpi[i]*(1-tpi[i]))-0.0000001

talpha50[i] = pbinom(L,m,tpi[i])+(1-pbinom(U,m,tpi[i]))

}

plot(tpi,talpha50,type="l",lty=1,xlab=expression(pi),
     ylab=expression(tilde(alpha)),mgp=c(2,1,0),xlim=c(0,1), 
     ylim=c(0,0.04), cex=0.8)
lines(tpi,rep(alpha,length(tpi)),lty=2,cex=0.8)
title(xlab="(b)",cex=0.9)


alpha = 0.001
tpi = seq(0.01,0.99,0.01)

m = 25
talpha25 = numeric()

for(i in 1:length(tpi)){

U=m*tpi[i]+qnorm(0.9995)*sqrt(m*tpi[i]*(1-tpi[i]))
L=m*tpi[i]-qnorm(0.9995)*sqrt(m*tpi[i]*(1-tpi[i]))-0.0000001

talpha25[i] = pbinom(L,m,tpi[i])+(1-pbinom(U,m,tpi[i]))

}

plot(tpi,talpha25,type="l",lty=1,xlab=expression(pi),
     ylab=expression(tilde(alpha)),mgp=c(2,1,0),xlim=c(0,1), 
     ylim=c(0,0.04), cex=0.8)
lines(tpi,rep(alpha,length(tpi)),lty=2,cex=0.8)
title(xlab="(c)",cex=0.9)

m = 50
talpha50 = numeric()

for(i in 1:length(tpi)){

U=m*tpi[i]+qnorm(0.9995)*sqrt(m*tpi[i]*(1-tpi[i]))
L=m*tpi[i]-qnorm(0.9995)*sqrt(m*tpi[i]*(1-tpi[i]))-0.0000001

talpha50[i] = pbinom(L,m,tpi[i])+(1-pbinom(U,m,tpi[i]))

}

plot(tpi,talpha50,type="l",lty=1,xlab=expression(pi),
     ylab=expression(tilde(alpha)),mgp=c(2,1,0),xlim=c(0,1), 
     ylim=c(0,0.04), cex=0.8)
lines(tpi,rep(alpha,length(tpi)),lty=2,cex=0.8)
title(xlab="(d)",cex=0.9)

graphics.off()