A lethal mitonuclear incompatibility in complex I of natural hybrids

The evolution of reproductive barriers is the first step in the formation of new species and can help us understand the diversification of life on Earth. These reproductive barriers often take the form of hybrid incompatibilities, in which alleles derived from two different species no longer interact properly in hybrids1–3. Theory predicts that hybrid incompatibilities may be more likely to arise at rapidly evolving genes4–6 and that incompatibilities involving multiple genes should be common7,8, but there has been sparse empirical data to evaluate these predictions. Here we describe a mitonuclear incompatibility involving three genes whose protein products are in physical contact within respiratory complex I of naturally hybridizing swordtail fish species. Individuals homozygous for mismatched protein combinations do not complete embryonic development or die as juveniles, whereas those heterozygous for the incompatibility have reduced complex I function and unbalanced representation of parental alleles in the mitochondrial proteome. We find that the effects of different genetic interactions on survival are non-additive, highlighting subtle complexity in the genetic architecture of hybrid incompatibilities. Finally, we document the evolutionary history of the genes involved, showing signals of accelerated evolution and evidence that an incompatibility has been transferred between species via hybridization.

The evolution of reproductive barriers is the first step in the formation of new species and can help us understand the diversification of life on Earth.These reproductive barriers often take the form of hybrid incompatibilities, in which alleles derived from two different species no longer interact properly in hybrids [1][2][3] .Theory predicts that hybrid incompatibilities may be more likely to arise at rapidly evolving genes [4][5][6] and that incompatibilities involving multiple genes should be common 7,8 , but there has been sparse empirical data to evaluate these predictions.Here we describe a mitonuclear incompatibility involving three genes whose protein products are in physical contact within respiratory complex I of naturally hybridizing swordtail fish species.Individuals homozygous for mismatched protein combinations do not complete embryonic development or die as juveniles, whereas those heterozygous for the incompatibility have reduced complex I function and unbalanced representation of parental alleles in the mitochondrial proteome.We find that the effects of different genetic interactions on survival are non-additive, highlighting subtle complexity in the genetic architecture of hybrid incompatibilities.Finally, we document the evolutionary history of the genes involved, showing signals of accelerated evolution and evidence that an incompatibility has been transferred between species via hybridization.
Biologists have long been fascinated by the question of how new species are formed and the mechanisms that maintain isolation between them.One key factor in the formation and maintenance of new species is the emergence of genetic incompatibilities that reduce viability or fertility in hybrids 1 .As originally described by the Dobzhansky-Muller model of hybrid incompatibility 2,3 , the unique sets of mutations that accumulate in diverging species may interact poorly when they are brought together in hybrids, given that they have never been tested against one another by selection.Owing to the technical challenges of identifying these interactions 4 , only around a dozen genes involved in hybrid incompatibilities have been precisely mapped 5 and exploration of the functional and evolutionary causes of hybrid incompatibilities has been limited to a small number of cases in model organisms 4 .
This knowledge gap leaves key predictions about the evolutionary processes that drive the emergence of hybrid incompatibilities untested.For one, theory suggests that incompatibilities should be more common within dense gene networks, both because genes involved in such interactions are expected to be tightly co-evolving and because the number of potentially incompatible genotypes explodes as the complexity of the genetic interaction increases 7,8 .Consistent with this prediction, mutagenesis experiments have highlighted the sensitivity of multi-protein interactions to changes in any of their components 8 .However, genetic interactions are notoriously difficult to detect empirically 9 , and this problem is exacerbated with complex genetic interactions 10 .Such technical challenges may explain the rarity of incompatibilities involving three or more genes in the empirical literature 8 (but see refs.9,11-14).
Another open question is the degree to which the genes that become involved in hybrid incompatibilities are predictable from their molecular or evolutionary properties.Researchers have proposed that rapid molecular evolution will increase the rate at which incompatibilities accumulate between species [4][5][6] .Although several known incompatibilities involve genes showing signatures of positive selection, it is unclear how unusual rates of protein evolution are in these genes relative to the genomic background 5,6 .The mitochondrial genome, in particular, has been proposed as a hotspot for the accumulation of Article genetic incompatibilities 15,16 , owing to substitution rates up to 25 times higher than the nuclear genome in many animals 17,18 and the potential for sexually antagonistic selection driven by its predominantly maternal inheritance 19,20 , among other factors 21 .At the same time, nuclear and mitochondrial proteins must interact with each other in key steps of ATP synthesis, increasing the likelihood of coevolution between these genomes 22,23 .These factors suggest that interactions between mitochondrial-and nuclear-encoded proteins could have an outsized role in the emergence of hybrid incompatibilities 15 , consistent with results from numerous species [24][25][26] .
As we begin to identify the individual genes underlying hybrid incompatibilities, the next frontier is evaluating the processes that drive their evolution.Over the past two decades, it has become abundantly clear that hybridization is exceptionally common in species groups where it was once thought to be rare 27,28 .As a result, it is now appreciated how frequently species derive genes from their relatives [29][30][31] .The effects of historical hybridization on the evolution of hybrid incompatibilities have been poorly investigated 32 , since the foundational theory in this area was developed before the ubiquity of hybridization was fully appreciated 7 .
Here we use an integrative approach to precisely map the genetic basis and physiological effects of a lethal mitonuclear hybrid incompatibility in swordtail fish and uncover its evolutionary history.The sister species Xiphophorus birchmanni and Xiphophorus malinche began hybridizing approximately 100 generations ago in multiple river systems 33 after premating barriers were disrupted by habitat disturbance 34 , and are a powerful system to study the emergence of hybrid incompatibilities in young species.Despite their recent divergence 35 (around 250,000 generations; 0.5% divergence per basepair), some hybrids between X. birchmanni and X. malinche experience strong selection against incompatibilities 35,36 .One incompatibility that causes melanoma has been previously mapped in this system and population genetic patterns suggest that dozens may be segregating in natural hybrid populations [35][36][37][38] .Moreover, the ability to generate controlled crosses 39,40 and the development of high-quality genomic resources 38,41 makes this system particularly tractable for studying hybrid incompatibilities in natural populations.Leveraging data from controlled laboratory crosses and natural hybrid populations, we pinpoint two nuclear-encoded X. birchmanni genes that are lethal when mismatched with the X. malinche mitochondria in hybrids, explore the developmental and physiological effects of this incompatibility, and trace its evolutionary history.

Mapping mitonuclear incompatibilities
To identify loci under selection in X. birchmanni × X. malinche hybrids, we generated approximately 1× low-coverage whole-genome sequence data for 943 individuals from an F 2 laboratory cross and 359 wild-caught hybrid adults, and applied a hidden Markov model to data at more than 600,000 ancestry-informative sites along the genome to infer local ancestry (approximately 1 informative site per kilobase 37,42 ; Methods and Supplementary Information 1.1.1-1.1.4).Using these results, we found evidence for a previously unknown incompatibility between the nuclear genome of X. birchmanni and the mitochondrial genome of X. malinche (Supplementary Information 1.1.5-1.1.10).Our first direct evidence for this incompatibility came from controlled laboratory crosses (Methods and Supplementary Information 1.1.1).Because the cross is largely unsuccessful in the opposite direction, all laboratory-bred hybrids were the offspring of F 1 hybrids generated between X. malinche females and X. birchmanni males and harboured a mitochondrial haplotype derived from the X. malinche parent species.Offspring of F 1 intercrosses are expected to derive on average 50% of their genome from each parent species.This expectation is satisfied genome wide and locally along most chromosomes in F 2 hybrids (on average 50.3%X. malinche ancestry; Supplementary Fig. 1).However, we detected six segregation distorters genome wide 40 , with the most extreme signals falling along a 6.5 Mb block of chromosome 13 and a 4.9 Mb block of chromosome 6 (Fig. 1a,d).
Closer examination of genotypes in the chromosome 13 region showed that almost none of the surviving individuals harboured homozygous X. birchmanni ancestry in a 3.75 Mb subregion (Fig. 1c and Supplementary Fig. 2; 0.1% observed versus 25% expected).This pattern is unexpected even in the case of a lethal incompatibility involving only nuclear loci (see Supplementary Information 1.1.1),but is consistent with a lethal mitonuclear incompatibility.Using approximate Bayesian computation (ABC) approaches we inferred the strength of selection against X. birchmanni ancestry in this region that was consistent with the observed genotypes and ancestry deviations.We estimated posterior distributions of selection and dominance coefficients and inferred that selection on this genotype in F 2 is largely recessive and essentially lethal (maximum a posteriori estimate h = 0.12 and s = 0.996, 95% credible interval h = 0.010-0.194and s = 0.986-0.999;Fig. 1b, Extended Data Fig. 1, Methods and Supplementary Information 1.2.1-1.2.2).
The degree of segregation distortion observed in F 2 individuals on chromosome 6 is also surprising (Fig. 1d).Only 3% of individuals harbour homozygous X. birchmanni ancestry in this region (compared with 0.1% in the chromosome 13 region and 25% on average at other loci across the genome; Fig. 1f), which is again lower than expected for a nuclear-nuclear hybrid incompatibility (Supplementary Information 1.1.1).ABC approaches indicate that selection on homozygous X. birchmanni ancestry on chromosome 6 is also severe (maximum a posteriori estimate h = 0.09 and s = 0.91, 95% credible interval interval h = 0.01-0.21and s = 0.87-0.94;Fig. 1e, Extended Data Fig. 1 and Supplementary Information 1.2.2).Thus, our F 2 data show that homozygous X. birchmanni ancestry in regions of either chromosome 13 or chromosome 6 is almost completely lethal in hybrids with X. malinche mitochondria (Fig. 1h).
To formally test for the presence of a mitonuclear incompatibility involving chromosome 13 and chromosome 6, or elsewhere in the genome, we leveraged data from natural hybrid populations.Most naturally occurring X. birchmanni × X. malinche hybrid populations are fixed for mitochondrial haplotypes from one parental species (Supplementary Information 1.1.2and 1.1.6).However, a few populations segregate for the mitochondrial genomes of both parental types, and we focused on one such population (the 'Calnali low' population, hereafter referred to as the admixture mapping population).Admixture mapping for associations between nuclear genotype and mitochondrial ancestry (after adjusting for expected covariance due to genome-wide ancestry 36 ) revealed two genome-wide significant peaks and one peak that approached genome-wide significance (Fig. 1g and Supplementary Tables 1-3).The strongest peak of association spanned approximately 77 kb and fell within the region of chromosome 13 identified using F 2 crosses (Fig. 1g).This peak was also replicated in another hybrid population (Methods, Supplementary Fig. 3 and Supplementary Information 1.1.5)and contains only three genes: the NADH dehydrogenase ubiquinone iron-sulfur protein 5 (ndufs5), E3 ubiquitin-protein ligase, and microtubule-actin cross-linking factor 1. Of these three genes, only ndufs5 forms a protein complex with mitochondrially encoded proteins, which along with other evidence implicates it as one of the interacting partners in the mitonuclear incompatibility (Fig. 1c, Extended Data Fig. 2 and Supplementary Fig. 4; see Supplementary Information 1.1.6-1.1.9).
We also identified a peak on chromosome 6 that approached genome-wide significance (Fig. 1g, Supplementary Fig. 5, Supplementary Table 2 and Supplementary Information 1.1.10)and fell precisely within the segregation distortion region previously mapped in F 2 hybrids (Fig. 1d and Supplementary Information 1.1.1).This peak contained 20 genes, including the mitochondrial complex I gene ndufa13 (Extended Data Fig. 2, Supplementary Fig. 5, Supplementary Information 1.1.10and Methods).Depletion of non-mitochondrial Nature | Vol 626 | 1 February 2024 | 121 parent ancestry at ndufa13 was unidirectional (Fig. 1f), consistent with selection acting only against the combination of the X. malinche mitochondria with homozygous X. birchmanni ancestry at ndufa13 (see Supplementary Information 1.2.3-1.2.4).Genomic analyses in natural hybrid populations confirmed this asymmetry (Extended Data Fig. 2).
Together, these results indicate that at least two X. birchmanni nuclear genes cause incompatibility when they are mismatched in ancestry with the X. malinche mitochondria (Fig. 1h and Supplementary Information 1.2.5).These genes, ndufs5 and ndufa13, belong to a group of proteins and assembly factors that form respiratory complex I (ref.43) (see Supplementary Table 1 for locations of the 51 annotated complex I genes in the Xiphophorus genome).Complex I is the first component of the mitochondrial electron transport chain that ultimately enables the cell to generate ATP.Both nuclear proteins interface with several mitochondrially derived proteins at the core of the complex I structure, hinting at the possibility that physical interactions could underlie this multi-gene mitonuclear incompatibility.

Interactions with X. birchmanni mitochondrial DNA
Admixture mapping analysis also identified a strong peak of mitonuclear association on chromosome 15, which we briefly discuss here and in Supplementary Information 1.1.10and 1.2.1.This peak was associated with X. birchmanni mitochondrial ancestry (Extended Data Fig. 3), indicating that it has a distinct genetic architecture from the incompatibility involving the X. malinche mitochondria and X. birchmanni ndufs5 and ndufa13.Specifically, analysis of genotypes at the admixture mapping peak indicates that X. birchmanni mitochondrial ancestry is incompatible with homozygous X. malinche ancestry on chromosome 15 (Fig. 1c and Extended Data Fig. 3).This region did not contain any members of complex I, but dozens of genes in this The region shaded blue shows the 99% quantiles of X. malinche ancestry at all ancestry informative sites genome wide.The dashed line represents expected X. malinche ancestry for this cross.The purple arrow points to the position of ndufs5.b, Results of ABC simulations estimating the strength of selection on X. malinche mitochondria combined with X. birchmanni ancestry at ndufs5.Shown is the posterior distribution from accepted simulations; the vertical line indicates the maximum a posteriori estimate for selection coefficient (s).c, Observed genotype frequencies of different genotype combinations of ndufs5 and mitochondrial haplotypes in admixture mapping population.d, Average ancestry of F 2 individuals on chromosome 6 reveals a large region of segregation distortion towards X. malinche ancestry.The region shaded blue shows the 99% ancestry quantiles and expected ancestry as in a, and the red arrow points to position of ndufa13.e, Results of ABC simulations estimating the strength of selection on X. malinche mitochondria combined with X. birchmanni ndufa13, as in b. f, Observed genotype frequencies of different genotype combinations of ndufa13 and mitochondrial haplotypes in the admixture mapping population.g, Admixture mapping of associations between nuclear ancestry and mitochondrial haplotype in natural hybrids using a partial correlation approach, controlling for genome-wide ancestry.The blue line indicates the 10% false-positive rate genome-wide significance threshold determined by simulations.The peak visible on chromosome 15 is driven by interactions with the X. birchmanni mitochondria and an unknown nuclear gene (Supplementary Information 1.1.10and 1.2.1).h, Schematic of identified interactions with the X. malinche mitochondrial genome from mapping data.The dashed line indicates a subtle interaction between ndufs5 and ndufa13 (see text, Fig. 2 and Supplementary Information 1.2.5).

Article
interval interact with known mitonuclear genes (see Supplementary Table 3 and Supplementary Information 1.1.10).The fact that we detect incompatible interactions with both the X. malinche mitochondria (at ndufs5 and ndufa13) and the X. birchmanni mitochondria (ndufs5 and chromosome 15) in our admixture mapping results supports the idea that mitonuclear interactions can act as 'hotspots' for the evolution of hybrid incompatibilities 15 .

Lethal effects in early development
The combination of X. birchmanni ndufs5 or ndufa13 with the X. malinche mitochondria appears to be lethal by the time individuals reach adulthood.To investigate the developmental timing of the incompatibility, we genotyped pregnant females from the admixture mapping population and recorded the developmental stages of their embryos 44 (swordtails are livebearing fish; Methods).We found a significant interaction between developmental stage and ndufs5 genotype, whereas ndufa13 genotype did not affect developmental stage (Fig. 2a,b, Supplementary Figs.6-9 and Supplementary Information 1.3.1).Genotyping results revealed that embryos with homozygous X. birchmanni ancestry at ndufs5 and X. malinche mitochondria are present at early developmental stages, but that these embryos did not develop beyond a phenotype typical of the first seven days of gestation (the full length of gestation is 21-28 days in Xiphophorus; Fig. 2a,b,d,e).Individuals with mismatched ancestry at ndufs5 whose siblings were fully developed still had a detectable heartbeat but had consumed less yolk than their siblings and remained morphologically underdeveloped (Fig. 2d, Extended Data Fig. 4 and Supplementary Figs.10-14).Unlike other species, in Xiphophorus this developmental lag could itself cause mortality, since embryos that do not complete embryonic development inside the mother do not survive more than a few days after birth (Supplementary Information 1.3.1 and Supplementary Table 4).Given that complex I inhibition lethally arrests development in zebrafish embryos 45,46 , we also tested the effects of complex I inhibition on X. birchmanni and X. malinche fry, and found a similar level of sensitivity (Supplementary Information 1.3.2).
In contrast to individuals with mismatched ancestry at ndufs5, those with ndufa13 mismatch survived embryonic development but suffered mortality in the early post-natal period (Fig. 2c).We tracked 74 F 2 fry from 24 h post birth to adulthood (Supplementary Information 1.3.3).We found that most fry with incompatible genotypes at ndufa13 had already suffered mortality by the time tracking began, with only 7 individuals found 24 h post birth that were homozygous X. birchmanni at ndufa13 (versus 19 expected; binomial P = 0.0005).No natural mortality was observed between 1 day and 3 months post birth (Supplementary Information 1.3.3).

Physiology and complex fitness effects
Our analysis of developing embryos indicates that individuals with the ndufs5 incompatibility exhibited abnormal embryonic development, whereas those with the ndufa13 incompatibility did not.This suggests that these genes may drive lethality through partially distinct mechanisms.Thus, we chose to further investigate the effects of ndufs5 and ndufa13, and their possible interactions.We sampled 235 F 2 embryos at a range of developmental stages and measured their overall rates of respiration (Supplementary Information 1.3.4-1.3.5).We also used imaging of these embryos to track cardiovascular phenotypes as these have been associated with ndufa13 defects in mammals 47 .We found that incompatible genotypes at ndufs5 and ndufa13 affected a range of phenotypes, including heart rate, length relative to compatible siblings, and length-corrected head size (Extended Data Figs. 4 and 5, Supplementary Figs.11-17 and Supplementary Tables 5-7).ndufa13 mismatch has a large effect on cardiovascular phenotypes, including heart rate and the size of the sinu-atrium (an embryo-specific heart chamber; Fig. 2g,h, Supplementary Figs. 14 and 16 and Supplementary Tables 6 and 8), whereas ndufs5 affects only heart rate (Supplementary Fig. 14 and Supplementary Table 6).We find initial evidence that cardiac defects persist into adulthood in surviving individuals with ndufa13 mismatch (Extended Data Fig. 5 and Supplementary Information 1.3.3).By contrast, ndufs5 mismatch has a major effect on rates of respiration and yolk consumption during development (Fig. 2f, Supplementary Figs.18-20 and Supplementary Tables 9-11).
Naively, the separable impacts of incompatible genotypes at ndufs5 and ndufa13 could indicate that even though these proteins are in physical contact in complex I (see below), they represent two distinct hybrid incompatibilities.We investigated this question by taking advantage of rare survivors of the ndufa13 incompatibility in an expanded dataset of 1,010 F 2 hybrids.Using this dataset, we were able to identify dozens of survivors of the ndufa13 incompatibility (3.4% of individuals) and found that genotypes at ndufa13 and ndufs5 were not independent (χ 2 association test P = 0.032; Supplementary Information 1.2.5).Upon further investigation, we found that the majority of survivors of the ndufa13 incompatibility had homozygous X. malinche ancestry at ndufs5, suggesting that harbouring even one X.birchmanni allele at ndufs5 may sensitize fry to the ndufa13 incompatibility.Indeed, we found that individuals that had heterozygous ancestry at ndufs5 were significantly under-represented among surviving ndufa13 incompatible individuals (Permutation test P = 0.015; Fig. 2i and Supplementary Information 1.2.5).These findings highlight a subtle but significant non-additive effect of ndufs5 and ndufa13 on survival.

Mitochondrial biology in heterozygotes
Because few individuals homozygous for incompatible genotypes at ndufs5 or ndufa13 survive past birth, our previous experiments focused on embryos.However, the small size of Xiphophorus embryos prevents us from using assays that directly target complex I.To further explore the effects of the hybrid incompatibility on complex I function in vivo, we turned to adult F 1 hybrids (Fig. 3a).Since F 1 hybrids that derive their mitochondria from X. malinche and are heterozygous for ancestry at ndufs5 and ndufa13 are fully viable, we tested whether there was evidence for compensatory regulation that might be protective in F 1 hybrids.We found no evidence for significant differences in expression of ndufs5 or ndufa13 (Supplementary Information 1.3.6 and Supplementary Figs.21 and 22) or in mitochondrial copy number (Fig. 3b and Supplementary Information 1.3.7) between F 1 hybrids and parental species.
With no indication of a compensatory regulatory response, we reasoned that we might be able to detect reduced mitochondrial complex I function in hybrids heterozygous for ancestry at ndufs5 and ndufa13.We quantified respiratory phenotypes in isolated mitochondria using an Oroboros O2K respirometer in adult hybrids and parental species (Methods, Supplementary Fig. 23 and Supplementary Information 1.3.8).We found that complex I efficiency was lower in hybrids (Fig. 3c and Supplementary Fig. 24, orthogonal contrast t = −2.53,P = 0.023, n = 7 per genotype), and that the time required for hybrids to reach maximum complex I-driven respiration was around 2.5 times longer (orthogonal contrast t = 4.303, P < 0.001; Fig. 3d and Supplementary Fig. 25).Conversely, overall levels of mitochondrial respiration were unaffected by genotype (Fig. 3e, orthogonal contrast t = 0.078, P = 0.94, n = 7 per genotype; Supplementary Information 1.3.8)as were other measures of mitochondrial integrity and function (Supplementary Figs. 26 and 27 and Supplementary Information 1.3.8 and 1.3.9).Together, these data point to reduced function of complex I without broader phenotypic consequences in individuals that are heterozygous for incompatible alleles 48 .
Given the physiological evidence for some reduction in complex I function in hybrids heterozygous at ndufs5 and ndufa13, we predicted that there might be an altered frequency of protein complexes incorporating both X. malinche mitochondrial proteins and X. birchmanni proteins at ndufs5 and ndufa13 in F 1 hybrids.To test this prediction, we took a mass spectrometry-based quantitative proteomics approach.We used stable isotope-labelled peptides to distinguish between the X. birchmanni and X. malinche ndufs5 and ndufa13 peptides in mitochondrial proteomes extracted from F 1 hybrids (n = 5; see Methods and Supplementary Information 1.4.1-1.4.4).Although endogenous  -9. f-h, Respiratory and morphometric measurements in F 2 embryos.To control for length, residuals of each variable regressed against length are plotted (Supplementary Information 1.3.4 and 1.3.5).Grey dots denote individual measurements, coloured points and bars show mean ± 2 × s.e.m. for each genotype, and brackets with asterisks denote significant differences from Tukey's honest significant difference test.f, Relationship between ndufs5 genotype and respiration rate in hybrids (n = 40 X. birchmanni, 102 heterozygotes, 47 X. malinche).Xb, X. birchmanni.g, Relationship between ndufa13 genotype and heart rate in hybrids (n = 39 X. birchmanni, 95 heterozygotes, 46 X. malinche).h, Relationship between ndufa13 genotype and width of sinu-atrium, a peristaltic canal between the yolk and embryonic atrium, in hybrids (n = 37 X. birchmanni, 82 heterozygotes, 42 X.malinche).i, Frequency of heterozygotes at ndufs5 among juveniles and adults of varying ndufa13 genotypes.Dots and lines represent observed frequency ± 2 × s.e.m.The dashed line represents expected frequency under additive selection (Supplementary Information 1.2.5).Article ndufa13 peptides were not observed frequently enough to quantify accurately, we found consistent deviations from the expected 50:50 ratio of X. birchmanni to X. malinche peptides for ndufs5 in F 1 hybrids, with a significant overrepresentation of matched ancestry at ndufs5 in the mitochondrial proteome (t = 3.96, P = 0.016; Fig. 3f, Supplementary Fig. 28 and Supplementary Information 1.4.5).Since we did not observe allele-specific expression of ndufs5 (Fig. 3f and Supplementary Information 1.3.6),this result is consistent with disproportionate degradation of X. birchmanni-derived ndufs5 peptides in the mitochondrial proteome or differences in translation of ndufs5 transcripts derived from the two species.

Mitonuclear substitutions in complex I
To begin to explore the possible mitochondrial partners of ndufs5 and ndufa13 among the 37 non-recombining genes in the swordtail mitochondrial genome, we turned to protein modelling, relying on high-quality cryo-electron microscopy (cryo-EM)-based structures [49][50][51] .Although these structures are only available for distant relatives of swordtails, the presence of the same set of supernumerary complex I subunits and high sequence similarity suggest that using these structures is appropriate (Supplementary Tables 12 and 13, Supplementary Figs.29-31 and Supplementary Information 1.4.6).
Barring a hybrid incompatibility generated by regulatory divergence (see Supplementary Information 1.3.6),our expectation is that hybrid incompatibilities will be driven by amino acid changes in interacting proteins 52 .We used the program RaptorX 53 to generate predicted structures of X. birchmanni and X. malinche Ndufs5, Ndufa13 and nearby complex I proteins encoded by mitochondrial and nuclear genes, which we aligned to a mouse cryo-EM complex I structure 49 (Fig. 4a, Supplementary Figs.29-31 and Methods).Using these structures, we visualized amino acid substitutions between X. birchmanni and X. malinche at the interfaces of Ndufs5, Ndufa13 and mitochondrial-encoded proteins (Extended Data Fig. 6 and Supplementary Figs.32 and 33).Whereas there are dozens of substitutions in the four mitochondrial-encoded proteins that are in close physical proximity to Ndufs5 or Ndufa13 (Supplementary Fig. 29; Nd2, Nd3, Nd4l and Nd6), there are only five cases where amino acid substitutions in either nuclear-encoded protein are predicted to be close enough to contact substitutions in any mitochondrial-encoded protein, all of which involve Nd2 or Nd6 (Fig. 4a and Extended Data Table 1; see Supplementary Fig. 33 for pairwise visualizations of interacting proteins).These paired substitutions in regions of close proximity between mitochondrial-and nuclear-encoded proteins suggest that nd2 and nd6 are the genes most likely to be involved in the mitochondrial component of the hybrid incompatibility (Fig. 4a,b Extended Data Fig. 6 and Supplementary Figs.33-35), and will be promising candidates for functional validation when such approaches become possible in swordtails.

Rapid evolution of complex I proteins
Theory predicts that hybrid incompatibilities are more likely to arise in rapidly evolving genes [4][5][6][7] .Consistent with this hypothesis, ndufs5 is among the most rapidly evolving genes genome-wide between X. birchmanni and X. malinche (Fig. 4c,d).Aligning the ndufs5 coding sequences of X. birchmanni, X. malinche and 12 other swordtail species revealed that all 4 amino acid substitutions that differentiate X. birchmanni and X. malinche at ndufs5 were derived on the X. birchmanni branch (Fig. 4c).Phylogenetic tests indicate that there has been accelerated evolution of ndufs5 on this branch (inferred ratio of non-synonymous substitutions per non-synonymous site to synonymous substitutions per synonymous site (d N /d S ) > 99, N = 4, S = 0, codeml branch test P = 0.005; Fig. 4c).Similar patterns of rapid evolution are observed at ndufa13, which also showed evidence for accelerated evolution in X. birchmanni (Fig. 4e; d N /d S = 1.2, N = 3, S = 1, codeml branch test P = 0.002).Although explicit tests for adaptive evolution at ndufs5 and ndufa13 could not exclude a scenario of relaxed selection (Extended Data Table 2 and Supplementary Information 1.5.1 and 1.5.2),our comparisons across phylogenetic scales highlight strong conservation in some regions of  3 | Physiology and proteomics of viable heterozygotes.a, Schematic of ancestry at loci of interest in X. birchmanni × malinche F 1 hybrids (the incompatibility is not lethal in F 1 ; Supplementary Information 1.2.2). b, Real-time quantitative PCR analysis of differences in mitochondrial copy number in liver tissue by genotype.The ratio is derived from the difference in C t between the single-copy nuclear gene Nup43 and the mitochondrial gene nd1 (n = 3 X.birchmanni, 6 F 1 , 6 X. malinche).c-e, Results of Oroboros O2K respirometer assay for adult X. birchmanni, X. malinche and hybrid individuals (n = 7 fish per genotype), with X. malinche mitochondria and heterozygous ancestry at ndufs5 and ndufa13.c, Complex I efficiency in F 1 hybrids and parental species.d, Time to reach the maximum rate of complex I-driven respiration after ADP addition in F 1 hybrids and parental species (see Supplementary Fig. 25 for example time-to-peak curves).For corresponding analysis of time to reach peak complex I-and complex II-driven respiration after addition of succinate, see Supplementary Fig. 26.e, Maximum respiration rates during full O2K protocol in F 1 hybrids and parental species.f, Allelic balance of ndufs5 in the F 1 hybrid transcriptome and proteome.Left, allele-specific expression of ndufs5 in adult F 1 hybrids (n = 3 fish).Right, results of quantitative mass spectrometry analysis of ndufs5 peptides in mitochondrial proteomes derived from adult F 1 hybrids (n = 5 fish).Dots show the proportion of area under the spectral curve contributed by the X. malinche allele in each individual.The data on the left are for endogenous peptides present in F 1 hybrids, and the data on the right are results for heavy isotope-labelled peptide controls.The peptides quantified for each species are shown on the graph.a-f, Coloured dots and bars show the mean ± 2 × s.e.m., and grey dots show individual data.the proteins and rapid turnover in others, complicating our interpretation of this test (Supplementary Fig. 36).
Rapid evolution of ndufs5 and ndufa13 could be driven by coevolution with mitochondrial substitutions, a mechanism that has been proposed to explain the outsized role of the mitochondria in hybrid incompatibilities 15,54 .Indeed, there is an excess of derived substitutions in the X. birchmanni mitochondrial protein Nd6, one of the proteins that physically contacts Ndufs5 and Ndufa13 (Extended Data Fig. 7 and Extended Data Table 2; codeml branch test P = 0.005).Moreover, several of the substitutions observed in both mitochondrial and nuclear genes are predicted to have functional consequences (Extended Data Table 3 and Supplementary Information 1.5.1),including ones predicted to be in contact between Ndufs5, Ndufa13, Nd2 and Nd6 (Fig. 4a,b and Extended Data Fig. 6).

Introgression of incompatibility genes
The presence of a mitonuclear incompatibility in Xiphophorus is especially intriguing, given previous reports that mitochondrial genomes may have introgressed between species 29 .While X. malinche and X. birchmanni are sister species based on the nuclear genome, they are mitochondrially divergent, with X. malinche and Xiphophorus cortezi grouped as sister species based on the mitochondrial phylogeny 29 (Fig. 5a,b).As we show, all X. cortezi mitochondria sequenced to date are nested within X. malinche mitochondrial diversity (Fig. 5b, Supplementary Fig. 37 and Supplementary Information 1.5.3 and 1.5.4).Simulations indicate that gene flow, rather than incomplete lineage sorting, drove replacement of the X. cortezi mitochondria with the X. malinche sequence (P < 0.002 by simulation; Fig. 5c and Supplementary Information 1.5.4).
The introgression of the mitochondrial genome from X. malinche into X.cortezi raises the possibility that other complex I genes may have co-introgressed 55 .Indeed, the nucleotide sequence for ndufs5 is identical between X. malinche and X. cortezi, and the sequence of ndufa13 differs by a single synonymous mutation (although conservation of both genes is high throughout Xiphophorus; Supplementary Figs.38  and 39).The identical amino acid sequences of the proteins suggest that hybrids between X. cortezi and X. birchmanni are likely to harbour the same mitonuclear incompatibility we observe between X. malinche and X. birchmanni, as a result of ancient introgression between X. malinche and X. cortezi (Fig. 5d and Supplementary Information 1.5.3-1.5.5).
This inference is supported by analysis of three contemporary X. birchmanni × X. cortezi hybrid populations 40 (Supplementary Fig. 40).We find that all known X. birchmanni × X. cortezi hybrid populations are fixed for the mitochondrial genome from X. cortezi (that originated in X. malinche) and show a striking depletion of X. birchmanni ancestry at ndufs5 and ndufa13 (Fig. 5e and Supplementary Fig. 41).This replicated depletion is not expected by chance (Fig. 5e and Supplementary Information 1.5.6,P = 0.0001) and instead indicates that selection has acted on these regions.These results suggest that the mitonuclear incompatibility observed in X. birchmanni × X. malinche is also active in hybridizing X. birchmanni × X. cortezi populations.This exciting finding shows that genes underlying hybrid incompatibilities can introgress together, transferring incompatibilities between related species.contacts between these proteins.Solid black lines highlight two areas predicted to be in close contact between interspecific substitutions (alpha carbon distance ≤ 10 Å for all models), and dashed lines show three additional areas in which there was weaker evidence for a predicted contact based on computational analyses (side chain distance ≤ 12 Å in at least one model).Asterisks denote residues with substitutions in X. birchmanni computationally predicted to affect protein function (Extended Data Table 3).b, Detailed view of the interface between Ndufs5, Nd2 and Nd6.Spheres highlight substitutions between X. birchmanni and X. malinche.For substitutions predicted to be in close proximity, residues are labelled with letters denoting the X. malinche allele, the residue number, and the X. birchmanni allele, respectively (Extended Data Table 3 and Supplementary Information 1.4.6).c, Gene tree for ndufs5 generated with RAxML, highlighting an excess of substitutions along the X. birchmanni branch.The scale bar represents the number of nucleotide substitutions per site.Derived non-synonymous substitutions are indicated by red ticks along the phylogeny; spacing between ticks is arbitrary.d, Distri bution of log 10 (d N /d S ) between X. birchmanni and X. malinche across all nuclear genes with values for ndufs5 and ndufa13 highlighted.e, Gene tree for ndufa13 (as in c), highlighting an excess of substitutions along the X. birchmanni branch.

Discussion
Here we investigate the genetic and evolutionary forces that drive the emergence of hybrid incompatibilities.Theory predicts that hybrid incompatibilities involving multiple genes should be common 7,8 , but with few exceptions 9,[11][12][13] , they remain almost uncharacterized at the genic level 8 .We have identified incompatible interactions in mitochondrial complex I that cause hybrid lethality in laboratory and wild populations.Our findings in naturally hybridizing species echo predictions from theory and studies in laboratory models 9,[11][12][13] suggesting that protein complexes may be a critical site of hybrid breakdown.
Researchers have proposed mitonuclear interactions as hotspots for the emergence of hybrid incompatibilities, given that mitochondrial genomes often experience higher substitution rates between species 17,18,56 , yet must intimately interact with nuclear proteins to perform essential cellular functions 22,23 .Our findings support this prediction, identifying incompatible interactions with both the X. malinche and X. birchmanni mitochondria.We also show that there has been exceptionally rapid evolution in both mitochondrial and interacting nuclear genes in X. birchmanni (Fig. 4).Whether driven by adaptation or some other mechanism, our findings support the hypothesis that the coevolution of mitochondrial and nuclear genes could drive the overrepresentation of mitonuclear interactions in hybrid incompatibilities 22,23,54 .More broadly, our results are consistent with predictions that rapidly evolving proteins are more likely to become involved in hybrid incompatibilities than their slowly evolving counterparts [4][5][6] .
Characterizing the incompatibility across multiple scales of organization enabled us to explore the mechanisms through which it acts [57][58][59] .Our results suggest that in the case of the X. malinche mitochondria hybrid lethality is mediated through arrested development in utero of individuals with mismatched ancestry at ndufs5, whereas individuals with ndufa13 mismatch have vascular defects and typically die shortly after birth.Intriguingly, individuals with ndufa13 mismatch that do survive are much less likely to harbour any X.birchmanni alleles at ndufs5 (Fig. 2i).Together, our results indicate that a subtle three-way interaction overlays two strong pairwise mitonuclear incompatibilities at ndufs5 and ndufa13.Evolutionary biologists have been fascinated by the idea that hybrid incompatibilities may commonly involve three or more genes following theoretical work by Orr 7 nearly 30 years ago, but this question has been challenging to address empirically.Our results highlight how the nuances of actual fitness landscapes may defy simplifying assumptions.
Finally, this mitonuclear incompatibility provides a new case in which the same genes are involved in incompatibilities across multiple species 30,38,60 .However, tracing the evolutionary history of the genes that underlie it adds further complexity to this phenomenon: we found that introgression has resulted in the transfer of genes underlying the incompatibility from X. malinche to X. cortezi, and evidence from X. birchmanni × X. cortezi hybrid populations indicates that the incompatibility is probably under selection in these populations as well.The possibility that hybridization could transfer incompatibilities between species has not been previously recognized, perhaps due to an underappreciation of the frequency of hybridization.The impact of past hybridization on the structure of present-day reproductive barriers between species is an exciting area for future inquiry.

Online content
Any methods, additional references, Nature Portfolio reporting summaries, source data, extended data, supplementary information, showing that X. birchmanni and X. malinche are sister species 29 .b, Phylogeny constructed from whole mitochondrial DNA sequences showing that X. cortezi mitochondria are nested within X. malinche mitochondrial diversity.c, Results of simulations modelling expected mitochondrial divergence between X. malinche and X. cortezi in a scenario with no gene flow.The first set of simulations used the average mitochondrial substitution rate observed between Xiphophorus species (red), and the second used the minimum mitochondrial substitution rate observed (teal).The dotted line shows observed divergence between mitochondrial haplotypes in X. malinche and X. cortezi, indicating that past mitochondrial introgression is more likely than incomplete lineage sorting.d, Clustal alignment of Ndufs5 sequences shows that X. malinche and X. cortezi (Xc) have identical amino acid sequences at Ndufs5, hinting at possible introgression of the nuclear ndufs5 gene, whereas X. birchmanni is separated from them by four substitutions.Similar patterns are observed for ndufa13.Colours indicate properties of the amino acid, and asterisks indicate locations where the amino acid sequences are identical.e, Non-mitochondrial parent ancestry is strongly depleted in natural X. cortezi × X. birchmanni hybrid populations fixed for the X. cortezi mitochondrial haplotype at ndufs5 and ndufa13.Upward arrowheads along x-axes show the locations of ndufa13 (red) and ndufs5 (purple) on chromosomes 6 and 13, respectively, which fall in minor parent ancestry deserts in both independently formed populations, as expected for strong hybrid incompatibilities.
Step curves show average X. birchmanni ancestry in 250-kb windows, and horizontal lines show the genome-wide average for each population.
Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made.The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material.If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.Left distribution shows the results of simulations modeling the X. malinche mitochondria -X.birchmanni ndufs5 component of the hybrid incompatibility, and right distribution shows the results of simulations modeling the X. malinche mitochondria -X.birchmanni ndufa13 component of the hybrid incompatibility.Each posterior distribution shows the accepted dominance coefficients in 500 simulations and the gray dashed line indicates the maximum a posteriori estimate for the dominance coefficients of interactions involving ndufs5 and ndufa13 (b) Posterior distributions of selection coefficients from ABC simulations fitting observed data in the Calnali Low hybrid population.Left distribution shows the results of simulations modeling the X. malinche mitochondria -X.birchmanni ndufs5 incompatibility and right distribution shows the results of simulations modeling the reverse X. birchmanni mitochondria -X.malinche ndufs5 incompatibility.In the main text we focus on inferred selection coefficients from F 2 hybrids for the incompatibility involving the X. malinche mitochondria (see Fig. 1) but present results from fitting the Calnali Low population data here and in Supplementary Information 1.2.2.Distribution shows the accepted selection coefficients in 500 simulations and the gray dashed line indicates the maximum a posteriori estimate.Note that we recover a much broader distribution and a lower maximum a posteriori estimate of the selection coefficient for the X. malinche mitochondria -X.birchmanni ndufs5 incompatibility here compared to F 2 simulations (Fig. 1).We interpret this result to be driven by the fact that simulations fitting the Calnali Low population data modeled >20 generations of admixture.Thus, a range of selection coefficients are consistent with the observation that no individuals have X.malinche mtDNA and homozygous X. birchmanni ancestry at ndufs5 in this population.(A) Non-mitochondrial parent ancestry is depleted around ndufs5 in all natural hybrid populations that are fixed for a particular mitochondrial haplotype.This pattern is strongly suggestive of a history of selection on this region in natural hybrid populations.In (B), ndufa13 is fixed for X. malinche ancestry only in the population with X. malinche mitochondria (Tlatemaco), mirroring the architecture of the genetic interaction between this gene and the mitochondria.Specifically, interactions with ndufa13 are only expected to be under selection in combination with the X. malinche mitochondrial haplotype (Fig. 1c).

Statistics
For all statistical analyses, confirm that the following items are present in the figure legend, table legend, main text, or Methods section.

n/a Confirmed
The exact sample size (n) for each experimental group/condition, given as a discrete number and unit of measurement A statement on whether measurements were taken from distinct samples or whether the same sample was measured repeatedly The statistical test(s) used AND whether they are one-or two-sided Only common tests should be described solely by name; describe more complex techniques in the Methods section.

A description of all covariates tested
A description of any assumptions or corrections, such as tests of normality and adjustment for multiple comparisons A full description of the statistical parameters including central tendency (e.g.means) or other basic estimates (e.g.regression coefficient) AND variation (e.g. standard deviation) or associated estimates of uncertainty (e.g.confidence intervals) For null hypothesis testing, the test statistic (e.g.F, t, r) with confidence intervals, effect sizes, degrees of freedom and P value noted Give P values as exact values whenever suitable.
For Bayesian analysis, information on the choice of priors and Markov chain Monte Carlo settings For hierarchical and complex designs, identification of the appropriate level for tests and full reporting of outcomes Estimates of effect sizes (e.g.Cohen's d, Pearson's r), indicating how they were calculated Our web collection on statistics for biologists contains articles on many of the points above.

nature portfolio | reporting summary
April 2023 birchmanni and 36 X. malinche wild-caught females, respectively.For analysis of embryo survival outside the womb, we used a total of 20 X. birchmanni fry from two females and 25 F1 fry from four females.In testing Xiphophorus fry sensitivity to Complex I inhibition, we used a total of 39 X. birchmanni fry from three broods, and 11 X. malinche fry from one brood.To construct de novo whole-mtDNA phylogenies, we used 5 X. birchmanni, 3 X.malinche, and 2 X. cortezi from Huichihuayan, San Luis Potosi, Mexico, and the Xiphophorus Stock Center at Texas State University.For mapping of the incompatibility between X. birchmanni and X. cortezi, we used 284 natural X. birchmanni x cortezi hybrids from Santa Cruz, Huextetitla, and Chapulhuacanito, San Luis Potosi, Mexico.

Sampling strategy
Sample sizes for the Powell et al.F2 cross were determined using a power analysis for a previous QTL mapping study, where sample size was chosen for an estimated 90% chance of detecting a QTL explaining 5% of phenotypic variance.Sample sizes for all individuals collected from natural populations were determined by sampling success during field seasons and the maximum permitted across all sampling sites by collection permits.Sample sizes for all analyses based on lab-born hybrids were either the total number available with the desired genotypes at the time of experimentation, or the maximum feasible while still preserving the viability of lab colonies (which can only periodically be refreshed with wild collections).

Data collection
Data collection procedures are listed in the main text and Supplement (being too extensive to list here).SMB, AED, and JJB performed wet lab work for genotyping, TRG and MS performed bioinformatics for genotyping, and BMM and QKL contributed to both.TRG performed Complex I pharmacological inhibition trials and all other measurements of fry survival.BMM, DLP, and SMB collected embryo stage data from wild-caught hybrids.BMM performed respirometry and morphometrics in live embryos, while ENKI and JCH performed Oroboros respirometry assays on isolated mitochondria.SMB performed qPCR for mitochondrial copy number, CYP performed differential expression and allele-specific expression analyses, and AM performed assays of mitochondrial membrane polarization.FL, RM, KS, and RDL performed mass spectrometry assays on mitochondrial isolates prepared by BMM and AED.RRS analyzed histology images of juvenile F2 hearts.BMM constructed all phylogenies and associated testing of evolutionary rates, performed structural modeling, and carried out coevolution tests.MS and BMM performed all simulation and statistical analyses.
Timing and spatial scale All natural hybrid samples and parental species were collected over three years of field trips to the CICHAZ field station in Calnali, Hidalgo, Mexico, from ten different sites within an area of Hidalgo and San Luis Potosi ~40 km in diameter.Sampling was roughly quarterly with the exception of the year 2020, when sampling was limited by the COVID-19 pandemic.
Data exclusions 18 individuals were excluded from whole embryo respirometry and morphometrics due to damage (punctured yolks, torn tissue, severe internal bleeding, etc) at the time of photography as identified by BMM; these damages could affect body measurements, and it was unclear when the damage occurred, such that a confounding effect on respirometry measurements could not be ruled out.Likewise, four heart compartment slides were excluded from histology analysis due to visible tearing of the tissue on the slide, which was deemed sufficient to affect cross-sectional area by BMM.These criteria were not pre-established, but the observer was unaware of the experimental treatment group when making the decision to exclude.

Reproducibility
Replicability was largely confirmed by comparing across biologically independent datasets and statistical approaches: the depletion of mismatched genotypes was repeated in segregation distortion in lab-born F2 hybrids, partial correlation analysis in two natural hybrid populations with both mitochondrial ancestries segregating, and permutation-based testing of depletion in three natural hybrid populations fixed for mitochondrial ancestry.Replicability of in silico protein modelling was tested by repeating structure prediction with five different random seeds and four different mammalian cryo-EM templates.Whole embryo respirometry was replicated across ten broods in two major sets separated by ~6 months.

Randomization
Randomization is not relevant to our study because all analyses either lacked experimental treatments or applied a single treatment to all individuals, with the genotypes of individual being the independent variable of interest.Where applicable, the covariate of genome-wide ancestry fraction was controlled for using partial correlation analysis, and any variation in the administration of respirometry methods were accounted for for by using date of testing as a blocking variable in downstream statistics.

Blinding
Blinding was not employed in this study.Many analyses featured only a single group, such that blinding was not relevant, and in others, the placement of individuals in experimental groups of interest (e.g.genotypes) was impossible until after all data collection and analysis had already occurred.In the case of Oroboros respirometry, treatment groups of individuals were visible to experimenters based on phenotype.

Did the study involve field work?
Field work, collection and transport

Field conditions
Collections were performed in September, November, January, February, March, May, and June in the Sierra y Huasteca region of Mexico, with sites ranging from ~10-30 C depending on the season, and rainfall varying from absent during the peak dry season (winter and spring) to frequent and substantial in the rainy season (summer and fall).

Fig. 1 |
Fig. 1 | Admixture mapping pinpoints mitonuclear incompatibility in Xiphophorus.a, Average ancestry of F 2 hybrids on chromosome (chr.)13 reveals a large region of segregation distortion towards X. malinche ancestry.The region shaded blue shows the 99% quantiles of X. malinche ancestry at all ancestry informative sites genome wide.The dashed line represents expected X. malinche ancestry for this cross.The purple arrow points to the position of ndufs5.b, Results of ABC simulations estimating the strength of selection on X. malinche mitochondria combined with X. birchmanni ancestry at ndufs5.Shown is the posterior distribution from accepted simulations; the vertical line indicates the maximum a posteriori estimate for selection coefficient (s).c, Observed genotype frequencies of different genotype combinations of ndufs5 and mitochondrial haplotypes in admixture mapping population.d, Average ancestry of F 2 individuals on chromosome 6 reveals a large region of segregation distortion towards X. malinche ancestry.The region shaded blue shows the 99% ancestry quantiles and expected ancestry as in a, and the red

Fig. 2 |
Fig. 2 | Effect of incompatibility on Xiphophorus hybrid embryos.a, Developmental stage and ndufs5 genotypes of hybrid embryos with X. malinche (Xm) mitochondria.b, Developmental lag of embryos with X. malinche mitochondria with varying ndufs5 genotypes compared with their most developed broodmate.c, Frequency of homozygous X. birchmanni ndufa13 ancestry over F 2 hybrid development.Dots and lines represent observed frequency ± 2 × s.e.m. (n = 208 embryos, 74 juveniles, 932 adults).d, F 2 siblings showing different phenotypes as a function of ndufs5 genotype: a heterozygote (top) and an ndufs5-incompatible sibling with matched scale (bottom left) and magnified (bottom right).e, Frequency of homozygous X. birchmanni ancestry along chromosome 13 in embryos with X. malinche mitochondrial DNA (mtDNA) that lagged their siblings by at least 1 developmental stage (red) versus those without developmental lag (blue) (Supplementary Information 1.1.7).Dashed line indicates ndufs5 location.Corresponding analyses of chromosomes 6 and 15 are shown in Supplementary Figs.7-9.f-h, Respiratory and morphometric Fig.3| Physiology and proteomics of viable heterozygotes.a, Schematic of ancestry at loci of interest in X. birchmanni × malinche F 1 hybrids (the incompatibility is not lethal in F 1 ; Supplementary Information 1.2.2). b, Real-time quantitative PCR analysis of differences in mitochondrial copy number in liver tissue by genotype.The ratio is derived from the difference in C t between the single-copy nuclear gene Nup43 and the mitochondrial gene nd1 (n = 3 X.birchmanni, 6 F 1 , 6 X. malinche).c-e, Results of Oroboros O2K respirometer assay for adult X. birchmanni, X. malinche and hybrid individuals (n = 7 fish per genotype), with X. malinche mitochondria and heterozygous ancestry at ndufs5 and ndufa13.c, Complex I efficiency in F 1 hybrids and parental species.d, Time to reach the maximum rate of complex I-driven respiration after ADP addition in F 1 hybrids and parental species (see Supplementary Fig.25for example time-to-peak curves).For corresponding analysis of time to reach peak complex I-and complex II-driven respiration after addition of succinate, see Supplementary Fig.26.e, Maximum respiration rates during full O2K protocol in F 1 hybrids and parental species.f, Allelic balance of ndufs5 in the F 1 hybrid transcriptome and proteome.Left, allele-specific expression of ndufs5 in adult F 1 hybrids (n = 3 fish).Right, results of quantitative mass spectrometry analysis of ndufs5 peptides in mitochondrial proteomes derived from adult F 1 hybrids (n = 5 fish).Dots show the proportion of area under the spectral curve contributed by the X. malinche allele in each individual.The data on the left are for endogenous peptides present in F 1 hybrids, and the data on the right are results for heavy isotope-labelled peptide controls.The peptides quantified for each species are shown on the graph.a-f, Coloured dots and bars show the mean ± 2 × s.e.m., and grey dots show individual data.

Fig. 4 |
Fig.4| Predicted structures of Xiphophorus respiratory complex I and evolutionary rates of incompatible alleles.a, Left, Xiphophorus respiratory complex I structures generated by RaptorX using alignment to a template mouse cryo-EM structure.Coloured protein structures include Ndufs5 and Ndufa13 and the four mitochondrial-encoded nd gene products in contact with Ndufs5 or Ndufa13.Right, expanded view showing the surface of predicted contacts between these proteins.Solid black lines highlight two areas predicted to be in close contact between interspecific substitutions (alpha carbon distance ≤ 10 Å for all models), and dashed lines show three additional areas in which there was weaker evidence for a predicted contact based on computational analyses (side chain distance ≤ 12 Å in at least one model).Asterisks denote residues with substitutions in X. birchmanni computationally predicted to affect protein function (Extended Data Table3).b, Detailed view of the interface

Fig. 5 |
Fig.5| Phylogenetic analysis and ancestry mapping suggest that genes underlying the mitonuclear incompatibility have introgressed from X. malinche into X.cortezi.a, Nuclear phylogeny of Xiphophorus species, showing that X. birchmanni and X. malinche are sister species29 .b, Phylogeny constructed from whole mitochondrial DNA sequences showing that X. cortezi mitochondria are nested within X. malinche mitochondrial diversity.c, Results of simulations modelling expected mitochondrial divergence between X. malinche and X. cortezi in a scenario with no gene flow.The first set of simulations used the average mitochondrial substitution rate observed between Xiphophorus species (red), and the second used the minimum mitochondrial substitution rate observed (teal).The dotted line shows observed divergence between mitochondrial haplotypes in X. malinche and X. cortezi, indicating that past mitochondrial introgression is more likely than incomplete lineage sorting.d, Clustal alignment of Ndufs5 sequences shows that X. malinche and

Fig. 1 |
ABC inference of additional selection parameters.(a) Posterior distributions of dominance coefficients from approximate Bayesian computation (ABC) simulations fitting observed data in F 2 hybrids.

Fig. 2 |
Ancestry depletion in natural hybrid populations.Average non-mitochondrial parent ancestry on (a) chromosome 13, (b) chromosome 6, and (c) chromosome 15 in natural hybrid populations fixed for X. birchmanni (Aguazarca, Acuapa, left column) or X. malinche (Tlatemaco, right column) mitochondrial haplotypes.Vertical dashed lines represent the position of ndufs5 and ndufa13 within the admixture mapping regions.
Data Fig. 3 | Chromosome 15 incompatibility.(a) Observed genotype frequencies at the peak associated marker (3.37 Mb) on chromosome 15 in the admixture mapping population.(b) Schematic of identified interactions with the X. birchmanni mitochondrial genome from our mapping data and strengthof selection underlying each interaction in hybrids (gray skull -moderate, black skull -near lethal).We discuss interactions with the X. birchmanni mitochondria in more detail in Supplementary Information 1.2.1-1.2.2.Extended Data Fig. 5 | Juvenile F 2 heart morphology by ndufa13 genotype.In all panels, colored points and bars show the mean ± 2 SE, and gray points show individual data.Individuals that were homozygous X. birchmanni at ndufa13 (the incompatible genotype) and homozygous X. malinche at ndufa13 (the compatible genotype) were raised in common laboratory conditions and sampling occurred at approximately 5 months of age (n = 5 juveniles per genotype).Measurements were taken from the sagittal section of largest cross-sectional area for the atrium (a-c) and ventricle (d-f), and area of occupancy for each cell class was calculated from the average of three quadrats (Supplementary Information 1.3.3).(g) Representative images of atria from incompatible (X.birchmanni) and compatible (X.malinche) individuals.Images are from the slide with maximum cross-sectional atrial area from each individual.(h) Example of histology analysis process in a representative atrium.Red square indicates randomly placed quadrat in which occupancy was calculated, and yellow borders represent areas which were manually annotated as cardiomyocytes.Note that the epicardium is visible at top left.All fish were raised and processed as one experimental group, with no independent attempts to test reproducibility at the experimental level.1 nature portfolio | reporting summaryApril 2023Corresponding author(s):Benjamin Moran (benmoran@stanford.edu)Molly Schumer (schumer@stanford.edu)Last updated by author(s): Nov 2, 2023 Reporting Summary Nature Portfolio wishes to improve the reproducibility of the work that we publish.This form provides structure for consistency and transparency in reporting.For further information on Nature Portfolio policies, see our Editorial Policies and the Editorial Policy Checklist.