Time-varying associations between an exposure history and a subsequent health outcome: a landmark approach to identify critical windows

Background Long-term behavioral and health risk factors constitute a primary focus of research on the etiology of chronic diseases. Yet, identifying critical time-windows during which risk factors have the strongest impact on disease risk is challenging. To assess the trajectory of association of an exposure history with an outcome, the weighted cumulative exposure index (WCIE) has been proposed, with weights reflecting the relative importance of exposures at different times. However, WCIE is restricted to a complete observed error-free exposure whereas exposures are often measured with intermittent missingness and error. Moreover, it rarely explores exposure history that is very distant from the outcome as usually sought in life-course epidemiology. Methods We extend the WCIE methodology to (i) exposures that are intermittently measured with error, and (ii) contexts where the exposure time-window precedes the outcome time-window using a landmark approach. First, the individual exposure history up to the landmark time is estimated using a mixed model that handles missing data and error in exposure measurement, and the predicted complete error-free exposure history is derived. Then the WCIE methodology is applied to assess the trajectory of association between the predicted exposure history and the health outcome collected after the landmark time. In our context, the health outcome is a longitudinal marker analyzed using a mixed model. Results A simulation study first demonstrates the correct inference obtained with this approach. Then, applied to the Nurses’ Health Study (19,415 women) to investigate the association between body mass index history (collected from midlife) and subsequent cognitive decline (evaluated after age 70), the method identified two major critical windows of association: long before the first cognitive evaluation (roughly 24 to 12 years), higher levels of BMI were associated with poorer cognition. In contrast, adjusted for the whole history, higher levels of BMI became associated with better cognition in the last years prior to the first cognitive interview, thus reflecting reverse causation (changes in exposure due to underlying disease). Conclusions This approach, easy to implement, provides a flexible tool for studying complex dynamic relationships and identifying critical time windows while accounting for exposure measurement errors. Supplementary Information The online version contains supplementary material available at (10.1186/s12874-021-01403-w).


Background
Long-term lifestyle, environmental or occupational exposures may have a major impact on the occurrence of chronic diseases. Yet, identifying critical time-windows during which risk factors have the strongest impact on disease risk is challenging. Chronic diseases often result from a long and accumulating pathological process that evolves over years before diagnosis [1]. In such context, the exposure time-windows close to the clinical event may be less meaningful in terms of etiology and, importantly, may be obscured by reverse causality (which occurs when behaviors or exposures change as the disease progresses in infra-clinic stages). For example, while obesity in midlife is a risk factor for subsequent dementia, it is consistently reported as a protective factor later in life; indeed advanced neuropathology, by altering olfactory function and taste, potentially leads to malnutrition and weight loss in late-life [2][3][4][5]. Thus, the only valuable approach to evaluate causal associations linking cumulative unhealthy body weight to cognitive aging might be to estimate the relationship between the entire history of exposure that precedes and begins well upstream of the period at risk of the event.
Exposure history metrics have been developed to estimate the cumulative effect of time-varying exposure on disease endpoints [6,7]. However, to be relevant in lifecourse epidemiology, such methodologies absolutely need to handle: (i) associations that can vary according to the age during exposure or the distance between the exposure and the disease endpoint, (ii) exposures that are measured only at sparse time points and with error, and (iii) exposure and disease endpoint time-windows that do not necessarily coincide. Finally, they should also apply to longitudinal disease outcomes (not only survival or binary) and their change over time such as cognitive decline that is a strong predictor of dementia. This work specifically aimed at extending the methodology of exposure history metrics to address all these important issues.
The most common exposure history metric is the cumulative index of exposure (CIE) [6,7]. Computed as the un-weighted sum of all past exposures, CIE assumes that past values of exposure are of equivalent importance (see Fig. 1, Scenario A) which may induce etiologically inadequate conclusions when the effect of the exposure on health outcomes likely varies according to the age or timing of exposures (see Fig. 1, Scenarios B and C). To address this challenge, Breslow et al. [8] and Thomas [9] introduced the concept of weighted CIE (WCIE), that combines information about duration, intensity and timing of the exposure. WCIE represents the time-weighted sum of past exposures, with weights reflecting the relative importance of exposures at different times. Therefore, unlike CIE, the estimated effects of the exposure history on health status are time-varying, a key aspect for the identification of critical windows of exposure in lifecourse epidemiology. Weights may be assigned a priori by using a known parametric weight function in the presence of biologically/clinically relevant assumptions [10][11][12] or when assessing very specific assumptions such as different life course hypotheses (e.g., progressive accumulation, critical period, recent period) [13,14]. Nevertheless, in most cases, the form of the weight function is not known and needs to be directly determined from the data. WCIE models are found in a wide range of applications with binary and time-to-event outcomes, such as concentration of occupational asbestos in relation to mesothelioma [15,16] or drug use linked to fall-related injuries [17]. In environmental health research, complex associations have been explored through distributed lag models (DLM) [18], which are very similar to WCIE models in the context where long-term exposure and outcome are collected at identical times across individuals.
A major limit of CIE, WCIE and even DLM methodologies is that they always require the exposure history to be complete. Yet, intermittent missing data on individual exposures are frequently reported in case-control and cohort studies, in particular when follow-up is long, thus preventing any analysis using these metrics. Second, although measurement error is inevitable, CIE and WCIE metrics also require the exposures to be measured without error and the violation of this assumption (which is likely in most epidemiological contexts) may lead to biased exposure-risk relationships [19]. To overcome these limitations, Mauff et al. [20] introduced a weighted cumulative association structure within a joint model to assess the impact of an endogenous time-varying exposure measured intermittently and with error on the risk of a time-to-event outcome; the exposure model accounted for independent centered Gaussian measurement errors through a standard linear mixed model. However, in this approach, exposure and outcome have to be evaluated during the same time-window, while as introduced previously, in some contexts the interest is in the cumulative effect of an exposure history that precedes evaluation of disease risk. In environmental health research, Chen et al. [21] proposed a reverse DLM, which models the exposure trajectory as a function of the health outcome rather than the inverse. By doing so, the reverse DLM naturally accounts for measurement errors and missing values in exposures. However, as a counterpart, this model conditions on the health outcome and thus considers it is measured without error. Finally, all these methodologies have been developed for survival or cross-sectional continuous or binary outcomes when we are interested in a longitudinal health outcome. To date, we are aware of only one pharmaco-epidemiological study applying the WCIE in association with a longitudinal outcome [22]. Nonetheless, it considered a complete error-free exposure history In this manuscript, we propose to extend the WCIE methodology to address life-course epidemiology research questions. We rely on a landmark approach which consists in limiting the analysis of the outcome of interest after a certain time of interest, called landmark, and considering as predictors the exposure history up to the landmark. Although mainly developed in dynamic prediction context [23,24], it also addresses the issue of non-concurrency between the exposure and the outcome time-windows central in life-long epidemiology.
In "Methods" section, we first introduce the proposed methodology, and describe the estimation procedure. In "Application to body mass index history starting at mid-life and subsequent cognitive trajectories after age 70" section, we apply the method in the Nurses' Health Study to investigate the association between BMI history starting in midlife and repeated cognition assessed after age 70. In "Simulation study" section, we demonstrate that the estimation procedure provides correct inference in a simulation study. Finally we discuss the results and methodology, and conclude.

General study framework
Within a prospective cohort study, we consider a sample of N individuals present in the cohort at a landmark time of interest from which the health outcome is to be studied, hereafter called time 0 (see Fig. 2). The temporal framework of the study consists in repeated measures of the exposure prior to the landmark time, that is from a negative time −S (e.g., time of entry in the cohort) up to time 0, as well as repeated measures of the outcome collected at and after the landmark time. Following the landmark framework [23], exposures data collected after time 0, if any, are not included. Both exposure and health outcome are collected at discrete and individualspecific occasions with error. The design of the Nurses' Health Study used in the application ("Application to body mass index history starting at mid-life and subsequent cognitive trajectories after age 70" section) naturally follows the landmark framework since the exposure was assessed several decades before the first assessment of outcome. In this manuscript, we will focus mainly on a longitudinal heath outcome and an exposure history constructed based on the time preceding the first outcome assessment. Nevertheless, the proposed approach also applies to other types of event, such as time-to-event and binary endpoints, and different timescales sur as age for the exposure and landmark age for the outcome (e.g., age at menopause).

Mixed model for the exposure
Let U il denote the exposure values collected for individual i (i = 1, ..., N) at the retrospective measurement time t Uil (l = 1, ..., m i ) before the landmark time so that t Uil ∈ [ −S, 0] (see Fig. 2). The times and the total number m i of repeated measures can differ from one individual to the other thus implicitly handling intermittent missing values. Although the methodology could apply beyond, we consider here a continuous exposure and describe its trajectory via a standard linear mixed effects model [25]: where U * i (t Uil ) is the true underlying exposure value at measurement time t Uil , X i (t Uil ) and Z i (t Uil ) are two vectors of covariates at time t Uil associated with the vectors of fixed effects β and of random effects b i , respectively; b i ∼ N (0, B). We consider random Gaussian measurement errors ε il with mean 0 and variance σ 2 ε ; b i and ε il are independent. The objective of this linear mixed model is to accurately model the underlying true exposure U * i (t U ) at any time t U in [ −S, 0] for inclusion in the WCIE. As such, both X i (t U ) and Z i (t U ) should include flexible functions of time, for instance a basis of splines [26] or fractional polynomials [27].

Weighted cumulative index of exposure
WCIE is defined as the weighted sum of the underlying true exposures during the period of interest. Without loss of generality, we consider the entire window [ −S, 0] although any sub-window could be considered instead: where w(t U ) is the weight function assigned to the history of the true exposure U * i (t U ) during the S + 1 time units (e.g., years) preceding the initial assessment of the outcome (see Fig. 1). Note that for ease of epidemiological interpretation, we chose to define WCIE as the sum of exposure levels at discrete times (e.g., every year) but the methodology also applies to WCIE defined as the integral of the exposure levels in [ −S, 0].

Model for the health outcome
Although transposable to cross-sectional or time-to-event outcomes, our primary interest was to estimate the timevarying effects of a past exposure history on a subsequent longitudinal health outcome (see Fig. 2). Thus, let's denote Y ij the measure for individual i (i = 1, ..., N) collected at time t ij ≥ 0 with j = 1, ..., n i . The times and the total number n i of repeated measures can differ from one individual to the other. Change over time of Y is modeled using a linear mixed model [25]. For ease of presentation, we consider here a linear trajectory over time, and thus introduce two time-varying effects of the exposure, one on the level at the landmark time (denoted WCIE I ) and one on the change over time of Y (denoted WCIE S ): whereX i (t ij ) is a vector of covariates at time t ij associated with the vectors of fixed effects α 1 and α 3 ; WCIE Ii and WCIE Si are the weighted cumulative exposure covariates associated to the initial level and to the slope of Y, respectively. They are defined according to Eq. 2 with different weights w I (t U ) and w S (t U ), and are associated with fixed effects γ I and γ S ; c 0i and c 1i are correlated individual random intercept and slope, respectively, with c i ∼ N (0,B); ε ij are the independent Gaussian measurement errors with mean zero and variance σ 2 ε .

Identification of the trajectories of association
The linear mixed model of the health outcome defined in Eq. 3 provides two trajectories of association (i.e., timevarying effects) of the history of past exposure on the subsequent health outcome, one for the initial level (noted γ * I (t U )), one for the slope (noted γ * S (t U )): They represent the mean difference of initial level (i.e. γ * I (t U )) or of change per unit of time (i.e. γ * S (t U )) when the exposure increases of 1-unit at time t U , adjusted for other covariates and when the exposure history is stricly similar at any other time in [ −S, 0].

Specification of weights
In some applications, the form of the weight function w . (t U ) (where subscript . refers either to I or S in Eqs. 3 and (4)) is known. However, most often, it has to be estimated directly from the data. As others previously [22,28,29], we chose to estimate the weight function by regression using a basis of splines [26], which are flexible enough to capture a variety of clinically plausible shapes [30]. Thus, the weight function can be written: where (B k ) k=1,...K refers to the K basis of splines functions and (θ .k ) k=1,...K the coefficients to estimate. For ease of calculation, we denote the intercept Although any type of splines could be considered, we favored natural cubic splines that limit erratic behavior at the boundary of the time window thanks to linearity constraints [31]. In that case, the basis of splines is built from K + 1 knots which are to be chosen inside the time window [ −S, 0]. The number of knots has to be carefully determined from the data. A large number of knots implies high flexibility but it may lead to overfitting. On the contrary, a small number of knots may result in an oversmooth estimate that is prone to under-fit bias [31]. In our work, we considered between one and five inner knots (i.e., K ∈[ 2, 6]) and relied on the Akaike information criterion [32] to select the final splines basis. In addition, in the absence of prior knowledge, we considered equidistant knots.
One additional advantage of approximating the weight function using splines is that the WCIE can be rewritten as a linear combination of K + 1 components: where H ki (for k=0,. . . , K) denote intermediate summary covariates of the exposure history [22,29,33].

Identifiability constraints
The time-varying effects of the exposure history defined in Eq. 4 is overparameterized; coefficients γ I and γ S , and the weights w I (t U ) and w S (t U ) cannot be simultaneously estimated from the data without further constraint.
Hauptmann et al. [28] proposed in a case-control study to estimate the total effect of the history of exposures γ .
S + 1. Following Sylvestre et al. [29] and Danieli et al. [22], we directly estimated the time-varying coefficients without constraining the weights by setting the coefficients γ I and γ S to 1. Therefore, the time-varying effects directly correspond to the weights and Eq. 4 becomes: With this parameterization, the overall mean effect of the history of exposure is 1 Finally, by considering an approximation of the weight functions by splines (Eq. 6) and by not constraining the weights (Eq. 7), the linear mixed model defined in Eq. 3 becomes: where θ Ik and θ Sk are K + 1 unconstrained parameters to be estimated, and H ki are the intermediate covariates defined from Eq. 6.

Maximum likelihood estimation using a two-stage procedure
By sharing the underlying true exposure level U * i (t U ), the sub-model for the exposure in Eq. 1 and the health outcome model in Eq. 3, define a shared parameter joint model, which would call for the simultaneous estimation of all parameters in order to avoid any bias [34]. However, given the peculiar temporal framework of the landmark approach with the inclusion of subjects present in the cohort at the landmark time and the distinct windows of observation for the two variables, a two-stage estimation procedure is unlikely to generate any bias -this assumption is confirmed by simulations in "Simulation study" section.

First stage
In the first stage, we classically compute the maximum likelihood estimators of φ U . Then, based on the estimated fixed effects parameters β and the best linear unbiased predictors of the random effects b i , the individual-specific true exposure can be predicted at any time t U preceding the initial assessment of the outcome:

Second stage
In the second stage, the true underlying exposure level . Therefore, in the case of an approximation by splines (as introduced in "Specification of weights" section), The model for the outcome thus becomes: and parameters φ Y can be again estimated by classical likelihood maximization.
In the end, the estimated time-varying effects of the prelandmark-time exposure on the level of the outcome at the landmark time (i.e. γ * I ) and on the change over time of the outcome after the landmark time (i.e. γ * S ) are:

Standard error estimation
One drawback of two-stage estimation approach is that the variances of the parameters obtained in the second stage do not account for the variability of the parameters estimated in the first stage. To properly take into account this variability of the first stage, we used parametric bootstrap [35]: instead of including in the second stage and computed the corresponding U * m (t U ) to be included in the second stage where φ Ym and V ( φ Ym ) were then obtained. The total variance V tot ( φ Y ) of the estimated parameters φ Y was obtained as the sum of the intra-and inter-replicate variances: The variance of γ .(t U ) was then easily derived using the Delta-method [36].

Implementation
The methodology can be implemented in any standard statistical software. We chose R programming (version 3.6.0) and function hlme of the R package lcmm (version 1.7.8) [37] for estimating the linear mixed models. The codes of the simulation studies are openly available in GitHub at https://github.com/MaudeWagner/WHistory.

Application to body mass index history starting at mid-life and subsequent cognitive trajectories after age 70
To emphasize the utility of our methodology, we investigated in a prospective cohort study the relationship of BMI history collected since mid-life with subsequent cognitive function and cognitive decline in older ages, where prior data indicate changing relations over time, with the possibility of reverse causation at older age [3][4][5].

The Nurses' Health Study
We relied on data from the Nurses' Health study (NHS). NHS began in 1976, when 121,700 female registered nurses, aged 30-55 years and residing in 11 US states, returned a mailed questionnaire about their lifestyle and health, including their weight and height [38]. Thereafter, the participants continued to complete biennial questionnaires, with updated data. BMI was computed as selfreported weight in kilograms divided by height in meters squared; self-reported weight was highly correlated with measured weight in a previous validation study in 184 participants [39]. From 1995 to 2001, all nurses who had reached age 70 or older with no history of stroke were invited to participate in a telephone-based study of cognitive function; the entry in this sub-study corresponds to our landmark time. Among eligible women, 19,415 (92%) completed the initial validated Telephone Interview for Cognitive Status (TICS, score range 0 to 41), a telephonebased modified version of the Mini-Mental State Examination [40]. Cognitive assessments were repeated at three occasions approximately every 2 years, with a high followup rate (>90% among those remaining alive at each followup point).
Following our landmark framework (see Fig. 2), we considered the history of BMI from the entry in the cohort through the first assessment of TICS (at time 0); more than 95% of participants had information on BMI up to 24 years before the initial cognitive interview. For the current analyses, among the 19,415 participants of the cognitive sub-study, we only excluded 34 women who did not have at least one measure of BMI between enrollment and the first cognitive interview or with missing data for educational level (an important potential confounder), leading to a study sample of 19,381 women.
At enrollment, the mean age of women was 50.5 yearsold (SD=2.5), 77.9% had the registered nurse diploma as highest educational level and 34.9% were overweight or obese whereas only 1.3% were underweight. At the landmark time of the initial cognitive interview, mean age was 73.3 years-old (SD=2.3), a majority of women was overweight or obese (53.4%), average TICS was 33.7 points (SD=2.8) and 11.2% of nurses had cognitive impairment. On average, 11.4 (SD=1.7) BMI measurements per person were collected over the 23.7 years (SD=1.2) of the window of exposure, followed by 3.2 (SD=1.1) TICS measurements collected subsequently over 4.9 years (SD=2.5). Figure 3 shows the observed individual trajectories of BMI in the 24 years of the window of exposure preceding the first cognitive interview and those of the TICS in the 8 years of the window of cognitive assessment in a subsample of 150 randomly selected women. In general, BMI tended to increase over time while TICS decreased.

Specification of the statistical models
We considered the following linear mixed effect models for the BMI and the subsequent TICS score: where age0 is the age of women at enrollment into the cognitive sub-study in years; education is the highest degree of diploma (binary variable: Registered nurse versus Bachelor's, Master's or Doctorate); V0 is an indicator for the first cognitive assessment which captures a possible firstpassing effect [41]; F(t Uil ) is a basis of natural cubic splines with 4 inner knots located at 20th, 40th, 60th and 80th percentiles of the elapsed time of exposure (the number of knots was determined by the AIC); and H ki are the intermediate covariate summarizing BMI history: . The weight function of the WCIE was approximated using natural cubic splines with K −1 = 2 inner knots located at 33th and 66th percentiles (the number of knots was determined by AIC). β, α, b i , c i , ε il , and ε ij have been previously defined in "Methods" section. In order to facilitate the interpretation of the results and as introduced in the methods, we assumed linear trajectories of cognitive decline. This choice seemed reasonable since cognition was evaluated in the NHS over a short follow-up. In addition, as expected, the association between BMI exposure history and cognition was non-constant, with shapes of associations on both the initial level and the slope of TICS suggesting opposite remote and recent effects. For the initial level (see Fig. 4, left panel), higher BMI from -24 to -20 years prior to cognitive assessment, corresponding to mid-life, were significantly associated with lower mean TICS at the first cognitive interview assessed after age 70. For example, for similar confounders (i.e., age and education level) and BMI trajectories at any times except at -24 years, a 1kg/m 2 increase of BMI was associated with a lower initial TICS level of 0.025 point on average (95% CI = -0.039; -0.009). Between -20 and -12 years, for similar confounders and BMI trajectories at any times except at the year evaluated, BMI levels were no longer significantly associated with initial TICS level whereas, from -12 to -5 years, higher BMI levels were again associated with worse cognitive function, consistently with that observed earlier in midlife. In contrast, starting around 5 years before cognitive assessment, higher BMI levels became significantly associated with higher initial level of TICS (likely reflecting reverse causation). For the slope of TICS (see Fig. 4, right panel), results showed no significant association with BMI levels from -24 to -13 years preceding the first cognitive interview. However, higher BMI levels between -12 and -5 years and lower BMI levels between -5 years and the first cognitive interview were associated with worse cognitive decline, similarly to what observed with the initial level of TICS (but to a lesser extent).

Results
To help interpret these complex trajectories of association, we examined hypothetical BMI trajectories (e.g., stable BMI of 25 kg/m 2 over the 24 years of exposure) in relation to subsequent changes in TICS (see Fig. 5). As observed with the trajectories of association, these predictions by profile of BMI generally showed that women with greater BMI over time had higher (better) initial mean levels of TICS compared to women who had a stable or a decreasing BMI with age, reflecting as expected the major influence of reverse causation. Overall, this application provides additional evidence that relationships between BMI and cognition are complex and largely depend on careful consideration of the window of exposure when BMI is assessed.
When considering piecewise constant weights instead of weights approximated by natural cubic splines, the trajectories of time-varying effects remained similar (see eFigure 1) suggesting that the splines function well reflected the data. Moreover, results remained generally the same after adjusting the linear mixed model of the outcome for additional potential confounders (i.e., chronic disease history, smoking status or postmenopausal hormone therapy) categorized as binary variables and collected over the same period as BMI history. Likewise, we observed similar trajectories of associations of BMI on cognitive function when we stratified analyses for these factors (results not shown).

Overview of the simulation design
To validate our methodology and its two-stage estimation procedure, we generated cohort data that mimicked the  . 95% confidence intervals were obtained by parametric bootstrap with 500 replicates. A negative estimate indicates that increased BMI is related to worse cognition/more cognitive decline and a positive estimate indicates better cognition/less cognitive decline NHS design and temporal framework of Fig. 2, according to a joint model with the shared true exposure as introduced in "A landmark approach for assessing the dynamics of association between an intermittently evaluated time-varying exposure and a subsequent health outcome" section. In the main simulation study, we generated individual-specific visit times for exposure, using a uniform distribution in [-6, 6] months around theoretical visits every 2 years from -24 years to landmark time 0. At each visit time, we simulated the observed value of the exposure U according to a mixed model as defined in Eq. 13 and considered a 10% proportion of missing visits completely at random (as observed in NHS). The trajectory of U was modelled according to natural cubic splines (with 4 inner knots at 20th, 40th, 60th and 80th percentiles of the observation times) at both the population and individual levels (with correlated individual random intercept and slope which modelled the within-subject correlation over time), and was adjusted for two covariates mimicking age at study entry (normal distribution: mean 51, standard deviation 3) and binary education level (0.25-probability Bernoulli distribution).
We also generated individual-specific visit times for the outcome with a uniform distribution in [-6, 6] months around theoretical visits every 2 years from 0 to 7 years. We then generated the values of the health outcome Y according to the mixed model as defined in Eq. 13 adjusted for age and education (on both the intercept and the slope) and the indicator for the first outcome assessment, and considering a 3% proportion of missing visits completely at random (as observed in NHS). In this model, the true underlying exposure level was considered.

Scenarios
Parameters in these generating models corresponded roughly to those obtained in the application data for BMI and TICS (parameters provided in eTables 1 and 2) with exception for the trajectory of association with BMI history for which we assumed three different relevant scenarios plotted in Fig. 1: • Scenario A: constant negative association over the time window of exposure (this corresponds to a standard CIE) with γ I (t U ) = −0.05 and γ S (t U ) = −0.01 for all t U ∈[ −24, 0]; • Scenario B: negative association far upstream of the outcome assessment and null association at the approach of the initial assessment using a truncated centred normal distribution: , 0]; • Scenario C: trajectories of association obtained in the application between BMI history and TICS trajectory (see Fig. 4).
In the main simulations, we considered a framework close to NHS with frequent exposure data (roughly every two years), a small proportion of missing visits (10% for the exposure and 3% for the outcome), and a small measurement error for the exposure (σ ε = 0.9). However, to further evaluate the methodology in less favorable situations, we also considered cases where (i) the exposure was measured every 4 years, (ii) the proportion of missing data was larger for the exposure (20%), and (iii) the error of measurement was larger (σ ε = 1.8).

Results
We systematically considered 500 replicates of samples of 1,000 participants each. For each scenario, we focused on the estimated trajectories of association with both the initial level and the slope of the repeated outcome, and on the coverage rate of its pointwise 95% confidence interval Fig. 6 Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) or the slope (top right panel) of the outcome of interest across 500 simulations of 1000 subjects each, and corresponding coverage rates of the 95% pointwise confidence interval (lower panels) for Scenario B (distant negative effect) obtained by parametric bootstrap. As shown in Fig. 6 for Scenario B, the two-stage estimation approach retrieved without bias the true generated trajectory of association, and the parametric Bootstrap provided satisfying coverage rate of the 95% confidence interval around the 95% nominal value. The same conclusions were drawn for Scenario A and C with other shapes of trajectories of association (see eFigure 2 and eFigure 3, respectively).
Furthermore, when considering less repeated measurement points for exposure (see eFigures 4, 5, 6), higher rate of missing data (see eFigures 7,8,9), higher measurement error (see eFigures 10,11,12), or higher number of inner knots of cubic splines in the definition of the BMI history (see eFigure 13), parameters were again well estimated with negligible bias and no departure from the expected 95% coverage rate of the 95% confidence interval. Overall the simulation study demonstrated that in this specific temporal landmark framework, the two-stage procedure combined with parametric bootstrap provided a correct inference.

Discussion
We proposed a flexible approach for estimating the trajectory of association between a time-varying exposure of any kind (e.g., related to lifestyle, occupation, environmental toxins), and a subsequent health outcome when exposure and outcome are assessed non-concurrently. We relied for that on the WCIE methodology [8,9] that we incorporated within a landmark approach to handle the non-concurrent windows for the exposure and outcomes. We followed previous works in WCIE and directly estimated the weights assigned to past exposures from the data using flexible natural cubic splines [31]. However, in contrast with most previous WCIE approaches that are restricted to complete error-free exposure, our approach considers that the exposure is endogenous, prone to measurement error, and assessed at discrete and individual-specific times; this naturally handles the intermittent missing visits usually encountered in cohort studies and limits the biases induced by neglecting the measurement error that is most likely present in many epidemiological contexts. Mauff et al. [20] had previously extended WCIE to handle endogenous exposure data. However, their method was limited to concomitant exposure and outcome assessments while in life-long epidemiology, exposures remote from the outcome risk-period are also of great interest. In addition, they had only considered time-to-event outcomes while we wanted to also consider repeated outcomes that are frequently encountered in longitudinal studies. Danieli et al. [22] have just proposed a WCIE methodology for repeated outcome data. However, their method is limited to exogenous complete error-free exposures, and evaluates the cumulative effect of the exposure on the current level of the outcome while in etiological studies such as ours, the interest is usually in the effect of the exposure on the rate of change [42].
Other methods were proposed in the literature to assess the association between two longitudinal processes. By modeling the association through correlated individual random effects, the bivariate mixed-effect model [43] can describe the correlation structure between the exposure and the outcome but usually fails at providing a proper evaluation of the association of the exposure on the outcome. Cross-lagged models (CLM) [44] and causal dynamic models (CDM) [45] circumvent this problem by directly analyzing the effect of the exposure assessed at a current time on the subsequent trajectory of the outcome. While CLM focuses on the outcome level and is intrinsically linked to the discrete visit process, CDM focuses on the change of the outcome in continuous time and is thus to be favored in etiological studies. However, to the best of our knowledge, none of them were extended to acount for the history of the exposure through a WCIE as we proposed, and both methods are limited to concomitant exposure and outcome assessments. Thanks to the separation between the assessment periods of the exposure and the outcome in the landmark approach, we were able to correctly estimate the time-varying effects of the exposure history using a two-stage approach, as demonstrated by simulations. This two-stage procedure makes our methodology easy to implement in standard statistical software, and easy to adapt to other contexts as long as exposure and outcome are not assessed simultaneously. First, any type of outcome could be considered such as time-to-event or binary endpoint as usually encountered in case-control studies. Note however that (i) in casecontrol studies (unlike prospective studies or case-control studies nested within a cohort), exposure data are collected retrospectively from the landmark time, when the groups were selected, which may result in non-differential measurement errors; (ii) with a time-to-event endpoint, our approach would still require that the exposure timewindow ends at the landmark time. Second, non-ignorable dropout for the repeated outcome could be easily taken into account by estimating the outcome parameters in the second stage using a shared random-effect model [46] rather than a standard mixed model. Third, with a repeatedly-measured outcome, any relevant nonlinear shape of trajectory over time could be considered instead of a linear trajectory. Fourth, the approach could be easily adapted to handle other types of exposures, such as categorical variables, by considering a generalized linear mixed model in the first stage. Fifth, the NHS women were highly homogeneous in age and we used the same time scale (i.e., elapsed time since entry in cognitive substudy) to (i) model the underlying true exposure and (ii) define the WCIE; however, in other applications, users could consider other timescales, and even two different timescales, such as age for the exposure trajectory and years before the landmark for the weights. Finally, for ease of epidemiological interpretation we defined the WCIE as the sum of weighted exposures at each time unit over the [ −S, 0] period but the integral over the continuous time window can also be considered, and leads to the same trajectory of association as illustrated in supplementary eFigure 14 with the trajectory of BMI association in the NHS. Simulation studies showed that the proposed twostage estimation model recovers various sets of shapes of the true weight function (i.e., constant over time or time-varying), and provides satisfactory estimates of the strength of the associations. In particular, splitting the estimation procedure into two parts does not seem to lead to bias when the exposure and outcome are nonconcurrent, even in the case of poorer exposure information (i.e., higher missing data or higher time-interval between two visits or higher measurement error).
The method presents however several limitations that warrant consideration. First, while our estimates are capable of capturing the shape of the weight function, the splines functions may not accurately reflect rapid and sudden changes in the relationship [29]. If these effects are expected, the user can improve the model by including one or more knot(s) in regions where a strong curvature is anticipated or by considering a piecewise constant trajectory of association instead of splines. Second, our methodology assumes a linear association between the exposure at each time and the outcome. Extensions to account for nonlinear associations both over time and over the exposure range is a direction for future research. Third, we considered a standard linear mixed model for the exposure which assumes zero-mean homoscedastic Gaussian independent measurement errors. Additional simulations showed that the method remained robust to errors with heavier tails of distribution when these errors were centered around 0 (see eFigure 15) but the method was likely to provide biased estimates in the case of errors not centered around zero (see eFigure 16). In epidemiological contexts where the distribution of the errors is expected to be non-standard (e.g., when an under-reporting is likely), the linear mixed model should thus be revised accordingly. Finally, as in previous WCIE methods, we only considered one time-varying exposure while it could be of interest in the future to simultaneously examine several exposure histories that are known to be interrelated (e.g., BMI and physical activity, or multi-pollutant/chemical mixtures).
Applied to the association between BMI history starting at midlife and cognition at older ages, our methodology confirmed the complex and age-dependent relationship. Indeed, with no a priori assumptions on the shape of relations, and after controlling for BMI history, age and educational level, we found that higher BMI levels at midlife but lower BMI levels at older ages were significantly associated with poorer cognition after the age 70 in women. This bidirectional relationship is consistent with previous studies, including our own work, that showed that obesity at middle age but low BMI late in life were associated with subsequent development of cognitive impairment and dementia [3-5, 47, 48]. Note however that our interpretations are based on a specific population of generally healthy women and that these results cannot be generalized to populations with other socioeconomic or health status. Although a nonlinear J-shaped association (i.e., underweight and overweight are associated with an increased risk of worst cognition) was sometimes reported for BMI in the literature [47], this is unlikely to be the case in our sample in which very few nurses were underweight during the whole exposure period. Finally, we note that the followup rate was high in the NHS and there was little or no attrition so that we did not need to account for death or possible informative dropout in our application.

Conclusions
This methodology offers great flexibility in capturing complex dynamic relationships over the life-course, but also allows application to the majority of epidemiological contexts when the shape of the effects is not known and even more importantly when the history of exposure is measured incompletely and with error. Overall, this method may significantly contribute to broadening the applications of WCIE for a variety epidemiological contexts.
information for exposure (i.e., measured every 4 years instead of every 2 years) across 500 simulations of 1,000 subjects each, and corresponding coverage rates (low panels) for Scenario C (effect mimicking the associations between BMI and TICS in the Nurses' Health Study).
Additional file 9: eFigure 7. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a larger proportion of missing data for the exposure (i.e., 20% instead of 10%) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario A (constant effect).
Additional file 10: eFigure 8. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a larger proportion of missing data for the exposure (i.e., 20% instead of 10%) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario B (distant negative effect).
Additional file 11: eFigure 9. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a larger proportion of missing data for the exposure (i.e., 20% instead of 10%) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario C (effect mimicking the associations between BMI and TICS in the Nurses' Health Study).
Additional file 12: eFigure 10. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a larger error of measurement (i.e., σ E = 1.8 instead of 0.9) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario A (constant effect).
Additional file 13: eFigure 11. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a larger error of measurement (i.e., σ E = 1.8 instead of 0.9) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario B (distant negative effect).
Additional file 14: eFigure 12. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a larger error of measurement (i.e., σ E = 1.8 instead of 0.9) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario C (effect mimicking the associations between BMI and TICS in the Nurses' Health Study).
Additional file 15: eFigure 13. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) and slope (top right panel) of the outcome of interest when considering a higher number of inner knots of cubic splines in the definition of the BMI history (i.e., 3 inner knots located at the 25th, 50th, and 75th percentiles instead of 2 located at the 33th and 66th percentiles) across 500 simulations of 1000 subjects each, and corresponding coverage rates (low panels) for Scenario C (effect mimicking the associations between BMI and TICS in the Nurses' Health Study).
Additional file 16: eFigure 14. Trajectories of associations between the body mass index history calculated every year (in black) or continuously (in blue) in the 24 years prior to the first cognitive interview on the initial level (left panel) or the slope (right panel) of the Telephone Interview for Cognitive Status (TICS) score approximated by natural cubic splines in the Nurses' Health Study (N=19,381), United States . 95% confidence intervals were obtained by parametric bootstrap with 500 replicates.
Additional file 17: eFigure 15. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) or the slope (top right panel) of the outcome of interest across 500 simulations of 1000 subjects each, and corresponding coverage rates of the 95% pointwise confidence interval (lower panels) for Scenario B (distant negative effect) in the context of asymmetric distribution of error measurements (Cauchy distribution with location 0 and scale 0.7).
Additional file 18: eFigure 16. Boxplots of the trajectory of association between the exposure history over the 24 years prior to the initial health outcome assessment on the initial level (top left panel) or the slope (top right panel) of the outcome of interest across 500 simulations of 1000 subjects each, and corresponding coverage rates of the 95% pointwise confidence interval (lower panels) for Scenario B (distant negative effect) in the context of asymmetric distribution of error measurements (Cauchy distribution with location -1 and scale 0.7).