Within-individual changes reveal increasing social selectivity with age in rhesus macaques

Significance The narrowing of social networks and prioritization of meaningful relationships with age is commonly observed in humans. Determining whether social selectivity is exhibited by other animals remains critical to furthering our understanding of the evolution of late-life changes in sociality. Here, we test key predictions from the social selectivity hypothesis and demonstrate that female rhesus macaques show within-individual changes in sociality with age that resemble the human social aging phenotype. Our use of longitudinal data to track how individuals change their social behavior within their lifetimes offers the most conclusive evidence to date that social selectivity is not a phenomenon unique to humans and therefore might have deeper evolutionary underpinnings.


Separating within-from between-subject effects
In these analyses we were specifically interested in how social behavior changes across an individual's lifetime, that is, in the within-subject effect of age. To separate the within-from the between-subject age components we used a within-subjects centering approach (as per 1, 2). In a standard random effects linear regression, the estimated effect of a fixed term in the model includes both the within-subject effect and the between-subject effect, and is given by the following regression equation: Equation 1 where xij is the x value of measurement i from subject j. β0 is the intercept of the regression equation and β1 is the slope showing the effect of predictor variable, x, on the response variable, y, for the i th measurement of subject j. The terms u0j and e0ij denote the random intercept and residual variance, respectively.
To isolate a within-individual effect, we separated chronological age (xij) into two terms (see Figure 1D). The between-subject effect was obtained by taking the average age x̄ of subject j, and was calculated by taking the mean age over all years that an individual was observed, resulting in one value of average age per individual x̄j ( Figure 1D). The within-subject effect (here called within-individual age) was calculated by subtracting average age from each age at which the individual was observed (xij -x̄j) resulting in many values per individual (one for each observation, centered around 0; Figure 1D). This within-subject centering model is given by the regression equation: yij = β0 + βW(xij -x̄j) + βB x̄j + u0j + e0ij Equation 2 Here βW is the slope for the within-individual effect of age (xij -x̄j) and βB is the slope for the between-individual effect of average age x̄j. This βB term therefore describes how the response variable (here some measure of social connectedness) is related to differences in age between individuals while βW describes how the response variable changes with age across an individual's lifetime, relative to their mean age.
We modified the equation above to account for the possibility that age-related changes in social behavior could be driven by selective disappearance of individuals with high or low measures of social connectedness. For instance, less social individuals may be more likely to die young because of poorer access to resources (3)(4)(5). To test for the presence of selective disappearance we transformed the equation above to test whether the slopes βB and βW differed significantly from one another (1, 2). The new formulation of the equation is as follows: yij = β0 + βWxij + (βB -βW) x̄j + u0j + e0ij Equation 3 Here βW still represents the within-individual effect of age xij and is equivalent to βW in Model 2 while the coefficient (βB -βW) of average age x̄j represents the difference between the betweenindividual and within-individual age terms. If this coefficient of average age is significantly positive it would indicate that individuals with low social connectedness (or with low values of the given response variable) disappear selectively from the population because the slopes of the within-versus between-subject age differ significantly (1, 2).

Analyzing who aging females directed their approaches toward and received approaches from
We conducted a post-hoc analysis to assess whether older females directed approaches to, and received approaches from important partners (i.e., kin, strong, and stable partners). To do so we subset our data to only include the last year of data that we had for all subjects. For all subjects and their partners present in this "last-year" dataset we calculated the mean dyadic sociality index (DSI) using all previous years of data for the subject and their partner. For all subjects present in the last-year dataset we also assessed whether their partners were "stable" social partners (recorded as a categorical variable no/yes). Dyads were considered to be stable partners if they had a DSI > 0 for at least two consecutive years. Additional details of this methodological approach can be found in the Methods section of the main manuscript (see Predictions 4b and 4c). We predicted that individuals who were strongly connected to the subject earlier in life and who were stable partners with the subject earlier in life would be more likely to approach and be approached by the aging female than those individuals who previously had a weak or unstable social relationship with the subject.
All statistical analyses were conducted at the level of the dyad. We fitted whether or not individuals approached or were approached by the aging female (coded as 0/1) as the response variables of interest in separate models and used a Bernoulli error distribution (logit link). For the two models looking at whether approaches were predicted by the strength of relationship between the female and her partner, the predictor variable of interest was the mean DSI for the focal individual and their partner -calculated based on all previous years of interaction. This was included in the model as a continuous fixed effect. For the two models looking at whether approaches were predicted by the stability of the relationship between the female and her partner, the predictor of interest was whether or not the focal individual and partner were stable social partners in previous years (no/yes). We included this in the model as a categorical fixed effect. We also included the partner's age (continuous), the partner's rank (categorical) and whether or not the partner was kin or non-kin (categorical) as fixed effects in all models (see Tables S17-S20). The inclusion of a fixed effect for kin or non-kin allowed us to assess whether females were preferentially approaching or being approached by related individuals. In each model we checked for an interaction between the predictor of interest (mean DSI and partner stability) and focal individual age to assess whether the likelihood of approaching or being approached by a strong or stable partner was dependent on a female's age. We removed the interaction term from the model when it was not significant. We included group and year as random effects to account for variation in approaches that might be due to differences between groups or years. Individual ID and partner ID were included as random effects in a multi-membership grouping term (6). This multi-membership grouping term accounts for the inherent multilevel structure of the data and allows each sample (dyad) to belong to more than one individual in a random effect at the same time.
We ran the model with the following weakly informative prior means and standard deviations (µ, σ): intercept (0, 1), mean DSI (0, 1), partner age (0, 1), partner rankL (0, 1), partner rankM (0, 1), partner relatedness_nonkin (0, 1), mean DSI:age (0, 1).  Table S12. Fixed and random effects from models looking at whether the stability between a female and her partner (i.e., being partners for at least two consecutive years) predicts the probability of those individuals being partners later in the female's life. Bolded terms indicate fixed effects where the 95% credible intervals did not overlap zero, providing evidence that those effects were significantly different from zero.
We ran the model with the following weakly informative prior means and standard deviations (µ, σ): intercept (0, 1), mean groom rate (0, 1), partner age (0, 1), partner rankL (0, 1), partner rankM (0, 1), partner relatedness_nonkin (0, 1), mean groom rate:age (0, 1).  Table S14. Fixed and random effects from models looking at whether the mean proximity rate (i.e., proximity strength) between a female and her partner predicts the probability of those individuals being proximity partners later in the female's life. Bolded terms indicate fixed effects where the 95% credible intervals did not overlap zero, providing evidence that those effects were significantly different from zero.
We ran the model with the following weakly informative prior means and standard deviations (µ, σ): intercept (0, 1), mean proximity rate (0, 1), partner age (0, 1), partner rankL (0, 1), partner rankM (0, 1), partner relatedness_nonkin (0, 1), mean proximity rate:age (0, 1).   Table S15. Fixed and random effects from models looking at whether the relationship stability between a female and her grooming partner (i.e., being grooming partners for at least two consecutive years) predicts the probability of those individuals being grooming partners later in the female's life. Bolded terms indicate fixed effects where the 95% credible intervals did not overlap zero, providing evidence that those effects were significantly different from zero.
We ran the model with the following weakly informative prior means and standard deviations (µ, σ): intercept (0, 0.5), mean DSI (0, 1), partner age (0, 1), partner rankL (0, 1), partner rankM (0, 1), partner relatedness_nonkin (0, 1).  Table S18. Fixed and random effects from a post-hoc analysis looking at whether, at their last time point in the dataset (i.e., oldest age), females were more likely to be approached by partners with whom they had previously been strongly connected (i.e., had a high mean DSI). Bolded terms indicate fixed effects where the 95% credible intervals did not overlap zero, providing evidence that those effects were significantly different from zero.
We ran the model with the following weakly informative prior means and standard deviations (µ, σ): intercept (0, 0.5), mean DSI (0, 1), partner age (0, 1), partner rankL (0, 1), partner rankM (0, 1), partner relatedness_nonkin (0, 1).  Table S19. Fixed and random effects from a post-hoc analysis looking at whether, at their last time point in the dataset (i.e., oldest age), females were more likely to approach partners with whom they had previously had a stable social relationship (i.e., been partners for at least two consecutive years). Bolded terms indicate fixed effects where the 95% credible intervals did not overlap zero, providing evidence that those effects were significantly different from zero.