- Research article
- Open Access
- Open Peer Review
Statistical methods for elimination of guarantee-time bias in cohort studies: a simulation study
BMC Medical Research Methodologyvolume 17, Article number: 126 (2017)
Aspirin has been considered to be beneficial in preventing cardiovascular diseases and cancer. Several pharmaco-epidemiology cohort studies have shown protective effects of aspirin on diseases using various statistical methods, with the Cox regression model being the most commonly used approach. However, there are some inherent limitations to the conventional Cox regression approach such as guarantee-time bias, resulting in an overestimation of the drug effect. To overcome such limitations, alternative approaches, such as the time-dependent Cox model and landmark methods have been proposed. This study aimed to compare the performance of three methods: Cox regression, time-dependent Cox model and landmark method with different landmark times in order to address the problem of guarantee-time bias.
Through statistical modeling and simulation studies, the performance of the above three methods were assessed in terms of type I error, bias, power, and mean squared error (MSE). In addition, the three statistical approaches were applied to a real data example from the Korean National Health Insurance Database. Effect of cumulative rosiglitazone dose on the risk of hepatocellular carcinoma was used as an example for illustration.
In the simulated data, time-dependent Cox regression outperformed the landmark method in terms of bias and mean squared error but the type I error rates were similar. The results from real-data example showed the same patterns as the simulation findings.
While both time-dependent Cox regression model and landmark analysis are useful in resolving the problem of guarantee-time bias, time-dependent Cox regression is the most appropriate method for analyzing cumulative dose effects in pharmaco-epidemiological studies.
Extensive studies have elaborated and documented protective effect of aspirin for the primary prevention of cardiovascular diseases, including myocardial infarction (MI), stroke, coronary artery disease (CAD) [1,2,3,4,5,6]. Distinct from other non-steroidal anti-inflammatory drugs (NSAIDs), aspirin has the capacity to irreversibly inhibit cyclooxygenase (COX) and suppress thromboxane production. Findings of experimental studies have suggested that COX-2 suppression [7,8,9] or anti-platelet effect of aspirin [10, 11] could interact with cancer cells and play a significant role in blocking tumor angiogenesis, invasiveness and metastatic potential. Accordingly, several observational studies have reported an inverse association between aspirin use and cancer incidence or mortality. However, the results have been inconsistent and controversy still remains regarding its anti-cancer effects [12,13,14,15], in part due to variations in study designs and statistical methods for analyzing the effectiveness of aspirin. For example, Fraser et al.  conducted a population-based study of breast cancer patients in Scotland, in which Cox’s proportional hazard models was used to assess the relationship between aspirin use and survival. The authors found significantly reduced all-cause (HR = 0.53, 95% CI = 0.45 to 0.63) and breast cancer-specific mortality (HR = 0.42, 95% CI = 0.31 to 0.55) among participants who consumed aspirin following a diagnosis of breast cancer. Similarly, using Cox proportional hazards modeling, Jacobs et al.  reported that total cancer incidence is significantly lower among regular long-term aspirin users (RR = 0.84, 95% CI = 0.76 to 0.93). In contrast, a study by McMenamin and colleagues  found a non-significant decrease in the risk of lung cancer-specific mortality (p = 0.5975) with aspirin use. In this study, the authors used a time-dependent Cox model whereby low-dose aspirin usage was treated as time-varying covariate, with users not considered to be exposed until after a lag of six months following initial prescription. Further, Cardwell et al. , in their conditional logistic regression analysis of nested case-control data, also found no evidence for protective effects of low-dose aspirin usage against colon cancer-specific mortality (OR = 1.06; 95% CI = 0.92 to 1.24).
The above-mentioned studies of aspirin use show that when evaluating drug effects, results are often confounded by time-related biases that tend to exaggerate the protective effects, even when the drug has no or little effect . Time-related biases can be further categorized into guarantee-time bias (also known as immortal time bias or time-dependent bias), time-window bias, or time-lag bias. For the purpose of this paper, we focus on guarantee-time bias, which is frequently encountered in time-to-event analyses evaluating drug effects. Guarantee-time bias refers to bias arising from incorrect handling of the time from the beginning of follow-up to the first treatment exposure . The bias occurs if treated patients are assumed to have already been treated from time of cohort entry. For example, in a study of patients with heart disease, heart transplant is a time-varying treatment. While some patients may receive a transplant at the start of follow-up, others may undergo transplant long after the beginning of follow-up or die before undergoing the transplant. Failure to account for time-varying feature of the treatment can result in biased estimation of the treatment effect, in this case, in favor of the treatment group .
This study aimed to evaluate the extent of guarantee-time bias in estimation of drug effects on disease outcomes. We conducted simulation studies in which time-fixed, time-dependent Cox model, and landmark analysis were compared for handling guarantee-time bias. For illustrative purposes, effect of cumulative rosiglitazone dose on the risk of hepatocellular carcinoma was investigated using the Korean National Health Insurance Data.
Guarantee-time bias in cohort studies can distort the results in favor of the treatment group, depending on the type of event (if it benefits the patient or not) [21, 23, 24]. In detail, consider this study as a cohort study and the situation where drug effect and event occurrence are independent. Event is defined as the disease of interest. Define the random variables W and T 0 as the time to initiation of drug usage and time to occurrence of event. The binary random variable Z is defined to be Z = I[W ≤ T 0], which represents drug usage. The random variables T 0 and T 1 are defined as the time to event conditional on Z = 0 or 1 respectively, i.e. T = (1 − Z)T 0 + ZT 1. As shown in Fig. 1, the person whose value of T 0 is smaller than W is allocated to the non-user group (person #1). Otherwise, individuals are allocated to the drug-user group, according to cumulative drug exposure (persons #2 to 4). In other words, to have received the treatment implies that the subject survived or was event-free, up to the initiation of drug use. As regards drug-user groups, the longer the survival times, the higher the probability of belonging to the high-dose group. This is the guarantee-time bias phenomenon, i.e. subjects must have lived long enough to receive the treatment and consume large amounts of drug.
Statistical modeling for Guarantee-Time bias
Our statistical modeling closely follows the paper by Nam and Zelen . For convenience, we consider the situation where the drug usage is binary. The probability density functions of W , T 0 , T 1 will be denoted by g(w) , q 0(t) and q 1(t) respectively. The survival functions will be denoted by G(w) = Pr[W > w] , Q 0(t) = Pr[T 0 > t] and Q 1(t) = Pr[T 1 > t]. Note that by definition Z = 1 if the time for drug usage is observed. If the drug usage does not cause a change in the survival function, the conditional density function of non-user group f(t| z = 0) and drug-user group f(t| z = 1) would be the same.
Under the hypothesis that q 0(t) = q 1(t), if we simply compare the survival functions of the two groups according to the drug usage, the two groups of conditional density functions are not equal to f(t| z = 0) ≠ f(t| z = 1), as it has been proven in the paper by Nam and Zelen . Therefore, it is not appropriate to simply compare the survival functions of two groups according to drug usage.
In the landmark method, selected landmark time τ 0 is set, and binary random variable Z(τ 0) is defined as Z(τ 0) = I(W < τ 0| T > τ 0). The conditional density functions for Z(τ 0) can be computed as follows :
where A(t) = 1 − G(t).
Hence if q 0(t) = q 1(t), it can be seen that the two functions are the same as
As a result, it is possible to make a valid comparison of two survival functions between the drug-user and non-user groups using the landmark method. In this study, landmark analysis was conducted using time-fixed Cox regression with drug exposure status defined prior to landmark time. One thing to mention in the landmark method is that it is very important to carefully select the landmark time.
Time-dependent Cox regression model
In a time-dependent Cox regression model, a time-dependent covariate in the model tracks whether the classifying event has occurred during the estimation process . This method eliminates guarantee-time bias by using drug usage as a time-dependent covariate; subjects are classified as unexposed until the start of drug usage and exposed thereafter [26,27,28]. The time-dependent Cox regression model has the advantage of using all study follow-up data since it starts analysis at the time of cohort entry. By including all data, this method has increased statistical power over the landmark method .
We conducted extensive simulations to compare the performance of our three methods, Cox regression, time-dependent Cox regression and landmark method, in terms of empirical type I error and power. We used a discrete-time survival model to describe the effect of cumulative drug dose on the outcome using three indicator variables Z 1(t) , Z 2(t) , Z 3(t) which represent low, moderate and high-drug dose groups. We set N as total sample size and d as the number of drug-users. Among N subjects, d subjects consumed drug for different durations and the dosage was fixed at a constant dose of 0.5 per unit of time. We assumed that drug is taken continuously without interruption until the end of the study or event.
We divided the total observation period (10-year) into ten intervals as needed for the discrete-time survival model. The overall simulation setting is represented in Fig. 2. We assumed that the drug user rate is 5%. First, we randomly assigned time to initiation of drug use between 0 and 10 for each drug-user. Using these initiation times, we calculated cumulative dose at each time t and categorized each subject into four groups: non-user group and low(Z 1(t) = 1), moderate(Z 2(t) = 1), high-dose(Z 3(t) = 1) groups. At each discrete time point, the probability of disease occurrence at time t is as follows:
If the event occurred within the observation period, survival time is the point of event occurrence. On the other hand, if the event did not occur, survival time is equal to 10 years and censored. For the landmark analysis, we considered subjects whose starting point of drug intake is before landmark time as drug-users, and otherwise non-user. We used 5 and 7 as landmark times, respectively. Subjects who died or who were censored before the landmark time were excluded from the analysis.
Simulations were performed under two different scenarios. In the first scenario, type I error rates for each of the statistical methods were evaluated under the null hypothesis, whereby β 1, β 2 and β 3 values were set to zero to represent a case of no significant protective effect of drug on disease outcome. For the power comparisons in scenario 2, we simulated data from the alternative hypothesis (i.e. cumulative drug dose is associated with disease outcome) for each of the three β 0 values denoting disease incidence: log(0.015), log(5∗0.015) , and log(10∗0.015). For each simulation setting, 1000 replication datasets were generated. All simulations were performed using the R statistical software .
The empirical type I errors of three methods under the null hypothesis of no association between cumulative drug dose and disease are illustrated in Table 1. Each value represents the total number of cases out of 1000 replications for which the null hypothesis is incorrectly rejected. As expected, Cox regression model with fixed covariate values tended to produce inflated type I error rates, especially in the moderate and high dose groups. Time-dependent Cox regression and landmark method, on the other hand, displayed well-controlled type I error rates; nearly all values remained close to the nominal level of significance of 0.05, regardless of the degree of disease incidence.
The power estimates for the alternative hypothesis are displayed in Table 2. Statistical power of the time-dependent Cox regression approach was higher than that of landmark analyses in the high-dose group. In the case of moderate and low-dose groups, the landmark method using landmark time 5 demonstrated higher statistical power compared to the time-dependent Cox regression. This is attributed to the fact that there was extremely small number of cases for the high-dose group, compared with the moderate group.
Table 3 summarizes the bias and mean squared errors (MSE) of hazard ratios for estimation of the association between cumulative drug dose and outcome. Time-dependent Cox regression had much lower bias than the landmark method and generally showed smaller MSE values compared with the other methods. For instance, for the low-dose group, when the β 0 value was log(0.015), the bias of the time-dependent Cox regression model was 0.0009, while the bias of the landmark analyses were 0.0673 and 0.0416 for τ = 5 and τ = 7, respectively.
Real data example
We applied our simulation findings to a real data from the Korean National Health Insurance Database (NHID). The sample was followed for 12 years from 2002 to 2013. The event of interest was incidence of hepatocellular carcinoma (HCC). Among 47,738 patients with incident diabetes, 203 hepatocellular carcinoma cases were identified. Full details of the study design have been described elsewhere [31, 32]. In the current study, total prescribed doses of rosiglitazone were calculated at each year and summed to produce cumulative doses of follow-up years. Discrete-time survival analysis was used to take into account of the intermittent administration of rosiglitazone and the varying dosage between patients. Based on cumulative dose of rosiglitazone use during the study period, drug users were categorized into low (<1350 mg), moderate (1350–4499 mg), high groups (≥4500 mg). In the landmark analysis, sixth year was chosen as a landmark time since the number of incident HCC was balanced at that time point. Data analysis was carried out using the SAS statistical software version 9.4 .
Results based on the NHID data are illustrated in Table 4. In the Cox regression, the risk of HCC incidence was lowest among subjects exposed to high cumulative doses of rosiglitazone (HR = 0.443, 95% CI = 0.218 to 0.899). However, the protective pattern of dose-response relationship and effect sizes were considerably attenuated in the time-dependent Cox regression and landmark analysis.
In this study, we examined the phenomenon of guarantee-time bias through statistical modeling and simulation study. Specifically, our simulation study assessed the performance of three methods, namely Cox regression, time-dependent Cox regression and landmark method based on time-fixed Cox regression. These methods depend upon the proportional hazard assumption. According to our simulation results, time-fixed Cox regression was shown to be vulnerable to guarantee-time bias [20, 22]. Pharmaco-epidemiological studies typically involve the use of time-varying exposures, such as cumulative dose. Hence, under such situations, applying the time-fixed Cox regression approach can induce bias due to model misspecification. Our results are in accordance with results previously reported by Mi et al. . However, unlike previous studies, we performed a simulation by including cumulative dose groups as exposures.
Landmark analysis has been suggested in previous studies as a simplified alternative to time-dependent Cox regression for elimination of guarantee-time bias in time-to-event data. In our simulation, landmark method was comparable to the time-dependent Cox regression method in terms of type I error. But landmark analysis tended to slightly increase MSE. One possible explanation for these results could be the sample size. For landmark analysis to produce efficient estimates, it is imperative that optimal landmark time is specified. As illustrated by our results, if the landmark point is too early, there is a greater possibility of imbalanced observations among treatment groups. On the other hand, if the landmark point is too late, a significant proportion of events may be omitted, giving rise to insufficient number of cases to achieve adequate power. Thus, it is important that landmark studies are designed with sufficient number of participants to maintain adequate statistical power. Additionally, the landmark method will produce estimates with minimal bias conditional upon the treatment being evenly distributed across the follow-up . Recently, the use of landmark super models, a pooled summary analysis of several landmarks, has been advocated to remedy the problem of low statistical power related to the landmark method [34, 35]. Further studies are warranted to compare the performance of landmark super models with the time-dependent Cox model.
Nevertheless, the landmark approach has several advantages over the time-dependent Cox regression model. Most notably, this method has the advantage of computational simplicity because drug use is defined as a time-fixed covariate by using landmark time. Thus, one can visualize the survival curve of drug users using the Kaplan-Meier method. But, unconditional Kaplan-Meier estimates of time-varying status of drug users are unstable because the number of drug users can be quite small at early time point .
In conclusion, to avoid guarantee-time bias in observational studies of drug effects, we recommend incorporating time-dependent exposure status in the analysis. While both time-dependent Cox regression model and landmark analysis were found to be useful in resolving the problem of guarantee-time bias, time-dependent Cox regression was the most appropriate method for analyzing cumulative and long-term drug exposure. We recommend the time-dependent Cox regression for estimating hazard ratios of cumulative doses. Alternatively, due to its graphical capabilities, the landmark method may be a suitable alternative for visualizing survival curves for treatment groups.
Mean squared error
National Health Insurance Database
Steering Committee of the Physicians' Health Study Research Group. Final report on the aspirin component of the ongoing physicians study. N Engl J Med. 1989;321:129–35.
Berger JS, Roncaglioni MC, Avanzini F, Pangrazzi I, Tognoni G, Brown DL. Aspirin for the primary prevention of cardiovascular events in women and men: a sex-specific meta-analysis of randomized controlled trials. JAMA. 2006;295(3):306–13.
Manson JE, Stampfer MJ, Colditz GA, Willett WC, Rosner B, Speizer FE, Hennekens CH. A prospective study of aspirin use and primary prevention of cardiovascular disease in women. JAMA. 1991;266(4):521–7.
Pignone M, Alberts MJ, Colwell JA, Cushman M, Inzucchi SE, Mukherjee D, Rosenson RS, Williams CD, Wilson PW, Kirkman MS. Aspirin for primary prevention of cardiovascular events in people with diabetes. Circulation. 2010;121(24):2694–701.
Ridker PM, Cook NR, Lee I-M, Gordon D, Gaziano JM, Manson JE, Hennekens CH, Buring JE. A randomized trial of low-dose aspirin in the primary prevention of cardiovascular disease in women. N Engl J Med. 2005;352(13):1293–304.
Seshasai SRK, Wijesuriya S, Sivakumaran R, Nethercott S, Erqou S, Sattar N, Ray KK. Effect of aspirin on vascular and nonvascular outcomes: meta-analysis of randomized controlled trials. Arch Intern Med. 2012;172(3):209–16.
Hwang I-C, Jeon J-Y, Kim Y, Kim HM, Yoon YE, Lee S-P, Kim H-K, Sohn D-W, Sung J, Kim Y-J. Association between aspirin therapy and clinical outcomes in patients with non-obstructive coronary artery disease: a cohort study. PLoS One. 2015;10(6):e0129584.
Ofman P, Petrone AB, Peralta A, Hoffmeister P, Albert CM, Djousse L, Gaziano JM, Rahilly-Tierney CR. Aspirin Use and Risk of Atrial Fibrillation in the Physicians' Health Study. J Am Heart Assoc. 2014;3(4):e000763.
Ekström N, Cederholm J, Zethelius B, Eliasson B, Fhärm E, Rolandsson O, Miftaraj M, Svensson A-M, Gudbjörnsdottir S. Aspirin treatment and risk of first incident cardiovascular diseases in patients with type 2 diabetes: an observational study from the Swedish National Diabetes Register. BMJ Open. 2013;3(4):e002688.
King A, Shipley M, Markus H, Investigators A. The effect of medical treatments on stroke risk in asymptomatic carotid stenosis. Stroke. 2013;44(2):542–6.
Azoulay L, Dell’Aniello S, Simon TA, Langleben D, Renoux C, Suissa S. A net clinical benefit analysis of warfarin and aspirin on stroke in patients with atrial fibrillation: a nested case–control study. BMC Cardiovasc Disord. 2012;12(1):49.
Cuzick J, Otto F, Baron JA, Brown PH, Burn J, Greenwald P, Jankowski J, La Vecchia C, Meyskens F, Senn HJ. Aspirin and non-steroidal anti-inflammatory drugs for cancer prevention: an international consensus statement. The lancet oncology. 2009;10(5):501–7.
Dubé C, Rostom A, Lewin G, Tsertsvadze A, Barrowman N, Code C, Sampson M, Moher D. The use of aspirin for primary prevention of colorectal cancer: a systematic review prepared for the US Preventive Services Task Force. Ann Intern Med. 2007;146(5):365–75.
Tsujii M, Kawano S, DuBois RN. Cyclooxygenase-2 expression in human colon cancer cells increases metastatic potential. Proc Natl Acad Sci. 1997;94(7):3336–40.
Tsujii M, Kawano S, Tsuji S, Sawaoka H, Hori M, RN DB. Cyclooxygenase regulates angiogenesis induced by colon cancer cells. Cell. 1998;93(5):705–16.
Fraser D, Sullivan F, Thompson A, McCowan C. Aspirin use and survival after the diagnosis of breast cancer: a population-based cohort study. Br J Cancer. 2014;111(3):623–7.
Jacobs EJ, Thun MJ, Bain EB, Rodriguez C, Henley SJ, Calle EE. A large cohort study of long-term daily use of adult-strength aspirin and cancer incidence. J Natl Cancer Inst. 2007;99(8):608–15.
Mc Menamin ÚC, Cardwell CR, Hughes CM, Murray LM. Low-dose aspirin and survival from lung cancer: a population-based cohort study. BMC Cancer. 2015;15(1):911.
Cardwell CR, Kunzmann AT, Cantwell MM, Hughes C, Baron JA, Powe DG, Murray LJ. Low-dose aspirin use after diagnosis of colorectal cancer does not increase survival: a case–control analysis of a population-based cohort. Gastroenterology. 2014;146(3):700–8. e702
Suissa S, Azoulay L. Metformin and the risk of cancer. Diabetes Care. 2012;35(12):2665–73.
Suissa S. Immortal time bias in pharmacoepidemiology. Am J Epidemiol. 2008;167(4):492–9.
Lévesque LE, Hanley JA, Kezouh A, Suissa S. Problem of immortal time bias in cohort studies: example using statins for preventing progression of diabetes. BMJ. 2010;340:b5087.
Giobbie-Hurder A, Gelber RD, Regan MM. Challenges of guarantee-time bias. J Clin Oncol. 2013;31(23):2963–9.
Suissa S. Immortal time bias in observational studies of drug effects. Pharmacoepidemiol Drug Saf. 2007;16(3):241–9.
Nam CM, Zelen M. Comparing the survival of two groups with an intermediate clinical event. Lifetime Data Anal. 2001;7(1):5–19.
Suissa S. Effectiveness of inhaled corticosteroids in chronic obstructive pulmonary disease: immortal time bias in observational studies. Am J Respir Crit Care Med. 2003;168(1):49–53.
Samet JM: Measuring the effectiveness of inhaled corticosteroids for COPD is not easy! In.: Am Thoracic Soc; 2003.
Suissa S, Ernst P. Bias in observational study of the effectiveness of nasal corticosteroids in asthma. J Allergy Clin Immunol. 2005;115(4):714–9.
Mi X, Hammill BG, Curtis LH, Lai ECC, Setoguchi S. Use of the landmark method to address immortal person-time bias in comparative effectiveness research: a simulation study. Stat Med. 2016;35(26):4824–36.
Ihaka R, Gentleman R. R: a language for data analysis and graphics. J Comput Graph Stat. 1996;5(3):299–314.
Han E, Jang S-Y, Kim G, Lee Y-h, Choe EY, Nam CM, Kang ES. Rosiglitazone use and the risk of bladder cancer in patients with type 2 diabetes. Medicine. 2016;95(6):e2786.
Kim G, Jang SY, Han E, Lee Y. Park Sy, Nam CM, Kang ES: Effect of statin on hepatocellular carcinoma in patients with type 2 diabetes: A nationwide nested case-control study. Int J Cancer. 2017;140(4):798–806.
SAS Institute Inc. Base SAS® 9.4 procedures guide: statistical procedures. Cary: SAS Institute Inc; 2014.
Van Houwelingen HC. Dynamic prediction by landmarking in event history analysis. Scand J Stat. 2007;34(1):70–85.
Nicolaie M, Houwelingen J, Witte T, Putter H. Dynamic prediction by landmarking in competing risks. Stat Med. 2013;32(12):2031–47.
Klein JP, Moeschberger ML: Survival analysis: techniques for censored and truncated data: Springer Science & Business Media; 2005.
This research received no specific grant from any funding agency in the public, commercial or not-for profit sectors.
Availability of data and materials
The data used in this study was simulated. The simulation program is in R and available from https://github.com/jihyeon1/yonsei_biostat.
Ethics approval and consent to participate
This study involved an analysis of publicly available data, which was used as an illustrative example for theoretical development and simulation study. Therefore, ethical approval was not obtained from the institutional review board.
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.