This program uses simulation to approximate the power of a test of a marker association with survival data assuming an exponential regression model. The calculations include a proportional hazards test using the model where the log-hazard ratio is linear in a continuous marker distribution. The choice of distributions include uniform, normal, and lognormal. All covariates are mean centered for the analysis. The magnitude of the effect size can be specified based on a unit of the prognostic factor or the interquartile range of the marker.

It is useful to note that if the underlying distribution is assumed to be uniform, and the proportional hazards model regression log-hazard ratio is linear in the marker, then effect size associated with the interquartile range will correspond to approximately same effect as a marker dichotomized at the median (or other quantile).

In addition to power calculations, the approximate confidence intervals of parameter effects and expected number of events are also presented.

The impact of the continuous prognostic model assumptions in terms of the expected events can be helpful in gaging assumptions on expected events relative to the overall trial design.

Input Items

The user is prompted to enter values for the following items. For items that have initial default values set, the default values are given in parentheses.

Output Items


This is a simulation-based calculation and therefore may take several seconds to calculate depending both on the number of simulations and sample size.

Statistical Code

The program is written in R.

View Code

function(accrual, follow, n, alpha, sided, lambda1, unit, HR1, dist, dist_param, nsim) {
    # Power for continous marker #
    # Note that covariates are mean centered in calculations #


    beta1 = log(HR1)

    # This is where the biomarker distibution is chosen
    if (dist == "unif") { dist.generate = runif }
    if (dist == "norm") { dist.generate = rnorm }
    if (dist == "lnorm") { dist.generate = rlnorm }

    testscorea = coefa = nevent = rep(NA, nsim)

    # Loop over simulations 
    for (j in 1:nsim) {

        # Generate Data #

        biomarker = dist.generate(n, 0, dist_param)

        biomarker0 = biomarker - mean(biomarker)
        if (unit == "IQR") {
            biomarker0 = biomarker0 / IQR(biomarker0)

        tt = (-log(runif(n))) / exp(log(lambda1) - beta1 * biomarker0) # generate true survival by taking exponential(1) / exp(alpha + beta * biomarker)
        cc = follow + runif(n) * accrual                               # generates true censoring times
        status = 1 * (tt < cc)                                         # generates 0,1 censoring indicator
        times = tt * status + cc * (1 - status)                        # observed failure times (times under observation)

        # Do Analyses #

        rega = coxph(Surv(times, status) ~ biomarker0)                                 # do cox regression continuous
        coefa[j] = rega$coef                                                           # cox reg coefficient (for continuous variable)
        testscorea[j] = rega$score                                                     # score test
        nevent[j] = sum(status)

    powercont = mean(testscorea > qchisq(1 - alpha * sided, df=1))
    coefcont = quantile(coefa, c(alpha / sided, .5, 1 - alpha / sided))
    events = mean(nevent)

    result = list(power_continuous = round(powercont, digits = 3),
                  expected_events = (floor(events * 10) / 10),
                  effect = round(exp(-coefcont[2]), digits = 3),
                  lower_ci = round(exp(-coefcont[3]), digits = 3),
                  upper_ci = round(exp(-coefcont[1]), digits = 3))
    return(jsonlite::toJSON(result, pretty = TRUE))