Natural Hybrid Origin of the Controversial “Species” Clematis × pinnata (Ranunculaceae) Based on Multidisciplinary Evidence

Interspecific hybridization is common and has often been viewed as a driving force of plant diversity. However, it raises taxonomic problems and thus impacts biodiversity estimation and biological conservation. Although previous molecular phylogenetic studies suggested that interspecific hybridization may be rather common in Clematis, and artificial hybridization has been widely applied to produce new Clematis cultivars for nearly two centuries, the issue of natural hybridization of Clematis has never been addressed in detail. In this study, we tested the hybrid origin of a mesophytic and cold-adapted vine species, Clematis pinnata, which is a rare and taxonomically controversial taxon endemic to northern China. Using field investigations, flow cytometry (FCM), phylogenomic analysis, morphological statistics, and niche modeling, we tested hybrid origin and species status of C. pinnata. The FCM results showed that all the tested species were homoploid (2n = 16). Phylonet and HyDe analyses based on transcriptome data showed the hybrid origins of C. × pinnata from either C. brevicaudata × C. heracleifolia or C. brevicaudata × C. tubulosa. The plastome phylogeny depicted that C. × pinnata in different sampling sites originated by different hybridization events. Morphological analysis showed intermediacy of C. × pinnata between its putative parental species in many qualitative and quantitative characters. Niche modeling results suggested that C. × pinnata had not been adapted to a novel ecological niche independent of its putative parents. These findings demonstrated that plants of C. × pinnata did not formed a self-evolved clade and should not be treated as a species. The present study also suggests that interspecific hybridization is a common mechanism in Clematis to generate diversity and variation, and it may play an important role in the evolution and diversification of this genus. Our study implies that morphological diversity caused by natural hybridization may overstate the real species diversity in Clematis.

Introduction

Natural hybridization between species is a long-lasting topic in evolutionary biology (Mallet, 2007; Soltis and Soltis, 2009; Yakimowski and Rieseberg, 2014; Goulet et al., 2017). It is well acknowledged that hybridization plays an important role promoting the diversification of plants (Soltis and Soltis, 2009). About 35% of vascular plant species were estimated to be the results of interbreeding of different species (Wood et al., 2009). By now, hybridization is recognized as an important evolutionary force and a remarkable portion of speciation, and study the process of hybridization can aid us to understand the origin of new adaptations and plant diversity (Li et al., 2017). However, natural hybridization presents great challenges for taxonomy, biodiversity estimation, and conservation because a reticulated phylogeny is forced into a hierarchical taxonomic system (Zhang et al., 2020; Draper et al., 2021).

Interbreeding of different species is only the first step for hybrid speciation (Zhang et al., 2020), and not all natural hybridizations finally produce new species. Hybridization may finally result in a new species only if the hybrid lineage could be established as viable progenies through vegetative (or clonal) propagation, or allopolyploidy events, or other homoploid speciation mechanisms (Comai, 2005; Sochor et al., 2015; White et al., 2018). However, studies have shown that many previously recognized plant species represented F1 hybrids, which cannot be accepted as a real species (Zha et al., 2010; Zhang et al., 2020; Li et al., 2017). Because F1 hybrids tend to have similar morphologies due to the complete combination of parental genomes, they have often been recognized as distinct species by taxonomists (Liao et al., 2021). This raises critical problems for morphological-based taxonomy, biodiversity estimations and biological conservation. Interspecific hybridization also have many evolutionary consequences, including the origin and transfer of adaptations, the blur of distinctive lineages, or the formation of maladaptive hybrids, that have great impact on biodiversity conservation (Draper et al., 2021). For example, Rosa pseudobanksiae, previous recognized as an endangered species (category: CR), has recently been tested to be mostly F1 hybrids (Zhang et al., 2020). For this reason, it is critical not only to clarify whether a plant has hybrid origins, but also to know whether the taxon holds species status or not.

Clematis L. (Ranunculaceae) is a horticulturally important genus in the buttercup family (Ranunculaceae) with a world-wide distribution. Plants of Clematis are herbaceous or woody vines and, rarely, erect shrubs, or perennial herbs (Tamura, 1987, 1995; Wang and Li, 2005). Clematis species are mostly diploid (2n = 16), with only a few polyploid species (Tamura, 1995). Taxonomy of Clematis has been considered to be very challenging. There are great differences in species estimation of this genus, ranging from 240 (Tamura, 1995), 297 (Grey-Wilson, 2000), 320 (Johnson, 1997), to 354 (Wang and Li, 2005).

Recent molecular phylogenetic analyses showed extensive incongruence between nuclear and organellar phylogenies, suggesting that natural hybridization may be common among Clematis species (Miikeda et al., 2006; Xie et al., 2011; He et al., 2021). A recent molecular dating analyses based on complete plastid genome sequences (He et al., 2021) showed that Clematis diverged from its sister genus Anemoclema in early Miocene and its major clades evolved in the end of Miocene and early Pliocene. Species diversification of Clematis was estimated to be rather late mainly in Quaternary era, indicating recent species radiation in this genus (Xie et al., 2011; He et al., 2021).

The importance of hybridization in Clematis has been well known to horticulturists. Artificial hybridization can generate novel genotypic and phenotypic variants as well as novel ecological adaptations of great horticultural value, and has been widely applied in Clematis cultivar breeding for almost two centuries (Johnson, 1997; Toomey and Leeds, 2001). Artificial crossing between closely related Clematis species, or even between morphologically diverged taxa, can be easily done to produce various cultivars. By now, hundreds of Clematis cultivars have been produced through hybridization (Toomey and Leeds, 2001). For natural Clematis, recurrent hybridization events may be one of the key factors causing a great deal of trouble in its taxonomy. On the other hand, Clematis may also provide a new study systems to investigate the role and effects of hybridization on its global diversification. However, natural hybridization has never been addressed in detail for this genus.

Clematis pinnata Maxim. is considered to be a rare and narrowly distributed species in Beijing and adjacent areas of northern China. Plants of C. pinnata are mesophytic, light-loving (when mature) and cold-adapted creeping vines. Taxonomic opinions differ (Xie et al., 2005; Wang and Xie, 2007) with some authors placing C. pinnata into sect. Clematis (Maximowicz, 1877; Handel-Mazzetti, 1939; Fang, 1980; Ting, 1980; Johnson, 1997), while others grouped it into sect. Tubulosae (Kuntze, 1885; Wang, 2001; Wang and Xie, 2007). In fact, C. pinnata shows morphological intermediacy between C. heracleifolia DC. (sect. Tubulosae) and C. brevicaudata DC. (sect. Clematis) or C. tubulosa Turcz. (sect. Tubulosae) and C. brevicaudata (Table 1 and Figure 1).

TABLE 1

Table 1. A comparison of morphological characters of Clematis pinnata and its putative progenitors.

FIGURE 1

Figure 1. Pictures of Clematis pinnata and its putative parents and their representative pollen grains under an electron microscope. C. brevicaudata (A) and its pollen (F); C. pinnata (B) and its pollen (E); C. heracleifolia (C) and its pollen (H); C. tubulosa (D) and its pollen (G).

Clematis heracleifolia and C. tubulosa have very similar morphological characters and were often treated as the same species under the name C. heracleifolia (Rehder and Wilson, 1913; Handel-Mazzetti, 1939; Ting, 1980). However, Wang and Xie (2007) argued that C. heracleifolia and C. tubulosa can be distinguished from one another by their pedicels, sepal shapes, and pollen types (Xie and Li, 2012; Figure 1). So, they separated them into two different species. Our field survey showed that although in a wider geographical context, C. heracleifolia and C. tubulosa share a large, overlapping distribution range, their populations are often well separated from each other. Therefore, we treated them as two species, and when referring to them as one of the putative parent species we use C. heracleifolia/C. tubulosa hereafter.

Like C. pinnata, C. brevicaudata, and C. heracleifolia/ C. tubulosa are also mesophytic and cold-adapted species, but with much wider distribution ranges in China (Fang, 1980; Ting, 1980). However, their morphological features are quite different (Table 1). Clematis brevicaudata is climbing, and C. heracleifolia/C. tubulosa is erect perennial herb. Whereas, plants of C. pinnata are long and creeping but often not climbing. Clematis pinnata’s leaves are usually pinnate with five leaflets, albeit ternate or bi-ternate leaves are occasionally present on the upper or basal parts of individuals (Xie et al., 2005). In contrast, C. brevicaudata often carries bi-ternate leaves and C. heracleifolia/C. tubulosa always bears ternate leaves. The young flowers of C. pinnata are bluish-white and bell-shaped (similar to C. heracleifolia/C. tubulosa). However, when the flower is in full bloom, the sepals spread and fade to become almost white. At this stage, C. pinnata was often misidentified as C. brevicaudata. Clematis pinnata’s stamens are nearly glabrous and have inconspicuous hairs (Wang, 2001), which was also a misleading character in previous taxonomic treatments (Fang, 1980; Ting, 1980; Wang and Xie, 2007). Fore the four taxa, chromosome number of C. brevicaudata and C. heracleifolia have been reported to be diploid (2n = 2x = 16) (Gong et al., 1985; Zhang and He, 1990; Chung et al., 2013; Wang et al., 2017).

Despite their great morphological differences, Clematis sect. Clematis and sect. Tubulosae were shown to be very closely related by molecular phylogenetic studies (Miikeda et al., 2006; Xie et al., 2011; Lehtonen et al., 2016; He et al., 2021), thus indicating hybrid possibility between them. A Japanese species, C. takedana, has been considered as a natural hybrid species between C. apiifolia (sect. Clematis, distributed in Japan and eastern China) and C. stans (sect. Tubulosae, endemic to Japan) since its publication (Makino, 1907), solely based on its intermediate morphological characters between the latter two coexisting species. This case is very similar to C. pinnata. Clematis brevicaudata and C. heracleifolia/C. tubulosa always coexist in the present site of the morphologically intermediate C. pinnata. Furthermore, all of the four taxa have the same flowering time in late July to early September. High proportion of both aborted pollen and achene of C. pinnata have been reported by Shi (2003). She also assumed that C. pinnata may be a hybrid taxon between C. brevicaudata and C. heracleifolia. However, her allozyme, ITS, and trnL-F analyses did not give sufficient results for testing this assumption. From all the above-mentioned evidence, we hypothesize that C. pinnata is a homoploid hybrid between C. brevicaudata and C. heracleifolia/C. tubulosa, and its species status needs to be tested.

When testing the hybrid origin of a plant, a clearly defined species concept must be adopted. Although morphological-based species concept has been widely used in plants because of its practicability for taxonomic purposes (Soltis and Soltis, 2009), it cannot be applied for testing hybridization hypothesis. The biological species concept, which emphasizes reproductive isolation, has played a major role in views of plant taxa (Soltis and Soltis, 2009). However, when testing hybrid hypothesis, this species concept is not operationally useful for empiricists (Feliner et al., 2017). In this study, we adopted broader conceptual frameworks for the definition of homoploid hybrid speciation in accordance with those proposed by Mallet (2007), Abbott et al. (2013), and Feliner et al. (2017). A homoploid hybrid species should have a clear genetic evidence of hybridization without polyploidization, give rise to an established self-evolving lineage, and be morphologically and ecologically distinct (especially in view that C. pinnata has the same phenological period with the putative parents) from its progenitors.

Therefore, using a multidisciplinary data including flow cytometry, phylogenomics, morphological statistics, and ecological niche modeling analysis, we test the possible homoploid hybrid origin and species status of C. pinnata. We focus on questions about the origin and process of the homoploid hybridization of C. pinnata. Flow cytometry determines ploidy level of C. pinnata compared to the related species. PhyloNet and HyDe analyses using transcriptome data give genetic evidence of hybridization, whereas the complete plastome analyses examine the maternal inheritance of C. pinnata. Using morphological statistics, we assess whether C. pinnata is intermediate between its putative parents in quantitative characters. We also model the niche overlap between C. pinnata and its putative parental species to test whether C. pinnata acquired new ecological niche independent of its parents. The aims of this study are (1) to test whether C. pinnata is of homoploid hybrid origin, (2) to determine which species are involved in hybridization events, (3) to determine the direction of hybridization, and (4) to test species status of C. pinnata.

Materials and Methods

Study Locations and Field Investigation

Distribution range of C. pinnata is concentrated in Beijing and adjacent areas in Tianjin and Hebei Province, only with a few records in Liaoning Province (explained in Taxonomic Treatment). From 2017 to 2019, we surveyed C. pinnata over its entire distribution range and found multiple distribution sites in Beijing and Liaoning Province (Table 1). In all its distribution sites, both of its putative parental species are present (Table 2 and Supplementary Figure S1). We recorded GPS coordinates and altitude of each population, and collected leaf tissues as well as specimens for flow cytometry experiments, phylogenomic studies, and morphological analyses. Individuals of C. pinnata are rare, often with a few plants in each site. So, we collected all the individuals in each site. In total, we gathered 27 C. pinnata individuals and one dubious plant from eight collecting sites in Beijing and Liaoning Province. A dubious individual was collected in Baihuashan (BHS), Beijing. It was not in flower when collected. The plant was creeping with pinnate leaves, which was similar to C. pinnata. However, the leaflets were smaller than typical C. pinnata and very similar to those of C. brevicaudata. We included this dubious individual in our phylogenomic study to determine its identity.

TABLE 2

Table 2. Locations of the population sites of Clematis pinnata and its putative progenitors collected in this study.

From all the eight collecting sites, we gathered 94 individuals of C. brevicaudata. We also collected 91 individuals of C. tubulosa from five of the eight sites, and 49 individuals of C. heracleifolia from the rest of the three sites. The individuals were randomly collected in each site, and the distance between individuals is more than 20 meters. Furthermore, five other Clematis species (C. hexapetala Pall., C. acerifolia Maxim., C. ochotensis (Pall.) Poir, C. intricata Bunge, and C. fruticosa Turcz.) commonly found in northern China were also included in the phylogenomic analysis as outgroups. Voucher specimens were deposited in the Herbarium of Beijing Forestry University (BJFC, Supplementary Table S1).

Ploidy Level Detection

We used a flow cytometry (FCM) method (Doležel et al., 2007; Bourge et al., 2018) to determine whether C. pinnata is homoploid with its putative parents. We used samples from Jiufeng forest park to check if C. pinnata is homoploid with the putative parents in the same site. Because only C. brevicaudata and C. tubulosa are present in Jiufeng population, we checked C. heracleifolia using samples from Yanqing district.

The FCM measurement was carried out in the Key Laboratory of Photobiology of Institute of Botany, the Chinese Academy of Sciences. Preparation of dried leaf samples of the four taxa basically followed a standard two-step method to isolate plant cell nuclei (Qu et al., 2018). For each sample, about 20–30 mg of dried leaf was chopped in 1 mL of LB01 lysis buffer (15 mM Tris, 2 mM EDTA-Na₂, 0.5 mM spermine tetrahydrochloride, 80 mM KCl, 20 mM NaCl, 0.1% Triton X-100, 15 mM β-mercaptoethanol, pH 7.5) to release cell contents. The resulting culture was gently pipetted and filtered through a 400-mesh screen to remove cell debris. The samples were then stained with 50 μL 1 mg⋅mL^–1 propidium iodide (PI) and 50 μL 1 mg⋅mL^–1 RNase A in an ice bath for 10 min before being analyzed using a MoFlo XDP flow cytometer (Beckman Coulter Inc., United States). We measured about 5000 nuclei from each sample (individual) and performed three replicates (three individuals) for each species. PI conjugated to Alexa Fluor 488 (Molecular Probes) was excited with a 488-nm argon-ion laser and the fluorescence was detected using a microscope equipped with a 625/26-nm HQ bandpass filter. Mean channel positions and coefficients of variation (CVs) of the G0/G1 peaks were calculated using Summit 5.2 software (Beckman Coulter Inc., United States). We used C. brevicaudata as the external standard reference for our flow cytometry measurements.

Transcriptome Analysis

To examine the possible hybrid origin of C. pinnata, we generated transcriptomes from samples of the eight collecting sites. Fresh leaf materials of the Clematis samples were collected from field and quickly deposited in liquid nitrogen. For C. pinnata, leaf tissues of all the individuals (from 1 to 10 individuals per site, Table 2) in each site were mixed together, respectively. Population transcriptome data of C. brevicaudata and C. heracleifolia/C. tubulosa were represented by five individuals for each species in each site. Five other Clematis species (C. ochotensis, C. intricata, C. hexapetala, C. fruticosa, and C. acerifolia, one individual for each species) were included for transcriptomes analysis, and Anemoclema glaucifolium was chosen as the outgroup.

We extracted total RNA from sampled leaves by using TRIzol Reagent (Invitrogen, Thermo Fisher Scientific, Shanghai, China) and sent the RNA samples to Biomarker Technologies¹ for cDNA library preparation and Illumina 2 × 150 bp paired-end sequencing. About 6 Gb of raw reads were obtained for each sampled species. The raw reads for the 30 newly generated transcriptomes were deposited in the NCBI Sequence Read Archive under BioProject PRJNA657443 (Supplementary Table S1). Then, they were filtered and trimmed using Trimmomatic v. 0.39 (Bolger et al., 2014). The clean RNA-Seq reads were de novo assembled using Trinity v. 2.5.1 (Grabherr et al., 2011) with default parameters. Because the assembly had a large number of redundant transcripts, we kept the longest isoforms of the related contigs by using the “get_longest_isoform_seq_per_trinity_gene.pl” utility in Trinity to let each unigene represent a collection of expressed sequences that apparently came from the same transcription locus. We used CD-HIT (Li and Godzik, 2006; Fu et al., 2012) to remove redundant sequences from the unigene and TransDecoder v 5.0 to predict protein-coding regions (Haas et al., 2013).

Single-copy orthologous nuclear genes (SCOGs) can be determined and screened out by reference-based methods (Zhao et al., 2013; Zeng et al., 2014, 2017; Feng et al., 2019) and all-against-all alignment methods without references (Sveinsson et al., 2014; Nishikawa et al., 2015; Majeed et al., 2019). The second strategy takes much more computational time, but can find more complete single-copy orthologous dataset than the former strategy, especially when there are no close genomic references available (Nishikawa et al., 2015). In this study, our method followed that of Nishikawa et al. (2015), which was based on all-against-all alignment with BLASTP as implemented in Proteinortho v. 6 (Lechner et al., 2011). The definition of an orthologous relationship between proteins was: e ≤ 1 × 10⁵, identity ≥ 25%, and alignment coverage ≥ 50%. The extracted SCOGs were aligned using MAFFT (Katoh et al., 2005), and each alignment was processed with a Python script² to remove ambiguous aligned positions.

For the phylogenetic analysis of the SCOGs data, both coalescent- and concatenation-based methods were applied. The coalescent-based method was conducted by the Accurate Species Tree Algorithm v. 4.4.4 (ASTRAL; Mirarab et al., 2014). Single-gene ML tree was reconstructed in RAxML v. 8.1.17 (Stamatakis, 2014) under the GTR + G model as suggested in the manual. Then, we used TreeShrink (Mai and Mirarab, 2018) to test whether long branch attraction caused incorrect phylogenetic relationship and then deleted all the trees with suspicious patterns of branch length. All the remaining ML trees were then used as input trees for ASTRAL. For the concatenation method, the aligned SCOGs sequences were linked using Geneious Prime 2019 (Kearse et al., 2012). Then, we used a ML method for phylogeny reconstruction as described above. In view that the presence of hybrid taxa may cause problems in phylogenetic inference (Xu, 2000), we exclude putative hybrid taxa (C. pinnata and C. ochotensis inferred by Phylonet analysis, see below) and re-do the above mentioned phylogenetic analysis to confirm the backbone of Clematis phylogeny.

Because complex evolutionary scenarios such as hybridization events do not proceed in tree-like manners (Huson and Bryant, 2006), we used the NeighborNet method implemented in SplitsTree 4.11.3 (Huson and Bryant, 2006) to reconstruct phylogenetic networks for the concatenated dataset. We excluded insertions/deletions (indels) and used the K2P model (Kimura, 1980) for distance analysis and support values at each node were estimated by running 1000 bootstrap replicates.

Then, we analyzed potential hybridization of C. pinnata using pseudolikelihood approach with PhyloNet (Wen et al., 2018), a software that detects phylogenetic networks based on the multispecies coalescent model, to directly infer the hybridization process. Because short gene alignments could cause random errors in phylogenetic inference in PhyloNet analysis (McLean et al., 2019), we selected gene alignments that were at least 500 bp in length for analysis. We divided SCOGs alignments into eight population datasets because PhyloNet cannot test too many hybridization events in a phylogenetic tree due to computational restrictions (Wen et al., 2018; Morales-Briones et al., 2021). Each dataset contains nine taxa: C. pinnata, both its putative parents, five other Clematis species, and Anemoclema as the outgroup. All the separate SCOGs phylogenies were inferred using an ML method as described above. We adopted the maximum pseudolikelihood method to model incomplete lineage sorting and gene flow using individual gene trees with the command InferNetwork_MPL (Yu and Nakhleh, 2015). Each Network search was allowed 0–4 reticulations and the log likelihood score for each network was also inferred. The best number of hybridization events was selected by plotting the likelihood scores. A sharp likelihood score increase is expected until it reaches the best number, after which the score increases slowly. For the dubious individual in BHS population, we also conducted another PhyloNet analysis by excluding another putative hybrid taxa C. ochotensis for comparison.

To further verify hybrid origin of C. pinnata and possible gene flow between the putative progenitors, we used an alternative method similar to D-statistic (or ABBA-BABA tests) analysis (Patterson et al., 2012) for comparing the results with the PhyloNet analysis. In this study, we used HyDe v. 0.4.1a (Blischak et al., 2018) to check for hybridization and gene flow. This method considered both hybridization and coalescence in a unified framework and can be used to fast detect both current and ancient hybridization events using genomic SNP data (Kubatko and Chifman, 2019). It assesses statistical significance of hybridization by testing for hybridization among all possible triplet combinations of the sampled species. Since HyDe uses single nucleotide polymorphism (SNP) data, the SCOG alignments were concatenated to construct a super SNPs (complete SNP data, hereafter) matrix using Phyutility v. 2.2 (Smith and Dunn, 2008) with removing the missing sites. We then used HyDe with a Python script, run_hyde.py (Blischak et al., 2018) to test possible hybrid origin of C. pinnata from each of the eight collection sites. The software assumed that each one of the three tested species (C. pinnata, C. brevicaudata, C. heracleifolia/C. tubulosa) is of a hybrid origin of the other two species. Thus, in total, 169 possible hybrid triplet combinations in each collecting site were tested.

Linkage disequilibrium (LD) may bias the results when analyzing genomic SNP data (Nielsen and Signorovitch, 2003; Malomane et al., 2018). Studies have shown that using LD-pruned data can effectively correct such biases (Hoeffding et al., 2017; Malomane et al., 2018; Davenport et al., 2020). However, accurately pruning LD sites from total SNP dataset needs high quality whole genomic reference (at the chromosome level). Because Clematis has huge genome size with the reported mean 1C value of 10.48 pg (Min: 6.90 pg and Max: 15.80 pg³), there are still no whole genome references published to date. For this reason, we used another approach to relieve the influence of LD. We randomly selected one SNP from each SCOG alignment using a python script⁴ to obtain a reduced SNP dataset, and then the above mentioned HyDe analysis was also conducted using this dataset.

A Z-statistical test was conducted and probabilities of genetic contributions of putative parents (γ and 1–γ) were calculated. According to Meng and Kubatko’s (2009) hybrid model, the hybrid taxa is either sister to “P1” with probability (γ) or sister to “P2” with probability (1–γ), The null hypothesis was that when hybridization was absent, the expected value of γ should be 0.

Chloroplast Genome Analysis

Because chloroplast (cp) genomes inherit maternally in Ranunculaceae (Corriveau and Coleman, 1988), we used genome skimming data to assemble complete cp genomes to determine the maternal parent of C. pinnata. For all the eight collecting sites, one individual of C. pinnata and its putative parental species was chosen for DNA extraction and library construction. Five other Clematis and one Anemoclema species (Zhang et al., 2015; Jiang et al., 2017) were also included in the analysis (Table 2).

We used genomic DNA extraction kits (Tiangen Biotech Co. Ltd., Beijing) to extract total genomic DNAs from silica-dried leaves. Extracted DNA was then sent to Biomarker Technologies (Beijing, China) for library construction and next-generation sequencing (NGS). Paired-end reads of 2 × 150 bp were generated on an Illumina Hiseq 4000 genome analyzer platform and raw reads were filtered using the FASTX-Toolkit⁵ to obtain high-quality data by deleting adaptors and low-quality reads.

We then used the Map to Reference option in Geneious Prime 2019 (Kearse et al., 2012) and reference sequences (MG675223.1 and MG675222) to filter out cp reads. Putative cp reads were used for de novo assembly using Geneious Prime 2019 (Kearse et al., 2012) with a low sensitivity setting to reconstruct the complete cp sequence. Gaps were bridged using 20 replicates in a FineTuning step in Geneious Prime 2019 (Kearse et al., 2012). When a contig containing a large single copy (LSC), a small single copy (SSC), and an inverted repeat (IR) region was assembled, the other IR region was determined and attached to the contig using the Repeat Finder function in Geneious Prime 2019 (Kearse et al., 2012) to construct the complete cp genome sequence. The assembled plastome sequences were annotated using Plastid Genome Annotator (Qu et al., 2019), and then checked manually in Geneious Prime 2019 (Kearse et al., 2012). We used the Organellar Genome DRAW tool to illustrate the newly sequenced cp genomes (Lohse et al., 2013) with accession numbers from MT796599 to MT706622 (Supplementary Table S2).

The complete cp genome dataset was aligned using MAFFT v 6.833 (Katoh et al., 2005). Ambiguous alignments and sites with more than 80% missing data were deleted automatically using a Python script⁶. We used the maximum likelihood (ML) and Bayesian inference (BI) methods for phylogeny reconstruction because these methods are less sensitive to long-branch attraction artifacts than the parsimony method (Bergsten, 2005). The ML trees were generated using RAxML v 8.1.17 (Stamatakis, 2014) under the GTR + G model, with bootstrap percentages computed after 100 replicates.

Bayesian inference was performed with MrBayes v3.2.3 (Ronquist and Huelsenbeck, 2003). The models were tested and set for BI analysis according to our previous study (He et al., 2019). Two Markov chain Monte Carlo (MCMC) chains were run independently. Each consisted of three hot chains and one cold chain for 5,000,000 generations. The trees were sampled every 100 generations. The convergence of the Markov chain was tested by calculating the standard deviation value of split frequencies (less than 0.01). The first 20% of trees were removed as burn-in. The remaining trees were used to construct the consensus tree.

Morphological Analysis

Both hybrids and hybrid species are expected to fall in between parental progenitors in morphology (Zhang et al., 2020). To test hybridization hypothesis of C. pinnata, 16 qualitative and ten quantitative morphological characters were chosen for analysis (Supplementary Tables S3, S4). In total, 278 specimens of the four tested taxa (C. pinnata: 29; C. brevicaudata: 124, C. heracleifolia: 56, and C. tubulosa: 69) from multiple herbaria were measured (information of all the specimens is presented in Supplementary Table S5). The selected qualitative characters are taxonomically important traits based on previous taxonomic study (Wang and Xie, 2007). We measured more than 90 leaflets, stems, and flowers of each species to retrieve the quantitative characters.

For the 10 quantitative characters, Kolmogorov–Smirnov test and F-test were used for normality test. Normal distribution data was used one-way ANOVA (analysis of variance) to analysis, and non-normal distribution data was used Kruskal–Wallis test to analysis via IBM SPSS statistics software v25 (SPSS Inc., Chicago, IL, United States), to test the significant differences among the four tested species. Significant differences between any two taxa were identified using a post hoc Tukey’s honest significant difference (HSD) test with false discovery rate (FDR) correction (Benjamini and Hochberg, 1995). Box charts and principal component analysis (PCA) were applied to visualize the differences between all the taxonomic units by R (R Core Team, 2018).

Ecological Niche Modeling

To determine whether C. pinnata has its own ecological niche independent of its parents, we used species distribution models (SDMs) with 278 distribution records (Supplementary Table S5) (C. pinnata: 29; C. brevicaudata: 124, C. heracleifolia: 56, and C. tubulosa: 69 from their entire distribution ranges) and 33 high-resolution environmental variables to project and predict the current distribution patterns of each taxon. The distribution records were obtained from our field investigation and the Chinese Virtual Herbarium database (CVH⁷). For CVH data, the identification of specimen was checked based on the previous taxonomic revision (Wang and Xie, 2007). When no exact GPS information was recorded for a specimen, the geographic coordinates were determined using Google Earth 7.0⁸. In order to avoid sampling bias, only one individual was retained in each 1.0 × 1.0 km square using the “spatially rarify occurrence data” tool in the SDM toolbox for the ArcGis software (Esri, Redlands, CA, United States).

Because many ecological factors, such as bioclimatic, sun light, vegetation, terrain, and soil, may influence the distribution of the four tested species, we obtained 33 environmental variables with 2.5′ spatial resolution (Supplementary Table S6) to conduct niche modeling. Those variables consist of 19 bioclimatic variables from WorldClim-Global Climate Data⁹ (Hijmans et al., 2005); six UV-B radiation variables from a global UV-B radiation dataset for macroecological studies¹⁰ (Beckmann et al., 2014); five vegetation and terrain variables from Harmonized World Soil Database v 1.2 of the Food and Agriculture Organization¹¹ (Food and Agriculture Organization [FAO] and International Institute for Applied Systems Analysis [IIASA], 2012), and three soil variables from the University of Wisconsin¹² (New et al., 1999).

Because collinearity between environmental variables can lead to wrong modeling results (Dormann et al., 2013), we applied the Pearson correlation coefficients using the “banding collection statistics” tool in ArcGis 10.2 and the Jackknife test implemented in MaxEnt to evaluate the contributions of each variable. All the environment variables were converted to ASCII format using an ArcGIS 10.2 conversion tool. Next, we obtained a map of Asia from the global administrative region database¹³. After removing overly correlated variables (r < 0.7) (Guisan and Thuiller, 2005), we performed predictive analysis for each species using MaxEnt 3.4.1¹⁴. We ran 10 bootstrap replicates, in which 25% was used for model testing and the other 75% of the presence data was randomly selected for model training to optimize the model. The model was parameterized with a maximum of 10,000 background points, a convergence threshold of 0.00001, and a maximum of 500 interactions (Lin et al., 2020).

We used area under the curve (AUC) values to evaluate model prediction accuracies. An AUC value is the area enclosed by the receiver operating characteristic curve and the abscissa, and the closer the value is to 1, the more predictively accurate the model is (Fielding and Bell, 1997; Babar et al., 2012). We used the jackknife method to evaluate the weight of each environmental factor on the distribution area, and based on these values, we classified potentially suitable habitats into four categories (which was widely applied by other niche modeling studies, such as Tang et al., 2017; Shitara et al., 2018; Lin et al., 2020, and many others) in a final predictive map of species: unsuitable (<0.2), barely suitable (0.2–0.4), moderately suitable (0.4–0.6), highly suitable (>0.6).

Then, the ecological niche overlap between tested species was checked following the methods by White et al. (2018). We tested niche overlap of species pairs using the ENMtools 1.4.4 (Warren et al., 2010) by calculating Schoener’s D (Schoener, 1968) and Warren’s I statistic values, which 0 indicates no overlap and 1 indicates full overlap (Warren et al., 2008). We then used a niche equivalence test species pairs by comparing statistics D and I using 100 pseudo-replicates (Warren et al., 2008) to test whether the species pairs have identical ecological niche modeling under the null hypothesis. When the statistic value is smaller than the pseudo-replicates value, the result indicates that the two species did not occupy the same ecological niche.

Next, to assess the extent of niche overlap between species pairs, we applied environmental PCA (PCA-env) method (Broennimann et al., 2012) as implemented in the R package of ecospat (Di Cola et al., 2017; R Core Team, 2018)¹⁵. We transformed the multidimensional space of environment variables into two-dimensional space by means of PCA. Following White et al. (2018), we set the resolution to 100 with each grid corresponding to a unique environmental space as suggested by Broennimann et al. (2012), and then we calculated the smoothness of species occurrence by using kernel density to project species onto a grid. We used Schoener’s D (Schoener, 1968) and Warren’s I statistics (Warren et al., 2008) to calculate niche overlap and we accepted the null hypothesis (two species’ niches were equivalent) when p < 0.05.

Results

Ploidy Level Detection

Our FCM results showed that the ratio of the mean (G₀/G₁) of C. × pinnata, C. tubulosa, and C. heracleifolia to that of the reference (C. brevicaudata) ranged between 0.95 and 1.22. These results showed that all tested species were of the same ploidy level, i.e., 2n = 16 (Supplementary Tables S7, S8 and Supplementary Figure S2).

Transcriptome Data Analysis

We generated 30 transcriptome datasets from the four tested species and other outgroup species in this study. Detailed information of the transcriptome data was presented in Supplementary Table S1. Because the sampled species are closely related, we filtered out 3198 SCOGs without missing data, using all-against-all alignment strategy. Among them, 28 SCOGs were tested to generate bad trees with suspicious branch length by TreeShrink. So, we kept 3170 SCOGs for phylogenetic reconstruction. Our results showed that the phylogenetic backbone without putative hybrid taxa (C. × pinnata and C. ochotensis, Supplementary Figure S3) is fully consistent with the results from phylogenetic analysis with inclusion of the hybrid taxa (Figure 2). Clematis × pinnata, C. brevicaudata, C. heracleifolia, and C. tubulosa formed a well-supported clade by both concatenated (ML BS = 100) and coalescent methods (Local posterior probabilities, LPP = 1, Figure 2). However, the resolution and support values within this clade are low, and samples of C. × pinnata were not tested to be monophyletic within the clade. The network analysis (Figure 3) showed that C. brevicaudata and C. heracleifolia/C. tubulosa are well separated. C. heracleifolia and C. tubulosa are closely related but formed two strains, respectively. Samples of C. × pinnata did not form a single strain in this analysis.

FIGURE 2

Figure 2. Phylogenetic relationships of eight populations of Clematis pinnata and its close relatives inferred by 3170 co-orthologous nuclear genes using concatenated and coalescent methods. (A) Concatenated phylogeny inferred using maximum likelihood method. Bootstrap percentages are indicated on the branches. * Shows that ML bootstrap values are 100, while – shows that support values are less than 50. (B) Coalescence based species tree, inferred by ASTRAL. Numbers at branches are local posterior probabilities. Local posterior probabilities with values equal to 1.00 were marked with * at branches. Population location name abbreviations are explained in Table 2.

FIGURE 3

Figure 3. Molecular phylogenetic network among Clematis pinnata and its putative parents from co-orthologous nuclear gene dataset (3170 SCOGs). Distance analysis was ran under K2P model (Kimura, 1980) and support values were estimated by running 1000 bootstrap replicates.

We selected 3000 SCOGs with aligned length longer than 500 bp for PhyloNet analysis. All the eight tested populations showed a gently increasing probability value when the maximum number of reticulations was set to two (Supplementary Figure S4). For this reason, we displayed and discussed the results based on the setting of two maximum hybridization events (Figure 4). The results showed that, except for the dubious sample in population BHS, samples of C. × pinnata from the other seven populations were tested to be hybrids between C. brevicaudata and C. heracleifolia/C. tubulosa. Another possible hybridization event occurred in C. ochotensis which belongs to C. sect. Atragene (Figure 4). For population BHS, the dubious sample was not tested to be a hybrid between C. brevicaudata and C. tubulosa but more closely related to C. brevicaudata in the PhyloNet analyses either with (Figure 4) or without C. ochotensis (Supplementary Figure S5).

FIGURE 4

Figure 4. Optimal species networks of all eight populations (A-H, population abbreviations are explained in Table 2) of C. pinnata and its putative parents inferred by 3000 co-orthologous nuclear genes (at least 500 bp in length) using the PhyloNet. The results are maximum pseudolikelihood trees with a maximum of two reticulations. Genetic contributions of the putative parents estimated by HyDe analysis (using reduced SNP data) were marked beside the curve lines.

In this study, we obtained 316,066 and 3170 SNPs for the complete and reduced SNP datasets, respectively, for HyDe analysis. In general, the two datasets generated similar results (Table 3 and Supplementary Table S9), and hybrid origin of C. × pinnata was confirmed by this analysis. The results from the reduced SNP data (as well as the complete SNP data) showed that the genetic contributions of the putative parents were nearly equal in six of the eight populations (Table 3). In WGY, the parental genetic contribution of C. brevicaudata was 0.33, while genetic contribution (1-γ) of C. heracleifolia was 0.67. For the dubious sample from population BHS, HyDe analysis for the complete SNP data detected weak possible hybrid origin of the dubious plant between C. brevicaudata and C. tubulosa (positive z-score = 3.74) with very low genetic contribution of C. tubulosa (1-γ = 0.01, Supplementary Table S9). However, analysis of the reduced SNP data did not support hybrid origin of the dubious plant (Table 3) which is consistent with the PhyloNet analysis.

TABLE 3

Table 3. HyDe analysis of Clematis pinnata and its putative parents in each population using reduced SNP dataset.

Chloroplast Genome Analysis

We obtained about 4 Gb of clean NGS data from each sample and filtered 132,798–443,873 cp reads for C. × pinnata samples and its close relatives by de novo assembly. Chloroplast genome sizes of C. × pinnata, C. brevicaudata, C. heracleifolia, and C. tubulosa ranged from 159,597 bp (C. tubulosa from YMS) to 159,667 bp (C. brevicaudata from YMS) and the overall GC content of all the four species were around 38%. All the acquired plastome sequences consisted of a pair of IRs (31,038–31,087 bp), separated by an LSC (79,392–79,419 bp) and an SSC (18,093–18,187 bp) regions (Supplementary Figure S6). The cp genomes of all the four Clematis species encoded an identical set of 112 genes, including 18 genes with introns, 78 protein-coding, 29 tRNA and four rRNA genes, and 25 genes are in IR region. Structural variation of the newly sampled cp genomes, such as gene inversion/translocation and IR expansion, was similar to that previously reported for other Clematis species (Liu et al., 2018).

Chloroplast phylogenomic analysis showed that all the samples of C. × pinnata and its possible progenitors formed a strongly supported clade (ML BS = 100, PP = 1, Figure 5). Within this clade, the branch lengths of subclades and terminal branches are very short, and the resolution as well as support values within the clade are relatively low. Individuals of C. × pinnata from different sites were separated from each other. Some C. × pinnata samples grouped with C. brevicaudata and clustered into a subclade (ML BS = 99, PP = 1), whereas others grouped with C. tubulosa and C. heracleifolia and were paraphyletic to the C. brevicaudata-C. × pinnata subclade. Unlike the nuclear phylogeny, the cp genomic phylogeny showed that both C. heracleifolia and C. tubulosa are not monophyletic.

FIGURE 5

Figure 5. Phylogeny of eight populations (one individual for each population) of Clematis pinnata and its close relatives inferred from complete plastome sequences using maximum likelihood and Bayesian inference methods. Bayesian phylogram is shown left below and the dendrogram with ML bootstrap values/posterior probability values are shown at right side. Population location name abbreviations are explained in Table 2. – shows that support values of ML bootstrap values are less than 50 and Bayesian PP values are less than 0.95.

Morphological Analyses

Qualitative characters suggested that C. × pinnata exhibited uniparental character states (e.g., hairs on both sides of lamina are the same with C. brevicaudata, hair of inside sepal is more similar to C. heracleifolia/C. tubulosa), intermediate character states (e.g., hair on stamens, sepal color, and spreading direction of the sepals), and new character states to its progenitors (e.g., leaf type) (Supplementary Table S3).

Eight of the ten quantitative characters showed that C. × pinnata measures fell in between C. brevicaudata and C. heracleifolia/C. tubulosa (Figures 6A–D,F,H–J). Whereas, the other two characters, pedicel and filament lengths, were larger in C. × pinnata than in its putative progenitors (Figures 6E,G). PCA results showed that C. × pinnata’s quantitative characters are intermediate to either C. brevicaudata and C. heracleifolia or to C. brevicaudata and C. tubulosa (Figures 6K,L).

FIGURE 6

Figure 6. Boxplots (A–J) and PCA plots (K,L) based on ten quantitative morphological characters of Clematis pinnata and its putative parents. (A–J) In each boxplot, groups sharing letters denote no significant difference between those groups, as identified by post hoc Tukey’s tests with false discovery rate corrections (Benjamini and Hochberg, 1995) for multiple comparisons. (K,L) Venn diagrams of PCA results show C. pinnata’s (b, 29 specimens) quantitative characters are intermediate to either C. brevicaudata (a, 124 specimens) and C. heracleifolia (c, 56 specimens) or to C. brevicaudata and C. tubulosa (d, 69 specimens).

Ecological Niche Modeling and Niche Overlapping Test

After removing overly correlated variables (r < 0.7) (Guisan and Thuiller, 2005), we reserved nine ecological variables for niche modeling (Supplementary Table S6). The AUC values of the four tested species ranged between 0.986 (C. brevicaudata) and 0.998 (C. × pinnata), indicating that the models had high levels of accuracy and validity. The predicted suitable habitats of the four species are shown in Supplementary Figure S7. Jackknife results (Supplementary Table S10) revealed that UVB4 and BIO18 were important for all the species, whereas BIO8 and BIO19 were only important for C. brevicaudata; UVB1 was important for C. × pinnata, C. heracleifolia and C. tubulosa; BIO15 was a key variable for C. × pinnata and C. tubulosa; and Soil pH was important for C. heracleifolia. The top four key environmental variables for C. × pinnata were the same as those for C. tubulosa.

Both niche equivalency test and PCA-env results rejected the null hypothesis that ecological niches for C. brevicaudata vs. C. heracleifolia and C. brevicaudata vs. C. tubulosa were overlapping (p < 0.05). Whereas, the hypothesis of niches overlapping of C. × pinnata vs. C. heracleifolia, C. × pinnata vs. C. tubulosa, C. × pinnata vs. C. brevicaudata, and C. heracleifolia vs. C. tubulosa, was accepted by the two statistical analyses (Table 4).

TABLE 4

Table 4. Niche overlap statistics between tested species pairs based on PCA-env and the equivalency test using Schoener’s D (Schoener, 1968) and Warren’s I (Warren et al., 2008).

Discussion

In this study, we present the first comprehensive study on the natural hybrid origin of Clematis species. Testing hybrid hypothesis and species status of C. × pinnata is the first important step for resolving its taxonomic problems, and hence obtains deeper insights into the diversification of the genus.

Phylogenomic Analysis Confirms Hybrid Origin of Clematis × pinnata

The application of genomic data has become an effective means to validate hybridization events. Many natural hybridization events have been proposed and validated in numerous organisms (Roberts and Roalson, 2018; Cao et al., 2019; Glémin et al., 2019; Wang et al., 2019; Zhang et al., 2019; Yang et al., 2020). In this study, our nuclear and cp genome analyses confirmed hybrid origin of C. × pinnata between C. brevicaudata and C. heracleifolia/C. tubulosa (Figures 2–5). Hybrid origin of C. × pinnata was further validated by PhyloNet and HyDe analyses, which have been designed for inferring hybridization (Blischak et al., 2018; Wen et al., 2018). Using populational transcriptome data, we not only confirmed hybrid origin of C. × pinnata but also obtained its parental genetic contributions in all the sampling sites except BHS (Table 3). From HyDe analysis, the genetic contributions of the progenitors of most populations were nearly equal, thus strongly suggesting that most plants of C. × pinnata represent the F1 of the parents, because segregation may occur in F2 or later generations causing significant biased parental contributions (Soltis and Soltis, 2009). Meanwhile, for the dubious sample from BHS, PhyloNet analysis failed to detect hybridization of this sample and that it clustered closely with C. brevicaudata. This result was further confirmed by HyDe analysis using the reduced SNP data (Table 3). The sample from BHS may represent a morphologically variable individual of C. brevicaudata

There are some methodological limitations in this study that should be concerned. The complete SNP data may generate biased results due to linkage disequilibrium (LD). In this study, our strategy (randomly select one SNP from each SCOG) greatly relieved the influence of LD. However, because there are no high-quality whole genome references for Clematis, and some of the 3170 SCOGs may be still physically closely related, it is possible that the effects of LD still exist. Furthermore, separating hybridization signal from other sources of incongruence, such as incomplete lineage sorting (ILS) and ancestral polymorphism in population, is difficult (Maddison, 1997). The ABBA-BABA tests implemented in HyDe analysis (or D-statistic test, Patterson et al., 2012) may fail to distinguish signal between hybridization and ancient population structure using SNP data, because this test do not test for hybridization per se, but for deviations from the Wright–Fisher model of random mating which can occur as a result from many different evolutionary processes (Huynh et al., 2019).

Another issue is concerned with our RNA-seq sampling. In the eight collecting sites, we collected multiple C. × pinnata individuals from five sites (Table 2). These multiple individuals from each site were mixed together for RNA-seq to represent genetic features in each population. However, in view that each C. × pinnata individual may have evolved independently, mixed samples may obtain inaccurate or incorrect result, especially when backcrossing occurred in some individuals. Further analysis using transcriptome data from separate individuals should be conducted to review the impact of this mixed sampling strategy.

Complete Chloroplast Genome Analysis Indicate Multiple Hybrid Directions

In recent years, complete cp genomes have been widely used in plant phylogenetic reconstructions (e.g., Moore et al., 2010; Li et al., 2019). Because cp genomes show uniparental inheritance in most angiosperm species, it is also a good molecular marker for inferring hybridization and introgression (Liu et al., 2020). In this study, the cp genome analysis demonstrated that multiple C. × pinnata samples do not form a single lineage (Figures 2, 3, 5). From cp genome phylogeny, we can clearly tell which one acted as the maternal parent of C. × pinnata (Figure 5). The result showed that all the three putative parental species may have contributed maternally to C. × pinnata. This results also depicted that different C. × pinnata individuals evolved from different hybrid events.

Causes and Consequences of Natural Hybridization Between Clematis Species

From our FCM and phylogenomic results, individuals of C. × pinnata are homoploid hybrids between C. brevicaudata and C. heracleifolia/C. tubulosa. Different hybridization events happened in different distribution sites, and plants of C. × pinnata may be predominantly F1 individuals. This conclusion is also supported by several lines of evidence. Our field investigations found that C. × pinnata’s individuals always occur in the places where both C. brevicaudata and C. heracleifolia/C. tubulosa are present. The flowering time of all four taxa overlap from July to September, thus creating the potential for cross-pollination. For this reason, if C. × pinnata is a species, it cannot easily survive from the parental introgression. In this condition, homoploid hybrid speciation is often accompanied by ecological isolation between daughter species and its parental species (Liu et al., 2014; Kadereit, 2015; White et al., 2018). However, the niche modeling and the niche equivalency test results demonstrated that C. × pinnata has not adapted to a new ecological niche independent of its parents’. From long-term observation, we found that C. × pinnata suffers from reduced fertility of a high proportion of pollen and fruit abortion, which was also reported by Shi (2003). For example, a C. × pinnata individual in Jiufeng, Beijing, noted only a few normal achenes per year, while most of them did not develop at all.

Except for several prerequisites for natural hybridization, e.g., close kinship, overlapping distribution, similar flowering period, shared pollinators, and same chromosome numbers (Ning et al., 2019), habitat disturbance has been often considered as one of the most important factors promoting hybridization (Arnold, 1997; Li et al., 2017; Liao et al., 2021). In the study area of Beijing and Liaoning Province, we found all the C. × pinnata individuals occurred along the mountain roadside. Whereas, C. brevicaudata and C. heracleifolia/C. tubulosa can occupy much larger distribution ranges than C. × pinnata. In autumn, the mountain roadside area will be cleared by forest workers for fire prevention in northern China. This provides opportunity for the parental species to contact each other, and open up habitat for the new hybrids (Li et al., 2017).

Recent studies have shown that F1 hybrids are common in angiosperms and that they can successfully impede gene flow and thus maintain species boundaries between parental species especially in areas where habitat disturbance is high (Liao et al., 2021). Morphologically, C. sect. Clematis and sect. Tubulosae have been considered distantly related in the genus due to their great morphological divergence. They were often placed into different subgenera due to their different floral characters (Tamura, 1995; Grey-Wilson, 2000). Only recently have molecular phylogenetic studies clarified their close relationship (Xie et al., 2011; Lehtonen et al., 2016; Yan et al., 2016; He et al., 2021). Representative species of these two sections, C. brevicaudata and C. heracleifolia/C. tubulosa shared part of their distribution areas in Beijing, Hebei, and Liaoning Provinces, which are exactly the distribution ranges of C. pinnata. The presence of C. × pinnata (F1 hybrid, with reduced fertility) can impede gene flow between parental species, and maintain species boundaries of C. brevicaudata and C. heracleifolia/C. tubulosa in their contact zones.

In Clematis, other species like C. ochotensis (sect. Atragene) also showed hybrid origin from other two different sections (sect. Montana and sect. Fruticella and/or sect. Meclatis, Figure 4). However, we neither have morphological evidence nor sufficient sampling to discuss about this hybridization event. This study demonstrates that interspecific hybridization between two morphologically highly diverged species can occur naturally in Clematis, and hybridization may play an important role in the evolution and diversification of the genus. Taxonomy of Clematis may have suffered from widely hybridization among morphologically diverged species.

Morphology, Species Status, and Taxonomy of Clematis × pinnata

This study provided opportunity to investigate how morphological characters of hybrids can be expressed in comparison with its morphologically highly diverged parents. Hybrids are often expected to be morphologically intermediate. However, morphological analyses of natural and artificial hybrids showed that characters of hybrids can be truly intermediate, or identical to those of either paternal or maternal parent, or even new traits (Soltis and Soltis, 2009). In this study, our morphological analysis showed that C. × pinnata exhibited all kinds of morphological outcomes, such as intermediate characters, uniparental characters, and new characters compared to its parental species (Supplementary Table S3 and Figure 6).

Clematis × pinnata shows extensive variation in morphological characters especially in its leaf types and shapes. Leaves of C. × pinnata are predominantly pinnate with five leaflets, but ternate, bi-ternate, or even simple leaves sometimes also occur in different development stages of the plants (Xie et al., 2005). Morphological characters often exhibit higher variability in hybrids than in hybrid species (Zhang et al., 2020). High morphological variation in C. × pinnata also supports the point that C. × pinnata represent early generation of hybrids rather than a hybrid species.

Our morphological analysis did not clearly distinguish C. heracleifolia from C. tubulosa. Except some qualitative characters, e.g., pedicel length, sepal shape, and pollen type (Wang and Xie, 2007), other morphological characters failed to distinguish the two species from one another. Furthermore, our genome analysis and niche equivalency test also did not clearly separate them. These results raise an interesting issue of species delimitation in C. sect. Tubulosae that need to be studied in the future.

All the previous taxonomic studies recognized C. × pinnata as a distinct species based solely on morphology (Fang, 1980; Xie et al., 2005; Wang and Xie, 2007). However, our study has clearly shown that plants of C. × pinnata has not formed a self-evolving lineage, and are generated by recurrent hybridization events between C. brevicaudata and C. heracleifolia/C. tubulosa in their overlapping zones. Clematis × pinnata cannot hold species status from our analysis, and we make a taxonomic treatment as below.

Taxonomic Treatment

Clematis × pinnata Maxim. (pro sp.) (= Clematis brevicaudata DC.: ♀or ♂ × C. heracleifolia DC.: ♀or ♂, C. brevicaudata DC.: ♀or ♂ × C. tubulosa Turcz.: ♀or ♂) in Bull. Acad. Imp. Sci. Saint-Pétersbourg, sér. 3. 22: 216. 1876 – Holotype: China. Near Beijing, ca. 1845, A.A. Tatarinov s. n. (LE!; isotype: PE!)

= C. tatarinowii Maxim., in l.c. Holotype: China. Beijing, ca. 1845, A.A. Tatarinov s. n. (LE; PE [photo!]; isotype: PE!)

= C. pinnata Maxim. var. tatarinowii (Maxim.) Kuntze in Verh. Bot. Vereins Prov. Brandenburg 26(2): 182. 1885.

= C. pinnata Maxim. var. ternatifolia W.T. Wang in Acta Phytotax. Sin. 39(4): 331. 2001. Holotype: China. Beijing, Pinggu, Mt. Nanjishan, 13. June. 1972, Pinggu Exped. 224. (PE!)

Distribution. – China. Beijing and adjacent areas of Hebei Province and Tianjin, and central Liaoning Province. Wang and Xie (2007) recorded that C. × pinnata occurred in Heilongjiang Province based on a single collection of E. Licent 9221 (collected in 20 August 1929) deposited in Tianjin Natural History Museum (TIE). However, we carefully checked the specimen in TIE and found that this specimen was collected in Yangjiaping, Zhuolu County, Hebei Province near the western border of Beijing. So, we don’t have evidence that C. × pinnata is distributed in Heilongjiang Province.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The name of the repository and accession number can be found below: National Center for Biotechnology Information (NCBI) BioProject, https://www.ncbi.nlm.nih.gov/bioproject/, PRJNA657443.

Author Contributions

RL, JH, YL, and LlL analyzed the data and prepared the draft. MY, RL, YL, SY, JC, LqL, and LX conducted the field surveys. JC, LP, LqL, and LX proposed the hybrid hypothesis and designed the study. JC and LX wrote and revised the manuscript. All the authors contributed to the article and approved the submitted version.

Funding

This study was supported by the Beijing Natural Science Foundation (grant no. 5182016) and the National Natural Science Foundation of China (grant no. 31670207).

Conflict of Interest

LP was employed by company Beijing Forestry University Forest Science Co. Ltd.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We appreciate the comments of Jenny Xiang from North Carolina State University and Yong Yang from Nanjing Forestry University for the early draft. We thank Xiao-Fei Xia from the Beijing Museum of Natural History for his help with SEM observation and photography; and Su-Hua Yang from the Key Lab of Photobiology CAS, for her help with flow cytometry measurement; and botanist Su Zhao from Liaoning Province for assisting with field collecting.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2021.745988/full#supplementary-material

Supplementary Figure S1 | Natural distribution (blue dash circles) and collecting sites (points) of Clematis pinnata.

Supplementary Figure S2 | Peak value image based on flow cytometry analysis of Clematis pinnata, C. tubulosa, and C. heracleifolia using C. brevicaudata as the external reference.

Supplementary Figure S3 | The backbone of Clematis phylogeny inferred from the concatenated and coalescent based method using 3170 co-orthologous genes with exclusion the two putative taxa, C. pinnata and C. ochotensis. ML bootstrap values (MLBS) of the concatenated analysis were marked on the branches with ^∗ indicating MLBS = 100.

Supplementary Figure S4 | Optimal species networks of the eight populations of Clematis pinnata and its putative parents, as well as of closely related species, inferred using PhyloNet. The results display maximum pseudolikelihood trees with maximum zero to four reticulations allowed. Bar charts show probabilities of different maximum numbers of reticulations. Population location name abbreviations are explained in Table 2.

Supplementary Figure S5 | Optimal species networks of the BHS populations of Clematis pinnata and its putative parents, as well as of closely related species, inferred using PhyloNet with exclusion of another putative hybrid species C. ochotensis. The results display maximum pseudolikelihood trees with maximum zero to four reticulations allowed. Bar charts show probabilities of different maximum numbers of reticulations. Population location name abbreviations are explained in Table 2.

Supplementary Figure S6 | Chloroplast genome maps for Clematis pinnata. Thick lines on the complete outer circle identify the inverted repeat regions (IRa and IRb). The innermost track of the plastome shows the G + C content. Genes on the outside of the map are transcribed in a clockwise direction and genes on the inside of the map are transcribed in a counter-clockwise direction. INV, inversion; TP, transposition; IR, inverted repeats; LSC, large single copy; SSC, small single copy; Pi, nucleotide variability.

Supplementary Figure S7 | Potential suitable distribution areas for Clematis pinnata and its putative parent species as predicted using MaxEnt modeling.

Supplementary Table S1 | Information of transcriptome data of Clematis pinnata and its relatives sampled in this study.

Supplementary Table S2 | Sampling information for the newly generated genome skimming data.

Supplementary Table S3 | Measurements and statistical summaries of 30 selected morphological characters of Clematis pinnata and its putative parents.

Supplementary Table S4 | Samples of Clematis pinnata and its putative parents for pollen morphology scanning using Scanning Electron Microscopy (SEM).

Supplementary Table S5 | Specimens’ information of Clematis brevicaudata, C. heracleifolia, C. tubulosa, C. pinnata used in morphological and niche modeling analysis.

Supplementary Table S6 | The 33 environmental variables used for the niche modeling.

Supplementary Table S7 | Flow cytometry analysis of Clematis tubulosa and C. pinnata using C. brevicaudata as an external standard reference.

Supplementary Table S8 | Flow cytometry analysis of Clematis heracleifolia using C. brevicaudata as an external standard reference.

Supplementary Table S9 | HyDe analysis of Clematis pinnata and its putative parents in each population using the complete SNP data. Sample size of C. pinnata: eight populations, 28 individuals.

Supplementary Table S10 | Percentage contributions and permutation importance of environmental variables included in MaxEnt modeling for Clematis pinnata and its putative parents.

Footnotes

1. ^ http://www.biomarker.com.cn/
2. ^ https://github.com/HeJian151004/get_homology
3. ^ https://cvalues.science.kew.org/
4. ^ https://github.com/Jhe1004/SNPRandomEextraction
5. ^ http://hannonlab.cshl.edu/fastx_toolkit
6. ^ https://github.com/HeJian151004/get_homology
7. ^ http://www.cvh.ac.cn
8. ^ https://www.google.com/earth
9. ^ http://www.worldclim.org
10. ^ https://www.ufz.de/gluv/index.php?en=32435
11. ^ http://www.fao.org/soils-portal/soil-survey/soil-maps-and-databases/harmonized-world-soil-database-v12/zh/
12. ^ http://nelson.wisc.edu/sage/data-and-models/atlas/maps.php
13. ^ https://www.gadm.org/maps.html
14. ^ http://biodiversityinformatics.amnh.org/open_source/maxent/
15. ^ https://www.r-project.org/

References

Abbott, R., Albach, D., Ansell, S., Arntzen, J. W., Baird, S. J. E., Bierne, N., et al. (2013). Hybridization and speciation. J. Evol. Biol. 26, 229–246. doi: 10.1111/j.1420-9101.2012.02599.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Arnold, M. L. (1997). Natural Hybridization and Evolution. NewYork, NY: Oxford University Press.

Google Scholar

Babar, S., Amarnath, G., Reddy, C. S., Jentsch, A., and Sudhakar, S. (2012). Species distribution models: ecological explanation and prediction of an endemic and endangered plant species (Pterocarpus santalinus L. f.). Curr. Sci. 102, 1157–1165. doi: 10.1371/journal.pone.0035867

PubMed Abstract | CrossRef Full Text | Google Scholar

Beckmann, M., Václavík, T., Manceur, A. M., Šprtová, L., Wehrden, H. V., Welk, E., et al. (2014). gIUV: a global UV-B radiation dataset for macroecological studies. Meth. Ecol. Evol. 5, 372–383. doi: 10.1111/2041-210X.12168

CrossRef Full Text | Google Scholar

Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. Roy. Statist. Soc. 57, 289–300. doi: 10.1111/j.2517-6161.1995.tb02031.x

CrossRef Full Text | Google Scholar

Bergsten, J. (2005). A review of long-branch attraction. Cladistics 21, 163–193. doi: 10.1111/j.1096-0031.2005.00059.x

CrossRef Full Text | Google Scholar

Blischak, P. D., Chifman, J., Wolfe, A. D., and Kubatko, L. S. (2018). HyDe: a Python package for genome-scale hybridization detection. Syst. Biol. 67, 821–829. doi: 10.1093/sysbio/syy023

PubMed Abstract | CrossRef Full Text | Google Scholar

Bolger, A. M., Lohse, M., and Usadel, B. (2014). Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30, 2114–2120. doi: 10.1093/bioinformatics/btu170

PubMed Abstract | CrossRef Full Text | Google Scholar

Bourge, M., Brown, S. C., and Siljak-Yakovlev, S. (2018). Flow cytometry as tool in plant sciences, with emphasis on genome size and ploidy level assessment. Genet. Applic. 2, 1–12. doi: 10.31383/ga.vol2iss2pp1-12

CrossRef Full Text | Google Scholar

Broennimann, O., Fitzpatrick, M. C., Pearman, P. B., Petitpierre, B., Pellissier, L., Yoccoz, N. G., et al. (2012). Measuring ecological niche overlap from occurrence and spatial environmental data. Global Ecol. Biogeogr. 21, 481–497. doi: 10.1111/j.1466-8238.2011.00698.x

CrossRef Full Text | Google Scholar

Cao, Z., Liu, X. H., Ogilvie, H. A., Yan, Z., and Nakhleh, L. (2019). Practical aspects of phylogenetic network analysis using PhyloNet. BioRxiv [preprint] 746362. doi: 10.1101/746362

CrossRef Full Text | Google Scholar

Chung, K. S., Oh, K. S., Park, M. S., Nam, B. M., and Chung, G. Y. (2013). Chromosome numbers of 28 taxa in 10 genera of the Ranunculaceae (buttercup family) from the Korean peninsula. Caryologia 66, 128–137. doi: 10.1080/00087114.2013.821831

CrossRef Full Text | Google Scholar

Comai, L. (2005). The advantages and disadvantages of being polyploid. Nat. Rev. Genet. 6, 836–846. doi: 10.1038/nrg1711

PubMed Abstract | CrossRef Full Text | Google Scholar

Corriveau, J. L., and Coleman, A. W. (1988). Rapid screening method to detect potential biparental inheritance of plastid DNA and results for over 200 angiosperm species. Amer. J. Bot. 75, 1443–1458. doi: 10.1002/j.1537-2197.1988.tb11219.x

CrossRef Full Text | Google Scholar

Davenport, K. M., Hiemke, C., McKay, S. D., Thorne, J. W., Lewis, R. M., Taylor, T., and Murdoch, B. M. (2020). Genetic structure and admixture in sheep from terminal breeds in the United States. Anim. Genet. 51, 284–291. doi: 10.1111/age.12905

PubMed Abstract | CrossRef Full Text | Google Scholar

Di Cola, V., Broennimann, O., Petitpierre, B., Breiner, F. T., D’Amen, M., Randin, C., et al. (2017). Ecospat: an R package to support spatial analyses and modeling of species niches and distributions. Ecography 40, 774–787. doi: 10.1111/ecog.02671

CrossRef Full Text | Google Scholar

Doležel, J., Greilhuber, J., and Suda, J. (2007). Estimation of nuclear DNA content in plants using flow cytometry. Nat. Protoc. 2, 2233–2244. doi: 10.1038/nprot.2007.310

PubMed Abstract | CrossRef Full Text | Google Scholar

Dormann, C. F., Elith, J., Bacher, S., Buchmann, C., Carl, G., Carré, G., Marquéz, J.R.G., et al. (2013). Collinearity: a review of methods to deal with it and a simulation study evaluating their performance. Ecography 36, 27–46. doi: 10.1111/j.1600-0587.2012.07348.x

CrossRef Full Text | Google Scholar

Draper, D., Laguna, E., and Marques, I. (2021). Demystifying negative connotations of hybridization for less biased conservation policies. Front. Ecol. Evol. 9:268. doi: 10.1111/10.3389/fevo.2021.637100

CrossRef Full Text | Google Scholar

Fang, M. Y. (1980). “Clematis sect. Viorna (Reichb.) Prantl subsect. Tubulosae Schneid,” in Flora Reipublicae Popularis Sinicae, Vol. 28. ed W.T. Wang (Beijing: Science Press), 93-96

Google Scholar

Feliner, G. N., Álvarez, I., Fuertes-Aguilar, J., Heuertz, M., Marques, I., Moharrek, F., et al. (2017). Is homoploid hybrid speciation that rare? An empiricist’s view. Heredity 118, 513–516. doi: 10.1038/hdy.2017.7

PubMed Abstract | CrossRef Full Text | Google Scholar

Feng, Y., Comes, H. P., Zhou, X. P., and Qiu, Y. X. (2019). Phylogenomics recovers monophyly and early Tertiary diversification of Dipteronia (Sapindaceae). Mol. Phylogenet. Evol. 130, 9–17. doi: 10.1016/j.ympev.2018.09.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Fielding, A. H., and Bell, J. F. (1997). A review of methods for the assessment of prediction errors in conservation presence/absence models. Environm. Conservation 24, 38–49. doi: 10.1017/S0376892997000088

CrossRef Full Text | Google Scholar

Food and Agriculture Organization [FAO] and International Institute for Applied Systems Analysis [IIASA] (2012). Global Agro-ecological Zones (GAEZ v3. 0). Rome: FAO.

Google Scholar

Fu, L., Niu, B., Zhu, Z., Wu, S., and Li, W. (2012). CD-HIT: accelerated for clustering the next-generation sequencing data. Bioinformatics 28, 3150–3152. doi: 10.1093/bioinformatics/bts565

PubMed Abstract | CrossRef Full Text | Google Scholar

Glémin, S., Scornavacca, C., Dainat, J., Burgarella, C., Viader, V., Ardisson, M., et al. (2019). Pervasive hybridizations in the history of wheat relatives. Sci. Adv. 5:eaav9188. doi: 10.1101/300848

CrossRef Full Text | Google Scholar

Gong, W., Long, Y., and Li, M. (1985). Karyotype studies on Clematis from Beijing China. Plant Sci. J. 3, 371–379.

Google Scholar

Goulet, B. E., Roda, F., and Hopkins, R. (2017). Hybridization in plants: old ideas, new techniques. Plant Physiol. 173, 65–78. doi: 10.1104/pp.16.01340

PubMed Abstract | CrossRef Full Text | Google Scholar

Grabherr, M. G., Haas, B. J., Yassour, M., Levin, J. Z., Thompson, D. A., Amit, I., et al. (2011). Full-length transcriptome assembly from RNA-seq data without a reference genome. Nat. Biotechnol. 29, 644–652. doi: 10.1038/nbt.1883

PubMed Abstract | CrossRef Full Text | Google Scholar

Grey-Wilson, C. (2000). Clematis, the Genus. Portland, OR: Timber Press.

Google Scholar

Guisan, A., and Thuiller, W. (2005). Predicting species distribution: offering more than simple habitat models. Ecol. Lett. 8, 993–1009. doi: 10.1111/j.1461-0248.2005.00792.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Haas, B. J., Papanicolaou, A., Yassour, M., Grabherr, M., Blood, P. D., Bowden, J., et al. (2013). De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nat. Protoc. 8, 1494–1512. doi: 10.1038/nprot.2013.084

PubMed Abstract | CrossRef Full Text | Google Scholar

Handel-Mazzetti, H. (1939). Plantae Sinenses a Dre. H. Smith annis 1921–1922, 1924 et 1934 lectae. 33. Ranunculaceae. Acta Horti Gotoburgensis 13, 37–219.

Google Scholar

He, J., Lyu, R. D., Luo, Y. K., Lin, L. L., Yao, M., Xiao, J. M., et al. (2021). An updated phylogenetic and biogeographic analysis based on genome skimming data reveals convergent evolution of shrubby habit in Clematis in the Pliocene and Pleistocene. Mol. Phylogenet. Evol. 164:107259. doi: 10.1016/j.ympev.2021.107259

PubMed Abstract | CrossRef Full Text | Google Scholar

He, J., Yao, M., Lyu, R. D., Lin, L. L., Liu, H. J., Pei, L. Y., et al. (2019). Structural variation of the complete chloroplast genome and plastid phylogenomics of the genus Asteropyrum (Ranunculaceae). Sci. Rep. 9, 1–13. doi: 10.1038/s41598-019-51601-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Hijmans, R. J., Cameron, S. E., Parra, J. L., Jones, P. G., and Jarvis, A. (2005). Very high resolution interpolated climate surfaces for global land areas. Int. J. Climatol. 25, 1965–1978. doi: 10.1002/joc.1276

CrossRef Full Text | Google Scholar

Hoeffding, L. K., Rosengren, A., Thygesen, J. H., Schmock, H., Werge, T., and Hansen, T. (2017). Evaluation of shared genetic susceptibility loci between autoimmune diseases and schizophrenia based on genome-wide association studies. Nord. J. Psych. 71, 20–25. doi: 10.1080/08039488.2016.1198420

PubMed Abstract | CrossRef Full Text | Google Scholar

Huson, D. H., and Bryant, D. (2006). Application of phylogenetic networks in evolutionary studies. Mol. Biol. Evol. 23, 254–267. doi: 10.1093/molbev/msj030

PubMed Abstract | CrossRef Full Text | Google Scholar

Huynh, S., Marcussen, T., Felber, F., and Parisod, C. (2019). Hybridization preceded radiation in diploid wheats. Mol. Phylogenet. Evol. 139:106554. doi: 10.1016/j.ympev.2019.106554

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiang, N., Zhou, Z., Yang, J. B., Zhang, S. D., Guan, K. Y., Tan, Y. H., et al. (2017). Phylogenetic reassessment of tribe Anemoneae (Ranunculaceae): nonmonophyly of Anemone s.l. revealed by plastid datasets. PLoS One 12: e0174792. doi: 10.1371/journal.pone.0174792

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, M. (1997). Släktet Klematis. Södertälje: Magnus Johnson Plantskola.

Google Scholar

Kadereit, J. W. (2015). The geography of hybrid speciation in plants. Taxon 64, 673–687. doi: 10.12705/644.1

CrossRef Full Text | Google Scholar

Katoh, K., Kuma, K., Toh, H., and Miyata, T. (2005). MAFFT version 5: improvement in accuracy of multiple sequence alignment. Nucl. Acids Res. 33, 511–518. doi: 10.1093/nar/gki198

PubMed Abstract | CrossRef Full Text | Google Scholar

Kearse, M., Moir, R., Wilson, A., Stones-Havas, S., Cheung, M., Sturrock, S., et al. (2012). Geneious Basic: an integrated and extendable desktop software platform for the organization and analysis of sequence data. Bioinformatics 28, 1647–1649. doi: 10.1093/bioinformatics/bts199

PubMed Abstract | CrossRef Full Text | Google Scholar

Kimura, M. (1980). A simple method for estimating evolutionary rate of base substitutions through comparative studies of nucleotide sequences. J. Mol. Evol. 16, 111–120.

Google Scholar

Kubatko, L. S., and Chifman, J. (2019). An invariants-based method for efficient identification of hybrid species from large-scale genomic data. BMC Evol. Biol. 19:112. doi: 10.1186/s12862-019-1439-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuntze, O. (1885). Monographie der Gattung Clematis. Verhandlungen des Botanischen Vereins der Prov. Brandenburg (Abhandlungen) 26, 83–202.

Google Scholar

Lechner, M., Findeiß, S., Steiner, L., Marz, M., Stadler, P. F., and Prohaska, S. J. (2011). Proteinortho: detection of (co-)orthologs in large-scale analysis. BMC Bioinf. 12:124. doi: 10.1186/1471-2105-12-124

PubMed Abstract | CrossRef Full Text | Google Scholar

Lehtonen, S., Christenhusz, M. J. M., and Falck, D. (2016). Sensitive phylogenetics of Clematis and its position in Ranunculaceae. Bot. J. Linn. Soc. 182, 825–867. doi: 10.1111/boj.12477

CrossRef Full Text | Google Scholar

Li, H. T., Yi, T. S., Gao, L. M., Ma, P. F., Zhang, T., Yang, J. B., et al. (2019). Origin of angiosperms and the puzzle of the Jurassic gap. Nat. Plants 5, 461–470. doi: 10.1038/s41477-019-0421-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, M., Chen, S., Zhou, R., Fan, Q., Li, F., and Liao, W. (2017). Molecular evidence for natural hybridization between Cotoneaster dielsianus and C. glaucophyllus. Front. Plant Sci. 8:704. doi: 10.3389/fpls.2017.00704

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, W. Z., and Godzik, A. (2006). CD-HIT: a fast program for clustering and comparing large sets of protein or nucleotide sequences. Bioinformatics 22, 1658–1659. doi: 10.1007/978-1-4899-7478-5_221

CrossRef Full Text | Google Scholar

Liao, R., Sun, W., and Ma, Y. (2021). Natural hybridization between two butterfly bushes in Tibet: dominance of F 1 hybrids promotes strong reproductive isolation. BMC Plant Biol. 21:133. doi: 10.1186/s12870-021-02909-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Lin, L. L., He, J., Xie, L., and Cui, G. F. (2020). Prediction of the suitable area of the chinese white pines (Pinus subsect. Strobus) under climate changes and implications for their conservation. Forests 11:996. doi: 10.3390/f11090996

CrossRef Full Text | Google Scholar

Liu, B., Abbott, R. J., Lu, Z., Tian, B., and Liu, J. (2014). Diploid hybrid origin of Ostryopsis intermedia (Betulaceae) in the Qinghai-Tibet Plateau triggered by Quaternary climate change. Mol. Ecol. 23, 3013–3027. doi: 10.1111/mec.12783

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, B. B., Campbell, C. S., Hong, D. Y., and Wen, J. (2020). Phylogenetic relationships and chloroplast capture in the Amelanchier-Malacomeles-Peraphyllum clade (Maleae, Rosaceae): evidence from chloroplast genome and nuclear ribosomal DNA data using genome skimming. Mol. Phylogenet. Evol. 147:106784. doi: 10.1016/j.ympev.2020.106784

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, H. J., He, J., Ding, C., Lyu, R. D., Pei, L. Y., Cheng, J., et al. (2018). Comparative analysis of complete chloroplast genomes of Anemoclema, Anemone, Pulsatilla, and Hepatica revealing structural variations among genera in tribe Anemoneae (Ranunculaceae). Front. Plant Sci. 9:1097. doi: 10.3389/fpls.2018.01097

PubMed Abstract | CrossRef Full Text | Google Scholar

Lohse, M., Drechsel, O., Kahlau, S., and Bock, R. (2013). OrganellarGenomeDRAW–a suite of tools for generating physical maps of plastid and mitochondrial genomes and visualizing expression data sets. Nucl. Acids Res. 41, W575–W581. doi: 10.1093/nar/gkt289

PubMed Abstract | CrossRef Full Text | Google Scholar

Maddison, W. P. (1997). Gene trees in species trees. Syst. Biol. 46, 523–536. doi: 10.1038/npg.els.0005137

CrossRef Full Text | Google Scholar

Mai, U., and Mirarab, S. (2018). TreeShrink: fast and accurate detection of outlier long branches in collections of phylogenetic trees. BMC Genom. 19:272. doi: 10.1186/s12864-018-4620-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Majeed, A., Singh, A., Choudhary, S., and Bhardwaj, P. (2019). RNAseq-based phylogenetic reconstruction of Taxaceae and Cephalotaxaceae. Cladistics 35, 461–468. doi: 10.1111/cla.12362

CrossRef Full Text | Google Scholar

Makino, T. (1907). Observations on the flora of Japan. Bot. Mag. (Tokyo) 21, 86–88.

Google Scholar

Mallet, J. (2007). Hybrid speciation. Nature 446, 279–283. doi: 10.1038/nature05706

PubMed Abstract | CrossRef Full Text | Google Scholar

Malomane, D. K., Reimer, C., Weigend, S., Weigend, A., Sharifi, A. R., and Simianer, H. (2018). Efficiency of different strategies to mitigate ascertainment bias when using SNP panels in diversity studies. BMC Genom. 19:22. doi: 10.1186/s12864-017-4416-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Maximowicz, C. J. (1877). Diagnoses plantarum olcano Japoniae et Mandshuriae. Bull. l’Academie impériale Sci. St. Pétersburg 22, 209–264.

Google Scholar

McLean, B. S., Bell, K. C., Allen, J. M., Helgen, K. M., and Cook, J. A. (2019). Impacts of inference method and data set filtering on phylogenomic resolution in a rapid radiation of ground squirrels (Xerinae: marmotini). Syst. Biol. 68, 298–316. doi: 10.1093/sysbio/syy064

PubMed Abstract | CrossRef Full Text | Google Scholar

Meng, C., and Kubatko, L. S. (2009). Detecting hybrid speciation in the presence of incomplete lineage sorting using gene tree incongruence: a model. Theor. Pop. Biol. 75, 5–45. doi: 10.1016/j.tpb.2008.10.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Miikeda, O., Kita, K., Handa, T., and Yukawa, T. (2006). Phylogenetic relationships of Clematis (Ranunculaceae) based on chloroplast and nuclear DNA sequences. Bot. J. Linn. Soc. 152, 153–168. doi: 10.1111/j.1095-8339.2006.00551.x

CrossRef Full Text | Google Scholar

Mirarab, S., Reaz, R., Bayzid, M.D. S., Zimmermann, T., Swenson, M. S., and Warnow, T. (2014). ASTRAL: genome-scale coalescent-based species tree estimation. Bioinformatics 30, i541–i548. doi: 10.1093/bioinformatics/btu462

PubMed Abstract | CrossRef Full Text | Google Scholar

Moore, M. J., Soltis, P. S., Bell, C. D., Burleigh, J. G., and Soltis, D. E. (2010). Phylogenetic analysis of 83 plastid genes further resolves the early diversifiation of eudicots. Proc. Natl. Acad. Sci. U.S.A. 107, 4623–4628. doi: 10.1073/pnas.0907801107

PubMed Abstract | CrossRef Full Text | Google Scholar

Morales-Briones, D. F., Kadereit, G., Tefarikis, D. T., Moore, M. J., Smith, S. A., Brockington, S. F., et al. (2021). Disentangling sources of gene tree discordance in phylogenomic data sets: testing ancient hybridizations in Amaranthaceae s.l. Syst. Biol. 70, 219–235. doi: 10.1093/sysbio/syaa066

PubMed Abstract | CrossRef Full Text | Google Scholar

New, M., Hulme, M., and Jones, P. (1999). Representing twentieth-century space-time climate variability. Part I: development of a 1961–1990 mean monthly terrestrial climatology. J. Clim. 12, 829–856. doi: 10.1175/1520-0442(1999)012<0829:RTCSTC>2.0.CO;2

CrossRef Full Text | Google Scholar

Nielsen, R., and Signorovitch, J. (2003). Correcting for ascertainment biases when analyzing SNP data: applications to the estimation of linkage disequilibrium. Theor. Popul. Biol. 63, 245–255. doi: 10.1016/S0040-5809(03)00005-4

CrossRef Full Text | Google Scholar

Ning, H., Pan, Y. Z., and Gong, X. (2019). Molecular evidence for natural hybridi zation between Ligularia nelumbifolia and Cremanthodium stenoglossum (Asteraceae, Senecioneae). Botany 97, 53–69. doi: 10.1139/cjb-2018-0022

PubMed Abstract | CrossRef Full Text | Google Scholar

Nishikawa, H., Iijima, T., Kajitani, R., Yamaguchi, J., Ando, T., Suzuki, Y., et al. (2015). A genetic mechanism for female-limited Batesian mimicry in Papiliobutterfly. Nat. Genet. 47, 405–409. doi: 10.1038/ng.3241

PubMed Abstract | CrossRef Full Text | Google Scholar

Patterson, N., Moorjani, P., Luo, Y., Mallick, S., Rohland, N., Zhan, Y., et al. (2012). Ancient admixture in human history. Genetics 192, 1065–1093. doi: 10.1534/genetics.112.145037

PubMed Abstract | CrossRef Full Text | Google Scholar

Qu, X., Yan, M., Zou, J., Jiang, M., Yang, K., and Le, J. (2018). A2-type cyclin is required for the asymmetric entry division in rice stomatal development. J. Exp. Bot. 69, 3587–3599. doi: 10.1093/jxb/ery158

PubMed Abstract | CrossRef Full Text | Google Scholar

Qu, X. J., Moore, M. J., Li, D. Z., and Yi, T. S. (2019). PGA: a software package for rapid, accurate, and flexible batch annotation of plastomes. Plant Meth. 15:50. doi: 10.1186/s13007-019-0435-7

PubMed Abstract | CrossRef Full Text | Google Scholar

R Core Team (2018). R: A Language and Environment for Statistical Computing. Vienna: R Foundation for Statistical Computing.

Google Scholar

Rehder, A., and Wilson, E. H. (1913). “Clematis,” in Plantae Wilsonianae, Vol. 1, ed <gnm>E. H. </gnm><snm>Wilson</snm> (Cambridge: Cambridge University Press), 319-343.

Google Scholar

Roberts, W. R., and Roalson, E. H. (2018). Phylogenomic analyses reveal extensive gene flow within the magic flowers (Achimenes). Amer. J. Bot. 105, 726–740. doi: 10.1002/ajb2.1058

PubMed Abstract | CrossRef Full Text | Google Scholar

Ronquist, F., and Huelsenbeck, J. P. (2003). MrBayes 3: bayesian phylogenetic inference under mixed models. Bioinformatics 19, 1572–1574. doi: 10.1093/bioinformatics/btg18

CrossRef Full Text | Google Scholar

Schoener, T. W. (1968). The Anolis lizards of Bimini: resource partitioning in a complex fauna. Ecology 49, 704–726. doi: 10.2307/1935534

CrossRef Full Text | Google Scholar

Shi, J. H. (2003). A Preliminary Study of the Hybrid Origin of Clematis pinnata Maxim. (Ranunculaceae). master’s thesis. Beijing: Institute of Botany, the Chinese Academy of Sciences.

Google Scholar

Shitara, T., Nakamura, Y., Matsui, T., Tsuyama, I., Ohashi, H., and Kamijo, T. (2018). Formation of disjunct plant distributions in Northeast Asia: a case study of Betula davurica using a species distribution model. Plant Ecol. 219, 1105–1115. doi: 10.1007/s11258-018-0862-y

CrossRef Full Text | Google Scholar

Smith, S. A., and Dunn, C. W. (2008). Phyutility: a phyloinformatics tool for trees, alignments and molecular data. Bioinformatics 24, 715–716. doi: 10.1093/bioinformatics/btm619

PubMed Abstract | CrossRef Full Text | Google Scholar

Sochor, M., Vašut, R. J., Sharbel, T. F., and Trávnícek, B. (2015). How just a few makes a lot: speciation via reticulation and apomixis on example of European brambles (Rubus subgen. Rubus, Rosaceae). Mol. Phylogenet Evol. 89, 13–27. doi: 10.1016/j.ympev.2015.04.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Soltis, P. S., and Soltis, D. E. (2009). The role of hybridization in plant speciation. Ann. Rev. Plant Biol. 60, 561–588. doi: 10.1146/anurev.arplant.043008.092039

CrossRef Full Text | Google Scholar

Stamatakis, A. (2014). RaxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics 30, 1312–1313. doi: 10.1093/bioinformatics/btu033

PubMed Abstract | CrossRef Full Text | Google Scholar

Sveinsson, S., McDill, J., Wong, G. K. S., Li, J., Li, X., Deyholos, M. K., et al. (2014). Phylogenetic pinpointing of a paleopolyploidy event within the flax genus (Linum) using transcriptomics. Ann. Bot. 113, 753–761. doi: 10.1093/aob/mct306

PubMed Abstract | CrossRef Full Text | Google Scholar

Tamura, M. (1987). A classification of genus Clematis. Acta Phytotax. Geobot. 38, 33–44.

Google Scholar

Tamura, M. (1995). “Clematis,” in Die Natürlichen Pflanzenfamilien, Vol. 17a(4). eds A. Engler, and K. Prantl (Leipzig: Engelmann), 372–373.

Google Scholar

Tang, C. Q., Dong, Y. F., Herrando-Moraira, S., Matsui, T., Ohashi, H., He, L. Y., et al. (2017). Potential effects of climate change on geographic distribution of the Tertiary relict tree species Davidia involucrata in China. Sci. Rep. 7, 1–18. doi: 10.1038/srep43822

PubMed Abstract | CrossRef Full Text | Google Scholar

Ting, C. T. (1980). “Clematis sect. Clematis subsect. Vitalbae Prantl,” in Flora Reipublicae Popularis Sinicae, Vol. 28. ed W. T. Wang (Beijing: Science Press), 178–199

Google Scholar

Toomey, M., and Leeds, E. (2001). An Illustrated Encyclopedia of Clematis. Portland, OR: Timber Press.

Google Scholar

Wang, N., Wang, K. L., Liu, Q. H., and Liu, Q. C. (2017). Karyotype analysis of seven wild Clematis species. Acta Pratacult. Sin. 26, 123–130. doi: 10.11686/cyxb2017025

CrossRef Full Text | Google Scholar

Wang, N., Yang, Y., Moore, M. J., Brockington, S. F., Walker, J. F., Brown, J. W., et al. (2019). Evolution of Portulacineae marked by gene tree conflict and gene family expansion associated with adaptation to harsh environments. Mol. Biol. Evol. 36, 112–126. doi: 10.1093/molbev/msy200

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, W. T. (2001). Notes on the genus Clematis (Ranunculaceae) (IV). Acta Phytotax. Sin. 39, 1–19. doi: 10.1360/aps040087

CrossRef Full Text | Google Scholar

Wang, W. T., and Li, L. Q. (2005). A new system of classification of the genus Clematis (Ranunculaceae). Acta Phytotax. Sin. 43, 431–488. doi: 10.1360/aps040091

CrossRef Full Text | Google Scholar

Wang, W. T., and Xie, L. (2007). A revision of Clematis sect. Tubulosae (Ranunculaceae). Acta Phytotax. Sin. 45, 425–457. doi: 10.1360/aps06114

CrossRef Full Text | Google Scholar

Warren, D. L., Glor, R. E., and Turelli, M. (2008). Environmental niche equivalency versus conservatism: quantitative approaches to niche evolution. Evolution 62, 2868–2883. doi: 10.1111/j.1558-5646.2008.00482.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Warren, D. L., Glor, R. E., and Turelli, M. (2010). ENMTools: a toolbox for comparative studies of environmental niche models. Ecography 33, 607–611. doi: 10.1111/j.1600-0587.2009.06142.x

CrossRef Full Text | Google Scholar

Wen, D., Yu, Y., Zhu, J., and Nakhleh, L. (2018). Inferring phylogenetic networks using PhyloNet. Syst. Biol. 67, 735–740. doi: 10.1093/sysbio/syy015

PubMed Abstract | CrossRef Full Text | Google Scholar

White, O. W., Reyes-Betancort, A., Chapman, M. A., and Carine, M. A. (2018). Independent homoploid hybrid speciation events in the Macaronesian endemic genus Argyranthemum. Molec. Ecol. 27, 4856–4874. doi: 10.1111/mec.14889

PubMed Abstract | CrossRef Full Text | Google Scholar

Wood, T. E., Takebayashi, N., Barker, M. S., Mayrose, I., Greenspoon, P. B., and Rieseberg, L. H. (2009). The frequency of polyploid speciation in vascular plants. Proc. Natl. Acad. Sci. U.S.A. 106, 13875–13879. doi: 10.1073/pnas.0811575106

PubMed Abstract | CrossRef Full Text | Google Scholar

Xie, L., and Li, L. Q. (2012). Variation of pollen morphology, and its implications in the phylogeny of Clematis (Ranunculaceae). Plant Syst. Evol. 298, 1437–1453. doi: 10.1007/s00606-012-0648-y

CrossRef Full Text | Google Scholar

Xie, L., Shi, J. H., and Li, L. Q. (2005). Identity of Clematis tatarinowii and C. pinnata var. ternatifolia (Ranunculaceae). Ann. Bot. Fenn. 42, 305–308.

Google Scholar

Xie, L., Wen, J., and Li, L. Q. (2011). Phylogenetic analyses of Clematis (Ranunculaceae) based on sequences of nuclear ribosomal ITS and three plastid regions. Syst. Bot. 36, 907–921. doi: 10.1600/036364411X604921

CrossRef Full Text | Google Scholar

Xu, S. (2000). Phylogenetic analysis under reticulate evolution. Mol. Biol. Evol. 17, 897–907. doi: 10.1093/oxfordjournals.molbev.a026370

PubMed Abstract | CrossRef Full Text | Google Scholar

Yakimowski, S. B., and Rieseberg, L. H. (2014). The role of homoploid hybridization in evolution: a century of studies synthesizing genetics and ecology. Amer. J. Bot. 101, 1247–1258. doi: 10.3732/ajb.1400201

PubMed Abstract | CrossRef Full Text | Google Scholar

Yan, S. X., Liu, H. J., Lin, L. L., Liao, S., Li, J. Y., Pei, L. Y., et al. (2016). Taxonomic status of Clematis acerifolia var. elobata, based on molecular evidence. Phytotaxa 268, 209–219. doi: 10.11646/phytotaxa.268.3.5

CrossRef Full Text | Google Scholar

Yang, Y., Sun, P., Lv, L., Wang, D., Ru, D., Li, Y., et al. (2020). Prickly waterlily and rigid hornwort genomes shed light on early angiosperm evolution. Nat. Plants 6, 215–222. doi: 10.1038/s41477-020-0594-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, Y., and Nakhleh, L. (2015). A maximum pseudo-likelihood approach for phylogenetic networks. BMC Genom. 16:S10. doi: 10.1186/1471-2164-16-S10-S10

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeng, L., Zhang, N., Zhang, Q., Endress, P. K., Huang, J., and Ma, H. (2017). Resolution of deep eudicot phylogeny and their temporal diversification using nuclear genes from transcriptomic and genomic datasets. New Phytol. 214, 1338–1354. doi: 10.1111/nph.14503

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeng, L., Zhang, Q., Sun, R., Kong, H., Zhang, N., and Ma, H. (2014). Resolution of deep angiosperm phylogeny using conserved nuclear genes and estimates of early divergence times. Nat. Commun. 5:4956. doi: 10.1038/ncomms5956

PubMed Abstract | CrossRef Full Text | Google Scholar

Zha, H. G., Milne, R. I., and Sun, H. (2010). Asymmetric hybridization in Rhododendron agastum: a hybrid taxon comprising mainly F1s in Yunnan, China. Ann.Bot. 105, 89–100. doi: 10.1093/aob/mcp267

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, B. W., Xu, L. L., Li, N., Yan, P. C., Jiang, X. H., Woeste, K. E., et al. (2019). Phylogenomics reveals an ancient hybrid origin of the Persian walnut. Mol. Biol. Evol. 36, 2451–2461. doi: 10.1093/molbev/msz112

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, C., Li, S. Q., Zhang, Y., Zhu, Z. M., Liu, J. Q., and Gao, X. F. (2020). Molecular and morphological evidence for hybrid origin and matroclinal inheritance of an endangered wild rose, Rosa× pseudobanksiae (Rosaceae) from China. Conserv. Genet. 21, 1–11. doi: 10.1007/s10592-019-01227-8

CrossRef Full Text | Google Scholar

Zhang, Y., Hong, Y., Ren, C., Tang, M., Hoot, S. B., and Yang, Q. E. (2015). Palynology, cytology, and molecular systematics of Anemone section Begoniifolia (Ranunculaceae). Plant Syst. Evol. 301, 411–424. doi: 10.1007/s00606-014-1082-0

CrossRef Full Text | Google Scholar

Zhang, Y. L., and He, S. Y. (1990). Chromosome studies on 6 species of Clematis in China. Plant Sci. J. 8, 114–121.

Google Scholar

Zhao, L., Zhang, N., Ma, P. F., Liu, Q., Li, D. Z., and Guo, Z. H. (2013). Phylogenomic analyses of nuclear genes reveal the evolutionary relationships within the BEP clade and the evidence of positive selection in Poaceae. PLoS One 8:e64642. doi: 10.1371/journal.pone.0064642.

PubMed Abstract | CrossRef Full Text | Google Scholar