Cognition-Mortality Associations Are More Pronounced When Estimated Jointly in Longitudinal and Time-to-Event Models

With aging populations worldwide, there is growing interest in links between cognitive decline and elevated mortality risk—and, by extension, analytic approaches to further clarify these associations. Toward this end, some researchers have compared cognitive trajectories of survivors vs. decedents while others have examined longitudinal changes in cognition as predictive of mortality risk. A two-stage modeling framework is typically used in this latter approach; however, several recent studies have used joint longitudinal-survival modeling (i.e., estimating longitudinal change in cognition conditionally on mortality risk, and vice versa). Methodological differences inherent to these approaches may influence estimates of cognitive decline and cognition-mortality associations. These effects may vary across cognitive domains insofar as changes in broad fluid and crystallized abilities are differentially sensitive to aging and mortality risk. We compared these analytic approaches as applied to data from a large-sample, repeated-measures study of older adults (N = 5,954; ages 50–87 years at assessment; 4,453 deceased at last census). Cognitive trajectories indicated worse performance in decedents and when estimated jointly with mortality risk, but this was attenuated after adjustment for health-related covariates. Better cognitive performance predicted lower mortality risk, and, importantly, cognition-mortality associations were more pronounced when estimated in joint models. Associations between mortality risk and crystallized abilities only emerged under joint estimation. This may have important implications for cognitive reserve, which posits that knowledge and skills considered well-preserved in later life (i.e., crystallized abilities) may compensate for declines in abilities more prone to neurodegeneration, such as recall memory and problem solving. Joint longitudinal-survival models thus appear to be important (and currently underutilized) for research in cognitive epidemiology.


INTRODUCTION
The boundaries which divide Life from Death are at best shadowy and vague. Who shall say where the one ends, and where the other begins?
-Edgar Allen Poe, The Premature Burial Cognitive abilities (e.g., verbal skill, abstract reasoning) decline at different rates across the adult lifespan and with differential sensitivity to biological and health-related influences (Baltes et al., 2006). Cognitive declines are also closely associated with elevated mortality risk, a relation that persists even after accounting for sociodemographic and health-related variables (Anstey et al., 2006;Aichele et al., 2016). Increased knowledge of cognition-mortality associations can inform strategies to support mental wellness in later adulthood, provide caregivers with insight about the scope and timeframe of end-of-life mental declines, and may be useful for diagnostic purposes: E.g., processing speed declines as indicative of elevated vascular risk (Batterham et al., 2012;Aichele et al., 2016). However, this knowledge hinges on methodological criteria.
To date, favored approaches for studying cognition-mortality associations have included comparing cognitive trajectories of survivors vs. decedents, modeling cognitive change as a function of time-to-death (i.e., terminal cognitive decline in samples of decedents), and the use of two-stage procedures for estimating survival times as a function of differences, and longitudinal changes, in cognitive performance. More recently, several studies have used joint longitudinal-survival models  to estimate, concurrently, intraindividual cognitive change conditionally on interindividual differences in mortality risk, and vice versa. This latter approach may provide increased statistical efficiency and lessen statistical bias (i.e., because longitudinal performance may inform death-related attrition, and mortality-related processes may influence longitudinal performance). In the current study, we compare outcomes from two-stage and joint procedures for modeling longitudinal cognitionsurvival associations. We first highlight considerations for estimating cognitive trajectories in the presence of death-related attrition. We then turn to cognitive ability as predictive of mortality risk (cognitive epidemiology), the focal application for the study.

Cognitive Trajectories, Accounting for Differences in Mortality Risk
Adult lifespan psychology studies have typically sought to describe cognitive changes that are well characterized by chronological age (Hertzog and Nesselroade, 2003;Sliwinski et al., 2006). Outcomes from these studies have prompted researchers to distinguish between broad "crystallized" and "fluid" abilities. Crystallized intelligence (Gc) reflects accumulated knowledge (e.g., education, experience, acculturation) and is generally thought to remain more-or-less stable during adulthood. Fluid intelligence (Gf), which is linked to basic information-processing efficiency, is more sensitive to biological and health-related influences and therefore declines more, and with greater variability, across the lifespan Cattell, 1966, 1967;McArdle et al., 2002;Baltes et al., 2006).
The above characterization of cognitive abilities as defined across functional and temporal dimensions also informs research on cognition-mortality associations. For example, White and Cunningham (1988) proposed that Gc, which is relatively unaffected by aging, is more markedly affected by end-of-life processes (i.e., evincing sharper terminal decline and stronger associations with death) than Gf, which diminishes steadily with age and which presumably has "less room" to decline in the final years of life. However, recent work suggests that elevated mortality risk is more strongly associated with Gf than Gc (Ghisletta et al., 2006;Aichele et al., 2015Aichele et al., , 2016 and that cognition-survival associations are quite robust across the adult lifespan (e.g., Bosworth et al., 1999;Thorvaldsson et al., 2006). However, if it is true that declines in crystallized abilities primarily manifest in close proximity to death, such declines may be especially prone to underestimation due to attrition (i.e., because individuals with sharp declines in Gc are more likely to drop out of a study due to illness and death). This suggests that understanding changes in different cognitive abilities during adulthood depends to some extent on how cognitive performance is modeled in relation to mortality risk.

Cognitive Performance of Survivors vs. Decedents
Most early studies of end-of-life changes in cognition assessed ability on one occasion, with follow-up recording of deaths, but no other events, at a single census point. Performance of survivors was then compared to that of decedents. Outcomes from these studies typically showed that decedents performed worse than survivors (e.g., Bosworth et al., 1999;Small et al., 2003). An important drawback of this methodology is that choice of census date influences group sample composition (i.e., at later census dates, there will be smaller percentages of survivors), and the sole dichotomous nature of survival status defined in this way ignores individual differences in length of survivorship (participants who died shortly before the census date are labeled similarly to those who died much earlier). This means that estimated cognition-mortality effects may either be attenuated or amplified contingent on choice of census date. Obtaining additional mortality-related information (e.g., changes in health status or illness-related dropout) has been shown to partially remedy these shortcomings (Rabbitt et al., 2005). However, an inherent limitation to this methodology is that it cannot accommodate information about individual differences in age at death or censoring due to survivorship.

Estimation of Cognitive Change in the Presence of Attrition
Longitudinal multilevel models (MLM) and structural equation models (SEM) are currently preferred by many researchers for estimating cognitive change. In studies of older adults, these models are usually applied to data with missing observations due to attrition (i.e., death or dropout). Not accounting for this missingness can contribute to underestimation of cognitive declines because higher functioning individuals remain in the study longer, giving an appearance of less decline and/or improved average functioning over time (a selection effect).
If attrition-related missingness is not directly contingent on differences in cognitive functioning (i.e., the outcome variable is not itself a source of attrition), then the data are considered missing at random (MAR), and, under maximum likelihood (ML) estimation, inclusion of covariates informative of death or dropout (e.g., smoking, general health status) can provide unbiased estimates of cognitive change. However, to the extent that missingness is directly contingent on differences in cognitive ability, then the data are considered to be missing not at random (MNAR), and unbiased estimates of longitudinal cognitive change can be obtained only if modeled conditionally on the corresponding basis for non-random attrition, e.g., death (Baraldi and Enders, 2012;Rouanet et al., 2017).
A related and highly debated issue in the application of longitudinal mixed models for estimating changes in health (and by extension cognitive ability) when attrition may be due to death is that such models assume an "immortal cohort" (Dufouil et al., 2004;Kurland et al., 2009;Wen and Seaman, 2018). That is, mixed models estimated under ML impute missing information following death (i.e., for decedents). Thus, it has been argued that the estimated average changes cannot be used for making valid inferences about the broader population (wherein all individuals remain alive) but rather only within the "immortal cohort" population assumed by the model. Attempts to remedy this problem initially took the form of weighting strategies and sensitivity analyses (Dufouil et al., 2004), but more recent innovation may now allow for estimating a conditional profile valid for "mortal" population inference, without separation of decedents' data, by jointly modeling the longitudinal outcome and semi-competing event times of dropout and death (Li and Su, 2018). In either case, person-specific effects partialed from the group average trajectory (i.e., random or intraindividual effects in multilevel models) may still be validly interpreted at the population level and thus useful for inferential and predictive purposes (Rouanet et al., 2017).

Joint Longitudinal and Time-to-Event Models
Joint longitudinal-survival models typically combine a MLM with a proportional hazards survival model, which are estimated concurrently within a single statistical framework (Henderson et al., 2000;Tsiatis and Davidian, 2004;Rizopoulos et al., 2014;Li and Su, 2018). This means that missingness in longitudinal scores is accounted for by differences in mortality-risk (i.e., survival status and survival time), and missing time of death information for survivors is accounted for by differences in longitudinal performance. In other words, joint models are a useful approach when longitudinal information is MNAR conditional on the event of interest (e.g., death). 1 Joint longitudinal-survival modeling has seen increasing use in biomedical research (e.g., Abdi et al., 2013;Thabut et al., 2013;Núñez et al., 2014;Li and Su, 2018), but it has very rarely been applied in research on adult developmental cognitive change (e.g., Ghisletta et al., 2006;Muniz-Terrera et al., 2011. Although the focus of these latter studies was to predict mortality risk contingent on differences and changes in cognitive abilities, it is informative to compare how longitudinal estimates of cognitive performance may differ when estimated in standalone MLM vs. when estimated in joint models using "real world" data, even if the group-level longitudinal parameter estimates (fixed effects) from these models arguably pertain to immortal cohorts.
Toward this end, we could identify only three prior studies that reported estimates of cognitive changes from both standalone longitudinal and joint longitudinal-survival modeling procedures. The first of these studies looked at Alzheimer's disease (AD) onset as the target event in a sample of 398 Swedish older adults. Intercepts and slopes of memory performance differed very little across frameworks, with slightly larger standard errors when estimated with the joint models (McArdle et al., 2005). In the remaining studies (both looking at death as the target event, and both using data from the Swedish Twin Registry; N = 618 and N = 551, respectively), estimates for intercepts of cognitive performance differed little, and slope parameters were non-significant, across frameworks. There was a very slight reduction in the standard error of the intercept in joint analysis from the earlier of these studies (Muniz-Terrera et al., 2011. Taken together, these results indicate at most slight differences in magnitude and accuracy of longitudinal parameter estimates in standalone vs. joint models. However, with so few comparative studies from which to draw conclusions (i.e., three studies, all based on samples of Swedish participants), it remains worth asking (a) whether different results might be found in similar comparative studies of different populations and (b) whether the inclusion of health-related covariates may have mitigated differences across frameworks (i.e., as might be expected under MAR).

Longitudinal Cognitive Performance as Predictive of Survival
Studies in cognitive epidemiology have shown that cognitive performance predicts differences in mortality risk in later adulthood (e.g., Anstey et al., 2001;Whalley and Deary, 2001;Bosworth and Siegler, 2002). Various interpretations have been given for cognition-mortality associations, and these are not necessarily mutually exclusive. For example, higher childhood IQ may influence later socioeconomic status and access to health care (and hence reduced mortality risk). Alternatively, associations between declines in cognitive function and differences in mortality risk may be mediated by an underlying health condition, such as cerebrovascular illness. In modeling cognition-mortality associations, it makes sense to place cognitive performance on the predictor side of the equation and survival status on the outcome side-if for no other reason than temporal precedence (death being final)-whether cognitive performance is viewed as a diagnostic indicator or as playing a causal role.
In previous work using data from the Manchester Longitudinal Study of Cognition (Rabbitt et al., 2004), we showed that lower baseline cognitive performance (intercepts) of Gc and Gf, and sharper decline (negative linear slope) in Gf, predicted increased mortality risk (Aichele et al., 2015). These results were based on a two-stage procedure in which longitudinal cognitive parameters (intercepts, slopes) were estimated prior to their inclusion as predictors in survival analysis. However, if Gc exhibits decline only in the terminal stage (White and Cunningham, 1988), a significant association between decline in Gc and mortality risk would more likely be found if longitudinal and survival processes were jointly estimated (i.e., to account for the effect of mortality-related dropout on estimated decline in Gc). This line of reasoning also applies to Gf to the extent that mortality-related declines in Gf do not overlap declines characterized by aging alone.
Testing these assumptions requires comparing two-stage vs. joint model estimates of cognition-mortality associations. We could identify only two prior longitudinal cognitive studies reporting results both from two-stage and joint longitudinal/time-to-event models. In the first of these studies (McArdle et al., 2005), joint modeling produced more accurate estimates (smaller standard errors) for the effects of cognitive ability on the target event (AD onset). Whereas slope-event estimates were non-significant across frameworks, there was sign flipping in the significant associations between intercepts of cognitive performance and disease onset across standalone (positive association) and joint (negative association) frameworks. In other words, the standalone model showed that better baseline memory predicted higher risk for AD onset, whereas the joint model showed that better baseline memory predicted lower risk for AD onset (a more theoretically admissible outcome). This suggests a selection effect such that AD-related dropout may have biased the estimate (and interpretation) of the memory-AD association in the standalone model (and hence a reason to prefer the joint framework).
The second study included eight cognitive measures (Ghisletta et al., 2006). In the two-stage analyses, none of the cognitive variables significantly predicted mortality risk. In the joint analyses, better baseline ability, and less decline in ability, predicted increased survival time across nearly all cognitive tasks, even when conditioned on age, sex, and socio-economic differences. Although these results highlight the joint modeling approach as potentially critical for accurate estimation of cognition-mortality associations, results from two-stage vs. joint models in this study cannot be directly compared because the former were based on models in which cognitive variables were analyzed simultaneously as predictors, whereas in the latter, the predictive effects of the cognitive variables were estimated independently.
Taken together, outcomes from these two studies indicate that estimated associations between cognitive trajectories and event incidence (whether AD onset or death) may be affected with respect to their accuracy (standard errors) and their magnitude and direction of effect. However, with scant extant research in this area, similar comparative studies are needed before firmer conclusions and recommendations for modeling strategies can be made.

Potential Implications for Public Health
Prevalence rates for cognitive impairment in later life are 15-20%, with more severe forms (i.e., dementias) affecting an estimated 5-7%, of adults age 60 years and over. Cost of care currently approaches $1 trillion (USD) globally (World Health Organization, 2017). Thus, it is critically important to identify cost-effective, non-pharmacological strategies for addressing mental declines in later life. At present, it is thought that well-preserved knowledge and skills (i.e., crystallized abilities) can to some extent compensate for other cognitive declines, such as in memory and abstract reasoning, that are more strongly influenced by mortality-related pathologies (Stern, 2009). A clearer understanding of the extent to which these abilities are differentially sensitive to mortality risk may thus allow for more accurately gauging the potential of crystallized abilities to serve as compensatory factors.

The Present Study
In a large-sample repeated-measures study of middle-aged and older adults, we compared two-stage vs. joint procedures for estimating longitudinal cognition-mortality associations. We estimated these models both with and without adjustment for survival-related covariates (smoking and self-rated health). We expected that differences across statistical frameworks in estimated average cognitive change and in cognition-mortality associations would be most evident in unadjusted models. We expected to observe sharper declines in cognitive abilities, stronger cognition-mortality associations, and improved estimation accuracy (narrower credible intervals) in results from joint analyses compared to results from two-stage analyses.

MATERIALS AND METHODS
Data for the analyses came from the Manchester Longitudinal Study of Cognition (MLSC). The broader aims and methods of the MLSC have previously been described at length (Rabbitt et al., 2004). We therefore briefly summarize participants and measures, whereas we describe statistical methods in greater detail. The MLSC was approved by the Department of Psychology Internal Ethics Committee, University of Manchester; by the University of Manchester Ethics Committee; and by the Greater Manchester National Institute for Health Trust Ethics Committee.

Participants
Participants were recruited by magazine, radio, or television advertisements. A first wave of participants entered the MLSC in 1983, with subsequent cohorts recruited yearly until 1993. Cognitive testing was carried out until 2003. Data were selected from participants who completed one or more cognitive assessments of both Gc and Gf tasks, who were assessed between the ages of 50 and 87 years, who had complete information for mortality status (i.e., dead vs. alive; and if dead, age at death) as of August, 2012, and who were free from dementia at time of assessment.
In the original sample, some individuals were observed prior to age 50 years and/or after age 87 years. Here, we set a minimum cutoff of 20 observations (n = 20) per time point to safeguard against leveraging effects from outliers, resulting in the current age range. Individuals with severe visual or auditory handicaps were excluded from the study. There were 4,197 women (70.5%) and 1,757 men (29.5%), giving a combined N = 5,954. Median age at study entry was 65.0 years [interquartile range (IQR) = 60.0-70.0 years]. In general, participants had higher income and education levels (38% completed one or more years of college) than U.K. averages, but there were substantial numbers of participants with lower incomes and less education than the U.K. population average (Rabbitt et al., 1993). Of participants, 1470 (24.7%) were tobacco smokers. Self-rated health (SRH) scores (scaled 1 = poor to 5 = excellent) showed that participants generally considered themselves to be healthy at study entry (Median = 4, IQR = 3-5).
Mortality information (dates and proximate causes for all deaths between 1983, when the study began, and August 2012, the most recent update) was obtained from a search of death certificates performed by Her Majesty's General Registry

Cognitive Measures
Cognitive data for the current analyses were obtained from measures of crystallized intelligence (Gc) and fluid intelligence (Gf). For each of these domains, three separate tasks were administered shortly after induction into the study and subsequently up to three times, on average at 4-year intervals (IQR = 3.0-5.0-year intervals). Thus, each individual completed up to 4 cognitive testing occasions over a period of approximately 12 years. Cognitive measures were selected on the basis that they were appropriate for assessment of cognitive change in samples of older adults according to lifespan developmental theory (Baltes et al., 2006), are well-known and documented in the empirical literature, have good psychometric properties, and could be administered by pencil-and-paper for effective assessment within a large sample.
Gc was measured by the Raven (1965) Mill Hill Vocabulary A and B (synonyms and word definitions) subtests and by the Wechsler Adult Intelligence Scale-Revised (WAIS-R) vocabulary scale (Wechsler, 1986). Gf was assessed by the Heim (1970) AH4-1 and AH4-2 tasks (logic, arithmetic, number series, verbal and visuospatial objects comparisons) and the Cattell and Cattell (1960) Culture Fair Intelligence Test (2A or 2B, at alternating occasions). Raw scores for each measure were rescaled to be centered at the first testing occasion, with mean = 50 and SD = 10 (T metric). We used longitudinal structural factor analyses with strict factorial invariance to aggregate data from individual tasks as factor scores (i.e., Gc, Gf) across measurement occasions, as described in Aichele et al. (2015). Factor scores were re-centered during estimation with mean = 0 at the first occasion. On average, participants completed 2 cognitive assessments (M = 2.2). Median time in study (across cognitive testing) was 4.0 years (IQR = 0.0-8.0 years). Longitudinal summary statistics based on measurement occasion (observations, age, and cognitive performance) are shown in Table 1.

Trajectories of Cognitive Change
In a first series of analyses, we estimated cognitive changes using multilevel models (MLM; Laird and Ware, 1982). Three MLM were specified, as described by the following latent change parameters: (a) intercept-only model (intercepts as fixed and random effects), (b) linear change model (intercept and linear change component as both fixed and random effects), and (c) quadratic change model (linear change model plus quadratic change component as a fixed effect only because few participants completed more than three assessments). Consistent with previous comparative studies noted in the introduction, we modeled cognitive change based on chronological age in years (centered at age 65 years, which was the median age at testing; for a similar approach, see Gerstorf et al., 2013).
We applied the MLM to data from all participants (i.e., inclusive of both survivors and decedents) and then to data only from decedents (to examine the potential selection effect of excluding survivors). Analyses were conducted separately by biological sex (men, women) because we previously demonstrated non-proportional mortality hazards for women and men (Aichele et al., 2015), and we wanted the standalone MLM to be consistent with the joint longitudinal-survival Age = chronological age in years at time of cognitive assessment. Age ranges at each measurement occasion are not shown as these were capped between 50 and 87 years for the current analyses. Number of deaths were approximated based on age at death (if known) relative to last observed age at assessment. If age at death was = 4 years from last assessment age, then absence at follow-up was attributed to death rather than dropout. Gc, crystallized intelligence, standardized factor score; Gf, fluid intelligence, standardized factor score. Cognitive factor scores were scaled based on the mean and standard deviation (based on combined scores of women and men) at first measurement occasion. Standard deviations for cognitive scores were near 1.0 across occasions and so are not reported. models. MLM models were first estimated using R statistical software (version 3.3.2; R Core Team, 2016) with package nlme (Pinheiro et al., 2016) with maximum likelihood estimation. We compared model fit using standard chi-square tests of changes in log-likelihood. Note that models selected based on these comparisons were used simply for providing starting values required by the joint modeling software, JMBayes, described below (Rizopoulos, 2012(Rizopoulos, , 2016. JMBayes uses Monte Carlo (MCMC) estimation with Gibbs sampling. Therefore, for comparative purposes, we re-analyzed the best-fitting MLM both in standalone and joint longitudinal-survival analytical frameworks using JMbayes MCMC estimation with noninformative priors.
Equation (1) shows the specification for the quadratic change model, under which intercept-only and linear change models were nested: In this equation, Y ti is the cognitive factor score (Gc, Gf) at age t for individual i; B 0 , B 1 , and B 2 correspond to fixed effects of intercept, linear change and quadratic change in chronological age in years; u 0i and u 1i are random effects for intercept and linear change, respectively, and e ti corresponds to residual component (random effects for quadratic change could not be estimated under ML, and we did not include them in models estimated with MCMC). To the best fitting models, we then added covariates (as fixed effects) that were later included in the survival models (see below). These were age at study entry (also centered at 65 years), the interaction of age at study entry with linear change in age, smoking (yes/no), and SRH at study entry. 2 Equation (2) thus adds to Equation (1): B 3 , the effect of age at study inception (StartAge i ); B 4 , the interaction between starting age and age-related linear change; and B 5 and B 6 , the effects of being a smoker and of SRH, respectively.

Cognition-Survival Associations
We estimated cognition-survival associations by two approaches. First, we used a two-stage method in which we estimated individuals' scores for intercepts and linear changes for each cognitive ability from the best-fitting MLM models (i.e., MCMCestimated individual random effects), and we then used those scores to predict mortality risk in a Cox proportional hazards (PH) modeling framework (Cox, 1972). The Cox PH model specifies a multiplicative association between covariates and the target event (in this case death), formulated as follows: where h i (t) is the cumulative hazard up to time t for individual i, h 0 (t) is the baseline hazard function at time t, and w i is a vector of baseline predictors with corresponding regression coefficients γ (described below). Survival models in the twostage procedure were estimated with a MCMC approach using R package spBayesSurv (Zhou et al., 2017) to ensure that both 2-stage and joint modeling were estimated with MCMC. MCMC estimation was run with 20,000 iterations per model with a burnin phase of 3000 iterations. Covariates included in the survival models are described in the following section.

Joint Longitudinal-Survival Models
We next conducted a series of joint longitudinal-survival analyses. Compared to a two-stage approach, joint modeling reduces bias in longitudinal parameters by concurrently accounting for sources of informative dropout (e.g., missingness related to mortality risk) during estimation. It does this by linking the longitudinal and survival statistical frameworks through shared parameterization (Henderson et al., 2000). The defining feature of the joint model is that longitudinal and survival data are modeled simultaneously with respect to a conditional joint density, rather than separately with two independent marginal densities . For fitting the joint models, we used the R package JMbayes (Rizopoulos, 2012(Rizopoulos, , 2016. JMbayes fits joint models under a Bayesian approach using MCMC estimation (with starting values obtained from standalone MLM and CoxPH models). We again specified MCMC estimation to run with 20,000 iterations per model and with a burn-in phase of 3000 iterations). We evaluated model convergence in reference to common MCMC diagnostic plots (trace, auto-correlation, and kernel density). We fit joint models using the JMbayes "shared random effects" association structure, which links the random effects (individual deviations from the trajectory parameters specified in the longitudinal sub-model) to the survival process. For this study, random effects were limited to individual intercept and linear change terms (specified in the longitudinal model of Equation 1), giving the following equation for the initial JMbayes survival sub-model: which expands on Equation (2) by adding α 1 and α 2 , the predictive effects corresponding respectively to the random intercept, u 0i , and random linear change component, u 1i , of the individual cognitive trajectory. Thus, α 1 and α 2 are the shared parameters for the joint longitudinal-survival models (Rizopoulos, 2016). For interested readers, we provide as Supplementary Materials working code (S1) and simulated data (S2) to demonstrate how to implement JMbayes. A longer tutorial is also available (Cekic et al., in press). Further discussion regarding implementation of two-stage vs. joint approaches can be found in Guo and Carlin (2004) and Rizopoulos (2012). For both two-stage and joint longitudinal-survival analyses, we first estimated mortality risk only as a function of the cognitive predictors (intercepts and linear changes for Gc and Gf). These were the shared parameters for the joint longitudinal-survival models, and we wanted to look at their predictive influence in isolation. We then expanded the analyses by further including covariates that from previous work (Aichele et al., 2015(Aichele et al., , 2016 we knew to be comparatively strongly predictive of mortality risk (age at study entry, smoking status, and SRH, thereby expanding Equation 3). Note that we did not include covariates related to socioeconomic status (i.e., education, occupational class) as in prior work we found these variables to have a minor influence on mortality risk in this sample. Survival analyses were run independently by sex because we previously showed non-proportional hazards across men and women. In total, there were 16 survival analyses using data from all participants, with permutations defined as: Statistical framework (2-stage/Joint) * sex(F/M) * cognitive domain (Gc/Gf) * covariates (no/yes).

RESULTS
The MLM estimated with ML all converged. Diagnostic plots for models re-estimated with MCMC indicated that 20,000 iterations were largely sufficient for obtaining stable parameter estimates.

Trajectories of Cognitive Change
Initial model fit comparisons based on changes in log-likelihood are reported as Supplementary Materials (S3). Fit statistics for standalone MLM with quadratic changes (the best-fitting models) estimated using MCMC are provided in Table 2. Explained variation for multilevel models is given as pseudo-R 2 , which was estimated as the proportional reduction in residual variance relative to that from an unconditional means model (i.e., fixed and random intercepts only), consistent with the QC, quadratic change model; QC.c, quadratic change with covariates for age at study entry, starting age * linear change interaction, smoking status, and selfrated health; Obs, observations; BIC, Bayesian information criterion. Pseudo-R 2 was estimated as proportional reduction in residual variance from an unconditional means model (Snijders and Bosker, 1999). method proposed by Snijders and Bosker (1999). Estimated trajectories of cognitive change are shown in Figure 1 (Gc) and Figure 2 (Gf). Estimates for fixed and random effects from the standalone MLM and from the longitudinal components of the joint analyses are provided in Tables 3-6. To facilitate comparison of estimated trajectories across analytic frameworks, we also report linear changes in cognition approximated as changes in standardized scores per decade for models with covariates ( Z/10y; Table 5, last column). In general, explained variation was higher for Gf than for Gc. Gf declined at a faster rate than Gc (the latter showed only slight declines in models adjusted for covariates). Men had higher levels of initial cognitive performance but subsequently showed steeper declines than women. Adjustment for covariates (age at study entry, smoking, and SRH) effectively removed downward curvature (i.e., accelerated declines in later life) that were visibly present in the unadjusted models. We previously reported similar outcomes from analyses of these data (Aichele et al., 2015), so we now turn our attention to the comparative analyses that are the basis for the current study.

Joint Longitudinal-Survival Models vs. Standalone MLM
To compare overall predictive accuracy across standalone vs. joint models, we first examined differences in residual (unexplained) variance (Tables 4, 6) rather than calculating pseudo-R 2 values because the latter are not applicable for the joint models. Note that for the joint models, residual variance here pertains only to the longitudinal sub-model (not the survival sub-model). Specifically, we examined whether the point estimate for residual variance from the joint analysis fell within the 95% CI for residual variance from the corresponding MLM analysis. For Gc, residual variances were slightly higher in joint analyses, but they mostly fell within the 95% CI from standalone MLM. In contrast, for Gf, residual variances were lower in joint models and fell outside the 95% CI from standalone MLM. Thus, for Gf, joint models appeared to confer slightly improved overall estimation accuracy for the longitudinal parameters.
With respect to the average trajectories (Figures 1, 2), in unadjusted models, cognitive changes estimated in joint analyses showed comparatively more downward curvature (significantly steeper quadratic declines) than those from corresponding standalone models. These differences were more pronounced at ends of the age range where observations were sparse (and CI were correspondingly wider due to increased uncertainty). The trajectories estimated in standalone MLM vs. joint models only differed significantly at far ends of the age range in models unadjusted for differences in starting age, smoking, and SRH.
Comparing the fixed longitudinal parameter estimates in models adjusted for covariates (Table 5), estimates for intercepts and for the effects of smoking and SRH on Gc and Gf scores differed dramatically across standalone vs. joint frameworks. However, the average trajectories for Gc and Gf differed litter across frameworks-and in fact less so after adjustment for these covariates (Figures 1, 2). This indicates that explained variance in Gc and Gf scores was partitioned differently across the two frameworks: I.e., estimated effects of health-related predictors were offset by changes in the estimated intercepts in standalone models; whereas in joint models, conditioning on survival status negated any effects of health covariates, and intercepts remained unaffected. With respect to random effects, within-person variation in intercepts of Gc and Gf was reduced, and within-person variation in linear changes in Gc and Gf was amplified, in joint vs. standalone models, whether unadjusted ( Table 4) or adjusted ( Table 6).

Sample Selection: All Participants vs. Decedents
Pseudo-R 2 values were 1%-3% higher for analyses restricted to data from decedents vs. those using data from all participants ( Table 2). Trajectories for decedents indicated worse overall cognitive performance (lower intercepts) compared to trajectories estimated with data from all participants, but the 95% credible intervals (CI) for the average trajectories of each sample largely overlapped (i.e., differences were non-significant), especially at far ends of the observed age range-with the exception of the unadjusted Gc trajectories for women ( Figure 1A) and the unadjusted Gf trajectories for women and men (Figures 2A, C). In other words, adjustment for differences in starting age and health-related covariates reduced differences in average cognitive declines across decedents vs. all participants.

Cognition-Survival Associations
Mortality risk estimates are summarized in Table 7. Significant proportional hazards for Gc predictors emerged only in joint analyses. Better baseline Gc (higher intercepts) predicted reduced mortality risk in unadjusted and adjusted joint models for women but only in adjusted joint models for men. Estimated proportional hazards for Gc intercepts were more negative (stronger effect) and had narrower 95% CI (greater accuracy) when estimated in joint vs. two-stage analyses. Linear changes in Gc significantly negatively predicted mortality risk (i.e., less Gc decline was associated with lower mortality risk) only in unadjusted joint models. Raw estimates for the proportional hazards of linear changes in Gc from the corresponding two-stage models were stronger (more negative) but less accurate (wider 95% CI which overlapped null values). For Gf, better baseline ability (more positive intercepts) and less decline (less negative linear changes) predicted reduced mortality risk in both two-stage and joint analyses, with and without adjustment for health-related covariates. For Gf intercepts, magnitudes and 95% CI of the estimated proportional hazards differed little across two-stage and joint analyses. In contrast, proportional hazards for Gf linear changes were stronger (more negative) and more accurate (narrower 95% CI) when estimated in joint vs. two-stage analyses.
For further comparison, we calculated hazard ratios for 1SD differences in the cognitive predictors. 3 For Gc, linear changes more strongly predicted mortality risk in unadjusted models, whereas Gc intercepts more strongly predicted mortality risk in adjusted models. For Gf, linear changes more strongly predicted mortality risk than did intercepts across analyses. With respect to other covariates, age at study inception was negatively associated with mortality risk-likely a selection effect wherein being older at study inception indicates increased survivability (Lindenberger et al., 2002). Smoking predicted higher mortality risk, and better SRH predicted lower mortality risk. Scaled hazard ratios for non-cognitive predictors differed little across two-stage and joint analyses.

DISCUSSION
We compared two-stage vs. joint longitudinal-survival estimation approaches using data from a large sample study of middle-aged and older participants. Differences in estimated 3 Percent reductions in mortality risk were calculated as (1 − HR Zp ) * 100, where HR Zp is the hazard ratio corresponding to + 1SD difference in predictor p (Zp). Zp was calculated as the square root of the estimated random variance for the given predictor. Thus, for women's intercept for fluid intelligence (in the joint model, with covariates): HR Zp = 1 −exp(−0.012 * √ 52.78) = 1 −0.916, or 8.4% lower risk per + 1SD difference. cognitive declines across standalone MLM vs. joint modeling frameworks were only evident at far ends of the age spectrum, and these differences were negligible after adjustment for differences in smoking and self-rated health (SRH). In contrast, cognition-mortality associations were conspicuously more pronounced when estimated using joint models vs. using a two-stage procedure. Importantly, longitudinal associations between crystallized intelligence and mortality risk emerged only when estimated jointly.

Cognitive Trajectories in Joint vs. Standalone Analyses
Compared to estimates from standalone multilevel models, cognitive trajectories estimated jointly with mortality risk showed sharper declines at ends of the age spectrum where participant data were fewer; however, these differences were negligible after adjustment for age at study entry, smoking, and SRH (Figures 1, 2), variables that in previous analyses we identified as comparatively strongly related to mortality risk in this participant sample. This outcome suggests that for research primarily focused on modeling longitudinal changes in cognitive performance, the use of standalone MLM with inclusion of appropriate health-and survival-related covariates may suffice for common use cases.
We could only identify two prior reports of cognitive declines across both two-stage and joint longitudinal-survival statistical frameworks. In these studies, memory and visuospatial reasoning showed non-significant declines in both frameworks (Muniz-Terrera et al., 2011. This differed from our results wherein performance in Gc and Gf declined (approximately −0.1 SD/decade and −0.6 SD/decade, respectively) in both frameworks. These differences may be because the studies estimated changes in different cognitive abilities or because of different sample selection criteria, including differences in sample size (and by extension, statistical power necessary to estimate reliable change).
We also examined the selection effects of excluding survivors (approximately 25% of participants) on cognitive trajectories estimated in standalone MLM. We found that decedents on average performed worse (lower baseline performance) than participants generally, which is consistent with previous studies (e.g., Bosworth et al., 1999;Small et al., 2003). However, these differences were greatly reduced in models adjusted for differences in starting age, smoking, and selfrated health. This is consistent with prior findings that adjustment for health-related covariates partially remedies the effects of death-related dropout on estimated performance (Rabbitt et al., 2005). Taken together, these results point to the importance of including survival-related covariates in models of adult developmental cognitive change estimated independently of mortality risk.

Cognition-Survival Associations in Joint vs. Two-Stage Methodologies
Gc was significantly associated with mortality risk only in joint analyses, whereas associations between Gf and  mortality risk were robust across both analytical frameworks. Scaled hazard ratios (i.e., risk reduction per + 1SD better cognitive performance) were 10-29% more pronounced when estimated using a joint vs. two-stage procedure. These outcomes were in general the result both of increased predictive strength and improved accuracy (narrower 95% CI) in cognition-mortality associations modeled jointly vs. using a two-stage procedure. These results align with those from the only other study we know of to report cognition-mortality associations for both two-stage and joint procedures, which showed that cognition-mortality associations were significant only in joint models (Ghisletta et al., 2006). However, in this earlier study, cognitive abilities were examined simultaneously in two-stage analyses but independently in joint analyses. Here, cognitive abilities were examined singularly across both methodologies. Thus, to our knowledge, this is the first study allowing for direct comparison of cognition-mortality associations across two-stage and joint longitudinal-survival frameworks.
That Gc was only significantly predictive of mortality risk in joint models 4 may have important substantive implications. White and Cunningham (1988) asserted that, compared to Gf, Gc would be less affected by aging and more markedly affected by end-of-life processes. To date there has been little support for this hypothesis (Bäckman and MacDonald, 2006), but the current results indicate that the longitudinal association between Gc and mortality risk may be underestimated when not modeled jointly (i.e., when mortality-related dropout remains unaccounted for). This is probably due to the increased statistical efficiency (Rizopoulos, 2012(Rizopoulos, , 2016 of joint models, where missingness in longitudinal scores is directly accounted for by observed differences in mortality risk, and missingness in mortality information (right censored time-of-death for survivors) is accounted for by observed differences in longitudinal performance. Additionally, the growing literature on cognitive reserve (Stern, 2009) suggests that well developed crystallized abilities (e.g., verbal intelligence, lifestyle skills, cultural and relational learning) may buffer against the deleterious impacts of neurodegenerative diseases. Just how well-preserved such crystallized abilities are in later life, and by extension their utility as compensatory factors for protecting against declines in abilities more prone to neurodegenerative conditions, may possibly be overestimated when not conditioned on differences 4 We previously showed that intercept for Gc was negatively associated with mortality risk in MLSC participants using a two-stage procedure (Aichele et al., 2015). Mortality information for the current study was updated to a more recent census and thus provided more information on decedents. This may explain why the association between Gc intercept and mortality risk was not observed in the current two-stage procedures. in mortality risk. This point bears consideration for future research in this domain.

Limitations
In joint models, longitudinal parameter estimates were adjusted for death-related dropout. However, study attrition can occur for other reasons: E.g., moving for a new job or retirement, or due to mobility problems unrelated to more serious health issues. We were not able to conduct a more extensive analysis of dropout-related missingness as we lacked information about participants' reasons for leaving. However, prior work has identified value in modeling other sources of informative dropout for estimating cognitive decline (e.g., Sliwinski et al., 2003). We encourage future researchers to consider such variables in addition to mortality risk. Further, we did not consider differences in cognitive declines with respect to specific cause of death, notably related to dementia, because these data were not available at the 2012 MLSC mortality census, which we used here. Dementia-related deaths reported at an earlier 2008 census constituted approximately 3% of all deaths in this sample (well beneath the U.K. population average in persons over 65 years of age).
On a technical level, at the time of analysis, JMBayes required output from standalone multilevel models estimated with maximum likelihood (ML) to obtain starting values for MCMC estimation. This approach encouraged us to test for the best model parameterizations (e.g., intercept-only vs. intercept plus linear slope models) in a standalone MLM framework. In hindsight, it would have been preferable to re-test model parameterization within the joint longitudinalsurvival framework (a) given that we observed some relative improvements in the accuracy of longitudinal parameters under joint estimation, (b) because a simpler subset of longitudinal parameters may have been preferred after conditioning on survival status, and (c) because model comparison under ML differs from that under MCMC. Although our primary substantive focus was in the area of cognitive epidemiology (i.e., cognition as predictive of mortality risk), we also compared effects on cognitive trajectories when estimated in standalone MLM vs. when estimated conditionally on mortality risk. It has been argued that these group-average trajectories pertain to "immortal cohorts" and should be interpreted accordingly. A more recently developed joint modeling approach (Li and Su, 2018) now appears to allow for closed-form estimation of longitudinal changes valid for inferences to mortal (i.e., non-decedent) populations. We look forward to using this approach in future research.

CONCLUSION
Cognition-mortality associations were more prominent when estimated with joint longitudinal-survival models vs. using a two-stage estimation procedure. This may have particular significance for detecting associations between crystallized abilities (e.g., verbal skills) and mortality-related processes. For applications wherein only the longitudinal component (e.g., cognitive change) is of substantive interest, inclusion of survival-related covariates in standalone multilevel models may provide estimates of declines similar to those obtained from a joint modeling approach.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the Department of Psychology Internal Ethics Committee, University of Manchester; the University of Manchester Ethics Committee; the Greater Manchester National Institute for Health Trust Ethics Committee. The patients/participants provided their written informed consent to participate in this study.

AUTHOR CONTRIBUTIONS
SA performed all analyses and was responsible for manuscript preparation. SC provided statistical guidance. PR was responsible for data collection and provided editorial feedback during manuscript preparation. PG assisted with statistical interpretation, provided editorial feedback during manuscript preparation, and was responsible for general oversight. All authors contributed to the article and approved the submitted version.