Email updates

Keep up to date with the latest news and content from BMC Medical Research Methodology and BioMed Central.

Open Access Research article

On reporting results from randomized controlled trials with recurrent events

Lisa Kuramoto1*, Boris G Sobolev2 and Meghan G Donaldson3

Author Affiliations

1 Centre for Clinical Epidemiology and Evaluation, Vancouver Coastal Health Research Institute, Vancouver, Canada

2 Department of Health Care and Epidemiology, University of British Columbia, Vancouver, Canada

3 SF Coordinating Center, San Francisco, USA

For all author emails, please log on.

BMC Medical Research Methodology 2008, 8:35  doi:10.1186/1471-2288-8-35

The electronic version of this article is the complete one and can be found online at: http://www.biomedcentral.com/1471-2288/8/35


Received:11 December 2007
Accepted:30 May 2008
Published:30 May 2008

© 2008 Kuramoto et al; licensee BioMed Central Ltd.

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

Background

Evidence-based medicine has been advanced by the use of standards for reporting the design and methodology of randomized controlled trials (RCT). Indeed, without this information it is difficult to assess the quality of evidence from an RCT. Although a variety of statistical methods are available for the analysis of recurrent events, reporting the effect of an intervention on outcomes that recur is an area that remains poorly understood in clinical research. The purpose of this paper is to outline guidelines for reporting results from RCTs where the outcome of interest is a recurrent event.

Methods

We used a simulation study to relate an event process and results from analyses of the gamma-Poisson, independent-increment, conditional, and marginal Cox models. We reviewed the utility of regression models for the rate of a recurrent event by articulating the associated study questions, preenting the risk sets, and interpreting the regression coefficients.

Results

Based on a single data set produced by simulation, we reported and contrasted results from statistical methods for evaluating treatment effect from an RCT with a recurrent outcome. We showed that each model has different study questions, assumptions, risk sets, and rate ratio interpretation, and so inferences should consider the appropriateness of the model for the RCT.

Conclusion

Our guidelines for reporting results from an RCT involving a recurrent event suggest that the study question and the objectives of the trial, such as assessing comparable groups and estimating effect size, should determine the statistical methods. The guidelines should allow clinical researchers to report appropriate measures from an RCT for understanding the effect of intervention on the occurrence of a recurrent event.

Background

Evidence-based medicine has been advanced by the use of standards for reporting the design and methodology of randomized controlled trials (RCT). Indeed, without this information it is difficult to assess the quality of evidence from an RCT. An increasing number of journals demand that submissions adhere to the Consolidated Standards for Reporting Trials (CONSORT) guidelines for improving report quality [1]. However, there are not yet available guidelines for reporting results from RCTs in which the subject may experience the same event multiple times during follow-up. Examples of recurrent events include falls, fractures, certain cancers, infections, chronic disease exacerbations, and hospitalizations [2-7]. Through a trial, clinical researchers attempt to determine whether the study outcome occurs more frequently in the intervention group than in the control group. In such trials clinicians are interested in a variety of questions, such as "How many events does the intervention prevent, on average, compared to the control?"; "Does the intervention decrease the event rate over the study period compared to the control?"; "What is the effect of intervention on the rate of subsequent event among those who experienced the preceding event?"; and "What is the protective effect of intervention on the rate of higher-order events compared to the control?"

Although a variety of statistical methods are available for the analysis of recurrent events, reporting the effect of an intervention on outcomes that recur is an area that remains poorly understood in clinical research [8,9]. Appropriate statistical techniques are not always used to analyze RCTs on recurrent falls [9]. Extensive work involving simulation studies based on varying event processes and case studies have compared recurrent event methods to illustrate their strengths and weaknesses [10-13]. Such methods include the gamma-Poisson model, and several extensions of the Cox proportional hazards model, including the independent-increment, marginal, and conditional models [14-20].

The purpose of this paper is to outline guidelines for reporting results from a trial of treatment that prevents a recurrent event. As an example, we are using the rationale of a randomized trial on falls prevention. Falls are the most common cause of injury among elderly people. One in three persons over the age of 65 falls at least once each year and this proportion increases to one in two people over the age of 80 [21,22]. Almost half of those who fall experience the event recurrently [23,24]. The goal of RCTs is to reduce the occurrence of falls with specific interventions strategies such as multi-factorial intervention, strength and balance retraining, medication rationalization and expedited cataract surgery.

In the Methods section we review the utility of regression models for the rate of a recurrent event by articulating the associated study questions, presenting the risk sets, and interpreting the regression coefficients. Based on a single data set produced by simulation, we report and contrast results from statistical methods for evaluating treatment effect from an RCT with a recurrent outcome in the Results section. Finally, we summarize our guidelines for reporting evidence from RCTs on recurrent events.

Methods

In this section, we relate study questions of interest in RCTs to methods for modelling recurrent event data. Recurrent event models were developed to account for potential dependence among observations within a subject. One approach allows for unobserved heterogeneity which is unmeasured, intraclass correlation where subjects have constant but unequal probabilities of experiencing the event [25]. Three other models, which were developed for the analysis of continuous time recurrent event data, are extensions of the Cox proportional hazards model. They first fit a Cox model that ignores dependence and then use the empirical sandwich estimator to adjust standard errors for the parameter estimates [17,18,20]. Several authors argued for a conditional approach that estimates the rate of kth event among those who have already experienced (k - 1) events [18,26]. This approach addresses the issue of constant susceptibility in a more natural way than marginal models [18,27]: while the association between event times remains unspecified, the event-specific rate functions condition on having had previous events.

There are substantial differences among the models described in this section, but all estimate the effect of factors on the occurrence and time to event while accounting for the dependence between observations. The methods that we review model the rate function, λ(t)-that is, the average intensity of a recurrent event at a certain time. We highlight differences in the model assumptions, risk sets, and rate ratio interpretation. The data structure required to fit each model is shown to illustrate the different risk sets, indicating which patients are considered to be at risk for events at certain times [25,28]. Examples of SAS code (SAS System version 9.1 for Windows, SAS Institute Incorporation, Cary, NC, USA) to fit each model are also presented.

Mean cumulative function

"How many events does the intervention prevent, on average, compared to the control?" is one study question in an RCT on recurrent events that could be addressed using the mean cumulative function (MCF). The MCF shows the population mean number of recurrent events by certain times [29]:

MCF(t) = E{N(t)}.

where N(t) is a random variable for the number of events that have occurred up to time t. The MCF curve changes as a function of time and its derivative gives the rate function, that is

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M1','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M1">View MathML</a>

The rate and intensity functions quantify different aspects of the recurrent event process: the intensity is the instantaneous risk of a recurrent event and the rate is the average intensity at time t [25,30]

λ(t)dt = E[dN(t)],

where dN(t) denotes the number of events in a small interval [t, t + dt).

We interpret the difference in MCFs between the intervention and control groups as an indicator of how many events the intervention would prevent, on average, by a certain time [31].

Gamma-Poisson model

A common study question for an RCT on recurrent events is "Does the intervention decrease the event rate over the study period compared to the control?", for which the gamma-Poisson model has been used. The gamma-Poisson model evaluates the relationship between the number of recurrent events and factors of interest when the data deviate from the Poisson model [15,16]. This model allows variation of the event rate among subjects in the same group according to an unobserved random variable, frailty, which defines how likely a subject is to experience the event compared to the average rate [16]. When the frailty follows a gamma distribution and a time homogeneous model is assumed then the marginal distribution of the total number of events is negative binomial [15].

Suppose Ni(t) counts the number of events that have occurred up to time t for subject i. Under the time-homogeneous, gamma-Poisson model, Ni(t) has a Poisson distribution with rate function

λi(t) = μi exp{α0 + βxi},(1)

where μi come from a gamma distribution with density function

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M2','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M2">View MathML</a>

(2)

In model 1, α0 is the logarithm of the baseline rate for the event, μi is the unobserved frailty for subject i, xi is a covariate value for subject i, β is the regression coefficient, and t represents the time from start of observation.

The expected value and variance of the frailty random variable is 1 and θ, respectively. Subjects with μi greater than 1 are considered more "frail" or more likely to experience the event at a higher rate; whereas, those with μi less than 1 are considered to experience the event at a lower rate [16].

Compared to the Poisson model which assumes the mean and variance for the number of events are equal, the gamma-Poisson model has an additional parameter which allows for over-dispersion. For a given set of covariates, this model assumes the expected number of events is t exp(α0 + βxi) and the variance is t exp(α0 + βxi) + θt2 exp(α0 + βxi)2 [32].

The rate function of any event for subject i averaged over the gamma-distribution is

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M3','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M3">View MathML</a>

(3)

Subjects are at risk of an event until they are censored. Suppose xi is a binary indicator of group membership, with value 0 if subject i belongs to the control group and 1 if the intervention group. Then, exp <a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M4">View MathML</a> from model 3 estimates the common rate ratio of event in the intervention group relative to the control. We interpret rate ratios less than 1 as indicating the overall rate of event, that is the rate of any event, in the intervention group is 100 [1 - exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M4">View MathML</a>]% lower than in the control.

The data structure for this model requires one record for each subject, regardless of the number events experienced. This record contains the total follow-up time and total number of events per subject. The data structure required for this model is illustrated through an example. Suppose subject 1 in the control group experiences a recurrent event at day 126, 216, and 314 from study start and is followed up for 365 days. In addition, subject 2 in the intervention group, who was followed for the same period of time, had events at day 42 and 350. Under the time-homogeneous gamma-Poisson model, the data for these subjects are represented as shown in Table 1. In this data set, pid is the subject identifier, time is the total follow-up time, nevent is the total number of events experienced, grp is the covariate for group membership, and logtime is the natural logarithm of time.

Table 1. Data structure for the time-homogeneous gamma-Poisson model

For these data, SAS can be used to fit a time-homogeneous gamma-Poisson model:

PROC GENMOD;

   MODEL nevent = grp/LINK = LOG DIST = NEGBIN OFFSET = logtime;

RUN;

A major limitation of the time-homogeneous gamma-Poisson model is it assumes that the recurrent event rate is constant over time, which is unlikely to hold in practice. Extensions to this model have been made to relax the independent increment assumption for recurrent events and the specification of the within subject correlation between recurrence times. For example, the general frailty model assumes that the counting process is a non-homogeneous Poisson process given the frailty and covariates, where the frailty is not restricted to follow a gamma distribution [33]. The proportional mean and rate model relaxed the non-homogeneous Poisson assumption for the counting process and directly models means and rates [17].

Independent-increment model

The study question "Does the intervention decrease the event rate over the study period compared to the control?" is also addressed by Lin's independent-increment model for the rate of recurrent events [17]. Originally this model was developed by Andersen and Gill to specify the intensity of a counting process with a Cox-type link function [14]. Lin et al. provided a rigorous formalization of the marginal rate model, which relaxes the assumption that the event history, Fi(t), can be completely described by time-dependent covariates, xi(t), that is, [17,30]

E[dNi(t)|Fi(t)] = E[dNi(t)|xi(t)].

In contrast to Cox's model where subjects are at risk of an event until its occurrence or they are censored, in the independent-increment model subjects still remain at risk after an event occurs. Unlike the gamma-Poisson model, the independent-increment model does not assume the recurrent event rate is constant over time. This model assumes that the number of events in disjoint time intervals are independent [27].

Under the independent-increment model, the rate function, λi(t), of any event for subject i is

λi(t) = Yi(t)λ0(t) exp{βxi(t)},(4)

where

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M5','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M5">View MathML</a>

In model 4, Yi is the at risk indicator of event for subject i, λ0(t) is the baseline rate function for the event, xi is a covariate value, which may be time-dependent but may not contain elements of the event history, for subject i, β is the regression coefficient, and t represents the time from start of observation.

From model 4 we observe that both the baseline rate functions, λ0, and regression parameters, β, are assumed to be common across events.

Subjects are at risk of the an event until they are censored. Suppose xi is a binary indicator of group membership, with value 0 if subject i belongs to a control group and 1 if an intervention group. Then exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M4">View MathML</a> estimates the common rate ratio of event for the intervention group relative to the control. The rate ratio is assumed to be constant over time and common across recurrent events. We interpret rate ratios less than 1 as indicating the overall rate of event in the intervention group is 100 [1 - exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M4">View MathML</a>]% lower than in the control. This model has a similar interpretation to the gamma-Poisson model except we no longer require the assumption of time-homogeneity or gamma distributed frailty.

Under the independent-increment model, the data for these subjects use the counting process format, where each subject is represented by a set of time intervals and event indicators. We illustrate these data in Table 2 using the example described in the Gamma-Poisson model subsection. In this data set, pid is the subject identifier, tstart is time of previous event or study start, tstop is time of event or censoring, status is an indicator of event, and grp is the covariate for group membership. Subject 1 experienced 3 events and then was censored at the end of follow-up, so there are 4 corresponding records for this subject. In contrast, subject 2 experienced 2 events before being censored, so there are only 3 records.

Table 2. Data structure for the independent-increment model

The corresponding SAS code to fit an independent-increment model is as follows:

PROC PHREG COVM COVS(AGGREGATE);

   WHERE (tstart < tstop);

   MODEL (tstart, tstop) * status(0) = grp/RISKLIMITS;

   ID pid;

RUN;

Conditional models

RCTs on recurrent events provide insight into the study question "What is the effect of intervention on the rate of subsequent event among those who experienced the preceding event?", which a condtional model can address. Pepe and Cai proposed the conditional model for the rate of recurrent events, where subjects are not considered to be at risk for event until all previous events have occurred [18].

Under the total, follow-up time conditional model, the rate function, λij(t), of the jth event for subject i is

λij(t) = Yij(t) λ0j(t) exp{βjxi(t)},(5)

where

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M6','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M6">View MathML</a>

From model 5 we observe that both the baseline rate functions, λ0j(t), and regression parameters, βj, can vary across events. The covariate xi may not contain elements of the event history.

In model 5, t represents the time from start of observation. The conditional model can also be formulated in terms of "gap time", the time from previous event:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M7','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M7">View MathML</a>

(6)

where

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M8','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M8">View MathML</a>

and <a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M9','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M9">View MathML</a> is the time of the event just prior to time t.

In contrast to the marginal model, subjects are considered at risk for an event at time t only if the previous event occurred before that time and they are still under observation. Suppose xi is a binary indicator of group membership, with value 0 if subject i belongs to a control group and 1 if an intervention group. Then, exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M10','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M10">View MathML</a> from model 5 estimates the event-specific rate ratio of the jth event from study start in the intervention group relative to the control, conditional on experiencing the previous events. The event-specific rate ratio for the jth event from model 6 represents the rate of the jth event from the time of the previous event in the intervention group relative to the control. We interpret rate ratios less than 1 as indicating that among those who experienced j - 1 events, the intervention reduces the rate of the jth event by 100[1 - exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M10','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M10">View MathML</a>]% compared to the control. While the conditional model using total follow-up time compares subjects who experienced the same number of events and have the same follow-up from study start, the gap-time conditional model compares subjects who have experienced the same number of events and have the same duration since their previous event.

Fitting these conditional models relies on creating the appropriate data sets. These data sets are illustrated through the example presented in Gamma-Poisson model subsection. Under the conditional model for total follow-up, the data set for these subjects follows the counting process format as shown in Table 3. Similar to the independent-increment model (equation 4), the number of records representing each subject depends on the number of events experienced. The data structure differs from that of the independent-increment model since we have a variable for the event number.

Table 3. Data structure for the conditional model for total follow-up time

Assuming that the most number of events observed per subject was seven, the corresponding SAS code for fitting a conditional, total follow-up time model is as follows:

PROC PHREG;

   MODEL (tstart, tstop) * status(0) = group1-group7/RISKLIMITS;

      group1 = grp * (event = 1);

      group2 = grp * (event = 2);

      group3 = grp * (event = 3);

      group4 = grp * (event = 4);

      group5 = grp * (event = 5);

      group6 = grp * (event = 6);

      group7 = grp * (event = 7);

   STRATA event;

RUN;

Under the conditional, gap time model, the data set for these subjects requires times between adjacent events, as shown in Table 4. Again, the number of records per subject depends on the number of events experienced. As opposed to time intervals, times between subsequent events are required.

Table 4. Data structure for the conditional model for gap time

Assuming that the most number of events observed per subject was seven, the corresponding SAS code for fitting a conditional, gap time model is as follows:

PROC PHREG;

   MODEL gaptime * status(0) = group1-group7/RISKLIMITS;

      group1 = grp * (event = 1);

      group2 = grp * (event = 2);

      group3 = grp * (event = 3);

      group4 = grp * (event = 4);

      group5 = grp * (event = 5);

      group6 = grp * (event = 6);

      group7 = grp * (event = 7);

   STRATA event;

RUN;

In these conditional model data sets, pid is the subject identifier, tstart is time of previous event or study start, tstop is time of event or censoring, gaptime is the time to event from previous event, event is the event number, status is an indicator of event, and grp is the covariate for group membership.

Marginal model

"What is the protective effect of intervention on the rate of higher-order events compared to the control?" is an important study question to help decide whether to start treatment. This question is addressed by the marginal model, proposed by Wei, Lin and Weissfeld, which allows for different effects on each subsequent event [20]. This model treats the ordered event like an unordered competing risk problem [27]. Estimates from the marginal model have a practically useful interpretation which allows comparison between groups at treatment onset [34].

Under the marginal model, the rate function, λij(t), of the jth event for subject i is

λij(t) = Yij(t)λ0j(t) exp{βjxi(t)},(7)

where

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M11','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M11">View MathML</a>

In model 7, Yij, is the at risk indicator of the jth event for subject i, λ0j(t) is the baseline rate function for the jth event, xi is a covariate value, which may be time-dependent, for subject i, βj is the regression coefficient for event j, and t represents the time from start of observation. From model 7 we observe that both the baseline rate functions, λ0j, and regression parameters, βj, can vary across events.

Subjects are at risk of the jth event until it occurs or they are censored. Furthermore, subjects are considered to be at risk for the jth event even if they did not yet experience the (j - 1)th event. Suppose xi is a binary indicator of group membership, with value 0 if subject i belongs to a control group and 1 if an intervention group. Then, exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M10','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M10">View MathML</a> estimates the average event-number-specific rate ratio of the jth event in the intervention group relative to the control. We interpret rate ratios less than 1 as indicating the transition rate from 0 to j events in the intervention group is 100 [1 - exp<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M10','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M10">View MathML</a>]% lower than in the control. The marginal event-number-specific rate ratios indicate whether subjects in the intervention group will have fewer higher-order events of a certain number from the time of treatment onset [34].

The data structure required for this model is illustrated through the example presented in the Gamma-Poisson model subsection. We would like to study the effect of intervention on the first four events. Under the marginal model, the data set for these subjects show times of event from study start for all events under study, as shown in Table 5. In this data set, pid is the subject identifier, tstart is time of study start, tstop is time of event or censoring, event is the event number, status is an indicator of event, and grp is the covariate for group membership. Both subjects are represented by the same number of records, namely four since we are interested in the first four events.

Table 5. Data structure for the marginal model

The corresponding SAS code to fit this marginal model is as follows:

PROC PHREG COVS(AGGREGATE);

   MODEL tstop*status(0)=group1-group4/RISKLIMITS;

      group1 = grp * (event = 1);

      group2 = grp * (event = 2);

      group3 = grp * (event = 3);

      group4 = grp * (event = 4);

   STRATA event;

   ID pid;

RUN;

Results

Using available statistical instruments for recurrent events, we report results from a simple simulation study of falls prevention to illustrate the utility of the methods. Although each of the models being compared has already been studied via simulation, we contrast reporting results in the context of an RCT based on a single data set. The measures discussed are the rate ratios from the recurrent event models described in the Methods section. These include the common rate ratio, which compares the average rate of event in the intervention group to the control, the conditional event-specific rate ratios, which summarize the effect of intervention on a specific event conditional on experiencing previous events, and the marginal event-number-specific rate ratios, which summarize the intervention effect on the transition rate of experiencing a certain number of events from study start. In addition, we report the event rate, a measure of the average number of event accrued per person-time, and the mean cumulative function (MCF), a measure of the average number of events experienced per subject within a certain time.

We simulated recurrent falls in two groups, control and intervention, using Matlab Version 7 software (see Additional file 1). Each group had 250 subjects, and all subjects were followed for 365 days. Fall rates were based on those observed in an RCT [35]. Times between falls were assumed to follow an exponential distribution with falls rates specified for each fall. In the control group the fall rates for all falls were held constant at 7.7 falls per 1000 person-days. In the intervention group the fall rate was 5.3 falls per 1000 person-days for the first fall, and changed to 3.3 for all subsequent falls. Dependence within subjects was modelled using a gamma frailty distribution with density function given in equation 2 and variance θ = 0.10. We report the effect of the first 4 falls only since higher-order event-specific estimates are unreliable when there are only a few subjects with a large number of falls [25,27].

Additional file 1. Matlab code to simulate recurrent falls data used in Results section

Format: TXT Size: 4KB Download fileOpen Data

Event rates

After 1 year, the control group had 675 falls, nearly double that of the intervention group with 373 falls. The total follow-up time in each group was 91,250 person-days. The average observed fall rates in the control and intervention groups were 7.4 (95%CI 6.8–8.0) and 4.1 (95%CI 3.7–4.5) falls per 1000 person-days, respectively. Compared to the control group, the rate of falls in the intervention was almost halved, a crude approximation of the anticipated effect size. This effect size can be used to design RCTs on recurrent events, specifically for determining the number of subjects.

Mean cumulative function, MCF

Figure 1 shows the MCF by group, estimated by a non-parametric estimator [36]:

thumbnailFigure 1. Estimated mean cumulative function (MCF) of falls by group (upper panel), their difference (lower panel), and 95% confidence intervals.

<a onClick="popup('http://www.biomedcentral.com/1471-2288/8/35/mathml/M12','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/8/35/mathml/M12">View MathML</a>

where ej is the number of events at time tj, nj-1 is the number of subjects at risk just beyond time tj-1, and j indexes the observed event times. A subject is at risk of event until the end of follow-up. At one year of follow-up, an average of 2.7 and 1.5 falls per subject were experienced in the control and intervention group, respectively. Both MCFs were approximately linear, which indicates that the rate of falls is relatively constant in each group [31,36]. The control group experienced more falls and had a higher fall rate than the intervention group. On average, the control group experienced 1 more additional fall by 301 days (Figure 1). From the MCF difference, we observed that 1.2 falls were prevented per year on average for each subject.

Common rate ratios

The time-homogeneous gamma-Poisson and independent-increment gave similar common rate ratio estimates of 0.55 (95% CI 0.48–0.63) and 0.55 (95% CI 0.48–0.62), respectively (Table 6). The gamma-Poisson and independent-increment models both infer that the rate of any fall in the intervention group is 45% lower in the intervention group than control. In practice the assumption of a constant recurrent event rate over time may not hold, so the independent-increment model is preferred over the time-homogeneous gamma-Poisson model. These common rate ratios indicate that the intervention had an impact on the risk of falls; however, it does not inform whether the effect changes for subsequent events.

Table 6. Effect of intervention on recurrent falls, as measured by common rate ratios and 95% confidence intervals

Conditional event-specific rate ratios

The majority of the control group experienced two falls within 1 year of follow-up: 228, 180, 122, and 77 subjects had fall 1, 2, 3, and 4, respectively. The number of falls in the intervention group was lower: 202, 104, 45, and 18 subjects had fall 1 to 4, respectively (Table 7). Higher-order events, up to 7 falls, were experienced by 38 subjects in the control group; whereas, in the intervention group, only 4 subjects had the highest-order event of 5 falls. In the conditional model, the risk set for a subsequent fall consisted of only subjects who experienced the previous falls, and total follow-up time decreased for later events. The crude rate ratios indicate a similar intervention effect on falls 2 and 3.

Table 7. Fall-specific characteristics for total events, number of subjects at risk, total follow-up in days, and crude rate ratios, as indicated by the marginal and conditional total time models

As expected, the rate ratios for the first fall from the conditional models give identical estimates, 0.68 (95% CI 0.57–0.83), since the total follow-up time and gap time to first falls refer to the same period (Table 8). For subsequent falls, the fall-specific rate ratios from the conditional models overlap and remain relatively constant ranging from 0.46 (95% CI 0.36–0.59) to 0.53 (95% CI 0.31–0.88). The rate ratio for fall 5, 0.38 (95% CI 0.13–1.07), may be unreliable due to the number at risk for this event, and effects could not be estimated for falls 6 or 7. Among subjects who experienced preceeding falls, the effect of intervention on the rate of the first four recurrent falls did not differ (Wald χ2 test = 6.6, df = 3, p = 0.08 for total follow-up time model, and Wald χ2 test = 6.7, df = 3, p = 0.08 for gap-time model).

Table 8. Effect of intervention on recurrent falls, as measured by fall-specific rate ratios and 95% confidence intervals

For recurrent falls, the rate ratios from the conditional, total follow-up time model indicate that conditional on experiencing the previous fall, the rate of second, third and fourth falls from study start are 54%, 47% and 50% lower in intervention than control. The rates of falls from the time of previous fall are 54%, 47%, and 47% lower in intervention than control, as estimated from the conditional, gap time model. The conditional models provide evidence of the constant difference in recurrent fall rates between the groups. The conditional fall-specific rate ratios evaluate how the intervention affected the rate of kth fall among those who experienced k - 1 falls.

For both the conditional total follow-up time model and conditional gap time model, subjects are considered to be at risk for an event only if the previous event occurred, so subjects at risk may not consist of all who were intially randomized. The number of subjects at risk for subsequent events should be reported to allow evaluation of how different the treatment groups are from the start of the study (Table 7).

Marginal event-number-specific rate ratios

In the marginal model, all subjects were considered to be at risk for the 1st, 2nd, 3rd, 4th, and higher-order falls regardless of experiencing previous events (Table 7). Subjects are at risk for a specific fall until its occurrence or censoring, so the total follow-up time accumulates over subsequent falls. The crude rate ratios decrease with fall events.

The fall-number-specific rate ratios decrease from 0.68 (95% CI 0.57–0.83) for fall 1 to 0.20 (95% CI 0.12–0.34) for fall 4 (Table 8). For higher-order events, the rate ratio for fall 5 was 0.10 (95% CI 0.03–0.27) and could not be estimated for falls 6 or 7. The marginal model indicated that there was a difference in the average effect of intervention on the first four falls (Wald χ2 test = 32.2, df = 3, p < 0.0001). Rate ratios based on the marginal model indicated that, on average, the transition rate from zero falls at the start of treatment to one, two, three and four falls were 32%, 58%, 70% and 80% lower, respectively, in the intervention group than the control. These rate ratios do not imply that the effect of intervention increased with recurrent falls. Rather, the marginal fall-number-specific rate ratios indicate that subjects in the intervention group will have fewer events overall.

Given an objective of an RCT is to compare groups which are similar in all aspects except for the treatment of interest, it is appropriate to use the marginal model since all subjects are considered to be at risk for each number-specific event from study start. In contrast, the groups being compared to evaluate the effect of subsequent events in the conditional models may not consist of all subjects initially randomized.

Discussion

Recurrent events arise in many contexts, such as falls in seniors considered in this paper. In evidence-based medicine there is increasing need for guidelines on what to report in the analysis of recurrent events [8]. In the Results section we have outlined briefly statistical methods for evaluation of treatment effect from an RCT with a recurrent outcome. These should allow clinical researchers to report appropriate measures from an RCT for understanding the effect of intervention on the occurrence of a recurrent event.

We used a simulation study to relate an event process and results from analyses of the gamma-Poisson, independent-increment, conditional, and marginal Cox models [15-18,20]. We showed that each model has different study questions, assumptions, risk sets, and rate ratio interpretation, and so inferences should consider the appropriateness of the model for the RCT. The gamma-Poisson and independent-increment models compare the common event rates between groups, with the assumption of independence of the number of events across time intervals being required in the latter, but not the former. The conditional model distinguishes between first and recurrent events, and conditions on having had previous events. In contrast, the marginal model treats the events as unordered, and all subjects are at risk for any event. In different trials the outcomes of interest and validity of assumptions will differ. Our guidelines for reporting results from an RCT involving a recurrent event suggest statistical methods which correspond to the objectives of the trial, such as addressing the study question of interest, assessing comparable groups and estimating effect size. First, the average event rate by intervention group is a measure of the average number of events accrued per person-time. These event rates serve an important role in determining sample size and follow-up time for the design of future RCTs involving recurrent events [37]. Second, the MCF by intervention group provides a measure of the average number of events experienced per subject within a certain time. The MCF allows us to determine how many events per subject the intervention would prevent, on average, compared to the control group [31]. Third, the common rate ratio, as measured by the gamma-Poisson and independent-increment models, quantifies the average rate of event in the intervention group relative to the control group. This rate ratio provides an estimate of the common effect size, thereby indicating whether the intervention had an impact on the event occurrence. Fourth, conditional event-specific rate ratios, which quantify the rate of the kth event in the intervention relative to the control, conditional on experiencing preceding events, should be reported. These rate ratios allow us to evaluate how the effect of intervention changes, if at all, on subsequent events. Lastly, we suggest reporting the marginal event-number-specific rate ratios, which represent the rate of transitioning to higher-order events from the start of treatment in the intervention group relative to the control group. These rate ratios allow us to evaluate the overall protective effect of intervention. For methods used in the assessment of goodness of fit for each model we refer the reader to the corresponding papers [17,27].

It has been argued that the average event rate might have little relevance in the context of recurrent events because this measure does not acknowledge dependence between events experienced by a subject [38]. However, by applying appropriate statistical methods for recurrent events we can make valid inferences on rates. Extensive simulation studies based on varying event processes and case studies have compared recurrent event methods to determine their strengths and weaknesses [10-13].

Regression methods for the analysis of recurrent events is not limited to modelling the rate of event. The mean number of recurrences can be modelled using semi-parametric Cox models and parametric models [17,39]. Proportional rates and proportional means models are equivalent when the rate only depends on covariates that do not directly impact the occurrence of event, namely external covariates [17,40]. Regression models for the intensity function, which condition on event history, are also available [14,19]. However, in RCTs treatment may affect event history, so conditioning on the event history may underestimate the treatment effect [41].

Conclusion

Our guidelines for reporting results from an RCT involving a recurrent event suggest that the study question and the objectives of the trial, such as assessing comparable groups and estimating effect size, should determine the statistical methods. Guidelines for reporting results from an RCT involving a recurrent event should allow clinical researchers to report appropriate measures for understanding the effect of intervention on the occurrence of a recurrent event.

Competing interests

The authors declare that they have no competing interests.

Authors' contributions

Study concept and design: BGS, LK. Analysis and interpretation: LK, BGS, MGD. Drafting of the manuscript: LK, BGS, MGD.

Acknowledgements

We are grateful to the reviewers for their insightful comments.

References

  1. Moher D, Schulz KF, Altman DG: The CONSORT statement: revised recommendations for improving the quality of reports of parallel-group randomised trials.

    Lancet 2001, 357:1191-1194. PubMed Abstract | Publisher Full Text OpenURL

  2. Close J, Ellis M, Hooper R, Glucksman E, Jackson S, Swift C: Prevention of falls in the elderly trial (PROFET): a randomised controlled trial.

    Lancet 1999, 353:93-97. PubMed Abstract | Publisher Full Text OpenURL

  3. Campbell AJ, Robertson MC, Gardner MM, Norton RN, Tilyard MW, Buchner DM: Randomised controlled trial of a general practice programme of home based exercise to prevent falls in elderly women.

    British Medical Journal 1997, 315:1065-1069. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  4. Hirte HW, Strychowsky JE, Oliver T, Fung-Kee-Fung M, Elit L, Oza AM: Chemotherapy for recurrent, metastatic, or persistent cervical cancer: a systematic review.

    Int J Gynecol Cancer 2007. PubMed Abstract | Publisher Full Text OpenURL

  5. Lavelle WF, Cheney R: Recurrent fracture after vertebral kyphoplasty.

    Spine J 2006, 6:488-493. PubMed Abstract | Publisher Full Text OpenURL

  6. Schokker S, Kooi EM, de Vries TW, Brand PL, Mulder PG, Duiverman EJ, Molen T: Inhaled corticosteroids for recurrent respiratory symptoms in preschool children in general practice: Randomized controlled trial.

    Pulm Pharmacol Ther 2007, 21:88-97. PubMed Abstract | Publisher Full Text OpenURL

  7. Metcalfe C, Thompson SG, Cowie MR, Sharples LD: The use of hospital admission data as a measure of outcome in clinical studies of heart failure.

    Eur Heart J 2003, 24:105-112. PubMed Abstract | Publisher Full Text OpenURL

  8. Robertson MC, Campbell AJ, Herbison P: Statistical analysis of efficacy in falls prevention trials.

    J Gerontol A Biol Sci Med Sci 2005, 60:530-534. PubMed Abstract | Publisher Full Text OpenURL

  9. Donaldson MG, Sobolev BG, Khan KM, Cook WL, Janssen PA: A systematic review of statistical methods reported in randomized controlled trials of falls prevention in older adults.

    Age & Ageing 2008.

    Under review

    OpenURL

  10. Boher J, Cook RJ: Implications of model misspecification in robust tests for recurrent events.

    Lifetime Data Analysis 2006, 12:69-95. PubMed Abstract | Publisher Full Text OpenURL

  11. Box-Steffensmeier JM, De Boef S: Repeated events survival models: The conditional frailty model.

    Statistics in Medicine 2006, 25:3518-3533. PubMed Abstract | Publisher Full Text OpenURL

  12. Metcalfe C, Thompson SG: The importance of varying the event generation process in simulation studies of statistical methods for recurrent events.

    Statistics in Medicine 2006, 25:165-179. PubMed Abstract | Publisher Full Text OpenURL

  13. Therneau TM, Hamilton SA: rhDNase as an example of recurrent event analysis.

    Statistics in Medicine 1997, 16:2029-2047. PubMed Abstract | Publisher Full Text OpenURL

  14. Andersen PK, Gill RD: Cox's regression model for counting processes: A large sample study.

    Annals of Statistics 1982, 10:1100-1120. Publisher Full Text OpenURL

  15. Cook RJ, Lawless JF: Analysis of repeated events.

    Statistical Methods in Medical Research 2002, 11:141-166. PubMed Abstract | Publisher Full Text OpenURL

  16. Klein JP, Moeschberger ML: Survival Analysis: Techniques for Censored and Truncated Data. 2nd edition. Springer; 2003. OpenURL

  17. Lin DY, Wei LJ, Yang I, Ying Z: Semiparametric regression for the rate and mean function of recurrent events.

    J Royal Statistical Society B 2000, 62:711-730. Publisher Full Text OpenURL

  18. Pepe MS, Cai J: Some graphical displays and marginal regression analyses for recurrent failure times and time dependent covariates.

    Journal of the American Statistical Association 1993, 88:811-820. Publisher Full Text OpenURL

  19. Prentice RL, Williams BJ, Peterson AV: On the regression analysis of multivariate failure time data.

    Biometrika 1981, 68:373-379. Publisher Full Text OpenURL

  20. Wei LJ, Lin DY, Weissfeld L: Regression analysis of multivariate incomplete failure time data by modeling marginal distributions.

    Journal of the American Statistical Association 1989, 84:1065-1073. Publisher Full Text OpenURL

  21. Tinetti ME, Speechley M, Ginter SF: Risk factors for falls among elderly persons living in the community.

    New England Journal of Medicine 1988, 319:1701-1707. PubMed Abstract OpenURL

  22. O'Loughlin JL, Robitaille Y, Boivin JF, Suissa S: Incidence of risk factors for falls and injurious falls among the community-dwelling elderly.

    American Journal of Epidemiology 1993, 137:342-354. PubMed Abstract | Publisher Full Text OpenURL

  23. Nevitt MC, Cummings SR, Kidd S, Black D: Risk factors for recurrent nonsyncopal falls. A prospective study.

    Journal of the American Medical Association 1989, 261:2663-2668. PubMed Abstract | Publisher Full Text OpenURL

  24. Campbell AJ, Borrie MJ, Spears GF: Risk factors for falls in a community-based prospective study of people 70 years and older.

    J Gerontol 1989, 44:M112-M117. PubMed Abstract OpenURL

  25. Ezell ME, Land KG, Cohen LE: Modeling multiple failure time data: A survey of variance-corrected proportional hazard models with empirical applications to arrest data.

    Sociological Methodology 2003, 33:111-167. Publisher Full Text OpenURL

  26. Klen JP, Goel PK, (Eds):

    Survival analysis: state of the art, Kluwer 1992 chap. Frailty models for multiple event times. OpenURL

  27. Therneau TM, Grambsch PM: Modeling Survival Data: Extending the Cox Model. Springer; 2000. OpenURL

  28. Kelly PJ, L-Y LL: Survival analysis for recurrent event data: an application to childhood infectious diseases.

    Statistics in Medicine 2000, 19:13-33. PubMed Abstract | Publisher Full Text OpenURL

  29. Nelson WB: Recurrent events data analysis for product repairs, disease recurrences, and other applications. 1st edition. ASA-SIAM; 2003. OpenURL

  30. Cai J, Schaubel DE:

    Handbook of Statistics, Advances in Survival Analysis, Elsevier 2004 chap. Analysis of Recurrent Event Data. 23 OpenURL

  31. Donaldson MG, Sobolev B, Kuramoto L, Cook WL, Khan KM, Janssen PA: Utility of the mean cumulative function in the analysis of fall events.

    J Gerontol A Biol Sci Med Sci 2007, 62:415-419. PubMed Abstract | Publisher Full Text OpenURL

  32. Lawless JF: Negative binomial and mixed Poisson regression.

    The Canadian Journal of Statistics 1987, 15:209-225. Publisher Full Text OpenURL

  33. Andersen PK, Borgan O, Gill RD, Keiding N: Statistical Models Based on Counting Processes. 1st edition. Springer-Verlag; 1993. OpenURL

  34. Metcalfe C, Thompson SG: Wei, Lin and Weissfeld's marginal analysis of multivariate failure time data: should it be applied to a recurrent events outcome?

    Statistical Methods in Medical Research 2007, 16:103-122. PubMed Abstract | Publisher Full Text OpenURL

  35. Donaldson MG: Falls risk in frail seniors: clinical and methodological studies. PhD thesis. University of British Columbia; 2007. OpenURL

  36. Nelson WB: Confidence limits for recurrence data: applied to cost or number of repairs.

    Technometrics 1995, 37:147-157. Publisher Full Text OpenURL

  37. Cook RJ: The design and analysis of randomized trials with recurrent events.

    Statistics in Medicine 1995, 14:2081-2098. PubMed Abstract | Publisher Full Text OpenURL

  38. Windeler J, Lange S: Events per person year-a dubious concept.

    British Medical Journal 1995, 310:454-456. PubMed Abstract | Publisher Full Text OpenURL

  39. Lawless J: Introductory overview lecture.

    2005.

  40. Kalbfleisch JD, Prentice RL: The statistical analysis of failure time data. 1st edition. John Wiley & Sons; 1980. OpenURL

  41. Schaubel DE, Zeng D, Cai J: A semiparametric additive rates model for recurrent event data.

    Lifetime Data Analysis 2006, 12:389-406. PubMed Abstract | Publisher Full Text OpenURL

Pre-publication history

The pre-publication history for this paper can be accessed here:

http://www.biomedcentral.com/1471-2288/8/35/prepub