Comparative Sperm Proteomics in Mouse Species with Divergent Mating Systems.

Sexual selection is the pervasive force underlying the dramatic divergence of sperm form and function. Although it has been demonstrated that testis gene expression evolves rapidly, exploration of the proteomic basis of sperm diversity is in its infancy. We have employed a whole-cell proteomics approach to characterize sperm divergence among closely related Mus species that experience different sperm competition regimes and exhibit pronounced variation in sperm energetics, motility and fertilization capacity. Interspecific comparisons revealed significant abundance differences amongst proteins involved in fertilization capacity, including those that govern sperm-zona pellucida interactions, axoneme components and metabolic proteins. Ancestral reconstruction of relative testis size suggests that the reduction of zona pellucida binding proteins and heavy-chain dyneins was associated with a relaxation in sperm competition in the M. musculus lineage. Additionally, the decreased reliance on ATP derived from glycolysis in high sperm competition species was reflected in abundance decreases in glycolytic proteins of the principle piece in M. spretus and M. spicilegus. Comparison of protein abundance and stage-specific testis expression revealed a significant correlation during spermatid development when dynamic morphological changes occur. Proteins underlying sperm diversification were also more likely to be subject to translational repression, suggesting that sperm composition is influenced by the evolution of translation control mechanisms. The identification of functionally coherent classes of proteins relating to sperm competition highlights the utility of evolutionary proteomic analyses and reveals that both intensified and relaxed sperm competition can have a pronounced impact on the molecular composition of the male gamete.

Introduction

The rapid evolution of male reproductive traits is widespread and particularly pronounced amongst components of the ejaculate, including sperm and seminal fluid proteins (Poiani 2006; Pitnick et al. 2008). This phenomenon can be attributed, at least in part, to sexual selection generated by competition between sperm of rival males to fertilize a limited supply of oocytes (Parker 1970). The evolution of larger testis, increased sperm production and sperm with enhanced motility and fertilization capacity is common in polyandrous species relative to less promiscuous taxa, where sperm competition is predicted to be less intense or absent (Gage 1994; Pizzari and Parker 2008; Soulsbury 2010; Firman et al. 2015). Importantly, the rapid evolution of ejaculate traits is mirrored by genotypic changes. Many genes expressed in the testis, including those encoding sperm proteins, and male accessory gland tissues evolve rapidly or display molecular signatures consistent with positive selection (Swanson and Vacquier 2002; Dean et al. 2008, 2009; Dorus et al. 2010; Vicens, Luke, et al. 2014). Although limited in number, reproductive genes have also been identified whose evolution correlates with predicted levels of postcopulatory sexual selection, an observation consistent with a role in processes governing competitive fertilization (Dorus et al. 2004; Herlyn and Zischler 2007; Ramm et al. 2009; Vicens et al. 2015). Despite these advances, mechanistic links between the molecular evolution of reproductive systems and specific phenotypes of relevance to competitive fertilization have yet to be widely established (but also see Greenspan and Clark 2011; Tollner et al. 2011; Yeh et al. 2012; Fisher et al. 2016).

The genetic basis of phenotypes relating to sperm competition, including gamete morphology, physiology and behavior, remain enigmatic in part because of their manifestation from a highly specialized and complex developmental process. Although analyses of the testis transcriptome should have the potential to link patterns of genome evolution to competitive fertilization phenotypes, an easily interpretable relationship has not been forthcoming. Nonetheless, the evolution of testis gene expression has been found to be rapid in relation to somatic tissues (Brawand et al. 2011) and correlated with predicted levels of sperm competition (Harrison et al. 2015). The testis transcriptome thus appears to be a responsive target for sexual selection, although the functional ramifications of testis transcriptome evolution remain largely unresolved. Proteomic analysis represents a powerful means of directly linking genomic and transcriptomic changes with cellular phenotypes. Despite the widespread application of proteomics to the characterization of sperm composition (Dorus et al. 2006, 2010; Martinez-Heredia et al. 2006; Baker et al. 2008a; 2008b; Dorus et al. 2010; Wasbrough et al. 2010; Chauvin et al. 2012; Skerget et al. 2013; Whittington et al. 2015; Bayram et al. 2016) limited progress has been made in understanding the proteomic basis of sperm diversification, in part because of an absence of targeted comparative analyses in species with divergent mating systems.

In this study, high-throughput proteomics were used to specifically study the sperm proteomes of three closely related mouse species (Mus musculus, M. spretus and M. spicilegus) that differ substantially in the strength of postcopulatory sexual selection (Gomendio et al. 2006). These species have been developed as a model for the study of sperm competition through detailed phenotypic analyses and exhibit pronounced differences in sperm production and quality (Gomez Montoto et al. 2011a). This includes higher levels of sperm viability, acrosome integrity, increased motility and ATP content in species experiencing enhanced sperm competition (Tourmente et al. 2015a). Importantly, enhanced motility in the high sperm competition species M. spretus and M. spicilegus is dependent on oxidative phosphorylation (OXPHOS), while species experiencing lower levels of sperm competition rely more heavily on ATP derived from glycolysis (Tourmente et al. 2015b). Additionally, asymmetries in fertilization success associated with levels of sperm competition have been observed in reciprocal heterospecific fertilization experiments: sperm from high sperm competition species have greater fertilization capacity than sperm from low sperm competition species, while ova from high sperm competition species exhibit stronger refractory properties to fertilization (Martin-Coello et al. 2009). The molecular basis of interspecific variation in sperm energetics and fertilization capacity has yet to be characterized.

The overarching goal of this study was to conduct the first targeted comparative analysis of sperm composition between closely related species with divergent mating systems. Based on existing functional and biochemical data, we predicted that divergence would be enriched in pathways associated with fertilization competency and energetics. We also explored patterns of molecular evolution associated with directional shifts in the intensity of sperm competition to determine if distinct changes might be related to the intensification or relaxation of sexual selection. Finally, our analyses, including the semi-quantitative comparison of over a thousand proteins, provides a more refined understanding of how gene regulation during spermatogenesis contributes to the composition and diversification of the fully differentiated male gamete.

Results

Mus Sperm Proteome Characterization

Tandem mass spectrometry (MS/MS) analysis of purified M. musculus sperm resulted in 50,904 high-quality spectral matches to a total of 1,244 proteins (supplementary table S1, Supplementary Material online). Of these, 99.4% were identified by multiple unique peptides, 70.3% were identified in more than one biological replicate and a survey of proteomic studies revealed that 68.3% have been previously identified as components of mouse sperm (supplementary table S2, Supplementary Material online) (Aitken et al. 2007; Baker et al. 2008a; 2008b; Dorus et al. 2010; Chauvin et al. 2012). In parallel, MS/MS characterization of sperm from M. spretus and M. spicilegus was conducted to investigate sperm proteome evolution across species experiencing different degrees of predicted sperm competition. Musmusculus has low levels of sperm competition, M. spretus has moderate-high levels of sperm competition, and M. spicilegus has high levels in relation to other Muroid rodents (Gomendio et al. 2006). MS/MS results for these species were comparable, if not slightly more comprehensive, to those for M. musculus with 54,436 and 62,783 high-quality spectral matches resulting in 1,272 and 1,357 identified proteins in M. spretus and M. spicilegus, respectively (supplementary tables S1 and S2, Supplementary Material online). 99.9% of these proteins were identified by multiple unique peptides and 75.0% were identified in more than one biological replicate. These searches were conducted against annotated M. musculus proteins in the absence of annotated genomes in the other species. Despite the low level of divergence between these taxa, even amongst sperm proteins (dN = 0.0039 between M. musculus and spretus), we sought to quantify the influence of divergence in interspecific MS/MS searches using a spectral matching algorithm that allows single amino acid substitutions. Only a single additional protein, ADAM2, was identified in M. spretus and no proteome differences were observed in M. spicilegus using this alternative approach. Divergence between these species thus has a minimal effect on protein identification.

A total of 1,012 proteins (63.1%) were identified in all three species with an additional 245 proteins (15.3%) identified in two of the three taxa analyzed (fig. 1). Over 10% of the proteins (109 out of 1012) common to all three species possess male fertility or sperm phenotypes, supporting the relevance of the proteome to sperm functionality (supplementary table S1, Supplementary Material online). Consistent with energetic demands of sperm motility, 67% of the most abundant proteins (top 5%) function in metabolism and energetics, including Gene Ontology Biological Process enrichments in the tricarboxylic acid cycle (P = 2.0E-13), acetyl-CoA biosynthesis (P =1.8E-7) and carbohydrate metabolism or glycolysis (P = 1.8E-5). Protein abundance estimates revealed highly comparable abundance distributions in each species (fig. 1), and a strong overall correspondence at the level of the protein (fig. 1). Lastly, hierarchical clustering analysis based on protein abundance revealed a closer relationship between the proteomes of M. spretus and M. spicilegus than either has with M. musculus (bootstrap probability = 0.948; Approximately Unbiased Tree Selection P value = 0.979; fig. 1). Thus, proteome-wide composition patterns reflect the established phylogenetic relationship between these species (Sarver et al. 2017).

F

Interspecific mouse sperm proteome comparisons. (A) Venn diagram displaying the overlap of sperm protein identification between species. Greater than 99% of proteins were identified by multiple unique peptides in each species and 78.3% of proteins were identified in at least two species. (B) Frequency distribution of normalized protein abundance estimates in each mouse species. (C) Ranked protein abundance estimates in Mus musculus and corresponding estimates in M. spretus and M. spicilegus. (D) Hierarchical clustering analysis conducted using protein abundance estimates in each species. Dendogram ordering, based on Euclidean distances, recapitulates the phylogenetic relationship of Sarver et al. (2017). Bootstrap probability is indicated and branch uncertainty was assessed with nonparametric resampling (10,000 iterations; bootstrap samples sizes ranged from 0.5 to 1.4). Heatmap ranges from low (yellow) to high (red) abundance.

Interspecific Protein Abundance Divergence

Protein abundance comparisons between species revealed significant correlations in all pairwise comparisons (R2 values range from 0.76 to 0.81; all P < 0.0001), with a marginally stronger correlation between the two more closely related high sperm competition species (supplementary fig. S3, Supplementary Material online). To systematically account for divergence between species we identified and corrected for all peptides absent in the M. spretus analysis that could be attributed to coding sequence divergence (see “Materials and Methods” section). This resulted in the removal of 3.2% of M. musculus peptides, corrected abundance estimates that were highly correlated with the uncorrected values (R2 = 0.99; P < 0.0001) and significant abundance changes for only 0.4% of the proteins analyzed (supplementary fig. S3, Supplementary Material online). Using these corrected values, a total of 125 (out of 3277 comparisons; 3.8%) interspecific abundance comparisons were found to be significant after Bonferroni correction for multiple testing, with a substantial number of proteins (33 out of 125; 26.4%) identified as significant in multiple pairwise species comparisons (supplemental table S2, Supplementary Material online). Several patterns are noteworthy. First, a significant enrichment in sperm phenotypes has been observed amongst mutants of these proteins (chi-square = 10.82; d.f. = 1; P = 0.001), including reduced sperm motility activation, asthenozoospermia and impaired acrosome reaction. Second, as would be expected based on the phylogenetic relationship between species, an increase (1.8-fold) in significant abundance differences was observed in comparisons between M. musculus and the two high sperm competition species relative to the comparison between these two more closely related species. Importantly, proteins identified as differentially abundant relative to M. musculus in both pairwise comparisons all exhibited a consistent directional change in abundance. Third, all seven identified sperm proteins involved in fertilization and zona pellucida binding were found to be significantly more abundant in one or both of the high sperm competition species relative to M. musculus. This included the sperm membrane proteins Zona Pellucida-Binding Protein 1-2 (ZPBP 1-2), Zonadhesin (ZAN), Zona Pellucida 3 Receptor (ZP3R) Cysteine-rich secretory GLI pathogenesis-related protein (GLIPR-like 2), Cd46 antigen and Parkinson-associated 7 (PARK 7). Finally, glycolytic proteins involved in sperm motility were observed amongst proteins of greater abundance in M. spretus relative to spicilegus, including Phosphoglycerate kinase (PGK2) and Glyceraldehyde-3-phosphate dehydrogenase spermatogenic (GAPDHS). This observation is consistent with the increased reliance of M. spicilegus on OXPHOS relative to M. spretus (Tourmente et al. 2015b).

Ancestral Reconstruction of Relative Testis Size

To contextualize our proteomic results within the evolutionary history of sperm competition in mice, we reconstructed ancestral body and testis sizes using available phenotypic data and estimated relative testis size across the clade. Relative testis size, the ratio of observed testis size to expected testis size based on body mass, is commonly used as a proxy for sperm competition in Muroid rodents (Firman and Simmons 2008; Gomez Montoto et al. 2011b). If the relative testis size of the most recent common ancestor (MRCA) of our study species was low then increases in predicted levels of sperm competition would be concentrated on the lineages leading to both M. spretus and M. spicilegus. Alternatively, if the relative testis size of the MRCA was intermediate in relation to our study species, then increases and decreases in sperm competition would be predicted to have occurred in M. spicilegus and M. musculus, respectively. Two ancestral state reconstruction approaches were employed (see “Materials and Methods” section), resulting in highly comparable testis size estimates (R2 = 0.995; P < 0.0001) (supplementary table S2, Supplementary Material online). This analysis yielded a relative testes size estimate in the MRCA (0.914) that was intermediate amongst our study species and comparable to that of M. spretus (1.07) (fig. 2). Thus, the lineages leading to M. spretus exhibit a relatively small increase in relative testis size (17.1%), while marked decreases (83.8%) and increases (55.5%) have occurred in the lineages leading to the low (M. musculus) and high (M. spicilegus) sperm competition species analyzed in this study. Our study species thus provide the opportunity to investigate sperm evolution in the context of both enhanced and relaxed regimes of sperm competition.

F

Relative testis size ancestral state reconstruction. Body and testis ancestral states were reconstructed independently to calculate ancestral relative testis size estimates. Lineages associated with the species analyzed in the current study are highlighted in red. Phylogenetic relationships and branch lengths (scaled to nucleotide substitutions) were based on Sarver et al. (2017).

Evolution of Sperm Competition and Proteome Divergence

An unsupervised soft clustering algorithm was implemented to identify emergent patterns in protein abundance variation and relate them to the evolutionary history of sperm competition. This analysis identified 6 distinct clusters that maximized average centroid distance between profiles (supplementary fig. S4, Supplementary Material online). Cluster membership deviated from a uniform distribution with an increase in protein membership in Clusters 1 (+48.0%) and 2 (+13.3%). These clusters include proteins that have relative abundance increases in either M. spicilegus (high sperm competition) or decreases in M. musculus (low sperm competition) (fig. 3). Although it would be premature to draw firm conclusions based on these patterns alone, more widespread changes in protein abundance do appear to be associated with recent shifts in sperm competition intensity identified by our ancestral state reconstruction. It is also noteworthy that clusters 5 (−16.7%) and 6 (−28.3%), which represent increased or decreased protein abundance specifically in M. spretus, were amongst the lowest in membership size (supplementary fig. S4, Supplementary Material online).

F

Protein abundance variation in relation to sperm competition. (A) Soft-clustering analyses were conducted to identify global patterns of protein abundance differences that are associated with interspecific variation in sperm competition. Protein abundance patterns relative to the cluster average are depicted for each protein (purple: high identity; green: lower identity). Deviation in cluster size relative to uniform membership across clusters is indicated; 60 proteins per cluster were expected based on a uniform distribution. (B) Fold enrichment of proteins belonging to significant Gene Ontology Categories relating to sperm head and acrosome, metabolism and flagellum. Significant enrichments are indicated by an asterisk for Clusters 1–4. Enrichment calculations were based on the expected frequency of proteins belonging to a given Gene Ontology category within the sperm proteome as a whole.

Gene Ontology analyses were performed to gain insights into the functional characteristics of the clustered proteins. This revealed distinct functional enrichments across Clusters 1–4, which are informative to predicted shifts in sperm competition. First, 11 proteins associated with impaired male fertility phenotypes were present in Cluster 1, including several which may contribute to the unique characteristics of sperm in this species. Three of these proteins are required for proper fertilization, including ATPase, class I, type 8B, member 3 (ATP8B3). ATP8B3 is an aminophospholid transporter required for lipid distribution in the acrosomal bilayer and for zona pellucida binding and penetration (Wang et al. 2004). Several proteins were also identified that are required for normal motile cilium physiology, included sperm associated antigen 17 (SPAG17), a protein indispensable to motile cilia function (Teves et al. 2013), armadillo repeat containing 4 (ARMC4) and ropporin 1-like (RPOPN1L). The latter two proteins are implicated in ciliary dyskinesia and Rpopn1l knockouts result in impaired motility, instability of the principle piece and infertility when combined in a double knock-out with ropporin 1 (Fiedler et al. 2013). Second, consistent with fertilization asymmetries between these species, a significant enrichment of acrosomal vesicle and zona pellucida complex proteins were observed in Cluster 2, including all 6 proteins involved in binding of sperm to the zona pellucida (fig. 3). This finding is consistent with our pairwise analyses between species, and based on our ancestral state reconstruction, this reduction occurred concurrently with the decrease of sperm competition in the lineage leading to M. musculus. Third, a unique enrichment in Cluster 1 and 2 was observed for inner dynein arm proteins, central axoneme apparatus and those contributing to microtubule motor activity. Those belonging to Cluster 1, and exhibiting a specific increase in abundance in M. spicilegus, included dynein light chain LC8 type 1 and 2 (DYNLL1 and 2), dynein light chain Tctex-type 1F (DYNLT1F), tubulin-beta 4b (TUBB4B) and sperm associated antigen 17 (SPAG17). In contrast to Cluster 1, the four belonging to Cluster 2 are exclusively dynein heavy chains, which possess the force-generating ATPase activity responsible for movement along the microtubule. This includes dynein heavy chains 1 (DYNC1H1), which is an axonemal dynein associated with asthenozoospermia and reduced ciliary beat frequency when mutated (Neesen et al. 2001). Thus, abundance changes in light and heavy chain dyneins have occurred on distinct lineages, corresponding to increases and decreases in sperm competition, respectively. Forth, an enrichment of principle piece and fibrous sheath proteins were observed in Clusters 3 and 4, which represents decreased abundance in M. spicilegus or both of the high sperm competition species relative to M. musculus (fig. 3). This observation is consistent with the results of our pairwise analyses and included a substantial number of glycolytic proteins, such as glyceraldehyde-3-phosphate dehydrogenase spermatogenic (GAPDHS), phosphofructokinase (PKM), aldolase A (ALDOA), pyruvate kinase (PK), dihydrolipoamide S-acetyltransferase (DLAT) and hexokinase 1 (HK1). These proteins also contribute to the overall enrichment of glycolytic proteins in Clusters 3 and 4 (fig. 3). This observation is generally consistent with the reduced reliance on glycolytically derived ATP in M. spretus and M.spicilegus (Tourmente et al. 2015b). However, broad patterns relating to OXPHOS dependency in high sperm competition species were not observed. The sole exception was two ATP synthase, H+ transporting mitochondrial F1 complex proteins (ATP5A) and valosin containing protein (VCP), a positive regulator of OXPHOS, which contributed to a significant enrichment of ATP metabolic proteins in Cluster 2. Finally, we note that with the exception of an increase in proteasome regulatory particles in M. spretus no distinct functional enrichments were observed in Clusters 5 or 6.

Post-Meiotic Testis Transcriptome Resembles Mature Sperm Proteome

Transcript abundance is one of the primary determinants of protein abundance (Vogel and Marcotte 2008) but estimates of its contribution to proteome composition vary substantially within and between analyses of various tissues and cell lines (Vogel et al. 2010; Schwanhausser et al. 2011; Bauernfeind et al. 2015; Uebbing et al. 2015) and it is poorly studied within developmental contexts. Spermatogenesis is also unique amongst cellular development because of nuclear condensation, global post-meiotic chromatin remodeling and the associated reduction in transcriptional activity. Importantly, this reduction is contemporaneous with a suite of dynamic morphological changes that occur as spermatids differentiate into elongated spermatids. To explore this in the context of spermatogenesis, protein abundance measures for all proteins identified in M. musculus were compared with gene expression estimates across four successive spermatogenesis stages (Soumillon et al. 2013). These stages included (1) spermatogonia, unspecialized germ cells undergoing mitosis to produce primary spermatocytes, (2) spermatocytes, gametocytes undergoing meiotic division, (3) spermatids, haploid male gametids maturing via individualization, flagellum elongation and acrosome formation and (4) mature spermatozoa, elongated gametes with condensed nuclei and chromatin structure. This analysis revealed a significant correlation between gene expression and protein abundance in spermatocyte (R2 = 0.061; P = 0.039) and spermatid (R2 = 0.074; P = 0.0125) developmental stages (fig. 4). Thus, significant abundance correlations between the proteome of mature sperm and the transcriptomes of early, post-meiotic spermatogenesis cells (but not other developmental time points) indicate that those cellular stages contribute most to the final sperm proteome. Importantly, this correspondence coincides with prominent late spermatogenic stages, including early spermiogenesis, that encompass axoneme elongation and acrosome development. As would therefore be predicted, genes specifically up-regulated in spermatids relative to other developmental time points share functional coherence with those independently identified amongst the most abundant sperm proteins (top 20%). This includes functional enrichments in GO reproductive processes, including reproduction (P = 5.60E-16), spermatogenesis (P = 5.85E-14), sperm motility (P = 1.6E-07), fertilization (P = 1.1E-06), sperm–egg recognition (P = 5.79E-6), spermatid development (P = 1.3E-04) and acrosome assembly (P = 9.1E-03). We also note that these correlations are low, albeit significant, particularly compared with transcriptome-to-proteome correspondence in steady-state cell culture systems (Schwanhausser, et al. 2011). This is consistent with the prominent role of translational control in spermatogenesis, particularly during late spermiogenesis when transcriptional activity has been curtailed. Finally, the absence of a relationship in spermatozoa confirms that RNA in mature sperm is quite distinct from those that contribute to the composition of the mature proteome (Fischer et al. 2012). As mature spermatozoa are transcriptionally inert (Hecht 1998), the transcripts identified at this stage represent RNA that may be actively retained during development to be delivered to the oocyte following fertilization (Sendler et al. 2013).

F

Spermatogenic gene expression and sperm proteome composition. (A) Linear regression analyses of spermatozoa protein abundance estimates in Mus musculus with gene expression levels across successive stages of spermatogenesis. Stage-specific expression values were obtained from Soumillon et al. (2013). All proteins are included which have corresponding gene expression data. (B) Percentage of genes subject to translational repression at the spermatogonia–spermatocyte and spermatocyte–spermatid transition based on Gan et al. (2013). Proteins exhibiting significant abundance differences between low (M. musculus) and high (either M. spretus or M. spicilegus) sperm competition species (orange), proteins exhibiting comparable abundance between low and high sperm competition species (grey) and testis expressed genes identified by Gan et al. (2013) that were not identified by MS/MS analysis in mature sperm (black) are displayed. (C) Stage-specific relative protein abundance relative to spermatogonia, based on Gan et al. (2013), for proteins exhibiting significant abundance differences between low (M. musculus) and high (either M. spretus or M. spicilegus) sperm competition species (orange), the remainder of the sperm proteome (grey) and other testis expressed genes not encoding protein products in mature sperm (black).

Translational Repression and Sperm Proteome Divergence

To ensure appropriate amounts of mRNA are available for the latter stages in spermatogenesis when progressive genome silencing occurs, many transcripts are produced at earlier stages and remain in a translationally repressed state (Kleene 1993; Lee et al. 1995). The enhanced transcriptome-to-proteome correspondence during morphological differentiation of elongated spermatids supports the hypothesis that genes encoding sperm components may be targets of translational repression. To assess the impact of translational repression on spermatogenesis and test this hypothesis, we compared the frequency of translational repression amongst sperm proteome genes (e.g. genes encoding integral sperm proteins) relative to other testis expressed genes using the data and criteria established by Gan et al. (2013). This revealed a significant enrichment of sperm proteome genes experiencing translational repression across the spermatogonia–spermatocyte transition (Fisher’s exact test; d.f. = 1; P < 0.0001) (fig. 4). Additionally, when this analysis was restricted to sperm proteome genes encoding proteins exhibiting interspecific abundance differences between low (M. musculus) and high (either M. spretus or spicilegus) sperm competition species, a significant enrichment was observed relative to the remainder of the proteome and other testis expressed genes across both the spermatogonia–spermatocyte (Fisher’s exact test; d.f. = 1; P = 0.031) and spermatocyte–spermatid (P = 0.045) transitions. Thus, post-transcriptional regulation disproportionally targets sperm genes and, amongst those, genes encoding proteins that exhibit abundance differences between species. An alternative indirect means of assessing the extent of translational repression is to examine the relative increase in protein abundance across spermatogenesis. As might be expected, testis expressed genes not contributing to the sperm proteome exhibit no increase over time, while sperm genes exhibit a substantial increase (fig. 4). Importantly, proteins exhibiting interspecific abundance variation experience an even greater temporal increase in protein abundance throughout spermatogenesis, indicative of the influence of translational repression. Although increased post-meiotic expression for some of these genes cannot be ruled out, this seems unlikely for sperm genes as a whole because of the well documented changes in chromatin across the genome and the evidence for translational repression presented in figure 4. Lastly, genes experiencing translational repression were found to be enriched for GO reproductive functions, including sperm motility (P = 2.59E-08), spermatid development/differentiation (P = 2.79E-04) and fertilization (P = 1.70E-02). This included Zona Pellucida Binding Protein (ZPBP), Zona Pellucida Receptor 3 (ZP3R) and the most abundant protein in the sperm proteome, Sperm Acrosome Associated 1 protein (SPACA1). Thus, translational repression is prominently involved in the regulation of genes encoding sperm components (relative to other testis-expressed genes), impacts genes mediating flagellum and acrosome development that contribute to sperm evolution and is therefore a regulatory mechanism that may be responsive to selection associated with sperm competition.

Discussion

To advance our understanding of sperm diversification and its relationship to regulatory mechanisms at work during spermatogenesis, we have conducted the first comparative analyses of sperm composition between closely related species experiencing different levels of sexual selection. The choice of species also allowed us to take advantage of the wealth of data relating to mouse spermatogenesis, including stage-specific gene expression data and fertility phenotypes. From a methodological standpoint, our study demonstrates that comparative proteomics can be conducted in a thorough and statistically robust manner between closely related species, as long as one of the species has a well-annotated genome. The sperm proteomes generated here were highly comparable and their composition was largely insensitive to searches allowing for amino acid substitutions, despite the fact that all searches were conducted against the M. musculus proteome. The availability of the M. spretus genome, albeit unannotated, also provided the opportunity to directly compare protein abundance estimates between the full set and filtered set of spectra after the removal of peptides without complete identity between M. spretus and M. musculus. Again, the impact of divergence was found to be minimal. Our observations complement an ambitious recent study by Bayram et al. (2016), which compared sperm proteomes across rodents and ungulates. Their analyses established that a “core” sperm proteome can be successfully identified across mammalian orders using a limited number of annotated genomes. As expected, this “core” proteome was largely comprised of highly conserved proteins, while rapidly evolving proteins were more commonly identified in an order- or species-specific manner. It is noteworthy that the total number of sperm proteins identified in their analysis of six rodent species was comparable to the number of proteins we have identified. Bayram et al. (2016) and the analyses presented here address methodological issues associated with opposite ends of the spectrum of molecular divergence present within mammals and, in conjunction, demonstrate the feasibility of and potential for extending comparative proteomics approaches beyond taxa with annotated genomes.

Although the species examined here have diverged substantially in their reproductive biology and experienced both marked increases and decreases in sperm competition, we nevertheless urge caution in interpreting our results as they are based on a limited number of species. Comparative analyses will be required in additional species, representing independent transitions in mating system and sperm competition regimes to assess the generality of our observations. Sperm competition is also only one of many possible influences on sperm evolution and additional neutral and selective processes will need to be considered to fully understand the proteomic variation delineated here. Nonetheless, our analyses revealed associations between the abundance of several functionally coherent groups of proteins and evolutionary transitions in the predicted strength of selection associated with sperm competition. Most importantly, many of these proteomic differences are directly relevant to well-characterized sperm competition phenotypes that distinguish these species, including gamete interactions, metabolism and motility (Martin-Coello et al. 2009; Gomez Montoto et al. 2011b; Tourmente et al. 2011, 2013; Tourmente et al. 2015a; Tourmente et al. 2015b).

Various lines of evidence support the influence of sperm competition on the evolution of gametic interactions, favoring sperm that are superior at binding and penetrating oocytes, and in turn, promoting the evolution of more effective blocks to polyspermy (Frank 2000; Levitan et al. 2007; Snook et al. 2011). Previous studies of the species examined here have revealed substantial variation in sperm competitive ability and heterospecific fertilization asymmetries that are consistent with interspecific differences in the intensity of sperm competition (Martin-Coello et al. 2009). Our analysis revealed a significant and systematic increase in the relative abundance of sperm membrane and acrosomal proteins, including those that govern sperm–egg interactions, in high sperm competition species relative to M. musculus. This includes proteins that are known to mediate sperm-zona pellucida interactions (ZAN, ZP3R, ZPBP1 and 2), acrosome biogenesis (MAN2B2), acrosomal dynamics (CD46 and ZPBP2) and sperm–egg fusion (ADAM3). As might be expected given the codiversification of male and female interacting partners (Clark et al. 2009; Claw et al. 2014; Vicens and Roldan 2014), several of these proteins (ZAN, ZPBP1 and 2) also exhibit patterns of evolution consistent with the influence of positive selection. It is therefore plausible that fertilization efficiency could be modulated through changes in expression and developmental incorporation of sperm components responsible for gamete interactions, in conjunction with adaptive molecular evolution of these proteins. It is also noteworthy that significant abundance differences were not observed between M.spretus and M. spicilegus, despite the dramatic increase in relative testis size in M. spicilegus. It will therefore be important to determine if differences in sperm–egg binding characteristics between these species involves the evolution of proteins governing these interactions (rather than their abundance) and conduct complementary comparative studies of the oocyte proteome to fully understand the molecular diversification of gamete interactions across these species.

As was expected based on interspecific variation in sperm motility and velocity, a suite of differences was observed in proteins localizing to the flagellum that, in some cases, are known to contribute to ciliary physiology and motility. Perhaps most striking were the compartmentalized changes in axonemal dynein constituents, including a decreased abundance of dynein heavy chain proteins in M. musculus and an increased abundance of dynein light chains in M. spicilegus. Although it would be premature to speculate about these distinct directional changes, knock-outs of both heavy and light dynein chains result in impaired flagellar motility (Neesen et al. 2001; Rashid et al. 2010). Additionally, the observed increase in numerous other proteins in M. spicilegus required for proper ciliary behavior is consistent with enhanced flagellar function and may ultimately relate to motility phenotypes specific to this species. Interestingly, principle piece proteins that were decreased in abundance in M. spicilegus and M. spretus were biased towards key glycolytic enzymes rather than structural proteins. This is generally consistent with variation in the relative metabolic reliance on glycolysis and OXPHOS amongst these species, including significant deficits in sperm performance following the inhibition of mitochondrial respiration in both high sperm competition species (Tourmente et al. 2015a).

A somewhat unexpected facet of our results was the observation that distinct protein repertoires experienced abundance changes on lineages associated with either an increase or decrease in sperm competition. While it is well-established in mice that sperm fertilization capacity increases rapidly in response to increased sperm competition (Firman et al. 2014) and that rates of sperm diversification are correlated with sperm competition intensity more generally (Rowe et al. 2015), the ramifications of relaxed sperm competition has received far less attention. The evolution of sperm in the absence of sperm competition has been discussed in the context of a degenerative process (van der Horst and Maree 2014) and it is certainly tempting to speculate that the decreased abundance of a coherent set of proteins critical to sperm function (e.g. egg binding and microtubule motors) may be related to deficits in M. musculus sperm quality, including decreases in acrosome integrity, proportion of motile sperm and velocity. However, in the absence of molecular links to precise structure–function relationships it is impossible to distinguish between degenerative outcomes and those that reflect selective trade-offs during spermatogenesis. Nonetheless, our results indicate that a recent relaxation of sexual selection has resulted in pronounced changes at the molecular level.

Spermatozoa morphogenesis is regulated by a remarkably large repertoire of testis-specific genes and complex patterns of gene expression (Ramskold et al. 2009). This complexity can be attributed to the dynamic processes the genome experiences during gamete differentiation, including meiosis, meiotic sex chromosome inactivation and the histone-to-protamine transition, and perhaps explains the paucity of knowledge about sperm genotype–phenotype relationships. Our quantitative analysis revealed the strongest correlation, albeit it still rather weak, with gene expression in spermatids, and to a lesser extent in spermatocytes. Thus, early, post-meiotic spermatogenesis cells make the largest contribution to the final sperm proteome and little, if any, relationship is observed at early time points or in transcriptionally quiescent mature spermatozoa. Although there is currently a tremendous amount of variation in estimates of how mRNA levels relate to protein abundance (Bauernfeind et al. 2015; Uebbing et al. 2015), steady-state analyses in cultured cells indicate that transcript abundance explains ∼25–40% of protein variation and that translational control plays a predominant role (Vogel et al. 2010; Schwanhausser et al. 2011). Our analysis, which represents one of the first in a developmental context, indicates a weaker relationship throughout spermatogenesis and is therefore consistent with a prominent contribution of post-transcriptional regulation to sperm development. However, we cannot rule out the possibility that our correlations might improve if transcriptome and proteomic analyses were conducted in a single male (or the same group of males); little is current known about variation in testis gene expression or sperm proteomes amongst fertile males. The influential nature of translational regulation in sperm development and evolutionary diversification (see below) is consistent with the observation that the testis is a highly permissive expressional environment (Kaessmann 2010) and that newly created and de novo genes frequently acquire robust testis expression (Kaessmann 2010; Zhao et al. 2014). As such, analysis of testis expression evolution alone may be insufficient to delineate functionally meaningful changes in relation to sperm form and function and future studies are likely to require integrated analyses of stage-specific transcriptomic and proteomic datasets.

Translational repression is a conserved feature of mammalian and insect spermatogenesis (Kleene 1996; Braun 1998) that, based on our analysis, has an enriched impact amongst sperm genes, and in particular those encoding proteins exhibiting evolutionary diversification. In light of the heightened importance of translational control in post-meiotic cells experiencing genome silencing and its direct relationship to developmental outcomes in spermiogenesis, translational repression may represent a regulatory mechanism that is particularly responsive to selection. The potential for translational repression to influence developmental variation by modulating the contribution of early transcription to transcript abundance at later stages, when the transcriptome is more directly aligned with proteome composition, suggests that further characterization of translational mediators in the testis and the evolution of their mRNA target repertoires may be particularly informative. Whereas the role of translational repression has been established for genes encoding structural proteins of the flagellum (Delbes et al. 2012), our analysis suggests that this mechanism also influences acrosomal proteins that are likely candidates involved in fertilization asymmetries between these species. More broadly, we speculate that translational repression is mechanistically associated with the need to compensate for the widespread transcriptional quiescence that results from genome repackaging and nuclear condensation during post-meiotic spermatogenesis. Little is currently known about the temporal pattern of protamine deposition across the genome, but one would expect that genes in regions experiencing early silencing may acquire higher levels of pre-meiotic expression accompanied by robust translational repression. Although firm conclusions will depend on a more comprehensive understanding of translational regulation in the context of histone-to-protamine repackaging, protamine incorporation does vary substantially across mammals (Corzett et al. 2002; Wykes and Krawetz 2003) and is, as of yet, a poorly understood parameter in the evolution of male germline gene regulation, spermatogenesis and sperm phenotypes.

Materials and Methods

Sperm Isolation and Purification

Adult males (20–25 weeks old) of three mouse species (Mus musculus, Mus spretus and Mus spicilegus) from wild-derived colonies were sacrificed by cervical dislocation. Mature sperm were purified from the caudal epididymis using a swim up approach in Dubecco phosphate-buffered saline (DPBS) medium at 37°C for 10 min. Isolated sperm were washed 3 times using 900 µl DPBS, followed by centrifugation at 3,000 g for 1 min. Sperm were then pooled from multiple randomly chosen individuals to obtain sufficient cell numbers for replicated (3×) protein isolation and mass spectrometry analyses in each species (see below).

Protein Quantitation and 1-D Gel Separation

Sperm pellets were lysed in Laemmli Buffer and quantified using Quick Star Bradford Protein Assay Kit 1 (Bio-Rad) according to the manufacturer’s protocol. Three replicate gels were run, each containing a single replicate from each of the three species. 25 µg of protein from each sample was separated on 4–12% Nupage Bis–Tris gels running on an XCell SureLock Mini-Cell PowerEase 200 system (Life-Technologies). Gels were then fixed in 45% methanol and 1.0% acetic acid for 1 h and immediately stained with Coomassie (0.1%w/v Coomassie, 34% methanol, 17%w/v ammonium sulfate, and 0.5% acetic acid) until protein bands were visible. Gels were transferred to a gel slicer where each lane was cut vertically and then horizontally into 16 slices. Gel slices were subjected to proteolytic digestion using an Automated Preparation Station (Perkin Elmer). Proteins were first reduced and alkylated and then successively treated with dithiothreitol and iodoacetamide followed by digestion with trypsin (500 ng/µl, Promega) at 37°C for 16 h.

Tandem Mass Spectrometry

All LC-MS/MS experiments were performed using a Dionex Ultimate 3000 RSLC nanoUPLC (Thermo Fisher Scientific Inc, Waltham, MA) system and a QExactive Orbitrap mass spectrometer (Thermo Fisher Scientific Inc, Waltham, MA). Separation of peptides was performed by reverse-phase chromatography at a flow rate of 300 nl/min and a Thermo Scientific reverse-phase nano Easy-spray column (Thermo Scientific PepMap C18, 2 µm particle size, 100 A pore size, 75 mm i.d.×50 cm length). Peptides were loaded onto a pre-column (Thermo Scientific PepMap 100 C18, 5 µm particle size, 100 A pore size, 300 mm i.d.×5 mm length) from the Ultimate 3000 autosampler with 0.1% formic acid for 3 min at a flow rate of 10 µl/min. After this period, the column valve was switched to allow elution of peptides from the pre-column onto the analytical column. Solvent A was water plus 0.1% formic acid and solvent B was 80% acetonitrile, 20% water plus 0.1% formic acid. The linear gradient employed was 2–40% B in 30 min. The LC eluant was sprayed into the mass spectrometer by means of an Easy-spray source (Thermo Fisher Scientific Inc.). All m/z values of eluting ions were measured in an Orbitrap mass analyzer, set at a resolution of 70000. Data dependent scans (Top 20) were employed to automatically isolate and generate fragment ions by higher energy collisional dissociation (HCD) in the quadrupole mass analyser and measurement of the resulting fragment ions was performed in the Orbitrap analyser, set at a resolution of 17500. Peptide ions with charge states of 2+ and above were selected for fragmentation. The mass spectrometry proteomics data have been deposited to the ProteomeXchange Consortium (http://proteomecentral.proteomexchange.org) via the PRIDE partner repository with the dataset identifier PXD003376.

Experimental Design and Statistical Rationale

To improve the power of detection and proteome coverage, sperm from randomly chosen males were pooled to produce three pooled replicates (25 µg each) for each species; 3 males were required for M. spretus and M. spicilegus per replicate, while as many as 15 males were required for Mus musculus. Sixteen gel slices were analyzed per replicate, resulting in 48 MS/MS runs per species. Peak lists generated using Proteome Discoverer (version 1.3.,ThermoFisher) and data from each MS/MS sequential run was analyzed by X!Tandem CYCLONE TPP (2011.12.01.1 – LabKey, Insilicos, ISB) (Craig and Beavis 2004) against the Global Proteome Machine mouse proteome database (GRCm38 refseq; 29,763 entries) with a fragment ion mass tolerance of 0.40 Da and a parent ion tolerance of 25.0 PPM. Iodoacetamide derivative of cysteine was specified as a fixed modification, whereas oxidation of methionine was specified as a variable modification. Two missed trypsin cleavages were allowed and nonspecific cleavages were excluded from the analysis. False Discovery Rates (FDRs) were estimated using a randomized decoy database. Peptide identifications were accepted if they could be established at >95.0% probability as specified by the Peptide Prophet algorithm and protein assignments were accepted if they could be established at >99.0% probability as specified by the Protein Prophet algorithm (Keller et al. 2005). Conservative estimates, based on the single highest protein and peptide FDRs, indicate that at most 0.016% of proteins and 2.3% of peptides would be incorrectly included as true identifications under these thresholds. This is, however, a substantial overestimate as the actual FDR distributions of the dataset are heavily skewed away from these lower inclusion thresholds. Proteins that contained identical peptides that could not be differentiated based on MS/MS analysis alone were grouped to satisfy the principles of parsimony. Additional statistical methodology and justification relating to protein quantitation are provided below.

Proteomic and Comparative Genomic Adjustments for Protein Divergence

Data obtained in each of the three species were searched against the M. musculus protein annotation in the absence of annotated genomes in the other two species. Two complementary approaches were utilized to account for potential divergence in M. spretus and M. spicilegus proteins, although we note that levels of protein divergence are extremely low between these closely related species (see “Results” section). First, X!Tandem was run using the parameters previously described with the addition of peptide identification with single substitutions, thus allowing the identification of peptides that differ from the M. musculus sequence by a single amino acid. This analysis led to the identification of only one additional protein (see “Results” section). Second, to systematically account for divergence between species we identified and corrected for all peptides absent in the M. spretus analysis that could be attributed to coding sequence divergence using the unannotated M. spretus genome sequence (Ref# ERS138732, http://www.ebi.ac.uk/). To this end, all peptides identified in M. musculus were queried against peptides identified in the M. spretus samples and those without matches were then queried against the M. spretus genome using a custom Perl script. Parent proteins of the unmatched peptides were queried against the M. spretus genome using tBLASTn with an E value cutoff of 10−3. Peptides without 100% identity and length matches were removed from PeptideProphet results, and quantitative proteomic analysis were completed as described below using only the set of identical peptides which could be theoretically identified in samples from both species. This resulted in the removal of 3.2% of M. musculus peptides. Corrected abundance estimates were highly correlated with the uncorrected values (r2 = 0.990) and resulted in significant changes in abundance in only 0.4% of the proteins analyzed (supplementary fig. S3, Supplementary Material online). Thus, both approaches revealed that divergence between these closely related species had a marginal effect on the overall composition of the proteome and protein abundance estimates. We note that the second correction approach could not be applied to M. spicilegus as the genome is not yet sequenced.

APEX Protein Quantitation and Analysis

The quantitation of proteins was conducted using the APEX Quantitative Proteomics Tool (Lu et al. 2007; Braisted et al. 2008). 50 proteins with the highest numbers of spectral counts, and peptide and protein identification probabilities in the M. musculus sperm were utilized for the training dataset (Skerget et al. 2015). 35 physicochemical properties available in the APEX tool were used for prediction of peptide detection/nondetection and protein detection probabilities (Oi) were computed using the Random Forest classifier algorithm. APEX protein abundances were calculated using a merged protXML file generated by the ProteinProphet algorithm. Modified peptides were not included and semi-tryptic and shared peptides were assigned to proteins using the default ProteinProphet and Apex parameters. Abundance estimates between replicates for each species were determined to be significantly correlated in all comparisons (Pearson’s r; all P values <0.000001). Down-sampling analyses revealed a substantially improved correlation in protein abundances between species using data pooled across replicates. Comparisons between species were therefore based on a single composite APEX value and assessed statistically using a Z-score with Bonferroni correction for multiple testing. However, we note that the majority of our results are not derived from statistical comparisons at the level of individual proteins but rather from emergent patterns across the full dataset. Soft clustering was performed using the Mfuzz R package (Kumar and Futschik 2007) for the 360 proteins identified in all three species and in the top quartile of absolute Z-scores in at least one of the pairwise comparisons between species. The optimality function was used to establish the most representative proximity of cluster components and cluster number was selected based on estimates from the minimum centroid distance and partition coefficient functions. Cluster membership was established based on the maximal membership values. Heatmaps were generated using heatmap.2 from the gplots R package using default parameters. Species and proteins were organized by dendrogram order, with the distance matrix computed by the dist R function using the euclidean method and the hclust R function using the complete linkage agglomeration method. Branching uncertainty was assessed with nonparametric bootstrap resampling using the pvclust R package (10,000 iterations of the Euclidean distance and complete linkage method with bootstrap samples sizes ranged from 0.5 to 1.4) (Suzuki and Shimodaira 2006). Confidence was assessed using Approximately Unbiased test P values (Shimodaira and Hasegawa 2001).

Ancestral State Reconstruction

Body and testes size were obtained from 9 Mus species, including M. musculus musculus, M. musculus domesticus, M. spretus, M. spicilegus, M. macedonicus, M. cookii and M. pahari (Gomez Montoto et al. 2011a), and M. caroli and M. minutoides (Vicens, Tourmente, et al. 2014). The Mus phylogeny and branch lengths of Sarver et al. (2017) were used. Ancestral body and testis size reconstructions were conducted independently using two approaches: the fastAnc function from the phytools R package (Revell and Graham Reynolds 2012) and the lognormal relaxed random walk (RRW) continuous trait substitution model in BEAST with a lognormal prior distribution (Drummond et al. 2012). The chain was run for 1000000 steps with sampling every 100. Convergence was assessed using the Tracer program, with ESS parameters confirmed to be >500. BEAST results were summarized using the TreeAnnotator program. Predicted testes mass (Y) was calculated using the formula of Kenagy and Trombulak (1986) for rodents: Y = 0.031X0.77, where X is the observed/reconstructed body mass. Relative testes size is then calculated as the ratio of reconstructed testes mass divided by the predicted testes mass.

Functional Classification, Gene Expression and Molecular Divergence Analysis

Mouse sperm proteome components were classified by their molecular and biological function using the PANTHER resource (Mi et al. 2013). Identification of enriched gene functions and estimation of significance levels was determined by the GOrilla tool (Eden et al. 2009). Information regarding mutant reproductive phenotypes associated with mouse sperm proteins was retrieved from Mouse Genome Informatics (http://www.informatics.jax.org) and additional functional annotation of significant proteins was conducted using the DAVID tool (http://david.abcc.ncifcrf.gov). Spermatogenesis stage-specific expression (Soumillon et al. 2013) and translational repression (Gan et al. 2013) data was integrated with sperm proteomic data and analyzed using chi-squared tests. Molecular divergence (dN and dS) of sperm proteins from Vicens, Luke et al. (2014) between M. musculus and M. spretus were estimated using the methodology of Nei and Gojobori (1986) as implemented in the PAML software package (Yang 2007).

Supplementary Material

Supplementary data are available at Molecular Biology and Evolution online.

Click here for additional data file.