Convergent evolution of conserved mitochondrial pathways underlies repeated adaptation to extreme environments.

Extreme environments test the limits of life; yet, some organisms thrive in harsh conditions. Extremophile lineages inspire questions about how organisms can tolerate physiochemical stressors and whether the repeated colonization of extreme environments is facilitated by predictable and repeatable evolutionary innovations. We identified the mechanistic basis underlying convergent evolution of tolerance to hydrogen sulfide (H2S)-a toxicant that impairs mitochondrial function-across evolutionarily independent lineages of a fish (Poecilia mexicana, Poeciliidae) from H2S-rich springs. Using comparative biochemical and physiological analyses, we found that mitochondrial function is maintained in the presence of H2S in sulfide spring P. mexicana but not ancestral lineages from nonsulfidic habitats due to convergent adaptations in the primary toxicity target and a major detoxification enzyme. Genome-wide local ancestry analyses indicated that convergent evolution of increased H2S tolerance in different populations is likely caused by a combination of selection on standing genetic variation and de novo mutations. On a macroevolutionary scale, H2S tolerance in 10 independent lineages of sulfide spring fishes across multiple genera of Poeciliidae is correlated with the convergent modification and expression changes in genes associated with H2S toxicity and detoxification. Our results demonstrate that the modification of highly conserved physiological pathways associated with essential mitochondrial processes mediates tolerance to physiochemical stress. In addition, the same pathways, genes, and-in some instances-codons are implicated in H2S adaptation in lineages that span 40 million years of evolution.

Stephen J. Gould made a strong case for the importance of contingency in evolution, famously quipping that replaying the “tape of life” would lead to different outcomes every time (1). However, despite the unpredictability of mutations, the effects of genetic drift, and other historical contingencies, convergent evolution of phenotypic traits and their underlying genes is common, indicating that natural selection sometimes finds repeatable and predictable solutions to shared evolutionary challenges (2, 3). A major challenge that remains is the identification of the ecological, genetic, and functional factors that might determine the repeatability and predictability of evolutionary outcomes (4).

Mitochondria and their genomes provide a fascinating model to ask questions about the predictability of evolution for two reasons: 1) Mitochondrial genomes were historically thought to be a prime example of contingency evolution because alternative genetic variants were assumed to be selectively neutral (5). This paradigm has been shifting though, with mounting evidence that mitochondria—and genes encoded in the mitochondrial genome—can play important roles in adaptation, especially in the context of physiochemical stress (6). 2) Mitochondria are critical for the cellular function of eukaryotes (7). Their function is dependent on the gene products from two genomes, the mitochondrial and the nuclear (8), which interact to ultimately shape whole-organism performance. Despite extensive characterization of allelic variation in mitochondrial genomes, it often remains unclear how variation in genes that contribute to mitochondrial function translates to variation in physiological and organismal function. Furthermore, it is not known whether exposure to similar selective regimes may cause convergent modifications of mitochondrial genomes and emergent biochemical and physiological functions in evolutionarily independent lineages.

Extreme environments that represent novel ecological niches are natural experiments to address questions about mechanisms underlying mitochondrial adaptations and illuminate the predictability of adaptive evolution of mitochondria. Among the most extreme freshwater ecosystems are springs with high levels of H2S, a potent respiratory toxicant lethal to metazoans due to its inhibition of mitochondrial ATP production (9). Multiple lineages of livebearing fishes (Poeciliidae) have colonized H2S-rich springs throughout the Americas and independently evolved tolerance to sustained H2S concentrations orders of magnitude higher than those encountered by ancestral lineages in nonsulfidic habitats (10). Here, we identify the molecular basis of an evolutionary innovation that facilitated the independent colonization of extreme environments (increased H2S tolerance) and ask if the underlying mechanisms have evolved in convergence in disparate lineages of livebearing fishes.

H2S toxicity and detoxification are associated with highly conserved physiological pathways in mitochondria (Fig. 1) (11, 12), providing a priori predictions about the potential molecular mechanisms underlying adaptation to this strong source of selection. Toxic effects of H2S result from binding to and inhibition of COX (cytochrome c oxidase, complex IV) in the oxidative phosphorylation (OxPhos) pathway, which contains subunits encoded in both the nuclear and the mitochondrial genomes (13). Animal cells can also detoxify low concentrations of endogenously produced H2S via the mitochondrial SQR (sulfide:quinone oxidoreductase) pathway, which is linked to OxPhos but entirely encoded in the nuclear genome (14). We have previously shown that genes associated with both pathways are under divergent selection and differentially expressed between fish populations in sulfidic and nonsulfidic habitats (10). These include nuclear and mitochondrial genes encoding subunits of the direct toxicity target (COX) and the nuclear gene encoding the enzyme mediating the first step of detoxification (SQR) (10). Tolerance to H2S may, therefore, be mediated by resistance (modification of toxicity targets that reduce the negative impact of H2S), regulation (modification of physiological pathways that maintain H2S homeostasis), or both (9).

Fig. 1.

(A) Physiological pathways associated with H2S toxicity and detoxification are located in the inner mitochondrial membrane. H2S inhibits OxPhos (orange enzymes, encoded by genes in the mitochondrial and nuclear genomes) by binding to cytochrome c oxidase (COX) (Complex IV). H2S can be detoxified through sulfide:quinone oxidoreductase (SQR) (green enzyme, encoded by a gene in the nuclear genome) and additional enzymes (indicated by the asterisk). (B) Relative activity of COX upon H2S exposure, which was primarily explained by an interaction between habitat type of origin and ambient H2S concentration (). (C) Activity of SQR as a function of H2S concentration, which was explained by an interaction between habitat type of origin and H2S concentration (). (D) Relative change in mitochondrial H2S concentrations in the liver of live fish exposed to different levels of environmental H2S. Variation in mitochondrial H2S levels were explained by habitat type of origin and exogenous H2S concentration (). (E) Relative spare respiratory capacity of isolated liver mitochondria at different levels of H2S. The interaction between habitat type of origin and drainage of origin best explained variation in spare respiratory capacity (). For all graphs, yellow colors denote P. mexicana from H2S-rich habitats, and blue denotes P. mexicana from nonsulfidic habitats. Symbols stand for populations from different river drainages (■: Tac; ▲: Puy; ●: Pich; see ).

Based on these previous results, we hypothesized that the repeated modification of enzymes in the OxPhos and SQR pathways in P. mexicana populations from sulfidic habitats leads to an increased ability to maintain mitochondrial function in the presence of H2S. In the present study, we used a series of in vivo and in vitro assays to identify the functional consequences of modifications to the OxPhos and SQR pathways in evolutionarily independent population pairs of P. mexicana from adjacent sulfidic and nonsulfidic habitats that are situated in different river drainages. In addition, we hypothesized that convergent molecular modifications in the same pathways underlie the convergent evolution of H2S tolerance across different lineages of poeciliid fishes. Hence, we also used phylogenetic comparative analyses of gene expression and analyses of molecular evolution to detect patterns of molecular convergence in 10 lineages of sulfide spring poeciliids and ancestors from nonsulfidic habitats.

Results and Discussion

Sulfide Spring P. Mexicana Exhibit a Resistant Toxicity Target.

If resistance is the primary mechanism of tolerance, we would predict that COX function is maintained in the presence of H2S in fish from sulfidic populations but not those from nonsulfidic populations. Quantification of COX function indicated that enzyme activity generally declined with increasing H2S concentrations, but this decline was reduced in populations from sulfidic habitats (Fig. 1; habitat × H2S: P < 0.001, ). Even though the drainage of origin was not retained as an explanatory variable in statistical models (), COX activity in one H2S-tolerant population [Tacotalpa (Tac)] declined just as in nonsulfidic populations (Fig. 1). The other two P. mexicana populations from sulfidic habitats [Puyacatengo (Puy) and Pichucalco (Pich)] maintained significant COX activity even at the highest H2S concentrations, which should reduce the negative impact of H2S on cellular respiration. These results are consistent with previous analyses (15) and indicate that resistance may contribute to H2S tolerance in some populations but cannot explain the repeated evolution of H2S tolerance by itself.

Sulfide Spring P. Mexicana Can Regulate Mitochondrial H2S through Increased Detoxification.

We also tested whether tolerant and intolerant populations differ in their ability to detoxify H2S by conducting enzyme activity assays of SQR. Activity of SQR was significantly higher in mitochondria from sulfidic populations at intermediate and high H2S concentrations (Fig. 1; habitat × H2S: P < 0.001 in ), likely helping fish from sulfidic habitats to maintain H2S homeostasis during environmental exposure. To test this prediction in vivo, we used a novel mitochondria-specific H2S-probe (MitoA) that allows for the monitoring of relative H2S levels inside the mitochondria of living organisms (16). We measured mitochondrial H2S concentrations in this manner using laboratory-reared fish that were exposed to varying levels of environmental H2S. Because laboratory-reared fish were not available for the population pair from Pich, only two population pairs were used for this analysis. Overall, mitochondrial H2S concentrations increased with environmental exposure (P = 0.001) and was higher in fish from nonsulfidic habitats (P < 0.001 in ). H2S concentrations in mitochondria isolated from livers (Fig. 1) and other organs () of fish from nonsulfidic habitats increased above control levels at all exposure concentrations. In contrast, mitochondrial H2S concentrations in isolates of fish from sulfidic populations did not usually exceed control levels and remained lower than levels in fish from nonsulfidic habitats. Together, these results indicate that populations of P. mexicana from sulfidic habitats can detoxify H2S at higher rates and, thus, regulate mitochondrial H2S upon environmental exposure.

Sulfide Spring P. Mexicana Can Maintain Mitochondrial Function in the Presence of H2S.

Modification of the OxPhos and SQR pathways in P. mexicana suggests that mitochondrial adaptations are key to the evolution of H2S tolerance. Therefore, mitochondrial function of fish from sulfidic habitats should be maintained upon exposure to H2S. We tested this hypothesis by quantifying different aspects of mitochondrial function (basal respiration, maximal respiration, and spare respiratory capacity) along a gradient of H2S concentrations using an ex vivo coupling assay. As expected, all aspects of mitochondrial function generally declined with increasing H2S (Fig. 1 and ). Comparison of mitochondrial function between adjacent populations in sulfidic and nonsulfidic habitats indicated no differences in basal respiration (). However, individuals from sulfidic populations were able to maintain maximal respiration and spare respiratory capacity at higher levels compared to individuals from nonsulfidic habitats of the same river drainage (Fig. 1), even though the magnitude of difference and the shape of response curves varied (; significant drainage × habitat interactions in ). These findings indicate that mitochondria of H2S-tolerant individuals continue to produce ATP in the presence of a potent inhibitor that reduces mitochondrial function in ancestral lineages.

Overall, our quantitative analyses indicate clear patterns of convergence in functional physiological traits associated with H2S tolerance. Nonetheless, further inspection of the results also reveals lineage-specific patterns (especially in H2S-dependent COX activity and mitochondrial respiration), indicating that evolutionary responses across lineages are similar but not necessarily identical. These idiosyncrasies are consistent with the results of previous comparative transcriptome analyses, which revealed a large number of genes that are under selection or differentially expressed in just a subset of lineages in addition to genes that are consistently differentially expressed and under selection across all lineages (17–19). Based on their functions, the OxPhos and SQR pathways undoubtedly include some major-effect genes influencing H2S tolerance in different populations of P. mexicana but tolerance—as an emergent physiological trait—is a complex trait impacted by other genes as well. In the future, quantitative genetic analyses will be required to understand how other loci contribute to tolerance within each lineage and how population-specific patterns of genetic differentiation might shape variation in functional physiology evident in our data.

Convergence among P. Mexicana Populations Is Shaped by Selection on De Novo Mutations and Standing Genetic Variation.

The convergent evolution of H2S tolerance in P. mexicana begs questions about the origin of adaptive alleles (20). On microevolutionary scales, convergence may be a consequence of the repeated assembly of related alleles into different genomic backgrounds either through selection on standing genetic variation or introgression (21, 22). However, the epitome of convergent evolution is, arguably, the independent origin of adaptive mutations at the same locus that lead to consistent functional outcomes (23). To identify convergence at a genomic level, we resequenced whole genomes of multiple P. mexicana individuals from sulfidic and nonsulfidic habitats. Analyzing phylogenetic relationships among P. mexicana populations (with Poecilia reticulata as an outgroup) using 13,390,303 single nucleotide polymorphisms (SNPs) distributed across the genome confirmed three independent colonization events of sulfide springs and distinct evolutionary trajectories for sulfide spring populations in different drainages (Fig. 2), as inferred by previous studies (24). If adaptive alleles arose separately through de novo mutation in each sulfide spring population, we would expect that putative adaptive alleles mirror these relationships as previously documented for H2S-resistance alleles in mitochondrial COX subunits (15). However, patterns of divergence () and local ancestry were highly variable across the genome. Classifying local patterns of genetic similarity using a hidden Markov model and a self-organizing map allowed us to identify genomic regions in which ancestry patterns deviate from the genome-wide consensus, including multiple regions with a strong signal of clustering by ecotype (sulfidic vs. nonsulfidic populations). Such clustering by ecotype occurred in less than 1% of the genome () but included genomic regions encoding key genes associated with H2S detoxification (e.g., SQR and ETHE1, Fig. 2 and Dataset S2). Clustering by ecotype indicates a monophyletic origin of putatively adaptive alleles at these loci that are shared across independent lineages of sulfide spring P. mexicana as a consequence of selection on standing genetic variation or introgression (25), although the latter scenario is less likely considering the geographic barriers and strong survival selection against migrants from sulfidic to nonsulfidic habitats (26). Consequently, multiple mechanisms—not just selection on de novo mutations (19)—played a role in the convergent evolution of H2S tolerance in P. mexicana.

Fig. 2.

(A) Phylogeny of different populations in the P. mexicana species complex (with P. reticulata as an outgroup) based on genome-wide SNPs. Colors indicate sulfidic (yellow) and nonsulfidic (blue) lineages. (B) Local ancestry patterns around genes encoding two enzymes involved in H2S detoxification, SQR and ETHE1. Gray bars represent the local ancestry pattern (cactus) associated with each region. Unrooted trees represent local ancestry relationships with sulfidic lineages colored in yellow and nonsulfidic lineages colored in blue. Cacti 10 and 19 show clear clustering by ecotype. In cacti 1, 5, and 12, four of five sulfidic individuals cluster together.

Convergent Modifications of Toxicity Targets and Detoxification Pathways Are Evident on Macroevolutionary Scales.

While selection on standing genetic variation and introgression can contribute to convergent evolution on microevolutionary scales, adaptive alleles are unlikely to be shared among lineages on macroevolutionary scales due to high phylogenetic and geographic distances separating gene pools (27). The absence of convergence in molecular mechanisms at broader phylogenetic scales might indicate the importance of contingency in evolution as asserted by Gould (3). In contrast, the presence of convergence would indicate that fundamental constraints limit the number of solutions for a functional problem (28).

We used phylogenetic comparative analyses of gene expression and analyses of molecular evolution to detect patterns of molecular convergence in 10 lineages of sulfide spring poeciliids and ancestors in nonsulfidic habitats (). This included members of five genera that span over 40 million years of divergence and occur in different biogeographic contexts (). We found evidence for convergence in both gene expression and sequence evolution. Variation in overall gene expression was strongly influenced by phylogenetic relationships (Fig. 3). However, 186 genes exhibited significant evidence for convergent expression shifts in sulfide spring fishes (Fig. 3 and Dataset S3), segregating lineages based on habitat type of origin, irrespective of phylogenetic relationships (Fig. 3). The only outlier was Limia sulphurophila, which clustered with nonsulfidic lineages despite significant expression differences with its sister, Limia perugiae. Functional annotation indicated that genes with convergent expression shifts were enriched for biological processes associated with H2S detoxification (SQR pathway, Fig. 3), the processing of sulfur compounds, and H2S toxicity targets in OxPhos ().

Fig. 3.

(A) Multidimensional scaling (MDS) plot of overall gene expression patterns across 20 lineages of poeciliid fishes. Black lines represent phylogenetic relationships among lineages; color represents habitat type of origin (yellow: sulfidic; blue: nonsulfidic). (B) Expression variation of 186 genes with evidence for convergent expression shifts (z-transformed fragments per kilobase of transcript per million mapped reads). Colors represent expression levels as indicated by the scale. The neighbor-joining tree on the left groups lineages based on expression similarity. The cladogram on the right shows the phylogenetic relationship among lineages. Pictures on the side are examples of sulfide spring fishes (from top to bottom): Xiphophorus hellerii, Pseudoxiphophorus bimaculatus, Gambusia holbrooki, G. sexradiata, G. eurystoma, Poecilia latipinna, P. sulphuraria (Pich), P. mexicana (Tac), P. mexicana (Puy), and Limia sulphurophila. (C) MDS plot of the expression of 186 genes with evidence for convergent expression shifts. (D) Boxplots with mean expression levels of different components of the SQR pathway across lineages from sulfidic (yellow) and nonsulfidic (blue) habitats.

We also identified 11 genes with elevated nonsynonymous to synonymous substitution rates across the phylogeny, including three mitochondrial genes that encode subunits of H2S’s toxicity target (COX1 and COX3) and OxPhos complex III (CYTB; Dataset S4). Most amino acid substitutions in COX1 and COX3 occurred in a lineage-specific fashion, but convergent substitutions across clades occurred at six codons in COX1 and two codons in COX3 (Fig. 4). These findings suggest that modifications of H2S toxicity targets and detoxification pathways are not only critical in the evolution of H2S tolerance in P. mexicana, but also they have evolved in convergence in other lineages that were exposed to the same source of selection.

Fig. 4.

Amino acid differences in COX1 and COX3 between lineages from sulfidic (yellow) and nonsulfidic (blue) habitats. Derived amino acids are shown in red. Bold letters indicate codons with convergent amino acid substitutions in different clades (separated by black horizontal lines) of sulfide spring fishes. Numbers above represent the position in the amino acid sequence.

Conclusions

We capitalized on past evolutionary genetics studies that compared P. mexicana populations from sulfidic and nonsulfidic environments (10) to test hypotheses about functional ramifications of genetic differences and their impact on organismal performance. As predicted, we found that the repeated evolution of H2S tolerance in independent P. mexicana populations is mediated both by modifications of a direct toxicity target (causing increased resistance to H2S) and a pathway involved in detoxification (causing an increased ability to regulate mitochondrial H2S). Similar modifications to COX and SQR have been hypothesized to mediate H2S adaptation in other groups of organisms (29–31), but the evolutionary context and the consequences for mitochondrial function in these cases remain unknown. Overall, our analyses indicated that closely related populations can exhibit substantial differences in what we assume to be highly conserved physiological pathways associated with the function of mitochondria. Modification of mitochondrial processes, consequently, can be critical in mediating adaptation to different environmental conditions on microevolutionary scales, underscoring the long overlooked role of mitochondria in adaptive evolution (6).

Our comparative transcriptome analyses across a broader sampling of sulfide spring fishes further indicated that colonization of novel niches with extreme environmental conditions can arise through the convergent modification of conserved physiological pathways. The convergent evolution of high H2S tolerance across species is the result of repeated and predicted modifications of the same physiological pathways, genes, and—in some instances—codons associated with mitochondrial function. This convergence at multiple levels of biological organization is likely a consequence of constraint because the explicit biochemical and physiological consequences of H2S limit the ways organisms can cope with its toxicity (32, 33). Due to these constraints, molecular convergence is not only evident on microevolutionary scales where selection can repeatedly assemble related alleles into different genomic backgrounds, but also on macroevolutionary scales including lineages separated by over 40 million years of evolution.

That said, there is an inordinate amount of genetic and gene expression variation that seemingly varies idiosyncratically across different lineages. In comparative analyses of highly quantitative traits (such as H2S tolerance), there is an inherent bias to emphasize the importance of shared modifications in adaptation, while we tend to dismiss lineage-specific patterns as noise. But how lineage-specific genetic and gene expression variation interacts with molecular mechanisms that have evolved in convergence remains largely unknown for most study systems. So, if we replayed the tape of life, the same characters may make an appearance in the same setting, but the overall plots may still unfold in very different ways when many characters are part of the story.

Methods

The following sections provide a synopsis of the procedures used in this study. Detailed materials and methods are provided in the .

Sampling.

Samples of P. mexicana for comparative biochemical and physiological analyses were collected from three population pairs from the Tac, Puy, and Pich drainages in Mexico, each including evolutionarily independent H2S-tolerant and ancestral, intolerant populations () (34). With the exception of measurements of mitochondrial H2S levels, which were conducted with common-garden-reared individuals, all assays were conducted with specimens collected in the field. For macroevolutionary analyses of convergence, we collected specimens from multiple species that have independently colonized H2S-rich springs in the United States, Mexico, and the island of Hispaniola, as well as from geographically and phylogenetically proximate lineages in nonsulfidic habitats. Our sampling included populations of P. latipinna and G. holbrooki in Florida; populations of the P. mexicana species complex (including P. mexicana, P. sulphuraria, and Poecilia limantouri), the G. sexradiata species complex (including G. sexradiata and G. eurystoma), P. bimaculatus, and X. hellerii in Mexico; as well as populations of the Limia perugiae species complex (L. perugiae and L. sulphurophila) in the Dominican Republic ().

Measurements of Enzyme Activity.

Measurements of COX and SQR activities were conducted using isolated mitochondria from three P. mexicana population pairs and followed previously established methods. COX activity was quantified by measuring a decrease in absorbance at 550 nm resulting from the protein’s oxidation of reduced cytochrome c (35, 36). Samples were assayed along a gradient of H2S concentrations between 0 and 0.9 µM. SQR activity was quantified by measuring the reduction of decyl-ubiquinone (37) along a gradient of H2S concentrations between 0 and 100 µM. Variation in enzyme activity was analyzed using linear mixed-effects models to estimate the effects of habitat type of origin (sulfidic vs. nonsulfidic), drainage, and H2S concentration. Alternative models were assessed by using Akaike information criteria with finite sample correction (AICC) (38).

Measurements of Mitochondrial H2S.

H2S exposure experiments and measurement of mitochondrial H2S were conducted using common-garden-raised F1 individuals from two population pairs (Tac and Puy). We used a mitochondria-specific H2S probe (MitoA) that can be injected into living organisms, where it accumulates within mitochondria due to its unique chemical structure (16). Inside mitochondria, MitoA reacts with H2S and is converted to MitoN. The ratio of MitoN/MitoA can be quantified using liquid chromatography with tandem mass spectrometry and serve as a metric of mitochondrial H2S. We used a MitoA protocol specifically validated for P. mexicana (39). Fish injected with MitoA were exposed to different environmental H2S concentrations for 5 h. Gill, liver, brain, and muscle tissues were sampled from fish and used for the quantification of MitoN/MitoA ratios. Standardized MitoN/MitoA ratios were analyzed using linear mixed-effects models with habitat type of origin, drainage, and ambient H2S concentration as predictor variables. Alternative models were assessed by using AICC.

Measurements of Mitochondrial Function.

Mitochondrial function was measured in mitochondria isolated from livers from three P. mexicana population pairs. Mitochondrial function was assayed using a Seahorse XFe96 Extracellular Flux Analyzer (Agilent Technologies, Santa Clara, CA), which allows for the quantification of oxygen consumption rates (OCR) of isolated mitochondria in 96-well plates (40). Measuring mitochondrial OCR in the presence of different substrates and inhibitors allows for the quantification of a variety of mitochondrial functions (41), and we measured basal respiration, maximal respiration, and spare respiratory capacity as indicators of mitochondrial function across a range of H2S concentrations between 0 and 90 µM. To compare responses in mitochondrial function to H2S between sulfidic and nonsulfidic populations, we employed a drug response analysis where separate n-parameter logistic regressions were fit for each mitochondrial isolate with the metrics of mitochondrial function (basal respiration, maximal respiration, and spare capacity) as dependent variables and H2S concentration as an independent variable (42). For each model, we quantified the area under the curve (AUC) based on Simpson’s rule, and higher AUC values represent an increased ability to maintain mitochondrial function in the presence of H2S. AUC values inferred for different mitochondrial isolates were then used as a dependent variable in linear models with habitat type of origin and drainage as predictor variables. Alternative models were assessed by using AICC.

Comparative Genomics and Local Ancestry Analysis.

To test hypotheses about the origin of putatively adaptive alleles, we reanalyzed data from Brown et al. (25), which included whole-genome sequences from all sites known to harbor H2S-tolerant populations of the P. mexicana species complex (n = 5) as well as adjacent nonsulfidic habitats (n = 5; ). Raw reads were mapped to the Xiphophorus maculatus reference genome (43), and we called SNPs using the Genome Analysis Toolkit (44). Local ancestry patterns were identified with Saguaro (45), which uses a combination of a hidden Markov model and a self-organizing map to build cacti (matrices of pairwise genetic distance between samples) that describe phylogenetic relationships among samples in specific genomic regions. Saguaro was executed for 29 iterations, which resulted in 30 cacti describing local ancestry patterns for segments of the genome.

Comparative Transcriptomics, Phylogenetic Comparative Analyses, and Molecular Evolution.

For comparative transcriptomic analyses, we isolated RNA from gills for 5 to 6 individuals each from 10 sulfidic and 10 nonsulfidic lineages distributed across multiple genera in the family Poeciliidae. Raw reads from transcriptome sequencing were mapped to the P. mexicana reference genome (46), and the number of RNA-sequencing reads mapped to each gene was determined for each individual using cufflinks (47). To analyze patterns of gene expression across lineages, we used a phylogenetic comparative approach for analyzing gene expression variation to explicitly account for the effects of evolutionary relationships. Specifically, we used individual-level expression data of the top 5,000 expressed genes as dependent variables in expression variance and evolution (EVE) models that identify genes exhibiting convergent shifts in expression upon the colonization of sulfidic habitats based on a phylogeny of focal taxa (48, 49). To identify genes with potential structural or functional changes in lineages from sulfidic habitats, we analyzed patterns of molecular evolution in protein coding genes included in the analysis of gene expression. For each gene, we generated a consensus sequence for all exons across individuals from the same populations and then used branch models implemented in the program codeml from the PAML package to test for evidence of positive selection (50).

Data Availability.

Data and code associated with biochemical and physiological analyses are available on GitHub (github.com/michitobler/convergent_h2s_evolution). All sequence data are available at the National Center for Biotechnology Information (NCBI) (https://www.ncbi.nlm.nih.gov under BioProject) (accession nos. PRJNA473350 and PRJNA608180).