Degree of anisogamy is unrelated to the intensity of sexual selection

Males and females often display different behaviours and, in the context of reproduction, these behaviours are labelled sex roles. The Darwin–Bateman paradigm argues that the root of these differences is anisogamy (i.e., differences in size and/or function of gametes between the sexes) that leads to biased sexual selection, and sex differences in parental care and body size. This evolutionary cascade, however, is contentious since some of the underpinning assumptions have been questioned. Here we investigate the relationships between anisogamy, sexual size dimorphism, sex difference in parental care and intensity of sexual selection using phylogenetic comparative analyses of 64 species from a wide range of animal taxa. The results question the first step of the Darwin–Bateman paradigm, as the extent of anisogamy does not appear to predict the intensity of sexual selection. The only significant predictor of sexual selection is the relative inputs of males and females into the care of offspring. We propose that ecological factors, life-history and demography have more substantial impacts on contemporary sex roles than the differences of gametic investments between the sexes.


Materials and methods
Sexual selection indices. We extracted data on the intensity of selection in 66 species from Janicke et al. 23 .
An extensive review of the literature published since then revealed only one additional relevant study 43 . Hermaphrodite species were excluded, resulting in a dataset of 64 species. For the details of the literature review see the supplementary material.
Following Janicke et al. 23 three indices were used to estimate the intensity of selection: opportunity for selection (I) that is the standardized variance in reproductive success; opportunity for sexual selection (I s ) that is the standardized variance in mating success; and the Bateman gradient (β ss ) that is the slope of an ordinary least-squares regression of reproductive success on mating success. For statistical calculations, the effect sizes of these three indices (ΔI, ΔI s , and Δβ ss ) were used following Janicke et al. These indices represent bias in the intensity of selection between males and females, with positive values indicating more intense selection in males. ΔI and ΔI s are the coefficients of variation ratio "lnCVR", defined as the natural logarithm of the ratio between the coefficients of variation from males and females, and Δβ ss was calculated as Hedges' d as described in 44,45 . Except for the species that was added by us (Lamprotornis superbus), all selection indices were taken from Janicke et al.

Sexual dimorphism.
Sexual size dimorphism (SSD) was calculated as log[male size / female size], which provides a statistically appropriate measure of dimorphism 46 . Thus, positive values indicate species where males are larger than females. Size dimorphism was calculated using body weight or body length; the latter are related to each other as weight = length 3 * c, where c is a constant which is the same or very similar in both sexes. Therefore, we assume that male weight / female weight = (male length) 3 /(female length) 3 . If body size of males and females were available in different dimensions (for example male body length and female body weight), one body size measure was converted to the other. For 20 and 13 species we were not able to find reliable body size data for females or males, respectively, and therefore, SSD was not calculated for these species. The rationale for using male size relative to female size relies on studies that show intense sexual selection and polygamous mating is often, but not always, associated with male-biased sexual size dimorphism [47][48][49][50][51] . Also, recent studies show that sexual size dimorphism predicts sexual selection in a wide range of taxa 52,53 . Sexual size dimorphism is thus an objective and accessible measurement of the extent of sexual dimorphism and it is often used as an indicator of sexual dimorphism 53 . Furthermore, we kept SSD in our dataset to be consistent with the previous analyses of sexual selection by Janicke et al. 23 .  23 . The bias in gametic investment between males and females is assumed to lead to different reproductive rates of males and females, and hence to intense sexual selection among males (pathway 1), female-biased parenting (pathway 2), and elaborated trait expression in males (pathway 3). The figure was created using Keynote https:// www. apple. com/ uk/ keyno te/. . Therefore, a positive value indicates a male-biased gamete size (or gametic investment), a positive value indicates a female-biased gamete size (or gametic investment). To estimate female gamete size, neonate weight or egg weight (i.e., fertilized ovum, including the weight of the shell and nutrients) were used, both in grams. For five species, egg volume was estimated from egg weight, and for one species it was estimated from dry egg mass (see supplementary material). To estimate male gamete size, sperm volume was used (for further details on calculating sperm volume see the supplementary material). For 12 species we did not find data on sperm sizes, and we used corresponding data form closely related species 54 (see supplementary material). Clutch size (number of eggs per breeding event) and testis weight (in grams) were also collected from published studies (see supplementary material; for one species testis weight was based on testes width and length). Importantly, our dataset contains a wide range of species that include species with femalebiased and with male-biased gametic investment bias, also approximately isogamous species. For more details, see supplementary material and the descriptive statistics ( Figures S3 and S4).

Phylogeny.
To represent the phylogenetic relationships between species, the most recent comprehensive phylogeny was used from timetree.org 55 that included all but seven species in the recent dataset; the latter species were added manually to the phylogeny (see Fig. 2 and supplementary material).
Statistical methods. Phylogenetic Generalized Least Squares (PGLSs) were used to test the relationships between life history traits and sexual selection 56,57 . Based on the Darwin-Bateman paradigm, we assessed the relationships between (1) gamete size bias (or gametic investment) and sexual selection, (2) sexual selection and parental care, and (3) sexual selection and sexual size dimorphism. In addition, multi-predictor models were constructed in which one of three indices of sexual selection was the dependent variable, and gamete size bias (or gametic investment), parental care bias and sexual size dimorphism were the predictors. Pagel's lambda was estimated using maximum likelihood to control for varying levels of phylogenetic dependence 57,58 . Note that this model accommodates a range of patterns of phylogenetic dependence, for instance, patterns that are similar to the Ornstein-Uhlenbeck (OU) model 59 . Diagnostic plots were used to check the assumptions of the models. All analyses were conducted in R using the package caper 60 .

Data imputation.
To use as many species for the analyses as possible, we followed previous phylogenetic analyses and imputed some of the life history trait values that were not available in the literature 41 . The following data were missing: 31% of female body mass, 20% of male body mass, 48% of testis mass, 36% of sperm length, 31% of female gamete mass, 27% of clutch size, 14% of opportunity for selection, 7% of opportunity for sexual selection and 23% of Bateman gradient (Fig. 2). In multi-predictor models, casewise deletion of species with missing data is expected to reduce the power of analysis and potentially bias the results 41 . To pre-empt these potential caveats, missing data were estimated using multiple imputation (see below) and the PGLSs were carried out using the imputed datasets as well. Imputation is a method by which missing data is replaced with estimated data. It performs better when the proportion of missing data is low 61 . Therefore, life history trait data of additional species were collected to increase the accuracy of imputation. We carried out substantial search in the literature and found relevant data from 259 and 12,042 additional species (see Table S2 in supplementary material). Imputation was carried out using the relevant data we had for the 64 species as well as additional species (see the list of species and their data in supplementary material).
Multiple imputation was performed using a Brownian model and including Pagel's lambda to allow for varying amounts of phylogenetic dependence. We used the Rphylopars R package 62,63 . The imputed means and variances were used to generate random imputed values and repeated ten times, resulting in 10 complete datasets. Each dataset was then analysed separately, and the means and standard deviations of fitted model parameters were calculated. A leave-one-out cross validation reliability check was performed to test the accuracy of the imputation by deleting the non-missing data points one by one and re-running the imputation. The correlation between the original and imputed data points was investigated (see Table 2 in supplementary material). Note that our imputation produced reliable results since the correlation coefficients are 0.8 (or above, see supplementary material).

Results
We found no support for pathway 1 of the Darwin-Bateman paradigm (Fig. 1) since the relationships between gamete size bias, gametic investment and sexual selection are non-significant except a weak association between gamete size bias and ΔI s using the original data (Table 1), although the latter association is no longer significant using the complete dataset that also includes the imputed data (Table 1). Importantly, neither gamete size bias nor gametic investment was a significant predictor of sexual selection indices in multi-predictor models (see below, Table 2). In addition, pathway 3 was also not supported since there was no relationship between sexual selection bias and sexual size dimorphism using the original dataset or the complete dataset (see Tables 1, 2 Fig. 3). Thus if parenting shifts toward female-biased care this is associated with more intense selection on males (Fig. 3). Taken together, among the predicted relationships of the Darwin-Bateman paradigm, only pathway 2 (parental care-sexual selection) has consistent support, while pathway 1 (gametic Table 1. Phylogenetically corrected linear relationships (PGLSs) between components of the Darwin-Bateman paradigm using the original and the complete dataset; the latter also include imputed values. PGLS discard any species that has a missing value, therefore the sample size of the models using the original and complete dataset may differ. Parental care bias was calculated only for those species that exhibit some level of care (37 species). N refers to the number of species. For the complete dataset multiple imputation using lambda model was used to estimate missing data, 10 completed datasets were generated and PGLSs were carried out in each dataset (see "Materials and methods"). The mean and the standard deviation of the statistics are shown. For each model the response variables are in bold and the predictor variables are listed below. ΔI, ΔI s and Δβ ss refer to the opportunity for selection, the opportunity for sexual selection and the Bateman gradient, respectively. Statistically significant associations are in italics. www.nature.com/scientificreports/ investment-sexual selection) has weak if any support and pathway 3 (sexual dimorphism-sexual selection) is not supported.
To investigate potential co-effects among pathways, multi-predictor models were constructed in which one of three indices of sexual selection was the dependent variable, and gametic investment, parental care bias and sexual size dimorphism were the predictors. Following the results of the bivariate PGLSs, among the life history traits, only parental care bias predicted sexual selection (see Table 2). Parental care bias has a negative association with ΔI and ΔI s in the original dataset, whereas in the complete dataset care bias was significantly associated with all three indices (Fig. 3). R 2 values are moderate (0.11-0.54) suggesting that further unexplored life history and/ or ecological traits may also influence the strengths of sexual selection acting on males vs females.

Discussion
Our study provided two major results. First, we found that the degree of anisogamy does not predict the intensity of sexual selection (Table 1). Second, male-biased sexual selection is associated with female-biased parenting (Fig. 3, Table 1) supporting pathway 2 so that parental care is related to sexual selection.
The lack of relationship between gametic investment and sexual selection is consistent across three indices of sexual selection and two estimates of gametic investment (Table 1), and the results remain robust in multipredictor models which control for the effects of other predictors ( Table 2). Note that the imputed datasets that Table 2. Phylogenetically corrected multiple predictor relationships (PGLSs) between the bias of sexual selection and the bias of parental care, gametic investment and sexual size dimorphism, using the original and the completed datasets. For each model the response variables are in bold and the predictor variables are listed below. ΔI, ΔI s and Δβ ss refer to the opportunity for selection, the opportunity for sexual selection and the Bateman gradient, respectively. The response variable is ΔI (Models 1 & 2), ΔI s (Models 3 & 4), and Δβ ss (Models 5 & 6). Parental care bias was only calculated for species that exhibit some level of care, thus the models of complete dataset used 37 species. Statistically significant results are in italic. See Table 1  www.nature.com/scientificreports/ use sperm length from 489 species and testis mass from 259 species provide arguable better model estimates than the non-imputed datasets, and in all models R 2 values are close to zero (Table 1). We suggest three explanations for the lack of association between anisogamy and sexual selection. First, anisogamy-the existence of sexes due to different gamete sizes of males versus females that has been suggested to evolve repeatedly from isogamic ancestors 64 -may only provide the starting condition for evolving different male vs female strategies, and the relationship between sexual selection and anisogamy may be restricted to the early stages of anisogamy evolution 29 . Once anisogamy produced different sexes in dioecious organisms, then ecology, life history and demography could all impact the development of sex-specific strategies that increase male and female reproductive success 5,65 . To study these scenarios further studies are warranted that will explore the associations between ecological, life-history and demographic traits, sexual selection and gametic biases (see also 53 ).
Second, anisogamy is related to sexual selection at all anisogamy ratios, although the impacts of ecology and life history possibly override the effects emerging from anisogamy. For example, the harsh environment can tone down behavioural differences between the sexes by selection for monogamy and biparental care 66 , or can enhance the difference between sexes in sexual selection 67 although it may not directly impact male vs female gametic investments 68 . A demographic property of populations, the adult sex ratio (ASR), also impacts sex roles in birds, fish and humans since biased sex ratios induce changes in pair bonding and parenting 7,[69][70][71][72][73][74] . Consistently with phylogenetic studies, recently a long-term study showed fluctuations in ASR that impacted female choice in wild populations of Darwin's finches 75 . Experimental studies showed that female-biased ASR facilitates males to develop traits improving their competitiveness for accessing mate, and even the presence or absence of the female could trigger trait development 76 . Operational sex ratio (OSR; the ratio of sexually active males to females) may also have a pronounced effect on sex roles 77 , although OSR has been proposed to emerge from sex roles, rather than vice versa 4,6,66 .
Third, the estimates of gametic investment may not be perfect, making it difficult to detect relationships between anisogamy and sexual selection. Although comparing the energy contents of female vs male gametes would be ideal, these data do not yet exist across a broad range of taxa. A potential drawback of using gamete mass as an indicator of gametic investment is that the quality of ejaculate and female gamete (and therefore their size) may change with age and health conditions [78][79][80][81][82] . Furthermore, testis size varies over the breeding season 83,84 and clutch size may respond to environmental factors 35,[85][86][87] ; so that these variations could potentially mask the predicted associations between gametic bias and sexual selection.
Our second major result is that our analyses supported only one element of the Darwin-Bateman paradigm, that parental care is related to sexual selection, since male-biased sexual selection is associated with female-biased parenting (Fig. 3, Table 1). Consistently in multi-predictor models, the only variable that predicts sex difference in sexual selection is parental care bias ( Table 2). The significant relationship between parental care and sexual selection supports previous phylogenetic analyses 35,[87][88][89][90][91] and overall, it is consistent with the results of Janicke et al. 23 . We acknowledge that parental care and sexual selection may have a more complex relationship than usually thought 37,92 : parental care may impact the intensity of sexual selection as Trivers 93 originally proposed, although sexually selected traits and mating opportunities may also drive parenting decisions [94][95][96] . Note that we prefer to use the term parental care as opposed to parental investment, since the latter assumes that caring is costly; an assumption that is often not tested thoroughly 97 .  Table 1). Shaded data points represent overlapping data. No imputed data were used in the figure.
The figure was prepared using RStudio version 4.0.0 https:// www.R-proje ct. org/. www.nature.com/scientificreports/ The association between sexual selection and care is complex and may emerge via the ASR or the OSR. Parental care may be costly for the caring sex that leads to higher mortality, therefore sex-biased care may produce biased ASR that has a knock-on effect on intrasexual competition; however if the mortality is higher during mating competition than during care, the common sex is predicted to perform more care while the rare sex is expected to spend more time in the mating pool 4,98 . Alternatively, changes in OSR could impact sexual selection 99 : as the caring sex becomes a limited source for the other sex as a mating partner, biased care leads to biased OSR that influence the intensity of competition 100,101 . Since OSR and ASR are expected to be associated (although empirical studies do not always support this association 102 ) to distinguish the impacts of ASR from OSR will require further investigations.
We propose three potential reasons for the lack of correlation between sexual size dimorphism and sexual selection. First, sexual size dimorphism could be due to various selective processes, not only sexual selection, including natural selection and fecundity selection 50,51,103 . Also, sexual selection may not have a monotonous relationship with sexual size dimorphism (as we assumed in the current study), because if sexual selection is manifested via the agility of males 104 than intense sexual selection produces small males relative to females 105 . Furthermore, body size itself may not be under sexual selection 34,106 . Second, a possible explanation for the difference between Janicke et al. 23 who reported an association between sexual dimorphism and sexual selection and our result, is that Janicke et al. 23 used a composite scoring system for dimorphism that included ornamentation and behaviour, whereas we used only one variable, body size, to represent dimorphism (or the lack of dimorphism). It is conceivable that ornaments and behaviour (such as courtship) are more directly related to some forms of sexual selection than body size dimorphism 104,106 . Therefore, the combined effects of ornamentation and behaviour could produce an association with sexual selection whereas body size dimorphism alone may not. Third, sexual size dimorphism could be reduced in species that do not have parental care behaviour 88 . By using a continuous parental care bias variable, this may have led to restricting the analyses to species that exhibit care. As species that do not have parenting behaviour play a major role to create a correlation between parental care and sexual dimorphism, excluding them may have reduced the ability to spot the association between dimorphism and sexual selection.
Nevertheless, a limitation of our study as well as those of Janicke et al. 23 is that the 64 species used in these studies may not fully represent the diversity of gametic traits nor sexual selection exhibited by the huge diversity of extant terrestrial and aquatic animals. Therefore, future analyses of the Darwin-Bateman paradigm are warranted by using more detailed analyses of gametic investment, sexual selection and parenting in an ecological and demographic framework.
Ideally, life-history traits from high number of taxonomically diverse species should be used to test the Darwin-Bateman paradigm. Even though we aimed at using the best available data, we had numerous missing data therefore imputation was used instead of removing species with missing data to keep as many species in the analyses. Imputation becomes a frequently used tool in ecology and evolutionary biology as an alternative to removing species with missing observations [107][108][109][110][111][112][113] . The use and reliability of imputation in ecology and evolutionary biology are discussed in several articles 40,45,61,114,115 , and it is suggested as a reliable method for estimating missing data in ecological datasets 116,117 . As the analyses of the original and the imputed dataset gave similar results nearly all cases, we argue our results are reliable.
In conclusion, our study confirms and extends the findings of Janicke et al. 23 by showing parental care bias predicts sexual selection. However, the results do not support a key assumption of the Darwin-Bateman paradigm since the extent of anisogamy is unrelated to the intensity of sexual selection. We argue that ecological, life history and demographic variables could influence sex roles, and we call for new studies that integrate these processes into the investigation of sex role evolution.