J CIbaEpUambl Vol. 45, No. IO, pp. 1131-l 136, 1992 Printed in Great Britain. All rights rwcrval

0895-4356/92 S5.00 + 0.00 Copyright 0 1992 Pcrgamon Press Ltd

SAMPLE SIZE CALCULATIONS FOR THE LOG RANK TEST: A GOMPERTZ MODEL APPROACH ALAN B. CANTOR* Departments of Pediatrics and Statistics, University of Florida, Gainesville, FL 32605, U.S.A. (Received in revised form 4 March 1992)

Abstract-Sample size calculations for clinical trials dealing with survivorship are often based on an exponential model. This model is inappropriate when a non-zero proportion of the population is expected to have indefinite survival. In such cases the Gompertx model offers a reasonable alternative. A method for calculating the required accrual time for a clinical trial in which the treatment arms have Gompertx survival distributions satisfying the proportion hazards assumption is developed. A computer program to perform this method is given, as well as an iterative method that can be used when a computer is not available. Survival analysis

Gompertx distribution

I. INTRODUCI’ION

Sample size

the hazards, error rates, accrual time, accrual rate, follow-up time, and loss rate, which can be In recent years, increased attention has been solved for any parameter in terms of the others. paid to the problem of determining the approSolution for accrual time requires an iterative priate sample size or accrual time needed for technique such as Newton’s method, by which studies in which the variable of interest is the convergence is rapid and insensitive to the initial length of time until the occurrence of a well estimate. defined endpoint, often death. The analysis of While the results of the log rank test will the data from such studies needs to take into be valid regardless of the underlying survival account the fact that, for some patients, the distributions, incorrect specification of those endpoint will not be observed. (Such obserdistributions leads to incorrect sample size calvations are said to be right censored.) Pasterculations. With the emphasis, in recent years, on nack and Gilbert [l] and Pastemack [2] present adequate sample size for clinical studies, it is results based on performing the likelihood ratio important to have methods of estimating the test for equality of hazards for two groups required sample size when the assumption of having exponential survival distributions. With exponentiality is clearly inappropriate. One the introduction of the log rank test [3] and the such non-exponential situation occurs when surdemonstration by Schoenfeld [4] of its asympvival curves exhibit non-zero horizontal asymptotic efficiency under equal censoring distributions, that non-parametric test became, for totes, thus indicating a proportion of indefinite survivors or “cures”. As shown by Cantor [6], many statisticians, the method of choice. attempts to estimate the required sample size by Rubenstein et al. [5] give an equation relating fitting exponential distributions can lead to gross underestimates in such cases. *All correspondence should be addressed to: Alan B. The literature on sample size calculation for Cantor, Pediatric Oncology Group Statistical Office, non-exponential models is relatively sparse. 4110 S.W. 34th Street. Suite 22. Gainesville, FL 3260% 6537, U.S.A. [Tel: (904) 392-5198; Fax: (904) 392-81621. Sposto and Sather [7j present general results 1131

ALANB.

1132

CANTOR

that assume proportional hazards and apply them to the model in which one group has a survival function of the form s,(t) = n + (1 - n)exp( -At) and the other survival function is S*(t) = [S,(t)]“. Note that while S,(t) expresses the fact that a proportion, n, is cured while the remainder had exponential survivorship, &(t) is not of this form. Shuster [8] takes a different approach. He considers a model for S, (t) with average hazard rl, for t < t,, and & for t 2 to. Letting 1r = 0 results in a non-zero proportion of “cures”. Again Sr(t) = [S, (t)lA. Extensive tables are provided. Finally, Halpern and Brown [9] present a computer program, which they make available, that estimates power using simulations. They allow models of the form S,(t)=JT,+(l -n,)exp(-l,t) and Their pro&(t) = nr+ (1 - n,)exp(-I,t). gram also permits the specification of sets of points to define arbitrary survival functions. II.

THE

GOMPERTZ

MODEL

survival function is s(t) = The Gompertz exp{(-O/y)[exp(yt) - 1]}, where 0 > 0. It can be derived from the assumption that the log hazard is linear. It was first described by Gompertz [lo], who observed that, in a life-table of aging males, the hazard seemed to increase exponentially. For positive values of y, the slope of the log hazard, the hazard increases and lim S(t) = S(co) = 0. If y < 0, then S(co) = exp(e/y), which will henceforth be denoted n. While S(t) is discontinuous at y = 0, that discontinuity is removable by defining s(t) = exp( -et) in that case. As a survival function describing a population with a non-zero proportion of cures, the Gompertz function has several attractive features: If S,(t) is Gompertz with parameters 8 and y, then [S, (t)lA is Gompertz with parameters BA and y. Thus the cure rate for [S, (t)lA is nA if the cure rate for S,(t) is n. With the discontinuity at y = 0 removed as described above, the exponential distribution becomes a special case in the interior of the parameter space. The two equations setting the partial derivatives of the log likelihood equal to zero can be reduced to one. That equation is easily solved for MLEs by Newton’s method. others [1 1] have noted By contrast, the difficulties in finding MLEs for S(t) = J7 + (1 - n)exp(-it).

0

2

I

Years LEGEND

6

10

6

12

Followed

-Kaplan-Meier -GomDertz

Fit

Fig. 1. Kaplan-Meier and Gompertz survival curves for one treatment arm in POG 7834, a Phase III study of acute lymphocytic leukemia.

This author has fit the Gompertz survival function to a large number of treatment groups studied by the Pediatric Oncology group and its predecessor, the Southwest Oncology Group. It was found that, for those treatment groups with apparent plateaus, good fits were always obtained and convergence to MLEs of the parameters posed no problem. As an example, Fig. 1 presents the Kaplan-Meier and Gompertz survival curves for one of the treatment groups in a study of acute lymphocytic leukemia. The Gompertz curve is based on the maximum likelihood estimates of the parameters. III. THE POWER

EQUATION

Suppose that two treatment groups indexed by 1 and 2 have survival functions S,(t) = [email protected](Vy)[exp(yQ

-

111

and S*(t) = [S,(t)]’

8, A > 0.

(1)

We wish to test H,:A = 1 against either a onesided or two-sided alternative at significance level a with power 1 - j3 if A = &. Accrual will be at rate r with a patient’s starting time being a uniform random variable over the interval [0, 7’j. Furthermore, follow-up will continue for an additional time r after accrual ends.

Sample Size Calculation

It is shown by Rubenstein ef al. [5j that the likelihood ratio test that the ratio of two exponential hazards is unity has error rates a and j? if

1133

by m the median survival time for those not cured, i.e. S(m) = (1 + n)/2.

(5)

This parameter is related to the others by m=iln where D, is the number of deaths in group i, and A,, is the alternative hazard ratio and @(Z,) = 1 -a, @(.) being the standard normal distribution function. It follows from the asymptotic efficiency of the log rank test that equation (2) must hold (asymptotically) for that test as well. Our approach will be to replace each Di by its expected value based on the above assumptions. Then the resultant equation can be solved for any parameter of interest given the others. The solution for T, accrual time, is of particular interest. Leaving the details of the derivation to Appendix A, the equation can be written rT

(In A)’

T (Z, + Z,)’

=i~211-ni+[F(xC!) (

- I;(+, )I~i/(YTIl)-‘, (3) where 41

= -In n, exp(yt)

Xi2= xii exp(V) Q = expV/r) nz=n:’ for i=l,2 and x’/(j *j!).

F(x) = i

(4)

j-l

The unique solution of equation (3) for T by Newton’s method is readily accomplished. A BASIC program is given in Appendix B. IV. AN ALTERNATIVE The

PAIR OF PARAMETERS

parameters y and 8, while describing the hazard associated with a Gompertz survival time variable, have no natural interpretation relative to the survival function itself. In this section an alternate pair of parameters, more readily interpreted, is introduced. If y < 0 then the limit of the survival function as t approaches infinity is exp(B/y). As indicated above, this value will be denoted R and is the proportion of indefinite survivors or “cures”. We will denote

l{

(6)

The parameters n and m describe the Gompertz distribution more naturally and are used in the program in Appendix B. V. AN EXAMPLE

Consider a clinical trial designed to test whether an experimental treatment provides superior survivorship to a control treatment. Suppose that the control treatment, based on earlier studies, is thought to have a cure rate of 30% and that half of the mortality occurs within 2 years. We wish to design a study with a one-sided significance level of 0.05 and a power of 80% if the experimental treatment has a cure rate of 50%. Patients will be followed for 2 years after the end of accrual, and an accrual rate of 40 patients/year is anticipated. Taking Z, = 1.645, Z, = 0.84 and A = ln(0.3)/ln(0.5), equation (3) can be solved for T obtaining an accrual time of 4.93 years or 198 patients. It is interesting to consider the effect of the follow-up time on the accrual needed. By eliminating the follow-up period completely (i.e. r = 0), the accrual time becomes 6.31 years. Thus by accruing 55 more patients over 17 months, the analysis can be done 5 months earlier. Reducing the follow-up period to 1 year requires only 0.56 year more or 22 more patients. By contrast, increasing the follow-up can reduce the necessary accrual. With 3 or 4 years of follow-up, 182 or 173 patients respectively are needed. As r increases without bound, the required sample size approaches 3.47 years for 139 patients--approximately the same as the number required for the comparison of the two cure rates using the normal approximation to the binomial. Note that larger values of m are associated with longer accrual times. Also, smaller values of m make follow-up time less worthwhile. In the above example, if m were 0.5 instead of 2, elimination of the 2 year follow-up period would require only 0.76 more years of accrual. of’ course it must be recognized that elimination or reduction of the foilow-up period could compromise our ability to determine the reasonableness of the Gompertz model.

ALANB. CANTOR

1134

Table 1. Comoarison of exoonential and Gomcertz models

Gompertz q = 0.25

Exponential

s, (2) 0.1 0.2

0.3 0.4

0.5

0.6

W)

Ne

Ng

Power

0.2 0.3 0.4 0.3 0.4 0.5 0.4 0.5 0.6 0.5 0.6 0.7 0.6 0.7 0.8 0.7 0.8 0.9

253 85 48 391 120 63 455 137 70 460 140 72 424 133 68 367 119 62

265 89 49 422 128 66 503 147 14 518 151 75 482 143 71 413 126 65

0.78 0.78 0.79 0.17 0.77 0.78 0.76 0.77 0.78 0.75 0.76 0.78 0.74 0.77 0.78 0.74 0.77 0.78

q = 0.50

q = 0.15

Ng

Power

Ng

Power

213 91 50 443 132 61 540 154 76 565 159 78 531 150 73 454 132 66

0.11 0.17 0.79 0.75 0.76 0.77 0.73 0.75 0.77 0.71 0.74 0.77 0.70 0.74 0.77 0.70 0.75 0.77

280 93 51 468 137 69 588 162 79 635 170 81 613 162 71 532 142 69

0.76 0.11 0.78 0.73 0.75 0.76 0.69 0.73 0.75 0.66 0.72 0.75 0.64 0.71 0.75 0.63 0.72 0.75

S,(2) and S,(2) are 2 year survival rates, q = I7,/S,(2) where II, is treatment 1 cure rate, a-sided a = 0.05, power = 0.8, follow-up = 2 years.

While it is natural to express the sample size goals in terms of the number of patients registered, it should be recognized that the actual requirement is based on the number of deaths which is a function of both the number of patients accrued and their time on study. Thus if accrual is slower than expected, fewer patients will be needed since they will be on study longer. Conversely, if accrual is faster, more patients will be needed. It may be worthwhile to monitor accrual and recalculate the accrual time and sample size needed if the rate differs from what was projected. VI. COMPARISON

TO EXPONENTIAL

MODEL

We next consider the effect of erroneously calculating the required sample size using an exponential model when a Gompertz model allowing for cures would be more appropriate. Suppose we intend to use a log rank test with 2-sided significance level 0.05 and a power of 0.80 is desired for a preassigned pair of 2 year survival rates S, (2) and S,(2). The expected accrual rate is 50 patients per year and we plan to follow patients for 2 years after accrual ends. The left side of Table 1 gives the required sample sizes, obtained by the method of Rubenstein [5] assuming exponentiality. The right side gives the required samples sizes under the Gompertz assumption with q representing the proportion of 2 year survivors who are assumed to be long term survivors. It can be seen that if most of the 2 year survivors eventually die (q = 0.25) then the error caused by the exponential assumption is not great. The sample

sizes are M-98% of what they should be and the actual power with the smaller sample usually exceeds 0.75. However, if most of the 2 year survivors are cured (q = 0.75), then samples sizes are often 70-85% of what they should be and the power, using the smaller sample size is generally between 0.70 and 0.75. VII. A COMPUTER-FREE

METHOD

Since the required accrual time is a function of six other parameters, a single table or graph that could be used to determine T is not feasible. Instead we can express T as T= f

(Z, + Z,r)2/(ln A)2,

(7)

where Pi is the proportion in arm i expected to die. A given value of Tallows calculation of the right side which in turn, leads to a new value of T. This iterative procedure can be performed as follows: (1) Choose a value for T. (2) Let Pi = 1 - ZIi + 17iaim/(TCi) for i = I,2 where Ci = ln[ 1 - ln(0.5 + ni/2)/(mni)], and ~i = G(Ui, r/m) - G(.lIi, (r + T)/m). The values of G can be read from Fig. 2. (3) Let T = :

(Z, + Z,)2/(ln A)2.

(4) Continue (if desired) by returning to step 2 with the new value of T.

Sample Size Calculation

0

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

Fig. 2. Graphs to facilitate accrual time calculation. Highest curve is for lZ = 0.1, next is for l7 = 0.2, etc.

Note that at any step the power can be obtained from

When this method is applied to the example of Section V with a starting value of T = 3.0 years, the first iteration produces a value of 4.64 years and a power of 0.77. The second iteration produces a value of 4.75 years and a power of 0.78. Although this method will, in theory, converge to the solution of equation (2) just as the program of Appendix B does, in practice its accuracy will be limited to the accuracy with which values of G can be read from Fig. 2.

For clinical trials with survival time as an endpoint, a reasonable survival model is needed for sample size calculations. When a non-zero proportion is expected to be cured, the Gompertz model presented here may be a good choice. Use of the program in Appendix B or the method of Section VII allow the accrual time needed to be easily determined based on this model. Acknowledgement-This research was supported by Grant CA29139 from the National Cancer Institute. REFERENCES

1 Pastemack ES, Gilbert HS. Planning the duration of ’ long-term survival time studies designed for accrual by cohorts. J Cluon Dis 1971: 24: 681-700. BS. Sample sixe for clinical trials designed 2. Pastemack for patient accrual by cohorts. J Chron L&I 1972; 25: 673-681. 3. Mantel N. Evaluation of survival data and two new

VHL CONCLUSIONS

Failure to give careful consideration to sample size when planning a clinical trial can have several unfavorable consequences. Important differences may be missed because of an inadequate sample size. Negative results might be naively interpreted as showing that there is no difference, while more sophisticated readers would realize that such results are uninformative. On the other hand, an excessive sample can be wasteful of clinical and financial resources.

rank order statistics arising in its consideration. Cancer 4%metl~ Rep 1966, SO: 163-170. 4. Schoenfeld D. The asymptotic properties of nonparametric tests for comparing survival distributions. B&m&&a 1968; 68: 316-319. 5. Rubenstein RV, Gail MH, Santner TJ. Planning the duration of a comparative clinical trial with loss to follow-up and a period of continued observation. J Cbntn Dis 1981; 34: 469-479. 6. Cantor AB. Power estimation for rank tests using censored data: conditional and unconditional. Contr 7 CBn Trials 1991; 12: 462-473. Sposto R, Sather HN. Determining the duration of ’ comparative clinical trials while allowing for cure: J Chmn Dis 1985; 38: 683-690.

ALAN B. CANTOR

1136 8.

Shuster JJ. Handbook of Sample Size Guidelines for Cutdeal Trials. Boca Raton, FL: CRC Press; 1990 9. Halpem J, Brown BW. Designing clinical trials with arbitrary specification of survival functions and for the log rank or generalized Wilcoxon test. Contr Clin Trials 1987; 8: 177-189. 10. Gompertz B. On the nature of the function expressive

of the law of human mortality, and on the new mode of determining the value of life contingencies. Phil Tram R Sot A 1825; 115: 513-580. 11. Goldman AI. Survivorship analysis when cure is a possibility: A Monte Carlo study. Stat Med 1984; 3: 153-163.

APPENDIX

A

In general, the expected number of deaths in each group is rT/2 times the probability of death for that group. A patient arriving at time I in (0, T) will die on study if his survival time is less than T + 7 - t. Thus if the survival function is S(t), and the arrival times are uniformly distributed in (0, T), the probability of death is 1 = T [1-S(T+s-t)]dt. (9) I0 For the Gompertz distribution with parameters y and 0, this becomes k

T{l-exp[(-O/y)exp(y(T+r

-t)-l)]}dr.

(10)

I The substitution u = exp[y(T + r -

r)]

;ields eXP(,r)

[email protected]/y)/u du. 111 s ex~lr(r+ The integral can be accomplished by writing the exponential function as a power series. The result is:

1+ exp(~/yYW)

1 -exp(e/y)

+ exp(6/y)/tiT)VI-68/y

exW)l-

(11)

FL-e/y ew(yV + 7))lI

where F(x) = f x’/(i .i!). (12) ,=I The above expression gives the value of P, Replacing tJ by 6A gives the value of P2 and multiplication of each by rT/2 gives E(Q) and I?(&). Result (3) follows. APPENDIX

B

Program to Calculate Accrual Time Based on Gompertz Model

DECLARE FUNCTION series# (x#) REM This program is called gompow.bas. Its purpose is to find a value of REM t, accrual time, that produces a given power for the logrank test when the REM gompertx model is assumed and the other parameters are specified. The REM program uses a Newton method iteration to solve the equation relating REM the narameters to newer fort. REM It &ems that convergence can sometimes fail if the initial value is too REM large, but it never seems to fail because the initial value is too small. REM Aninitial value of 1 has always worked. REM REM input parameters INPUT “cure rates in group 0 and group 1”; pio#, pil# INPUT “median survival time for non-cures in group 0 “; m# INPUT “followup time and accrual rate “; tau#, r# INPUT “zalpha and zbeta “; zalpha#, zbeta# INPUT “initial value of accrual time “; t# betag = LOG(l# - LOG((l# + piO#) /2#) ! LOG(pio#)) I m# t%Yt=s$es(-EXP(beta# * tau#) * LOG(piO#)): i2# = senes(-EXP(beta# * tau#) * Lffi(pil#)) 10 & = seties(-EXP(heta# * (t# + tau#)) * LOG(piO#)) f4# = series(-EXP(beta# * (t# + tau#)) * LOG(pi I#)) nint = nint + 1 a#=t#-t#*piO#-piO#/beta#*(f3#-fl#):b#=t#-t#*pil#-pil#/betll#*(f4#-fl#) c# = r# I 2 * (LOG([email protected]#) I LOG(pio#))) A2 I (zalpha# + zbeta#) A 2 d# = pio# h(l-EXP(be.ta#*(t# + tau#)))-1 e# = il# h(l-EXP(beta#*(t# + tau#)))-1 newt R =t#-(l/a#+l/b#-c#)/(d#/a#A2+e#/b#A2) IF ABS(t# - [email protected] c .Ol THEN PRINT ‘It= “; INT(CSNG(newt#) * 1000 + .05) I IO(K);“after “; nint; ” iterations”: END ENDIF IF nint = 30 THEN PRINT “failed to converge”: END IF ABS(t# - newt#) > .ftOOOOlTHEN t# = newt* GOT0 10

END‘



FUNCTION series# (x#) s#=O

se#=l n#=O fact# = 1 DO WHILE (EXP(x#) - se#) I (n# + 1) > .OOOO I n#=n#+ 1 fact# = fact# * n# s#=s#+x#An#ln#lfact# s&se# + x# h n# I fact# .series# = s# END FUNCTION

Sample size calculations for the log rank test: a Gompertz model approach.

Sample size calculations for clinical trials dealing with survivorship are often based on an exponential model. This model is inappropriate when a non...
585KB Sizes 0 Downloads 0 Views