Multiple time scales in modeling the incidence of infections acquired in intensive care units

Abstract

Background

When patients are admitted to an intensive care unit (ICU) their risk of getting an infection will be highly depend on the length of stay at-risk in the ICU. In addition, risk of infection is likely to vary over calendar time as a result of fluctuations in the prevalence of the pathogen on the ward. Hence risk of infection is expected to depend on two time scales (time in ICU and calendar time) as well as competing events (discharge or death) and their spatial location. The purpose of this paper is to develop and apply appropriate statistical models for the risk of ICU-acquired infection accounting for multiple time scales, competing risks and the spatial clustering of the data.

Methods

A multi-center data base from a Spanish surveillance network was used to study the occurrence of an infection due to Methicillin-resistant Staphylococcus aureus (MRSA). The analysis included 84,843 patient admissions between January 2006 and December 2011 from 81 ICUs. Stratified Cox models were used to study multiple time scales while accounting for spatial clustering of the data (patients within ICUs) and for death or discharge as competing events for MRSA infection.

Results

Both time scales, time in ICU and calendar time, are highly associated with the MRSA hazard rate and cumulative risk. When using only one basic time scale, the interpretation and magnitude of several patient-individual risk factors differed. Risk factors concerning the severity of illness were more pronounced when using only calendar time. These differences disappeared when using both time scales simultaneously.

Conclusions

The time-dependent dynamics of infections is complex and should be studied with models allowing for multiple time scales. For patient individual risk-factors we recommend stratified Cox regression models for competing events with ICU time as the basic time scale and calendar time as a covariate. The inclusion of calendar time and stratification by ICU allow to indirectly account for ICU-level effects such as local outbreaks or prevention interventions.

Background

The individual patients’ time at-risk (the length of stay at-risk in ICU) is a key determinant for ICU-acquired infections [1]. For infections which are based on rather exogenous acquisition routes [2], the risk might also depend on calendar time due to local outbreaks or fluctuations in the prevalence of other infectious patients, contaminated health-care worker in the same ICU or prevention interventions on the ICU-level [3]. In other words, the use of calendar time (in combination with stratification by ICU-level) allows to indirectly account for these transmission-associated effects. Thus, there are two time scales to be addressed in a risk factor analysis of ICU-acquired infections.

In time-to-event analysis, one basic time scale has to be selected. In general, possible choices for the basic time scale may include age, time since enrollment in a study, calendar time, or time since an event such as disease diagnosis. In many cases there may be no single time scale that is clearly more appropriate than others. The choice of time scale is, however, crucial: it affects the interpretation of the model and how risks and rates are assumed to vary over time; and in some cases different choices can even lead to apparently contradictory results [4–6]. The time scale therefore has to be chosen with care and the choice taken into account when interpreting results.

In the Cox proportional hazards model, the specification of the basic time scale plays a major role since time effects are not explicitly modeled as they are absorbed into the unspecified baseline hazard [4]. Thus, the underlying time scale provides the most flexible way to control for time effects. For instance, if the primary interest is in how a factor (such as a drug treatment) affects mortality for a disease, then a time scale of age might be most appropriate if the mortality rate is highly age-dependent but relatively unaffected by time since diagnosis, while a time scale of time since diagnosis might be more appropriate if the converse is true.

Time since enrollment (or time on study) is one of the most frequently used time scales, though there has been considerable debate about whether this choice has always been the most appropriate [4–6]. The use of calendar time as the basic time scale for longitudinal observational data (the so-called real-time approach) treats the population of interest as a dynamic population rather than a closed cohort [7, 8]. An application is, for instance, the impact of environmental exposure on pregnancy outcomes [6].

For ICU-acquired infections caused by transmissible pathogens such as MRSA and Vancomycin-resistant Enterococci, calendar time is often a natural choice when studying the effect of interventions on the ICU-level [9]. This is because hazards of acquiring the infection are likely to vary over calendar time as a result of fluctuations in the prevalence of the pathogen on the ward [3]. These fluctuations are typically unobserved as they result from asymptomatic carriage, making direct adjustment for the ICU-level prevalence impossible.

In addition to patient-individual characteristics, the risk of acquiring a MRSA infection in an ICU might also depend on spatio-temporal factors, i.e., when (calendar time) and where (which ICU) a patient requires intensive care.

The choice of calendar time as the basic time scale also controls for time-varying factors acting on the ICU-level such as changes in medical management, hygiene practices, patterns of antibiotic usage, staffing levels, and seasonal factors [3, 9–11].

The occurrence of ICU-acquired infection also depends on a second time scale, the patients’ individual time at-risk (i.e. time since patient admission to the ICU), with longer stays creating more opportunity for infection. This ICU exposure time is one the most important determinants for ICU-acquired infections and is frequently used for studying patient individual risk-factors such as age, morbidity, patient-individual antibiotic treatment or invasive devices [12]. Here, we discuss these two time scales which are crucial for the incidence of ICU-acquired Methicillin-Resistant Staphylococcus aureus (MRSA) infections.

The patients’ individual time at-risk ends with the occurrence of MRSA infection, ICU discharge or death in ICU since after the two latter events the risk of ICU-acquired infection is zero. Therefore, ICU discharge and death in ICU are competing events for ICU-acquired MRSA infections which should be considered in a risk factor analysis [12–14]. Ignoring these competing events can easily lead to heavily biased risk estimates [15] and wrong conclusions about the impact of risk factors [16]. Due to the presence of competing events, there are two metrics (the rate and the risk metric) in a risk factor analysis [17, 18]. Thus, for a complete analysis, it is necessary to perform event-specific hazard rate analyses (for MRSA infection, discharge and death) as well as a summary analysis for the cumulative risk of MRSA infection [19]. Furthermore, to account for the patients’ environment or geographical space, multi-level techniques are necessary [13].

The major aim of this paper is to find an appropriate model to study the incidence of MRSA infections by accounting for multiple time scales, competing risks and the hierarchical nature of the data. To do this, we explore, compare and combine the aforementioned time scales in a real ICU data setting. We calculate hazard rates with respect to the corresponding time scale and perform analyses based on the stratified Cox proportional hazards model to study patient-individual risk factors in a competing-risk framework.

Methods

Spanish ICU data

We used a multi-center data base from the Spanish surveillance network HELICS-ENVIN (http://hws.vhebron.net/envin-helics/), embedded in the HELICS project (Hospitals in Europe Link for Infection Control through Surveillance) [20]. We included ICUs which contributed to the registry between January 2006 and December 2011 and we included only patients who stayed at least two days in an ICU due to the definition of hospital-acquired infections. We excluded ICUs which contributed fewer than 500 patient admissions to the cohort to ensure a sufficient amount of patient time at risk for each ICU. The study population contains 81 intensive care units with 84,843 admissions (693,180 admission-days).

Statistical methods

In the Additional file 1 is a Lexis diagram [21] of individual patient data from one selected ICU over 100 days in calendar time. It demonstrates how the data depend on the two time scales. In the following, we compare the two time scales in several steps. For the ICU time scale, the time origin is the time of admission. For the calendar time scale, patient admissions entered the model with staggered or delayed entry with left-truncation occurring at the time of admission.

Overall hazard rates

We used a penalized likelihood approach [22] to estimate the overall hazard rates λ^k(t) separately for each event k: for ICU-acquired MRSA infection (the event of interest), death and discharge without MRSA. The overall hazard rates depend both on ICU or calendar time. The variation of the overall hazards due to different ICUs was accounted by using a shared frailty model [22]. More formally, let c represent the calendar time and c₀ the truncation time, i.e., the admission time in calendar time scale. Thus, for each competing event k (MRSA, death without MRSA, discharge without MRSA) and the i-th ICU, the event-specific hazard rate with a shared frailty term ${Z^{k}_{i}}$ is defined as

$$\text{ICU time: } \lambda^{k}(c-c_{0}|{Z^{k}_{i}})={Z^{k}_{i}} {\lambda_{0}^{k}}(c-c_{0}) $$

$$\text{calendar time: } \tilde{\lambda}^{k}(c|c_{0},\tilde{Z}^{k}_{i})=\tilde{Z}^{k}_{i} \tilde{\lambda}_{0}^{k}(c|c_{0}) $$

with baseline hazard ${\lambda _{0}^{k}}(.)$ (or $\tilde {\lambda }_{0}^{k}(.)$) for event k; the term |c⁰ denotes the left-truncation time. The frailty term ${Z^{k}_{i}}$ (or $\tilde {Z}^{k}_{i}$) is a random effect which varies across ICUs (patients within ICU share the same frailty) and is assumed to be Gamma distributed with shape parameter 1/θ _k and inverse scale parameter 1/θ _k , thus E(Z^k)=1 and Var(Z^k)=θ _k . Large values of θ _k signify a closer positive relationship between patients within ICU and greater heterogeneity across ICUs.

Patient-level risk factors

We used event-specific Cox proportional hazards models (rate metric) and a Fine & Gray model [23] (risk metric) to explore covariate effects of vector X (gender, age, type of diagnosis, antibiotic treatment 48 h before and/or after ICU admission, trauma, days in hospital before ICU admission, APACHE II (Acute Physiology And Chronic Health Evaluation) score) comparing results between models where the timescale was calendar time and ICU time. The assumption of proportional hazards was checked via the inspection of the Schoenfeld residuals [24]; note that proportionality due to the rate metric does not lead to proportionality due to the risk metric but even if proportionality is not fulfilled the hazard ratio has the meaningful interpretation of an time-averaged effect [25]. We then model both times together by including the second time scale as a covariate. We stratified for ICU in order to allow the hazard to be different across ICUs, and hence we did not use the frailty terms.

Models with one time scale

Model 1a: ICU time as basic time scale

$$\lambda^{k}(c-c^{0}|X)=\lambda_{0i}^{k}\left(c-c^{0}\right)\text{exp}\left(\sum_{j} {\beta^{k}_{j}} X_{j}\right) $$

This is an event-specific Cox model with ICU time as the basic time scale. The exponential of the regression coefficients ${\beta ^{k}_{j}}$ are corresponding hazard ratios of variable X _j and event k.

Model 2a: calendar time as basic time scale (the real-time approach)

$$\tilde{\lambda}^{k}(c|c^{0},X)=\tilde{\lambda}_{0i}^{k}(c|c^{0})\text{exp}\left(\sum_{j} \tilde{\beta}^{k}_{j} X_{j}\right) $$

This is an event-specific Cox model with calendar time as the basic time scale and staggered or delayed entry with left-truncation occurring at the time of admission.

Model 3a: subdistribution for MRSA, basic time scale is ICU time

In a competing risks setting, the cumulative incidence function for event k (CIF^k(t)) depends on all event-specific hazards [26]. This can be seen with following formula with t=c−c⁰:

$$\begin{aligned} \text{CIF}^{k}(t)={\int_{0}^{t}} \left(\text{exp}\left(-\sum_{\text{all events }i} {\int_{0}^{u}} \lambda^{i} (v) dv \right)\right) \times \lambda^{k} (u) du \end{aligned} $$

We are basically interested in MRSA infection as our event of interest. In simple words, the formula above is the product of the time-dependent probability of staying alive at-risk on ICU and the conditional probability of acquiring a MRSA infection. Note that the probability of staying alive at-risk on ICU depends on the competing events hazards (discharge or dying without MRSA) in addition to the MRSA hazard.

Fine & Gray [23] defined the subdistribution hazard which is in our setting the probability of the MRSA infection given that a patient has stayed in ICU up to time t without a MRSA infection or has had the competing event (death ore discharge) prior to time t [27]. Thus, the risk sets for the subdistribution hazard are unnatural (discharged and died patients remain technically at-risk). However, unlike the event-specific hazard, the subdistribution hazard is directly linked to the corresponding cumulative incidence function of MRSA infection. Based on the subdistribution hazard, Fine & Gray proposed a proportional hazards model to study risk factors on the risk metric. The resulting subdistribution hazard ratios of an exposure can be interpreted as effects which can be seen when plotting cumulative incidence functions, grouped by exposure categories.

Models with two time scales

Model 1b: model 1a plus year of admission as a covariate

$$\lambda^{k}(c-c^{0}|X)=\lambda_{0i}^{k}\left(c-c^{0}\right)\text{exp}\left(\sum_{j} {\beta^{k}_{j}} X_{j}+\gamma f\left(c^{0}\right)\right) $$

with f(c⁰) as the calendar year of admission. If necessary, other more detailed functions (which, for instance, includes also the calendar month of admission) can be chosen.

Model 2b: model 2a plus length of stay at-risk as a time-dependent covariate

$$\tilde{\lambda}^{k}(c|c^{0},X)=\tilde{\lambda}_{0i}^{k}\left(c|c^{0}\right)\text{exp}\left(\sum_{j} \tilde{\beta}^{k}_{j} X_{j}+\gamma g\left(c-c^{0}\right)\right) $$

with g(c−c⁰) as a function to categorize ICU time in 0-4, 5-9, 10-14 and 15+ days. If necessary, other categorizations can be chosen.

Model 3b: model 3a plus year of admission as covariate.

Results

Overall hazards

The overall hazard rates are displayed in Fig. 1. The bottom row is clearly showing what is happening (with respect to event risk) over a long time period (perhaps reflecting more what is going on at the ICU level) whereas the top row is showing what is happening at the patient level in the short term since admission to an ICU. With the ICU time scale, the hazard of MRSA acquisition increases from 0 to 0.0015 within the first 15 days from admission (perhaps due to the fact that an MRSA infection often follows a MRSA colonisation which itself also requires some exposure time in ICU). After 15 days, the MRSA hazard slightly decreases. The hazard reflects the (unadjusted) instantaneous risk for an ICU patient to acquire a ICU-acquired MRSA infection during his or her ICU stay. This hazard has to be interpreted jointly with the corresponding hazards of the competing events since the cumulative risk of MRSA infection depends also on the discharge and death hazards. The hazards of death or discharge without MRSA are also increasing within the first week after admission (see Fig. 1: up to 0.025 at day 7 for death and up to 0.25 at day 7 for discharge). After 10-15 days, the hazard of death without MRSA is about 0.02 and the hazard of discharge without MRSA is much lower (about 0.05) in contrast to earlier days.

Fig. 1

Event-specific hazard rates depending on ICU time (upper panels) or calendar time (lower panels)

The MRSA hazard rate with respect to calendar time shows how the hazard steadily decreases from the January 2006 to December 2010 on the ICU level (Fig. 1 bottom left), potentially due to the implementation of prevention strategies [28]. During the years 2006–2010, the death hazard rate without MRSA slightly increased from 0.014 to 0.017. The discharge hazard rate remains more or less constant (about 0.1) over calendar time. Details regarding the random effects of ICU are displayed in the Additional file 1.

Effects of patient individual risk-factors

The covariates of interest are listed in Table 1. Descriptive statistics indicate the following variables are associated with calendar time. The mean number of days in hospital before ICU admission steadily increased from about 9.5 days in year 2006 up to 18.5 days in year 2010+. The proportion of patient admissions with trauma steadily decreased from 11.5 % in 2006 to 7.1 % in 2010+. The proportion of admissions receiving antibiotic treatment 48 h before and/or after ICU admission steadily decreased from about 30 to 19 %.

Table 1 Description of study population

One time scale

In Table 2, the effects of several risk factors are listed from Models 1a and 2a. The hazard ratios do not significantly differ much (estimates for one model are generally within the 95 % confidence interval of the corresponding estimates in the other model) but they have different interpretations. If ICU time is used as the time scale, then the effects of interest are directly adjusted for ICU time effects (such as the patient-individual intensity of the underlying morbidity – generally the more critical care is needed the longer the ICU stay). By using calendar time, the effects of interest are directly adjusted for calendar time effects such as local outbreaks, changing hygiene practices or intervention strategies. Thus, the combination of using the calendar time scale and the additional stratification by ICU addresses the spatio-temporal dynamics of ICU-acquired MRSA infections. That means that the risk factors are studied beyond this spatio-temporal pattern.

Table 2 Results from multivariate analysis (one time scale). Event-specific regression analyses based on Cox proportional regression models as described in main text. Proportionality could be assumed after the inspection of Schoenfelds’ residuals. Event-specific hazard ratios are displayed with 95 %-confidence interval in brackets

The estimated effects in Table 2 are different according to the time scale used. This is especially true for the APACHE II score, type of diagnosis, antibiotic treatment and trauma where the effects are more pronounced when calendar time is taken as the basic time scale. One explanation is that these factors are highly associated with the severity of patients’ illness and the ICU time captured already part of this severity.

Two time scales

If calendar year is introduced in the model as an additional covariate (Table 3), the results for the patient-individual factors barely change (comparing effects between models 1a and b) even though calendar year is associated with all outcomes, MRSA as well as the competing events death or discharge. This consistency might be due to the noncorrelation of calendar and ICU time. Comparing models 2a and b, one can see that including the covariate length of stay at-risk in the model changes the hazard ratios; they are now very similar to those from model using ICU time scale (model 1b). Table 3 shows the results of the two models 1b and 2b accounting for both time scales simultaneously. Both models yield very similar results.

Table 3 Results from multivariate analysis (two time scales). Event-specific regression analyses based on Cox proportional regression models as described in main text. Proportionality could be assumed after the inspection of Schoenfelds’ residuals. Eventspecific hazard ratios are displayed with 95 %-confidence interval in brackets

Subdistribution approach

Since most covariates are also associated with the competing events for MRSA, it is also required to study the covariate effects on the cumulative incidence function (CIF) [19, 23]. The CIF depends on the underlying time scale. The CIF is a very relevant quantity when the ICU time scale is used since it reflects how the individual patient risk accumulates with time spent in the ICU (Fig. 2). In contrast, displaying the cumulative incidence function with the calendar time scale is not useful. One way to account for calendar time is adjusting for the admission year as a covariate (model 3b). Table 4 shows how the covariates are associated with the subdistribution hazard of MRSA as our event of interest. The importance of competing events in such a risk factor analysis can be seen when comparing the results of model 3b with those from model 1b. Depending on the metric used (rate or risk), the effects can be very different. For instance, APACHE II score is moderately associated with an increased hazard rate of MRSA infection (e.g., HR=1.52 (95 %-CI: 1.00-2.32)) for APACHE II score >31 vs. 0-10) but highly associated with an increased risk of MRSA infection (subdistribution HR=5.79 (95 %-CI: 3.82-8.76)) for APACHE II score >31 vs. 0-10). This phenomena can be explained by considering the effects on the competing events: patients with higher APACHE II scores are associated with an increased death hazard (they die faster) but also with an decreased discharge hazard (they stay longer at ICU). Since most patients are discharged, the actual at-risk time is prolonged for patients with higher APACHE II scores. Therefore, more MRSA infections occur in this patient group since longer stays create more opportunity for infection. The difference between the two metrics can also go in the other direction: patients older than 80 years have a (non-significantly) lower hazard to acquire a MRSA infection than patients aged between 61-80 years (HR=0.73 (95 %-CI: 0.47-1.11)). However, the cumulative risk for MRSA infection is much lower (subdistribution HR=0.47 (95 %-CI: 0.31-0.71)) since the older patient group die faster without MRSA (HR=1.60 (95 %-CI:1.52-1.70)) and is discharged faster (HR=1.12 (95 %-CI: 1.09-1.15)) meaning that the at-risk time for MRSA infection is reduced.

Fig. 2

Cumulative incidence function of ICU-acquired MRSA depending on ICU time scale and stratified per calendar year on admission

Table 4 Results from multivariate analysis (subdistribution analysis). Regression analyses based on adapted Cox proportional regression model (Fine & Gray model) as described in main text. Proportionality could be assumed after the inspection of Schoenfelds’ residuals. Subdistribution hazard ratios are displayed with 95 %-confidence interval in brackets

The inclusion of the calendar year of admission in the subdistribution analysis has no impact on the patient-level covariates. Even though the effect of admission year is strong, the estimates from model 3a and 3b are almost the same.

Discussion

In this paper, we compared two important time scales (ICU and calendar time) for modeling the incidence of ICU-acquired MRSA infections. Both time scales have a strong influence on the hazard rate and cumulative risk of MRSA infections in ICUs but also of the competing events (death and discharge). We showed that hazard ratios of patient individual risk-factors of MRSA infections can differ depending on the underlying time scale in Cox regression models. This difference can be overcome by using both time scales simultaneously.

One strength of this study was the application of advanced statistical methods on a large data base which was necessary to model calendar time and patient-individual characteristics simultaneously. We used competing risks models based on the semi-parametric Cox proportional hazards model and shared frailty models even though other parametric frailty models could be applied in such settings. In a recent review about analyses of hospital-acquired infection risk factors, Brown and colleagues also emphasized the need to adjust for the at-risk time and period effects (calendar time) [1]. In addition to the time scales, we gave emphasis on competing risks since the interpretation and conclusions of the results might be very different [18].

This study has limitations. First, the data have been collected from volunteer ICUs, thus data can be subject to reporting, information or selection bias. Second, part of the ICUs contributed the whole study period from January 2006 and December 2011 whereas others contributed only 3 months per year (April to June) and some ICUs started their continuous contribution at some time between 2006 and 2011. This might affect our findings, particular for calendar time. Third, it can be assumed that transmission dynamics play the main role in MRSA colonisation (which might later lead to MRSA infection). In our setting, we were not able to study MRSA colonisation rates because such an evaluation requires regular swabs from patients to detect asymptomatic carriage. Instead, we focused on MRSA infections because this is both the clinically important outcome and the outcome for which the best data are available. Fourth, we studied only time-independent risk factors in the subdistribution approach. Time-dependent risk factors (such as antibiotic treatment or invasive devices during ICU) can be introduced in the event-specific analyses but there are challenges in interpreting results from risk metric approaches in presence of time-dependent risk factors [18].

The choice of the basic time scale has been discussed in the statistical literature. Andersen and Keiding [29] stated that the underlying time scale should be chosen ’for which the variation in the hazard is unknown or is expected to be dramatic and a parametric description is less important’. According to Pencina et al. [30], one should ask the question: ’Which approach seems to better capture the nature of the data that are to be modeled, and which is more suitable to answer the research questions?’.

Given these general statements and based on our findings regarding MRSA infections, ICU time should be chosen as the basic time scale when studying patient individual risk-factors. To indirectly account for exogenous factors (such as local outbreaks or implementation of prevention strategies on the ICU-level), we recommend to use calendar time as a covariate in Cox regression models and stratify by ICU. A further advantage of this time scale choice is that the cumulative incidence functions of MRSA infections can be studied in a competing event framework.

Conclusions

Risk factor analyses of general ICU-acquired infections are already complex due to the presence of competing risks during the time in ICU as well as unit-level effects. Moreover, the analysis of ICU-acquired infections requires the involvement of calendar time. By accounting for both time scales, we believe that this approach provided deeper insights into the disease-risk association since the additional use of calendar time allows to indirectly account for transmission-associated effects.

Abbreviations

HR, hazard ratio; ICU, intensive care unit; MRSA, methicillin-resistant staphylococcus aureus; sHR, subdistribution hazard ratio