Optimal Versus Realized Trajectories of Physiological Dysregulation in Aging and Their Relation to Sex-Specific Mortality Risk

While longitudinal changes in biomarker levels and their impact on health have been characterized for individual markers, little is known about how overall marker profiles may change during aging and affect mortality risk. We implemented the recently developed measure of physiological dysregulation based on the statistical distance of biomarker profiles in the framework of the stochastic process model of aging, using data on blood pressure, heart rate, cholesterol, glucose, hematocrit, body mass index, and mortality in the Framingham original cohort. This allowed us to evaluate how physiological dysregulation is related to different aging-related characteristics such as decline in stress resistance and adaptive capacity (which typically are not observed in the data and thus can be analyzed only indirectly), and, ultimately, to estimate how such dynamic relationships increase mortality risk with age. We found that physiological dysregulation increases with age; that increased dysregulation is associated with increased mortality, and increasingly so with age; and that, in most but not all cases, there is a decreasing ability to return quickly to baseline physiological state with age. We also revealed substantial sex differences in these processes, with women becoming dysregulated more quickly but with men showing a much greater sensitivity to dysregulation in terms of mortality risk.

Introduction

The aging process involves many physiological changes affecting the homeostatic state of the organism (1). Studies of biomarkers of aging thus have the potential to simultaneously shed light on the underlying aging process and to provide measurement tools for understanding how aging proceeds within an individual. While a large number of studies have examined changes in biomarkers during the aging process, almost all of these have considered markers one at a time and without reference to how optimal levels of the markers change with age [see, e.g., Ref. (2)]. However, it is clear on the one hand that optimal levels can and do change with age (3, 4), and on the other that information contained in biomarkers is highly redundant, non-linear, and complex (5–7). While some studies have succeeded in incorporating one or the other of these complexities, none have yet attempted both simultaneously. The need for new methodological developments that would analytically integrate the number of biomarkers and complex interrelationships among them is recognized in the literature (2).

Recent advances in biodemographic models of aging resulted in development of mathematical models of aging that are based on biological theory, incorporate several essential mechanisms of aging-related changes in organisms, and work with relevant socio-demographic and other information collected in longitudinal data on aging. These models, known as the stochastic process models (SPM) or the quadratic hazard models (4, 8–10), allow us to evaluate “hidden components” of aging-related changes including adaptive capacity, resistance to stresses, physiological norms, and effects of allostatic adaptation – all variables that play important roles in the processes of aging (11–15); their inclusion in the model is therefore important to understand regulatory mechanisms driving observed aging-related changes in physiological variables and their influence on risks of death or developing a chronic disease, as well as for evaluating a genetic component in such processes. The advantage of the SPM-based approaches is that “hidden components of aging” can be estimated from the data although corresponding variables are not directly observed in the study.

Recently, Cohen et al. (5) presented a novel approach for measuring physiological dysregulation via the joint distribution of multiple biomarkers based on calculation of a multivariate distance called the Mahalanobis distance (which we denote D_M) (16, 17). Such a measure represents the deviation of a current physiological state of an organism from the “normal” physiological state. In contrast to previous approaches, such as those used to calculate allostatic load, D_M does not depend on a choice of biomarkers known to be associated with poor health, and thus avoids a potential for circularity (18). The interpretation of D_M as physiological dysregulation has been validated by showing that it (a) increases with age within individuals, (b) predicts mortality controlling for age, (c) gains predictive power as more variables are included, and (d) is not particularly sensitive to which markers are included, nor to their individual correlations with age (5). It has even been shown to function in a wild bird species (19). While individual trajectories of D_M with age have been examined in relation to health outcomes (20), the modeling approach used lacked the advantages of SPMs listed above, notably the ability to quantify adaptive capacity and optimized versus realized trajectories. Differences between men and women have also not been characterized.

The primary goal of this paper is to understand individual trajectories of dysregulation during aging in terms of adaptive capacity, mortality risk, and male–female differences. To do this, we will implement the measure of multivariate distance (D_M) in the SPM framework. Joining the two approaches (SPM, D_M) in the same model can further increase our capacity to uncover important aging-related processes. First, D_M can be used as a quantitative measurement of how an individual departs from the physiological norm and is thus directly informative about parameters to be estimated by SPM. Second, although successfully applied to different data in different settings in one- or two-dimensional cases (3, 4, 21–27), the SPM may face computational difficulties in applications to multidimensional data because of a large number of parameters to be estimated and a computationally intensive likelihood maximization procedure. D_M is useful in application of the SPM because it allows us to work with multiple physiological variables in a one-dimensional model while still allowing us to estimate and interpret all components of this one-dimensional SPM in the same way as in the original model. The SPM originates from the Woodbury–Manton model (28) and complexities associated with multidimensional versions of such models are recognized, see, e.g., Ref. (29), so the current approach is an important step forward in this research area. This paper is the first one addressing the combined effect of multiple physiological variables on mortality in the SPM (8) framework. We carry this out by applying the model to the Framingham Original Cohort (FHS) data. We estimate the model separately for female and male FHS participants to investigate sex differences in the dynamic behavior of the measure of physiological dysregulation and its possible relationships to the observed sex-specific differences in mortality risks.

Data and Methods

Measure of Statistical Distance and Physiological Dysregulation

For a given set of physiological variables (biomarkers) represented by a column vector X measured in an individual at age t, X(t), D_M is defined as (5):

$D_M(X(t))=\sqrt{{(X(t)-\overline{X})}^T{S}^{-1}(X(t)-\overline{X})}$ (1)

where $\overline{X}$ is a vector of means, S is the variance–covariance matrix for respective variables calculated from some population defining the “normal” state (which can be the same population or some other “reference” population), and T denotes transposition.

Mathematical Model of Age Dynamics of Physiological Variables, Aging, and Mortality

Recent advances in mathematical modeling of aging allow us to link the age dynamics of physiological variables and mortality risks and evaluate “hidden components” of aging-related changes incorporating several major concepts of aging in the structure of the model. These models, originating from the Woodbury–Manton’s random walk model (28), became known as the SPM (8), see Section “Introduction.” Theoretical background of the models with survival functions induced by stochastic covariates is presented in Ref. (28, 30, 31). These results are directly applicable to the SPM version by Yashin et al. (8). The random evolution of covariates with age or time (which is a natural assumption, for example, for the evolution of biomarkers affecting survival) is not explicitly modeled in standard approaches, such as logistic or Cox models, and it is captured in such stochastic models. This is a fundamental feature of such models, which is an appealing property for practical applications, but the presence of stochastic covariates requires computations of marginal hazards and survival functions from the conditional hazards. When the conditional hazard is a quadratic function of such random covariates [as in Ref. (8)], such expressions can be computed as shown in Ref. (32, 33). Below, we briefly describe the SPM (8) in which we implement the measure of statistical distance (Eq. 1).

Briefly, the SPM consists of two equations representing the dynamics of chosen variables as a function of age t and conditional mortality rate at age t given the vector of the variables measured at respective ages. Let Y(t) be a K-dimensional stochastic process representing a (K × 1) vector of variables at age t. The conditional hazard of death is specified as:

$\mu (t,Y(t))={\mu }_0(t)+{(Y(t)-f_0(t))}^{T}Q(t)(Y(t)-f_0(t))$ (2)

Here, μ₀(t) is the “residual” or “baseline” hazard that represents the mortality risk, which would remain if the variables in Y(t) followed the “optimal trajectory” represented by a (K × 1) vector function f₀(t). This baseline hazard μ₀(t) models the effect of other risk factors not included in Y(t) that impact the risk of death. The age trajectory of variables f₀(t) that minimizes the risk of death is referred to as the “physiological norm.” The matrix Q(t) is a (K × K) positive-definite symmetric matrix.

The dynamics of the stochastic process Y(t) is given by the following equation:

$dY(t)=a(t)(Y(t)-f_1(t))dt+b(t)dW(t)$ (3)

with initial condition Y(t₀). Here, W(t) is a vector Wiener process with independent components [W(t) is assumed to be independent of the initial vector Y(t₀)], which defines random paths of the variables in Y(t). A matrix of diffusion coefficients b(t) regulates how this “randomness” propagates to variability of trajectories of Y(t).

The vector function f₁(t) (referred to as the “mean allostatic state”) describes the effect of allostatic adaptation (13). This state characterizes the allostatically modified set point for homeostatic regulation, that is, the physiological state that organisms are forced to maintain by the process of homeostatic adaptation, which represents averaged effects of complex interplay among factors controlled by the ontogenetic program, senescence, and long-term stresses exceeding the limits of the organism’s homeostatic regulation. This dynamic behavior [i.e., that Y(t) adapts to changes in the function f₁(t)] in the model is possible due to the presence of the negative feedback mechanism in Eq. 3 with coefficients of homeostatic regulation given by a matrix a(t). Age-related changes in these coefficients characterize changes in adaptive capacity with age. Specifically, the elements of this matrix characterize the rate of the adaptive response for any deviation of the variables from the state f₁(t). This allows us to test hypotheses on decline in adaptive capacity with age [see also more discussion and examples in Ref. (8, 10)]. Note that the general model with age-dependent components in Eqs 2 and 3 can be simplified to the version with age-independent components. This may be helpful in some applications where testing the hypotheses on age-dependence of respective characteristics is not important but in our applications age-dependence of all components [except f₀(t) and b(t)] is essential as it allows for addressing the hypotheses formulated below in Section “Application to Framingham Heart Study Data.”

The SPM described above can be applied to any set of variables for which longitudinal measurements are available in the data. The measure of physiological dysregulation (D_M) represents an example of such a variable: it can be calculated for each individual at each observation point (exam) so that each individual has a trajectory of measurements of D_M at different ages.

Mathematical Model of Age Trajectories of Physiological Dysregulation, Aging, and Mortality

Let Y(t) = D_M[X(t)] as given by Eq. 1. This is now a one-dimensional process represented by Eq. 3 where respective components [W(t), f₁(t), a(t), b(t)] are scalar counterparts of corresponding vectors/matrices. The one-dimensional version of Eq. 2 is:

$\mu (t,Y(t))={\mu }_0(t)+Q(t){(Y(t)-f_0(t))}^2$ (4)

where f₀(t) and Q(t) are also scalars. As D_M measures physiological dysregulation in an organism represented as deviations from the “normal state” (as specified in its calculation), it is natural to assume that f₀(t) = 0 in Eq. 4. That is, in this implementation we assume that the “normal” values of physiological variables from the definition of D_M (Eq. 1) minimize the risk of death at all ages [see also Section “Results and Discussion” regarding the specification of the “normal state” in D_M calculations and alternative implementations of the “norm” f₀(t)]. Note that the definition of D_M involves calculations of a vector of means and the variance–covariance matrix for respective variables from some “reference” population defining the “normal” state. The question of which “reference” population should be used for these purposes can be dictated by data availability and research paradigms in specific studies, but the results are generally not overly sensitive to specifics of the choice (34). Reference populations that differ markedly from study populations across many demographic aspects (e.g., age, race, and sex) perform poorly, but differences in any single aspects have mostly minor effects. In our application, we used the population of individuals from the original Framingham cohort aged 40 years and younger at the baseline exam, as described below in Section “Application to Framingham Heart Study Data.”

Application to Framingham Heart Study Data

The Framingham Original Cohort was launched in 1948 and has continued with biennial examinations to the present. The FHS Original Cohort included 5,209 respondents (nearly all subjects were Caucasians) aged 28–62 years at baseline residing in Framingham, MA, USA, between 1948 and 1951. Data on 5,079 individuals (2,785 females, 2,294 males) from the Original Cohort were available for this study. Data relevant to our analyses include ages at exams, sex, ages at death/censoring (available data contain information on the number of days since the first exam until the event/censoring from which we calculated respective ages dividing the number of days by 365.25), and various physiological variables measured at exams. For our analyses, we selected physiological variables measured at a sufficient number of exams and whose dynamic characteristics have been shown to be related to the risk of death or onset of aging-related diseases in participants of the FHS Original Cohort in our earlier studies including those which used the original SPM approach [see, e.g., Ref. (3, 4, 24, 26, 27, 35)]. The list includes blood glucose (BG), body mass index (BMI), total cholesterol (CH), diastolic blood pressure (DBP), hematocrit (HC), systolic blood pressure (SBP), pulse pressure (PP), and ventricular rate (VR). Their relevance to research on aging has been discussed in the literature, see, e.g., an earlier publication by Manton et al. (36).

The original data on physiological measurements contain missing values (either by design when at some exam the respective variable was not measured or intermittent missing values when some individuals missed an exam). Following Engels and Diehr (37), we imputed intermittent missing longitudinal values for an individual using available data for that individual. Intermittent missing values were imputed using a linear approximation of available observations in the adjacent exams. For missing data at initial exams, two methods were used: the “next observation carried backward” method and the average of measurements at the next two exams. We repeated all computations using the data generated by these two methods and the results were not sensitive to the imputation method so only the results that used the latter method are reported here.

The data on original variables (BG, BMI, CH, DBP, HC, PP, SBP, and VR) were each transformed using the Box–Cox transformation and then standardized to a zero mean and a unit variance so that all transformed variables are on the same scale in the analyses. We used these standardized transformed variables in calculations of D_M, as defined in Eq. 1. We calculated the means and the variance–covariance matrix in Eq. 1 using baseline measurements for individuals aged 40 years and younger (altogether, there were 2,012 such individuals in the Original Cohort, 1,104 females and 908 males). All analyses were performed separately for females and males (to avoid any assumptions on relationships between parameters of the model in two sexes) using the sets of physiological variables included in calculations of D_M variants, as described below.

The values of D_M depend on the specific set of variables included in its definition. Therefore, to test if the results are sensitive to the choice of variables in D_M, it is necessary to perform analyses with different variants of D_M based on different subsets of variables. The most comprehensive approach would be to use all possible combinations of variables as in Cohen et al. (5). However, as the present study applies the SPM that involves intensive computations, the analyses of all such variants are not feasible. Therefore, we selected a “basic” variant of D_M and performed sensitivity analyses for an additional (limited) set of variants to test whether the results are sensitive to the choice of the variables in the construction of D_M. Our “basic” variant is “the most complete” D_M when all variables are included in its construction. Note that, as PP = SBP − DBP, we cannot use these three variables simultaneously in analyses. Therefore, we considered three “basic” variants each including seven variables: BG, BMI, CH, HC, VR, and two variables out of DBP, PP, and SBP. We also computed several additional variants of D_M to perform sensitivity analyses to test: (a) whether the removal of one variable from such “basic” variants substantially changes the results; and (b) whether the removal of HC and VR (which, unlike the other variables, were not measured until exam 4) from such “basic” variants substantially changes the results. Altogether, we analyzed 24 different D_M variants.

In our applications of SPM (Eqs 3 and 4) to these D_M variants we used the following functional forms of the model’s components: (1) a linear function of age for the feedback coefficient: a(t) = a_Y + b_Y(t − t_min), where a_Y < 0, b_Y ≥ 0, and t_min = 28 (which is the minimal age at the first exam in the FHS Original Cohort); (2) a linear function of age for the mean allostatic trajectory: $f_1(t)=a_{f_1}+b_{f_1}(t-t_{\min })$; (3) a linear function of age for the multiplier in the quadratic hazard term [named “vulnerability index” in Arbeev et al. (4) as it characterizes the “robustness,” or “vulnerability,” a component of stress resistance]: Q(t) = a_Q + b_Qt; (4) the Gompertz baseline hazard: $\ln {\mu }_0(t)=\ln a_{{\mu }_0}+b_{{\mu }_0}(t-t_{\min })$; (5) a constant diffusion coefficient: b(t) = σ₁; and (6) a normally distributed initial value $Y(t_0)~N(f_1(t_0),{\sigma }_0^2)$.

This is a parsimonious specification of the model that still allows us to test several hypotheses on the impact of systemic dysregulation on mortality risk and its relation to different aging-related processes represented by components of the model:

(1) H0: Q(t) = 0 [i.e., a_Q = 0 and b_Q = 0, so that there is no quadratic term in the hazard rate (Eq. 4) and mortality is described by the baseline Gompertz rate μ₀(t)]. This is the most important hypothesis because the failure to reject this hypothesis would not allow us to claim that D_M is a legitimate characteristic for capturing physiological dysregulation affecting mortality risk in this sample.

(2) H0: Q(t) = a_Q (i.e., b_Q = 0). Rejection of this hypothesis would mean that the effect of physiological dysregulation on mortality risk is age dependent. This hypothesis allows us to make inferences on aging-related change in stress resistance associated with the variables used in D_M. For example, if Q(t) increases with age then the J-shape of the hazard rate (as a function of D_M at a fixed age t) narrows with age indicating that the same value of D_M induces a larger increase in mortality risk at old ages compared to younger ages, which can be interpreted as the aging-related decline in resistance to stresses associated with the dynamics of the respective variables, see Yashin et al. (8) and Arbeev et al. (4).

(3) H0: f₁(t) = 0 (i.e., $a_{f_1}=0$ and $b_{f_1}=0$). This hypothesis allows us to make inferences on the presence of a systemic dysregulation for the specific set of variables used in D_M. Rejection of this hypothesis would mean that there is a systemic dysregulation in an organism that forces the trajectories of physiological variables to deviate from their “normal values” specified in calculation of D_M [or, equivalently, forces Y(t) to deviate from f₀(t)].

(4) H0: $f_1(t)=a_{f_1}$ (i.e., $b_{f_1}=0$). This hypothesis allows us to make inferences on the changes in the level of systemic dysregulation with age for the specific set of variables used in D_M. Rejection of this hypothesis would mean that the level which the trajectory of Y(t) is forced to follow changes with age. For example, increasing f₁(t) with age would mean an increasing level of systemic physiological dysregulation with age (so that the trajectories of physiological variables deviate further with age from their “norms” defined at younger ages), which can be associated with the manifestation of the aging process.

(5) H0: a(t) = a_Y (i.e., b_Y = 0). This hypothesis allows us to make inferences on the decline in adaptive capacity with age for the specific set of variables used in D_M. Rejection of this hypothesis would mean that the value of the feedback coefficient in Eq. 3 changes with age. As noted, this feedback coefficient in the model is associated with the adaptive capacity of an organism [i.e., the rate of the adaptive response for any deviation of Y(t) from f₁(t)]. In general, the larger the absolute value of a(t), the faster Y(t) tends to f₁(t). This means that if the absolute value of a(t) declines with age then more time is needed for the trajectory of Y(t) to go back to f₁(t) at old ages compared to younger ages, which is the manifestation of the decline in adaptive capacity with age. Note also that Y(t) can deviate from f₁(t) in both directions, i.e., closer to zero as well as toward larger values, and if the absolute value of a(t) declines with age then such periods when the trajectory of D_M departs farther from f₁(t) tend to be longer at old ages. This may mean that aging-related changes observed in the aging human body may also manifest effects of compensatory adaptation which tend to reduce dysregulation effects when f₁(t) becomes large enough. These considerations stimulate further development of methods of dynamic modeling of aging-related changes and their connections with health and survival outcomes. See also further discussion on relationships between adaptive capacity, dysregulation and different components of aging-related changes in Section “Results and Discussion.”

Testing the above hypotheses involves fitting the restricted models (with restrictions on the parameters specified as in the parentheses above) along with the original (unrestricted) model using the same likelihood estimation procedure (8). As we deal with nested models in all cases, the likelihood ratio test can be used to make statistical inference (see Table S1 in Supplementary Material). The likelihood maximization in the SPM is performed using the constrained optimization procedure (implemented, e.g., in the MATLAB’s Optimization toolbox). Constrained optimization is needed because constraints on parameter values are necessary in respective components of the model both for the mortality risk as well as for the risk factor Y(t). The details can be found elsewhere [e.g., Ref. (8, 9)]. Also note that Y(t) can always be transformed (e.g., using the logarithm or the Box–Cox transformation) similar to the original variables X(t), if necessary. The original scale was used in the paper for the sake of interpretability of the resulting trajectories.

Results and Discussion

Estimates of Parameters and Results of Testing Hypotheses

Estimates of parameters of the SPM applied to different variants of D_M calculated for individuals from the Framingham Original Cohort as described in Section “Application to Framingham Heart Study Data” are given in Table 1 (for females) and Table 2 (for males). The tables show that, although the parameter estimates differ for different variants of D_M, they still follow a common pattern so that generally the components of the model look similar in applications to considered D_M variants (a few exceptions are described below). The tables also reveal systematic differences between estimates of parameters in females and males. This is addressed in more detail in Section “Sex-Specific Differences in Estimates in the Model.”

TABLE 1

Table 1. Estimates of parameters of the stochastic process model applied to different variants of D_M calculated for females from the Framingham Original Cohort.

TABLE 2

Table 2. Estimates of parameters of the stochastic process model applied to different variants of D_M calculated for males from the Framingham Original Cohort.

Tables 1 and 2 also present the results of the hypothesis tests specified in Section “Application to Framingham Heart Study Data” (see notes in the tables for definitions of symbols or the absence of those used in different columns of the tables). Most importantly, the tables show that the null hypotheses H0: Q(t) = 0 is rejected for all variants of D_M (p < 0.0001 for all D_M in males and most D_M in females, with maximal p = 0.0007 for D_M with BG, BMI, CH, HC, SBP, and VR). This allows us to conclude that there is a quadratic term in the hazard rate (Eq. 4) so that mortality is not entirely captured by the baseline Gompertz rate μ₀(t), regardless of the variant of D_M. This means that D_M captures effects of deviations of physiological variables on mortality in this application. All considered D_M affect mortality risk in this sample so that non-zero values of D_M result in a higher mortality rate compared to the baseline mortality rates at respective ages. Note also that, although we have non-zero estimates of the multiplier Q(t) and the quadratic term is present in Eq. 4, we also have non-zero estimates of the baseline mortality rate μ₀(t). There are many more risk factors that affect mortality risk and which are not included in definitions of D_M and, correspondingly, in the quadratic part of the hazard. This “residual” mortality rate μ₀(t) summarizes the effects of such unspecified factors affecting mortality risk.

The null hypothesis H0: Q(t) = a_Q is also rejected in all cases, with p ≤ 0.0001 for males and p-values ranging between 0.0005 and 0.01 for females (except D_M with BG, BMI, CH, DBP, and SBP for which p = 0.038). Tables 1 and 2 also show that in all cases the parameter b_Q is positive. These observations indicate that the effect of physiological dysregulation (represented by D_M) on mortality risk is age dependent and that the J-shape of the mortality rate (as a function of D_M at any fixed age) narrows with age. This narrowing pattern with age has the effect that the same level of physiological dysregulation (i.e., the same value of D_M) induces a larger increase in mortality risk at old ages than it does at younger ages (compared to the baseline mortality rate for respective age). This can be interpreted as the aging-related decline in resistance to stresses associated with deviant dynamics of respective physiological variables [see also Yashin et al. (8) and Arbeev et al. (4)].

The null hypotheses H0: f₁(t) = 0 and H0: $a_{f_1}$ are rejected in all cases (all p < 0.0001). This result shows the presence of a systemic dysregulation in an organism that forces the trajectories of physiological variables to deviate from their “normal values” specified in calculation of D_M (equivalently, forces D_M to deviate from zero). The effect of this dysregulation is that, on average, an individual tends to have the values of physiological variables that deviate from the “normal” state (i.e., he/she has non-zero D_M) so that the resulting mortality risk is higher than it could be if that person managed to keep his/her physiological variables equal to the “norm” (i.e., to have zero D_M) in which case the mortality risk would be equal to the baseline mortality at the respective age. Note also that in applications to all variants of D_M parameter $a_{f_1}=0$ is positive. This indicates an increasing level of systemic physiological dysregulation with age that manifests itself in trajectories of physiological variables drifting further away from their “norms” as age progresses.

The results for the null hypothesis H0: a(t) = a_Y are mixed. For many variants of D_M (16 for females and 6 for males), the decline in the absolute value of the feedback coefficient a(t) with age (associated with the decline in adaptive capacity of an organism) is significant with p < 0.0001, and for some other variants, the decline is not significant, with some estimates b_Y = 0. This means that the conclusions about the decline in adaptive capacity depend on the variables used in calculations of D_M and that different physiological variables have different behaviors in terms of the change with age in the ability of an organism to push them back to the trajectories specified by the mean allostatic trajectories. For those variables that indicated the decline in the absolute value of the feedback coefficient a(t) with age, this ability worsens with age and more time is needed for the trajectories of physiological variables to go back to their mean allostatic trajectories at old ages compared to the time needed at younger ages. However, for some physiological variables, this ability of an organism seems not to be affected by age. These results confirm our earlier observations on mixed patterns of age dynamics of adaptive capacity in analyses of different physiological variables and different outcomes (4, 26). A strong negative correlation between intercept and slope for the feedback coefficient a(t) (−0.86 for males, −0.92 for females) indicates that for the D_M variants with worse initial adaptive capacity (i.e., smaller absolute values of a_Y) the subsequent decline tends to be slower than that for the D_M variants with better initial adaptive capacity (i.e., larger absolute values of a_Y).

Sex-Specific Differences in Estimates in the Model

Tables 1 and 2 reveal systematic differences between estimates of parameters in females and males. These differences can be better understood from Figures 1 and 2. Figure 1 compares the estimates of female and male patterns of different components of the model applied to different variants of D_M. It shows that the baseline mortality rate μ₀(t) is higher in males. Males also have higher levels of the multiplier in the quadratic hazard term Q(t), which also increases faster with age than that in females. This means that males have a narrower J-shape of the mortality risk as a function of D_M at each age and, moreover, this J-shape narrows even further at a faster rate in males than the J-shape for females. This confirms our earlier observations (4) and suggests that males have generally lower resistance to stresses associated with deviant dynamics of the respective physiological variables than females: at each age males are more vulnerable to deviations from the “normal” state (i.e., the same value of D_M results in a larger increase in the mortality risk in males compared to females). Increasing patterns of Q(t) in both sexes show that the same value of deviations from the norm (i.e., the same value of D_M) causes a larger increase in mortality rate (compared to the baseline rate at respective age) at old ages than it does at younger ages. This is true for both sexes but a faster increase of Q(t) with age in males implies that this “additional price” which an organism has to pay in terms of increasing mortality rate for deviations of physiological parameters from the “norm” increases faster with age in males than in females. This, along with a higher baseline mortality rate, results in the higher mortality rate for males that is observed in human populations.

FIGURE 1

Figure 1. Estimates of different components of stochastic process model applied to different variants of D_M. (A) baseline mortality rate μ₀(t); (B) multiplier in quadratic hazard term Q(t); (C) mean allostatic trajectory f₁(t); (D) absolute value of feedback coefficient a(t). Different lines correspond to estimates of respective components in specific D_M variants for females (solid lines) and males (dashed lines).

FIGURE 2

Figure 2. Box plots of parameter estimates in different variants of D_M applied to data on females (F) and males (M). Note: points are shown as outliers (“+”) if they are larger than q₃ + w(q₃ − q₁) or smaller than q₁ − w(q₃ − q₁), where q₁ and q₃ are the 25th and 75th percentiles, respectively, and w = 1.5.

Mean allostatic trajectories f₁(t) show the opposite pattern – they are higher and increase faster in females. This indicates that the physiological variables summarized by the respective D_M tend to deviate farther from the “norm” in females and with age this gap widens faster in females. As such deviations generally correspond to abnormal values of physiological variables that are indicators of diseases/conditions (e.g., diabetes or hypertension), this can contribute to the observed higher prevalence of aging-related diseases and conditions in females. These results on sex differences agree with previous studies showing higher mortality in males but greater propensity for clinical frailty in females (38–40). However, the effect of this on mortality is attenuated by the observed lower values of the multiplier in the quadratic hazard term Q(t) in females compared to males. Nevertheless, for some compositions of D_M, we can see no major differences in mean allostatic trajectories for females and males. This observation implies that for some specific indicators of health the male–female mortality–morbidity paradox (41, 42) may not hold and may even be reversed (43).

The results on the age dynamics of the feedback coefficient a(t) are mixed, although females have a tendency to have worse adaptive capacity (both initial adaptive capacity at younger ages and adaptive capacity at old ages). But generally, it confirms our earlier observations that there is no universal behavior of the decline in adaptive capacity for different physiological variables in females and males (4, 26).

Figure 2 displays box plots of parameter estimates in different variants of D_M applied to female and male data. This figure illustrates the sex differences in the components of the model described above from the perspective of model’s parameters. It also shows the estimates of parameters σ₀ and σ₁ indicating that females and males also differ in terms of variability of D_M (both at baseline and dynamic variability over age), and, respectively, in terms of variability of the underlying physiological variables, with females having higher variability than males.

Advantages, Limitations, and Further Perspectives

The advantages of the presented mathematical model of age trajectories of physiological dysregulation, aging, and mortality combine those of the underlying model (SPM) and the measure of physiological dysregulation (D_M), i.e., a robust measure of physiological dysregulation with a biologically explicit and easily interpretable model of age-related changes in physiology in relation to mortality (4, 5, 8, 10, 26, 27). An additional advantage is that the use of this measure of physiological dysregulation (D_M) allows us to perform analysis of a one-dimensional SPM instead of using its multidimensional version.

Physiological variables can be analyzed within the SPM approach by using them explicitly in the hazard or in the summary measures such as the distance D_M. The first approach is more flexible in the sense that we have the capability to specify parameters of the model for each variable separately and, thus, to investigate their effects on the time-to-event outcome in detail, as well as to make inference on their dynamic properties from respective parameters in Eq. 3. However, this flexibility has its price as we deal with a multidimensional model and the computational workload in the likelihood estimation procedure is essential. Therefore, if one is interested in the combined effect of the variables on the outcome of interest (and also is not investigating the dynamic properties of each individual variable) then the approach implementing the summary measures such as D_M becomes beneficial. It substantially reduces the computational burden as it works in a one-dimensional setting while it still has the same components, which can be interpreted in the context of aging-related processes. In the case of a large number of physiological variables this may result in a considerable computational advantage because the calculations of D_M are much faster than the maximization of the likelihood function of the multidimensional SPM for the original set of physiological variables.

From a theoretical standpoint, use of a global measure of physiological dysregulation is likely to be more robust and less noisy than individual biomarkers, which may vary for many reasons unrelated to long-term aging processes. The interpretation of D_M as a robust measure is enhanced by the fact that our results here are highly concordant with previous studies on the Women’s Health and Aging Study (WHAS), the Baltimore Longitudinal Study on Aging (BLSA), and InCHIANTI (5, 20, 44), despite major differences in (a) statistical approach; (b) biomarker choice (only three of the biomarkers here – CH, glucose, and hematocrit – were among the 40+ used in those studies); (c) participant characteristics (e.g., WHAS contained exclusively older women from the Baltimore area; InCHIANTI was in Italy); and (d) follow-up time (~5–10 years for WHAS, BLSA, and InCHIANTI versus ~55 years here). This replication of results confirms the idea that D_M represents an underlying system-level property of physiological dysregulation that can be detected and interpreted relatively robustly under different contexts and with different markers.

The proposed model provides an approach to measure the impact of systemic dysregulation in an organism on mortality risk, which is an alternative to that using a cumulative measure of health deterioration such as an index of cumulative deficits or deficits index (DI) (45, 46), also an efficient approach to investigate aging-related processes of health deterioration. The DI, which takes into account the cumulative contribution of different variables (including those with possibly minor effects of individual variables on the risk) on mortality risk, can also be implemented in the SPM as in our earlier work (22). Remarkably, both approaches revealed similar sex-specific differences in dynamics of model’s components, see Figure 1 here and Figure 1 in Yashin et al. (22): baseline mortality is higher in males; U-shape (J-shape) of mortality as a function of DI (D_M) is narrower and it narrows faster with age in males; and mean allostatic trajectories are higher and increase faster in females (with respective trajectories estimated at zero for males in the case of DI). Such similarities are observed despite, again, major differences in (a) the measure used in the model (D_M and DI); (b) variables used in construction of respective measures (biomarkers here versus questions from the survey questionnaires in the DI studies); (c) participant characteristics, sample size, and follow-up time (the model with DI was applied to the National Long Term Care Survey (NLTCS) data, a nationally representative survey of more than 49,000 Medicare enrollees); and (d) specifications of the SPM [gamma-Gompertz (logistic) baseline hazards, non-symmetric U-shapes, non-zero “norms,” quadratic functions for mean allostatic trajectories, and constant feedback coefficients in the SPM applied to the NLTCS data]. This similarity in the results obtained using these two different approaches can thus indicate manifestation of different aspects of the same general process of aging-related deterioration in an organism that cause sexual dimorphism in the dynamics of aging-related processes.

There are two limitations in the implementation of the model described in this paper. First, we note that, as a single imputation method was used, the uncertainty and SEs could in general be underestimated. Another limitation is that it defines the “normal” state of physiological variables from which deviations are evaluated by D_M using the mean values at the baseline exam in individuals 40 years and younger. That is, it is assumed that these values minimize the risk of death (as a function of these physiological variables) at all ages. Note that, if we insert the formula of D_M (Eq. 1) in Eq. 4 then, essentially, the vector of means $b_{f_1}=0$ in D_M represents the “physiological norm” or “optimum” for respective physiological variables in the corresponding multidimensional SPM for the original variables included in the definition of D_M. Thus, we actually assume here that this “physiological optimum” is age independent. While this assumption can be true for some variables, for other variables it is likely to be incorrect. Many hormone-related variables experience changes in accordance with the ontogenetic program, and these changes are likely to modify the optimal value for these variables [e.g., this could be manifested by the menopause; see also discussion on age-dependent physiological norms in Yashin et al. (3)]. Thus, one possibility to enhance the results presented here is to assume age-dependent “optima” in the definition of D_M. However, this can generally be a challenging problem because one needs to decide which values to use to define the “optima” for different ages before performing analyses of SPM rather than estimate these “physiological optima” from the SPM itself as in our earlier applications (3, 4, 21–24). A possible remedy could be to use a “surrogate” definition of the “physiological optimum” suggested in Yashin et al. (26) taking the average age trajectories of physiological variables for long-lived individuals (say, 90+) on the premise that the long-lived individuals are those who, on the average, managed to keep the age trajectories of physiological variables close to the optimal ones that minimize the risk of death at respective ages. However, there are still computational challenges in this approach as one will need a large enough number of long-lived individuals to reliably estimate the means and the variance–covariance matrices in respective age groups as needed for D_M. Alternatively, we can relax the assumption on zero “optimum” f₀(t) (and this is also a testable hypothesis in this approach). Although in the current implementation of the model a non-zero “optimum” for D_M could mean different things depending on the patterns of changes in the underlying physiological variables, this approach could provide insights on whether the values of biomarkers that are “optimal” for younger individuals are still “optimal” at older ages.

Another direction for applications of the approach presented in this paper is to investigate genetic effects of different candidate genes or single-nucleotide polymorphisms (SNPs) on the measures of physiological dysregulation and on respective aging-related characteristics such as mean allostatic trajectory, resistance to stresses or adaptive capacity which can be evaluated in SPM. For example, one can evaluate carriers of which alleles or genotypes have higher levels of physiological dysregulation, a faster increase in systemic dysregulation with age, a faster decline in resistance to stresses or adaptive capacity calculating respective estimates of SPM components for carriers of different alleles or genotypes. This can be done using either a stratified analysis by genotype or allele using the original SPM (8) as we did here stratifying by sex or implementing D_M in the genetic SPM (9), which has an advantage of using additional information on non-genotyped individuals to increase the power. This can also be done either for individual genes/SNPs or for cumulative “genetic doses,” that is the variable calculating the number of some alleles in an individual’s genome that were pre-selected using some criteria or analyses [e.g., “longevity alleles” as in Yashin et al. (27)].

The desire to have useful biological interpretation of the results of analyses stimulates the need for further development of concepts and models describing the aging-related changes developing in the human body. Although the “dysregulation” variable (D_M) is useful because it allows us to quantify aging-related changes showing how the particular individual differs from the “ideal” standard (e.g., observed in healthy young adults), its applications to studying mechanisms of aging-related changes requires its further development. From the point of view of this paper, the dysregulation is the result of everything that causes the constructed measure (D_M) to deviate from zero where the zero value characterizes the healthy young adults. At the same time, such measure could be too general if one would like to understand causes and mechanisms of aging-related changes, for example, to produce useful recommendations for improving person’s health. One approach might be to examine D_M by specific physiological system (47). Alternatively, in Ukraintseva and Yashin (48) and Arbeev et al. (49), the importance of distinguishing among three components of aging-related changes was emphasized. These include basal (senescent), ontogenetic, and exposure related components. The approach described above has a potential to include these components into the model and use them to address questions about various mechanisms driving dynamics of aging-related changes during the life course. For example, allostatic load can be linked with the exposure-related component. Changes in physiological norm can be linked with the ontogenetic component and changes in stress resistance and adaptive capacity with the basal (senescent) component. All three components contribute to dysregulation in physiological state, although their interaction may sometimes have side effects that might slow down or even reduce D_M. For example, the reduction of adaptive capacity may slow down accumulation of D_M due to weakening negative feedback, especially when f₁(t) continues to increase. This process could partly be due to compensatory adaptation. The effects, however, may not last long because in case of weak feedback regulation the “noise” component in Eq. 3 will make stronger effects on Y(t). Note that different components of aging-related changes can be evaluated from available longitudinal data (24, 27). More detailed analyses of dysregulation components will be performed elsewhere.

Conclusion

Implementation of the measure of physiological dysregulation (Mahalanobis distance) in the framework of the SPM of aging allows one to investigate physiological dysregulation in relation to different hidden mechanisms of aging-related changes, and, ultimately, to yield estimates of how such dynamic relationships can produce an increase in the risk of death with age and how it may be related to the observed sex-specific differences in mortality risks. Results of application of the method to individuals from the Framingham original cohort indicate that physiological dysregulation increases with age; that increased dysregulation results in increased mortality, and increasingly so with age; and that, in most but not all cases, there is a decreasing ability to return to baseline physiological state with age. We also show substantial sex differences in these processes, with women becoming dysregulated more quickly but with men showing a much greater sensitivity to dysregulation in terms of mortality risk.

Author Contributions

Designed the study: KA and AC; performed statistical analyses and data preparation: KA and LA; contributed software for computation of statistical distance: AC and EM; contributed to methodological work with stochastic process model implementing statistical distance: KA, ES, IA, and AY; wrote the paper: KA, AC, and EM; contributed to the final version and interpretation of results: ES, AK, SU, KC, and AY.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Funding

The Framingham Heart Study (FHS) is conducted and supported by the National Heart, Lung, and Blood Institute (NHLBI) in collaboration with Boston University (Contract No. N01-HC-25195). This manuscript was not prepared in collaboration with investigators of the Framingham Heart Study and does not necessarily reflect the opinions or views of the Framingham Heart Study, Boston University, or NHLBI. Research reported in this publication was partly supported by the National Institute on Aging of the National Institutes of Health (NIA/NIH) under Award Numbers R01AG030612, R21AG045245, R01AG046860, P01AG043352, and P30AG034424 (KA, LA, ES, AK, IA, SU, and AY). The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health. This work was also partly supported by the NIA/NIH grant 5U01AG023712-10 “Long Life Family Study (LLFS): Duke University/University of Southern Denmark Field Center” (KA, LA, ES, AK, SU, KC, and AY). This paper was not prepared in collaboration with investigators from other LLFS Field Centers, and thus does not necessarily reflect the opinions or views of the entire LLFS group. AC is a member of the FRQ-S-supported Centre de recherche sur le vieillissement and Centre de recherche du CHUS and is a funded Research Scholar of the FRQ-S. This research was additionally supported by CIHR grant #s 110789, 120305, 119485.

Supplementary Material

The Supplementary Material for this article can be found online at http://journal.frontiersin.org/article/10.3389/fpubh.2016.00003

References

1. Fried LP, Xue Q-L, Cappola AR, Ferrucci L, Chaves P, Varadhan R, et al. Nonlinear multisystem physiological dysregulation associated with frailty in older women: implications for etiology and treatment. J Gerontol A Biol Sci Med Sci (2009) 64:1049–57. doi: 10.1093/gerona/glp076

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Crimmins E, Vasunilashorn S, Kim JK, Alley D. Biomarkers related to aging in human populations. Adv Clin Chem (2008) 46:161–216. doi:10.1016/s0065-2423(08)00405-8

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Yashin AI, Ukraintseva SV, Arbeev KG, Akushevich I, Arbeeva LS, Kulminski AM. Maintaining physiological state for exceptional survival: what is the normal level of blood glucose and does it change with age? Mech Ageing Dev (2009) 130:611–8. doi:10.1016/j.mad.2009.07.004

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Arbeev KG, Ukraintseva SV, Akushevich I, Kulminski AM, Arbeeva LS, Akushevich L, et al. Age trajectories of physiological indices in relation to healthy life course. Mech Ageing Dev (2011) 132:93–102. doi:10.1016/j.mad.2011.01.001

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Cohen AA, Milot E, Yong J, Seplaki CL, Fueloep T, Bandeen-Roche K, et al. A novel statistical approach shows evidence for multi-system physiological dysregulation during aging. Mech Ageing Dev (2013) 134:110–7. doi:10.1016/j.mad.2013.01.004

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Cohen AA. Complex systems dynamics in aging: new evidence, continuing questions. Biogerontology (2015). doi:10.1007/s10522-015-9584-x

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Cohen AA, Milot E, Li Q, Bergeron P, Poirier R, Dusseault-Belanger F, et al. Detection of a novel, integrative aging process suggests complex physiological integration. PLoS One (2015) 10(3):e0116489. doi:10.1371/journal.pone.0116489

CrossRef Full Text | Google Scholar

8. Yashin AI, Arbeev KG, Akushevich I, Kulminski A, Akushevich L, Ukraintseva SV. Stochastic model for analysis of longitudinal data on aging and mortality. Math Biosci (2007) 208:538–51. doi:10.1016/j.mbs.2006.11.006

CrossRef Full Text | Google Scholar

9. Arbeev KG, Akushevich I, Kulminski AM, Arbeeva LS, Akushevich L, Ukraintseva SV, et al. Genetic model for longitudinal studies of aging, health, and longevity and its potential application to incomplete data. J Theor Biol (2009) 258:103–11. doi:10.1016/j.jtbi.2009.01.023

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Yashin AI, Arbeev KG, Akushevich I, Kulminski A, Ukraintseva SV, Stallard E, et al. The quadratic hazard model for analyzing longitudinal data on aging, health, and the life span. Phys Life Rev (2012) 9:177–88. doi:10.1016/j.plrev.2012.05.002

CrossRef Full Text | Google Scholar

11. Troncale JA. The aging process: physiologic changes and pharmacologic implications. Postgrad Med (1996) 99(111–114):120–2.

PubMed Abstract | Google Scholar

12. Hall DM, Xu L, Drake VJ, Oberley LW, Oberley TD, Moseley PL, et al. Aging reduces adaptive capacity and stress protein expression in the liver after heat stress. J Appl Physiol (2000) 89:749–59.

PubMed Abstract | Google Scholar

13. McEwen BS, Wingfield JC. The concept of allostasis in biology and biomedicine. Horm Behav (2003) 43:2–15. doi:10.1016/s0018-506x(02)00024-7

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Rankin MM, Kushner JA. Adaptive beta-cell proliferation is severely restricted with advanced age. Diabetes (2009) 58:1365–72. doi:10.2337/db08-1198

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Robb EL, Page MM, Stuart JA. Mitochondria, cellular stress resistance, somatic cell depletion and lifespan. Curr Aging Sci (2009) 2:12–27. doi:10.2174/1874609810902010012

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Mahalanobis PC. Mahalanobis distance. Proc Natl Inst Sci India (1936) 49:234–56.

Google Scholar

17. De Maesschalck R, Jouan-Rimbaud D, Massart DL. The Mahalanobis distance. Chemometr Intell Lab Syst (2000) 50:1–18. doi:10.1016/s0169-7439(99)00047-7

CrossRef Full Text | Google Scholar

18. Singer BH, Ryff CD, Seeman T. Operationalizing allostatic load. In: Schulkin J, editor. Allostasis, Homeostasis, and the Costs of Physiological Adaptation. Cambridge: Cambridge University Press (2004). p. 113–49.

Google Scholar

19. Milot E, Cohen AA, Vezina F, Buehler DM, Matson KD, Piersma T. A novel integrative method for measuring body condition in ecological studies based on physiological dysregulation. Methods Ecol Evol (2014) 5:146–55. doi:10.1111/2041-210X.12145

CrossRef Full Text | Google Scholar

20. Milot E, Morissette-Thomas V, Li Q, Fried LP, Ferrucci L, Cohen AA. Trajectories of physiological dysregulation predicts mortality and health outcomes in a consistent manner across three populations. Mech Ageing Dev (2014) 141:56–63. doi:10.1016/j.mad.2014.10.001

CrossRef Full Text | Google Scholar

21. Yashin AI, Arbeev KG, Kulminski A, Akushevich I, Akushevich L, Ukraintseva SV. Cumulative index of elderly disorders and its dynamic contribution to mortality and longevity. Rejuvenation Res (2007) 10:75–86. doi:10.1089/rej.2006.0500

CrossRef Full Text | Google Scholar

22. Yashin AI, Arbeev KG, Kulminski A, Akushevich I, Akushevich L, Ukraintseva SV. Health decline, aging and mortality: how are they related? Biogerontology (2007) 8:291–302. doi:10.1007/s10522-006-9073-3

CrossRef Full Text | Google Scholar

23. Yashin AI, Arbeev KG, Kulminski A, Akushevich I, Akushevich L, Ukraintseva SV. What age trajectories of cumulative deficits and medical costs tell us about individual aging and mortality risk: findings from the NLTCS-medicare data. Mech Ageing Dev (2008) 129:191–200. doi:10.1016/j.mad.2007.12.005

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Yashin AI, Arbeev KG, Akushevich I, Ukraintseva SV, Kulminski A, Arbeeva LS, et al. Exceptional survivors have lower age trajectories of blood glucose: lessons from longitudinal data. Biogerontology (2010) 11:257–65. doi:10.1007/s10522-009-9243-1

CrossRef Full Text | Google Scholar

25. Arbeev KG, Ukraintseva SV, Kulminski AM, Akushevich I, Arbeeva LS, Culminskaya IV, et al. Effect of the APOE polymorphism and age trajectories of physiological variables on mortality: application of genetic stochastic process model of aging. Scientifica (Cairo) (2012) 2012:568628. doi:10.6064/2012/568628

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Yashin AI, Arbeev KG, Ukraintseva SV, Akushevich I, Kulminski A. Patterns of aging related changes on the way to 100: an approach to studying aging, mortality, and longevity from longitudinal data. N Am Actuar J (2012) 16:403–33. doi:10.1080/10920277.2012.10597640

CrossRef Full Text | Google Scholar

27. Yashin AI, Arbeev KG, Wu D, Arbeeva LS, Kulminski A, Akushevich I, et al. How lifespan associated genes modulate aging changes: lessons from analysis of longitudinal data. Front Genetics (2013) 4:3. doi:10.3389/fgene.2013.00003

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Woodbury MA, Manton KG. A random-walk model of human mortality and aging. Theor Popul Biol (1977) 11:37–48. doi:10.1016/0040-5809(77)90005-3

CrossRef Full Text | Google Scholar

29. Manton KG, Stallard E, Singer B. Projecting the future size and health status of the United States elderly population. Int J Forecast (1992) 8:433–58. doi:10.1016/0169-2070(92)90057-g

CrossRef Full Text | Google Scholar

30. Myers LE. Survival functions induced by stochastic covariate processes. J Appl Probab (1981) 18:523–9. doi:10.2307/3213300

CrossRef Full Text | Google Scholar

31. Yashin AI, Manton KG, Stallard E. The propagation of uncertainty in human mortality processes operating in stochastic environments. Theor Popul Biol (1989) 35:119–41. doi:10.1016/0040-5809(89)90013-0

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Yashin AI. Dynamics in survival analysis: conditional Gaussian property vs. Cameron-Martin formula. In: Krylov NV, Lipster RS, Novikov AA, editors. Statistics and Control of Stochastic Processes. New York: Springer (1985). p. 446–75.

Google Scholar

33. Yashin AI. An extension of the Cameron-Martin result. J Appl Probab (1993) 30:247–51. doi:10.2307/3214636

CrossRef Full Text | Google Scholar

34. Cohen AA, Li Q, Milot E, Leroux M, Faucher S, Morissette-Thomas V, et al. Statistical distance as a measure of physiological dysregulation is largely robust to variation in its biomarker composition. PLoS One (2015) 10(4):e0122541. doi:10.1371/journal.pone.0122541

CrossRef Full Text | Google Scholar

35. Yashin AI, Arbeev KG, Akushevich I, Arbeeva L, Kravchenko J, Il’yasova D, et al. Dynamic determinants of longevity and exceptional health. Curr Gerontol Geriatr Res (2010) 2010:381637. doi:10.1155/2010/381637

CrossRef Full Text | Google Scholar

36. Manton KG, Stallard E, Woodbury MA, Dowd JE. Time-varying covariates in models of human mortality and aging: multidimensional generalizations of the Gompertz. J Gerontol (1994) 49:B169–90. doi:10.1093/geronj/49.4.B169

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Engels JM, Diehr P. Imputation of missing longitudinal data: a comparison of methods. J Clin Epidemiol (2003) 56:968–76. doi:10.1016/s0895-4356(03)00170-7

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Mitnitski A, Song X, Skoog I, Broe GA, Cox JL, Grunfeld E, et al. Relative fitness and frailty of elderly men and women in developed countries and their relationship with mortality. J Am Geriatr Soc (2005) 53:2184–9. doi:10.1111/j.1532-5415.2005.00506.x

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Puts MT, Lips P, Deeg DJ. Sex differences in the risk of frailty for mortality independent of disability and chronic diseases. J Am Geriatr Soc (2005) 53:40–7. doi:10.1111/j.1532-5415.2005.53008.x

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Clegg A, Young J, Iliffe S, Rikkert MO, Rockwood K. Frailty in elderly people. Lancet (2013) 381:752–62. doi:10.1016/s0140-6736(12)62167-9

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Verbrugge LM. Sex differentials in health. Public Health Rep (1982) 97:417–37.

PubMed Abstract | Google Scholar

42. Verbrugge LM, Wingard DL. Sex differentials in health and mortality. Women Health (1987) 12:103–45. doi:10.1300/J013v12n02_07

CrossRef Full Text | Google Scholar

43. Kulminski AM, Culminskaya IV, Ukraintseva SV, Arbeev KG, Land KC, Yashin AI. Sex-specific health deterioration and mortality: the morbidity-mortality paradox over age and time. Exp Gerontol (2008) 43:1052–7. doi:10.1016/j.exger.2008.09.007

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Cohen AA, Milot E, Li Q, Legault V, Fried LP, Ferrucci L. Cross-population validation of statistical distance as a measure of physiological dysregulation during aging. Exp Gerontol (2014) 57:203–10. doi:10.1016/j.exger.2014.04.016

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Mitnitski AB, Mogilner AJ, Rockwood K. Accumulation of deficits as a proxy measure of aging. ScientificWorldJournal (2001) 1:323–36. doi:10.1100/tsw.2001.58

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Kulminski A, Yashin A, Arbeev K, Akushevich I, Ukraintseva S, Land K, et al. Cumulative index of health disorders as an indicator of aging-associated processes in the elderly: results from analyses of the national long term care survey. Mech Ageing Dev (2007) 128:250–8. doi:10.1016/j.mad.2006.12.004

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Li Q, Wang S, Milot E, Bergeron P, Ferrucci L, Fried LP, et al. Homeostatic dysregulation proceeds in parallel in multiple physiological systems. Aging Cell (2015) 14:1103–12. doi:10.1111/acel.12402

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Ukraintseva SV, Yashin AI. How individual age-associated changes may influence human morbidity and mortality patterns. Mech Ageing Dev (2001) 122:1447–60. doi:10.1016/S0047-6374(01)00277-9

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Arbeev KG, Ukraintseva SV, Arbeeva LS, Yashin AI. Mathematical models for human cancer incidence rates. Demogr Res (2005) 12:237–71. doi:10.4054/DemRes.2005.12.10

CrossRef Full Text | Google Scholar