Population genomics and epigenomics provide insights into the evolution of facultative asexuality in plants

doi:10.21203/rs.3.rs-3111723/v1

Download PDF

Article

Population genomics and epigenomics provide insights into the evolution of facultative asexuality in plants

https://doi.org/10.21203/rs.3.rs-3111723/v1

This work is licensed under a CC BY 4.0 License

Journal Publication

published 16 May, 2024

Read the published version in Communications Biology →

Version 1

posted

You are reading this latest preprint version

Many organisms, particularly plants, can reproduce both sexually and asexually but with varying frequencies. Because variations in the reproductive system profoundly affect fitness and genomic diversity, natural selection should act on changes between sexual and asexual types of propagation. However, the genomic evidence of natural selection driving switches between reproductive systems remains scarce, and the evolutionary consequences of such switches on the genomic and epigenomic diversity remain controversial. Here, by analyzing population genomes and epigenomes of Spirodela polyrhiza (L.) Schleid., a facultatively asexual angiosperm, we show that while the genomic diversity and DNA methylation levels are very low, they vary among populations. These variations are associated with both the demographic history and the frequency of sexual reproduction in each population. A genome-wide scan revealed that multiple genes involved in flowering and embryogenesis were under positive selection, consistent with the hypothesis that natural selection drove the evolution of asexuality in this plant. Together, these results provide new insights into the evolution of plant reproductive systems and suggest that natural selection can drive the evolutionary changes of asexuality, which in turn alters the levels of genomic and epigenomic diversity.

Biological sciences/Evolution/Population genetics

Biological sciences/Genetics/Epigenomics

Biological sciences/Evolution/Evolutionary genetics

Biological sciences/Genetics/Evolutionary biology

Biological sciences/Plant sciences/Plant reproduction

Understanding the evolution of sexual reproduction has long been at the centre of evolutionary biology. Theories suggest that asexual reproduction is beneficial for the short term but costly for the long term, mainly due to accumulations of deleterious mutations and low effective population size^1–5. Facultative asexuality, where organisms can reproduce both sexually and asexually depending on environmental conditions, should be optimal for one individual’s lifespan^6,7. While rather few animals such as aphids (aphidoidea)⁸, water fleas (cladocerans)⁹, and rotifers¹⁰ reproduce facultatively asexually, up to ~ 80% of the flowering plants, including important crops and keystone species, can reproduce both sexually and asexually¹¹. Asexual reproduction in plants involves different types of vegetative reproduction (e.g. runners, tubers, bulbs, corms, suckers, plantlets), as well as apomixis, the formation of seeds without fertilization¹². Because changes between sexual and asexual reproduction profoundly affect fitness, natural selection might act on the genetic controls of sexual and asexual reproduction in facultative asexual organisms, which in turn can alter the levels of genomic diversity, heterozygosity and effectiveness of selection in the population^2,5,13−15. However, direct evidence supporting this prediction remains scarce, mainly due to the lack of a suitable facultative asexually reproducing system in which the signature of selection can be detected at genomic levels.

Evolutionary changes in sexual and asexual reproduction might also affect the maintenance and dynamics of epigenetic markers, i.e., DNA methylations. In plants, cytosine methylation can occur in three sequence contexts: CpG, CHG, and CHH (H = A, T, or C), which are controlled by different mechanisms and have different dynamics during reproduction¹⁶. Typically, CpG and CHG methylation are maintained by methyltransferases1 (MET1) and CHROMOMETHYLASE3 (CMT3), whereas CHH methylation is mostly maintained by CMT2 or the RNA-directed DNA methylation (RdDM) pathway¹⁷. However, CMT3 was found to also play a role in the maintenance of CHH methylation, as well as de novo methylation¹⁸. During sexual reproduction, DNA methylations are highly dynamic¹⁹. In both male and female gametogenesis, the megaspore mother cell and its male analogue experience dramatic chromatin changes during cell specification, such as heterochromatin decondensation and an enlarged nuclear volume^20,21. During male gametogenesis, sperm DNA is highly methylated in the CpG and CHG context but has low CHH methylation in retrotransposons^19,22,23. During female gametogenesis, CpG and CHH methylation remains largely steady²⁴. After fertilization, CHH methylation increases during embryogenesis and can approach 100% at individual cytosines, which then decreases likely through a passive mechanism after germination^25–27. In contrast, during vegetative reproduction, DNA methylation is likely steady since meiosis and embryogenesis are lacking^28–30. Although a recent study comparing DNA methylations among 34 angiosperm species suggested that clonally propagated species often have low CHH methylation³¹, the extent to which asexual reproduction affects genome-wide methylation levels remains unclear.

Here, we investigated the population genome and epigenome of a facultatively asexual flowering plant, Spirodela polyrhiza (the giant duckweed; Lemnaceae), which has a global distribution. This species, like other duckweeds from the genera Spirodela, Landoltia and Lemna, is characterized by leaf-like fronds derived from fused stems³² and, with multiple roots on each frond³³ and with a highly reduced vascular system³⁴. Spirodela polyrhiza reproduces vegetatively via budding under normal conditions but rarely switches to sexual reproduction under unfavourable conditions^35,36. Recent studies showed that despite its global distribution in diverse habitats, the genomic diversity, spontaneous mutation rates and DNA methylation levels in S. polyrhiza are very low^37–42, which might be associated with its overall low frequency of sexual reproduction. DNA methylation profiling of two genotypes suggests that DNA methylation in S. polyrhiza, which is substantially lower than in other plants, likely due to the loss of the RdDM pathway, varied between genotypes⁴¹. Further insights into the evolutionary origin and consequences of asexuality on genomic and epigenomic variations in S. polyrhiza are required to understand the demographic history and to identify the footprint of selection on the genome.

Low genomic variations in S. polyrhiza

We sequenced the genomes of 131 globally distributed S. polyrhiza genotypes with an average of ~ 25 X coverage. Together with previously published samples^1,2, we analysed the genomic diversity of 228 S. polyrhiza individuals across five continents (Supplementary Data 1). We identified 1,241,981 high-quality biallelic single-nucleotide polymorphisms (SNPs) and 166,075 short insertions and deletions (INDELs, less than 50 bp of length). Based on an updated genome annotation of S. polyrhiza (see Section 1.1 of the Supplementary Information), we found that most of the SNPs (70.3%) are in the intergenic regions (Supplementary Fig. 1). Of all the SNPs located in the protein-coding regions, 61,039 were identified as nonsynonymous and 44,287 as synonymous (Supplementary Data 2). Consistent with our previous study, the genome-wide nucleotide diversity at neutral sites is 0.0013 (Supplementary Table 1), the lowest among multicellular eukaryotes. The species-wide efficacy of selection (π_N/π_S ratio) is 0.39, the highest among studied organisms⁴³, indicating relatively weak purifying selection in S. polyrhiza, despite its large effective population size^1,2.

In addition to SNPs and small INDELs, we also characterize the genome-wide structural variations (SVs, ≥ 50 bp in length) in S. polyrhiza (see Section 1.2 of the Supplementary Information). We identified 3,205 high-quality SVs, including 2,089 deletions, 291 duplications, and 825 insertions. The abundance of the SVs is substantially lower than all reported multicellular organisms^44,45. Among all identified SVs, 155 duplications and 169 deletions affected protein-coding genes (Supplementary Table 2). Interestingly, in comparison to the genome-wide average, the MADS-box genes (see Section 1.3 of the Supplementary Information) showed a higher tendency to overlap with SVs (permutation tests; P = 0.063). In addition, we also observed a significant overrepresentation of small INDELs in the MADS-box genes (permutation tests; P = 0.008). These results indicate that MADS-box genes might be under positive or relaxed purifying selection in S. polyrhiza.

Population structure and demographic history of S. polyrhiza

Because S. polyrhiza is facultatively asexual, genotypes collected from the geographic proximity can be derived from the same clonal family. Using a previously established grouping threshold that was developed in S. polyrhiza², we identified 159 likely clonal families (Supplementary Data 3).

Population structure and principal component analyses revealed four populations in S. polyrhiza (Fig. 1a and 1b). Consistent with our previous study, the four populations are largely concordant with their geographic origins, namely America, South East Asia (SE-Asia), Europe, and India (Supplementary Fig. 2), with a few exceptions that can be due to recent migration events or artefacts during long-term duckweed maintenance.

We inferred the population history with a Maximum Likelihood (ML) phylogeny combined with Approximate Bayesian Computation (ABC). Using Colocasia esculenta (from the Araceae family) as an outgroup, the maximum likelihood phylogeny of all 228 genotypes indicated an early split of the American population from the other populations and subsequent splits of the Indian and European populations from SE-Asia (Fig. 1c). The European population constitutes the most recent split (Fig. 1c-d). Here, genotypes collected from the transcontinental region (e.g. Russia) showed intermediate features of SE-Asian and European populations, suggesting this as a likely migration route. Furthermore, the Indian population possibly originated via Thailand and Vietnam, as genotypes from these countries show intermediate features between Indian and SE-Asian populations.

We modelled the demographic history using an ABC modelling approach to further validate the evolutionary history of the four populations in S. polyrhiza. Three plausible demographic scenarios were simulated, allowing for either the SE-Asian, American or an additional putative population to function as the ancestral population (Supplementary Fig. 3). We found that the scenario, in which the American population and Asian population were derived from an additional putative ancestral population, constituted the most supported model (Fig. 1d). While the American population was separated from other populations around 1 million generations ago, the European population was derived from the SE-Asian population only 12,000 generations ago (see Section 1.4 of the Supplementary Information).

Determinants of genomic diversity among populations

Among the four populations, nucleotide diversity (π) and the efficacy of selection (π_N/π_S ratio) varied among populations (Fig. 2b). While the SE-Asian population has the highest π and lowest π_N/π_S ratio, the American population has the lowest π and highest π_N/π_S ratio. Interestingly, while the European population has a much smaller π compared to the SE-Asian population, the π_N/π_S ratio of the European population remains similar to the latter, likely due to its recent split from the SE-Asian population.

Using genome-wide SNPs, we found that linkage disequilibrium (LD) is comparable to Arabidopsis thaliana⁴⁶, suggesting considerable historical sexual reproduction in S. polyrhiza. However, the extent of LD decay varied substantially among populations (Fig. 2b and Supplementary Fig. 4). While the Asian population showed the most rapid LD decay (about 12 kb at r² = 0.2), the European population had very long LD blocks (> 100 kb). The Indian and American populations had intermediate LD decay. This pattern is further supported by the recombination rates calculated for all four populations, where the Asian population had the highest recombination rate compared to the other three (Fig. 2b). Different LDs and recombination rates found among populations indicate that the frequencies of sexual reproduction varied among populations. In addition, we found that the variations of heterozygosity in S. polyrhiza showed a similar pattern with the genomic diversity and recombination rate among four populations (Fig. 2 and Supplementary Fig. 5).

Interestingly, the changes in genomic diversity and levels of heterozygosity are associated with two SVs involving MADS-box genes that are involved in sexual reproduction. One SV is an 84 bp insertion at the last coding sequence (CDS) of gene SpGA2022_005278, a homolog of AGL62 from the Mα subclade of MADS-box genes (Supplementary Fig. 6). In A. thaliana, AGL62 is a transcription factor that suppresses endosperm cellularization by activating the expression of a putative invertase inhibitor, InvINH1, in the micropylar region of the endosperm^47,48 (Supplementary Fig. 7). The insertion likely disrupted the function of the AGL62-like gene, thus reducing the suppression of endosperm development, which is likely required for sexual reproduction (Fig. 2a). Consistently, we found the insertion was at a higher abundance in the SE-Asian population (87.5%) than in other populations (Fig. 2c and 2d). In addition, the insertion positively correlates with heterozygosity within the European population (Supplementary Table 3 and Supplementary Fig. 8). Another SV is a 69 bp deletion at 1.8 kb upstream of SOC1 (SpGA2022_007306, Supplementary Fig. 5), the ortholog of which is a positive regulator of the flowering process in A. thaliana⁴⁹ (Fig. 2a). The deletion was exclusively found in the Indian population with the alternate allele frequency of 73% (Fig. 2c). It is plausible that the deletion, due to its disruption potential at the cis-regulatory region, reduces the ability of the SOC1-like gene to respond to the upstream floral activators (e.g. CO) in S. polyrhiza, thus reducing the frequency of sexual reproduction in the Indian population (Fig. 2d). Consistently, this deletion negatively correlates with heterozygosity in the Indian population (Supplementary Table 3, Supplementary Fig. 9).

Population epigenomic diversity in S. polyrhiza

As changes in sexual reproduction can also alter epigenomic dynamics, we further investigated the patterns of population epigenomic diversity in S. polyrhiza. We selected five individuals from each population and quantified their shoot DNA methylation levels at single-base resolution using whole genome bisulphite sequencing. Similar to a recent study⁴⁰, we found that only 1.5% of cytosines are methylated in S. polyrhiza (7.2% of CpG, 2.2% of CHG, and 0.1% of CHH; see Supplementary Table 4), and the average species-wide methylation levels are much lower than all other investigated angiosperm species (Supplementary Fig. 10)^31,50. The hierarchical clustering of 20 methylomes in CpG, CHG, and CHH contexts show high consistency with their genetic similarity (Supplementary Fig. 11–13). Few discrepancies were mostly found within the same population or between the recently diverged SE-Asian and European populations, which indicates that population demography can strongly affect epigenomic variations in S. polyrhiza.

We then compared the genome-wide weighted methylation level (wML) among populations. For CpG methylation, the four populations did not differ at genome-wide levels (Fig. 3a). However, CpG methylation level in gene bodies (but not in TEs) was higher in the American population than the other three populations (Fig. 3d, g). For CHG, the Indian population had the lowest genome-wide methylation level among all four populations (Fig. 3b, 3e and 3h). Interestingly, for CHH, the American population had the lowest genome-wide methylation level (P < 0.01, pairwise Wilcoxon test; Fig. 3c), while the European and Indian populations showed a gradual reduction of methylation in comparison to the SE-Asian population. The pattern was the same for both gene bodies and TEs (P < 0.01, pairwise Wilcoxon test; Fig. 3f and 3i). The genome-wide reduction of CHH methylation is consistent with the hypothesis that clonal reproduction reduces CHH methylation, and the effects gradually accumulate over clonal generations⁵¹.

The footprint of selection on the genome

To identify the genomic signature of selection at the species level we performed genome-wide scan. To reduce false positives, we used the µ-statistics from RAiSD⁵², the composite likelihood ratio CLR statistic from SweeD⁵³, and the T statistic from LASSI⁵⁴. We found 69 genes showed strong signature of selection using all three methods (Supplementary Fig. 14). Manual inspection indicated that several orthologs of these genes are related to gametogenesis (e.g., NOTCHLESS) and embryogenesis (e.g., NUP214, CPSF, CDK, AGP, and ACR4) in Arabidopsis thaliana^55–59. Further enrichment analysis indeed showed that embryo lethal genes were enriched in these 69 genes (P = 0.016, \({\chi }^{2}\) test). In addition, the A. thaliana orthologs of several genes under selection are also associated with controlling sexual reproduction, including floral development (DRMY1 and ACR4)^55,60, flowering time (NF-Y AT2G27470, NF-YAT1G72830, and CPSF), pollen development (EFOP3, ELMOD, and CLC) ^57,61−65, seed development (NUP214, NF-Y AT2G27470 and NF-YAT1G72830, and Transducin/WD40 )^56,61,66. Furthermore, among these 69 genes, we also found several genes involved leaf development and vascularity (SECA2, RbgA, PHABULOSA/PHAVOLUTA)^67–69, light signaling (NF-Y, CCR4-NOT, and PPP)^61,70,71, root development (GEND1, WAVY, and ACR4)^55,72,73 and DNA damage repair (ATM and Xrcc3)^74,75 (Supplementary Data 4).

To further understand the selection that drove the evolution within individual populations, we identified the signature of positive selection in a three-population tree using patterns of linked allele frequency differentiation and calculating the corresponding composite-likelihood ratio (CLR, see Methods). In total, we found 1,883 genes on the SE-Asian branch, 593 genes on the Indian branch and 401 genes on the European branch (Fig. 4a; see Section 1.5 of the Supplementary Information) which showed strong signatures of selection (top 1% of CLR values).

We found that genes under positive selection in the European branch are enriched with reproduction and development-related GO terms (Supplementary Fig. 15). Among these, SpGA2022_013448, in chromosome 9, is an ortholog of FLOWERING LOCUS KH DOMAIN (FLK) that delays flowering by up-regulating FLC family members in A. thaliana⁷⁶. This gene showed a strong signature of selection in the European branch but not in other branches (Fig. 4c and d). Similarly, SpGA2022_006111, on chromosome 3, is an ortholog of the A. thaliana BIG BROTHER (BB) that negatively regulates floral organ size and is also under selection in Europe⁷⁷ (Fig. 4c).

In the SE-Asian population, we found that gene SpGA2022_051517, a CHROMOMETHYLASE3 (CMT3) homolog in A. thaliana, is involved in maintaining CpG, particularly CHG methylation¹⁷. This is consistent with the higher CpG and CHG methylation levels observed in the SE-Asian population when compared to the European and Indian populations (Fig. 3a and 3b). Within the Indian population, we found that five MADS-box genes have been under selection exclusively along this branch. Given that there are 43 MADS-box genes in the genome, the fact that five of them have been targeted by selection, constitutes a significant enrichment of such genes under selection (P = 0.0075, Fisher’s Exact test). For example, SpGA2022_013078 is an ortholog of AGAMOUS-LIKE6 (AGL6), which is involved in flower and meristem identity specification in rice⁷⁸; SpGA2022_052274, a homolog of APETALA3 (AP3), is involved in the petal and stamen specification in A. thaliana⁷⁹; and SpGA2022_006905 belongs to the SHORT VEGETATIVE PHASE (SVP-group) which controls the time of flowering and meristem identity⁸⁰.

We found 77 genes under positive selection (top 1% CLR values) in both the European and Indian populations (Fig. 4b), significantly more genes than expected by chance (P < 2.2e-16, Fisher’s Exact Test). Among these, gene SpGA2022_055195, a homolog to CYP78A9 of cytochrome P450 monooxygenases in A. thaliana, belongs to a highly conserved gene family CYP78A. Previous studies in A. thaliana and other species found that CYP78A9 plays a critical role in promoting cell proliferation during flower development and further impacts seed size^81–83. Overall, these data consistently suggest that genes involved in reproduction and development were under selection in Indian and European populations, which might have led to reduced sexual reproduction in these two populations.

Here, we characterized the genomic and epigenomic diversity, as well as the demographic history of a facultative asexual flowering plant, S. polyrhiza. We found that among populations of S. polyrhiza, demographic history and reproductive system jointly determine the population’s genomic and epigenomic diversity. Analyses on the footprint of selection suggest that natural selection drove the reduced vascular system and increased asexuality in S. polyrhiza.

Theory predicts that asexual reproduction reduces genomic diversity and the efficiency of purifying selection⁸⁴. Consistent with this prediction, at the species level, we found that S. polyrhiza has very low genomic diversity and reduced purifying selection (seen as an increased π_N/π_S ratio), when compared to a wide range of spermatophyte plants⁴³. Within species, the SE-Asian population, which has the highest frequency of sexual reproduction based on the estimated recombination rate, has the highest genomic diversity, the lowest π_N/π_S ratio and the highest heterozygosity (Fig. 2b), supporting the theoretical prediction^2,5,13−15. The low π_N/π_S ratio found in the European population, which has the lowest sexual reproduction and genomic diversity, is most likely due to its migration history. The demographic model suggested that the European population derived from the SE-Asian population very recently (Fig. 1d). It is likely that the π_N/π_S ratio in the European population remained the same as its ancestral population and has not reached an equilibrium level yet.

While there are fewer genome-wide SVs in S. polyrhiza compared to other species^44,45, we found these variants and small INDELs are in tendency enriched in MADS-box genes. This indicates that the loss-of-function of MADS-box genes, many of them involved in flower development and sexual reproduction, is under natural selection. The results are consistent with the observation that the number of functional MADS-box genes was dramatically reduced in S. polyrhiza⁸⁵.

Single-base resolution methylomes of 20 individuals showed that the overall CpG, CHG and CHH methylation levels in S. polyrhiza shoots are very low, consistent with previous studies^40,86. The low levels of DNA methylation might be associated with reduced sexual reproduction: while CpG and CHG methylations in plants are important for controlling cross-overs during meiosis⁸⁷ and are increased during male gametogenesis, CHH methylation is highly accumulated during embryogenesis^19,25−27. In facultative asexual plants, due to reduced sexual reproduction and meiosis, the selection of genetic mechanisms maintaining or increasing the CpG, CHG and CHH methylation is reduced or absent, which might have led to the reduced CpG, CHG and CHH methylation levels. Consistently, a recent study suggests that S. polyrhiza has lost several genes in the RdDM pathway⁸⁶. Interestingly, within species, the methylation profile of the 20 individuals correlates with their genetic distance (Supplementary Fig. 11–13), indicating a gradual evolution of DNA methylomes in S. polyrhiza. For example, the Indian and European populations, which diverged from SE-Asian populations around 51,000 and 12,000 generations ago, gradually decreased their DNA methylations (Fig. 3a-c).

At the species level, using a genome-wide scan approach, we found a strong signature of natural selection on genes involved in flower and seed development. indicating that the evolution of reproduction, likely, an increased clonal propagation in S. polyrhiza, was driven by natural selection. This is consistent with the pattern that many aquatic organisms reproduce clonally⁸⁸. In addition, several genes related to vascularity, root development and DNA damage repair were also under strong selection, suggesting the reduced root and vascular development and low mutation rate in S. polyrhiza were likely also driven by natural selection.

Among populations, we found strong positive selection on genes involved in sexual reproduction and development in Indian and European populations, two recently evolved populations that showed reduced genomic recombination. These results suggest that natural selection favours clonal reproduction in S. polyrhiza during the recent colonization process, a pattern was frequently found in many invasive species^89,90. However, despite strong selection favouring clonal reproduction, substantial recombination in the S. polyrhiza genome, mostly in the SE-Asian population, remained, reflecting that sexual reproduction is essential to overcome the costs involved in clonal reproduction in long term.

Taken together, population genomes and epigenomes of S. polyrhiza suggest that natural selection acting on the mating system and organ development can shape genome-wide genomic and epigenomic variations.

DNA sample preparation and sequencing

We sequenced 131 genotypes that were primarily collected from Asia and Europe (Supplementary Data 1). These samples were cultivated in N-medium⁹¹ until DNA isolation using a CTAB method. Library preparations were carried out following the protocol described in Xu et al. (2019)³⁷. All libraries were sequenced either on Illumina HiSeq X Ten or Illumina Hiseq 4000 platforms for paired-end sequencing with a read size of 150 bp. Low-quality reads and adapter sequences were trimmed with AdapterRemoval (v2.033)⁹². On average, 33.8 million reads per genotype were obtained. The clean reads were aligned to the S. polyrhiza reference genome^93,94 using BWA-MEM (https://github.com/lh3/bwa) with default parameters. Reads without alignment hits or with multiple alignment positions were removed. SAMtools “rmdup” function was used to remove PCR duplicates⁹⁵.

Single nucleotide variant calling

We filtered SNPs using GATK⁹⁶ (v4.1.4.1, Java 11) with options: "QD < 2.0 | QUAL < 30.0 | SOR > 3.0 | FS > 60.0 | MQ < 40.0 | MQRankSum < -12.5 | ReadPosRankSum < -8.0". Then, VCFtools (v0.1.13)⁹⁷ was used to remove SNPs that have the following features: 1) non-biallelic; 2) minor allele frequency (MAF) less than 0.01; 3) mean sequencing depth greater than eight or less than 41; 4) missing genotypes greater than 20%; 5) located in small SNP clusters (more than two SNPs in a ten base-pair window). SnpEff (version 5.0c)⁹⁸ was used to annotate SNPs and INDELs.

We estimated genome-wide nucleotide diversity and genome-wide π_N/π_S ratios using SNPgenie (v2019.10.31)⁹⁹. The genome-wide heterozygosity for each individual was calculated using VCFtools (v0.1.13)⁹⁷. We estimated the genetic associations between heterozygosity and the SVs of AGL62 and SOC1 using RVTESTS¹⁰⁰ with single variant Wald test.

Population structure and linkage disequilibrium (LD)

We grouped genetically similar genotypes by defining clonal genotype pairs that have no more than 0.01% different homozygous sites and no more than 2% different heterozygous sites. These thresholds were previously adopted by Ho et al. (2019)³⁸.

Prior to the population structure analysis, we removed SNPs that 1) deviated from Hardy-Weinberg Equilibrium (Fisher exact test, P < 0.01) or 2) linked loci (each pair of SNP have correlation coefficient r² > 0.33 in a sliding window with a size of 50 SNPs and step of 5 SNPs), using VCFtools (v0.1.13)⁹⁷ and Plink (v1.9)¹⁰¹.

Principal component analysis (PCA) and population structure analysis were carried out using Plink (v1.9)¹⁰¹ and fastStructure (v1.0)¹⁰², respectively. The simple mode (as default) from fastStructure was used for the population structure analysis. The K value was estimated using a heuristic function in fastStructure.

For each of the 159 clonal families, we selected the least missingness genotype (i.e. the genotype with the highest sequencing coverage of that clonal family) as the representative genotype. SNP information from all 159 representative genotypes was used to estimate the linkage disequilibrium decay for each of the four populations. PopLDdecay (v3.41)¹⁰³ was used to measure LD decay. For each population, we used the following filters: SNP of missing allele > 20% and MAF < 0.05. The allele frequency correlation (denoted as r²) of pairwise SNPs within 100 kb physical distance was calculated.

Phylogenetic tree reconstruction

We used BLAST + version 2.9.0¹⁰⁴ to identify orthologous fragments between the genomes of S. polyrhiza and Colocasia esculenta (Araceae). For each SNP from the core set, the reference allele and its flanking 300 bp (upstream 150 bp and downstream 150 bp, respectively) sequences were extracted from the S. polyrhiza genome and then aligned to the C. esculenta reference genome¹⁰⁵. The hit thresholds were set as 1) alignment identity > 70%; 2) e-value > 1e − 6; 3) minimum aligned sequence length ≥ 50 (the aligned sequence must cover SNP position); 4) keep the best hit; and 5) ignore short deletions from C. esculenta. The orthologous alleles from C. esculenta were used as the outgroup genotype. We identified only 13,120 SNPs that have orthologous fragments in the C. esculenta genome. Those data were further used to infer the maximum-likelihood (ML) phylogenetic tree using RAxML-ng (v1.0.1)¹⁰⁶. The best hit model was estimated to be ‘TVM + G4’ using Modeltest-ng (0.1.6)^107,108. The bootstrapping converged after 700 iterations of the ML tree search. ITOL v5¹⁰⁹ and the Python package ETE2¹¹⁰ were used for tree visualization.

Selection analysis

Genome-wide scans of selection were performed on all 20 chromosomes of all sampled populations. Selective sweeps were inferred by three programmes: RAiSD⁵², SweeD⁵³ and LASSI⁵⁴. RAiSD uses the µ statistic, which provides information on the SFS, LD, and genomic diversity to evaluate the presence of positive selection⁵². SweeD calculates the traditional composite likelihood ratio (CLR) to infer loci under selection⁵³. LASSI employs the T statistic, which uses a likelihood model based on the haplotype frequency spectrum to detect hard and soft sweeps⁵⁴. As recommended by the authors of LASSI, we selected the top 5% T scores as candidates for selection. For RAiSD and SweeD we selected the top 1% scores. After finding the common genes under selection according to all three programmes, we reported the genes that have orthologs in A. thaliana. The embryo lethal genes from A. thaliana¹¹¹ were used for the enrichment analysis.

To test for population/branch-specific signals of selection, we ran a composite likelihood ratio (CLR) approach as implemented in 3P-CLR¹¹². Briefly, this method uses three-population trees coupled with genomic data as input, from which patterns of linked allele frequency differentiation are calculated. By doing this, this algorithm can tell apart signals of selection that happened in either branch of the tree or in the ancestral lineage, as well as outputting the loci with the highest CLR¹¹². In our case, we used either a North America-Asia-Europe, or a North America-Asia-India population tree as input, and 3P-CLR output the CLR across windows along each chromosome in the S. polyrhiza genome. We then selected the top 1% windows for each branch of the input tree and reported the genes that are present in each window. To further validate the evidence of positive selection on the regions with the highest CLR, we ran scans of Tajima’s D and genomic diversity along the same windows and contrasted them with the same signal along the other population branches. We expect negative Tajima’s D and low genomic diversity values along the populations with high CLR values.

DNA methylation in S. polyrhiza

We selected five genotypes from each of the four populations (America, India, SE-Asia, Europe) for single-base whole-genome bisulphite sequencing (WGBS. The genotypes originated from distinct clonal families, except for two European genotypes that came from the same clonal family. (Supplementary Table 5).

FastQC (v0.11.5, https://www.bioinformatics.babraham.ac.uk/projects/fastqc/) was used to summarize statistics of the sequencing reads. Trimmomatic (Version 0.36)¹¹³ was used to filter out low-quality reads with the parameters “SLIDINGWINDOW: 4:15, LEADING:3, TRAILING:3, ILLUMINACLIP: adapter.fa: 2: 30: 10, MINLEN:36”. Bismark (version 0.16.3)¹¹⁴ was used to align bisulphite-treated reads to the S. polyrhiza reference genome^93,94. Identical reads aligned to the same genomic regions were deemed as duplicated reads and thus were removed. Cytosines covered by less than five sequencing reads were excluded from the study. Only after applying these filters the sequencing depth and coverage were then summarized. The sodium bisulphite non-conversion rate was calculated as the percentage of non-converted cytosines to all cytosines in the reads that mapped to the chloroplast genome¹¹⁵ (GenBank: JN160603.2). For each cytosine site, a binomial test was performed to determine if the cytosine was methylated. If the methylation frequency at the site was lower than the background, which was estimated as the non-conversion rate, then the site was considered unmethylated, and the reads supporting methylation at this site were excluded¹¹⁶.

We calculated two different methylation parameters: the proportion of methylated cytosines (mC methylation) and weighted methylation level (wML)¹¹⁶. For both parameters, only cytosines covered by more than four sequencing reads were involved in the calculation. Those cytosines with low reads supporting methylation but not passing the binomial test were considered as un-methylated cytosines. The mC proportion was calculated by dividing the number of methylated cytosines by the total number of cytosines. Genomic regional wML was calculated using the methylKit (v1.17.5)¹¹⁷ and the regioneR (v1.28.0)¹¹⁸, with input based on the cytosine report generated with the Bismark pipeline. Line plots that show the wML patterns across the gene body and transposable elements, as well as their 2 kb flanking regions, were generated using ViewBS (v0.1.11)¹¹⁹. The hierarchical clustering, based on the methylation profiles’ similarity, was done using methylKit. The comparison between the genetic phylogenetic tree and hierarchical clustering based on the methylome was made using the R packages ggtree (v3.4.4)¹²⁰, treeio (v1.20.2)¹²¹, ape (v5.6.2)¹²², and phytools (v1.2.0)¹²³.

Data availability

The raw genomic and bisulfite sequencing reads involved in this study can be retrieved from NCBI under Accession numbers Bioproject PRJNA701543 and Bioproject PRJNA934173. The scripts for the data analyses are deposited in https://github.com/Xu-lab-Evolution/Great_duckweed_popg. The authors declare that the data and corresponding computational codes supporting the conclusions of this study are available within the article, its Supplementary Information files, or upon request.

Acknowledgements

We thank Dr Martin Schäfer, Marie Sárazová and Laura Böttner for supporting plant sample maintenance and relevant wet-lab experiments. This project is supported by German Research Foundation (427577435 and 438887884 to S. X and 422213951 to M.H), the Center for Adaptation to a Changing Environment (ACE) at ETH Zurich (to S. X. and A. W.), the Swiss National Science Foundation (P400PB_186770 to M.H.), the Volkswagen Foundation (97236 to M.H.) and through career development measures of the University of Münster (to M.H.) The project was inspired by discussions with the members of the CRC TRR 212 (NC3) – Project number316099922, and Research Training Group 2526 (GenEvo) – Project number 407023052. Parts of this research were conducted using the supercomputer Mogon and/or advisory services offered by the Johannes Gutenberg University Mainz (hpc.uni-mainz.de), which is a member of the AHRP (Alliance for High Performance Computing in Rhineland Palatinate, www.ahrp.info) and the Gauss Alliance e.V. The authors gratefully acknowledge the computing time granted on the supercomputer Mogon at the Johannes Gutenberg University Mainz (hpc.uni-mainz.de) and PALMA-II at the University of Münster.

Author contributions

Y. W., P. D. and S. X. performed data analysis. A. C. performed the experiments. K. J. A., H. Z., A. W., K. S. S., and S.X. contributed to the giant duckweed collections and resources. S. X. and M.H. conceived and supervised the project. S. X., Y. W., and P. D. wrote the manuscript. All authors contributed to the final version of the manuscript.

Competing interests:

The authors declare no competing interests.

Kondrashov, A. S. Deleterious Mutations and the Evolution of Sexual Reproduction. Nature 336, 435–440 (1988).
Muller, H. J. The relation of recombination to mutational advance. Mutat Res 1, 2–9 (1964).
Case, T. J. & Taper, M. L. On the Coexistence and Coevolution of Asexual and Sexual Competitors. Evolution 40, 366–387 (1986).
Doncaster, C. P., Pound, G. E. & Cox, S. J. The ecological cost of sex. Nature 404, 281–285 (2000).
Hartfield, M. Evolutionary genetic consequences of facultative sex and outcrossing. J Evolution Biol 29, 5–22 (2016).
Green, R. F. & Noakes, D. L. G. Is a Little Bit of Sex as Good as a Lot. J Theor Biol 174, 87–96 (1995).
Lynch, M. & Gabriel, W. Phenotypic Evolution and Parthenogenesis. Am Nat 122, 745–764 (1983).
Simon, J. C., Rispe, C. & Sunnucks, P. Ecology and evolution of sex in aphids. Trends Ecol Evol 17, 34–39 (2002).
Hebert, P. D. N. Population Biology of Daphnia (Crustacea, Daphnidae). Biol Rev 53, 387–426 (1978).
Wallace, R. L. Rotifers: Exquisite metazoans. Integr Comp Biol 42, 660–667 (2002).
Klimeš, L., Klimešová, J., Hendriks, R. & van Groenendael, J. Clonal plant architecture: a comparative analysis of form and function. The Ecology and Evolution of Clonal Plants, 1–29 (1997).
de Meeus, T., Prugnolle, F. & Agnew, P. Asexual reproduction: Genetics and evolutionary aspects. Cell Mol Life Sci 64, 1355–1372 (2007).
Keightley, P. D. & Otto, S. P. Interference among deleterious mutations favours sex and recombination in finite populations. Nature 443, 89–92 (2006).
Jaron, K. S. et al. Convergent consequences of parthenogenesis on stick insect genomes. Sci Adv 8, eabg3842 (2022).
Tucker, A. E., Ackerman, M. S., Eads, B. D., Xu, S. & Lynch, M. Population-genomic insights into the evolutionary origin and fate of obligately asexual Daphnia pulex. Proc Natl Acad Sci U S A 110, 15740–15745 (2013).
Niederhuth, C. E. & Schmitz, R. J. Covering your bases: inheritance of DNA methylation in plant genomes. Mol Plant 7, 472–480 (2014).
Matzke, M. A. & Mosher, R. A. RNA-directed DNA methylation: an epigenetic pathway of increasing complexity. Nat Rev Genet 15, 394–408 (2014).
Nozawa, K. et al. DNA methyltransferase CHROMOMETHYLASE3 prevents ONSEN transposon silencing under heat stress. PLoS Genet 17, e1009710 (2021).
Gehring, M. Epigenetic dynamics during flowering plant reproduction: evidence for reprogramming? New Phytol. 224, 91–96 (2019).
She, W. et al. Chromatin reprogramming during the somatic-to-reproductive cell fate transition in plants. Development 140, 4008–4019 (2013).
She, W. J. & Baroux, C. Chromatin dynamics in pollen mother cells underpin a common scenario at the somatic-to-reproductive fate transition of both the male and female lineages in Arabidopsis. Front Plant Sci 6 (2015).
Slotkin, R. K. et al. Epigenetic reprogramming and small RNA silencing of transposable elements in pollen. Cell 136, 461–472 (2009).
Calarco, J. P. et al. Reprogramming of DNA methylation in pollen guides epigenetic inheritance via small RNA. Cell 151, 194–205 (2012).
Ingouff, M. et al. Live-cell analysis of DNA methylation during sexual reproduction in Arabidopsis reveals context and sex-specific dynamics controlled by noncanonical RdDM. Genes Dev. 31, 72–83 (2017).
Bouyer, D. et al. DNA methylation dynamics during early plant life. Genome Biol 18, 179 (2017).
Kawakatsu, T., Nery, J. R., Castanon, R. & Ecker, J. R. Dynamic DNA methylation reconfiguration during seed development and germination. Genome Biol 18, 171 (2017).
Narsai, R. et al. Extensive transcriptomic and epigenomic remodelling occurs during Arabidopsis thaliana germination. Genome Biol 18 (2017).
Verhoeven, K. J. F., Jansen, J. J., van Dijk, P. J. & Biere, A. Stress-induced DNA methylation changes and their heritability in asexual dandelions. New Phytologist 185, 1108–1118 (2010).
Verhoeven, K. J. & Preite, V. Epigenetic variation in asexually reproducing organisms. Evolution 68, 644–655 (2014).
Van Antro, M. et al. DNA methylation in clonal duckweed (Lemna minor L.) lineages reflects current and historical environmental exposures. Mol Ecol 32, 428–443 (2023).
Niederhuth, C. E. et al. Widespread natural variation of DNA methylation within angiosperms. Genome Biol 17 (2016).
Landolt, E., Jäger-Zürn, I. & Schnell, R. Extreme adaptations in angiospermous hydrophytes. (Gebrüder Borntraeger, 1998).
Bog, M., Appenroth, K. J. & Sree, K. S. Key to the determination of taxa of Lemnaceae: an update. Nord J Bot 38 (2020).
Kim, I. Structural Differentiation of the Connective Stalk in Spirodela polyrhiza (L.) Schleiden. Applied Microscopy 46, 83–88 (2016).
Hicks, L. E. Flower production in the Lemnaceae. Ohio J Sci 32, 115–132 (1932).
Fourounjian, P., Slovin, J. & Messing, J. Flowering and Seed Production across the Lemnaceae. INT J MOL SCI 22 (2021).
Xu, S. et al. Low genetic variation is associated with low mutation rate in the giant duckweed. Nat Commun 10, 1243 (2019).
Ho, E. K. H., Bartkowska, M., Wright, S. I. & Agrawal, A. F. Population genomics of the facultatively asexual duckweed Spirodela polyrhiza. New Phytol 224, 1361–1371 (2019).
Sandler, G., Bartkowska, M., Agrawal, A. F. & Wright, S. I. Estimation of the SNP mutation rate in two vegetatively propagating species of duckweed. G3-Genes Genom Genet 10, 4191–4200 (2020).
Michael, T. P. et al. Comprehensive definition of genome features in Spirodela polyrhiza by high-depth physical mapping and short-read DNA sequencing strategies. Plant J. 89, 617–635 (2017).
Harkess, A. et al. Unusual predominance of maintenance DNA methylation in Spirodela polyrhiza. bioRxiv, 2020.2012. 2003.410332 (2020).
Bog, M. et al. Strategies for Intraspecific Genotyping of Duckweed: Comparison of Five Orthogonal Methods Applied to the Giant Duckweed Spirodela polyrhiza. Plants (Basel) 11 (2022).
Chen, J., Glemin, S. & Lascoux, M. Genetic Diversity and the Efficacy of Purifying Selection across Plant and Animal Species. Mol Biol Evol 34, 1417–1428 (2017).
Zhou, Y. F. et al. The population genetics of structural variants in grapevine domestication. Nat Plants 5, 965–979 (2019).
Guan, J. et al. Genome structure variation analyses of peach reveal population dynamics and a 1.67 Mb causal inversion for fruit shape. Genome Biol 22 (2021).
Cao, J. et al. Whole-genome sequencing of multiple Arabidopsis thaliana populations. Nat Genet 43, 956–963 (2011).
Kang, I. H., Steffen, J. G., Portereiko, M. F., Lloyd, A. & Drews, G. N. The AGL62 MADS domain protein regulates cellularization during endosperm development in Arabidopsis. Plant Cell 20, 635–647 (2008).
Hoffmann, T. et al. The identification of type I MADS box genes as the upstream activators of an endosperm-specific invertase inhibitor in Arabidopsis. Bmc Plant Biol 22 (2022).
Lee, J. & Lee, I. Regulation and function of SOC1, a flowering pathway integrator. J Exp Bot 61, 2247–2254 (2010).
Cokus, S. J. et al. Shotgun bisulphite sequencing of the Arabidopsis genome reveals DNA methylation patterning. Nature 452, 215–219 (2008).
Ibanez, V. N. & Quadrana, L. Shaping inheritance: how distinct reproductive strategies influence DNA methylation memory in plants. Curr Opin Genet Dev 78 (2023).
Alachiotis, N. & Pavlidis, P. RAiSD detects positive selection based on multiple signatures of a selective sweep and SNP vectors. Commun Biol 1 (2018).
Pavlidis, P., Zivkovic, D., Stamatakis, A. & Alachiotis, N. SweeD: Likelihood-Based Detection of Selective Sweeps in Thousands of Genomes. Mol Biol Evol 30, 2224–2234 (2013).
Harris, A. M. & DeGiorgio, M. A Likelihood Approach for Uncovering Selective Sweep Signatures from Haplotype Data. Mol Biol Evol 37, 3023–3046 (2020).
Demko, V., Ako, E., Perroud, P. F., Quatrano, R. & Olsen, O. A. The phenotype of the CRINKLY4 deletion mutant of Physcomitrella patens suggests a broad role in developmental regulation in early land plants. Planta 244, 275–284 (2016).
Braud, C., Zheng, W. & Xiao, W. Identification and analysis of LNO1-like and AtGLE1-like nucleoporins in plants. Plant Signal Behav 8, e27376 (2013).
Zhao, H., Xing, D. & Li, Q. Q. Unique features of plant cleavage and polyadenylation specificity factor revealed by proteomic studies. Plant Physiol 151, 1546–1556 (2009).
Takatsuka, H., Umeda-Hara, C. & Umeda, M. Cyclin-dependent kinase-activating kinases CDKD;1 and CDKD;3 are essential for preserving mitotic activity in Arabidopsis thaliana. Plant J 82, 1004–1017 (2015).
Johnson, K. L., Kibble, N. A., Bacic, A. & Schultz, C. J. A fasciclin-like arabinogalactan-protein (FLA) mutant of Arabidopsis thaliana, fla1, shows defects in shoot regeneration. Plos One 6, e25154 (2011).
Zhu, M. et al. Robust organ size requires robust timing of initiation orchestrated by focused auxin and cytokinin signalling. Nat Plants 6, 686–698 (2020).
Zhao, H. et al. The Arabidopsis thaliana Nuclear Factor Y Transcription Factors. Front Plant Sci 7, 2045 (2016).
Chantha, S. C., Gray-Mitsumune, M., Houde, J. & Matton, D. P. The MIDASIN and NOTCHLESS genes are essential for female gametophyte development in Arabidopsis thaliana. Physiol Mol Biol Plants 16, 3–18 (2010).
Chen, X. et al. Full-length EFOP3 and EFOP4 proteins are essential for pollen intine development in Arabidopsis thaliana. Plant J (2023).
Zhou, Y. et al. Members of the ELMOD protein family specify formation of distinct aperture domains on the Arabidopsis pollen surface. Elife 10 (2021).
Jossier, M. et al. The Arabidopsis vacuolar anion transporter, AtCLCc, is involved in the regulation of stomatal movements and contributes to salt tolerance. Plant J 64, 563–576 (2010).
Gachomo, E. W., Jimenez-Lopez, J. C., Baptiste, L. J. & Kotchoni, S. O. GIGANTUS1 (GTS1), a member of Transducin/WD40 protein superfamily, controls seed germination, growth and biomass accumulation through ribosome-biogenesis protein interactions in Arabidopsis thaliana. Bmc Plant Biol 14, 37 (2014).
Skalitzky, C. A. et al. Plastids Contain a Second Sec Translocase System with Essential Functions. Plant Physiol 155, 354–369 (2011).
Jeon, Y., Ahn, H. K., Kang, Y. W. & Pai, H. S. Functional characterization of chloroplast-targeted RbgA GTPase in higher plants. Plant Mol Biol 95, 463–479 (2017).
McConnell, J. R. et al. Role of PHABULOSA and PHAVOLUTA in determining radial patterning in shoots. Nature 411, 709–713 (2001).
Schwenk, P. et al. Uncovering a novel function of the CCR4-NOT complex in phytochrome A-mediated light signalling in plants. Elife 10 (2021).
Farkas, I., Dombradi, V., Miskei, M., Szabados, L. & Koncz, C. Arabidopsis PPP family of serine/threonine phosphatases. Trends Plant Sci 12, 169–176 (2007).
Guo, Z. F., Wang, X. Y., Hu, Z. B., Wu, C. Y. & Shen, Z. G. The pentatricopeptide repeat protein GEND1 is required for root development and high temperature tolerance in Arabidopsis thaliana. Biochem Bioph Res Co 578, 63–69 (2021).
Mochizuki, S. et al. The Arabidopsis WAVY GROWTH 2 protein modulates root bending in response to environmental stimuli. Plant Cell 17, 537–547 (2005).
Liu, C. H. et al. Repair of DNA Damage Induced by the Cytidine Analog Zebularine Requires ATR and ATM in Arabidopsis. Plant Cell 27, 1788–1800 (2015).
Bleuyard, J. Y. & White, C. I. The Arabidopsis homologue of Xrcc3 plays an essential role in meiosis. Embo J 23, 439–449 (2004).
Lim, M. H. et al. A new Arabidopsis gene, FLK, encodes an RNA binding protein with K homology motifs and regulates flowering time via FLOWERING LOCUS C. Plant Cell 16, 731–740 (2004).
Disch, S. et al. The E3 ubiquitin ligase BIG BROTHER controls arabidopsis organ size in a dosage-dependent manner. Curr Biol 16, 272–279 (2006).
Li, H. F. et al. The AGL6-like gene OsMADS6 regulates floral organ and meristem identities in rice. Cell Res 20, 299–313 (2010).
Krizek, B. A. & Meyerowitz, E. M. The Arabidopsis homeotic genes APETALA3 and PISTILLATA are sufficient to provide the B class organ identity function. Development 122, 11–22 (1996).
Lee, S., Choi, S. C. & An, G. Rice SVP-group MADS-box proteins, OsMADS22 and OsMADS55, are negative regulators of brassinosteroid responses. Plant J 54, 93–105 (2008).
Fang, W. J., Wang, Z. B., Cui, R. F., Li, J. & Li, Y. H. Maternal control of seed size by EOD3/CYP78A6 in Arabidopsis thaliana. Plant J 70, 929–939 (2012).
Sotelo-Silveira, M. et al. Cytochrome P450 CYP78A9 Is Involved in Arabidopsis Reproductive Development. Plant Physiol 162, 779–799 (2013).
Qi, X. L., Liu, C. L., Song, L. L., Li, Y. H. & Li, M. PaCYP78A9, a Cytochrome P450, Regulates Fruit Size in Sweet Cherry (Prunus avium L.). Front Plant Sci 8 (2017).
Ellegren, H. & Galtier, N. Determinants of genetic diversity. Nat Rev Genet 17, 422–433 (2016).
Theissen, G. & Theissen, G. Stranger than Fiction: Loss of MADS-Box Genes During Evolutionary Miniaturization of the Duckweed Body Plan Loss of MADS-Box Genes in Duckweeds. Compend Pl Genome, 91–101 (2020).
et al. Unusual predominance of maintenance DNA methylation in <em>Spirodela polyrhiza</em&gt.bioRxiv, 2020.2012.2003.410332 (2020).
Underwood, C. J. et al. Epigenetic activation of meiotic recombination near Arabidopsis thaliana centromeres via loss of H3K9me2 and non-CG DNA methylation. Genome Res. 28, 519–531 (2018).
Santamaria, L. Why are most aquatic plants widely distributed? Dispersal, clonal growth and small-scale heterogeneity in a stressful environment. Acta Oecol 23, 137–154 (2002).
Wang, Y. J. et al. Invasive alien plants benefit more from clonal integration in heterogeneous environments than natives. New Phytol 216, 1072–1078 (2017).
Gutekunst, J. et al. Clonal genome evolution and rapid invasive spread of the marbled crayfish. Nat Ecol Evol 2, 567–573 (2018).
Appenroth, K.-J. & Sree, K. in ISCDRA Duckweed Forum; Rutgers Duckweed Stock Cooperative: New Brunswick, NJ, USA. 180–186.
Schubert, M., Lindgreen, S. & Orlando, L. AdapterRemoval v2: rapid adapter trimming, identification, and read merging. BMC Res Notes 9, 88 (2016).
Wang, W. et al. The Spirodela polyrhiza genome reveals insights into its neotenous reduction fast growth and aquatic lifestyle. Nat Commun 5 (2014).
Cao, H. X. et al. The map-based genome sequence of Spirodela polyrhiza aligned with its chromosomes, a reference for karyotype evolution. New Phytol 209, 354–363 (2016).
Li, H. et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078–2079 (2009).
McKenna, A. et al. The Genome Analysis Toolkit: A MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res 20, 1297–1303 (2010).
Danecek, P. et al. The variant call format and VCFtools. Bioinformatics 27, 2156–2158 (2011).
Cingolani, P. et al. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome of Drosophila melanogaster strain w(1118); iso-2; iso-3. Fly 6, 80–92 (2012).
Nelson, C. W., Moncla, L. H. & Hughes, A. L. SNPGenie: estimating evolutionary parameters to detect natural selection using pooled next-generation sequencing data. Bioinformatics 31, 3709–3711 (2015).
Zhan, X., Hu, Y., Li, B., Abecasis, G. R. & Liu, D. J. RVTESTS: an efficient and comprehensive tool for rare variant association analysis using sequence data. Bioinformatics 32, 1423–1426 (2016).
Purcell, S. et al. PLINK: A tool set for whole-genome association and population-based linkage analyses. Am J Hum Genet 81, 559–575 (2007).
Raj, A., Stephens, M. & Pritchard, J. K. fastSTRUCTURE: Variational Inference of Population Structure in Large SNP Data Sets. Genetics 197, 573-U207 (2014).
Zhang, C., Dong, S. S., Xu, J. Y., He, W. M. & Yang, T. L. PopLDdecay: a fast and effective tool for linkage disequilibrium decay analysis based on variant call format files. Bioinformatics 35, 1786–1788 (2019).
Camacho, C. et al. BLAST plus: architecture and applications. Bmc Bioinformatics 10 (2009).
Yin, J. M. et al. A high-quality genome of taro (Colocasia esculenta(L.) Schott), one of the world's oldest crops. Mol Ecol Resour 21, 68–77 (2021).
Kozlov, A. M., Darriba, D., Flouri, T., Morel, B. & Stamatakis, A. RAxML-NG: a fast, scalable and user-friendly tool for maximum likelihood phylogenetic inference. Bioinformatics 35, 4453–4455 (2019).
Flouri, T. et al. The Phylogenetic Likelihood Library. Syst Biol 64, 356–362 (2015).
Darriba, D. et al. ModelTest-NG: A New and Scalable Tool for the Selection of DNA and Protein Evolutionary Models. Mol Biol Evol 37, 291–294 (2020).
Letunic, I. & Bork, P. Interactive Tree Of Life (iTOL) v4: recent updates and new developments. Nucleic Acids Res 47, W256-W259 (2019).
Huerta-Cepas, J., Dopazo, J. & Gabaldon, T. ETE: a python Environment for Tree Exploration. Bmc Bioinformatics 11 (2010).
Meinke, D. W. Genome-wide identification of EMBRYO-DEFECTIVE (EMB) genes required for growth and development in Arabidopsis. New Phytol 226, 306–325 (2020).
Racimo, F. Testing for Ancient Selection Using Cross-population Allele Frequency Differentiation. Genetics 202, 733–750 (2016).
Bolger, A. M., Lohse, M. & Usadel, B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30, 2114–2120 (2014).
Krueger, F. & Andrews, S. R. Bismark: a flexible aligner and methylation caller for Bisulfite-Seq applications. Bioinformatics 27, 1571–1572 (2011).
Wang, W. Q. & Messing, J. High-Throughput Sequencing of Three Lemnoideae (Duckweeds) Chloroplast Genomes from Total DNA. Plos One 6 (2011).
Schultz, M. D., Schmitz, R. J. & Ecker, J. R. 'Leveling' the playing field for analyses of single-base resolution DNA methylomes. Trends Genet 28, 583–585 (2012).
Akalin, A. et al. methylKit: a comprehensive R package for the analysis of genome-wide DNA methylation profiles. Genome Biol 13 (2012).
Gel, B. et al. regioneR: an R/Bioconductor package for the association analysis of genomic regions based on permutation tests. Bioinformatics 32, 289–291 (2016).
Huang, X. S., Zhang, S. L., Li, K. Q., Thimmapuram, J. & Xie, S. J. ViewBS: a powerful toolkit for visualization of high-throughput bisulfite sequencing data. Bioinformatics 34, 708–709 (2018).
Yu, G. C., Lam, T. T. Y., Zhu, H. C. & Guan, Y. Two Methods for Mapping and Visualizing Associated Data on Phylogeny Using Ggtree. Mol Biol Evol 35, 3041–3043 (2018).
Wang, L. G. et al. Treeio: An R Package for Phylogenetic Tree Input and Output with Richly Annotated and Associated Data. Mol Biol Evol 37, 599–603 (2020).
Paradis, E., Claude, J. & Strimmer, K. APE: Analyses of Phylogenetics and Evolution in R language. Bioinformatics 20, 289–290 (2004).
Revell, L. J. phytools: an R package for phylogenetic comparative biology (and other things). Methods Ecol Evol 3, 217–223 (2012).

There is NO Competing Interest.

Supplementarydata.zip
Supplementary Data 1-6
Supplementaryinfov6.docx
Supplementary Information of Population genomics and epigenomics provide insights into the evolution of facultative asexuality in plants

Download PDF

Journal Publication

published 16 May, 2024

Read the published version in Communications Biology →

Version 1

posted

You are reading this latest preprint version

Population genomics and epigenomics provide insights into the evolution of facultative asexuality in plants

Status:

Journal Publication

Version 1

Abstract

Figures

Introduction

Results

Determinants of genomic diversity among populations

The footprint of selection on the genome

Discussions

Methods

DNA sample preparation and sequencing

Single nucleotide variant calling

Population structure and linkage disequilibrium (LD)

Phylogenetic tree reconstruction

Selection analysis

Declarations

References

Additional Declarations

Supplementary Files

Status:

Journal Publication

Version 1