Asynchrony of senescence among phenotypic traits in a wild mammal population.

The degree to which changes in lifespan are coupled to changes in senescence in different physiological systems and phenotypic traits is a central question in biogerontology. It is underpinned by deeper biological questions about whether or not senescence is a synchronised process, or whether levels of synchrony depend on species or environmental context. Understanding how natural selection shapes patterns of synchrony in senescence across physiological systems and phenotypic traits demands the longitudinal study of many phenotypes under natural conditions. Here, we examine the patterns of age-related variation in late adulthood in a wild population of Soay sheep (Ovis aries) that have been the subject of individual-based monitoring for thirty years. We examined twenty different phenotypic traits in both males and females, encompassing vital rates (survival and fecundity), maternal reproductive performance (offspring birth weight, birth date and survival), male rutting behaviour, home range measures, parasite burdens, and body mass. We initially quantified age-related variation in each trait having controlled for annual variation in the environment, among-individual variation and selective disappearance effects. We then standardised our age-specific trait means and tested whether age trajectories could be meaningfully grouped according to sex or the type of trait. Whilst most traits showed age-related declines in later life, we found striking levels of asynchrony both within and between the sexes. Of particular note, female fecundity and reproductive performance declined with age, but male annual reproductive success did not. We also discovered that whilst home range size and quality decline with age in females, home range size increases with age in males. Our findings highlight the complexity of phenotypic ageing under natural conditions and, along with emerging data from other wild populations and laboratory models, suggest that the long-standing hypothesis within evolutionary biology that fitness-related traits should senesce in a synchronous manner is seriously flawed.

Introduction

The identification of genetic and environmental manipulations that extend the lifespan of laboratory model organisms has revolutionised our understanding of the ageing process and is central to modern biogerontology (Partridge, 2010). It is becoming increasingly clear that senescent declines in health and function may begin well in advance of eventual mortality in both humans and laboratory organisms (Herndon et al., 2002, Papadopoulos et al., 2002, Christensen et al., 2009, Bansal et al., 2015). A question of growing importance, especially given the continued increase in human life expectancy, is whether interventions that extend lifespan in the laboratory also extend so-called ‘healthspan’, or instead leave individuals in a frail state for longer (Christensen et al., 2009, Bansal et al., 2015). A perhaps more fundamental biological question is to what degree senescence is synchronous across physiological systems and phenotypic traits in a given species (Promislow et al., 2006, Martin et al., 2007). Whilst influential evolutionary biologists have hypothesised that natural selection should shape senescence to be synchronous across physiological systems (Maynard-Smith, 1962, Williams, 1999), empirical data from humans and laboratory model organisms suggests that asynchrony is commonplace and that health- and life-span are readily uncoupled (Herndon et al., 2002, Martin et al., 2007, Christensen et al., 2009, Bansal et al., 2015). However, the benign and protected conditions experienced by laboratory model organisms and modern human societies are associated with lifespans hugely in excess of those observed under natural conditions. To understand how patterns observed in the laboratory generalise to more challenging environments and how natural selection has shaped the multifaceted process of senescence, we require studies that investigate patterns, causes and consequences of synchrony of senescence in the wild.

Age-related declines in survival probabilities and reproductive performance are widely observed in wild vertebrates (Nussey et al., 2013). Investigation of the causes of the dramatic variation in ageing observed among species, populations and individuals in nature could offer important insights into the biology of ageing (Nussey et al., 2013, Jones et al., 2014). To date, the predominant focus of studies of senescence in the wild has been on those traits most proximate to fitness (i.e. survival and fecundity), but efforts to measure other salient phenotypic traits (e.g. body mass, secondary sexual characters, parental investment, ranging behaviour) and markers of relevant physiological processes (e.g. endocrine function, sarcopenia, oxidative stress, telomere length) are rapidly increasing (Nussey et al., 2013). Within this literature, there is mounting evidence for both differences in ageing rates between the sexes and asynchrony among phenotypic traits within sexes in the way they change with age in later life (Nussey et al., 2013). The evidence includes: reproductive cessation well before the end of life in some female mammals (Packer et al., 1998, Ward et al., 2009); evidence from a range of vertebrates for asynchrony of senescence among maternal traits associated with successful reproduction (Nussey et al., 2009, Evans et al., 2011, Massot et al., 2011, Hayward et al., 2013); asynchronous senescence among male secondary sexual traits and male reproductive performance (Nussey et al., 2009, Evans et al., 2011, Kervinen et al., 2015). Furthermore, some studies have observed so-called ‘terminal declines’ in traits associated with fitness, which are to some degree age-independent and occur over the period immediately preceding death (Martin and Festa‐Bianchet, 2011, Nussey et al., 2011, Hammers et al., 2012). Conceivably, phenotypic traits or physiological measures could follow a progressive decline with chronological age, an age-independent and more sudden decline, or a combination of the two. Were declines with chronological age to predict phenotypic variation in later adulthood better than age-independent declines in relation to time to death, it suggests a physiological system selected to be maintained for a set period of time, with weaker selection beyond this. However, if age-independent declines provide a better fit to data than chronological age, it would suggest individual variation in the onset and rate of physiological deterioration. It should be noted, however, that age and time to death are typically confounded and their statistical separation can be challenging (Martin and Festa‐Bianchet, 2011, Nussey et al., 2011). To date, few studies in the wild have directly compared the synchrony of senescence patterns among more than a handful of (typically reproductive) traits and thus the evolutionary basis of apparent asynchrony in trait senescence remains poorly understood.

Long-term, individual-based ecological studies provide detailed longitudinal data, commonly encompassing the entire lifespan of individuals, from relatively long-lived species. Survival, reproductive, behavioural, genetic, biometric and, increasingly, physiological data are routinely collected in an increasing number of studies (Nussey et al., 2013). Natural systems are obviously a great deal more variable than their laboratory counterparts, and environmental sources of mortality (e.g. predation, parasites, starvation) may mean that only a relatively small number of individuals survive to experience senescence. Furthermore, the so-called ‘selective disappearance’ of lower quality phenotypes can readily mask within-individual changes with age in studies in the wild (Nussey et al., 2008). Recent studies of ageing in the wild have sought to account for environmental variability and selective disappearance by statistical means or by decomposing change with age into component processes (van de Pol and Verhulst, 2006, Rebke et al., 2010). Here, we apply the former statistical approach to examine the synchrony of senescence among a wide range of phenotypic traits in an unmanaged population of Soay sheep (Ovis aries).

The long-term study of Soay sheep resident in the Village Bay area on the island of Hirta, St. Kilda, is one of the most detailed individual-based studies of a wild vertebrate population anywhere in the world. For three decades (1985 to present day), individuals in the population have been marked and followed from birth to death, with regular re-capture of individuals, producing a wealth of information on age-specific survival, fecundity, maternal reproductive performance, male reproductive behaviour, infection with parasites, ranging behaviour, and morphology (Clutton-Brock and Pemberton, 2004). The Soay sheep on St. Kilda are unmanaged and entirely free from predation. However, they experience many environmental challenges characteristic of temperate wild vertebrate systems including over-winter food limitation, thermoregulatory challenges associated with winter weather, and infection with parasites (Clutton-Brock and Pemberton, 2004). Most mortality occurs over-winter associated with interactions among these environmental pressures, and a characteristic feature of the population is that it experiences occasional severe over-winter mortality, during which the youngest and eldest appear most susceptible (Clutton-Brock and Pemberton, 2004). Although first winter mortality can be high among lambs, individuals that survive to maturity can be long-lived: females can survive up to sixteen years (mean = 5.31 years; median = 5 years) and males to eleven years (mean = 2.67 years; median = 2 years). Phenotypic traits including body mass, parasite burdens, horn size (an important secondary sexual trait in males), and home range size and quality are all known to be associated with over-winter survival or lifetime reproductive success and are thus under natural selection in this population (Hayward et al., 2011, Morrissey et al., 2012, Johnston et al., 2013, Regan et al., in review). Separate studies of the population have documented senescence, typically from around five or six years onwards, in traits including annual survival and fecundity, body mass, parasite burden, and maternal reproductive performance traits (Robinson et al., 2006, Hayward et al., 2009, Nussey et al., 2011, Colchero and Clark, 2012, Hayward et al., 2013), but no direct comparison of senescence patterns among sexes or traits has been made. Furthermore, whilst a previous study of the population found evidence that female body mass follows a pattern of terminal decline rather than a progressive decline with age (Nussey et al., 2011), broader comparisons of such patterns among other traits have not yet been made.

Here, we undertake analyses of age-related variation in twenty different traits measured in males and females during later adulthood. Our principal aims are: (1) to determine whether age-related variation is best-explained by changes with chronological age or by changes associated with time-to-death, and (2) to statistically compare ageing trajectories among functionally-linked groups of traits in order to determine the extent to which ageing rates are synchronous among traits. Available evolutionary theory predicts synchrony of senescence in fitness-related traits (Maynard-Smith, 1962), whilst empirical data on wild populations published to date suggests that some degree of synchrony is usually observed (Nussey et al., 2009, Evans et al., 2011, Massot et al., 2011, Hayward et al., 2013).

Materials and methods

Study population and data collection

Soay sheep are descendants of domestic sheep that were present throughout northwest Europe during the Bronze age, and probably reached the St. Kilda archipelago 3000–4000 years ago (Clutton-Brock and Pemberton, 2004). A population has lived on the island of Soay since their arrival on the archipelago, but there is compelling evidence that they interbred with the now extinct dunface breed (a precursor of modern blackface sheep) sometime in the mid-nineteenth century (Feulner et al., 2013). The largest island of the archipelago, Hirta (638 ha), was evacuated of its indigenous human population and their modern domestic stock in 1930. In 1932, 107 Soay sheep were reintroduced to the island from the neighbouring island of Soay, following which they increased to carrying capacity and have since remained unmanaged. Monitoring of the population began in 1959. Since then, there have been two periods of intensive study in the Village Bay area: 1959 to 1968 (Jewell et al., 1974) and 1985 to the present (Clutton-Brock and Pemberton, 2004). Our study uses field data collected during this latter period of intensive, individual-based study of the Village Bay population, which has followed a standardised annual routine of data and sample collection.

The Village Bay study area is 170 ha, and contains approximately 30% of the total island sheep population. Three times per year (February–March, July–August, and October–December), ten censuses of the study area are carried out, during which the locations of all individuals are recorded (Clutton-Brock and Pemberton, 2004). Females give birth in spring (March–April) and on average around 13% of litters are twins, the remainder being singletons. Soay lambs develop rapidly and are typically weaned by mid-June, although they may continue to be suckled throughout the summer. Around 90% of lambs born in the Village Bay study area are caught within a few days of birth, tagged for future identification, weighed and blood and tissue sampled for genetic analysis. Daily monitoring of the population during the spring means that the day of birth of these lambs is known precisely. Each August, as many sheep from the study population as possible (usually 50–60%) are rounded up in a series of temporary traps, caught and processed over the course of a two week period. Previously unknown individuals are tagged and captured individuals have a variety of measures taken including body weight, hind and foreleg lengths, horn lengths, circumference and growth and testes circumference in males. Faecal samples are taken whenever possible and a variety of parasite egg counts are made (see below for details). During the rut (October to November) males compete for access to females as they come into oestrus and the study area is monitored daily and the identities of each ram and ewe in consort are noted. Most mortality occurs during the winter months, and regular censuses and mortality searches during the winter months mean that the carcasses of most over-winter mortalities are found and death dates can be assigned with confidence to most animals (Clutton-Brock and Pemberton, 2004).

Phenotypic traits

From the current Soay sheep database we identified twenty phenotypic traits for further investigation (Table 1). The traits were divided into six classes: vital rates (annual survival, male annual reproductive success and two traits underpinning female annual reproductive success); biometric measures (7 traits); parasitological measures (3 traits); home range measures (2 traits); maternal performance (3 traits); rut behaviour (1 trait). Each of these traits was measured longitudinally in every individual in each year where available. The traits are defined in detail within these classes below.

Table 1

Schematic of Soay sheep annual cycle and timing of data collection, with respect to an individual's age.

We calculated age based on an approximate ‘sheep year’ running from May through to April rather than using a calendar year. Thus, a sheep born in year t (typically in March or April of that year) was assigned to age zero for all traits measured between May in year t and the end of April in year t + 1, including all reproductive traits associated with its reproduction in the spring of year t + 1. This is described below, and details of the timing of trait collection are also presented in Table 1. We restricted our data selection to individuals that were known to have died of natural causes between 1985 and 2014 inclusive, and which had a known birth year. Trait measures were also restricted to the time period of 1985–2014 inclusive, although not all measures were available from 1985 (see Table 2 for more details).

Table 2

List of phenotypic traits for which GLMMs were run, with units, sample sizes and details of error distribution used. See Figure S1 for plotted trait distributions.

Trait	Units	Season of measurement	Sex	Years available	Sample size (individuals)	Error distribution
Vital rates
Annual survivala	Yes/no	Annual	Both	1985–2014	3223 (821)	Binomial
Female fecundity	Yes/no	Spring	F	1985–2014	2662 (593)	Binomial
Female twinning rate	Yes/no	Spring	F	1985–2013	1778 (479)	Binomial
Male ARS	No. of offspring	Spring	M	1987–2014	411 (176)	Poisson
Biometric
Weight	kg	August	Both	1985–2014	1394 (559)	Normal
Foreleg length	mm	August	Both	1988–2014	1238 (486)	Normal
Hind leg length	mm	August	Both	1988–2014	1297 (494)	Normal
Testes circumference	mm	August	M	1985–2014	135 (78)	Normal
Horn length	mm	August	M	1985–2014	133 (73)	Normal
Horn circumference	mm	August	M	1989–2014	128 (67)	Normal
Horn growth	mm	August	M	1988–2014	118 (64)	Normal
Parasitology
Ked count	Number	August	Both	1988–2014	1294 (493)	Poisson
FECb	Eggs/g	August	Both	1988–2014	1230 (485)	Poisson
FOCb	Oocysts/g	August	Both	1993–2014	1089 (442)	Poisson
Maternal traits
Offspring survivalc	Yes/no	Annual	F	1985–2014	1772 (468)	Binomial
Date of birthc	Days since 1st March	Spring	F	1985–2014	1912 (471)	Normal
Offspring weight at birthd	kg	Spring	F	1985–2014	1704 (459)	Normal
Male reproductive behaviour
Number of consorts	Rut census observations	Rut	M	1987–2014	411 (176)	Normal
Home range behaviour
Home range size	Ha	Annual	Both	1985–2014	2646 (658)	Normal
Home range quality	% cover Holcus	Annual	Both	1985–2014	2646 (658)	Normal

Individual identity and longevity not included in model.

Observation included as random effect in model.

Offspring sex and twin status included as fixed effects in model.

Capture age in days, offspring sex and twin status included as fixed effects in model.

Previously studies suggest that the Soay sheep have largely stopped growing and showing signs of improvement in fitness-related traits by the age of three or four, and there is little evidence for senescence beginning before five or six years of age in any trait (Colchero and Clark, 2012, Hayward et al., 2013). We therefore decided to investigate senescence patterns from age four onwards. Since very few females or males survived or were measured beyond thirteen and nine years of age, respectively, we grouped measurements from these ages onwards into a single final age class (‘13’ for females and ‘9’ for males, consisting of < 10% of individuals in both sexes). We also calculated ‘years until death’ for each trait (as the difference in years between an animal's longevity and its current age) and found that, when considering only measures from ages four or more, there were very few observations in females or males beyond eight and four years from death, respectively. We therefore grouped observations at or beyond these years until death into a final class. The final numbers of measurements and individuals measured for each trait available are presented in Table 2.

Vital rates

Survival: available data on individual birth and death date was used to construct an individual life-history for each animal. Thus, each animal was represented by as many data records as the number of years in which they were observed to be alive. Survival was monitored for an individual alive in April of year t to May 1st of year t + 1, with survival of that period coded 1 and death coded 0.

Female fecundity and twinning: Female fecundity was coded 1 if a female gave birth in that sheep year and 0 if she did not (see Table 1); female twinning rate was only scored for females which gave birth and was coded 1 if the female had twins and 0 if she had a singleton.

Male annual reproductive success (ARS): the number of paternities assigned to a male in a given year. Only males observed during rut censuses in a given year (observed at least once in a census between October and November) were included as potentially having any ARS at a given age. Paternities were assigned based on genotyping individuals at 315 SNPs in the program MasterBayes (Hadfield et al., 2006) and using 14–18 microsatellite loci (Overall et al., 2005). This method was able to infer paternity for 4593 individuals in the study population (for more details see Bérénos et al., 2014).

Biometric measures

Body weight (kg): measured to the nearest 0.1 kg at capture in August.

Leg measurements (mm): foreleg length was measured to the nearest mm as the length of the metacarpal with hoof and knee joint bent away from it; hind leg length was measured in mm from the tubercalcis of the fibular tarsal bone to the distal end of the metatarsus (Milner et al., 1999).

Horn measurements (mm): horn length was measured from the base of the horn, along the outer curvature to the tip; horn circumference was measured around the base of the horn at the point closest to the skull; horn growth was measured from the base of the skull, along the outer curvature to the first growth increment (Johnston et al., 2013). Only data from normal-horned males were included.

Testicular circumference (mm) was measured at the widest point of the scrotum as a proxy for testes mass (Preston et al., 2012).

Parasitological measures

Strongyle faecal egg count (FEC): estimated as the number of strongyle nematode eggs present per gramme of faeces collected in the August of year t using a modified McMaster egg counting technique (described in Craig et al. (2006)). Faecal egg counts were grouped into bins of 100 and counts > 1000 were collapsed into the highest bin (around 1% of samples). This is a combined count for five strongyle species which have eggs indistinguishable by eye (Gulland and Fox, 1992).

Coccidian faecal oocyst count (FOC): was estimated as the number of oocysts present per gramme of faeces, using the modified McMaster method. Coccidia are protozoan parasites, consisting of 13 species which occur at high prevalence in the population, but with oocytes of indistinguishable morphology which are grouped into one count of FOC (Craig et al., 2007). Faecal oocyst counts were grouped into bins of 100 and counts > 2000 were collapsed into the highest bin (around 3% of samples).

Ked count: the number of live keds (Melophagus ovina) counted during a one-minute search of the wool on a sheep's belly performed at capture in August. Keds are wingless ectoparasitic flies that live on the wool of the sheep and are blood feeders and can be observed in large numbers on young animals in some years (Clutton-Brock and Pemberton, 2004).

Faecal egg counts for other parasites, including the nematodes Nematodirus spp., Capillaria longipes and Trichuris ovis, and the cestode Moniezia expansa were also undertaken (Craig et al., 2008), but lacked sufficient prevalence or abundance in adults to permit meaningful analyses.

Home range measures

Home range size (hectares) and home range quality were calculated as follows. Briefly, census locations for each sheep at each age (May in year t to April in year t + 1) were collated and animals with fewer than ten census observations at a given age were excluded from further analyses. 70% kernel home ranges were calculated using the package ‘adehabitatHR’ (Calenge, 2006), and home range size was calculated in hectares. Based on the Ordnance Survey Grid, the study area was divided into 160 one hectare squares (100 m × 100 m). Between the years of 2008 and 2012, complete species lists were compiled for the vascular plants present in each hectare, and the percentage cover of each species was scored by eye (to the nearest 5%). Mean percentage covers of plant species were calculated for each individual by determining the hectares contained within the home range, extracting the vegetation cover in each hectare and weighting by the proportion of the hectare that was within the home range. Holcus lanatus was selected to represent high quality plant communities, based on previous studies of vegetation preferences and a separate principal components analysis (PCA) which included the fourteen most abundant plant species (Regan et al. in review). Home range quality was calculated as the proportion of H. lanatus cover within an individual's home range.

Rut behaviour

Number of rut consorts: Soay sheep have a highly promiscuous mating system in which females are observed to mate with many males during a 1–4 day oestrus period (Clutton-Brock and Pemberton, 2004). Male Soays exhibit mate-guarding behaviour (consorting) in order to block access to a female in oestrus by other males. Consorts were defined as a close spatial relationship between a male and female (typically within 5 m), with frequent courtship and attempted defence of the ewe by the ram (Preston et al., 2001). Daily observations were made during the rut period and the identities of males and females in consort were noted. We calculated the total number of times a male was observed in consort over the course of a rut from these observations. Note that whilst the number of consorts has been shown to be correlated with the number of paternities a male is subsequently assigned through genotyping, consort censuses have rather poor predictive power for male ARS (Coltman et al., 1999a).

Maternal performance

Date of lamb birth: the number of days after March 1st in year t + 1 on which a female gave birth. As with all maternal performance traits, the unit of analysis was the lamb, rather than the mother: thus twin lambs both received a separate data record.

Birth weight (kg): the weight, measured to the nearest 0.05 kg, of the lamb born in year t + 1 upon its first capture, which was within a week of birth in 95% of cases.

Offspring survival: whether the lamb born in year t + 1 (contributing to maternal performance at age t) survived to May 1st in year t + 2.

Statistical analyses

All analyses were performed using R version 3.2.0 (R Core Team 2012).

Generalised linear mixed-effects models: mapping age-related changes in phenotype

We began by examining the distributions of all phenotypic traits in both sexes (Fig. S1). Home range area and the number of rut consorts were both right-skewed and so were log transformed and then assumed to be normally distributed in further analyses. We initially sought to calculate age-specific predictions for each trait in both males and females (where appropriate), whilst accounting for individual and annual variation and selective disappearance effects. We followed the approach used in a previous study of senescence in different reproductive traits in wild red deer (Nussey et al., 2009). We constructed generalised linear mixed-effects models (GLMMs) for each trait using the package ‘lme4’ (Bates et al., 2014). All models included individual identity and year of measurement as random effects to account for variation among individuals and years, except survival for which individual identity was excluded (as death is observed only once per individual). In models of faecal egg count (FEC) and faecal oocyst count (FOC), we fitted GLMMs with a Poisson error distribution which also included an observation-level random effect to account for over-dispersion (following Elston et al., 2001). In all models, we included linear and quadratic terms for individual longevity as covariates, to account for covariance between lifespan and trait value (‘selective disappearance’, van de Pol and Verhulst, 2006). Additional fixed effects terms were included in maternal performance traits models. Sex of offspring and whether the offspring were twins or singletons were included as categorical fixed effects in all models, since both of these variables are known to affect lamb birth date, birth weight and survival (Wilson et al., 2005) whilst the capture age in days was included as a continuous covariate in the model of offspring birth weight, since lambs are not all captured on the same day after birth and they grow rapidly. Full details of the sample sizes and error structures used for each trait are presented in Table 2.

For each trait, we fitted and compared up to seven different GLMMs to determine whether variation in the trait was better explained by chronological age or years remaining until death and, where appropriate, whether rates of change in relation to these two variables differed between the sexes. Our ‘null’ model with respect to senescence included all fixed and random terms described above but excluded age and years remaining until death. We then fitted further models describing age-specific change in different ways. Firstly, with individual age as a categorical fixed effect, allowing the trait mean to vary independently across all ages. Then, in traits measured in both sexes (see Table 2) we also included a categorical fixed effect of sex (trait means differ between sexes but age trajectories are identical) and an age-by-sex interaction (age trajectories can differ among the sexes). We also fitted models with years-to-death as a fixed factor (instead of age), to test for variation in the trait associated with age-independent loss-of-function. Again, where the trait was measured in both sexes, we included sex as a fixed categorical variable and a sex-by-age interaction. Age and years-to-death were fitted as categorical fixed effects rather than continuous covariates because our aim was to generate age-specific estimates and standard errors for each trait in order to compare them; fitting age and years-to-death as covariates would have required using a specific functional form (e.g. linear, quadratic) and not allowed us to generate age-specific standard errors. We compared the fit of the models for each trait using Akaike Information Criterion (AIC). We then extracted the predicted means and standard errors for each age and year-before-death, after controlling for the other fixed and random effects, using the package ‘lsmeans’ version 2.17 (Lenth and Hervé, 2014).

Years remaining until death was not fitted in survival models because it is perfectly correlated with survival. It was also not fitted in models of female breeding success, twinning or maternal traits. Due to the way we structured our ‘sheep year’ (Table 1), females almost invariably fail to reproduce in their final year of life, simply because the vast majority die over the winter preceding their ‘last’ lambing season. As a result, we could not construct models including years to death that were directly comparable to other traits for these maternal reproductive performance measures.

Generalised additive models: testing for synchrony of senescence among traits

Variation in most traits was best described by age, rather than years to death (see Results) so we focussed examination of synchrony of senescence on predictions from our GLMMs including age as a categorical fixed effect. Where the trait was measured in both sexes, we used predictions from models with a sex-by-age interaction to predict age-dependent variation in both sexes. We only tested for asynchrony among traits which showed a significant change with age in the GLMMs (omitting foreleg length, ked count, FOC, female twinning and male consorts; see Results). We applied generalised additive models (GAMs) to the age-specific predictions of the remaining traits from our GLMMs to test whether groups of traits identified a priori, showed similar or different senescence patterns (following Nussey et al., 2009). GAMs fit non-parametric smoothing functions relating a trait to a covariate (in this case, age) and allow comparison of trajectories among groups of variables without the need to assume that traits all follow a particular function with respect to senescence (e.g. quadratic or higher order polynomial function).

Lower birth dates (when expressed in Julian days since 1st March) and strongyle faecal egg counts are associated with higher fitness (Wilson et al., 2005, Hayward et al., 2011), unlike other traits which would be expected to show positive associations with fitness. To ensure that all age-related changes were in the same direction with respect to fitness and function, we made the predicted values of birth dates and egg counts negative prior to further analysis. The age-specific predicted means of all traits were then standardised to be on comparable scales by subtracting the age-specific predicted mean at age 4 for each sex from the age-specific predicted mean at each age and then dividing this by the range of age-specific predicted means from age 4 to 13 in females and 4 to 9 in males. We then fitted GAMs through these standardised trait values: the value at age 4 was zero for each trait for each sex. Our standardized age-specific predictions were weighted by the inverse of the age-specific standard error from the appropriate trait GLMMs to ensure greater weight to the (younger) ages for which more data were available.

Our aim was to determine whether all traits showed synchrony of senescence, whether each trait in each sex showed distinct senescence trajectories (‘total asynchrony’), or whether the senescence trajectories among traits could be grouped based on physiological or ecological similarities or between the sexes (various forms of ‘partial asynchrony’). We hypothesised nine potential heterochrony groupings among the 15 remaining traits, and also fitted a ‘null’ age model which did not include any age term and a model which fitted a single age function to all traits (‘total synchrony’). The details of the trait groupings for the nine heterochrony models are presented in Table 4. The response variable in these comprised all age-specific predictions from all traits, with explanatory variables coding either the age or the group of that trait in each model. For instance, in the model of ‘total synchrony’ all traits belonged to the same group and an identical age function was fitted to all traits whilst in the ‘total asynchrony’ model each trait belonged to a different group and thus had a potentially distinct age function fitted. For each of the heterochrony models, we fitted separate models where: (1) the intercept associated with the grouping varied (only the mean at age 4 differed among groups, the senescence trajectory was the same); (2) the senescence trajectory associated with the grouping varied (senescence trajectories varied among groups but the mean at age 4 did not); (3) both the intercept and senescence trajectory associated with the grouping varied (both mean at age 4 and senescence trajectories varied among groups). The best-fitting model was determined by selecting the model with the lowest AIC value.

Table 4

A comparison of 11 different models of the senescence trajectories of phenotypic traits in wild Soay sheep. The first two models are the null model with respect to senescence and the model of synchrony of senescence, under which all traits follow the same trajectory. The remaining 9 models encompass different scenarios of heterochrony. For each model AIC values are presented, and for each heterochrony model three variant models were fitted, for which only the intercept or the age trajectory alone varied among groups, or for which both were allowed to vary. The model with the lowest AIC is the best fitting and is denoted in bold with grey fill.

Results

Most traits showed strong age-dependence in later adulthood (Fig. 1), and in the majority of traits chronological age explained trait variation better than the number of years remaining until death (Table 3), indicating that trait variation was largely better-explained by current age rather than remaining lifespan. Five of our 20 phenotypic traits did not show any evidence of age-dependence from four years of age onwards: foreleg length, ked count, FOC, female twinning rate and male rut consorts (Table 3). Many of the ageing patterns observed recapitulate those observed in previous studies of this system, including: (i) declining survival probability in both sexes from 5–6 years onwards (Colchero and Clark, 2012); (ii) declines in female fecundity, offspring birth weight and offspring survival with maternal age but not in female twinning rate (Hayward et al., 2013); (iii) strongyle FEC increased with age in later life in both males and females (Hayward et al., 2009); and (iv) August weight in females declined in old age (Nussey et al., 2011). That said, the GLMMs also revealed many novel aspects of the ageing process in the Soay sheep. Male ARS showed no sign at all of declining at later ages, instead increasing up to around seven years of age and remaining at around 2 offspring per year on average until the eldest age class (Fig. 1). This increasing trend was accompanied by an evident decrease in testes circumference and body mass in males (Fig. 1). Horn morphology in males showed a complex pattern of variation with age: total horn length increased, which is not surprising given that each year a new horn increment is grown; but horn circumference decreased whilst our measure of annual horn growth (length of the most recently grown horn increment) varied with age without showing an obvious trend (Fig. 1).

Fig. 1

Age-dependent variation in twenty traits measured in wild Soay sheep. Points and bars are means and standard errors estimated from generalised linear mixed models (see text and Table 2 for details), with females represented by solid lines and symbols, and males represented by broken lines and open symbols. In the “Fem. Fecundity” plot, the probability of a female giving birth at a given age is plotted in black, and the probability of her twinning (given that she reproduced) is plotted in grey. Note that male and female weight predictions are plotted separately to allow age trends to be clearly visualised given that males are much heavier than females. Units are given in the main text.

Table 3

Table of AIC values comparing generalised linear mixed-effects models of each trait which included no age or years to death term (“Null”); either age as a fixed factor or years until death (YTD) as a fixed factor, plus sex as a fixed factor where the trait was measured in both sexes; and an interaction between either age or years until death and sex, where the trait was measured in both sexes. The model best explaining variation in the trait (i.e. the lowest AIC value, unless a simpler model has an AIC value which is less than 2 higher than that of the best model) is highlighted in bold italics.

Trait	AIC value
Null	Age	Age + sex	Age ∗ sex	YTD	YTD + sex	YTD ∗ sex
Survival	2863.117	2715.22	2592.656	2599.343	NA	NA	NA
Weight	6074.058	6816.51	6039.261	6031.031	6814.792	6028.522	6018.787
Foreleg length	7151.253	7284.47	7162.22	7157.373	7279.416	7157.146	7155.941
Hind leg length	7602.844	7719.81	7614.598	7605.329	7718.115	7612.775	7598.902
Testes circumference	1173.119	1151.13	NA	NA	1156.124	NA	NA
Horn length	1398.316	1328.52	NA	NA	1354.046	NA	NA
Horn circumference	926.045	916.88	NA	NA	919.1126	NA	NA
Horn growth	1114.24	1091.05	NA	NA	1100.489	NA	NA
Keds	1093.488	1114.93	1105.941	1113.032	1109.181	1101.437	1103.797
Faecal strongyle egg count	3426.779	3504.26	3394.422	3400.567	3510.888	3402.721	3404.63
Faecal Coccidia oocyst count	5637.975	5658.64	5641.37	5650.53	5652.877	5636.516	5639.961
Home range area	4809.752	4825.97	4779.724	4767.779	4833.229	4786.812	4791.489
Home range quality	16955.24	16939.14	16938.923	16934.192	16964.77	16964.578	16960.905
Female fecundity	2950.807	2296.43	NA	NA	NA	NA	NA
Female twinning rate	1654.162	1660.27	NA	NA	NA	NA	NA
Offspring date of birth	11800.82	11790.19	NA	NA	NA	NA	NA
Offspring birth weight	1977.585	1966.58	NA	NA	NA	NA	NA
Offspring first year survival	1764.19	1756.36	NA	NA	NA	NA	NA
Male annual reproductive success	1702.948	1663.91	NA	NA	1682.089	NA	NA
Male rut consorts	1030.405	1049.155	NA	NA	1043.943	NA	NA

Of the nine traits which were measured in both sexes, three showed no change with age in either sex (FOC, ked count, foreleg length). There was no statistical support for an interaction between sex and age in survival or FEC, suggesting similar rates of decline and increase, respectively, in males and females (Table 3, Fig. 1). However, there was strong support for models including sex-by-age interactions for weight, home range size and home range quality (Table 3). Male weight increased with age until six years old and then declined thereafter, whilst female weights remained stable until around 10 years of age before declining (Fig. 1). Home range size declined with age in females but actually increased from seven years onwards in males, whilst home range quality appeared to decline in very late life (from 11 years onwards) in females but increased slightly with age in males (Fig. 1).

Only two traits were better fitted by years until death than chronological age: weight and hind leg length (Table 3, Fig. 2). There was a dramatic decline in weight of around 1 kg across the two years preceding death in females, a finding that has been documented before in this system (Nussey et al., 2011). However, males also showed this pattern, losing around 2 kg on average across the two years prior to death (Fig. 2). Hind leg length showed little evidence of change with respect to either age or years to death in females, but in males slight decreases were evident over the years preceding death (Fig. 2). Since only these two traits showed compelling evidence for age-independent declines in later life, we proceeded to focus on testing whether or not particular groups of traits followed different trajectories with chronological age in later life.

Fig. 2

Two morphometric traits (weight and hind leg length) for which years remaining until death explained more variation in our models than chronological age, plotted against years to death. Females are represented with solid lines and symbols and males are represented by broken lines and open symbols (with the sexes plotted separately for weight as in Fig. 1). Points and error bars are predicted means and standard errors from GLMMs including years to death as a factor along with its interaction with sex. Units are given in the main text.

Further analyses supported the presence of considerable and complex asynchrony of senescence both between sexes and among traits within sexes (Table 4, Fig. 3). The best supported model of senescence was the most complex, involving separate age functions for every trait in each sex (Table 4), and this model outperformed the next best by a considerable margin (AIC difference of 83.59). Models with either a single age function for all traits (synchrony of senescence) or separate functions for each sex performed extremely poorly (respective ΔAIC = 232.18 and 233.01 compared to the total asynchrony model).

Fig. 3

Predictions from the best-fitting generalised additive model (see Table 4) which supported total synchrony with all traits measured followed a different ageing trajectory in males and females. Solid lines show predicted trait values across age; broken lines show predicted standard errors. Black lines represent traits measured in females and grey lines represent traits measured in males. Abbreviations are as follows: ARS — annual reproductive success; FEC — faecal egg count; HR size — home range size; HR quality — home range quality; Lamb BWT — lamb birth weight; Lamb DOB — lamb date of birth.

Discussion

Our analyses constitute the broadest assessment of longitudinal phenotypic changes in later life so far conducted in a wild organism. The results highlight both the complexity and asynchrony of changes occurring during the senescent phase of life in both males and females in our study population. Our models clearly show that senescence trajectories across different traits are highly divergent and cannot be readily simplified or grouped (Table 4, Fig. 3). Previous studies of natural populations have documented apparent differences in patterns of senescence among phenotypic traits in birds, mammals and reptiles (e.g. Nussey et al., 2009, Lecomte et al., 2010, Evans et al., 2011, Massot et al., 2011). Studies of model organisms, notably nematode worms, also point to variability in ageing rates among different physiological systems and health measures (e.g. Herndon et al., 2002, Bansal et al., 2015). Although the present study and this previous research supports the idea that ageing is asynchronous across physiological systems and phenotypes, understanding both the generality and the evolutionary causes of asynchrony in senescence remains an important challenge within both evolutionary ecology and bio-gerontology.

Asynchrony and evolutionary theories of ageing

The existing body of work on the evolutionary theory of ageing currently offers little to help our understanding of among trait asynchrony in ageing rates within populations. Available theory that specifically addresses this issue appears to consist of a verbal model that argues for synchronous senescence arising as a consequence of natural selection. This model – first suggested by Maynard-Smith (1962) and reiterated by Williams (1999) – imagines that each trait has some critical lower value below which death occurs. Under this model, the intensity of selection for a trait increases as its mean value across the population decreases. As a consequence, a trait with a mean value that is closest to its critical threshold relative to other traits will experience the strongest selection for improvement. Conversely, trait means that are furthest above the threshold are under the weakest selection to improve. Selection against mortality declines with age, and it follows that age-specific trait values will decline with age as well. However, Maynard Smith's model argues that natural selection favours a situation in which all traits evolve to the same trait values relative to their specific critical thresholds. As a result, Maynard Smith argued that all survival-related traits should senesce at the same rate.

As discussed above, available empirical evidence suggests that Maynard Smith's prediction, and therefore the model's assumptions about how natural selection acts on senescence across traits, are wrong. Although never formalised mathematically, the model appears to rest on a form of threshold selection, in which individuals with values below some critical point have one fitness value and those above that threshold have another (Falconer and Mackay, 1996). Although such threshold relationships may exist, we hypothesise that few continuously varying physiological processes or phenotypic traits actually relate to age-specific survival in this way. It seems more likely to us that a unit of trait decline will confer some proportional increase in risk of mortality or decrease in reproductive performance. Furthermore, Maynard Smith's model considers only trait relationships with survival, neglecting equally important and potentially even more complex relationships with fecundity and reproductive performance. To develop more nuanced evolutionary predictions regarding trait synchrony, we require a better understanding of whether and how age-specific selection gradients differ among phenotypic traits. These relationships are readily obtained using classical phenotypic selection and demographic methods (e.g., Lande and Arnold, 1983, Moorad, 2014), and estimating variation in age-dependent selection across traits represents an important next step for studies of ageing in the wild. It will be similarly important to develop our understanding of the quantitative genetic relationships among traits at different ages in wild animals. Genetic correlations among traits in late life, or indeed in earlier life when selection is much stronger, may constrain or facilitate the evolution of asynchrony in senescence rates in response to prevailing selection pressures. Although age-specific genetic covariance matrices among traits have been estimated in laboratory model systems (e.g. Tatar et al., 1996), we are not aware of any such estimates from the wild (Charmantier et al., 2014).

Although Maynard Smith's verbal model predicting synchrony of senescence appears over-simplistic and is contradicted by available empirical evidence, evolutionary genetic theory can provide explanations and predictions regarding observed patterns of asynchrony in ageing rates. The challenge is to obtain the required estimates of age-specific selection gradients and age-dependent genetic covariance matrices among traits in order to parameterise evolutionary models. This is no small challenge, as these require both detailed information on the relatedness structure of the population and very large sample sizes. However, the Soay sheep system, with its SNP-based pedigree and detailed longitudinal records for a wide range of phenotypes and fitness estimates (Bérénos et al., 2014, Bérénos et al., 2015) is one of several wild systems where this may be possible.

Sex differences in senescence

Our analyses have revealed hitherto unappreciated differences in senescence patterns between the sexes in the Soay sheep of St. Kilda. Sex differences in senescence are predicted by evolutionary theory for polygynous species in which males tend to have shorter life expectancies than females (Williams, 1957, Bonduriansky et al., 2008). Specifically, males are predicted to senesce earlier and faster in such systems either because selection weakens more strongly with age or because costs of reproduction are greater in males (Bonduriansky et al., 2008). Comparative studies offer some support for this prediction (Clutton-Brock and Isvaran, 2007) and studies of other species of ungulate suggest that males show earlier or stronger declines in survival probability, body mass and annual fecundity (Mysterud et al., 2001, Toigo et al., 2007, Nussey et al., 2009). However, although males had lower average survival probability than females in later adulthood, we found no evidence for differences in the ageing pattern in later life, a finding corroborated by a previous Bayesian mark-recapture analysis in this population (Colchero and Clark, 2012). More surprisingly, whilst female fecundity and key maternal reproductive performance parameters (offspring birth weight and survival) declined approximately linearly with age, male annual reproductive success increased and then plateaued, but did not decline with age. This result is at odds with previous studies, which suggested declines in male ARS from 6 or 7 years onwards (Coltman et al., 1999b, Robinson et al., 2006). However, these studies used less well resolved, microsatellite-based pedigrees and at least a decade's less data. In both studies, the apparent decline in ARS in very old males is based on very small samples and was not directly shown to be statistically significant. Our analyses, based on a larger sample size and a more complete pedigree, suggests that there is no evidence for a decline in male ARS in very old age in this population.

Our analyses do offer evidence that males are senescing, however. Testes size declined with age, a result mirrored in previous study suggesting that testosterone levels during the annual rut also decline from around 4 years onwards (Preston et al., 2012). More strikingly, males lost on average 2 kg over the two summers before they died (around 5% of average body mass at four years of age), almost double the comparable loss documented here and previously in females (Nussey et al., 2011; note that females are around 35% lighter than males at 4 years). Our results may reflect alterations in the behaviour and life history tactics of elderly males to maintain reproductive fitness in the face of competition from healthier, younger males. Two non-exclusive possibilities suggest themselves: first, elderly males could be showing some form of ‘terminal investment’ (Clutton-Brock, 1984), in which all available resources are ploughed into the remaining reproductive attempt following the onset of physiological senescence, and second, their greater accrued experience in the rut and knowledge of the habitat could be allowing them to compensate for loss of function. Further analyses of available behavioural and life history data would allow us to test these possibilities and gain a better understanding of how sexual selection and senescence interact in this population.

Age-related changes in ranging behaviour

The pronounced and sex-dependent ageing patterns in home range size and quality are rare examples of behavioural changes in late life from a wild mammal. Evidence from seabirds suggests that ranging behaviour associated with foraging during reproduction may be altered in older individuals, and this effect may be sex dependent. In grey-headed albatrosses, older males were found to take longer foraging trips during incubation and show lower mass gain over the trip compared to middle-aged males (Catry et al., 2006). One study of wandering albatrosses found striking evidence that incubating older males – but not females – tend to forage in more southerly Antarctic waters and spend longer away from the nest and more time flying between spells on the water (Lecomte et al., 2010). However, no evidence of age-related changes in similar foraging metrics was found in either sex in a different population of the same species (Froy et al., 2015). Evidence for changes in space use in late life in mammals has thus far been limited to a study of the locations of moose carcasses killed by wolves on Isle Royale (Montgomery et al., 2014). This study showed differences in the locations of wolf-killed ‘senescent’ moose (identified based on degree of osteoarthritis or periodontal disease rather than age) compared to ‘non-senescent’ animals. Senescent moose were more likely to be found in habitats associated with lower predation risk, suggesting habitat selection changes associated with age-related pathologies (Montgomery et al., 2014). Our finding that female Soay sheep have smaller home ranges as they age could be associated with reduced mobility. Whilst sarcopenia-like changes in muscle structure and osteoarthritis have been detected in wild mammals (Hindle et al., 2009a, Hindle et al., 2009b, Peterson et al., 2010, Arthur et al., 2015), we think this an unlikely explanation. Home ranges and foraging routes are rather limited in Soay sheep; certainly compared to the thousands of miles travelled by foraging albatrosses. Furthermore, we would expect to see similar changes in both sexes but instead observe an increase in home range size with age in males. It may be that older females are competitively excluded from some areas of higher quality grazing by younger conspecifics. However, the striking differences between the sexes point to marked sexual differences in behavioural changes in later life. The larger average home range sizes in males are largely due to their behaviour during the rut when they roam widely in search of females, and it seems likely that variation in male behaviour at this time of year drives age-related variation in home range size and quality. On the other hand, female behaviour may be expected to be more consistent throughout the year. Further work that explores seasonal differences in ranging behaviour is clearly warranted to better understand the patterns observed here, and the possibility that age-related changes in adult ranging behaviour and habitat use could underpin sex differences in age-related fitness declines in the wild certainly deserves further investigation.

Conclusions

The assumption that mutations, drugs and environmental interventions that extend lifespan will also extend healthspan rests, at least to some degree, on synchrony of senescence among different organs, systems and phenotypic traits. Yet asynchrony of senescence has been documented in humans, laboratory models and wild animals (Walker and Herndon, 2010, Nussey et al., 2013, Bansal et al., 2015). The present study represents the most striking evidence to date of asynchrony of senescence among phenotypic traits from a natural population. Laboratory models offer exquisite insights into effects of single genes, pathways and interventions in a single environment on lifespan and the maintenance of physiological function. However, the use of non-model systems, including those in the wild, offer important and complementary insights: from the identification of novel pathways and mechanisms that might regulate ageing and repair physiological damage (Austad, 2010) to a more general understanding of the cause of variation in ageing within genetically heterogeneous populations of long-lived species experiencing challenging environments (Nussey et al., 2013). Importantly, in the context of questions about the synchrony of senescence, evolutionary theory and studies of wild populations can help explain how and why natural selection under variable environments may couple or uncouple senescence across physiological systems and phenotypes. Studies like the present one are a descriptive first step in this process. Subsequent research determining the genetic basis of asynchrony of senescence, and estimating age-dependent selection on different traits in natural populations should help illuminate the evolutionary origins of asynchrony of senescence.

The following is the supplementary data related to this article.

Fig. S1

Histograms of all non-binary traits used in the present analyses, with females in red and males in blue. Note that sex distributions are overlaid on one another and not stacked. Only data from sheep aged 4 and older are included. On the x-axes, number of rut consorts, faecal strongyle egg count, faecal coccidia oocyst count and home range area are natural log-transformed.