Interactions between mechanisms of reproductive isolation

Interactions between mechanisms of reproductive isolation | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results Interactions between mechanisms of reproductive isolation View ORCID Profile Alexandre Blanckaert , View ORCID Profile Vitor C. Sousa doi: https://doi.org/10.1101/2025.02.20.639255 Alexandre Blanckaert a cE3c - Centre for Ecology, Evolution and Environmental Changes & CHANGE - Global Change and Sustainability Institute, Faculdade de Ciências, Universidade de Lis-boa , Campo Grande, Portugal Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Alexandre Blanckaert For correspondence: blanckaert.a{at}gmail.com Vitor C. Sousa a cE3c - Centre for Ecology, Evolution and Environmental Changes & CHANGE - Global Change and Sustainability Institute, Faculdade de Ciências, Universidade de Lis-boa , Campo Grande, Portugal Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Vitor C. Sousa Abstract Full Text Info/History Metrics Supplementary material Data/Code Preview PDF Abstract Speciation is responsible for the diversity of species observed today and corresponds to the build-up of reproductive isolation between populations. Reproductive isolation can be generated by different mechanisms that have been extensively characterized, yet how their interactions affect speciation remains largely unknown. Here, we explicitly model the interaction of three key mechanisms (local adaptation, mate choice and genetic hybrid incompatibilities) quantifying their relative contribution to the evolution of reproductive isolation. We modeled two populations exchanging migrants using Fisher’s Geometric Model for local adaptation, phenotype matching for mate choice, and multiple pairs of Bateson-Dobzhansky-Muller Incompatibilities (DMI). All three mechanisms were determined by the same set of loci, creating conditions for interactions between barriers both at the genetic and population levels. We found very few cases where the three barriers evolved together. Instead, two barriers could evolve depending on the migration rate: either local adaptation and genetic incompatibilities for limited migration, or local adaptation and mate choice for higher migration. Our results showed that local adaptation due to ecological differentiation was the first to evolve and by far the most effective reproductive barrier. Finally, we demonstrated that in a polygenic model, populations could become locally adapted and evolve strict mate choice, yet they would not accumulate incompatibilities provided that there was sufficient gene flow. Speciation has been a central theme in evolutionary biology since its early origin ( Darwin, 1859 ). Yet, despite a truly large body of work ( Kirkpatrick and Ravigné, 2002 ), general rules of how speciation works and which mechanisms are involved remain to a large extent a mystery. So far, most findings provide a case by case set of rules rather than universal ones (see Table 1 in Sobel et al., 2010 ). Two key questions have been under active research and debate. First, what should be measured to capture progress in the speciation process ( Westram et al., 2022 ; Mallet and Mullen, 2022 )? Second, what is the end of the speciation process? This was historically defined as complete reproductive isolation (RI) between diverging populations ( Coyne and Orr, 2004 ). However, recent advances in genomics reveal that gene flow between well-defined species is far from uncommon ( Payseur and Rieseberg, 2016 ), making it difficult to draw the line of when speciation completes ( Hey and Pinho, 2012 ; Dopman et al., 2024 ). Furthermore, it has been recently argued that partial RI may not be a transient state on the way to speciation but that it may also be “a final equilibrium” state ( Servedio and Hermisson, 2020 ; Barraclough, 2024 ) further blurring the notion of species. View this table: View inline View popup Download powerpoint Table 1: List of parameters and their default values. Status indicates whether the parameter was always constant across all simulations (Fixed), could change between simulations (Variable) or within simulations (Evolvable). Different mechanisms of RI can drive the speciation process. They have been classified depending on the timing (pre- or postzygotic) and nature of selection acting upon the mechanism of RI (extrinsic or intrinsic) ( Coyne and Orr, 2004 ). Disentangling which RI mechanism(es) have the largest impact during speciation, and whether their interactions are synergistic or antagonistic remained a central question of speciation research. In other words, do interactions lead to coupling (“coincidences of barrier effects, resulting in a stronger overall barrier to gene flow”, sensu Butlin and Smadja, 2018 , synergistic interactions) or to a reduction of the overall strength of the RI barrier (antagonistic)? There have been many attempts to characterize the joint evolution of different mechanisms of RI (e.g., mate choice, hybrid load), or the evolution of one conditioned on the other, with assortative mating models via “magic trait” or reinforcement being known examples. Indeed, the coupling of RI barriers is considered an important step in the speciation process ( Barton, 1983 ; Butlin and Smadja, 2018 ; Dopman et al., 2024 ). Yet, it may lead to overall weaker rather than stronger RI ( Aubier et al., 2024 ), and it remains unclear under which conditions multiple barriers lead to synergistic or antagonistic outcomes. Theoretical studies have shown that multiple RI barriers do not always act synergistically (e.g., Maisonneuve et al., 2024 ). There has been a recent push to fully characterize the different RI barriers acting between pairs of nascent species. For example, male mate choice and (almost complete) hybrid sterility have been reported between Tetranychus cinnabarinnus and T. urticae spider mites ( Cruz et al., 2024 ). Similarly, natural populations of Xiphophorus birchmanni and X. malinche swordtail fish adapted to different thermal environments have strong assortative mating (female choice). Their hybrids display “widespread misregulation” of thermotolerance genes, suggesting an overlap between genes under local adaptation and the ones responsible for genetic incompatibilities ( Payne et al., 2024 ). However, quantifying the role of these mechanisms can be challenging, even for species amenable to laboratory manipulation (e.g., Aguillon et al., 2025 ). In addition, there are no general rules emerging. For example, Castillo et al. (2015) conducted experimental evolution of Caenorhabditis remanei on different substrates, measuring the evolution of local adaptation, mate choice and genetic incompatibilities, concluding that local adaptation was not relevant for the evolution of RI. This contrasts with other studies highlighting the importance of local adaptation in natural populations (reviewed in Table 1 of Rundle and Nosil, 2005 ), with many empirical examples of ecological speciation: Darwin finches ( Grant and Grant, 2024 ), ciclids ( McGee et al., 2020 ), freshwater stickleback ( Hatfield and Schluter, 1999 ; Marques et al., 2019 ), whitefish ( Rogers and Bernatchez, 2007 ) Heliconius butterfly ( Mallet and Barton, 1989 ), and pea aphids ( Via et al., 2000 ). While the spatial context of speciation (i.e., allopatric, parapatric and sympatric) has historically been strongly debated, several authors proposed to consider whether speciation occurs with or without gene flow, irrespective of the spatial context ( Butlin et al., 2008 ; Stankowski and Ravinet, 2021 ; Bolnick et al., 2023 ). When there is gene flow between populations (non sterile F 1 s), genetic incompatibilities create a hybrid load that may act as a sieve, leading to removal of genetic incompatibilities from both populations. Thus, provided that the mutations involved in incompatibilities do not confer a local selective advantage, or are not tightly linked to a mutation that does, we expect genetic incompatibilities to be lost when there is gene flow ( Agrawal et al., 2011 ; Bank et al., 2012 ; Blanckaert and Hermisson, 2018 ). In addition, migration between two populations in contact with each other creates divergent selection on mate choice traits ( McPeek and Gavrilets, 2006 ; Yamaguchi and Iwasa, 2013 ) to avoid the production of unfit hybrids (whether the cause is intrinsic or extrinsic). Overall, this makes migration and the strength of gene flow a key parameter to any study relative to speciation. Here, we characterized the joint evolution of three specific mechanisms of RI, determining whether there were synergistic or antagonistic interactions between the different RI barriers. We investigated whether the evolution of one mechanism impacted the evolution of the others, and quantified their interaction by assessing their joint effect on RI, compared to their individual contributions. The three mechanisms of RI were chosen to represent each a different type: mate choice (here an intrinsic prezygotic mechanism), accumulation of genetic incompatibilities (intrinsic postzygotic) and local adaptation (extrinsic pre- and postzygotic). We decided to investigate this question, using a model where all three RI mechanisms interact with each other at the genetic level. By considering a common underlying genetic architecture, we expected this setup to favor the joint evolution of RI barriers, increasing the likelihood of observing interactions between barriers. In details, we assumed a polygenic basis for adaptation (e.g., adaptation to host plants of T. urticae ( Villacis-Perez et al., 2024 ); see also ( Mackay, 2001 )), with the phenotype of an individual defined by the contribution of many loci. A subset of these loci formed intrinsic genetic incompatibilities, as such shared genetic basis should be common ( Schluter and Conte, 2009 ; Kulmuni and Westram, 2017 ) and such pattern has been identified in Mimulus guttatus (copper tolerance, Macnair and Christie, 1983 ) or Arabidopsis thaliana (resistance against biotrophic pathogens, Świadek et al., 2017 ). Recently, Frayer et al. (2025) compiled a list of known genetic incompatibilities - among them the most likely explanation for 18 of DMIs was host pathogen conflicts and adaptation for 4.2%. Moreover, the fitness of an individual is defined as the compound effect of extrinsic and intrinsic fitness; the first term is determined by the phenotype (i.e., how well adapted it is to the environment) and the second by the intrinsic genetic load due to genetic incompatibilities. Finally, (intrinsic) mate choice is based on the phenotype distance between potential mates (phenotype matching). This last assumption corresponds to a “classical magic trait” (according to Servedio et al., 2011 ) and has been reported in various organisms: wing pattern in Heliconius butterfly ( Mallet and Barton, 1989 ; Jiggins et al., 2001 , 2005 ) or body size in stickleback fishes ( Snyder and Dingle, 1989 ; McKinnon et al., 2004 ); for a more complete list of biological putative examples, see Table 1 of Servedio et al. (2011) . As argued by Thibert-Plante and Gavrilets (2013) , these traits might be more common than expected, as selection favors their formation when initially absent. Overall, our approach included all five major, important features of speciation models, as identified by Kirkpatrick and Ravigné (2002) . Under this model, we showed that local adaptation was the first mechanism to evolve and also the key mechanism responsible for the strength of reproductive isolation. Typically, depending on the level of migration between populations, either mate choice or the accumulation of genetic incompatibilities evolved in tandem with local adaptation. Interestingly, we found that even populations that evolve local adaptation and strict mate choice, might not accumulate genetic incompatibilities, provided that there is sufficient gene flow. Model To understand the evolution of different types of RI barriers and their interactions, we considered a model with a common architecture. We investigated a polygenic model where loci determining the phenotype, which affected local adaptation and mate choice, could also be involved in genetic incompatibilities, as illustrated in Figure 1 . Due to the complex polygenic genetic basis of the phenotype and genetic incompatibilities, we used individual-based simulations. To investigate the long term evolution, we followed > 10 N generations of two populations of N diploid dioecious individuals exchanging migrants at rate m (two-island-island model). To model local adaptation, we considered that the two islands represent two environments (A and B) with different phenotypic optima, and for simplicity, we referred to the population residing in each environment as population A and B, respectively. We let the phenotype, genetic incompatibilities and mate choice evolve, exploring the impact of key processes in the build-up of RI. As detailed below, we varied: (1) the extent of migration, (2) the genetic architecture of phenotypes and intrinsic genetic incompatibilities, (3) life history of species (juvenile vs adult dispersal stage, as they can affect evolutionary outcomes ( Johst and Brandl, 1997 ; Débarre and Gandon, 2011 )), and (4) initial levels of divergence between populations (ancestral vs secondary contact), keeping the other parameters fixed (see Table 1 for the complete list of parameters). Download figure Open in new tab Figure 1: Schematic illustration of the genetic basis of the different mechanisms of reproductive isolation (RI). The diploid genome (panel A) consists of L loci separated in n LB linkage blocks depicted by the blue squares. Their effect on each phenotypic trait are summed to generated the phenotype, depicted in yellow. In addition, a modifier locus (in dark orange) controls the strength of choosiness ( P C ), which determines how steep the probability of accepting a mating is (panel C, see eq. (3) ). The phenotype determines the extrinsic fitness through a quadratic function (panel D, see eq. (1) ). Finally, genetic incompatibilities are depicted by the red arrows between blue squares (panel B for the full map for the “default” scenario), and their number determined the intrinsic fitness (panel E, see eq. (2) ). The overall fitness (before mate choice) is obtained by combining extrinsic and intrinsic fitness, as depicted in panel F. The life cycle of individuals included the following sequence of events at each generation: (i) selection (based on intrinsic and extrinsic fitness components), (ii) migration, (iii) mating and (iv) zygote formation (including recombination and mutation). We considered two possible life cycles differing in whether selection or migration occurred first, i.e., differing in the dispersal stage of the organism. With “adult dispersal”, selection occurred before migration, whereas with “juvenile dispersal” selection occurred after migration. Since fitness of an individual was computed when it reached adulthood, the difference between the two life cycles was caused by the dispersal stage and whether immigrating individuals had their fitness determined by the environment they were born in (“adult dispersal”) or by the environment they inhabited at the time of reproduction (“juvenile dispersal”). Migration was modeled as a Poisson process, of parameter N×m (per default, m = 0.05 and N = 10, 000). We investigated migration ranging from 5×10 −5 to 0.2, corresponding to 2 Nm values between 1 to 4,000, covering most of the reported values in Hey and Pinho (2012) (we did not considered lower values as m = 5 × 10 −5 and m = 0 generated the same patterns in our simulations). Each population produced N offspring with a 50:50 sex ratio. Genotype and phenotype The genome of each diploid individual was composed of two different (autosomal) elements: a modifier locus determining mate choice and one large genomic element determining the phenotype and genetic incompatibilities, composed of L biallelic loci (per default, L = 500) distributed equally over n LB linkage blocks (fixed to n LB = 5). These L biallelic loci determined the n t quantitative traits (per default, n t = 1). Mutation was bidirectional, symmetric and happened at rate µ (per default, µ = 10 −6 per generation per locus). All mutations were pleiotropic and isotropic (affecting all n t phenotypic traits equally), with phenotypic effects for the n t traits drawn from a multivariate centered normal distribution with the covariance matrix , where is the identity matrix and the variance of mutational effects (fixed to σ µ = 0.05). Mutational effects were additive at the phenotypic level, hence the phenotype of a (diploid) individual was given by , for a given phenotypic trait j , with X i taking the value 0 if the individual was homozygous for ancestral allele, 1 if heterozygous and 2 if homozygous for the derived allele at locus i . For simplicity, we will refer to the phenotype that is at the optimum in the focal environment as the “optimal phenotype”, and to the phenotype at the optimum in the other environment as the “alternative phenotype”. Given that we only explored a low number of phenotypic traits (1 to 3), we did not correct the mutational effect size by the size of the phenotypic space. Recombination between adjacent loci (within a block) happened at rate ρ (fixed to ρ = 10 −3 per individual and per generation), while recombination between linkage blocks happened at rate ρ LB (fixed to ρ LB = 10 −2 ). For mate choice, we modeled choosiness as a single modifier locus ( Haldane, 1941 ; Nei, 1967 ; Feldman et al., 1996 ) with multiple alleles, rather than a quantitative trait with many underlying loci. Mutation happened at a fixed rate 10 µ at the choosiness locus and the new allele differed from the current allelic value by a deviation drawn from a normal centered distribution with variance (fixed to σ p = σ µ = 0.05). The effect of alleles was additive and allelic values exceeding the boundaries ( 0.5) were set to the corresponding boundary value. The choosiness locus was situated at the beginning of the main genomic element, at a genetic distance ρ MC from the first position of the first linkage block (fixed to ρ MC = 0.1). Sex was not genetically determined, and was randomly assigned to each individual at birth according to the expected sex ratio (50:50). Fitness The fitness w of an individual before mate choice was obtained by considering both extrinsic (due to the environment, here A) and intrinsic (due to the genetic incompatibilities) w int components: , and corresponded to the relative probability of an individual to initiate a mating. Therefore, since w does not account for mate choice, it represents a component of fitness. Local adaptation (extrinsic fitness) The extrinsic component of fitness was modeled using Fisher’s Geometric Model ( Fisher, 1930 ), as it captures the distribution of fitness effect of mutations across different organisms (e.g., vesicular stomatitis virus, Escherichia coli, Drosophila melanogaster or Triticum durum , Martin and Lenormand, 2006 ) or Arabidopsis thaliana ( Stearns and Fenster, 2016 ). The (absolute) fitness of each individual depended on the (Euclidean) distance between their phenotype, given by the vector over the n t phenotypic traits, and the environmental optimum. The relation between phenotype and fitness was obtained from Gros et al. (2009) with a slightly different choice of parametrization, to explicitly account for fitness in the alternative environment (see SI section M1). The fitness of an individual in environment A was given by: where Z j is the phenotypic value of trait j , is the optimal phenotypic value for trait j in environment A, Δ AB is a constant corresponding to the squared (Euclidean) distance between the phenotype of two individuals at the optima in environment A and environment B ; per default, Δ AB = 100), Q controls the shape of the fitness landscape ( Q ∈ [0, ∞[; per default, Q = 2), and Ω B is a constant corresponding to the fitness of an individual in environment A with a phenotype perfectly adapted to environment B . In this fitness landscape ( Gros et al., 2009 ) when Q < 2, epistasis between pairs of mutations at the optimum was positive on average, otherwise it was negative (with null epistasis for Q = 2). The fitness in environment B was obtained by substituting by , and Ω B by Ω A (per default, Ω A = Ω B = 0.1). The shape of the landscape determined the fitness of the F 1 individuals with regards to the parents (see SI section M2 and equation (S3)). Under the default parameters ( Table 1 ), weak hybrid vigor was the default behavior ( was larger than 0.55, the mean fitness of parents). Genetic incompatibilities (Intrinsic fitness) We modeled intrinsic fitness, as pairs of Bateson-Dobzhansky-Muller incompatibilities (DMIs, Bateson, 1909 ; Dobzhansky, 1936 ; Muller, 1942 ). Importantly, as mentioned above, the loci involved in incompatibilities also affected the phenotype and therefore could be involved in local adaptation. We considered n DMI pairs of incompatibilities (per default, n DMI = 50), with the loci chosen among the L loci that affected the individual phenotype. Per default, each locus could form at most a single DMI, with the two interacting loci chosen independently of their position along the genome. Each DMI was codominant and their interaction multiplicative, resulting in the following formula for the intrinsic fitness component: where X k and Y j are the counts of the number of derived alleles (0, 1 or 2) at the two loci forming the DMI, and ϵ is the strength of epistasis (i.e, the reduction in fitness due to a single DMI; per default, ϵ = 0.1). As extensions, we explored some alternative genetic architectures (e.g., DMIs between loci not affecting the phenotype), detailed in the Appendix. Mate choice We modeled mate choice using a phenotype matching approach following Irwin (2020) , where the probability of accepting a mate depended on the phenotypic distance between the two individuals, and was given by a Gaussian function. Rather than defining a variance, the shape of the Gaussian function was controlled by one parameter, P C (Fig. S1), that corresponded to the probability of an individual with the “optimal phenotype” accepting a mate with the “alternative phenotype” (see SI section M3). This model could applied to either male or female choice as no aspect is sex-specific, and both choosing and chosen individuals can mate multiple times. Here, motivated by the spider mite example ( Potter et al., 1976 ), we assumed that males were the choosy individuals. While mate choice has often been modeled as female choice (e.g., Servedio and Boughman, 2017 ) and male choice is rarer and possible under stricter conditions ( Fitzpatrick and Servedio, 2018 ), it is still frequently found in nature (reviewed in Barry and Kokko, 2010 ). Mate choice proceeded as followed: a male and χ potential females (sampling with replacement) were chosen (fixed to χ = 100), according to their respective fitness. Second, the male considered mating with the first potential female, accepting said mating according to his choosiness ( P C ) and how phenotypically similar the two potential mates were, with probability where and are the phenotypes at trait j of males and females, respectively. If the mating was accepted, a single offspring was produced and a new set of one male and χ females were selected. If the first potential female was rejected, the male would be considered for mating with the second female. This process continued until a mating was accepted or all χ females were rejected. If a male failed to accept any of the χ potential females, the male was discarded (but available to be chosen to initiate a mating again based on his fitness). This process continued until N offspring were produced. The number of potential females considered by a male for a mating ( χ ) therefore corresponds to a cost of choosiness, ranging from costly ( χ = 1) to cost-free ( χ → ∞). As mentioned above, choosiness ( P C ) was determined genetically by the modifier locus, with values varying between 0 (strict assortative mating) and 1 (no choosiness) at the phenotypic level, and therefore between 0 and 0.5 at the allelic level. Starting conditions and population metrics We let RI barriers evolve for 50 N generations, considering three possible starting conditions: an “Ancestral” state, a “Secondary contact” (SC) state and a “Secondary contact with initial mate choice” (SC + initial mate choice) (as detailed in the Appendix). Regardless, populations were always in a monomorphic state at the start of the simulations (i.e., all individuals were identical in each population). We used the following genetic map as “default”. It consisted in L = 500 loci affecting a single phenotypic trait ( n t = 1). Among them, 100 loci formed n DMI = 50 DMIs pairs, with each locus involved in at most a single DMI pair. This is the map illustrated in Figure 1B . A complete list of the genetic maps investigated are described in Appendix. For each barrier, we used a specific metric to capture its evolution through time. For local adaptation, we measured the within population mean phenotype, for mate choice, the within population mean choosiness and for the accumulation of DMIs, we used two metrics: the between population intrinsic hybrid load and the number of equivalent fixed DMIs needed to generate the observed hybrid load (detailed in Appendix). Finally, we measured RI as the reduction in gene flow at a neutral marker, following Westram et al. (2022) . We used the sojourn time of a neutral allele (number of generations until loss) introduced by transferring one individual from one population to the other. This was measured in isolated populations; further details are given in Appendix. Simulations were written using C++ (2011), with the gsl (v2.7, Galassi et al. (2002) ) and boost (v1.80) libraries. Analysis was conducted using Mathematica (v12, Wolfram Research, Inc. (2019) ), RStudio (v2022.07.2+576) and R (v4.2.0, R Core Team (2022) ), with figures generated using the ggplot2 package (v3.5.1, Wickham (2016) ). Results Evolution of strict mate choice or accumulation of genetic incompatibilities depends on migration rate We identified two possible regimes in the evolution of the three mechanisms of reproductive isolation (RI; Fig. 2 ). Local adaptation was always the first mechanism to establish (reaching 90% of the mean final value took between 0.2 N and 0.3 N generations for m = 0, and 1.4 N to 11 N for m = 0.2), but whether it was followed by the build-up of intrinsic postzygotic genetic barriers (accumulation of genetic incompatibilities) or the evolution of prezygotic barriers (strict conspecific mate choice) depended on the migration rate. For convenience, we referred to migration values leading to the first regime as weak and and the second as strong. Although the mean population phenotypic values moved further away from the environment optimum (here z A = 5) as migration rate increased, migration did not prevent local adaptation, as a large proportion of individuals had phenotypes close to the optimum and similar to the ones of isolated individuals (Fig S2). We noticed that local adaptation stopped occurring in some replicates when migration became close to extrinsic hybrid load (Fig. S3). We did not find any case where the three barriers evolved simultaneously consistently across all 100 replicates ( Fig. 3A ). There were only 27 cases out of 890 where local adaptation, strict mate choice, and nonnegligible hybrid load evolved between populations (using the minimum mean hybrid load or choosiness observed in isolation as thresholds - i.e., hybrid load > 0.194 and choosiness < 0.034). Download figure Open in new tab Figure 2: Evolution (in population A) of the phenotype (panel A), mate choice (B; on a log scale) and corrected hybrid load (C) for different migration rates (given in the header) and different dispersal stages over 50 N generations. The solid lines correspond to the mean over the 100 replicates, and the dashed lines to the 10 th and 90 th percentiles. Relation between the evolution of mate choice and corrected hybrid load are displayed in Figures 3 and S7. Stronger RI barrier corresponded to convergence of mean phenotype to 5, mean choosiness to 0 and mean corrected hybrid load to 1. For reference, the migration rates given in the headers correspond to 2 Nm values of 0, 10, 1000 and 4000, respectively. The genetic map used was the “default”, populations started in the “ancestral” state and other parameters correspond to the default values given in Table 1 . Download figure Open in new tab Figure 3: Evolution (in population A) of strict mate choice and accumulation of DMIs are (mostly) exclusive and depends on the migration rate. A/ Relation between (-log) choosiness (strict mate choice corresponds to large values) and mean corrected hybrid load (higher values indicate DMI accumulation) at the end of the 50 N generations. Color indicates the migration rate (on a log-scale). A detailed figure of this relation per migration rate is given in Figure S7. B/ Relative contrast between scaled log choosiness (defined between 0 (no mate choice) and 1) and scaled hybrid load (defined between 0 (no hybrid load) and 1) as a function of migration rate ( m + 10 −6 on a log-scale) and dispersal stage. Color indicates whether RI mechanisms evolved (beyond local adaptation that always evolved). We defined evolution of the RI mechanisms using the minimum values observed in isolation; i.e., “No barrier” corresponded to a hybrid load ≤ 0.194 and choosiness > 0.034, “Strict mate choice” to a hybrid load ≤ 0.194 and choosiness ≤ 0.034, “F1 low fitness” a hybrid load > 0.194 and choosiness > 0.034 and “Strict mate choice & F1 low fitness” a hybrid load > 0.194 and choosiness ≤ 0.034). Other parameters correspond to the default values and are given in Table 1 . In the absence of migration, populations evolved toward their phenotypic optimum. In the weak migration regime, once local adaptation has evolved, it was followed by the accumulation of DMIs, but no evolution of mate choice. This regime of local adaptation, (almost) random mating and DMIs accumulation across populations extended to weak migration rates, with the threshold depending on the dispersal stage (up to m < 5×10 −5 for “adult dispersal” and m < 5 × 10 −4 for “juvenile dispersal”, Fig. 3 , S4). Mate choice behaved neutrally in this context, but was displaced to values lower than P C = 1 by mutational bias. Indeed, with weak migration, evolving towards strict assortative mating ( P C → 0) was not advantageous since the population was relatively homogeneous. However, there was no selective pressure to be less choosy P C → 1; i.e., random mating), as males were allowed to consider multiple females when reproducing (i.e., for χ > 2 the cost of being choosy was lower than the effect of drift, Fig. S5). In the moderate to strong migration regime, once populations became locally adapted, strict mate choice evolved, but surprisingly, despite evolution of prezygotic barriers, there was almost no accumulation of DMIs. This corresponded to a regime of local adaptation, strict assortative mating and few or no DMIs. The evolution of strict mate choice ( P C < 0.034) can be explained by matings between individuals, adapted to the two environments, producing individuals with both unfit phenotypes (low extrinsic fitness component), and with potential expression of DMIs (low intrinsic fitness component). Thus, there was a strong selective pressure towards the evolution of assortative mating to avoid such unfit matings. Conversely, under these conditions, gene flow was still prevalent enough to prevent the accumulation of DMIs. That is, we found that despite strong extrinsic and prezygotic barriers to gene flow, both in terms of local adaptation (the optimum phenotype of environment A was ∼ 10x fitter than an immigrant with the optimum phenotype of environment B, and 1.78x fitter than the F 1 phenotype), and mate choice (a male at optimum phenotype of A will accept mating with a immigrant female with the optimum phenotype of B in ∼ 1% cases), in most cases the between-population hybrid load remained close to 0 when migration rate was higher than 0.005. We expected that incompatible alleles expressed in F 1 s would create an intrinsic hybrid load that became stronger with increasing migration rates, due to a higher proportion of F 1 individuals. However, provided that there was sufficient gene flow (i.e., through F 2 s and further backcrosses), selection against incompatible alleles involved in DMIs should prevent their accumulation (seen for m ≥ 0.005). Alternatively, very strict mate choice could result in a low mean hybrid load if immigrants behave as an isolated population in the alternative environment. This was not the case, as the presence of gene flow was further confirmed by the multimodal distribution of phenotypes after 50 N generations (Fig. S2), with peaks corresponding to F 1 individuals (phenotypic value close to 0), and backcrosses with intermediate phenotypes (F 2 s around ± 2.5). The absence of incompatible alleles at DMI loci was further confirmed by checking the genomes of individuals (after 50 N generations; Fig. S6). The dispersal stage influenced the migration rate at which we saw a transition between the two regimes ( Fig. 3B ), with a lower migration rate threshold with “adult-dispersal”. This is expected because immigrating individuals with maladapted phenotypes were more likely to mate in the “adult dispersal” regime (their fitness depends on their birth environment), leading to stronger selective pressure to avoid them. Interestingly, we did not find a parameter space where both strict mate choice and DMI evolved simultaneously (and consistently) at significant strength ( Fig. 3B , no values of m with a majority of purple points). As migration increases, two selective pressures, absent in isolation, emerged: selection for stricter mate choice and selection against DMIs. Our results suggested that selection against incompatibilities was more effective than selection for stricter mate choice, resulting in a Goldilocks zone, where neither strict mate choice evolved nor DMI accumulated. Yet, within a given set of parameters, we observed a correlation between stricter mate choice and higher hybrid load, suggesting a trend for the two mechanisms to coevolve, in particular for moderate to strong migration rate (Fig. S7). This was likely a consequence of evolution of mate choice: once it became strict enough ( P C < 10 −3 ), gene flow was sufficiently reduced such that the two populations behaved as isolated, and DMIs would start accumulating ( Fig. 3 , S8). As a consequence, we expected reinforcement (i.e., the evolution of strict mating due to genetic incompatibilities) to be unlikely under our model. To further confirm the absence of reinforcement, we investigated secondary contact cases, with or without initial strict assortative mating. Under all considered scenarios, the postzygotic barrier due to DMIs eroded over time, either until it was mostly gone or towards some hybrid load equilibrium ( Fig. 4 ). When secondary contact occurred between populations without initial mate choice (“SC” starting with P C = 1), we found that the evolution of mean choosiness towards strict mate choice was similar to when two populations diverged from an ancestral population (“Ancestral”), and the hybrid load due to DMIs converged to similarly low values (mean hybrid load ≤ 0.10) in both scenarios ( Fig. 4 left vs mid panels). In contrast, when secondary contact occurred between populations already with DMIs and strict mate choice (“SC + initial mate choice” starting with P C = 10 −4 ), we found a distinct dynamics, with maintenance of very strict mate choice for scenarios with strong migration ( m > 0.05), and evolution of hybrid load towards high equilibrium values (mean hybrid load ∼ 0.35), irrespective of the migration rate. Given that the erosion of the postzygotic barrier due to DMIs happened relatively slower than the evolution of choosiness, this supports effective purging of DMIs. Download figure Open in new tab Figure 4: Evolution (in population A and over 100 N generations) of mate choice (panel A) and corrected hybrid load (B) for different dispersal stages and different starting conditions: the “Ancestral” (default) scenario, the secondary contact (“SC”) scenario and the secondary contact with initial mate choice scenario (“SC + initial mate choice”). Different migration rate are given indicated by different colors, with each thin line corresponding to a different replicate and the thick line to the mean over the 30 replicates. The final states and link between mate choice and hybrid load are displayed in Figure S9. Other parameters correspond to the default values and are given in Table 1 Overall, this indicates that reinforcement (as defined above) is not possible under the investigated scenarios, due to the effective removal of DMIs when there is enough gene flow ( m ≥ 0.005), consistent with previous theoretical work on continent-island models ( Agrawal et al., 2011 ; Bank et al., 2012 ). We considered alternative genetic architectures to better understand how the shared genetic architecture between local adaptation and a second RI mechanism shaped their interactions. First, in the absence of local adaptation (Fig. S10), mate choice and DMIs alone could not maintain the population differentiation, either due to a collapse of the RI barrier (at weak migration rates) or ecological exclusion of one of the populations (at strong migration rates, SI section R2). Second, having the DMI loci involved in local adaptation (as presented above) or not (“neutral map”) led to a similar evolution of each of the three RI barrier investigated ( Fig. 5 , SI section R4). Finally, when mate choice was based upon assortative mating at a neutral phenotypic trait, mate choice was unable to evolve from an ancestral state nor from a secondary scenario (Fig. S11), but could persist if already present, and if F 1 s intrinsic fitness was sufficiently low (SI section R3). Download figure Open in new tab Figure 5: Evolution (in population A and over 50 N generations) of the phenotype (panel A) and number of equivalent fixed DMIs (B; computed according to equation (A1)) for different DMIs architecture (given in the header) and different migration rates. Each thin line corresponds to a different replicate, with the thick line corresponding to the mean over the 30 replicates. Color corresponds to different migration rates: m = 0 in blue, m = 0.005 in purple, m = 0.05 in pink and m = 0.2 in orange. For the B panel, the black dashed line corresponds to the weak mutation strong selection approximation for the proportion of DMI fixed between populations, . For the “neutral” architecture, DMIs are formed between loci not affecting the phenotype, while for the “default” architecture, the DMIs were formed between loci affecting the phenotype. Other parameters correspond to the default values and are given in Table 1 No accumulation of DMIs with migration irrespective of the network of DMI interactions Across migration rates, the intrinsic hybrid load was either absent or evolved towards some non-zero equilibrium value ( Figure 2C ). In the later case, the hybrid load matched what was observed in the absence of gene flow ( m = 0), and could be obtained starting from an ancestral population or a secondary contact scenario. Such equilibrium hybrid load was equivalent to having five DMI pairs differently fixed between the two populations. This indicated that in scenarios where DMIs did accumulate, only about 10-15% of all possible DMIs were differently fixed between the two populations (in terms of equivalent number of DMIs; Fig. S12). This could be partially explained by the design of the genetic architecture and mutation process, as forward and backward mutation rates (at each biallelic locus) were considered equal. To understand why only 10 to 15% of the DMIs were differentially fixed, we considered a simple deterministic model describing how DMIs would accumulate between populations in the limit of weak mutation and strong selection (i.e., strong epistasis, ϵ → ∞). In each population, at equilibrium and assuming a low mutation rate, about half of the loci would be fixed for the derived allele (i.e., A i or B i at loci A and B of DMI pair i ) and the other half for the ancestral allele (i.e., a i or b i ). In addition, for each DMI pair, the derived alleles ( A i or B i ) could only arise in the ancestral background, whereas the ancestral alleles could arise in any background. Under these conditions, the population had three possible states: fixed for a i b i , A i b i or a i B i . The expected probability of having A i fixed in the first population at a given time-point was (see SI section M5 for details). Therefore, the probability of having a DMI between the two isolated populations was , while in the absence of back mutation, this probability rose from to . We observed a lower number of DMIs in the simulations than predicted by the weak mutation strong selection approximation. This was the case irrespective of the epistasis strength ( ϵ = 0.9 or ϵ = 0.01, Fig. S13), or the genetic map (“default”, “neutral”,”network” or “optimized” maps). Yet, we found that the scaled mutation rate (i.e., θ = 4 Nµ ) influenced the proportion of equivalently fixed DMIs (Fig. S14). These results suggested that, in large populations or when mutation rate is high, multiple mutations might segregate at the same time within each population. If both A i and B i were polymorphic simultaneously in a given population, selection against those incompatible alleles would prevent their fixation. Indeed, in isolated populations at generation 50 N , almost 40% of all pairs of DMIs had both alleles A i and B i segregating in the population, yet with the rarer among the two alleles never exceeding a frequency of 0.015 (Table S1). Local adaptation as the main barrier to gene flow in scenarios with migration To quantify the barriers to gene flow we computed the mean sojourn times of a neutral allele from a focal immigrant individual. Regardless of the migration rate at which populations evolved, the mean sojourn time of an unlinked neutral marker was usually lower or similar to the expected sojourn times under neutrality, indicating that some RI barriers effectively reduced gene flow (e.g., local adaptation) where other did not or did not evolve (e.g., mate choice at weak migration rates, and genetic incompatibilities at strong migration rates, Fig. 6 ). The combination of barriers always decreased the mean sojourn time compared to the effect of each barrier individually. This result was qualitatively the same if we considered a linked neutral marker instead (Fig. S15). This was the case regardless of the sex of the focal individual, although the choosy sex (males) had longer sojourn times when mate choice was involved (“MC”, “LC+MC”, “MC+DMI”, “all”), consistent with higher gene flow when immigrants were of the choosy sex. Indeed, we demonstrated that as long as choosiness was not too costly, the choosy sex was responsible for most of the gene flow between populations (see Appendix and Fig. A1). Download figure Open in new tab Figure 6: Impact of reproductive isolation (RI) barriers on the relative sojourn time of neutral alleles at a freely recombining marker introduced by a focal immigrant individual. The sojourn time corresponds to the number of generations until the introduced marker was lost and used as a proxy for the barrier strength, with lower values indicating stronger barriers against gene flow. Sojourn times for specific RI barrier were obtained by turning off the other RI barriers, and the x-axis indicates which barriers were kept active: local adaptation (LA), mate choice (MC), postzygotic isolation (DMI), a combination of two (e.g., LA+MC) or all three barriers together (all). Depending on the migration rate, some RI barrier may not have evolved; in these instances the RI barrier was highlighted in gray. Results shown for a male (orange) and female (purple) focal immigrant individual. Sojourn times were measured relative to the mean sojourn time under neutrality across all replicates, i.e., neutrality corresponds to a value of 1 (dashed horizontal line). Lower values indicate less introgression, i.e., stronger barriers. Panel A corresponds to the “adult dispersal” stage and panel B to the “juvenile dispersal” one. Additional migration rates are displayed in Fig. S16. The corresponding figures for the absolute times (Fig. S17), and for the relative sojourn time of a linked marker (Fig. S15) are given in Supplement. Finally, Fig. S19 displays the relationship between individual barrier strength and sojourn time. Other parameters correspond to the default values. Across all combinations of barriers and dispersal stages, the shortest sojourn times were similar irrespective of the migration rate, indicating that RI barriers evolved to similar levels with and without migration. Irrespective of the dispersal stage and migration rates, cases where all barriers were active (“all”) showed shorter sojourn times, i.e., stronger barriers. Yet, cases with local adaptation (“LA”), local adaptation and mate choice (“LA+MC”) and local adaptation and genetic incompatibilities (“LA+DMI”) often led to similar sojourn times to those found by combining all barriers (“all”), although this depended on the migration rates and dispersal stage. In agreement with purging of DMIs, when migration m ≥ 0.05, the sojourn times when only the DMI barrier was active (“DMI”) were similar to neutral. In contrast, without migration, DMIs evolved successfully and acted as a RI barrier reducing sojourn times to ≈ 0.40 of the values observed under neutrality. Interestingly, across most scenarios we found differences depending on the dispersal stage in the reduction of gene flow ( Fig. 6 ). First, for scenarios with local adaptation (“LA”, “LA+MC”, “LA+DMI”, “all”), the reduction in sojourn time was stronger for “juvenile dispersal” than “adult dispersal”, i.e., local adaptation was a much stronger barrier to gene flow with “juvenile dispersal”. Second, for “juvenile-dispersal”, when all mechanisms were active (“all”), the reduction in sojourn time was mostly independent of the migration rate and of the sex of the focal immigrant individual, since focal individuals failed to reproduce in the first generation. For “adult-dispersal”, this was no longer the case, as sojourn times depended on migration rate and sex of focal individual (strongest joint effect for “all” barriers observed for m = 0.005, Fig. S4). The evolution of stronger barriers at intermediate migration rates with “adult-dispersal” was due to a balance between two mechanisms: increasing migration triggered the evolution of mate choice, but also ensured that at each generation there was a proportion of individuals with the alternative phenotype in the resident population. As migration further increased, this proportion individuals with the alternative phenotype became sufficiently large, that the focal immigrant was likely to successfully reproduce. Therefore, under the “adult dispersal” life cycle, we predict that strong barriers to gene flow leading to speciation would only evolve for scenarios with the “right” amount of migration. A high variance (overdispersion) for the mean sojourn time under some scenarios was observed, mainly for barriers involving mate choice (“MC”, “LA+MC”, “MC+DMI”, “all”) and for intermediate migration rates ( m = 0.0005 and m = 0.05). This overdispersion reflects the stochasticity across replicates, since in some replicates, populations evolved strict mate choice, while not in others. Finally, we investigated the sojourn time of neutral alleles in scenarios of secondary contact with preexisting choosiness and fixed levels of the strength of each barrier, thus removing the stochastic variation in evolution of barriers across runs. With “juvenile dispersal”, local adaption remained the strongest barrier to gene flow, i.e., shorter sojourn times for “LA” (Fig. S18), with addition of mate choice (“LA+MC”) further reducing the sojourn times (but only in non-choosy sex, and when each barrier was relatively weak). Genetic incompatibilities (“DMI”) were the weakest of all barriers, providing shorter reductions in sojourn time than the other mechanisms, but equivalent to mate choice (“MC”) in the choosy sex. With “adult dispersal”, local adaptation (“LA”) was one of the weakest barriers. Surprisingly, although with “juvenile dispersal” genetic incompatibilities (“DMI”) were one of the weakest barriers, with “adult dispersal” “DMI” provided the strongest barrier (i.e., shorter sojourn times), even stronger than local adaptation and without sex-specific effects. Discussion We investigated how the interaction of different mechanisms of reproductive isolation (RI) impacts divergence and gene flow between populations, under a polygenic model of local adaptation and genetic incompatibilities with a shared genetic architecture. For the evolution of RI, Table 2 provides a summary of general trends observed. In scenarios with migration, local adaptation was unsurprisingly always the first barrier to evolve, with either accumulation of genetic incompatibilities (DMIs) or the evolution of strict assortative mating evolving afterwards. The evolution of the last two barriers were somewhat mutually exclusive: whether strict mate choice or accumulation of DMIs evolved was dictated by migration rate, and to a lesser extent by the life cycle. Surprisingly, even in cases where local adaptation and strict mate choice evolved, DMIs acted as a very weak barrier, consistent with purging of genetic incompatibilities provided there was sufficient gene flow. Thus, reduction in gene flow (strength of RI) was mostly a consequence of the action of local adaptation, with mate choice being also effective for cases with sufficient migration. Importantly, we modeled here mate choice as male choice but the model is fully symmetric and all genetic elements were autosomal. Therefore, all conclusions apply to both male and female mate choice. View this table: View inline View popup Download powerpoint Table 2: Summary table of our results, indicating the general trend, within the parameter ranges explored. Each row corresponds to a different combination of genetic architecture, different initial conditions, and (possibly other parameter changes). We summarized the general trend in the evolution of phenotypic divergence (in most cases a proxy for local adaptation), evolution of strict mate choice, and accumulation of DMIs. Figures supporting the results are indicated below each statement, and not indicated when results were not shown. Local adaption at the center of speciation with gene flow Under our model with polygenic architectures, local adaptation played a major role in the evolution of RI as it was always the first barrier to evolve. The key role of local adaptation, in models of speciation with gene flow, had already been pointed out previously using two-locus models ( Agrawal et al., 2011 ; Bank et al., 2012 ). In addition to being the first RI barrier to evolve, with “juvenile dispersal”, local adaptation was the main mechanism responsible for the reduction of sojourn times of neutral introgressed alleles, a proxy for the level of gene flow between populations. The importance of ecological speciation has received a lot of attention (as detailed in Nosil, 2012 - but see Anderson and Weir, 2022 ). Our results are in agreement with the role of local adaptation in ecological and speciation scenarios with gene flow, but bring an important nuance. Indeed, when considering the “adult dispersal” life cycle, local adaptation drives the evolution of phenotypic differentiation between the two populations, and mate choice or accumulation of DMIs only evolve once populations become locally adapted. However, provided that populations evolved strict mate choice or accumulated DMIs, the strength of these two barriers becomes similar to local adaptation, as we found for secondary contact models from fixed population states (Fig. S18). For both the evolution from an ancestral population or from secondary contact from fixed population states, the relative strength of local adaptation depended on the life cycle, being stronger with “juvenile dispersal”. Since the only difference between the two life cycles was whether the fitness of an immigrant individual was affected or not by the new environment, this indicates that local adaptation is mostly effective upon arrival of the immigrant individual, but becomes ineffective to reduce gene flow when the focal individual successfully reproduces. This result highlights the brittle nature of local adaptation/ecological differentiation as a genetic barrier to gene flow ( Blanckaert et al., 2020 ). Evolution of mate choice when there is sufficient migration In our model, mate choice was never able to function as an independent barrier. It evolved only in the presence of local adaptation and non negligible migration. Therefore it always evolved as a consequence of another barrier, and corresponds to a form of coupling of multiple RI barriers (sensu Butlin and Smadja, 2018 ).This delay in the evolution of mate choice versus local adaptation has been previously reported across different mate choice models ( Thibert-Plante and Gavrilets, 2013 ) when investigating prezygotic barriers exclusively. Mate choice was also effective at limiting gene flow between species, but only for the non-choosy sex and when a very strict choosiness evolved. This is in line with previous simulations studies, highlighting the relatively ineffective role of assortative mating to maintain a narrow hybrid zone between populations ( Irwin, 2020 ). As in Irwin (2020) , our results partly stem from modeling choosiness using a Gaussian function and a matching phenotype approach. In particular, this implies that mating of F 1 s with identical phenotypes is very likely, and F 1 s are also likely to be accepted as potential partners by resident individuals. Our approach to model mate choice creates a “now or never” dichotomy for reducing gene flow, which is seen in differences between the choosy and non-choosy sex. If the focal immigrating individual manages to reproduce, then most of the RI barrier is gone and the introgressed marker behave somewhat similarly to a neutral mutation. Having the choosiness function more akin to a step function instead of a Gaussian could mitigate the effects described here and generate stronger barriers to gene flow. The suitability of each approach depends on whether the trait used for mate choice is qualitative or quantitative. In the later case, a continuous smooth function is a more likely representation of the biological process involved in mate choice (especially since we are not considering the environmental variance that could contribute to the phenotype and would flatten the choosiness function). Importantly, we never observed reinforcement, (i.e., the evolution of mate choice due to the presence of genetic incompatibilities), as the two mechanisms seem to exist in different parameter spaces. This is likely because the selection against DMIs, caused by the hybrid load, is more sensitive to migration than selection for stricter mate choice, resulting in a disjoint parameter space: combination of parameters where DMIs could accumulate did not overlap with the ones where stricter mate choice could evolve. While surprising given that occurrence of the two RI barriers is common in nature, we observed rare cases where both strict mate choice evolved and DMIs accumulated. In those rare cases, DMIs accumulated because of the evolution of mate choice, rather than the opposite. In addition, mate choice was not able to evolve if the trait used for assortative matings was not under local adaptation (i.e., if it was not a magic trait), but could persist if F 1 s individuals were sufficiently unfit. Finally, we chose to model mate choice as assortative mating, with both males and females being able to engage in multiple matings. This makes the results obtained here valid for both males and females. Genetic incompatibilities do not accumulate in locally adapted populations with strict mate choice when there is sufficient gene flow Two types of incompatibilities were modeled in our approach: extrinsic and intrinsic ones. Indeed, Fisher’s Geometric Model naturally generates epistasis between mutations due to the quadratic relationship between phenotype and fitness. Within this framework, Fraïsse et al. (2016) showed that many empirical patterns of speciation (e.g., Haldane rule, heterosis, hybrid breakdown) can be mimicked, making it a useful model to study RI. However, such approach makes it extremely challenging to characterize the effect of genetic incompatibilities, as they depend both on the environment and the genetic background. Therefore, our model also included intrinsic incompatibilities using the classical Bateson-Dobzhansky-Muller model ( Bateson, 1909 ; Dobzhansky, 1936 ; Muller, 1942 ). Assuming that the same locus affecting local adaptation can generate genetic incompatibilities may seem to be a strong assumption, yet examples of local adaptation generating DMIs have been reported ( Macnair and Christie, 1983 ; Świadek et al., 2017 ; Frayer et al., 2025 ). In particular, adaptation of Saccharomyces cerevisae to low glucose medium ( Anderson et al., 2010 ; Kvitek and Sherlock, 2011 ) or nyastin ( Ono et al., 2017 ) reveals the existence of genetic incompatibilities between beneficial mutations. The same pattern was also observed in Methylobacterium extorquens ( Chou et al., 2014 ). Regardless, we obtained the same results when the DMIs affected neutral loci, or when the architecture was “optimized” to favor DMI accumulation (Fig. S20), suggesting that local adaptation did not favor their accumulation. Generally, under our model, genetic incompatibilities seem to arise more as a consequence of isolation (either due to weak migration rates or strong reproductive barriers), and not as an active mechanism to prevent gene flow, matching the predictions from a two-locus model ( Bank et al., 2012 ). Our results also matched with the outcome of experimental evolution in Saccharomyces cerevisiae, where Dettman et al. (2007) observed the evolution of DMIs as a by-product of divergent selection between isolated populations. The same pattern was also observed in experimental evolution of isolated populations of Neuspora ( Dettman et al., 2008 ). Whether such accumulation of DMIs as a by-product could happen with ongoing gene flow remain largely unexplored ( Kulmuni and Westram, 2017 ). Our study reveals a duality in the effect of these intrinsic genetic incompatibilities in the speciation process. Accumulation of DMIs is relatively constrained, only happening when gene flow was absent or limited. This effect was further amplified when forward and backward mutation rates were symmetric, as the within-population mutational load generated by the possible DMIs tend to select against the fixation of any of the incompatible alleles (Fig S14), even in the absence of migration. Despite this, multiple DMIs have been reported across many organisms ( Presgraves, 2003 ; Kao et al., 2010 ; Corbett-Detig et al., 2013 ; Powell et al., 2020 ; Frayer et al., 2025 ), including a characterization of the underlying genetic mechanisms (reviewed in Table 1 of Kitano and Okude, 2024 ). Potentially, all of these DMIs may have formed in isolation, as described in the conceptual definition of the DMI model ( Coyne and Orr, 2004 ) whether due to temporal, geographical or prezygotic barriers. However, isolation is rather unlikely, as gene flow has been reported for some of the cases mentioned above (e.g., between Xiphophorus malinche and X. birchmanni Powell et al., 2020 ). Alternatively, it is possible that DMIs accumulation happens at a much slower rate, as discussed for cichlids ( Stelkens et al., 2010 ), that was not captured with our design here. Finally, DMIs are often thought as fixed differences between species. Yet as argued by Cutter (2012) and demonstrated in Drosophila melanogaster Corbett-Detig et al. (2013) , this may not be the case, and many genetic incompatibilities may be segregating within populations - a feature that we recovered in our simulations (Table S1). Under our model, as soon as migration was non negligible, genetic incompatibilities were unlikely to accumulate between populations, despite the evolution of phenotypic differentiation and strict mate choice. Interestingly, this pattern has been observed in nature in cases of ecological divergence. For example, cichlid fish species maintain their integrity despite coexisting with gene flow in at least parapatry, and in the presence of mate choice, yet hybrid breakdown remained limited ( Stelkens et al., 2015 ) and DMIs accumulation happened at a much slower rate than speciation ( Stelkens et al., 2010 ). Such pattern also matches the example of Anopheles gambiae and A. coluzzi . As summarized in Mallet and Mullen (2022) , the two species display ecological differences and strong assortative mating, but are sympatric and apart from a few islands of differentiation, 99% of the genome is homogenized, and therefore not propitious to the accumulation of DMIs ( Turner et al., 2005 ; Lee et al., 2013 ). The lack of characterized DMIs may be partly due to the challenge associated with their detection ( Nouhaud et al., 2020 ; Blanckaert and Payseur, 2021 ), but our model here also highlights that low levels of gene flow between populations are sufficient to purge DMIs. Finally, it is possible that some of the assumption of our model may limit the accumulation of DMIs. First, we are assuming a constant population size and therefore soft selection ( Wallace, 1968 , 1975 ) - low fitness phenotypes therefore have a higher chance of leaving offspring when frequent enough compared to a hard selection scenario. While extinction of either or both populations is a non negligible risk, demographics and evolutionary mechanisms have been shown to act in opposite direction with regards to the hybrid speciation process ( Blanckaert et al., 2023 ). Therefore, DMIs may play a more important role when considering hard selection with both demographic and evolutionary mechanisms - such as allowing the populations to be affected by changes in effective population size. Second, mutations were limited to point mutations and structural variants (e.g., inversions, translocation) were not included. These large-effect mutations can be central to local adaptation and ecotype formation (e.g., 22 and 12 putative inversion have respectively been reported in Littorina saxatalis and L. fabalis ( Johannesson et al., 2024 ), at least 4 putative inversions in Coelopa frigida ( Mérot et al., 2021 ), and 3 inversions in Heliconius numata ( Jay et al., 2021 )), especially in the early stages of speciation (Kirk-patrick and Barton, 2006). Third, we focused on pairwise epistasis. However, higher order epistasis (Ayala-López and Bank, 2025) may play an important role, as it can partially alleviated the within-population mutational load, making the accumulation of DMIs more likely. Conclusions We explicitly modeled the interaction between local adaptation, genetic incompatibilities (DMIs) and mate choice, by considering that the same loci could be involved in the phenotype (under selection in the two environments) and associated with DMIs, and that mate choice depended on phenotype and that choosiness was due to a modifier locus linked to sites determining the phenotype. Surprising, despite a common genetic basis that we expected to favor the joint evolution of all three RI barriers, we found that long term evolution from an ancestral population did not lead to this joint build up, i.e., locally adapted populations with strict mate choice and accumulation of DMIs. Instead, we found that populations evolved towards two regimes: either local adaptation and accumulation of DMIs with weak mate choice, or local adaptation with strict mate choice and no DMIs. This suggests two possible paths of speciation, one where ecological and genetic isolation dominate, and another where ecological and behavioral isolation dominate. The migration rate determined which regime was possible, with stronger migration rates favoring the latter. It is thus not surprising that we recover results similar to the ones that have been derived for when two of the three RI mechanisms interact. This implies that under the genetic architecture with polygenic adaptation and multiple DMIs considered here, it is unlikely that populations evolve towards build up of strong pre- and postzygotic barriers when there is gene flow. How loci involved in local adaptation and DMIs are distributed along the genome, and whether these different regimes lead to genetic signatures that can be detected based on polymorphism population genomics data remains to be further studied. Acknowledgments We thank João Carvalho, Inês Fragata and the Evolutionary Genomics and Bioinformatics lab for discussion of the project and Claudia Bank for comments on the manuscript. Footnotes Conflict of interest: The authors declare no conflict of interest. Data availability: The source code and analysis are available on github ( https://github.com/alex-blanckaert/GeneticBarrier ). Funding: This work was funded by the strategic project UIDB/00329/2025 granted to CE3C from the Portuguese National Science Foundation — Fundação para a Ciência e a Tecnologia (FCT). This study was funded by HFSP Young Investigator Grant RGY0081/2020 to V.S, FCT Grant 2022.03475.PTDC to A.B. and V.S. and FCT Grant HYBRIDOMICS to V.S. ( https://doi.org/10.54499/PTDC/BIA-EVL/4345/2021 ). We thank the National Network for Advanced Computing (RNCA) and INCD ( https://incd.pt/ ) for access and use of their computing infrastructure, funded by FCT to A.B. and V.S. (2022.38032.CPCA.A0 and 2022.15858.CPCA.A2). This version of the manuscript has been revised to streamline the manuscript. https://github.com/alex-blanckaert/GeneticBarrier References ↵ A. Agrawal , J. Feder , and P. Nosil . Ecological divergence and the origins of intrinsic postmating isolation with gene flow . International Journal of Ecology , 2011 ( 1 ): 435357s , 2011 . OpenUrl ↵ S. Aguillon , S. K Haase Cox , Q. Langdon , T. Gunn , J. Baczenas , S. Banerjee , A. Donny , B. Moran , P. Fascinetto-Zago , C. Gutiérrez-Rodríguez , et al. Pervasive gene flow despite strong and varied reproductive barriers in swordtails . Nature Ecology & Evolution , pages 1 – 12 , 2025 . ↵ J. Anderson , J. Funt , D. Thompson , S. Prabhu , A. Socha , C. Sirjusingh , J. Dettman , L. Parreiras , D. S Guttman , A. Regev , et al. Determinants of divergent adaptation and dobzhansky-muller interaction in experimental yeast populations . Current Biology , 20 ( 15 ): 1383 – 1388 , 2010 . OpenUrl CrossRef PubMed Web of Science ↵ S. Anderson and J. Weir . The role of divergent ecological adaptation during allopatric speciation in vertebrates . Science , 378 ( 6625 ): 1214 – 1218 , 2022 . OpenUrl CrossRef PubMed ↵ T. Aubier , M. Kopp , I. Linn , O. Puebla , M. Rafajlović , and M. Servedio . Negative coupling: the coincidence of premating isolating barriers can reduce reproductive isolation . Cold Spring Harbor Perspectives in Biology , page a041435 , 2024 . J. Ayala-Lopez and C. Bank . What can we gain from modelling complex hybrid incompatibilities? Evolutionary Journal of the Linnean Society , 4 ( 1 ): kzae034 , 2025 . OpenUrl ↵ C. Bank , R. Bürger , and J. Hermisson . The limits to parapatric speciation: Dobzhansky– Muller incompatibilities in a continent–island model . Genetics , 191 ( 3 ): 845 – 863 , 2012 . OpenUrl Abstract / FREE Full Text ↵ T. Barraclough . Does selection favour the maintenance of porous species boundaries? Journal of Evolutionary Biology , page voae030 , 2024 . ↵ K. Barry and H. Kokko . Male mate choice: why sequential choice can make its evolution difficult . Animal Behaviour , 80 ( 1 ): 163 – 169 , 2010 . OpenUrl CrossRef Web of Science ↵ N. Barton . Multilocus Clines . Evolution , 37 ( 3 ): 454 , 1983 . OpenUrl CrossRef PubMed Web of Science ↵ W. Bateson . Heredity and variation in modern lights . Darwin and modern science , pages 85 – 101 , 1909 . ↵ A. Blanckaert and J. Hermisson . The limits to parapatric speciation II: Strengthening a preexisting genetic barrier to gene flow in parapatry . Genetics , 209 ( 1 ): 241 – 254 , 2018 . OpenUrl Abstract / FREE Full Text ↵ A. Blanckaert and B. Payseur . Finding hybrid incompatibilities using genome sequences from hybrid populations . Molecular biology and evolution , 38 ( 10 ): 4616 – 4627 , 2021 . OpenUrl CrossRef PubMed ↵ A. Blanckaert , C. Bank , and J. Hermisson . The limits to parapatric speciation 3: Evolution of strong reproductive isolation in presence of gene flow despite limited ecological differentiation . Philosophical Transactions of the Royal Society of London. Series B, Biological Sciences , 375 ( 1806 ), 2020 . ↵ A. Blanckaert , V. Sriram , and C. Bank . In search of the goldilocks zone for hybrid speciation ii: hard times for hybrid speciation? Evolution , 77 ( 10 ): 2162 – 2172 , 2023 . OpenUrl CrossRef PubMed ↵ D. Bolnick , A. Hund , P. Nosil , F. Peng , M. Ravinet , S. Stankowski , S. Subramanian , J. Wolf , and R. Yukilevich . A multivariate view of the speciation continuum . Evolution , 77 ( 1 ): 318 – 328 , 2023 . OpenUrl CrossRef PubMed ↵ R. Butlin and C. Smadja . Coupling, reinforcement, and speciation . The American Naturalist , 191 ( 2 ): 155 – 172 , 2018 . OpenUrl CrossRef PubMed ↵ R. Butlin , J. Galindo , and J. Grahame . Sympatric, parapatric or allopatric: the most important way to classify speciation? Philosophical Transactions of the Royal Society B: Biological Sciences , 363 ( 1506 ): 2997 – 3007 , 2008 . OpenUrl CrossRef PubMed ↵ D. Castillo , M. Burger , C. Lively , and L. Delph . Experimental evolution: assortative mating and sexual selection, independent of local adaptation, lead to reproductive isolation in the nematode caenorhabditis remanei . Evolution , 69 ( 12 ): 3141 – 3155 , 2015 . OpenUrl CrossRef PubMed ↵ H. Chou , N. Delaney , J. Draghi , and C. Marx . Mapping the fitness landscape of gene expression uncovers the cause of antagonism and sign epistasis between adaptive mutations . PLoS genetics , 10 ( 2 ): e1004149 , 2014 . OpenUrl ↵ R. Corbett-Detig , J. Zhou , A. Clark , D. Hartl , and J. Ayroles . Genetic incompatibilities are widespread within species . Nature , 504 ( 7478 ): 135 – 137 , 2013 . OpenUrl CrossRef PubMed Web of Science ↵ J. Coyne and H. Orr . Speciation . Sinauer Associates Sunderland, MA , 2004 . ↵ M. Cruz , S. Magalhães , M. Bakirdöven , and F. Zélé . Wolbachia strengthens the match between pre-mating and early post-mating isolation in spider mites . bioRxiv , pages 2024 – 05 , 2024 . ↵ A. Cutter . The polymorphic prelude to Bateson–Dobzhansky–Muller incompatibilities . Trends in ecology & evolution , 27 ( 4 ): 209 – 218 , 2012 . OpenUrl PubMed ↵ C. Darwin . The origin of species . London: Murray , 1859 . ↵ F. Débarre and S. Gandon . Evolution in heterogeneous environments: between soft and hard selection . The American Naturalist , 177 ( 3 ): E84 – E97 , 2011 . OpenUrl CrossRef PubMed Web of Science ↵ J. Dettman , C. Sirjusingh , L. Kohn , and J. Anderson . Incipient speciation by divergent adaptation and antagonistic epistasis in yeast . Nature , 447 ( 7144 ): 585 – 588 , 2007 . OpenUrl CrossRef PubMed Web of Science ↵ J. Dettman , J. Anderson , and L. Kohn . Divergent adaptation promotes reproductive isolation among experimental populations of the filamentous fungus Neurospora . BMC evolutionary biology , 8 ( 1 ): 35 , 2008 . OpenUrl PubMed ↵ T. Dobzhansky . Studies on hybrid sterility. II. Localization of sterility factors in Drosophila pseudoobscura hybrids . Genetics , 21 ( 2 ): 113 , 1936 . OpenUrl FREE Full Text ↵ E. Dopman , K. Shaw , M. Servedio , R. Butlin , and C. Smadja . Coupling of barriers to gene exchange: causes and consequences . Cold Spring Harbor Perspectives in Biology , page a041432 , 2024 . ↵ M. Feldman , S. Otto , and F. Christiansen . Population genetic perspectives on the evolution of recombination . Annual review of genetics , 30 ( 1 ): 261 – 295 , 1996 . OpenUrl CrossRef PubMed Web of Science ↵ R. A. Fisher . The genetical theory of natural selection . Oxford University Press , Oxford , 1930 . ↵ C. Fitzpatrick and M. Servedio . The evolution of male mate choice and female ornamentation: a review of mathematical models . Current zoology , 64 ( 3 ): 323 – 333 , 2018 . OpenUrl PubMed ↵ M. Frayer , N. Robles , M. Rodríguez-Barrera , J. Coughlan , and M. Schumer . The molecular evolutionary basis of species formation revisited . Trends in Genetics , 2025 . ↵ C. Fraïsse , P. A. Gunnarsson , D. Roze , N. Bierne , and J. J. Welch . The genetics of speciation: Insights from fisher’s geometric model . Evolution , 2016 . ↵ Mark Galassi , Jim Davies , James Theiler , Brian Gough , Gerard Jungman , Patrick Alken , Michael Booth , Fabrice Rossi , and Rhys Ulerich . GNU scientific library. Network Theory Limited Godalming , 2002 . ↵ P. Grant and R. Grant . From microcosm to macrocosm: adaptive radiation of darwin’s finches . Evolutionary Journal of the Linnean Society , 3 ( 1 ): kzae006 , 2024 . OpenUrl ↵ P.A. Gros , H. Le Nagard , and O. Tenaillon . The evolution of epistasis and its links with genetic robustness, complexity and drift in a phenotypic model of adaptation . Genetics , 182 ( 1 ): 277 – 293 , 2009 . OpenUrl Abstract / FREE Full Text ↵ J. Haldane . The relative importance of principal and modifying genes in determining some human diseases . Journal of Genetics , 41 : 149 – 157 , 1941 . OpenUrl CrossRef Web of Science ↵ T. Hatfield and D. Schluter . Ecological speciation in sticklebacks: environment-dependent hybrid fitness . Evolution , 53 ( 3 ): 866 – 873 , 1999 . OpenUrl CrossRef PubMed Web of Science ↵ J. Hey and C. Pinho . Population genetics and objectivity in species diagnosis . Evolution , 66 ( 5 ): 1413 – 1429 , 2012 . OpenUrl CrossRef PubMed Web of Science ↵ D. Irwin . Assortative mating in hybrid zones is remarkably ineffective in promoting speciation . The American Naturalist , 195 ( 6 ): E150 – E167 , 2020 . OpenUrl CrossRef PubMed ↵ P. Jay , M. Chouteau , A. Whibley , H. Bastide , H. Parrinello , V. Llaurens , and M. Joron . Mutation load at a mimicry supergene sheds new light on the evolution of inversion polymorphisms . Nature genetics , 53 ( 3 ): 288 – 293 , 2021 . OpenUrl CrossRef PubMed ↵ C. Jiggins , I. Emelianov , and J. Mallet . Assortative mating and speciation as pleiotropic effects of ecological adaptation: Examples in moths and butterflies . 2005 . ↵ C. D Jiggins , R. Naisbit , R. Coe , and J. Mallet . Reproductive isolation caused by colour pattern mimicry . Nature , 411 ( 6835 ): 302 – 305 , 2001 . OpenUrl CrossRef PubMed Web of Science ↵ K. Johannesson , R. Faria , A. Le Moan , M. Rafajlović , A. Westram , R. Butlin , and S. Stankowski . Diverse pathways to speciation revealed by marine snails . Trends in Genetics , 40 ( 4 ): 337 – 351 , 2024 . OpenUrl CrossRef PubMed ↵ K. Johst and R. Brandl . Evolution of dispersal: the importance of the temporal order of reproduction and dispersal . Proceedings of the Royal Society of London. Series B: Biological Sciences , 264 ( 1378 ): 23 – 30 , 1997 . OpenUrl CrossRef Web of Science ↵ K. Kao , K. Schwartz , and G. Sherlock . A genome-wide analysis reveals no nuclear Dobzhansky-Muller pairs of determinants of speciation between S. cerevisiae and S. paradoxus, but suggests more complex incompatibilities . PLoS genetics , 6 ( 7 ): e1001038 , 2010 . OpenUrl PubMed M. Kirkpatrick and N. Barton . Chromosome inversions, local adaptation and speciation . Genetics , 173 ( 1 ): 419 – 434 , 2006 . OpenUrl Abstract / FREE Full Text ↵ M. Kirkpatrick and V. Ravigné . Speciation by natural and sexual selection: models and experiments . the american naturalist , 159 ( S3 ): S22 – S35 , 2002 . OpenUrl CrossRef PubMed Web of Science ↵ J. Kitano and G. Okude . Causative genes of intrinsic hybrid incompatibility in animals and plants: what we have learned about speciation from the molecular perspective . Evolutionary Journal of the Linnean Society , 3 ( 1 ): kzae022 , 2024 . OpenUrl ↵ J. Kulmuni and A. Westram . Intrinsic incompatibilities evolving as a by-product of divergent ecological selection: Considering them in empirical studies on divergence with gene flow . Molecular Ecology , 26 ( 12 ): 3093 – 3103 , 2017 . OpenUrl CrossRef PubMed ↵ D. Kvitek and G. Sherlock . Reciprocal sign epistasis between frequently experimentally evolved adaptive mutations causes a rugged fitness landscape . PLoS genetics , 7 ( 4 ): e1002056 , 2011 . OpenUrl PubMed ↵ Y. Lee , C. Marsden , L. Norris , T. Collier , B. Main , A. Fofana , A. Cornel , and G. Lanzaro . Spatiotemporal dynamics of gene flow and hybrid fitness between the m and s forms of the malaria mosquito, anopheles gambiae . Proceedings of the National Academy of Sciences , 110 ( 49 ): 19854 – 19859 , 2013 . OpenUrl Abstract / FREE Full Text ↵ T. Mackay . The genetic architecture of quantitative traits . Annual review of genetics , 35 ( 1 ): 303 – 339 , 2001 . OpenUrl CrossRef PubMed Web of Science ↵ M. Macnair and P. Christie . Reproductive isolation as a pleiotropic effect of copper tolerance in mimulus guttatus? Heredity , 50 ( 3 ): 295 – 302 , 1983 . OpenUrl CrossRef Web of Science ↵ L. Maisonneuve , C. Smadi , and V. Llaurens . Which cues are sexy? the evolution of mate preference in sympatric species reveals the contrasted effect of adaptation and reproductive interference . Evolution Letters , 8 ( 2 ): 283 – 294 , 2024 . OpenUrl PubMed ↵ J. Mallet and N. Barton . Strong natural selection in a warning-color hybrid zone . Evolution , 43 ( 2 ): 421 – 431 , 1989 . OpenUrl CrossRef PubMed Web of Science ↵ J. Mallet and S. Mullen . Reproductive isolation is a heuristic, not a measure: a commentary on westram et al., 2022 . Journal of Evolutionary Biology , 35 ( 9 ): 1175 – 1182 , 2022 . OpenUrl CrossRef PubMed ↵ D. Marques , K. Lucek , V. Sousa , L. Excoffier , and O. Seehausen . Admixture between old lineages facilitated contemporary ecological speciation in lake constance stickleback . Nature communications , 10 ( 1 ): 4240 , 2019 . OpenUrl PubMed ↵ G. Martin and T. Lenormand . A general multivariate extension of fisher’s geometrical model and the distribution of mutation fitness effects across species . Evolution , 60 ( 5 ): 893 – 907 , 2006 . OpenUrl CrossRef PubMed Web of Science ↵ M. McGee , S. Borstein , J. Meier , D. Marques , S. Mwaiko , A. Taabu , M. Kishe , B. O’Meara , R. Bruggmann , L. Excoffier , et al. The ecological and genomic basis of explosive adaptive radiation . Nature , 586 ( 7827 ): 75 – 79 , 2020 . OpenUrl CrossRef PubMed ↵ J. McKinnon , S. Mori , B. Blackman , L. David , D. Kingsley , L. Jamieson , J. Chou , and D. Schluter . Evidence for ecology’s role in speciation . Nature , 429 ( 6989 ): 294 – 298 , 2004 . OpenUrl CrossRef PubMed Web of Science ↵ M. McPeek and S. Gavrilets . The evolution of female mating preferences: differentiation from species with promiscuous males can promote speciation . Evolution , 60 ( 10 ): 1967 – 1980 , 2006 . OpenUrl CrossRef PubMed Web of Science ↵ C. Mérot , E. Berdan , H. Cayuela , H. Djambazian , A. Ferchaud , M. Laporte , E. Normandeau , J. Ragoussis , M. Wellenreuther , and L. Bernatchez . Locally adaptive inversions modulate genetic variation at different geographic scales in a seaweed fly . Molecular biology and evolution , 38 ( 9 ): 3953 – 3971 , 2021 . OpenUrl CrossRef PubMed ↵ H. Muller . Isolating mechanisms, evolution and temperature . In Biol. Symp , volume 6 , pages 71 – 125 , 1942 . OpenUrl ↵ M. Nei . Modification of linkage intensity by natural selection . Genetics , 57 ( 3 ): 625 , 1967 . OpenUrl FREE Full Text ↵ P. Nosil . Ecological speciation . Oxford University Press , 2012 . ↵ P. Nouhaud , A. Blanckaert , C. Bank , and J. Kulmuni . Understanding admixture: haplodiploidy to the rescue . Trends in Ecology & Evolution , 35 ( 1 ): 34 – 42 , 2020 . OpenUrl PubMed ↵ J. Ono , A. Gerstein , and S. Otto . Widespread Genetic Incompatibilities between First-Step Mutations during Parallel Adaptation of Saccharomyces cerevisiae to a Common Environment . PLoS biology , 15 ( 1 ): e1002591 , 2017 . OpenUrl CrossRef PubMed ↵ C. Payne , R. Bovio , D. Powell , T. Gunn , S. Banerjee , V. Grant , G. Rosenthal , and M. Schumer . Genomic insights into variation in thermotolerance between hybridizing swordtail fishes . Molecular Ecology , 33 ( 23 ): e16489 , 2024 . OpenUrl ↵ B. Payseur and L. Rieseberg . A genomic perspective on hybridization and speciation . Molecular ecology , 25 ( 11 ): 2337 – 2360 , 2016 . OpenUrl CrossRef PubMed ↵ D. Potter , D. Wrensch , and D. Johnston . Aggression and mating success in male spider mites . Science , 193 ( 4248 ): 160 – 161 , 1976 . OpenUrl Abstract / FREE Full Text ↵ D. Powell , M. García-Olazábal , M. Keegan , P. Reilly , K. Du , A. Díaz-Loyo , S. Banerjee , D. Blakkan , D. Reich , P. Andolfatto , et al. Natural hybridization reveals incompatible alleles that cause melanoma in swordtail fish . Science , 368 ( 6492 ): 731 – 736 , 2020 . OpenUrl Abstract / FREE Full Text ↵ D. Presgraves . A fine-scale genetic analysis of hybrid incompatibilities in Drosophila . Genetics , 163 ( 3 ): 955 – 972 , 2003 . OpenUrl Abstract / FREE Full Text ↵ R Core Team . R: A Language and Environment for Statistical Computing . R Foundation for Statistical Computing , Vienna, Austria , 2022 . URL https://www.R-project.org/ . ↵ S. Rogers and L. Bernatchez . The genetic architecture of ecological speciation and the association with signatures of selection in natural lake whitefish (coregonus sp. salmonidae) species pairs . Molecular biology and evolution , 24 ( 6 ): 1423 – 1438 , 2007 . OpenUrl CrossRef PubMed Web of Science ↵ H. Rundle and P. Nosil . Ecological speciation . Ecology letters , 8 ( 3 ): 336 – 352 , 2005 . OpenUrl CrossRef Web of Science ↵ D. Schluter and G. Conte . Genetics and ecological speciation . Proceedings of the National Academy of Sciences , 106 ( Supplement 1 ): 9955 – 9962 , 2009 . OpenUrl CrossRef PubMed ↵ M. Servedio and J. Boughman . The role of sexual selection in local adaptation and speciation . Annual Review of Ecology, Evolution, and Systematics , 48 ( 1 ): 85 – 109 , 2017 . OpenUrl CrossRef ↵ M. Servedio and J. Hermisson . The evolution of partial reproductive isolation as an adaptive optimum . Evolution , 74 ( 1 ): 4 – 14 , 2020 . OpenUrl CrossRef PubMed ↵ M. Servedio , GS. Van Doorn , M. Kopp , A. Frame , and P. Nosil . Magic traits in speciation:’magic’but not rare? Trends in ecology & evolution , 26 ( 8 ): 389 – 397 , 2011 . OpenUrl PubMed ↵ R. Snyder and H. Dingle . Adaptive, genetically based differences in life history between estuary and freshwater threespine sticklebacks (gasterosteus aculeatus l .). Canadian Journal of Zoology , 67 ( 10 ): 2448 – 2454 , 1989 . OpenUrl CrossRef ↵ J. Sobel , G. Chen , L. Watt , and D. Schemske . The biology of speciation . Evolution , 64 ( 2 ): 295 – 315 , 2010 . OpenUrl CrossRef PubMed Web of Science ↵ S. Stankowski and M. Ravinet . Defining the speciation continuum . Evolution , 75 ( 6 ): 1256 – 1273 , 2021 . OpenUrl CrossRef PubMed ↵ F. Stearns and C. Fenster . Fisher’s geometric model predicts the effects of random mutations when tested in the wild . Evolution , 70 ( 2 ): 495 – 501 , 2016 . OpenUrl PubMed ↵ R. Stelkens , K. Young , and O. Seehausen . The accumulation of reproductive incompatibilities in african cichlid fish . Evolution: International Journal of Organic Evolution , 64 ( 3 ): 617 – 633 , 2010 . OpenUrl CrossRef PubMed Web of Science ↵ R. Stelkens , C. Schmid , and O. Seehausen . Hybrid breakdown in cichlid fish . PloS one , 10 ( 5 ): e0127207 , 2015 . OpenUrl CrossRef PubMed ↵ M. Świadek , S. Proost , D. Sieh , J. Yu , M. Todesco , C. Jorzig , A. Rodriguez Cubillos , B. Plötner , Z. Nikoloski , E. Chae , et al. Novel allelic variants in acd6 cause hybrid necrosis in local collection of arabidopsis thaliana . New Phytologist , 213 ( 2 ): 900 – 915 , 2017 . OpenUrl CrossRef PubMed ↵ X. Thibert-Plante and S. Gavrilets . Evolution of mate choice and the so-called magic traits in ecological speciation . Ecology letters , 16 ( 8 ): 1004 – 1013 , 2013 . OpenUrl CrossRef PubMed ↵ T. Turner , M. Hahn , and S. Nuzhdin . Genomic islands of speciation in anopheles gambiae . PLoS Biol , 3 ( 9 ): e285 , 2005 . OpenUrl CrossRef PubMed ↵ S. Via , A. Bouck , and S. Skillman . Reproductive isolation between divergent races of pea aphids on two hosts. ii. selection against migrants and hybrids in the parental environments . Evolution , 54 ( 5 ): 1626 – 1637 , 2000 . OpenUrl CrossRef PubMed Web of Science ↵ E. Villacis-Perez , F. De Graeve , B. De Beer , S. Ali Alshami , R. De Jong , T. De Meyer , and T. Van Leeuwen . Independent genetic mapping experiments identify diverse molecular determinants of host adaptation in a generalist herbivore . Molecular Ecology , page e17618 , 2024 . ↵ B. Wallace . Polymorphism, population size, and genetic load . Population biology and evolution , pages 87 – 108 , 1968 . ↵ B. Wallace . Hard and soft selection revisited . Evolution , pages 465 – 473 , 1975 . ↵ A. Westram , S. Stankowski , P. Surendranadh , and N. Barton . What is reproductive isolation? Journal of evolutionary biology , 35 ( 9 ): 1143 – 1164 , 2022 . OpenUrl CrossRef PubMed ↵ Hadley Wickham . ggplot2: Elegant Graphics for Data Analysis . Springer-Verlag New York , 2016 . ISBN 978-3-319-24277-4 . URL https://ggplot2.tidyverse.org . ↵ Wolfram Research, Inc . Mathematica v. 12. Champaign, Illinois, USA ., 2019 . https://www.wolfram.com . ↵ R. Yamaguchi and Y. Iwasa . Reproductive character displacement by the evolution of female mate choice . Evolutionary Ecology Research , 15 ( 1 ): 25 – 41 , 2013 . OpenUrl View the discussion thread. Back to top Previous Next Posted October 29, 2025. Download PDF Supplementary Material Data/Code Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following Interactions between mechanisms of reproductive isolation Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share Interactions between mechanisms of reproductive isolation Alexandre Blanckaert , Vitor C. Sousa bioRxiv 2025.02.20.639255; doi: https://doi.org/10.1101/2025.02.20.639255 Share This Article: Copy Citation Tools Interactions between mechanisms of reproductive isolation Alexandre Blanckaert , Vitor C. Sousa bioRxiv 2025.02.20.639255; doi: https://doi.org/10.1101/2025.02.20.639255 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Evolutionary Biology Subject Areas All Articles Animal Behavior and Cognition (7635) Biochemistry (17691) Bioengineering (13892) Bioinformatics (41937) Biophysics (21452) Cancer Biology (18588) Cell Biology (25504) Clinical Trials (138) Developmental Biology (13378) Ecology (19899) Epidemiology (2067) Evolutionary Biology (24320) Genetics (15609) Genomics (22506) Immunology (17736) Microbiology (40394) Molecular Biology (17181) Neuroscience (88605) Paleontology (666) Pathology (2832) Pharmacology and Toxicology (4824) Physiology (7641) Plant Biology (15156) Scientific Communication and Education (2045) Synthetic Biology (4294) Systems Biology (9825) Zoology (2271)