No evidence for accumulation of deleterious mutations and fitness degradation in clonal fish hybrids: Abandoning sex without regrets

Abstract Despite its inherent costs, sexual reproduction is ubiquitous in nature, and the mechanisms to protect it from a competitive displacement by asexuality remain unclear. Popular mutation‐based explanations, like the Muller's ratchet and the Kondrashov's hatchet, assume that purifying selection may not halt the accumulation of deleterious mutations in the nonrecombining genomes, ultimately leading to their degeneration. However, empirical evidence is scarce and it remains particularly unclear whether mutational degradation proceeds fast enough to ensure the decay of clonal organisms and to prevent them from outcompeting their sexual counterparts. To test this hypothesis, we jointly analysed the exome sequences and the fitness‐related phenotypic traits of the sexually reproducing fish species and their clonal hybrids, whose evolutionary ages ranged from F1 generations to 300 ky. As expected, mutations tended to accumulate in the clonal genomes in a time‐dependent manner. However, contrary to the predictions, we found no trend towards increased nonsynonymity of mutations acquired by clones, nor higher radicality of their amino acid substitutions. Moreover, there was no evidence for fitness degeneration in the old clones compared with that in the younger ones. In summary, although an efficacy of purifying selection may still be reduced in the asexual genomes, our data indicate that its efficiency is not drastically decreased. Even the oldest investigated clone was found to be too young to suffer fitness consequences from a mutation accumulation. This suggests that mechanisms other than mutation accumulation may be needed to explain the competitive advantage of sex in the short term.


| INTRODUC TI ON
The evolutionary success of sex has been one of the greatest puzzles evolutionary biologists have pondered for over a century (Weismann, 1889). Although, theoretically, sexually reproducing organisms should be outcompeted by the asexual organisms due to inherent costs of sex, they remain predominant in the world of eukaryotes. Numerous hypotheses have been proposed to explain this paradox (see Shcherbakov, 2010, for one recent hypothesis, reviewed in, e.g., Hartfield & Keightley, 2012). Some are linked to inherent ecological properties of reproductive modes, such as the penalization of asexuals for their inability to rapidly generate novel variants, thus increasing their intraspecific competition (Doncaster, Pound, & Cox, 2000;McDonald, Rice, & Desai, 2016;Peck, Yearsley, & Waxman, 1998) or vulnerability to pathogens (Hamilton, 1980). Another popular class of theories emphasizes an increased extinction risk of the clones (ultimately decreasing their diversification), or the gradual deterioration of their genomes (hereafter referred to as mechanisms of "clonal decay"). These are based on the hypothesis that every asexual lineage will acquire harmful mutations, and selection will not be efficient enough to maintain the fittest genome (Kondrashov, 1982) because individual mutations may not escape their genomic background due to a whole-genome linkage, leading to accumulation of mutations in a stochastic (Muller's ratchet) or a deterministic (Kondrashov's hatchet) manner (Kondrashov, 1993). Consistent with this hypothesis, most asexual lineages appear to be evolutionarily short-lived, and occupy terminal positions on the tree of life (reviewed in Janko, Drozd, Flegr, & Pannell, 2008;Schwander & Crespi, 2009). Several genomic studies indeed have found evidence of higher rates of accumulation of nonsynonymous mutations in asexuals (e.g., Henry, Schwander, & Crespi, 2012;Hollister et al., 2015;Howe & Denver, 2008;Johnson & Howard, 2007;Neiman, Hehman, Miller, Logsdon, & Taylor, 2010;Paland & Lynch, 2006).
Despite broad consensus that accumulation of mutation should ultimately lead to the demise of asexual lineages, it remains unclear whether the process proceeds fast enough to exterminate the asexuals before they become competitive threats to the sexual populations. Some argue that higher rates of mutation accumulation can considerably harm the clonal populations even over a short period of time (e.g., Neiman et al., 2010), especially in synergy with other mechanisms (Pound, Cox, & Doncaster, 2004;West, Lively, & Read, 1999). However, a non-negligible number of studies have failed to confirm higher loads of mutations in asexual organisms (e.g., Allen, Light, Perotti, Braig, & Reed, 2009;Brandt et al., 2017;Naito & Pawlowska, 2016;Pellino et al., 2013;Warren et al., 2018), and several studies have also demonstrated comparable or higher net diversification rates in asexual clades than in the sexual ones, suggesting these lineages may not inherently suffer from a decreased diversification, an increased extinction or both (Fontaneto, Tang, Obertegger, Leasi, & Barraclough, 2012;Johnson, FitzJohn, Smith, Rausher, & Otto, 2011;Liu et al., 2012).
Such contradictions do not necessarily invalidate detrimental effects of mutation accumulation as a factor, but they indicate that the full explanation is more complex. For example, asexuals may have some mechanisms that temporarily delay the process of clonal decay, such as "minimal sex," ameiotic recombination, gene conversions, increasing ploidy with consequent genome refreshment and masking mutations with heterozygous states (e.g., Flot et al., 2013;Halligan & Keightley, 2003;Loewe & Lamatsch, 2008;Sémon & Wolfe, 2007). It is also possible that reduction in efficacy of selection stemming from asexuality is not drastic enough to universally cause the expected trends in mutation accumulation patterns (Allen et al., 2009). It thus remains a major challenge in evolutionary biology to identify the short-term-acting mechanisms that protect the sexuals from competitive displacement by the asexuals.
The so-called sexual-asexual complexes, that is taxa in which related sexual and asexual forms coexist, play prominent roles in studies of the persistence of sex since both types of organisms enter the same evolutionary arena for mutation, selection and drift (Birky & Barraclough, 2009;Janko, Drozd, & Eisner, 2011). In some complexes, asexual lineages may emerge from the sexual populations very dynamically , which poses an additional question regarding the role of mutation accumulation process in saving sex: even if individual clonal lineages become debilitated in a relatively short term, whole asexual populations, composed of multiple frequently emerging clones, may escape extinction and pose a serious short-term or long-term threat to the sexual counterparts (Paland, Colbourne, & Lynch, 2005). Interestingly, Janko et al. (2008) showed that continuous influx of new clones into an asexual population causes inevitable replacement of the old clones by drift without invoking an age-dependent decay of the clonal lineages. This implies that ages attained by the individual clonal lineages may not correspond to their theoretical maxima determined by an age-dependent decay, but rather by the rate of influx of the new clones into an asexual community of a finite size. This drift-like clonal turnover process proved sufficient to explain some phylogenetic data (Janko, 2014), suggesting that in some cases clones may vanish from a population before mutation accumulation or other such processes compromise their fitness.
Even less is known about the direct impact of mutation accumulation on clone fitness. Indeed, if increased mutation accumulation in asexuals is relevant to the maintenance of sex, it should have phenotypic consequences when relatively old asexual lineages "… suffer in comparison to younger asexual lineages (and sexuals) with regard to important traits such as mitochondrial function and the rate of offspring production and population growth" (Neiman et al., 2010). Unfortunately, studies investigating fitness-related traits are scarce, and their results are contradictory; age-dependent fitness deterioration has been demonstrated in some lineages of viruses, bacteria and yeasts (Andersson & Hughes, 1996;Chao, 1990;Goddard, Godfray, & Burt, 2005), but not in water frogs with genomes transmitted clonally for over 25 ky (Guex, Hotz, & Semlitsch, 2002).
In this study, we focused on the Cobitis taenia hybrid complex of the European bottom-dwelling loaches, and examined whether the genomes of sexual, young asexual and old asexual lineages differed in accumulation of nonsynonymous mutations, and whether the old and the young clones differed in traits that impacted their fitness.
Our model taxon consists of several parapatrically distributed sexual species that inhabit Europe and coexist over much of the continent with their hybrids that regularly establish successful, persistent clonal lineages with a gynogenetic reproductive mode (i.e., females lay unreduced eggs that require activation by sperm from sexual males to trigger development). Production of clonal eggs is achieved via "premeiotic endomitosis" when oogonial chromosomes are duplicated, and subsequent meiotic divisions do not yield any variability since crossovers occur among identical copies of chromosomes (Dedukh et al., 2019;Juchno, Arai, Boroń, & Kujawa, 2016). Phylogenetic and crossing experiments (Choleva et al., 2012;Janko, Culling, Rab, & Kotlik, 2005;Janko et al., 2012;Janko, Vasil'ev, et al., 2005) showed that range shifts related to the Pleistocene climatic oscillations have repeatedly placed the Cobitis parental species into a secondary reproductive contact, provoking the dynamic formation of hybrid clones. Their diversity is created in a two-step process (Figure 1 insert). First, diploid clonal F1 females form through crosses of the parental species C. elongatoides with either C. taenia or C. tanaitica (the hybrid males are sterile; Choleva et al., 2012;Janko et al., 2018); for simplicity, the haploid genomes of those species have been hereafter denoted as E, T and N, respectively, and hybrid forms with particular combinations of the parental genome (so-called genomotypes) have been denoted by these letters; for example, ETT indicates a triploid hybrid with one elongatoides and two taenia genomes. In the second step, mating with sexual males may occasionally result in the incorporation of sperm genomes leading to the formation of secondary triploid clones.
Tetraploids are occasionally produced in an analogous way, but generally die before maturity Juchno et al., 2014). Contemporary Cobitis populations are composed of a mixture of clones of quite different evolutionary ages. Some F I G U R E 1 Sampling sites of Cobitis study materials. E, T, N-genotypes of C. elongatoides, C. taenia and C. tanaitica, respectively; number of specimens in parentheses. Red area = range of C. taenia; yellow area = range of C. elongatoides; blue area = range of C. tanaitica. Dashed line delimitates range of the oldest asexual lineage, the "hybrid clade I" clones are very recent. Other clones originated during the early Holocene with establishment of the secondary contact zones between the parental species. Molecular dating using mtDNA data indicates that the origin of the oldest known lineage dates to ~0.3 Mya (Janko, Culling, et al., 2005;Majtánová et al., 2016).
This lineage, referred to as the "hybrid clade I," originated in the Ponto-Caspian region as a diploid C. elongatoides-tanaitica (EN) hybrid, and during its subsequent expansion, it has incorporated additional genomes on multiple occasions, resulting in a monophyletic clonal assemblage of the EN, the EEN, the ENN and the ETN clonal genomotypes. Owing to a lack of evidence on the observation of asex → sex transitions, the emergence of asexuality in the Cobitis appears to be a unidirectional process that has been more or less continuously occurring for at least hundreds of thousands of years, which should be enough to observe the putative effects of clonal decay (Loewe & Lamatsch, 2008). Mechanisms ensuring such long-term coexistence with sexual forms have not been properly understood, although Kotusz et al. (2014) and Bartoš et al. (2019) documented niche segregation between the Cobitis forms, suggesting some role of ecological factors.
Dynamic emergence of asexuality, the presence of diploid and polyploid clones, and a relatively well-resolved evolutionary history and ecology all make the Cobitis an attractive model for investigation of genomic and phenotypic consequences of asexual reproduction, allowing, among other objectives, the disentanglement of polyploidy effects from asexuality and/or hybridization. Here, we simultaneously analysed the genetic variability in nuclear and mitochondrial genes, together with fitness data, to test predictions of mutation accumulation hypotheses. Specifically, we investigated whether asexual genomes carried traces of reduced efficacy of purifying selection, which may be manifested by d N /d S ratios closer to 1 along the "asexual" branches of the phylogenetic trees (Wertheim et al., 2015), by generally higher rates of accumulation of the nonsynonymous mutations or by more radical amino acid substitutions in clonal Note: Upper-case letters represent a haploid genome of these species: E = C. elongatoides; L = C. paludica; N = C. tanaitica; T = C. taenia; for example, EET stands for triploid with two genomes of C. elongatoides and one of C. taenia. Categories were used to divide data sets for use in linear models (LME/GLM/GLME), where sexuals were also used for interspecific differences between pairs of species. F1 denotes F1 laboratory hybrids. Clone ID denotes the membership of given individual to particular multilocus lineage inferred from SNP data, which represent independent clonal lineage, with number of individuals per lineage in parentheses.  (Pellino et al., 2013;Sharbrough, Luse, Boore, Logsdon, & Neiman, 2018). We further compared mutation frequency spectra to test whether the spread of novel nonsynonymous variants differed between the sexual and the asexual lineages (Hollister et al., 2015).
Finally, we extended our analyses to the phenotypic data from additional individuals and tested whether traits related to the body condition and fecundity tended to decay in the ageing clones.

| MATERIAL S AND ME THODS
Origins and numbers of analysed specimens are listed in Figure 1 and Tables 1 and 2. The taxonomical identities of all fish used in this study were determined with a set of previously published and validated species diagnostic markers including allozyme and microsatellite loci and flow cytometric determination of ploidy Janko et al., 2012).

| Genomic data
Exome data were acquired using the exome capture method described in Janko et al. (2018) and briefly below. Using probes for targeted enrichment of gDNA loci, we collected sequence data from nuclear exomic and mtDNA segments of 22 females belonging to three hybridizing sexual species, 26 hybrids of different ploidy, and the genome composition and one specimen of Cobitis paludica as an outgroup (Table 1; Figure 1). Samples of parental species represent all major phylogroups, sensu Janko, Culling, et al. (2005), to capture maximum coverage of intraspecific variability in our data. Hybrid specimens include two or more individuals from several independent clonal lineages characterized previously by Janko et al. (2012), to compare genetic variability within and between clonal lineages.

| Sequencing and SNP calling
Isolated gDNA was sheared with Bioruptor, tagged by indices, pooled and hybridized to custom-designed sequence-capture probes.
Captured fragments were sequenced on an Illumina NextSeq in 75 bp paired-end (PE) mode. Fastq files were trimmed by the fqtrim tool (Pertea, 2015) with minimum read lengths of 20 bp, and 3′ end trimming when average base quality within a sliding window drops below 15, which is a commonly recommended threshold (e.g., Suren et al., 2016).
To further process the reads, we used a published Cobitis reference transcriptome (GGJF00000000; Janko et al., 2018), from which Janko et al. (2018) discarded potentially paralogous contigs, as described in Gayral et al. (2013). In short, we considered as potentially paralogous contigs those that possessed spurious heterozygosity patterns when the same heterozygotic positions occurred across distantly related species (including the outgroup). Such patterns are unlikely to originate from meaningful biological processes, and probably result from mapping of paralogous reads onto one contig present in the reference. Loci with such properties were removed from the ref- We also compared our NGS data with Sanger sequences of several loci obtained by Choleva et al. (2014) to assess the accuracy of our method.

| Tests of relaxed selection in mitochondrial genomes
Because haploid mtDNA data are not biased by the problem of phasing, we first analysed mitochondrial genomes by mapping reads to published C. elongatoides mitochondrion (Accession no. NC_023947.1) and applying the relax software (Wertheim et al., 2015) to individual haplotypes to test the hypothesis that selection is relaxed along asexual branches, resulting in d N /d S ratios closer to 1 than along sexual branches. Since the software requires phased sequences representing single ORF with no heterozygous sites or stop codons, we discarded all noncoding sequences (D-loop and tRNAs), and the reverse-oriented ND6 gene from our mtDNA reference. We also masked several observed heterozygous positions and stop codons with "N" using the hyphy package (Pond, Frost, & Muse, 2005), since these may represent sequencing errors or mutated premature termination codons. In case of overlapping ORFs in several mitochondrial genes, we split the sequence at the start codon of the second ORF, and padded the 3′ end of the first ORF with "N" to keep the correct reading frame of the whole alignment.  Table S1).
We then used both the alignment and Newick tree file in hyphy 3.8 to perform relaxed selection analysis in relax with the assumption of strict sex → asex transitions, which is justified based on known pathways of asexual hybridization in Cobitis (Choleva et al., 2012).
The software requires two types of branches to be specified in each tree. As test branches, we selected those leading exclusively to asexual individuals (i.e., asexual branches). As reference branches, we considered the intraspecific tree branches within all three sexual species (i.e., sexual branches; Figure 2). We left the interspecific branches unclassified, since sequence evolution at the interspecific level is likely to bear traces of stronger negative selection, as many mildly deleterious mutations segregating within species would not become fixed between species. Choleva et al. (2014) showed that C. tanaitica is fixed for introgressed elongatoides-derived mitochondrion and appears paraphyletic to C. elongatoides, potentially indicating repeated introgression events in the past. Therefore, we also left one internal branch within C. tanaitica unclassified, as we could not determine whether it represents an intra-or interspecific branch (see Figure 2). We also left unclassified the branches leading to laboratory progeny used in F1 crosses, as they represent non-natural experimental strains with unknown reproductive modes in the F1 generation.
We compared the fit with data of the two models allowing site-specific d N /d S ratios, that is the null model assuming the same ratios for asexual (test) and sexual (reference) branches of the ML phylogenetic tree (see Figure 2), and the alternative model assuming that d N /d S ratios differ between types of branches.

| Tests of reduced efficacy of selection in nuclear DNA
While mitochondrial data allowed the application of sophisticated methods (i.e., relax), we were also interested in trends in nuclear DNA, which is highly heterozygous in hybrids and thus had to be F I G U R E 2 Maximum likelihoodbased reconstruction of phylogenetic relationships of mtDNA haplotypes as input for relax analysis. Asexual branches were tested (blue) against intraspecific sexual branches, which served as reference branches (red) in the relax model. Blue box represents lineages of the oldest asexual lineage, the "Hybrid clade I." Letter codes after each sample indicate its genomotype analysed differently. Specifically, we inspected whether ratios of synonymous and nonsynonymous mutations differ between sexual and asexual populations by comparing the frequency spectra in each population, and testing whether radicality of nonsynonymous We also tested whether radicality of acquired nonsynonymous mutations could be higher in asexuals because lack of segregation in asexuals should reduce the efficacy of selection against mutations with strong effects, especially when the clone possesses another functional copy of a given gene (Hollister et al., 2015;Sharbrough et al., 2018). We performed two tests of this hypothesis. We first calculated for all categories of SNPs the proportions of polymorphic sites carrying premature stop codons, which likely represent strong-effect mutations. Second, we scored the "radicality" of nonsynonymous mutations with BLOSUM90 and PAM100 substitution matrices, and compared score distributions among sexual and asexual data sets. We have chosen these matrices due to low numbers of mutation per gene, which may leave some amino acid substitutions unobserved.

| Estimating fitness effects in old versus recent clones
Inspection of SNPs, albeit very instructive, may only provide an indirect proxy for inferences about fitness decay in clones. The next logical step, albeit rarely undertaken, is to test for age-dependent fitness deterioration of clones. To make such a comparison, we examined additional samples from 7 sites in the Central European hybrid zone, where asexuals belonging to recent and ancient lineages coexist. All specimens used for phenotypic analyses were genetically assessed by standard markers to determine their genomotype .
Because diploid and triploid Cobitis hybrids differ in many traits (Bobyrev, Burmensky, Vasilev, Kriksunov, & Lebedeva, 2003;Juchno & Boron, 2010;Kotusz, 2008;Maciak et al., 2011), we restricted analyses to triploid females, and investigated traits related to body condition, growth and fecundity, predicting more prominent fitness decay in older clones compared with recent ones (Neiman et al., 2010). For this part of the study, we captured and analysed 48 EEN females belonging to the ancient (hybrid clade I) and 52 females belonging to recent clonal lineages (EET and ETT genomotypes; Table 2; Figure 1). To address seasonal variation, we performed sampling before and after the reproductive season (April and September, respectively). We caught the fish by electrofishing, sacrificed them by overdose of 2-di-phenoxyethanol and fixed them in 4% formaldehyde after removing a piece of fin tissue for genotyping.
Each specimen was measured (SL; 1 mm accuracy), weighed (Wt; 1 g) and dissected for internal organs: heart, gonad, liver and spleen, which were weighed. The bodies were cooked, then remea- oocytes. We note that these latter measures were performed on a subset of fish of similar size and age to avoid biases due to age differences (Table 2).
To further test for fecundity differences, we allowed four and three hybrid females belonging to ancient (EEN) and recent (ETT) clones, respectively, to mate at will with C. elongatoides males under seminatural conditions. We kept each triploid hybrid female in an identical tank, together with a single male, with spawning substratum following the design of Bohlen (1999). After each spawning event, we counted eggs from each clutch and let them develop under standardized conditions until they hatched and reached the third fry stage, when we evaluated their survival rate.

| Statistical analyses of data
We used R packages nlme (Pinheiro, Bates, DebRoy, & Sarkar, 2016) and lme4 (Bates, Maechler, Bolker, & Walker, 2015) to analyse differences among fish types in clutch sizes, survival rates, SNP data, site frequency spectra, amino acid radicality and d N /d S ratios. To address specifics of each data set, we used several linear models including linear mixed effects (LME) models, generalized linear models (GLM) and generalized linear mixed-effects (GLME) models. We

| Test of relaxed selection on mtDNA
Since we detected some stop codons and heterozygous positions in our putative mtDNA markers, we first attempted to test whether such contigs may not constitute "nuclear mitochondrial DNA" loci Finally, we noted that putative mtDNA contigs were sequenced with above-average coverage in comparison with other loci (Wilcoxon rank sum test with continuity correction, p-value < 2.2e−16). Such observations seem to contradict the expectations of NUMTs, thus supporting our assumption that these contigs represent mitochondrial DNA.
Smart model selection (SMS), used by the phyml server, selected GTR + G as the best model of sequence evolution (Table S1); the resulting tree is depicted in Figure 2.
Comparison of the two models in relax by likelihood-ratio test  Figure 3a). This difference between SNP types was highly significant (GLME with binomial family of error distribution, LRT against null model, p-value = 7.874699e−05). Moreover, the site frequency spectra within each species showed that mutations inducing synonymous changes were significantly more frequent than mutations changing the amino acid (LME, LRT, p = 1.61937e−107; Figure S1).
We further noted 8,840 of the so-called private asexual SNPs.
The hypothesis that such sites most likely accumulated after the hybrid origin of clonal lineages is corroborated by the fact that their proportion in hybrids' genomes was tightly correlated with the mtDNA distance of any clone from its nearest sexual haplotype   (Figure 4).

Mutation accumulation and radicality in clones and polyploids
Private asexual SNPs had a significantly higher probability of being nonsynonymous than SNPs fixed among parental species (GLME with binomial family of error distribution, LRT, p = .0008). Both intraspecific sexual polymorphisms and private asexual SNPs had similar proportions of nonsynonymous SNPs (binomial GLME, LRT, p = .369; Figure 3a). Nevertheless, we noticed that old clones had a significantly lower proportion of nonsynonymous SNPs than young clones (binomial GLME, LRT, p = 0.00055) or than sexual species when measured on intraspecific sexual polymorphisms (binomial GLME, LRT, p = .12). Altogether, these data indicate no tendency towards a higher nonsynonymous mutation rate in clones, and even suggests the possibility of lower rates in the old clonal lineage.
To evaluate the impact of polyploidization on mutation accumulation, we compared diploid and triploid asexuals within both major types of hybrids (i.e., ET vs. EET and ETT in elongatoides-taenia and EN vs. EEN in elongatoides-tanaitica hybrid types). We noticed that triploids possessed significantly higher proportions of private asex- We also scored the "radicality" of nonsynonymous mutations using the PAM100 substitution matrix, and compared score distributions among sexual and asexual data sets (Figure 5a). Again, fixed interspecific nonsynonymous substitutions were significantly less radical than intraspecific segregating mutations (LME, LRT, p = .0009) or private asexual SNPs (LME, LRT, p = .004). However, the radicality of mutations segregating within sexual species (intraspecific polymorphisms) did not differ from that of nonsynonymous private asexual SNPs (LME, LRT, p = .339). We observed no significant differences between old and young clones (LME, LRT, p = .084).
Analogous results were obtained with the BLOSUM90 matrix (data not shown).

No evidence for reduced efficacy of selection from d N /d S ratios and mutation frequency spectra
A more formal test of selection was based on the per-gene distributions of d N /d S ratios in all genomotypes (Figure 5b). A subset of 1,993 ORFs passed a threshold of 50 or more codons fully resolved in all investigated individuals. To address the hybrid origin of asexuals' variability, we compared the d N /d S ratio of asexuals with that estimated from elongatoides-taenia or elongatoides-tanaitica sequence pairs, as well as estimates from laboratory F1 hybrid offspring. The usefulness of this approach was corroborated by the fact that elongatoides-taenia pairwise d N /d S values did not differ from those in F1generation ET hybrids (LME, LRT, p = .43).

F I G U R E 4
Correlation between accumulation of private mutations and loss of heterozygosity (LOH) events. Green and blue lines represent second-order polynomial models fitted to diploid data: g(x) corresponds to the best fit, while f(x) corresponds to the model with the quadratic coefficient forced to positive values.
Note that the second model fits convex parabola, consistent with the assumption that accumulating LOH events erase increasing proportion of private SNPs. Letters correspond to haploid genomes as follows: Most genes in all comparisons had values below 1, consistent with the influence of negative selection. Within-species pairwise comparisons yielded significantly higher values than interspecific comparisons, and values observed within asexual genomotypes (LME, LRT, p < 1e−8). However, the d N /d S ratios calculated from interspecific pairs of sequences were significantly higher than values from old clones (LME, LRT, p =0.01), young clones (LME, LRT, p = 0.02) or asexuals in general (LME, LRT, p = 0.004), which was at odds with predictions based on mutation accumulation hypotheses.
Moreover, we found no differences between old and young clones (LME, LRT, p = 0.35).
Looking at the distribution of private asexual SNPs, we found that nonsynonymous substitutions often appeared as singletons, while synonymous substitutions were significantly more likely to be shared by  two or more members of the same clone (Figure 5c; binomial GLME, LRT, p = 5.7e−6). Since the ancient clonal group (hybrid clade I) had the highest sample size, we also looked at its mutation frequency spectra and noticed that nonsynonymous substitutions are again significantly less frequent than synonymous ones (Poisson GLM, LRT, p = 1.1e−6).

| No evidence for fitness deterioration in old clones
We found a significant correlation between fish length and weight (ANCOVA, F 1, 56 = 115.13, p = 3.346e−15) and also a significant effect of season, with spring samples having higher SSI (t 95 = 8.34, p «0 .000), HSI (t 97 = 3.75, p = 0.000) and GSI (t 97 = 4.740, p « 0.000) values, and gonads with higher counts of generally smaller oocytes compared with autumn samples (Figure 6a). Some effect of location was also present as fishes from the Złotnica River site had smaller relative SL/Wt ratios (F 1, 56 = 33.56, p = 3.291e−07) and slower growth rates (F 1, 251 = 103.3101, p < 2.2e−16) than those from the remaining sites.
Taking the aforementioned effects into account, we found no differences between old and young clones in condition-related mea- A reproductive experiment corroborated the aforementioned differences, since females of the oldest clone had significantly larger clutches compared with those from the younger clones (average of 394 and 176 eggs per clutch, respectively; GLME with female ID as a random factor, Z = 4.364, p = 1.28e−05). We noted generally low survival rates (on average ~16% of eggs per clutch hatched and survived till third fry stage), consistent with previous experiments (Bohlen pers.comm.;Juchno et al., 2014), but we found no significant differences between old and recent clones in survival (zero-inflated GLME, Z = −.001, p = .999).

| D ISCUSS I ON
Various mechanisms have been proposed to explain why sexual reproduction evolved and became a dominant force in nature, but no clear answer has been obtained. We tested whether the mutational deterioration of clones, one of the most cited candidate processes, may serve as an efficient short-term mechanism explaining the stable coexistence of sexual species with clonal lineages, some reaching ages of up to 300 ky.

| No evidence for accumulation of deleterious mutations and fitness decay in asexuals
The observed negative selection prevailed among >3,000 studied nuclear and mitochondrial loci, as indicated by higher proportions of nonsynonymous and more radical mutations, by higher d N /d S values observed within the parental species rather than between the parental species, as well as by generally lower intrapopulation frequencies of nonsynonymous mutations compared with the synonymous ones. This is in line with the hypothesis that nonsynonymous or radical intraspecies polymorphisms may be often mildly deleterious, and may be consequently removed by negative selection before reaching fixation.
We further found that clonal genomes accumulated novel mutations in a time-dependent manner, so that the older clonal lin- The fitness-related data are also at odds with the prediction of clonal decay since the ancient "hybrid clade I" outperformed the young clones in some important traits, such as GSI and fecundity.
Given that proportions of asexuals in mixed sexual-asexual populations often correlate with their competitive success (Hellriegel & Reyer, 2000;Vrijenhoek & Parker, 2009), the evolutionary success of the oldest Cobitis clonal lineage is further indicated by its occupation of the largest distribution range of all Cobitis clones (Janko, Culling, et al., 2005) and its ability to achieve higher densities than co-occurring younger clones at many sites . We have to note that old and recent clones used in our study have different genomic compositions (elongatoides-tanaitica vs. elongatoides-taenia genome combinations, respectively), suggesting their dissimilarity may not necessarily reflect asexuality-specific processes but also reflect differences between parental species that were inherited in clonal hybrids. Alternatively, it may be also explained by a longer period of interclonal selection acting upon the older clones than the younger clones (Wetherington, Kotora, & Vrijenhoek, 1987;Wetherington, Weeks, Kotora, & Vrijenhoek, 1989). Nevertheless, the absence of a detectable effect on fitness continues to indicate that even the oldest Cobitis clones successfully coexist and compete for resources with recent clones and sexual species.

| How do asexuals escape mutational deterioration?
Altogether, the acquired data provide no evidence for reduced efficacy of the selection in asexuals, suggesting that mutation accumulation may not notably hamper the demographic and the competitive performances of the Cobitis clones even after 300 ky of asexuality.
What are the reasons behind these observations?
First, the observed patterns were unlikely to be methodological artefacts, since our study covered most of the Cobitis distribution range and provided a relatively high number of SNPs that were validated in the experimental F1 hybrids and previously published sequences, and estimated values of the phenotypic traits did not deviate from previous reports of the Cobitis genomotypes, including the polyploids (Halačka et al., 2000;Juchno & Boron, 2010). On the other hand, we note that the inherent biological nature of the studied organisms may have affected our data since comparisons between asexuals and their sexual counterparts (or simulated hybrids) conflate the variability frozen by hybridization in the distant past with that segregating within the contemporary sexual populations. Additionally, comparisons between the old and the recent clones potentially suffer owing to the presence of different species in their ancestries. Ideal asexual systems to study such processes would include both young and old clonal lineages, as well as sexual individuals with the same genomic compositions. However, even in the case of nonhybrids, natural asexuals often differ from sexuals in other ways than just the mode of reproduction. For instance, asexuality is often connected with polyploidy (Lundmark, 2006), which may affect many biological processes to a greater extent than hybrid origin itself (e.g., Maciak et al., 2011).
Carefully addressing these pitfalls, we clearly documented that mutations did accumulate in time, as predicted by the clonal decay model, but ~300 ky of asexuality did not leave detectable traces of the deleterious mutation accumulation or decreased clonal fitness.
This study adds to the growing list of publications that have failed to detect traces of clonal decay in asexual organisms, including hybrids (e.g., Pellino et al., 2013;Warren et al., 2018).
What does the growing amount of negative evidence tell us about the validity of mutation accumulation theories for evolution of sex and clonality? Commonly, negative results have been attributed to special mechanisms enabling investigated asexuals to avoid mutation accumulation, such as deviations from a strict clonality including minimal sex, ameiotic recombination, gene conversions, beneficial effects of polyploidization or an increased efficiency of DNA repair (e.g., Maciver, 2016;Roach & Heitman, 2014;Schön & Martens, 2003;Warren et al., 2018). However, given the available data, it is unlikely that any such mechanism may have efficiently slowed down the clicks of the ratchet in the Cobitis.
Refuting cryptic sexuality in predominantly clonal organisms is notoriously difficult (e.g., Birky, 2010), but available evidence argues against its role in the Cobitis. Extensive crossing experiments have never documented allelic segregation (Choleva et al., 2012;Janko et al., 2018), and Majtánová et al. (2016) reported remarkable stability of the hybrids' karyotypes with no large-scale chromosomal restructuring or recombination. The current study also shows that the vast majority of the private hybrid Polyploidization is another mechanism that may refresh a clonal genome by temporarily masking deleterious mutations, although it also frees gene copies to mutate, ultimately increasing the per-genome mutation rate (Otto & Whitton, 2000). Consistent with this expectation, we found higher numbers of private asexual SNPs in triploid loaches than in diploids, but the differences were relatively small. This may suggest that the contemporary triploids may have spent considerable parts of their history as diploid clones before polyploidization, leaving relatively little time for accumulation of extra mutations. Importantly, we found no evidence for the increased nonsynonymous mutation rates in triploids. Hence, although polyploidization likely plays a very important role in clonal evolution, at least in the Cobitis, its major benefits do not seem to stem from a deleterious mutation masking; after all, the diploid clones have persisted at the Balkan together with their triploid derivatives at least since the last interglacial period (Janko, Culling, et al., 2005). Major benefits of triploidy may instead relate to other effects, such as metabolic changes invoked by the modified cell architecture, the altered dosage of parental genomes or the modified crosstalk between the allopolyploids' genomes due to different stoichiometric relations between transcription factors and their binding sites (Bartoš et al., 2019;Beukeboom & Vrijenhoek, 1998;Maciak et al., 2011).

| The null model of clone persistence assumes delayed effects of accumulated mutations
In summary, although we may not rule out the existence of some mechanisms counteracting the effects of clonality, the aforementioned candidate mechanisms do not appear strong enough to stop or considerably slow the ratchet in the Cobitis. Instead, we propose that the apparent absence of signs of clonal decay, as in the Cobitis and other similar cases, may have a much simpler reason, which we propose as a null hypothesis for any tests of clonal decay; namely, the simplest plausible explanation is that the reduction in efficacy of a purifying selection caused by asexuality is not strong enough to allow sufficient mutation accumulation that will compromise the performance of clones over several hundreds of thousands of generations.
Consistent with this explanation, we found similar frequency spectra of non-/synonymous mutations in the sexual and the asexual populations, indicating that selection efficiently restricts the spread of putatively deleterious mutations even in the clones. The counterintuitive observations of lower nonsynonymous mutation loads and the incidence of premature stop codons in old clones may thus be explained by the classical population genetic theory, since a selection-based removal of the deleterious mutations requires varying amounts of time depending on the selection coefficient, the population size and the genetic background. It follows that recent clones, despite having acquired lower absolute numbers of mutations, will possess a higher fraction of nonsynonymous private asexual mutations, or more radical SNPs.
Different mutation loads in the old and the young clones may thus reflect a time-lag necessary to remove the slightly deleterious mutations (Johnson & Howard, 2007) rather than some special mechanisms that alleviate the mutation accumulation. Paradoxically, the delay in fitness decay may be caused directly by clonal reproduction because the maintenance of heterozygosity may mask the negative effects of acquired deleterious recessives (Halligan & Keightley, 2003) until new mutations accumulate in high numbers (Otto, 2007) or become exposed to selection in the homozygous states (Guex et al., 2002;Leslie & Vrijenhoek, 1978, Leslie & Vrijenhoek, 1980. This was in line with recent reanalysis of published sequences of asexual animals , which showed that significant deviations from neutrality, indicative of selection against ageing clones, could only be detected in asexual complexes when clones achieved substantial ages beyond 1 million generations. Hence, even if a clonal lineage were to accumulate deleterious mutations since its birth, it may take considerable time before they start to negatively affect the evolutionary lifespan of that clone. In younger asexual complexes, the replacement of clonal lineages may be dominated by a more neutral drift-like process of the clonal turnover, which does not depend on the effects of accumulated mutations (Janko et al., 2008). By analogy to a genetic drift, the clonal turnover predicts that a major portion of clonal diversity may be formed by the young clonal lineages, which did not have enough time to spread, while ancient lineages will be less diverse but more geographically widespread. Interestingly, both predictions have been met in the Cobitis and in some other taxa (Janko, 2014;Janko et al., 2012;Quattro, Avise, & Vrijenhoek, 1992).

| CON CLUS IONS
To date, the existence of successful clones in nature has usually been explained by their special abilities to avoid accumulation of the deleterious mutations. Our study shows that even if clones clearly accumulate mutations as they age, the reduction in efficacy of selection due to asexuality may not be sufficient to cause their fitness to deteriorate, even after hundreds of thousands of generations of asexuality. Thus, although mutation accumulation may ultimately drive any clone to extinction, individual clones may be replaced by more recent clones in a drift-like process of the clonal turnover if the influx rates of new clones are high (Janko et al., 2008). Clones in such complexes may simply not be given an opportunity to become old enough for accumulated mutations to start to manifest their effects.
Successful and widespread clones may therefore often represent merely temporary winners in the ongoing turnover race.
This may seem trivial, but it has serious implications for our understanding of the persistence of sex: Mutation accumulation may be relevant to the survival of clones after they reach some substantial age, but before they do so, the persistence of sexual species in mixed sexual-asexual complexes must rely on different mechanisms that offer short-term advantages to sex.

ACK N OWLED G EM ENTS
We would like to thank Jan Šipoš and Petr Pyszko for valuable advice on linear mixed-effect models, Jacek Stefaniak for preparing maps of samples, and Jan Eisner for help with the fitting of second-order polynomial models to our data. This research was