Quantitative characterization of biological age and frailty based on locomotor activity records

Quantification of human locomotor activity

For this study, we used two large-scale repositories of wearable accelerometer track records made available by the 2003−2006 National Health and Nutrition Examination Survey (NHANES, 12053 subjects, age range 5−85 years old) and the UK Biobank (UKB, 95609 subjects, age range 45−75 years old). For both NHANES and UKB, a continuous, 7-day long activity track was collected for each subject, as well as data for a comprehensive set of clinical variables and death records up to nine years following the activity monitoring. Human physical activity is usually collected in the form of a series of direct sensor readouts, such a 3D accelerations, sampled at a specified frequency of time. However, the NHANES database provides sequences of transformed variables such as the number of steps and the activity counts per minute. Fig. 1A shows plots of two representative 2-day long activity tracks from a middle age (age 43) individual and an older (age 65) individual, who displayed the same level of overall activity. However, their patterns of activity were qualitatively different; the transitions between the different levels of physical activity appeared to be random. Therefore, instead of trying to determine the precise shape of the activity time series, we chose to apply a Markov chain approximation, which is a simple yet powerful probabilistic model from stochastic processes theory (see [9] for a review of its applications, including the stochastic modelling of biological systems).

Figure 1. Quantitative description of human locomotor activity tracks. (A) Individuals with the same daily average level of activity can yet differ by their chronological age, health status and activity distribution during the day. Representative 2-day long locomotor activity tracks of two NHANES 2003−2006 cohort participants aged 43 (upper) and 65 (lower) illustrate how movement patterns can be visually different while having the same level of daily average activity. We quantify individual sample by dividing activity levels into 8 bins (left panel, histograms) and then counting the probabilities W_ij of random jumps from each discrete activity state j to every other state i per unit time (right panel, color corresponds to intensity of transitions with respect to the population average). (B) The eigenfrequencies of the Markov chain transition matrices are calculated for same two middle-aged and old individuals and represented by vertical bars (note the difference in the positions of the bars). The distribution of the eigenfrequencies in the relevant age-cohorts of 35-45 y.o. 65-75 y.o. are illustrated by overlaid transparent histograms (the light green and dark blue, respectively). Power Spectral Densities (PSD) reconstructed for Markov chain transition matrices (see Appendix A for details) reproduces the approximately a scale-invariant segment of the true PSD of the signal on time-scales up to tens of minutes. This characteristic shift of the cross-over frequency with age has been reported in numerous studies of human and animal locomotor activity (see text).

A statistical description of the participants’ activity was based on the concept that any future state of a Markov chain is determined only by its current state and the probabilities of transitioning between different states. We discretized the physical activity measurements over time into eight bins representing activity states (numbered from 1 to 8 and corresponding to increasing activity levels; see histograms to the left of the activity tracks in Fig. 1A). We counted the transitions between every consecutive pair of activity states along the track. For every pair of states i and j, the number of transitions from state j to i was then normalized to the number of times that state j was encountered along the entire activity record. This calculation yielded the kinetic transition rate, i.e. the probability of a stochastic “jump” from state j to state i per unit time. We then combined these transition rates into the transition matrix (TM) elements (shown as bins in heatmaps to the right of the activity tracks in Fig. 1A). The TM elements represent a complete description of the underlying Markov chain model (see Materials and Methods and the Fig. 1A description for additional details).

On a more technical level, the TM element values have the meaning of transition rates and hence can be related to the time scales characterizing the organism’s responses to external perturbations. To make this connection, we checked explicitly that the TM elements satisfied a detailed balance condition [10] and hence the TM eigenvalues represent inverse equilibration times. Using the relation between the autocorrelation function of the time series and the Markov chain TM from Appendix A, we plotted a reconstructed a Power Spectrum Density (PSD) in Fig. 1B for the physical activity signals corresponding to the same two study participants shown in Fig. 1A. Fig. 1B also shows the discrete sets of TM eigenvalues (the TM spectra) for the same individuals. The cross-over frequency on the PSD plots coincides with the lowest TM eigenvalues, corresponding to a time scale in the range of tens of minutes. The time scale corresponding to these eigenvalues is considerably longer than any period associated with body motion and, therefore, should reflect the organism’s physiological state. The observed decrease of the limiting time scale with age (see the density of the TM eigenfrequencies distributions for cohorts of 35−45 y.o. and 65−75 y.o. individuals in Fig. 1B) signifies a reduction of temporal correlations of physical activity in older subjects.

A transition matrix is a conceptually simple and intuitive aggregate characteristics of physical activity time series. TM elements are kinetic transition rates and the spectral properties of TM are directly related to the organism’s responses at physiological time scales. Therefore, TM elements calculated for each sample are a set of useful descriptors characterizing the physiological state of an organism and will be referred to as the physiological variables or, collectively, a physiological state vector or the organism state representation.

Human locomotor activity reveals aging trajectory

To reveal the intrinsic structure of the physical activity data for the entire NHANES study population, we used Principal Components Analysis (PCA), which is commonly used for multivariate data analysis and visualization [11]). PCA is an unsupervised method and can be employed without prior assumptions regarding the functional dependence of the biologically relevant variables on age. Fig. 2A shows the distribution of the transformed data along the first three PCs, PC₁ vs. PC₂ vs. PC₃. Each point represents the average activity profile representing the age-matched cohorts of men and women. The physiological state vector changes in the course of human lifespan, meaning the aging trajectory is continuous, and can be subdivided into distinct phases that are recognizable as the subsequent human development phases. We used one of commonly accepted systems of age classification [12] and divided the trajectory into four segments: (I) childhood and adolescence (younger than 16 years old); (II) early adulthood (16−35 y.o.); (III) middle ages (35−65 y.o.); and (IV) older ages (older than 65 y.o.). Although there was a significant difference between the trajectories of male and female participants, their overall shape and direction were relatively similar.

Figure 2. Principle Component Analysis (PCA) reveals low-dimensional aging trajectory. (A) The graphical representation of the PCA for 5−85-year-old NHANES 2003−2006 participants follows a winding aging trajectory. Samples were plotted in the first three PCs in 3D space along with 2D projections. To simplify the visualization, the PC scores are shown for the age-matched averages for men (squares) and women (diamonds) and color-coded by age. The Roman numerals and corresponding arrows illustrate the approximately linear dynamics of PC scores over sequential stages of human life: I) age<16; II) age 16−35; III) age 35−65; and IV) age >65. (B) Age-dependence of PCA scores along chronological age for NHANES 2003-2006 cohort aged 35+ is shown by age-cohort average values. Human physiological state dynamics has a low intrinsic dimensionality: only the principal component score, PC₁, which corresponds to the largest variance in data, showed a notable correlation with age (Pearson's r = 0.62 for PC₁ and r < 0.2 for other PCs) and therefore could be used as a natural biomarker of age. Shaded regions illustrate the spread corresponding to one standard deviation in each age-matched cohort for PC₁. The inset shows the increase of variance in biological age (PC₁) in the age- and sex-matched cohorts along the chronological age.

According to the Gompertz law, the risk of mortality in human populations increases exponentially in mid-life, starting around age 40 [8,13]. We observed the relevant turning point, a significant shift in character of the physiological state dynamics, between the aging trajectory of segments II and III exactly at this age (Fig. 2A), corresponding to the transition from early adulthood to middle age. In addition, we found another cross-over at approximately 65 years old, corresponding to the boundary between middle age and older age (segments III and IV in Fig. 2A), and occurring in the vicinity of the average human healthspan, defined as the survival free of chronic age-related diseases. According to a recent World Health Organization report [14], the average healthspan is approximately 63 years old. Since aging is the focus of this study, we limited the subsequent analysis to participants older than 40 years old.

In this restricted dataset, aging manifested itself as the approximately linear evolution of the participants' physiological state along the PC₁ direction, which by definition is the direction of the most variance in the data. Only PC₁ scores in this group of participants were strongly associated with chronological age (Pearson’s correlation coefficient r=0.62; see Fig. 2B). Based on our observations, we propose that the first PC score, PC₁, represents a natural definition of biological age, or bioage, which is a quantitative measure of the aging process in the most relevant age range. It increases linearly with chronological age for participants older than 40, as shown in Fig. 2B. In addition, its correlation with age persists (r=0.47) even in the cohort of the most healthy individuals (according to an implementation of the Frailty Index adopted for NHANES in [2]), suggesting that the association cannot only be attributed to the development of illness. The non-linearity in the bioage dynamics with age is weak and is not sufficient to explain the exponential growth of mortality risks with age. The exponential fit of our data yields the doubling rate of approximately 0.02 per year, which is far less than the doubling rate of 0.085 per year according to the Gompertz law. Hereafter, we refer to the biological age dynamics (the dependence of the biological age on the chronological age) and the associated variation of the physiological variables with age as "aging drift".

To characterize the effect of diseases on the dynamics of biological age, we assessed the effect of type 2 diabetes mellitus (T2DM), a common age-related disease, on biological age as defined by PC₁. We compared the mean and standard deviation of biological age in age-matched cohorts of T2DM patients and healthy subjects (Fig. 2B). Generally, the T2DM patients appeared to be older according to their biological age (PC₁) when compared to their chronological age-matched peers. The biological age difference between the healthy and T2DM groups did not significantly change with chronological age.

Biological age acceleration predicts mortality and the incidence of chronic diseases (morbidity)

The biological age acceleration (BAA) is commonly defined as the difference between the biological and the chronological age of an individual. It is a natural measurement of a person’s aging process relative to that of their peers, and can be associated with their lifespan and the presence of chronic age-related diseases (see [15–17]). We propose a more general and robust definition of aging acceleration associated with an arbitrary variable: namely, the residual from the average of the same variable in a cohort of age- and gender-matched individuals (see a recent example in [5]). In this way, aging acceleration can be calculated for any measurement, not simply those expressed in years of life, but also for more sophisticated measures of aging progress, including “biological age” PC₁.

First, we tested the hypothesis that the BAA of “biological age” PC₁ is associated with all-cause mortality. We used the death records available for 4612 NHANES participants aged 40+, including 550 observed death events during the followup of up to 9 years, and obtained the BAA by adjusting the PC₁ score for gender and chronological age. A Cox proportional hazard regression using the BAA value as a co-variate yielded the hazard ratio HR=1.58, 95% CI=1.54−1.62, see Table 1. Notably, we observed a prominent correlation between PC₁ and the average level of physical activity, specifically its log-scaled value (Fig. 3B). As expected, the BAA of the measurement of total activity was also significantly associated with mortality in the same dataset.

Figure 3. Hazards ratio models show high correlation with each other and are strongly associated with average level of physical activity and the largest variance in physiological measurements (PC₁). (A) Scatter plots of estimated mortality hazard ratio (log-scale) vs PC₁ scores and log-hazard ratio estimated by “LogMort” model trained on NHANES survival follow-up data shows high correlation (see text). (B) Different models for hazard ratio of mortality and morbidity show high correlation between each other and the PC₁ "biological age" in NHANES samples. Models for mortality and morbidity were built using Cox proportional hazards method based on either NHANES or UKB death follow-up data and UKB follow-up on diagnosis. All values were adjusted by age and gender and thus represent the corresponding Biological Age Acceleration (BAA) values.

We confirmed the association of “biological age” PC₁ with mortality risks using the independent UKB dataset, which encompasses another 93597 individuals aged 43−78 years (mean of 62.4±7.8 yr) with 285 recorded deaths during a 3 year follow-up. We estimated the "biological age" PC₁ scores for UKB participants using the same Principal Component vector obtained in the NHANES dataset and found a highly significant contribution in the Cox proportional model HR=1.81, 95% CI=1.76−1.86. This suggests that the biological age signature is relatively robust and can be applied to different datasets.

BAA also turned out to be a significant risk factor for the prospective incidence of chronic related diseases. Following [18], we observed that there is a large cluster of age-related diseases, such as cardiovascular diseases (coronary heart disease, heart attack, heart failure, stroke), diabetes, hypertension, and cancer. All the diseases from the list are characterized by an exponentially increasing probability of incidence, with a doubling time similar to that of the mortality rate of eight years identified in the Gompertz mortality law in human population. Assuming the single underlying risk factor, i.e. the aging itself, we defined the healthspan as the age marked by the first diagnosis of any of the aforementioned diseases. To test the association between the BAA and the healthspan, we selected 43533 UKB participants without any age-related disease at the moment of locomotor activity assessment (1331 disease events during 3 year follow-up) and tested the association with the prospective first incidence of chronic disease using the Cox model. The observed hazard ratio (HR=1.16, 95% CI=1.13−1.19) demonstrated a highly significant effect (p=2.4E−7; see Table 1). These data suggest that BAA is a risk factor that marks the increased probability of the prospective incidence of age-related diseases in the healthy individuals.

Supervised proportional hazards models and the biological age

The prevalence of mortality and/or the incidence of major diseases in a population can be inferred from the death records or clinical data and represent the ultimate objective measure of an individual’s resilience to disease. In this section, we introduce and characterize another natural biological age measure–the log-hazard ratio of a risk model–fitted directly to the experimentally observed occurrence of death or the incidence of chronic diseases. We started by confirming that the empirical mortality in the NHANES dataset follows the well-established Gompertz law. To do this, we fit the age-at-death statistics using a parametric Cox-Gompertz model, which is a version of the Cox-proportional hazard model with an explicit Gompertz assumption on the follow-up mortality [19]. We obtained a mortality rate doubling constant of 0.08 per year, which is close to the empirical value of Γ≈0.085 per year. This constant corresponds to a mortality rate doubling time of 8 years [20] and an average life expectancy of 75 years, which is close to 79, the reported value for the United States population (see [21]).

Having established that the expected pattern of exponential mortality exists in the NHANES dataset, we used a Cox proportional hazards model [22], with gender and all of the physiological variables obtained from the locomotor activity as covariates. The mortality risks model was trained on data for NHANES participants aged 40 and older and was then used to estimate the logarithm-scaled risk of mortality for participants from both NHANES and UKB datasets. For simplicity, we will denote these predicted risks as “logMort” to distinguish them from the “biological age” and other models utilized in this study. The risk of death was found to increase exponentially as a function of biological age (Fig. 3A; the determination coefficient of the corresponding log-linear model R²=0.81). This further supports our conjecture that the PC₁ score is a quantitative measure of the biological aging progress. These data suggest that the logarithm of the risk of death may serve as a viable but essentially equivalent approach to evaluate biological age.

The estimated log-hazard ratio robustly predicted the risk of mortality and chronic diseases across the sex- and age-adjusted NHANES and UKB cohorts (see Table 1 for the “logMort” model summary). For every calculation, we used the log-hazard ratio predictor as a covariate and adjusted for chronological age and gender. The resulting BAA estimates were strongly associated with mortality risks in NHANES population (HR=1.76, 95% CI=1.72−1.80). The observations were confirmed in the independent dataset of UKB participants, with a hazard ratio similar to that obtained for the “biological age” PC₁ after adjusting for age and gender (HR=1.81, 95% CI=1.76−1.86). The log-hazard ratio of the mortality model was also significantly predictive of the prospective morbidity (p=1.7E−6), with a similar result observed for the unsupervised “biological age” PC₁ in Table 1.

The apparent similarity between the ability of the supervised (“LogMort”) and the unsupervised (PC₁) biological age models to predict the incidence of death and age-related diseases is a consequence of the high degree of correlation characteristic for biological data and hence the redundancy in the physiological variables measurements. The concordance between the biological age predictors also hints at another practical opportunity to train bioage models. In most cases, the death records are scarce even in very large studies, since the mortality rate is small and requires a very long participant followup to gather sufficient information on the lifespan of the population. The incidence of disease is much higher than the mortality, since the healthspan is shorter than lifespan, and hence any dataset with aging human subjects would contain a considerable fraction of people suffering from age-related diseases. Given this information, we could build a hazards model of the incidence of disease in the UKB dataset, “LogMorb (UKB)”, using the individuals who were considered healthy at the time of the physical activity assessment. As expected, all of the presented mortality and morbidity risks model’s predictors, including the mortality ratio built in UKB data, “LogMort (UKB)”, and the unsupervised biological age PC₁, demonstrated remarkably high correlations across samples after adjustment for age and gender in both NHANES (Fig. 3B) and UKB (data not shown). Our findings indicate a substantial overlap between the signatures of the mortality and morbidity risks that can be used to define human lifespan and healthspan, respectively.

Biological age acceleration and frailty

Finally, we checked the association of BAA with the Frailty Index (FI), a classical measurement of general health using functional and clinical information rather than physiological characteristics. In general, FI is proportional to the overall number of health deficits or diseases and is significantly associated with mortality, see, e.g [1,2]. We had already observed that the BAA computed using biological age is higher in groups diagnosed with T2DM, a disease whose presence increases the FI of the affected individual (see Fig. 2B). We made the same observation using the supervised biological age models (the “LogMort” model), in which the predicted mortality log-hazard ratio was significantly higher in the sex- and age-adjusted cohorts in the NHANES dataset that constituted the "frail" and "most frail" individuals (Fig. 4C), stratified according to a version FI tailored for NHANES dataset [2].

Figure 4. Hazards ratio model distinguished low and high-risk populations and hazardous lifestyles. The effect of unhealthy lifestyle such as smoking caused reversible effect on estimated hazards ratio in the NHANES (A) population and the UK Biobank (B) datasets; (C) Distribution of logarithm of estimated hazards ratio in frailty cohorts shown by median ± standard error of mean (S.E.M.). “Frail” and “most frail” cohorts are stratified on the basis of the respective Frailty Index (FI) values computed according to [2] and are characterized by significant difference in the predicted log-mortality.

The BAA is a continuous measurement with a broader significance beyond predicting the degree of frailty. For the individuals from the healthiest cohorts (i.e. excluding individuals with the "frail" and "most frail" designations), the BAA was identified as a signature of having an elevated risk of chronic diseases incidence (see Table 1) and death associated with a hazardous lifestyle. To demonstrate this association, we compared the biological age of current smokers, those who never smoked, and those who stopped smoking. The calculation revealed a significant difference of BAA between the individuals who currently smoke and those who never smoked (see Fig. 4A). Interestingly, the BAA level of those who stopped smoking was reversed to that of never-smokers. We validated this effect in the UKB population using the model trained only on the NHANES dataset, which produced the same arrangement of BAA differences (see Fig. 4B).

NHANES dataset

Locomotor activity records and questionnaire/laboratory data from the National Health and Nutrition Examination Survey (NHANES) 2003-2004 and 2005-2006 cohorts were downloaded from [www.cdc.gov/nchs/nhanes/index.html]. NHANES provides locomotor activity in the form of 7-day long continuous tracks of “activity counts” sampled at 1 min^-1 frequency and recorded by a physical activity monitor (ActiGraph AM-7164 single-axis piezoelectric accelerometer) worn on the hip. Of 14,631 study participants (7176 in the 2003-2004 cohort and 7455 in the 2005-2006 cohort), we filtered out samples with abnormally low (average activity count <50) or high (>5000) physical activity. We also excluded participants aged 85 and older since the NHANES age data field is top coded at 85 years of age and we desired precise age information for our study. The mortality data for NHANES participants is obtained from the National Center for Health Statistics public resources (4017 in the 2003-2004 cohort and 3985 in the 2005-2006 cohort).

To calculate a statistical descriptor of each participant’s locomotor activity, we first converted activity counts into discrete states with bin edges b_k, k=1..K. Activity level states 1...K−1 were then defined as half-open intervals b_k≤ak+1, state 0 as a1 and state K as aK, where a is the activity count value. In this study, we defined 8 activity states with bin edges b_k=e^k-1, k=1…7. Thus, each sample was converted into a track of activity states and a transition matrix (TM) was then calculated for each participant (see below). To ensure that our analysis dealt only with days on which a participant actually performed some physical activity, we applied an additional filter. We excluded days with less than 200 minutes corresponding to activity states >0. Only participants with 4 or more days that passed this additional filter were retained, yielding a total of 11839 samples (age, years: 35±23, range 6−84; women: 51%). For PCA and Survival analysis, the only samples used were those for participants aged 40 and older with known follow-up on survival/mortality outcome (age, years: 60±13, range 40−84; women: 50%). Once PCA loading vectors were identified, we plotted all NHANES samples’ scores in Fig. 2A, including those for which survival/mortality data were not available.

Transition matrices (TM) T_ij, i=1…8, j=1…8 were calculated as a set of transition rates from each state j to each other state i (the diagonal elements correspond to the probability of remaining in the same activity state). TM elements were calculated as T_ij = N(j→i)/N(j), where N(j) is the number of minutes corresponding to state j and N(j→i) is the number of times the state j was immediately followed by state i (in the consecutive minute along the sample record). We next converted the TM from a discrete point map to continuous notation: W_ij = T_ij -I_ij, where I_ij is the identity matrix. W_ij is the proper TM for which the apparatus of the Markov chains can be used. We used this property to calculate Power Spectral Densities (PSD) and eigenfrequencies (shown in Fig. 1B) based on the assumption that the Markov chain model can be an approximation of observed activity records.

We flattened 8×8 TM of each sample into a 64-dimensional descriptor vector for Principal Component Analysis (PCA) and Survival analysis. Additionally, we converted the flattened descriptor to log-scale to ensure approximately normal distribution for elements of the locomotor descriptor (a useful property for the stability of the linear models that we applied in PCA and Survival analysis). All near-zero elements (<10^-3, which corresponds to less than 10 transitions during a week) were imputed by the value of 10^-3 before log-scaling.

UKB dataset

We accessed data from UK Biobank (UKB) under the approved research projject 21988 (formerly 9086). At the time the present study was conducted (2015-2017), locomotor activity data were available for 103710 UKB participants. Physical activity was measured using Axivity AX3 tri-axial accelerometers worn on the wrist for 7 consecutive days. The data were recorded in the low-level format as continuous tracks of 3D acceleration values sampled at 100Hz. Some tracks indicated that hardware errors occurred during the monitoring period. Participants with more than 10 such hardware errors in their track were excluded from our analysis, leaving 102914 participants. To make it possible to apply the PCA and Survival analysis models established using NHANES data to the UKB data, we downsampled the original UKB tracks to 1 min^-1 (as used in NHANES). For this purpose, individual acceleration records were split into 1-minute slices, and for each slice, the natural logarithm of the sum over the power spectral density (PSD) of the signal within that slice was calculated. Each of these PSDs was calculated from the absolute values of acceleration using the Welch method with 512 points Hann window function and 50% window overlap.

The downsampled UKB tracks represent the level of physical activity per minute but are quantitatively different from the NHANES activity counts. We used a quantile normalization procedure to re-scale the UKB values to the range of discrete activity states of NHANES. We selected NHANES participants in the age range 45-75 and dropped 16 of participants with the lowest and highest average activities. The combined tracks from the remaining 2398 NHANES participants were used to calculate the occupancy fractions p_k = N(k)/N for each NHANES activity state (here N(k) is the number of times the state k was seen and N is the total number of minutes in all tracks). Then we randomly selected 5000 UKB participants from the same age range and similarly dropped 16 of participants with the lowest and highest average activities; this resulted in selection of 3288 UKB participants. Using the combined UKB tracks from selected participants, we found UKB bin edges b′_k such, that the occupancy fractions for the corresponding activity states, were equal to the occupancy fractions in NHANES. Note that such quantile normalization automatically accounts for shift, linear and monotonic non-linear scaling of values, and so the resulting UKB activity states are roughly equivalent to the ones from NHANES. Once bin edges for UKB were obtained, the downsampled UKB tracks were processed exactly as described above for NHANES. TMs and corresponding descriptors were obtained for 95609 UKB participants (age, years: average 61±8, median 62, range 43−78; women: 56%).

Survival analysis

We estimated hazards ratio using Cox proportional hazards model fit to NHANES 2003−2006 linked mortality data. The covariates used in the model included gender label and locomotor activity variables in the form of natural logarithm of transition matrix elements. The total number of covariates was 65 (64 elements of transition matrix and one gender label), so we used regularization parameter λ=0.01. Once fit, the model (“LogMort”) then was applied to produce hazard ratio estimations for NHANES and UKB participants. The model did not explicitly include age of participants. Hazard ratio models for mortality (“LogMort (UKB)”) and morbidity (“LogMorb (UKB)”) were trained in similar way using UKB 3-year follow-up data on death and diagnosis ICD10 codes, respectively. Only UKB participants without any of cardiovascular, diabetes, hypertension, cancer diagnoses at the moment of loocomotor activity measurements were used to train the “LogMorb (UKB)” model (43533 UKB participants with 1331 diseased during 3-year follow-up).

The resulting hazard ratio score was further tested for significance of association with mortality risks again using Cox proportional hazards approach. Now, chronological age and gender were explicitly used as covariates along with the hazard ratio and, optionally, PC₁ score, the latter being an approximation to biological age (see below). Both hazard ratio and PC₁ were linearly detrended by chronological age and gender. This was done to ensure that the obtained significance parameters reflect the contribution of the age- and gender-adjusted part of hazard ratio or PC₁. All procedures were performed in the same way for NHANES and UKB. PC₁ scores for NHANES were obtained using PCA. To obtain PC₁ for UKB participants we calculated projections of UKB variables onto corresponding first eigen vector of NHANES data covariance matrix.

Empirical mortality (i.e. incidence rate depending on age only) was estimated using NHANES death register follow-up data to check consistency with Gompertz law of mortality using parametric Cox-Gompertz proportional hazard model in the form of maximal likelihood optimization adapted from [19] with M₀ and Γ the parameters of Gompertz mortality law, t_i, Δt_i, and δ_i the age, follow-up time and death event indicator of participant i, respectively.

$$logL(M_0,\Gamma )={\displaystyle \sum _{i=1}^N\frac{M_0}{\Gamma }}{e}^{(\Gamma t_i)}(1-e^{(\Gamma \Delta t_i)})+{\displaystyle \sum _{i=1}^N\delta }i(log{M}_0+\Gamma t_i+\Gamma \Delta t_i)$$$\mathrm{ }$

All analyses were conducted using a set of in-house scripts developed in Python [www.python.org] and R [www.r-project.org].