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