Maintenance of deleterious mutations and recombination suppression near mating-type loci under selfing

I reviewed the paper entitled “The fate of recessive or overdominant mutations near mating-type loci under partial selfing”, by Tezenas et al. In this interesting manuscript, the authors theoretically investigate how mutations are purged or can accumulate near a mating-type locus, with a two-locus two-alleles model. The dominance of the introduced mutation, and the mating system of the populations are the main variables affecting the purging of the mutation. The main findings are that the presence of a mating-type locus always decreased the purging probability and increased the purging time of the mutations under selfing. The authors also showed that deleterious mutations can accumulate near the mating-type locus over evolutionary time scales. 

I found the introduction to be very interesting, presenting the previous models used to study similar questions, what are their limitations and what will be improved in the manuscript. It also presents the mathematical tools previously used and that will be used in the manuscript.

I liked the fact that the authors put a lot of details around the mathematical tools they used in M&M and Results sections, even if it can make the model hard to follow for non-theoretical people. 

The discussion is interesting, and proposes to link theoretical outcomes with empirical data, notably in fungi. The limitation section is also nice and discusses some potential unexplored mechanisms that could modify the outcomes of the presented model.

Overall, I enjoyed reading the paper, and I only have minor comments related to the clarity of some sections of the manuscript. I hope you will find them helpful. 

Comments:

L141: “Each individual reproduces at rate 1:” -> Consider reformulating the sentence, the beginning is confusing to me.

Section 2.1: I would explain somewhere what is the difference between intra- and inter-tetrad selfing events somewhere, as it is probably not obvious for most readers, notably regarding the production of heterozygotes offspring (it is explained in lines 386-387, but it could be described sooner). 

Lines 147-149: In the case of inter-tetrad selfing, when referring to a “pairs of gametes” among all the possible viable possibilities, you only sample one gamete per tetrad right? As written it is not clear.

L240: Remove the “in English”.

L280-290: By simulations, do you mean numerical iterations, or individual-based simulations? If I understood well this is numerical iterations, but I am not 100% sure. Add a little bit more detail here about the simulations performed.

L322: “the deleterious mutation is purged from the population before it reaches a significant frequency with probability one” -> “the deleterious mutation is purged from the population with probability one before it reaches a significant frequency” ?

L321-329: Give numerical values for the purging time. You are not referring to any figures or tables, so numerical values would help the reader to understand what “… purged … take a much longer time…” means.

L330-348: I think this section can be shortened and moved toward the Discussion section. It is long, and even if interesting, it can be removed from here as the Results section is already long enough to avoid comparisons with previous models. 

L367-368: Again, a numerical value would help to understand what “long” means.

Fig.4: The mean is not very informative here, and is driven by rare, extreme times of purging. I would remove the “mean” lines to lighten the figure. In addition, you are not discussing this metric in the main text.

L516: “Deleterious population” -> You mean deleterious mutation? If not, the sentence is unclear.

L523: Abu Awad & Roze is not the only model investigating multi-locus dynamics of interacting deleterious mutations with different mating systems, indeed Lande and Lande & Porcher paved the way (Lande, 1977; Lande & Porcher, 2015), more recently some people added dominance to these models also (Clo & Opedal, 2021). Nevertheless, those models are a bit different because in these multi-locus models, mutation accumulation is rather positive for fitness in selfing populations due to compensatory effects. 

Literature cited:

Clo, J. & Opedal, Ø.H. 2021. Genetics of quantitative traits with dominance under stabilizing and directional selection in partially selfing species. Evolution 75: 1920–1935. John Wiley & Sons, Ltd.

Lande, R. 1977. The influence of the mating system on the maintenance of genetic variability in polygenic characters. Genetics 86: 485–498. Oxford University Press.

Lande, R. & Porcher, E. 2015. Maintenance of quantitative genetic variance under partial self-fertilization, with implications for evolution of selfing. Genetics 200: 891–906. Oxford University Press.

The paper analyses how linkage to a mating type locus (which is necessarily heterozygous) influences the persistence/extinction time of a newly arisen mutation at a selected locus, focusing on scenarios where this mutation is either deleterious and (partially) recessive, or overdominant. The paper allows for arbitrary linkage between the two loci, partial selfing which may be within or between tetrads produced by a individual. The analysis focuses on the initial stochastic dynamics of the mutation and models this as a multi-type branching process (where types correspond to different two-locus genotypes carrying the newly arisen muttaion). The analysis focuses on the dominant eigenvalue of the rate matrix C, which governs the long-term evolution of the expected numbers of each type, as a proxy for extinction times.

I think the paper is carefully written and well-organised. Nevertheless I have some substative comments/queries:

- In general, I wonder if it might be better to show 1/rho or 1/|rho| instead of rho in the figures, as 1/|rho|  has a more intuitive biological interpretation (as being rough proxies for extinction times, when rho is negative). This would perhaps also give a better visualisation of the sheltering effect: for example, rho=0.01 and 0.001 correspond to extinction times (and sheltering effects) that are different by an order of magnitude but appear as essentially the same color (e.g. in figure 2). By contrast,  rho=0.1 and 0.01 (which again correspond to a difference in factor of 10) appear as very different colors. Thus, in my opinion, heat maps of rho do not give a sense of the magnitude of the sheltering effect, especially in parameter regimes where it may be most substantial.  

- Lines 311-316: With a recombination rate r0 per base pair and d base pairs, the assumption that the total recombination rate r=r0*d breaks down beyond r~0.1-0.2. By the same logic, it might be better to have recombination rate rather than base pairs along the x-axis of fig 5.

- If the goal is to understand whether the sheltering effect is stronger in the vicinity of a mating type locus, would it not be more meaningful to look at the derivative of rho with respect to r as r→ 0, rather than r→ 0.5. Evaluating the derivative at r=0.5 essentially corresponds to a comparing the effect of two locations that are both quite far from the mating type locus. For example, a chromosome of 100 cM (respectively 150 cM) corresponds to a recombination fraction r~0.432 (resp. r~0.475) between the two ends of the chromosome (assuming Haldane’s map function). Even if one allows for crossover interference, r=0.5 would still correspond to rather large map scales, that are not so relevant to questions considered here (about the formation of strata and sheltering effects in the vicinity of the mating type locus)?

- Line 341: Based on the notation above, w23 should perhaps be w32; same comment for line 369.

- Overdominance scenario (line 349 onwards): It may be interesting to discuss the effect of the mating type locus on unlinked (r=0.5) loci from the point of view of it being a baseline “genomewide” effect. I was also curious whether one can come up with a clean biological interpretation for rho to be negative in this r=1/2 limit (for 0<f<1): based on App 7.3.2. Based on Fig 2, rho appears to switch sign at some threshold value of s4: can we understand this threshold in a more intuitive way for the r=1/2 case?

- Lines 359-362: This does not really correspond to a special case of over-dominance, but just describes a completely recessive b allele.

- Figure 2: In the f=1 panel (bottom row), there appears to be a discontinuous change in the blue shades at a critical value of r, independent of s4: (e.g., around r=0.1 for pin=0, around r=0.2 for pin=1 etc.). Is this an artefact of the color scale or a real effect- and if it is the latter, then what does it correspond to?

- Section 3.2: I think in general it may be good to distinguish more carefully between how the presence of a mating type locus affects sheltering vs. how (tight) linkage to a mating type locus affects sheltering. The first question is related to how the presence of a mating type locus changes dynamics at an unlinked (r=0.5) locus in comparison to a situation where there is no mating-type locus. For instance, how different is the largest eigenvalue (as calculated here) or extinction time of a deleterious allele from that in a single-locus model with only the load locus (and no mating-type locus).

- Figure 3: It wasn’t clear to me what is meant by the values of the derivative- “scaled by their minimum”- what does the minimum refer to here? I was also wondering whether it’d be more natural to show the derivative (at r=0.5) scaled by the value of rho at r=0.5, since this gives an indication of the relative change in rho. However, more importantly: see my comment above on evaluating the derivative of rho w.r.t. r at r=0.5 vs. r=0.

- Figure 4: Perhaps it may be more insightful to show this as a semilog plot.. I understand that the distribution is not predicted by the Gumbel law, but what sort of distribution is this? Is it a power law (at least in some regime of t), for example, as one would expect for a critical BP? Also what is the value of rho in this figure: i.e., how close is it to zero?

This preprint outlines branching-process based results to demonstrate how recessive and overdominant deleterious alleles are sheltered under mating-type loci under different types and levels of self-fertilisation. Results are presented using both analytical formulae and simulations, which uncover novel behaviour regarding how deleterious mutations could be sheltered for long periods of time. For example, the model is used to investigate the trade-off between the purging of deleterious mutations under selfing, compared to the maintenance of heterozygosity due to linkage to the mating-type locus.

This is a great preprint that clearly outlines the model derivations, along with the key results while using numerical computations to illustrate important outcomes of how recombination and selfing interact to maintain deleterious variation around mating-type loci. It is also clearly written, with the theoretical assumptions neatly laid out and derivations well-explained – I found it a genuine pleasure to read.

My only major concern is that I found some of the biological justification of the model to be lacking. It looks at deleterious mutation dynamics around mating-type alleles under different levels of selfing, including complete selfing. However, I have trouble understanding how there can be a combination of complete self-fertilisation and mating type loci in a population, since the presence of mating-type loci should restrict selfing by definition? Perhaps the authors could provide some further biological justification as to how the two cases could co-exist, or whether these results hold under other types of balancing-selection type cases under selfing.

Minor comments (L = line number):

• L2: change “having extended with time” to “having extended over time”.
• L199: Explain what N^3 is for non-mathematical biology readers.
• L206: change ‘descendance’ to ‘descendant’
• Figure 4 : it’s hard to see the different distributions of time for each recombination rate as they overlap. Perhaps they could be made transparent in some way so it’s clearer to see how they compare?
• Appendix 6, Tables 2-4: In these table entries, I propose removing the ‘X’ to denote multiplication. For in-text formulae it is not too intrusive but given the amount of algebra on display in these tables, then including such multiplication symbols makes each cell harder to parse. On a positive note, I found the inclusion of these tables a great addition to the manuscript, and they helped aid me in understanding the model derivations.
• Page 34, set of equations just under L828: on the 5th line of equations, I think there’s a typo: the last entry is meant to end with s^u and not s^v.