[R] Conditional power, predictive power

[francogrex]

Or more elegantly the function below where a and b are the parameters of the
beta prior, xa and xb are the current number of events in group A and B
respectively; na and nb are the current total number of subjects in group A
and B respectively; Na and Nb are the final total number of subject in group
A and B respectively; and alpha is the alpha level (from 0 to 1).
Don't forget to empty your folder: "DATA/predictivepower" (that you create
in advance of course) after each computation.

-----------------------------------------------------------------------------
##Predictive Probability Interim Analysis, John Cook, 2006
predictive.power=function (a,b,xa,xb,na,nb,Na,Nb,alpha){
setwd("C:/R-2.4.1/DATA/predictivepower")
prepower=function(a,b,s0,f0,s1,f1){
factorial(s1+f1)/ (factorial(s1)*factorial(f1))*
beta(s0+s1+a,f0+f1+b)/beta(s0+a,f0+b)}

Nta=Na-na
for(i in 0:Nta) {sink("pw.power1.txt",append=TRUE);
v=prepower(a,b,xa,(na-xa),i,Nta-i)
dput(v)
sink()
}
Ntb=Nb-nb
for(i in 0:Ntb) {sink("pw.power2.txt",append=TRUE);
v=prepower(a,b,xb,(nb-xb),i,Ntb-i)
dput(v)
sink()
}

x=scan(file="pw.power1.txt")
y=scan(file="pw.power2.txt")

zth=function(xa,ya,Na,xb,yb,Nb){
tha=(xa+ya)/Na
thb=(xb+yb)/Nb
th=(xa+xb+ya+yb)/(Na+Nb)
z=(tha-thb)/sqrt(th*(1-th)*((1/Na)+(1/Nb)))
z
}

for(i in 0:Nta){sink("pw.zval.txt",append=TRUE)
zval=(zth(xa,i,Na,xb,0:Ntb,Nb))
cat(zval,sep="\n")
sink()
}
zz=scan(file="pw.zval.txt")

for(i in 1:(Nta+1)){sink("pw.predict.txt",append=TRUE)
pw=((x[i]*y[1:(Ntb+1)]))
cat(pw,sep="\n")
sink()
}
p=scan(file="pw.predict.txt")
z=-qnorm(alpha)
mdf <- cbind(zz, p) 
d1=subset(mdf,apply(mdf>(z),1,any))
d2=subset(mdf,apply(mdf<(-z),1,any))
s1=sum(d1[,2])
s2=sum(d2[,2])
s3=1-(s1+s2)
print(c("p A>B"=s1,"p B>A"=s2, "p A=B"=s3))
}
-- 
View this message in context: http://www.nabble.com/Conditional-power%2C-predictive-power-tf3603396.html#a10094068
Sent from the R help mailing list archive at Nabble.com.