Timing of mortality during development alters the evolution of aging

TL;DR

The timing of when survival benefits occur during development strongly influences how quickly or slowly aging evolves in animals.

Contribution

A new theoretical framework showing that early-life survival benefits accelerate aging, while late-life benefits slow it down.

Findings

Early-life survival benefits increase reproductive success and accelerate aging.

Late-acting survival benefits reduce the rate of aging evolution.

Mortality timing during development shapes evolutionary forces on aging.

Abstract

For most animals, intrinsic senescence-induced mortality increases with age, while deaths from extrinsic threats, such as predation or accidents, decline during development as individuals grow and mature. Age-dependent modulation of extrinsic mortality is known to influence the evolution of aging. Yet, how the precise timing of changes in mortality shapes evolutionary forces, i.e., when they enhance and when they slow-down aging, remains poorly understood.

To address this, we developed two complementary mathematical models that integrate survival benefits arising during development with the progressive increase in mortality associated with senescence. Agent-based simulations and deterministic analyses revealed a strong and consistent influence of the timing of developmental survival benefits on their evolutionary impact: early-life survival benefits accelerate the rate of ageing, while late-acting benefits slow it down. This difference arises because early-life benefits more strongly enhance individual reproductive success. Early-life benefits thereby increase the overall growth potential of the population more strongly than late-life benefits, thereby diminishing the relative evolutionary advantage of increased longevity.

Our results underscore the importance of mortality timing in the evolution of aging and provide a theoretical framework for connecting developmental trajectories to aging dynamics.

The online version contains supplementary material available at 10.1186/s12862-025-02464-x.

Figures

Early survival benefits accelerate, and late benefits decelerate the aging rate in the DA model. (a) Evolution of the mean aging rate during agent-based simulations. Colors represent populations with timings of survival benefits early (red, *d* = 0) or late (green, *d* = 35) in life, compared to a population without developmentally controlled survival benefits (blue). The filled area on the plot shows 95% confidence interval of the mean. (b) Difference in aging rate *Δb* within the first 1’000 simulation time steps between simulations starting at the same aging rate *b*_*0*_ with and without s

Deterministic calculations confirm the dependence of aging evolution on the developmental timing of survival benefits. Phase diagram distinguishing accelerating or decelerating effects of developmental survival benefits in the DA model calculated by the Euler-Lotka equation. The change in the intrinsic growth rate *Δr* after a 5% change in the aging rate *Δb* was calculated for life histories with and without survival benefits occurring at age *d* for different starting points of the aging rate *b*_*0*_. Blue area marks parameters in which the survival benefit decelerated the evolution of slow

hDA model fits human mortality data. a. Fit of hDA (red), GMM (green) and Siler (purple), Heligman (blue) models to mortality data (black circles) from United Kindgom 1953 b-d. comparison of R^2^, RMSE and AIC calculated for DA, GMM, Siler, and Heligman models fitted to population data 1948–2022. Boxplot show 25th − 75th percentile. Line: median value. dots: outliers (>1.5 times the interquartile range away from the 25th or the 75th percentile). e-g. fitted values of $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{

Early survival benefits accelerate the aging rate, and late survival benefits decelerate the aging rate in the hDA model. a. left: Evolution of aging rate in the simulation of sexually reproducing individuals for two values of $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document}$ _*max*_. Starting aging rate *b*_*0*_ = 0.15 y^− 1^. Fertility function used: Brass polynomial, based on ref [[31](#CR

Fitness landscape distinguishes regions where developmental processes accelerate or constrain the evolution of slower aging. (a) fitness landscape for combination of different aging rates (*b*) and late developmental processes ( $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document}$ ). (b) same as a, but for $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \use

Differential impact on population growth and lifespan gain explains opposite effect of early and late survival benefits for slow aging. a. Example of mortality rates for the GMM for two different aging rates (dark and light red), and for developmental survival benefits (cyan, right panel only) as a function of age. Net mortality is the difference between GMM and survival benefits, shown by vertical arrows. Relative change in mortality upon reduced aging rate is larger in the presence of survival benefits (right). b. Change in lifespan after reduction of aging rate for scenario with (right) and

Full text

Background

Differences in developmental programs produce a wide range of age-specific mortality patterns among species [1]. In many animals, the risk of death from extrinsic causes, such as predation or accidents, declines during development as individuals grow and mature [2–6]. Conversely, in many species [1], the risk of death caused by intrinsic failure increases due to progressive physiological decline as they age [7, 8]. The interaction of mortality modulation through development and through aging is central to understanding the evolution of aging rates [9, 10].

Classic theory by George Williams [11] emphasized a role for extrinsic mortality in shaping the evolution of aging [12, 13]. Williams argued that high extrinsic mortality weakens the selection for slower aging because most individuals die before reaching old age [11]. Building on this verbal argument, subsequent mathematical models showed that this effect only occurs under certain population density-dependences [9, 12, 14–18] or when extrinsic mortality is age-dependent [15, 17, 18]. These models also showed that age-dependent modulation of mortality can, in principle, either strengthen or weaken the selection for slower aging, but the relation between the time-dependence of mortality and the evolutionary consequence on aging remains incompletely understood.

Age-dependence of extrinsic mortality was already implicitly considered by Williams, who distinguished between adult and juvenile mortality [11]. Importantly, however, extrinsic mortality can also vary after reaching adulthood. In humans, for example, mortality from accidents sharply declines in young adults due to behavioral shifts and reduced risk-taking [19, 20]. Other adaptive traits, such as strengthened immunity, can similarly reduce extrinsic mortality with age. These age-dependent modulations of mortality are a key aspect of diverse developmental strategies and are intimately connected to the evolution of aging rates.

Notably, the distinction of extrinsic and intrinsic mortality is not always self-evident [21], as mortality from the same external threat is often shaped by both developmental processes and age-related decline. For example, the ability to evade predators can improve during development due to muscle growth, but deteriorate in old age due to aging-related sarcopenia. In this study, we define changes in extrinsic mortality as age-related modulation of mortality risks that occurs independent of physiological senescence [14]. By focusing on the processes that modulate the mortality risk, rather than on the ultimate cause of death, this definition resolves the above-mentioned ambiguity and allows us to explore how developmental processes influence the evolution of aging rates. Under this definition, extrinsic mortality is indeed rarely age-independent: instead, it commonly declines during development as animals grow, mature, and adjust their behavior to environmental challenges [6, 12, 21–23].

To explore how developmental changes in survival shape the evolution of aging, we developed two complementary mathematical models, which we analyzed by agent-based simulations and by deterministic calculations. Both models showed that the developmental timing of survival benefits is critical for their evolutionary impact: while early survival benefits accelerate the aging rate, late survival gains slow it down. This pattern arises because early and late survival benefits differentially affect population growth, thereby altering the strength of selection for slow aging. Our results underscore the evolutionary significance of time-resolved, developmentally driven changes in mortality and offer a powerful framework for understanding the co-evolution of development and aging.

Results

Early survival benefits accelerate, and late benefits decelerate the aging rate in a minimal model

To investigate the interactions between aging and developmentally occurring survival benefits, we created a parsimonious model, called the Development-Aging (DA) model. The model assumes that mortality due to aging increases exponentially with time according to the Gompertz law [24].

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${h}_{G}\:({t})={a}^*{e}^{bt}$$\end{document}

where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:a$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:b$$\end{document} are constants, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:t$$\end{document} is time (age), and hG (t) is the hazard (mortality) rate. The model further assumes that the deviations from the Gompertz law are due to developmentally occurring survival benefits that modulate this mortality. For this minimal model, we describe survival benefits \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\beta\:$$\end{document} as a step function that decreases mortality hG from the Gompertz model by 20% from age \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:d$$\end{document} onwards, yielding an overall mortality.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:h\left(t\right)={h}_{G}\left(1-\beta\:\right),\:\text{w}\text{i}\text{t}\text{h}\:\beta\:\left(t\right)=\left\{\begin{array}{c}0\:if\:t<d\\\:0.2\:if\:t\ge\:d\end{array}\right.$$\end{document}

We call the parameter \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:b$$\end{document} the aging rate, i.e., the exponential rate of increase in the hazard rate if there were no developmentally occurring survival benefits.

To ask how the timing of developmental survival benefits affects the evolution of aging, we used agent-based simulations, where individual agents underwent repeated rounds of death, reproduction, and mutation. Simulations were carried out with a population of 10’000 individuals. The mortality rate of each individual was determined by the DA model, and each individual was given a 2% chance of mutation to alter the rate of aging (i.e., parameter b). Reproduction was sexual with constant fertility of one offspring per time step between ages 15 to 100 and zero fertility outside of this interval. Population density-dependence was implemented by filling freed spaces by random selection of offsprings created in the previous reproduction step. In this minimal model, we did not consider tradeoffs between aging and other life-history traits.

As expected, the aging rate rapidly declined during simulations (Fig. 1a) and eventually plateaued when the benefit of any additional slow-down of aging was counteracted by genetic drift (Supplemental Fig. 1). To ask how developmentally controlled survival benefits impacted the evolution of aging, we systematically varied the age \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:d$$\end{document} at which the reduction in mortality occurred. We then determined how this timing affected the rate of decline in the aging rate during simulations. For small \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:d$$\end{document} , i.e., an early onset, survival benefits accelerated aging rate. However, with increasing \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:d$$\end{document} this effect became weaker and eventually shifted to slow down the aging rate (Fig. 1a, b). The effect size on the rate of evolution scaled with the magnitude of the survival benefit (parameter β), but the difference between early and late timing was qualitatively unaltered by this parameter (Fig. 1c, d).

Fig. 1. Early survival benefits accelerate, and late benefits decelerate the aging rate in the DA model. (a) Evolution of the mean aging rate during agent-based simulations. Colors represent populations with timings of survival benefits early (red, d = 0) or late (green, d = 35) in life, compared to a population without developmentally controlled survival benefits (blue). The filled area on the plot shows 95% confidence interval of the mean. (b) Difference in aging rate Δb within the first 1’000 simulation time steps between simulations starting at the same aging rate b0 with and without survival benefit (20%) occurring at age d. (Δb = (b0 – b1000, no benefit) - (b0 – b1000, benefit at age d) = b1000, benefit at age d – b1000, no benefit)). Positive values (blue) indicate deceleration of the evolution of slow aging by the survival benefit (i.e., faster aging). Red indicates accelerated evolution (i.e., slower aging). Each value represents the mean of 50 simulations of 10’000 individuals over 1’000 time points. (c) same as b, but with a fixed value of d = 5 and different strengths of survival benefits. (d) same as c but with d = 40

Deterministic calculations confirm the conclusions from agent-based simulations

To validate the results from agent-based simulations, we computed the selective fitness as the intrinsic rate of population increase r by the Euler-Lotka Equation [25] for a wide range of b0 and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:d$$\end{document} of the DA model. These calculations closely matched results from agent-based simulations (Fig. 2), showing that the effect of developmental survival benefits on aging is not restricted to specific implementations used in the simulations, such as population density-dependence or population size.

Fig. 2. Deterministic calculations confirm the dependence of aging evolution on the developmental timing of survival benefits. Phase diagram distinguishing accelerating or decelerating effects of developmental survival benefits in the DA model calculated by the Euler-Lotka equation. The change in the intrinsic growth rate Δr after a 5% change in the aging rate Δb was calculated for life histories with and without survival benefits occurring at age d for different starting points of the aging rate b0. Blue area marks parameters in which the survival benefit decelerated the evolution of slower aging (Δrbenefit at age d > Δr no benefit). Red area: (Δrbenefit at age d < Δr no benefit). White area: 1^− 13^ >(Δr benefit at age d - Δr no benefit) >-1^− 13^

We conclude that, in the DA model, the timing of developmental survival benefits qualitatively changes their impact on the evolution of aging. The age d where the effect of survival benefits changed sign usually occurred around the age of first reproduction but also depended on the initial aging rate b0. This additional dependence suggests a more complex relationship between development and aging than a simple distinction between juvenile and adult mortality [26].

An extended DA model describes deviations of human mortality from the Gompertz model by early and late survival benefits

To assess if these findings also apply to more complex mortality dynamics, we extended the DA model to capture realistic human mortality dynamics (human DA model, hDA). We chose human mortality dynamics for this analysis because of their complex shape, and for the unrivaled quality and availability of human mortality data. To capture the complex age-dependent mortality in humans, the hDA model uses the Gompertz-Makeham model (GMM), hGMM (t) = c + ae*^bt^, and subtracts developmental survival benefits with either exponentially decelerating or logistic age dependence:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{h}_{hDA}\:\left(t\right)\:=\:c\:+\:a*{e}^{b\:t}-\:\gamma\:\left(t\right)\:-\:\lambda\:\left(t\right)$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}\text{where}\, \gamma\, \left(t\right)&={\gamma}_{max}\left(1-\frac{1}{1+{t}^{h}}\right)\, \text{and}\cr \lambda\,\left(t\right)&={\lambda}_{max}\left(1-\frac{1}{1+{e}^{n\left(t-k\right)}}\right)\end{aligned}$$\end{document}

The terms \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\gamma\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} are scaling factors, corresponding to the maximal magnitude of the respective survival benefits. As for the DA model, we call parameter b of the Gompertz term the aging rate.

The hDA model provided good fits to human mortality data across juvenile and adult age classes for the population data from eight countries in the period from 1948 to 2020. Quality of fit was similar to other multi-parameter models of human mortality, such as Heligman’s model [27] (Fig. 3a-d). Specifically, the hDA model captured the rapid decline in mortality at a young age, and a plateau in the hazard rate during early maturity, as well as the exponential increase in mortality at ages above 30 years (Fig. 3a). As expected, the hDA model provided better fits than the GMM, which is not designed to capture all-age mortality, and the Siler model [28], which does not capture the mortality plateau in early adulthood (Fig. 3a).

The half-way points of λ(t), determined by fits, ranged between 28 and 80, depending on year and country (mean = 36.56, sd = 7.22). The halfway point of γ(t) is inherently fixed to 1 year by the model. Thus, γ(t) represents early survival benefits, and λ(t) represents those occurring later in development. Although the hDA model is agnostic to the exact biological nature of the processes providing the benefit, these processes may involve developmental growth and developmental maturation for γ(t), and adaptive immunity, behavioral changes, or learning for λ(t).

Fitted parameters of the hDA model followed consistent trends across countries and time (Fig. 3e-g, Supplemental Fig. 2) that matched intuitive expectations. For example, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max continuously declined between 1948 and 2020 (Fig. 3e, f). This trend is expected, given the reduction of childhood mortality by advances in health care [29], which reduces the survival benefit of having passed early childhood. Also, parameters that did not follow a monotonous trend, such as parameter k (Fig. 3g), showed high consistency across countries. Finally, fitted parameters matched well with commonly used measures of population health: \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max was highly correlated to childhood mortality under 5 years (Fig. 3h), Parameter c was correlated with mortality from chronic disease (Fig. 3i), and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max also closely followed life expectancy at birth (Fig. 3j).

In summary, the hDA model describes human mortality by an increasing hazard rate according to the classic GMM, corrected by survival benefits starting either early in life or during adulthood. This extended model fits human mortality data well and offers interpretable parameters.

Fig. 3hDA model fits human mortality data. a. Fit of hDA (red), GMM (green) and Siler (purple), Heligman (blue) models to mortality data (black circles) from United Kindgom 1953 b-d. comparison of R^2^, RMSE and AIC calculated for DA, GMM, Siler, and Heligman models fitted to population data 1948–2022. Boxplot show 25th − 75th percentile. Line: median value. dots: outliers (>1.5 times the interquartile range away from the 25th or the 75th percentile). e-g. fitted values of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\gamma\:_{max}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\lambda \:_{max}$$\end{document} , and k of hDA model over time for indicated countries. h. correlation of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\gamma\:_{max}$$\end{document} and child mortality [30] for 1972–2020. i. correlation of c with mortality from cardiovascular disease (CVD), cancer, diabetes or chronic respiratory disease (CRD) [30] for 2000–2019. j. correlation of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\gamma\:_{max}$$\end{document} and life expectancy from 1948–2015. e-j. color represents countries, as indicated at the bottom of the figure

Early survival benefits decelerate, and late benefits accelerate the evolution of slower aging in the hDA model

To investigate how the magnitude of early and late survival benefits affected the evolution of aging, we ran agent-based simulations of the hDA model with initial parameters obtained from fits to human mortality data. We used the same simulation design as for the minimal DA model (Fig. 1), except that sexual reproduction was age-dependent, following a Brass polynomial, matching human reproduction dynamics [31].

Similar to the minimal DA model, the aging rate rapidly declined at the start of the simulations in the extended hDA model (Fig. 4), followed by deceleration and an eventual plateau. Increasing \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max or \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} elevated this plateau. This behavior is expected, given that the reproduction span in the hDA model is finite. Therefore, reduced extrinsic mortality elevates the threshold of the aging rate below which mortality during the reproductive period reaches near-zero, diminishing the benefit of further slow-down of aging.

However, unlike the plateau, the rate of initial decline in the aging rate was differently impacted by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max than by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} . Higher \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max, representing early survival benefits, accelerated aging, whereas higher \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} , representing late survival benefits, slowed down aging (Fig. 4c, d). This behavior of the hDA model is consistent with that of the minimal DA model (Fig. 1). Equivalent results were also obtained in simulations using asexual and age-independent, constant, reproduction, except that the effect on the plateau was less pronounced (Supplemental Fig. 3). This consistency suggests qualitative robustness of our conclusions to using other fertility functions or asexual reproduction.

In summary, simulations of the hDA model confirm a crucial impact of the timing of mortality during development on the evolution of aging that is consistent with conclusions from the minimal DA model.

Fig. 4. Early survival benefits accelerate the aging rate, and late survival benefits decelerate the aging rate in the hDA model. a. left: Evolution of aging rate in the simulation of sexually reproducing individuals for two values of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max. Starting aging rate b0 = 0.15 y^− 1^. Fertility function used: Brass polynomial, based on ref [31]. Other life history parameters: fits of the hDA model to mortality data from United Kingdom 1953. right: mean difference in aging rate after 600’000 timesteps (n = 5 independent simulations, with m = 10’000 individuals) as a function of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max. b. as a, but for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} . c. left: as a but shorter timescale. right: mean difference in aging rate after 1’000 timesteps (n = 50 independent simulations, with m = 10 000 individuals). d. as c, but for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} . Conclusions were not affected by assortative mating, different reproduction functions for males and females, and for asexual simulations (Supplemental Fig. 3). Simulations with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max < 0.055 y^− 1^ lead to population extinction and are not shown in a and c

Low mortality during reproduction weakens the selection of slower aging

To validate results from agent-based simulations, we computed the intrinsic rate of population increase r by the Euler-Lotka Equation [25] over a wide range of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max. These calculations yielded a fitness landscape that plateaued at low aging rate b, and at high \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} or \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max (Fig. 5a, b) for all fertility functions tested (Supplemental Figs. 4 and 5). Consistent with agent-based simulations, the slowest rate of aging that could evolve, marked by the edge of the plateau, increased for larger survival benefits \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max (Fig. 5a, b; Supplemental Figs. 4 and 5).

From this landscape, we estimated the selection strength for slower aging, as the change in fitness σ = Δr/Δb resulting from a 5% reduction in the rate of aging b (Supplemental Fig. 6) and computed how this selection pressure changed with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma_{max}$$\end{document} by computing \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta\lambda_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta\sigma$$\end{document} / \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta\gamma_{max}$$\end{document} (Fig. 5c, d, Supplemental Figs. 4 and 5). Across a wide range of parameters outside of the fitness plateau, increasing \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} increased the selection pressure for slower aging (red area in Fig. 5c). Conversely, increasing \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max decreased the selection pressure for most values outside of the fitness plateau (blue area in Fig. 5d). The only parameter regions where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max did not affect σ, or had opposite effects, were those with very low total lifetime reproduction or very low mortality during the reproductive period.

In summary, deterministic calculations confirm agent-based simulations for the biologically relevant parameter space, showing that the survival benefits of early and late developmental processes have opposing effects on the evolution of aging over a wide range of parameters outside of the fitness plateaus. Early survival accelerates aging (i.e., decreases selection for slower aging), whereas late-acting benefits slow down aging (i.e., increase selection for slower aging).

Fig. 5. Fitness landscape distinguishes regions where developmental processes accelerate or constrain the evolution of slower aging. (a) fitness landscape for combination of different aging rates (b) and late developmental processes ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} ). (b) same as a, but for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max. (c) numerical approximation of derivative of selection pressure for slower aging (Δr/Δb) with respect to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} . The red and blue color scale indicates normalization to the maximal positive or negative value, respectively. (d) same as c, but for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max

Discussion

While early verbal arguments by Williams suggested that extrinsic mortality leads to the evolution of faster aging [11], subsequent mathematical analysis revealed a more nuanced relationship that is influenced by population density [10, 12] and age-dependence [9, 15, 18]. The role of density-dependence has been extensively studied [10, 12]. However, while age-dependent extrinsic mortality has been acknowledged to be essential for the evolution of aging in the absence of density effects [9], the specific impact of timing, relative to the onset of reproduction, has not been systematically investigated. Here, we revealed a consistent pattern: reduced mortality early in life accelerates the aging rate, and reduced mortality later in life slows the aging rate.

The impact of survival benefits on the evolution of aging can be understood by separating two distinct effects [10]: The first effect, following Williams’ reasoning, is that improved survival of extrinsic threats increases the fraction of individuals that live long enough to be affected by aging. As illustrated in Fig. 6a-b, the first effect amplifies the relative change in mortality and consequently the lifespan gains from a given reduction in the aging rate. This effect thereby increases the selective advantage of a slower aging rate. The second, opposing, effect is that improved survival accelerates population growth, which in turn diminishes the relative fitness value of individuals reaching an old age [12]. This effect reduces the selective benefit of slower aging, thereby accelerating aging.

Fig. 6. Differential impact on population growth and lifespan gain explains opposite effect of early and late survival benefits for slow aging. a. Example of mortality rates for the GMM for two different aging rates (dark and light red), and for developmental survival benefits (cyan, right panel only) as a function of age. Net mortality is the difference between GMM and survival benefits, shown by vertical arrows. Relative change in mortality upon reduced aging rate is larger in the presence of survival benefits (right). b. Change in lifespan after reduction of aging rate for scenario with (right) and without (left) survival benefits. fast aging: b = 0.085 y^− 1^, slow aging: b = 0.065 y^− 1^. c. Summary of differential impact of early and late survival benefits on aging evolution

When extrinsic mortality is age-independent and in the absence of density effects, these two evolutionary effects precisely cancel each other out, resulting in no net selection on aging [10, 12, 15]. However, when extrinsic mortality is modulated with age, one of the two effects can dominate, generating a net selection pressure [12, 21, 22, 32]. Both early and late survival benefits increase the fraction of individuals surviving to an old age, which slows down aging via the first effect. However, because late survival benefits influence only a fraction of the reproductive window and many individuals have already died by that age, they have a weaker effect on population growth. Late life benefits thus have a smaller contribution to the second effect and consequently lead to a net increase in the selection for slower aging. Conversely, the stronger impact of early benefits on population growth tip the balance in the opposite direction (Fig. 6c).

The hDA model does not specify biological mechanisms responsible for the survival benefits and instead uses phenomenological functions to capture age-dependent mortality based on empirical human data. Several known age-dependent processes could underlie such survival benefits. Early benefits are likely driven by postnatal growth and development, while late survival benefits may stem from learning, adaptive immunity, or, in some species, complex social behaviors. The wide variety of such processes across taxa highlights the importance of considering developmentally controlled survival benefits to understand the diversity of aging rates.

Our findings on the timing of mortality-modulating factors, such as extrinsic mortality and developmental processes, offer a new perspective on experimental evolution studies that have yielded unexpected results. For example, a landmark study in guppies found that populations under high predation pressure evolved slower aging, rather than faster aging [33], contrasting predictions from classic theory [11]. Previous explanations invoked density dependence [13, 34], environmental and genetic variation [35], or changes in food abundance [36]. Our model provides an alternative interpretation: because the predators used in the study predominantly affected young individuals [37], the increase in mortality was likely concentrated early in life. According to our model, such early-life mortality reduces the selection for slow aging, aligning with the experimental observations. This matching between experiment and theory underlines the applicable significance of our model.

An important question in human biology is how ongoing changes in birth and mortality rates through cultural change will impact the evolution of aging rates. Our model suggests that the evolutionary outcome is sensitive to the balanced effects of fertility rates, the current aging rate, and extrinsic mortality. The hDA model may help to extract these key parameters from demographic data.

In addition to distinguishing the effects on the evolution of the aging rate (i.e., the slope of the aging rate over time), our model also makes predictions on the slowest rate of aging that could evolve in the absence of tradeoffs (i.e., the aging plateau). Our interpretation for the emergence of the plateau is that when survival benefits reduce mortality so strongly that aging has no impact on the reproductive output, then the selection pressure for slower aging diminishes. The model thus predicts that evolution of negligible senescence is favored when developmental survival benefits are weak. Consistently, testudines species with near-negligible senescence in captivity were found to have increasing mortality with age in the wild, suggesting that they do not express strong developmental survival benefits [38].

Conclusions

Our study shows that the timing of changes in extrinsic mortality plays a crucial role in shaping the evolutionary pressures on aging. Specifically, reducing mortality early in life diminishes selection for slower aging, whereas reducing mortality later in life strengthens it. This time-dependence arises from the interplay between two opposing forces: on one hand, lower mortality increases the number of individuals reaching old age, thereby intensifying selection for longevity; on the other hand, it also accelerates population growth, which weakens this selection pressure. Early-life mortality reductions have a stronger impact on the population growth rate and thus reduce the evolutionary advantage of slower aging. In contrast, late-life mortality reductions primarily increase the proportion of long-lived individuals, enhancing selection for slow aging. Recognizing this nuanced effect of mortality timing helps explain previously puzzling outcomes in studies of extrinsic mortality and aging evolution.

Methods

Minimal DA model

Mortality increases exponentially with time according to the Gompertz law.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${h}_{G}\: {({t})}={a}^*{e}^{bt}$$\end{document}

Survival benefits \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\beta\:$$\end{document} were a step function that decreases mortality by 20% from age \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:d$$\end{document} onwards, yielding an overall mortality

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:h\left(t\right)={h}_{G}\left(1-\beta\:\right),\:with\:\beta\:\left(t\right)=\left\{\begin{array}{c}0\:if\:t<d\\\:0.2\:if\:t\ge\:d\end{array}\right.$$\end{document}

Reproduction was sexual with constant fertility of one offspring per time step between ages 15 to 100 and zero fertility outside of this interval.

Details of the hDA model

The hDA model describes the hazard rate \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:h\left(t\right)$$\end{document} as the sum of effects of aging (α), early life mortality benefits (γ) and later life mortality benefits (λ) on survival, each of which has a distinct time dependence:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:h\left(t\right)=\alpha\:\left(t\right)-\lambda\:\left(t\right)-\gamma\:\left(t\right)$$\end{document}

aging follows the Gompertz-Makeham function [24], a commonly used model of aging:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\alpha\:\left(t\right)=c+a{e}^{bt}$$\end{document}

c, a, and b are constant parameters, and t represents time. The parameter b is called the aging rate, the exponential rate at which mortality due to aging increases with time (i.e., in the absence of developmental survival benefits).

Later life mortality benefits reduce mortality with age by a logistic function which plateaus to a maximal benefit \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} :

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\lambda\:\left(t\right)={{\uplambda\:}}_{max}\left(1-\frac{1}{1+{e}^{n\left(t-k\right)}}\right)$$\end{document}

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:k$$\end{document} is the age at which the later life benefits reach half their maximum. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:n$$\end{document} determines the steepness of the curve.

We modeled the early life benefits by a decelerating function.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:\gamma\:\left(t\right)={{\upgamma\:}}_{max}\left(1-\frac{1}{1+{t}^{h}}\right)$$\end{document}

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\gamma\:}_{max}$$\end{document} corresponds to the maximal early survival benefits and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:h$$\end{document} relates to the steepness of the curve.

Mortality data

Data on human mortality rates were acquired from the human mortality database [39] on 20 July 2023. We have analyzed the population data (cross-sectional). Analysis was restricted to mortality data after 1948. We focused on the 8 countries with a population of over 10 million in 1948. Population data was taken from 1948 to 2020 for the same countries.

Fitting the hDA model

The Least-Squares method (implemented in the scipy.optimize Python package) was used to fit the parameters of the realistic DA model to log-transformed mortality data.

Computation of fitness functions and the benefit of slower aging

We used the intrinsic rate of population increase \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:r$$\end{document} as a measure of Darwinian fitness, defined by the Euler-Lotka Equation ([25]

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\int\:}_{x=0}^{+\infty\:}{e}^{-rx}{l}_{x}{m}_{x}=1$$\end{document}

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{l}_{x}$$\end{document} is the likelihood of survival to age x. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{m}_{x}$$\end{document} corresponds to the number of offspring produced by individuals of age x. Based on the hDA model, survivorship \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{l}_{x}$$\end{document} was calculated as:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{l}_{x}=exp\left(-H\left(t\right)\right)$$\end{document}

With \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:H$$\end{document} corresponding to the cumulative hazard rate.

For mx, we used a Brass Polynomial defined as

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&\:{m}_{x}\left(v,p,w\right)\cr&=\left\{\begin{array}{c}{v\left(x-p\right)\left(p+w-x\right)}^{2}\\\:0\end{array}\right.\begin{array}{c}if\left(p+w\right)>x>p\\\:otherwise\end{array}\end{aligned}$$\end{document}

For models shown in Supplemental data, we used constant (mx = m0), increasing ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{m}_{x}={m}_{0}+0.1\:x$$\end{document} ), or decreasing ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{m}_{x}={m}_{0}-0.1x$$\end{document} ) fertility functions. For the decreasing function, fertility was not allowed to go below 0.

Fitness landscapes in Fig. 4 were based on parameters from fits to United Kingdom population mortality data from 1953. Fitness landscapes produced using parameters from other countries and years were very similar. Landscapes were generated using Rust. Numerical integration was performed using the peroxide Rust crate.

Agent-based simulations

For all agent-based simulations, time was discrete, with each time step corresponding to the age increment of 1 and with all relevant rates and probabilities defined on a per-time-step basis. Each simulation was run in 5 replicates unless stated otherwise in the main text. Simulations were carried out with a population of 10’000 individuals.

In all simulations, individual agents underwent repeated rounds of death, reproduction, and mutation. During the death part of the simulation step, the cumulative change of dying during the current step was calculated for each individual, and a Bernoulli trial determined survival. During the reproduction step, offspring were created, with a 2% chance of a mutation to alter the rate of aging (i.e., parameter b). As shown in Supplementary Fig. 7, with lower mutation rates, the qualitative dynamics of the simulations remain unchanged, except that the rate of change is significantly slower. Out of the created offspring, a random subset of offspring was selected to fill free spaces up to the population limit.

Details of agent-based simulation of the hDA model

Parameters were taken from fits to United Kingdom population mortality from 1953. Conclusions were robust to using parameters from fits to cohort or population form data from other countries and years.

Each agent was represented by four attributes: its sex, its age x, its aging rate b, and its maximum survival benefits ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma$$\end{document} max). When initializing the simulation, each agent was assigned x and b drawn from a normal distribution (mean(x) = 20, SD(x) = 10; mean (b) = 0.14, SD (b) = 0.005). \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} was initialized as a constant for all individuals. The population size was set to 10’000.

For every round of the simulation, cumulative probability of death was computed based on the hDA model for each agent. Survivorship was determined by a Bernoulli trial, followed by either increasing the age by 1 or removing the agent from the simulation.

Reproduction was modeled by the Brass Polynomial as described above using parameters obtained from human fertility data of men and women [31]. To optimize computations, simulations were discontinued for agents whose age was above the fertile period ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:w+p\:$$\end{document} from the Brass polynomial). Simulations, where post-reproductive individuals were not removed, led to qualitatively identical conclusions. However, in such a setting, the eventual existence of non-dying but no longer reproducing individuals in the control condition ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{\lambda\:}_{max}$$\end{document} = 0) led to a sudden drop of the mean b to zero (Supplemental Fig. 8). The probability p of an agent to reproduce in each round was calculated from the brass polynomial by:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\:{p}_{x}\left(v,p,w\right)=\frac{{m}_{x}(v,p,w)}{\text{m}\text{a}\text{x}{\{m}_{x}\left(v,p,w\right)\::\:x\mathbb{\:}\in\:\mathbb{\:}\mathbb{N}\}}$$\end{document}

All agents were then randomly paired as male-female couples. If both agents of a couple were able to reproduce, they produced a new agent of random sex. Simulations with constant fertility rate and assortative mating with regard to the age of the agent led to identical conclusions (Supplemental Fig. 3). The aging rate b of the new agent was assigned the mean b of its parents. 2% of all newly created agents were then randomly selected to mutate and assigned an aging rate b drawn from a normal distribution with mean being its pre-mutation value and SD = 0.012.

Parameters used in the simulations are summarized in Supplemental Table 1. Simulations were all run using Rust.

Supplementary Information

Below is the link to the electronic supplementary material.

Supplementary Material 1

Supplementary Material 2