[R] calculating power of log rank test
raymond chiruka
rtchiruka at yahoo.com
Fri Oct 19 11:31:48 CEST 2007
hie
Im trying to calculate the power of the logrank test for different values of rho .I was just wandering whether the following programme would do it. any suggestions are welcome
s=50
number=1
count1=0;count2=0;count3=0;count4=0;count5=0;count6=0;count7=0;count7=0;
count8=0;count9=0
while(s!=0){
n=20
n1=n/2
n2=n/4
a11=1 ;a12=1.4 ;a21=16 ;a22=a12 * a21
t1<-array(1,c(n1))
t2<-array(2,c(n1))
treatgrp=matrix(c(t1,t2))
st11<-array(1,c(n2))
st12<-array(2,c(n2))
st21<-array(1,c(n2))
st22<-array(2,c(n2))
strata=matrix(c(st11,st12,st21,st22))
time11=array(rweibull(n2,a11,1))
time12=array(rweibull(n2,a12,1))
time21=array(rweibull(n2,a21,1))
time22=array(rweibull(n2,a22,1))
time=matrix(c(time11, time12, time21, time22))
censorTime=runif(n,0,1)
m=cbind(treatgrp,strata,time,censorTime)
colnames(m)=c(treat,strata,time,censorTime)
act.surv.time<-pmin(m[,"time"],m[,"censorTime"])
m<-cbind(m,act.surv.time)
m<-cbind(m,0)
m[m[,3]>m[,4],6]<-1
colnames(m)[6]<-"censoring"
c1= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=-2)
c2= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=-.15)
c3= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=-1)
c4= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=-0.5)
c5= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=0)
c6= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=0.5)
c7= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=1)
c8= survdiff(Surv(act.surv.time,censoring)~treatgrp ,data=b,rho=1.5)
c9= survdiff(Surv(act.surv.time,censoring)~treatgrp,data=b,rho=2)
d1=(sqrt(c1$chisq))
d2=(sqrt(c2$chisq))
d3=(sqrt(c3$chisq))
d4=(sqrt(c4$chisq))
d5=(sqrt(c5$chisq))
d6=(sqrt(c6$chisq))
d7=(sqrt(c7$chisq))
d8=(sqrt(c8$chisq))
d9=(sqrt(c9$chisq))
if(d1>1.96)count1=count1+1
if(d2>1.96)count2=count2+1
if(d3>1.96)count3=count3+1
if(d4>1.96)count4=count4+1
if(d5>1.96)count5=count5+1
if(d6>1.96)count6=count6+1
if(d7>1.96)count7=count7+1
if(d8>1.96)count8=count8+1
if(d9>1.96)count9=count9+1
number=number+1
s=s-1
}
power1=count1/number
power2=count2/number
power3=count3/number
power4=count4/number
power5=count5/number
power6=count6/number
power7=count7/number
power8=count8/number
power9=count9/number
print(power1)
print(power2)
print(power3)
print(power4)
print(power5)
print(power6)
print(power7)
print(power8)
print(power9)
I know the programming is not the best but I'm still new at this.
__________________________________________________
Do You Yahoo!?
Tired of spam? Yahoo! Mail has the best spam protection around
http://mail.yahoo.com
More information about the R-help
mailing list