Research Article

Natural hybridization reveals incompatible alleles that cause melanoma in swordtail fish

See allHide authors and affiliations

Science  15 May 2020:
Vol. 368, Issue 6492, pp. 731-736
DOI: 10.1126/science.aba5216

Mapping vertebrate incompatibility alleles

Deleterious gene interactions may underlie the observed hybrid incompatibilities. However, few genes underlying hybrid incompatibilities have been identified, and most of these involve species that do not hybridize in natural conditions. Powell et al. used genome sequencing to map genes likely responsible for incompatibilities that reduce fitness in naturally occurring hybrid swordtail fish. These gene combinations result in malignant melanoma, which is found in naturally hybridizing populations but is not present in the parental populations (see the Perspective by Dagilis and Matute). Using genome and population resequencing, the authors performed a genome-wide association study to identify potentially causative mutations. Using an admixture mapping approach that assessed introgression between multiple swordtail fish species, the authors suggest that lineages carry different genes that interact with the same candidate gene, resulting in the observed melanomas and providing insight into convergent hybrid incompatibles that arise between species.

Science, this issue p. 731; see also p. 710


The establishment of reproductive barriers between populations can fuel the evolution of new species. A genetic framework for this process posits that “incompatible” interactions between genes can evolve that result in reduced survival or reproduction in hybrids. However, progress has been slow in identifying individual genes that underlie hybrid incompatibilities. We used a combination of approaches to map the genes that drive the development of an incompatibility that causes melanoma in swordtail fish hybrids. One of the genes involved in this incompatibility also causes melanoma in hybrids between distantly related species. Moreover, this melanoma reduces survival in the wild, likely because of progressive degradation of the fin. This work identifies genes underlying a vertebrate hybrid incompatibility and provides a glimpse into the action of these genes in natural hybrid populations.

The emergence of reproductive barriers between populations is the first step in the process of speciation and drives Earth’s biological diversity, yet surprisingly little is known about how it occurs at the genetic level. The Dobzhansky-Muller model of hybrid incompatibility (13) posits that new mutations arising in diverging species can interact negatively in hybrids, generating lower hybrid viability or causing hybrid sterility. Although empirical work provides support for the general predictions of this model (4), progress in this area has been limited by a lack of knowledge about which genes interact to generate hybrid incompatibilities. Despite the effort devoted to this problem, only a dozen incompatible interactions have been mapped to the single-gene level [reviewed in (5)]. With so few known cases, it has been difficult to evaluate whether common genetic and evolutionary mechanisms underlie the emergence of incompatibilities (68).

With an increasing appreciation that hybridization is common across the tree of life (913), there has been renewed interest in identifying hybrid incompatibilities and understanding how these genes act as barriers in nature. Of hybrid incompatibilities that have been mapped to the single-gene level, most have been identified with crosses between model species that no longer naturally hybridize (4, 5). As a result, it is unclear whether these mapped incompatibilities were important in the initial divergence between species or arose after these lineages had stopped exchanging genes.

One such example is the melanoma receptor tyrosine-protein kinase (xmrk) gene in swordtail fish (genus Xiphophorus). xmrk is one of two identified genes in vertebrates that drive hybrid incompatibility (the other being the regulator of mammalian recombination hotspots, prdm9) (14) and was one of the earliest described hybrid incompatibilities (15). In crosses between Xiphophorus maculatus and Xiphophorus hellerii, a malignant melanoma develops in a subset of F2 hybrids, emanating from natural pigmentation spots on the body and fins. This hybrid incompatibility is the result of an interaction between the xmrk gene derived from X. maculatus and an unknown locus derived from X. hellerii (16).

Despite work that demonstrated the role of xmrk in the development of hybrid melanomas, its importance as a barrier between species has been debated (17). This is because X. maculatus and X. hellerii diverged ~3 million years ago and do not naturally hybridize (18). Moreover, because melanoma in X. maculatus x X. hellerii laboratory hybrids develops later in life, it is unclear whether it affects survival and reproduction (17).

Melanoma occurs in hybrids between recently diverged swordtail species

We identified a phenotypically similar melanoma in natural hybrids formed between the swordtail fish species Xiphophorus birchmanni and Xiphophorus malinche. X. birchmanni and X. malinche are closely related and hybridize in the wild (~0.5% differences per base pair and ~250,000 generations diverged) (19). Although a subset of hybrids are viable and fertile, there is evidence of selection against hybrid incompatibilities (2022). In some populations, hybrids develop melanoma early in life (13 ± 4% of males develop melanoma before sexual maturity) (fig. S1).

Melanoma in X. birchmanni x X. malinche hybrids develops from a phenotype derived from X. birchmanni called the “spotted caudal,” which is a dark blotch on the caudal fin generated by clusters of macromelanocyte cells (Fig. 1A and fig. S2) (23). The spotted caudal trait occurs at intermediate frequencies in X. birchmanni but is absent from X. malinche populations (Fig. 1B). Some hybrid populations have a high frequency of the trait and exhibit phenotypes that extend beyond the range of those observed in X. birchmanni (Fig. 1, B and C, and figs. S3 and S4), including invasion of macromelanocyte cells into the body, where they are normally absent. Tracking of hybrids in the laboratory documents the progression of the trait from a phenotype typical of the X. birchmanni spot to the expanded trait found in some hybrids (Fig. 1D and fig. S1) (24). Histological sections from hybrid individuals revealed penetration of melanocytes into the musculature and invasion of surrounding tissues, which is indicative of a malignant melanoma (Fig. 1E and fig. S5).

Fig. 1 Hybridization generates a high incidence of melanoma.

(A) Naturally hybridizing species X. malinche (top) and X. birchmanni (middle) differ in morphological traits, including the presence of a melanin pigment spot that is polymorphic in X. birchmanni. In hybrids, this spotting phenotype can transform into a melanoma (bottom). (B) Whereas X. birchmanni populations segregate for the presence of this spot, the trait is absent in X. malinche populations; hybrid populations have high frequencies of this trait. (C) The trait is at higher frequencies in hybrid populations and covers more of the body. Shown here is invasion area, or the melanized body surface area outside of the caudal fin (normalized for body size). Hybrid phenotypes are shown from three populations on the Río Calnali (fig. S3). AGCZ, Aguazarca; CALL, Calnali low; CHAF, Chahuaco falls. (D) Spots expand more over a 6-month period in hybrids than in X. birchmanni individuals. (E) A cross section of the caudal peduncle from a Chahuaco falls hybrid (10× magnification). Melanoma cells invading the body and muscle bundles are visually evident (indicated with blue stars). (F) Gene ontology categories enriched in melanoma tissue compared with normal caudal tissue (24, 45). The size of the dots reflects the number of genes identified, and the color corresponds to the P value. Categories with undefined odds ratios (not plotted) are listed in table S1. In (B) and (D), the plot shows the mean, and whiskers indicate two standard errors of the mean. Individual points show the raw data.

We performed mRNA-sequencing of hybrid individuals that varied in the degree of expansion of their spot (figs. S2 and S6). Functional enrichment analysis indicated changes in the regulation of a number of melanoma-associated gene categories, such as pigment cell differentiation and regulation of cytoskeletal organization, including several implicated in melanoma in other fish species (Fig. 1F and table S1) (24, 25).

Melanoma is extremely rare in nonhybrid individuals (Fig. 1C) (6, 26), and we have not identified a single wild-caught X. birchmanni male with melanoma (1296 males collected from 2017 to 2019, 0 with melanoma). Laboratory-reared individuals indicate that environmental levels of ultraviolet irradiance or other natural carcinogens do not underlie differences in the frequency of melanoma between hybrid and parental populations (24). The presence of melanoma in hybrids, but not the parental species, suggests that this melanoma is a hybrid incompatibility generated by interactions between alleles in the X. birchmanni and X. malinche genomes.

Loci associated with the spotting phenotype in X. birchmanni

We mapped the genetic basis of the spotted caudal in X. birchmanni and the genetic basis of melanoma in interspecific hybrids. We generated de novo assemblies for both species using a 10X-based linked read approach followed by assembly into chromosomes with Hi-C data and annotated the resulting assemblies with RNA-sequencing (RNA-seq) data (24).

We collected low-coverage whole-genome sequence data for 392 adult male X. birchmanni individuals from a single collection site and performed a genome-wide association study (GWAS), scanning for allele frequency differences between spotted cases (n = 159) and unspotted controls (n = 233), evaluating the impact of population structure and low-coverage data (24). We identified a strong association between the spotting pattern and allele frequency differences on chromosome 21 at an estimated false positive rate of 5% (by permutation) (Fig. 2A and figs. S7 and S8) (24). Two distinct signals are evident on this chromosome, ~5 Mb apart (Fig. 2A) (24). The first peak is centered on the xmrk gene (fig. S9) (27), which arose through duplication of a gene homologous to the mammalian epidermal growth factor receptor ~3 million years ago (11, 28, 29). xmrk controls pigmentation patterns and drives hybrid melanomas in relatives of X. birchmanni and X. malinche (16). The signal at the second peak on chromosome 21 contains a single gene, the melanosome transporter gene myosin VIIA and Rab interacting protein (myrip) (Fig. 2A) (24, 30).

Fig. 2 Combined genome-wide association and admixture mapping approaches identify the genetic basis of the melanoma hybrid incompatibility.

(A) Results of genome-wide association scan for allele frequency differences between spotted cases and unspotted controls. (Top) Results can be seen for all chromosomes, and the red line indicates the genome-wide significance threshold, determined by permutation (24). (Bottom) Results from chromosome 21, where two distinct regions are strongly associated with spotting. (B) Admixture mapping in hybrids identifies associations between X. birchmanni ancestry on chromosome 21 and spot presence. Plotted here are log likelihood differences between models with and without ancestry at the focal site included as a covariate. The red line indicates the genome-wide significance threshold, determined by permutation (24). (C) When we treated melanocyte invasion as the focal trait and mapped associations with ancestry, we again identified associations with X. birchmanni ancestry on chromosome 21 but also identified a second region on chromosome 5 associated with X. malinche ancestry.

Genetic architecture of the melanoma incompatibility

For the melanoma phenotype, we used an admixture mapping approach, focusing our efforts on a hybrid population with high incidence of melanoma (19 ± 3% of adult males, the Chahuaco falls population). To infer local ancestry, we generated ~1X low-coverage whole-genome sequence data for 209 adult males from this population and applied a hidden Markov model to 680,291 ancestry informative sites genome-wide (20, 22, 31) [approximately one ancestry informative site per kilobase, (24)]. Simulations and analyses of laboratory-generated crosses indicate that we should have high accuracy in local ancestry inference (figs. S10 to S13) (24).

Using these data, we performed admixture mapping for spot presence and melanocyte invasion of the body (59% of individuals had spots, and 19% of individuals had melanoma). Admixture mapping for the presence of the spotted caudal revealed one strongly associated region on chromosome 21 (log likelihood difference of linear models, 23) (24) where spotting correlated with X. birchmanni ancestry (Fig. 2B). Because of the lower resolution of admixture, mapping this peak is broad, but the signal occurs in the same regions identified by our GWAS scan (24), confirming that the genetic basis of the spot is the same in X. birchmanni and hybrids. Simulations accounting for effect size inflations owing to the “winner’s curse” (32) suggest that X. birchmanni ancestry in this region explains ~75% of the variation in spot presence or absence (24).

Admixture mapping for melanoma identified an additional significant region on chromosome 5. In this case, melanocyte invasion was associated with X. malinche ancestry (Fig. 2C). A contingency test indicated a nonrandom association between X. birchmanni ancestry at the chromosome 21 peak and X. malinche ancestry at the chromosome 5 peak, with the prevalence of melanoma (Fisher’s exact test, P = 0.0005) (Fig. 3A). Individuals heterozygous for X. malinche ancestry at this region on chromosome 5 appear to have a lower risk of melanoma (Fig. 3A). Moreover, regardless of melanoma phenotype, spotted individuals that were heterozygous at the chromosome 5 peak had smaller spots than individuals that were homozygous (Student’s t test on log-normalized area P = 0.007) (fig. S14).

Fig. 3 Interactions between chromosomes 5 and 21 are associated with melanoma in hybrids.

(A) Proportion of individuals with melanoma as a function of ancestry at the associated regions on chromosome 5 and chromosome 21. The blue dashed line indicates the expected proportion of cases if melanoma risk were equally distributed among individuals with at least one birchmanni allele at chromosome 21. We only had one observation for the bir-bir and bir-het genotypes. (B) The xmrk sequence in X. birchmanni harbors two mutations (G364R and C582S) that transform xmrk to a constitutively active state (33, 46). The schematic compares the ancestral form of the protein (egfrb) to the predicted structure of xmrk in X. birchmanni. Proteins are shown in red, and the cell membrane is shown in gray. In xmrk, residues R364 and S582 promote intramolecular disulfide bonds that cause protein dimerization and phosphorylation (blue circles) (33, 46). (Single-letter abbreviations for the amino acid residues are as follows: C, Cys; G, Gly; R, Arg; and S, Ser. In xmrk, amino acids were substituted at certain locations; for example, G364R indicates that glycine at position 364 was replaced by arginine.) (Inset) A partial clustal alignment of X. birchmanni egfrb and xmrk with these substitutions highlighted. Colors indicate properties of the amino acid, and asterisks indicate locations where the amino acid sequences are identical. (C) Clustal alignment showing the N terminus of cd97 in X. birchmanni and X. malinche. We observed a substitution in a conserved epidermal growth factor–binding domain (gray rectangles). (Inset) The substitution found in X. birchmanni is not present in closely related species. (D) Expression of cd97 based on RNA-seq data in melanoma, spotted, and unspotted tissue from Chahuaco falls hybrids (four biological replicates per group). (E) Real-time quantitative PCR of cd97 from caudal fin tissue from X. malinche, X. birchmanni, and natural and F1 hybrids (four to nine biological replicates per group). In (D) and (E), large solid dots indicate the mean, and whiskers indicate two standard errors of the mean. Individual points show the raw data.

Linking molecular changes to hybrid melanoma

Our GWAS identified associations between the spotted caudal and both xmrk and myrip (Fig. 2A), making it initially unclear whether either or both of these genes interacts with malinche ancestry on chromosome 5 to produce melanoma. Although both are associated with the spotting pattern that precedes melanoma, myrip is not expressed in an RNA-seq dataset of adult caudal tissue, nor is it expressed in the melanoma itself (fig. S15). By contrast, xmrk is expressed in caudal tissue and has higher expression in spotted than unspotted tissues (fig. S15). In addition, functional studies have linked xmrk to the development of melanoma. We identified two amino acid substitutions in X. birchmanni that fall within the extracellular domain of xmrk known to drive the oncogenic properties of xmrk in vitro (Fig. 3B) (33), and transgenic studies have demonstrated that overexpression of xmrk causes the formation of tumors (25, 34, 35). Although myrip may not be directly involved in the development of melanoma, past work in other swordtail species has suggested the presence of “patterning” loci linked to xmrk [reviewed in (23)]. Given myrip’s role in melanosome transport, we speculate that it could play a role in pigmentation patterning, which occurs in the first several weeks of life.

The region on chromosome 5 associated with X. malinche ancestry and melanoma contains only two genes, a gene called cd97 and a fatty acid transporter gene (Figs. 2C and 3C). Although this is unusually high resolution given our admixture-mapping approach, subsampling the data indicated that this scale of resolution is likely the result of high recombination rates in this region (24). We therefore sought to better characterize the two genes in this region.

The ortholog of cd97 in mammals plays a role in epithelial metastasis and is associated with tumor invasiveness in cancers (3638). Accordingly, we found that cd97 is up-regulated in RNA-seq data from melanotic tissue in hybrids, whereas the fatty acid transporter gene is not (Fig. 3D); nor is this pattern of up-regulation observed in any other gene within 100 kb of this region (24). In addition, of five amino acid changes between X. birchmanni and X. malinche in cd97, one occurs in a conserved epidermal growth factor–like calcium-binding domain (Fig. 3C).

We further investigated differences in expression of cd97 using a targeted quantitative polymerase chain reaction (qPCR) approach. We found that cd97 was expressed at low levels in the caudal fin tissue of X. birchmanni, regardless of spotting phenotype, but at similarly high levels in X. malinche and in natural and artificial hybrids [analysis of variance (ANOVA) P = 1−5, Tukey post hoc; all groups different from X. birchmanni at P < 0.005) (Fig. 3E). Higher expression of cd97 in X. malinche and hybrids is not tissue-specific and surprisingly does not appear to be driven by cis-regulatory differences (fig. S16) (24). We do not know whether the link between X. malinche ancestry at cd97 and melanoma is driven by coding or regulatory differences (24). However, in mammals, overexpression of cd97 has been linked to tumor metastasis; a similar mechanism could be involved here, given that high expression of cd97 coincides with invasion of other tissues with melanoma cells.

Independent evolution of a melanoma incompatibility

Although the role of xmrk in the maculatus-hellerii hybrid incompatibility has been known for 30 years (39), the identity of the interacting gene is not known. Laboratory crosses have narrowed to a ~5 Mb region on chromosome 5 but have not yet identified the underlying gene, although candidates have been proposed (16). We identified a distinct region on chromosome 5 (Fig. 2C), more than 7 Mb away from the region identified in the hellerii-maculatus cross. Alignments of chromosome 5 confirm that cd97 is at the same location in all four species (fig. S17), and linkage disequilibrium between cd97 and the region identified in hellerii-maculatus decays to background levels in hybrids (24). Using simulations, we ruled out a lack of power to detect an association between this region and melanoma, assuming a similar effect size to that seen in the hellerii-maculatus cross (fig. S18). These results indicate that the incompatibility has a partially distinct genetic basis in the two crosses generating hybrid melanoma. However, we may not have mapped all components of the melanoma incompatibility, particularly if other genes have subtle impacts on melanoma risk (24).

These mapping results are surprising because they suggest that a melanoma incompatibility involving xmrk emerged independently in two distinct lineages. Despite the evolutionary distance between these species (fig. S19), it is possible that the melanoma incompatibility arose through similar evolutionary paths in both cases. X. hellerii and its relatives lack xmrk (39), either because the lineage leading to this clade diverged before xmrk arose (fig. S19) or because of an ancient loss of xmrk. By contrast, many species in the lineage leading to X. birchmanni and X. malinche retained xmrk (fig. S19) (24), but we found that xmrk has been deleted in X. malinche since its divergence from X. birchmanni (fig. S20) (24). We speculate that the loss of xmrk in X. malinche could have changed the level of constraint on interacting genes in this lineage, and if so, similar evolutionary mechanisms could be at play in X. hellerii.

Selection on melanoma in natural populations

Although the melanoma that forms in birchmanni x malinche hybrids appears to be deleterious from its development early in life (fig. S1) and its malignancy (Fig. 1E and fig. S5), we wanted to evaluate its impact in natural hybrid populations. Over several years, we observed shifts in the frequency of the spotted caudal trait between juvenile and adult males (24). Specifically, in hybrid populations with high incidences of melanoma, juvenile males had a significantly higher frequency of the spotted caudal trait than that of adult males (two-sample z test, both P < 0.02) (Fig. 4A). By contrast, this pattern was not observed in the X. birchmanni parental population or in a hybrid population with a low incidence of melanoma (Fig. 4A). Phenotype tracking of laboratory-raised individuals shows that once it appears, the spotted area always expands over time, indicating that we do not expect reversal of spotting due to some form of phenotypic plasticity (Fig. 1D) (24). We also did not find evidence for systematic shifts in ancestry genome-wide between the juvenile and adult male life stages that could explain this pattern (24). However, we do see a shift toward X. birchmanni ancestry at the melanoma risk locus (in the top 1% of changes genome-wide) (fig. S21) (24).

Fig. 4 Impact of the spotted caudal melanoma in natural hybrid populations.

(A) Frequency of spotting in juvenile and adult males across populations with high (circles, Calnali low and Chahuaco falls) or low (squares, Aguazarca and X. birchmanni) melanoma incidence. Asterisks indicate significant differences by age class (*P < 0.05, **P < 0.01; ns indicates nonsignificant differences in a two-sample z test). Gray points indicate the raw data, black points indicate the mean, and error bars indicate one standard error of the mean. (B) Results of approximate Bayesian computation simulations indicate that the change in frequency of the spotting phenotype between juvenile and adult males is consistent with strong viability selection (24). Shown here are posterior distributions of viability selection coefficients consistent with the observed frequency change data in (left) Chahuaco falls and (right) Calnali low. (C) Because of where the melanoma develops, it can cause (top) the degradation of a fin essential in swimming or (bottom) the growth of tumors on the fin (overhead and side view of the same individual). (D) Visualization of the difference in fast-start response between individuals with low and high melanoma invasion (upper and lower 25% quantiles shown here). This representation is for visualization only; the statistical analysis comes from a linear model.

The observed shifts in spotting phenotype and ancestry at the melanoma risk locus, combined with an absence of substantial genome-wide shifts in ancestry, suggest that viability selection acts against spotted hybrids during maturation (24). Using an approximate Bayesian approach, we inferred that the strength of viability selection required to generate observed phenotypic shifts was extremely high and consistent across the two hybrid populations where melanoma is common (maximum a posteriori estimate of s ~ 0.2; 95% confidence intervals 0.05 to 0.44 and 0.04 to 0.38) (Fig. 4B).

Histology showed degradation of the muscle tissue that connects to the caudal fin in advanced melanomas (Figs. 1E and 4C, fig. S5, and movie S1). We thus measured its impact on swimming performance using two approaches. We did not find differences between phenotypes in ability to swim against a current (24). However, individuals with three-dimensional melanoma had slower escape responses when startled (linear model, t = –2.6, p = 0.014) (Fig. 4D) (24). This result is intriguing because fish with expanded spotting are likely more visible (movie S2), which could affect detection by avian and piscine predators.

Given the evidence for reduced survival of spotted individuals in populations with high rates of melanoma, it is surprising that this trait is still segregating in some hybrid populations (Fig. 4, A and B, and fig. S22) (24). Simulations suggest that high levels of gene flow from X. birchmanni would be required to maintain spotting at observed frequencies in hybrid populations (fig. S22) (24). However, because our inferences are based on viability rather than direct measures of fitness, we stress that there may be weaker effects of melanoma on overall fitness. Alternatively, other factors, such as mating advantages in individuals with large spots (40, 41), may explain its maintenance.


The involvement of xmrk in a melanoma hybrid incompatibility in two distantly related swordtail species pairs raises the question of whether certain genetic interactions are particularly prone to breakdown in hybrids. Genes that interact with many other genes or those that are involved in evolutionary arms races may be especially likely to generate hybrid incompatibilities (such as observed in Arabidopsis) (42, 43). Indeed, the only other known hybrid incompatibility in vertebrates, the recombination hotspot regulator prdm9, causes hybrid sterility in multiple crosses in mice (44). Whether unifying molecular or evolutionary forces drive the evolution of hybrid incompatibilities will become clearer as more incompatibilities are mapped to the single-gene level.

Supplementary Materials

Materials and Methods

Figs. S1 to S41

Tables S1 to S3

Appendix S1

References (51108)

Movies S1 and S2

References and Notes

  1. Materials and methods are available as supplementary materials.
Acknowledgments: We thank E. Calfee; H. Fraser; B. Kim; M. Lipson; P. Moorjani; M. Przeworski; and members of the Reich, Rosenthal, and Schumer laboratories for feedback on this work. We appreciate the help of A. Moyaho-Martinez in the design of swim trials, J. Lim and members of the Rosenthal laboratory in collecting samples, and O. Juárez-Mora for help conducting swim trials. We thank the federal government of Mexico for permission to collect fish. Stanford University and the Stanford Research Computing Center provided computational support for this project. Funding: This work was supported by NSF LTREB 1354172 to G.G.R.; funding from the Hagler Institute for advanced study to M.Scha.; and a Hanna H. Gray fellowship, L’Oreal for Women in Science grant, and NIH 1R35GM133774 grant to M.Schu. Author contributions: D.L.P., M.G.-O., and M.Schu. designed the project; D.L.P., M.G.-O., M.K., A.P.D.-L., S.B., D.B., M.Schu., and M.Scha. collected data; D.L.P., M.G.-O., M.K., K.D., P.R., and M.Schu. performed analyses; D.R., P.A., M.Scha., and G.G.R. provided expertise and technical support. Competing interests: The authors declare no competing interests. Data and materials availability: Data are available through NCBI (BioProject PRJNA610049; SRA accessions: SRX7847847-SRX7847871, SRX7866838-SRX7867011, SRX7861514-SRX7861761, and SRX7860174-SRX7860180) and Dryad (47). Code is available on github ( and archived at Zenodo (4850).

Stay Connected to Science

Navigate This Article