Socioecological factors and partner preferences in the play behaviour of wild vervet monkeys, Chlorocebus pygerythrus

.

Play has been recorded in every observed species of primate (Burghardt, 2005) and is likely to serve multiple, overlapping adaptive purposes (Bateson, 1981).Sex differences in social play have been classically proposed to reflect the diverging needs of adults.For example, the observation that the juvenile males of many species often engage in higher levels of social play compared to females is hypothesized to prepare them for the higher levels of aggression, fighting, defence and dominance-related behaviours they participate in during adulthood (Marley et al., 2022;Paukner & Suomi, 2008).However, recent research by Marley et al. (2022) challenges this classical perspective, as they revealed variation in social play is not strongly phylogenetically patterned, and there are fewer sex differences in social play across many mammalian species than expected.Marley et al. (2022) suggested that any sex differences that do arise can result from social and ecological factors beyond preparing animals for adulthood.For example, ecological factors such as resource distribution or predation risk can influence group composition and population size, which in turn impact group social structure, mating strategies and dominance hierarchies (Marley et al., 2022).Such social structures and mating strategies are highly variable across mammalian species, and intrasexual competition in males and females may be one potential socioecological predictor in variation in sex differences in social play (Marley et al., 2022).
There are two overlapping hypotheses in the social play literature related to intrasexual competition.The first, the motor-training hypothesis, states that the development or acceleration of motor skill acquisition depends on the rate and intensity of play, which could later translate into fight or flight competence (Bergh€ anel et al., 2015;Byers & Walker, 1995).Support has been offered for this hypothesis in a longitudinal study of Assamese macaques, Macaca assamensis, where male juveniles play more and acquire motor skills more quickly compared to female juveniles, at the expense of growth (Bergh€ anel et al., 2015).Typically, adult male mammals engage in higher intrasexual competition in order to secure mates (Kokko & Jennions, 2008), and motor competence may increase fitness by assisting in predator avoidance, preventing injury during fights or increasing future dominance ranks (Bekoff, 1988;Blumstein et al., 2013).In species where both females and males experience high levels of intrasexual competition (e.g.access to ecological resources such as food or access to coalition partners; Stockley & Campbell, 2013;Clutton-Brock, 2009) and both engage in physical combat, we may expect equally high levels of social play during adolescence as both sexes would value from practising fighting skills (Marley et al., 2022).Such competence in motor training may also lead to behavioural flexibility in physically and cognitively preparing for future unexpected events (Bergh€ anel et al., 2015;Graham & Burghardt, 2010;Spinka et al., 2001).The second, the social skills hypotheses, argue that play may strengthen social bonds or assist in assessing dominance relationships, which may promote the development of social skills and the integration of juveniles within a group (Graham & Burghardt, 2010;Shimada & Sueur, 2018;Smith, 1978).Social play may be a means of evaluating the physical and competitive capabilities of others and may establish or uphold dominance relationships without the risks associated with overt aggression (Cafazzo et al., 2018).Determining dominance relationships in the moment may therefore immediately benefit (or disadvantage) the individual (Pellis & Pellis, 2009) and could potentially influence future dominance hierarchies (Blumstein et al., 2013;Cafazzo et al., 2018).In species where sexes differ in patterns of philopatry and dispersal, the dispersing sex may engage in social play behaviours more than the philopatric sex in order to practise fighting skills needed in combat and develop social skills needed for integration in new social groups (Burghardt, 2005;Marley et al., 2022;Mitani et al., 2012).Furthermore, if play is used as a tool to assess such social relationships, increase affiliation or forge alliances, we may expect to see a positive relationship between the amount of social play and factors that demonstrate or maintain affiliation (e.g.kinship, grooming, spatial proximity; Mancini & Palagi, 2009).
Here, we test both the motor-training and social skills hypotheses of social play in wild vervet monkeys, a species for which juveniles require a sophisticated set of social skills needed to integrate into complex adult societies and dominance hierarchies (Young et al., 2017), and adults have both sex-specific roles (e.g.patterns of dispersal; Isbell et al., 1991) and shared behaviours (e.g.high intrasexual competition in both sexes).At the same time, male vervets, who are the dispersing sex, are more likely to face increased risk of predation (Bonte et al., 2012) and require developed fighting skills needed for integration (Marley et al., 2022), while philopatric female vervets also require developed fighting skills to participate in intertroop encounters (Isbell et al., 1990).Despite a male-bias in body mass, female vervets are co-dominant and have the potential to be useful allies during aggressive interactions (Young et al., 2017).Consequently, both vervet males and females should benefit from an acceleration of motor skill development if these skills translate into fighting competence.If social play provides some form of motor training, allowing individuals to gain skills or develop motor control they may need during combat, (Prediction 1) we would expect high frequencies of adolescent play in both male and female vervets given the high levels of aggression within and across both sexes (McFarland et al., 2021;Vilette et al., 2024;Young et al., 2017).
Social play may also afford benefits related to dominance rank, both during adolescence and adulthood.Vervets usually have linear dominance hierarchies, which are relatively stable in females, but which may nevertheless vary across sites in their steepness, suggesting local differences in competition for resources (Henzi et al., 2013).Male dominance ranks are more variable over time and point to the importance of physical condition in establishing dominance (Bramblett et al., 1982).If social play serves an immediate social function, by testing or establishing dominance relationships (Blumstein et al., 2013;Seyfarth & Cheney, 2003), we predicted (Prediction 2) that individuals would play with more closely ranked individuals, as these closely ranked relationships would need to be tested more frequently for an individual to determine its dominance rank in relation to others.One could theoretically play with every individual in the group to accurately place oneself within the dominance hierarchy, but this is unnecessary if other social information is available (Lutz et al., 2019) and, in primates, social information diffuses quite rapidly, allowing for the judgement of dominance relationships indirectly (Seyfarth & Cheney, 2003).
Finally, if play is used to assess social relationships or used as a tool to increase affiliation between individuals, we expected (Prediction 3) a positive relationship between play frequency and behaviours/conditions that support affiliation such as grooming and close spatial proximity, and we would also expect higher rates of social play to occur between kin rather than nonkin.

Study Site and Animals
We collected these data on the Samara Private Game Reserve, Eastern Cape Province, South Africa (32 22 0 S, 24 52 0 E) on three troops of vervet monkeys (RBM, RST, PT) who reside within a dense strip of vegetation along the Melk River.This area supports a high population density and uncharacteristically large vervet monkey groups, with a mean of 40 animals and an upper limit exceeding 70 (Pasternak et al., 2013) due to reduced dispersion possibilities away from the river.Troop territories can overlap extensively and vary from approximately 64 ha to 176 ha (Pasternak et al., 2013).In this semiarid environment, troops are vulnerable to the effects of periodic and severe drought, and group sizes therefore fluctuate markedly over time (C.M. Nord, 2021).

Baseline Data Collection
Data for these analyses were collected between June 2014 and June 2017.Troops were followed for 10 h a day, 5 days a week, by at least one observer on each study day.All troops were fully habituated to human observers, with RBM and RST having been studied since 2008 and PT since 2012.All troop members were individually identifiable from natural markings (Pasternak et al., 2013) and maternal sibling relationships were largely known.For this study, 'kin' were defined as having the same mother and 'kinship' refers solely to maternal sibling status.Study population structure across the study period is provided in Table 1.
We used electronic data loggers to collect scan samples from all troop members visible within a 15 min period every 30 min (Altmann, 1974).We recorded each animal's behavioural state (social, resting, moving, foraging), the identity of and distance to (in metres) their nearest adult female, adult male and juvenile, and, for animals engaged in social behaviour (e.g.grooming, play), the identities of their social partner.As the field site was characterized by high visibility and all animals were fully habituated, we are confident that there was no systematic bias in the number of scan samples collected from different individuals (C.M. Nord, 2021).
Over the 3-year period, we collected a total of 120 529 scan samples of juveniles, from which we derived 5145 individual play dyads.We constructed an activity budget to reflect the annual allocation of time by the adults and juveniles of all three troops to six activity categories (Table 2), and given that adult play was rare, our analysis is directed at juvenileejuvenile dyadic play.Any play bout that contained more than two individuals was excluded to avoid repeat data (Shimada & Sueur, 2018;Wang et al., 2021), as were bouts with unidentified partners, leaving us with 4208 play bouts across 949 unique dyads.Descriptive statistics of play bouts across ageesex categories are presented in the Appendix, Table A1.
Play has been defined as behaviour that is (1) not entirely functional, (2) spontaneous and voluntary, (3) distinct from serious behaviour in terms of structure or context, (4) recurrent in similar manners and (5) executed in the absence of stress (Burghardt, 2005).We considered social play to be an interaction between two or more individuals that could include play fighting, play chasing or teasing, as defined by Petr u et al. (2009).Following Shimada and Sueur (2018), we did not differentiate different types of play in these analyses.Juveniles were defined as those who had not yet reached sexual maturity (from birth to less than ~3.5 years for females, less than ~5 years for males; Jarrett et al., 2018).Births in this population are seasonal and all juveniles were allocated to a birth cohort (Blersch et al., 2023).As the dates of birth of juveniles were known, we calculated age differences (days) between members of each play dyad.

Dominance Interactions
We recorded all observed dominance interactions and identified the aggressor (the individual that initiated the aggressive interaction), the victim (the individual that received the aggressive interaction) and the outcome of the interaction from the perspective of the aggressor (win, lose, draw, unknown).Dominance interactions included both aggressive interactions (i.e.directed threat displays, chasing or physical contact) and submissive interactions (i.e.displacing or supplanting another from a resource; see Nord et al., 2022).We recorded 32 398 dominance interactions and their outcomes from adults and juveniles over the study period.
Of these dominance interactions, 3387 (10.5%) included physical contact, where adult females initiated the highest proportion of the interactions (49.2%), followed by juvenile females (20.3%), adult males (18.5%) and juvenile males (12.2%).Adult females were also the receivers of the highest proportion of physical contact (36.8%), followed by juvenile females (27%) and juvenile males (16.5%), with adult males receiving the lowest proportion (16.5%).Injuries obtained were documented during daily censuses.Instances of injuries resulting from predation attempts or accidents were notably infrequent, and it was assumed that injuries primarily stemmed from conspecific aggression (McFarland et al., 2021).

Construction of Social Networks and Simple Ratio Index
We constructed social networks using grooming and proximity data collected during behavioural scans of juveniles.We collected a total of 7542 directed grooming dyads and 114 693 spatial proximity dyads.
Familiar social partners, which we defined as social partners outside of play, were identified through grooming interactions and spatial proximity, using the simple association ratio index (SRI).The SRI uses observational data to estimate the probability of observing two individuals together, given that one of the pair has been seen.The SRI was calculated for each dyadic pair using the formula below developed by Cairns and Schwager (1986) where X is the number of times individuals A and B were observed in the same group, as some males immigrated between study groups during the study period.Y AB is the number of times both individuals were observed in different groups, Y A is the number of times individual A was observed without individual B, and Y B is the number of times individual B was observed without individual A.
We then used the 'netTS' package with the 'create.a.network' command and specified 'SRI ¼ True' (Bonnell & Vilette, 2020) to generate SRI estimates for the grooming and nearest-neighbour spatial associations (proximity SRI) for each dyad.

Calculating Dominance Rank
Based on the actual temporal sequence of dominance interactions, we calculated sequential estimations of individual dominance strengths as Elo-ratings (Albers & de Vries, 2001).We estimated Elo-ratings using the package 'EloRating' (Neumann et al., 2011).Higher Elo-scores correspond to more dominant individuals (winning more outcomes) and scores decrease as individuals are more subordinate (losing more outcomes).Individual Elo-ratings can be used to track rank changes and estimate rank order within a group at any specified time, making this method advantageous when frequent changes in group composition occur (Vilette et al., 2020).The number of points allotted to an individual after an outcome is determined by the k value, which was set to the default value of 100 (Neumann et al., 2011).

Statistical Methods
All data were analysed in a Bayesian framework, using R version 1.2.1335 (R Core Team, 2018).To determine which demographic variables predicted play frequency, we generated a hierarchical model using the 'brms' package (Bürkner, 2017), with frequency of play dyads as the response variable.Because the response variable was dyadic, so too were the fixed and random effects.We entered the sex, age difference in days, rank difference, kinship, grooming network SRI and proximity network SRI between the dyads as predictor variables.Focal subject and partner identities were entered as crossed random effects.We set weakly informative priors centred on zero (i.e.normal (0,1)) for our model.All numeric predictor variables were centred and scaled.We set mixed-sex dyads as the reference level for participant sex.We accounted for variation in the number of scan samples by setting the total number of scan samples of the combined dyad over the study period as an offset variable.We also accounted for difference in year of birth and thus any effect of juvenile cohort by including it in the model as a fixed effect.We ran four chains for 3500 iterations and confirmed convergence of the chains (i.e.all b R ¼ 1.00).Bulk effective sample size estimates (Bulk-ESS) and tail effective sample size estimates (Tail-ESS) confirmed that the posterior means and medians were reliable.
We used the posterior predictive check ('pp_check') function to assess the adequacy of model performance (Gabry et al., 2019) and the 'testdispersion' function of the 'DHARMa' package to test for dispersion issues (Hartig, 2017).Given that our response variable of interest (play bouts) involved count data, we initially used a Poisson distribution; however, posterior predictive checks indicated dispersion issues, which are common in Poisson models (Vilette et al., 2022).We changed to a negative binomial distribution, and when this did not resolve the dispersion issue, we fitted the model with a hurdle negative binomial distribution (Hilbe, 2017;Vilette et al., 2022).The hurdle negative binomial model had no dispersion issues and was qualitatively similar to the Poisson model (see Supplementary Material), and we report the findings from this model in the main text.Posterior density distributions and other graphical outputs were generated with the 'ggplot2' package (Wickham, 2016) and the 'bayestestR' package (Makowski et al., 2019).
Credible intervals (CI) were set to 95% as these bounds, while not prescriptive (McElreath, 2020) are interpretatively familiar (Henzi et al., 2021;McElreath, 2020).Similarly, we also estimated the probability of direction ('pd') for each effect using the 'bayes-testR' package (Makowski et al., 2019).The pd estimate ranges from 50% to 100% and uses the posterior distributions to determine the certainty of the direction of an effect (either positive or negative).A pd of ~97.5%, ~99.5% and ~99.95% may be considered to indicate weak, moderate and strong evidence for an effect, respectively (Colquhoun, 2014).We extracted conditional and marginal R 2 values from each model using the 'bayes_R2' function from the 'brms' package (Bürkner, 2017;Gelman et al., 2019).We used the 'pairs' function in R to visually assess the presence of collinearity between variables in our model and found none (see Supplementary Material).

Ethical Note
Data collected here were purely observational and adhered to the laws and guidelines of South Africa and Canada.Procedures were approved by the University of Lethbridge Animal Welfare Committee (Protocol 1505).This study also adheres to the ASAB/ ABS Guidelines for the use of animals in research.

RESULTS
We recorded play behaviours from 176 of 243 potential femaleefemale juvenile dyads, 160 of 217 potential maleemale juvenile dyads and 306 of the potential 489 mixed-sex juvenile dyads that had the opportunity to play during the study period.Model main effects for posterior estimates, standard errors (SE), upper and lower 95% CIs, bulk and tail effective sample size estimates and probability of direction estimates are presented in Table 3.
Compared to mixed-sex dyads, we found that maleemale dyads played more frequently (Fig. 1).Play frequencies of femaleefemale dyads were no different from mixed-sex dyads (pd ¼ 80.03%).Rank difference had an effect on play bout frequency across all dyad sex categories (Fig. 2), where the smaller the difference in dominance rank between a dyad, the higher the play frequency between them.There was a strong effect of kinship, with maternal siblings being more likely to play than nonsiblings.This was consistent across all dyad sex categories (Fig. 3).After controlling for cohort identity, there was still strong evidence of an effect of age difference on play frequency, where juveniles closer in age played more (Fig. 4).We found no effect between play frequency with grooming network SRI (pd ¼ 93.77%) or proximity network SRI (pd ¼ 61.91%).Our model also revealed variation with respect to group level effects (Table 3).There was variation across individual focal identities and across individual partner identities.

DISCUSSION
Our results indicate that the Samara population of vervet monkeys preferentially chose social play partners.Contrary to our first prediction, we found that maleemale dyads had a higher frequency of play compared to mixed-sex and femaleefemale dyads.As predicted, we found a preference to play with individuals of similar dominance rank, suggesting play may be used as a tool to determine one's dominance rank in relation to others.Additionally, we found that vervets were more likely to play with kin rather than nonkin, indicating that social play may be used to assess social relationships or maintain affiliation between kin.However, we found no evidence that vervets play preferentially with familiar social partners or individuals of close proximity.Taken together within the context of social skills hypotheses, our findings indicate that play in both male and female vervets may serve social functions, such as testing or establishing dominance relationships, navigating rank trajectory and/or affiliation maintenance among kin.
We found that maleemale dyads had higher frequencies of play compared to mixed-sex and femaleefemale dyads.This was unexpected due to high levels of physical aggression occurring between and across the sexes and because females are co-dominant (see Supplementary Material, Fig. S8), emphasizing the need for both sexes to be physically competent.With respect to sex differences in play behaviour, primates exhibit considerable variation, with many species displaying mixed outcomes or lack of any discernible sex differences (Marley et al., 2022).However, playrelated sex differences appear to be variable over time in some species and may be influenced by the availability of same-sex and opposite-sex peers.For instance, male rhesus macaques, Macaca mulatta, raised in mixed-sex peer groups displayed higher frequencies of play than male infants reared in social groups with only male peers (Goldfoot & Wallen, 1978).Fluctuations in maleefemale ratios in birth cohorts may therefore influence play frequencies from year to year in our study population.How social play is measured can also influence what conclusions are reached regarding sex differences.Play can be measured through various metrics, including rates, initiations and frequencies, and studies may choose to focus on specific types of play or movements.Studies can also vary widely based on factors such as context, sample size and the statistical methods employed (Marley et al., 2022).
Our findings related to dominance rank corroborate other evidence that mammals initiate and play more with individuals that they can potentially dominate during play (Biben, 1986;Owens, 1975).Experiencing subordinate or disadvantaged positions during play may be of value as it allows for the individual to practise defensive strategies that could be used in an aggressive interaction (Ward et al., 2008) and may help an individual prepare emotionally for future unexpected or stressful events ( Spinka et al., 2001).However, determining the extent to which play behaviour has benefits relating to dominance rank is complicated, given that the dominance hierarchies in males and females in this population are extensively intertwined (Young et al., 2017) and that males and females have different strategies (outside of play behaviour) that influence rank.For instance, Young et al. (2017) found that males that interacted more frequently with female spatial associates and that engaged in grooming sessions with well-connected females were more likely to enhance their dominance rank.As such, we

Male-male
Male-female Female-female 0 5 Play rate Low rank difference High rank difference  3).

Male-male
Male-female    3).may also expect there to be sex-specific benefits related to play and dominance rank.
Previous data primarily derived from captive vervet populations (Bramblett et al., 1982;Fairbanks, 1980) and those with small troop sizes (Horrocks, 1986;Horrocks & Hunte, 1983), have shown a strong matrilineal influence on the relative rank of all individuals in the group.Juveniles enter the hierarchy low in relative rank but climb to positions appropriate to their matriline at approximately 4 years of age (Bramblett et al., 1982), suggesting that play behaviour may serve an immediate benefit to juveniles in establishing and testing dominance ranks in the period before the effects of the matriline are clear.However, within the Samara population, dominance ranks are shallower compared to other vervet populations and rank and grooming effects typically expected by a matrilineal hierarchy are absent (Henzi et al., 2013;Vilette et al., 2024).Samara vervets showed no preference to groom up the dominance hierarchy, showing no evidence that females attempt to maintain coalitions or grooming relationships with high-ranking females (Henzi et al., 2013).This might indicate that the Samara vervets do not sustain matrilineal kin association (Henzi et al., 2013;Vilette et al., 2024), or that the process of adult rank acquisition has changed (Horrocks & Hunte, 1983).Nevertheless, it appears that our study population does not follow typical matrilineal rank acquisition previously documented in vervets and, as such, play may be useful in continuously navigating its rank trajectory, as a juvenile cannot guarantee it will fall into rank under its mother.
Many species of nonhuman primates exhibit a similar preference to play more with siblings than with peers (Cheney, 1978;Fedigan, 1972;Owens, 1975;Southwick, 1965;van Lawick-Goodall, 1968).Most nonhuman primates grow up with older siblings and will typically acquire younger siblings at some point throughout their life (Suomi, 2014).The interactions and relationships between kin can comprise a major component of sociality and may serve fitness-related functions (Silk et al., 2006).For instance, the enduring relationships between kin may offer protection against conspecific competitors (Vilette et al., 2022), and directing play towards kin may serve as a tool to increase affiliation between kin.Although we found no evidence that vervets play preferentially with familiar social partners or individuals of close proximity, the low pd score for proximity (pd ¼ 61.91%) compared to grooming (pd ¼ 93.77%) may be related to grooming networks in our vervet population.Female juveniles are more integrated within grooming networks compared to males, and therefore females allocate more of their time to grooming, which may also partly explain why they are playing less frequently than males.
Additionally, these results raise several questions that may be used to guide future research on wild vervet monkey populations.It may be useful to first determine how wild vervets living in large groups acquire rank, as these large groups do not appear to sustain matrilineal kin associations given a lack of rank-related female relationships in the face of demographic stress (Henzi et al., 2013;Vilette et al., 2024).We could then compare the play partner preferences of the Samara population to a separate wild population of vervets that follows rank acquisition patterns typical of matrilineal vervet societies.Additionally, given that social structure (e.g.sex ratio, sibling availability, group size) and ecological conditions (e.g.drought) fluctuate temporally, quantifying such changes in analyses would prove useful.Multilayer network techniques (Finn et al., 2019) have become increasingly desirable, as they allow us to quantify how different interaction levels may influence one another.Bonnell et al. (2020) developed a dynamic multilayer approach that allows for the comparison of interaction layers over time.Such an analysis would be useful in play research, as it would allow us to quantitatively compare how different social layers (e.g.grooming interactions, dominance interactions, play bouts, kinship, spatial proximity, etc.) relate to one another over time and across varying levels of ecological stress.

Data Availability
Data supporting the findings of this study are available within the article and its supplementary materials.

Declaration of Interest
None.

Figure 1 .
Figure 1.Rates of play predicted by the model for each dyad sex category.Density plots present the range of play rates predicted by the model, with the height of the density curve indicating the probability of the predicted rate and the spread of the curve indicating its uncertainty (see Table3).

Figure 2 .
Figure 2. Rates of play predicted by the model for dyads with a low and high rank difference.Density plots present the range of play rates predicted by the model, with the height of the density curve indicating the probability of the predicted rate and the spread of the curve indicating its uncertainty (see Table3).

Figure 3 .
Figure 3.Rates of play predicted by the model for kin and nonkin.Density plots present the range of play rates predicted by the model, with the height of the density curve indicating the probability of the predicted rate and the spread of the curve indicating its uncertainty (see Table3).

Figure 4 .
Figure 4. Rates of play predicted by the model for dyads with low and high age difference.Density plots present the range of play rates predicted by the model, with the height of the density curve indicating the probability of the predicted rate and the spread of the curve indicating its uncertainty (see Table3).

Table 2
The annual percentage of scan samples allocated to each behaviour category for adults and juveniles of the Samara troops

Table 1
Composition of the Samara study troops over the study period