When one phenotype is not enough: divergent evolutionary trajectories govern venom variation in a widespread rattlesnake species

Abstract

Understanding the origin and maintenance of phenotypic variation, particularly across a continuous spatial distribution, represents a key challenge in evolutionary biology. For this, animal venoms represent ideal study systems: they are complex, variable, yet easily quantifiable molecular phenotypes with a clear function. Rattlesnakes display tremendous variation in their venom composition, mostly through strongly dichotomous venom strategies, which may even coexist within a single species. Here, through dense, widespread population-level sampling of the Mojave rattlesnake, Crotalus scutulatus, we show that genomic structural variation at multiple loci underlies extreme geographical variation in venom composition, which is maintained despite extensive gene flow. Unexpectedly, neither diet composition nor neutral population structure explain venom variation. Instead, venom divergence is strongly correlated with environmental conditions. Individual toxin genes correlate with distinct environmental factors, suggesting that different selective pressures can act on individual loci independently of their co-expression patterns or genomic proximity. Our results challenge common assumptions about diet composition as the key selective driver of snake venom evolution and emphasize how the interplay between genomic architecture and local-scale spatial heterogeneity in selective pressures may facilitate the retention of adaptive functional polymorphisms across a continuous space.

Keywords: adaptive trait; diet; phenotypic variation; population structure; structural polymorphism; venom.

Figure 1.

Geographical variation in venom and diet of adult C. scutulatus. (a) Distribution of samples for which the major venom types were identified based on toxin genotypes; stars represent the sampling locations of the two representative individuals used for the genome-transcriptome-proteome analyses. (b) Two-ring pie-charts showing the proportion of mammals and reptiles from stomach contents (inner charts) and venom types (outer ring) for each locality.

Figure 2.

Toxin genotype and niche modelling. (a) Presence–absence matrix of toxin genes and admixture plot (TESS) with K = 2. (b) Niche models and sample distribution of the Mojave-Sonoran clade of Crotalus scutulatus with individuals represented by proportion of genetic clusters. Grey lines delineate ecoregion boundaries.

Figure 3.

Association between venom phenotypic variation, neutral genetic differentiation and environment. (a) Non-metric multidimensional scaling (NMDS) analysis of venom profiles shows great overall variation. Variation along NMDS1 is strongly correlated with the marked east–west environmental cline across Arizona (electronic supplementary material, table S12 and figures S9 and S10), whereas environmental associations along NMDS2, broadly separating the A–B transition, are weaker because global-scale variation hinders the detection of local-scale patterns. (b–g) Local-scale analysis along two transects (b) reveals sharp clines in various temperature (c–e) and precipitation (f,g) variables (see electronic supplementary material, table S11 for bioclimatic variable description) across the venom A–B transition zone.