[R] Log rank test power calculations

[Wittner, Ben, Ph.D.]

I suspect that someone will write saying that what you want is in one or a few
packages, but you might also find the code below useful.
(Note, the first two functions require the third.)
-Ben

##
##  Do power calculation of Freedman for Cox PH as set forth on page 733
##  of 5th edition of Rosner, "Fundamentals of Biostatistics"
##
##            k = ratio of (# subjects in exposed group) to
##                         (# subjects in control)
##            t = max time of follow-up
##           RR = (hazard rate for exposed group) /
##                (hazard rate for control group)
##    lambda[j] = Rosner's lambda_(j-1)
##     delta[j] = Rosner's delta_(j-1)
##
survivalUtil.power <- function(n.1, n.2, t, RR, lambda, delta, alpha=0.05) {
  pf <- survivalUtil.probFail(t, RR, lambda, delta)
  k <- n.1/n.2
  m <- n.1*pf$p.E + n.2*pf$p.C
  pnorm(sqrt(k*m)*abs(RR - 1)/(k*RR + 1) - qnorm(1 - alpha/2))
}


##
##  Do sample-size calculation of Freedman for Cox PH as set forth
##  on page 735
##  of 5th edition of Rosner, "Fundamentals of Biostatistics"
##
##            k = ratio of (# subjects in exposed group) to
##                         (# subjects in control)
##            t = max time of follow-up
##           RR = (hazard rate for exposed group) /
##                (hazard rate for control group)
##    lambda[j] = Rosner's lambda_(j-1)
##     delta[j] = Rosner's delta_(j-1)
##
##  Value: n1 - number of people needed in exposed group
##         n2 - number of people needed in control group
##
survivalUtil.sampSize <- function(k, t, RR, lambda, delta,
                                  alpha=0.05, beta=.2) {
  m <- ((1/k)*(((k*RR + 1)/(RR - 1))^2)*
        ((qnorm(1 - (alpha/2)) + qnorm(1 - beta))^2))
  pf <- survivalUtil.probFail(t, RR, lambda, delta)
  n.1 <- (m*k)/((k*pf$p.E) + pf$p.C)
  n.2 <- m/((k*pf$p.E) + pf$p.C)
  list(n.1=n.1, n.2=n.2)
}


##
##  Compute probabilties of failure used by survivalUtil.power() and
##  survivalUtil.sampSize()
##
survivalUtil.probFail <- function(t, RR, lambda, delta) {
  A <- cumprod(1 - lambda)[1:t]
  B <- cumprod(1 - (RR*lambda))[1:t]
  C <- cumprod(1 - delta)[1:t]
  D <- lambda[2:(t + 1)]*A*C
  E <- RR*lambda[2:(t + 1)]*B*C
  list(p.C=sum(D), p.E=sum(E))
}

> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
> On Behalf Of Daniel Brewer
> Sent: Thursday, January 31, 2008 10:53 AM
> To: r-help at stat.math.ethz.ch
> Subject: [R] Log rank test power calculations
> 
> Does anyone have any ideas how I could do a power calculation for a log
> rank test.  I would like to know what the suggested sample sizes would
> be to pick a difference when the control to active are in a ratio of 80%
> to 20%.
> 
> Thanks
> 
> Dan
> 
> --
> **************************************************************
> Daniel Brewer, Ph.D.
> Institute of Cancer Research
> Email: daniel.brewer at icr.ac.uk
> **************************************************************
> 
> 
> The Institute of Cancer Research: Royal Cancer Hospital, a charitable
> Company Limited by Guarantee, Registered in England under Company No.
> 534147 with its Registered Office at 123 Old Brompton Road, London SW7
> 3RP.
> 
> This e-mail message is confidential and for use by the a...{{dropped:2}}
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

The information transmitted in this electronic communication is intended only
for the person or entity to whom it is addressed and may contain confidential
and/or privileged material. Any review, retransmission, dissemination or other
use of or taking of any action in reliance upon this information by persons or
entities other than the intended recipient is prohibited. If you received this
information in error, please contact the Compliance HelpLine at 800-856-1983 and
properly dispose of this information.