Impact Factor 4.151
2017 JCR, Clarivate Analytics 2018

Frontiers journals are at the top of citation and impact metrics

Original Research ARTICLE

Front. Genet., 17 May 2018 |

Phylogeography of Parasyncalathium souliei (Asteraceae) and Its Potential Application in Delimiting Phylogeoregions in the Qinghai-Tibet Plateau (QTP)-Hengduan Mountains (HDM) Hotspot

Nan Lin1,2,3, Tao Deng3, Michael J. Moore4, Yanxia Sun1, Xianhan Huang3, Wenguang Sun3, Dong Luo3, Hengchang Wang1*, Jianwen Zhang3* and Hang Sun3*
  • 1Key Laboratory of Plant Germplasm Enhancement and Specialty Agriculture, Wuhan Botanical Garden, Chinese Academy of Sciences, Wuhan, China
  • 2University of Chinese Academy of Sciences, Beijing, China
  • 3Key Laboratory for Plant Diversity and Biogeography of East Asia, Kunming Institute of Botany, Chinese Academy of Sciences, Kunming, China
  • 4Department of Biology, Oberlin College, Oberlin, OH, United States

Biogeographic regionalization can help to better understand diversity in biogeography, conservation, and macroecology. Historical regionalization schemes typically focus on species distributions, often rarely considering the rich context that phylogeographic information can provide. We investigated whether phylogeographic data could help to delineate floristic regions in the Qinghai-Tibet Plateau (QTP)-Hengduan Mountains (HDM) region by analyzing phylogeographic structure in the herb Parasyncalathium souliei (Asteraceae). We sequenced the plastid psbA-trnH and trnL-rpl32 spacer regions for 417 individuals in 36 populations across the geographic range of the species. To estimate the phylogeographic history of this species, a series of population genetic, phylogenetic, molecular dating, and haplotype network analyses were conducted, as were tested for historical demographic expansions. Using occurrence data, species distribution modeling was used to estimate geographic distributions at three time points: the present, the Mid-Holocene and the Last Glacial Maximum. Significant phylogeographic structure was evident (NST> GST; P < 0.05) among the 37 haplotypes detected. Four major haplogroups were identified based on phylogenetic analyses. Private haplotypes were restricted to geographically distinct regions that generally corresponded to previously identified biogeographic subregions within the QTP-HDM region. Our results imply Pliocene-Pleistocene diversification of P. souliei and suggest that the species may have been geographically widespread early in its history. This study may provide valuable evidence for phylogeographic regionalization using chloroplast genetic data in a common, widespread endemic species from the QTP-HDM.


The Himalayas, especially their core regions, i.e., the Qinghai-Tibet Plateau (QTP) and the Hengduan Mountains (HDM), comprise one of the key high-altitude biodiversity hotspots in the world (Myers et al., 2000). Geological and climatic differences have created profound ecological heterogeneity in the QTP and HDM: eastern regions have deep valleys characterized mainly by a warm and wet climate, whereas the high-elevation central and western Himalayas are characterized by a cold and dry climate (Wu and Wu, 1996; Favre et al., 2015). In an influential study, Wu and Wu (1996) divided these regions into two floristic subkingdoms, Sino-Himalayan and the Tibetan Plateau, based on the distribution of taxa endemic to these areas (Figure 1A). The former subkingdom includes most parts of the Yunnan Plateau, HDM, and the East Himalayan region, while the latter subkingdom is composed of most regions of the QTP including the Tangut the Tibet–Pamir–Kunlun, and the Western Himalayan region. Wu and Wu (1996) further divided the center of the HDM region itself into three subregions: the northern, southern, and Three river-gorges regions (Wu and Wu, 1996; Figure 1A). In the past two decades phylogeographic studies have been conducted in the QTP-HDM region (see introduction of Ren et al., 2017). Almost all of these studies have demonstrated that past climatic changes, including at least some related to geological events, have affected the demographic history of the organisms and have further suggested that multiple plant refugia probably existed in the Himalayas (Sun et al., 2017; Chen et al., 2018). Nevertheless, the results of these studies have rarely been incorporated into analyses of phylogeographic regionalization, in contrast to other areas such as South America (Amarilla et al., 2015).


Figure 1. Biogeographic regions proposed for the flora of the QTP-HDM: (A) Previous regions proposed and modified from Wu and Wu (1996); (B) Simplified phylogeoregions proposed in this study based on private haplotypes: phylogeoregions I, II, and III represent the north and south Hengduan Mountains phylogeoregions, and the Eastern Himalayan phylogeoregion, respectively. The dots in (B) correspond to the 36 natural populations of Parasyncalathium souliei in this study.

Parasyncalathium souliei J. W. Zhang, Boufford and H. Sun (Asteraceae: Cichorieae) (Zhang et al., 2011a) is a widespread species of the QTP-HDM that represents an excellent opportunity to explore the phylogeography of this region. As the only species in the recently described genus Parasyncalathium, P. souliei is known from five of the QTP-HDM subregions of Wu and Wu (1996): the North Hengduan Mountains (hereafter called NHDM), the South Hengduan Mountains (hereafter called SHDM), Three river-gorges, Southeast Xizang and Southern wing of the East Himalayan (Figure 1A; Wu and Wu, 1996). Previous biogeographic studies have suggested that this species may have originated relatively recently (Zhang et al., 2011b; Kilian et al., 2017). P. souliei is small, biennial rosette-forming herb up to 19 cm tall that occurs in well-delimited patches within open bare ground (Zhang et al., 2011b). Heads are congested and are surrounded by a ring of basal leaves, and fruits are likely wind-dispersed (Zhang et al., 2011b). Individuals are relatively small and have a short life span, with relatively high reproductive and mutation rates (Robert et al., 2013). P. souliei is a pioneer species in the QTP-HDM, and thus particularly sensitive to environmental/habitat change (Zhang et al., 2011b).

In this study, we employed an integrative approach to identify phylogeographic patterns within P. souliei and to determine which factors are responsible for its current geographic distribution. We also propose the concept “phylogeoregion” to delimit regions with phylogeographically concordant patterns. More specifically, our goals were to explore: (1) how topography/geography and Pleistocene climatic oscillations have shaped the present-day distributions of P. souliei; (2) whether significant divergence and genetic structure exists among P. souliei populations across different floristic regions/subregions; and (3) whether phylogeographic information may be useful in characterizing potential phylogeoregions in the QTP-HDM region.

Materials and Methods

DNA Samples, PCR, and Sequencing

All the Material was collected and deposited in the Herbarium of the Kunming Institute of Botany, China. Otherwise, no specific permissions were required for material collection of the nature reserves about sampling and investigate in field of Provincial Forestry Department; the locations of Parasyncalathium souliei are not privately-owned and are endangered or protected in the field. A total of 417 individuals representing 36 populations of P. souliei were collected from across the geographic range of the species. Within each population, 3–19 individuals growing at least 40 m apart were collected to maximize sampling of genetic diversity within populations. For most populations, all individuals found were sampled to obtain maximum coverage of P. souliei. Fresh leaves of all individuals were dried in the field using silica gel. Total genomic DNA was extracted using a QIAGEN DNeasy Plant Mini Kit (QIAGEN, Beijing, China). Two plastid DNA spacer regions (psbA-trnH and trnL-rpl32) were amplified via PCR in 80 μl reactions, each containing 30–40 ng of template DNAs, 50 mM Tris-HCl, 1.5 mM MgCl2, 0.5 mM dNTPs, 2 mM of each primer and 0.75 U of Taq polymerase (TaKaRa, Dalian, China). The following PCR program was used: one cycle at 94°C for 4 min; 32 cycles at 94°C for 1 min, 55°C for 1 min, 72°C for 1 min; and one cycle at 72°C for 7 min. PCR products were purified with a QIAquick PCR Purification Kit (BioTeke, Beijing, China), and were sequenced using an ABI 3730XL automated DNA sequencer (Applied Bio-systems, Foster City, California, U.S.A.). The sequences were deposited in GenBank (accession numbers MH023542-MH023958).

Phylogenetic Analyses and Haplotype Network

Sequences were aligned using MAFFT 6.864 with default parameters (Katoh et al., 2002), and were then manually adjusted using BioEdit 7.0.9 (Hall, 1999). All the sequences were assigned to different haplotypes using DnaSP 5.0 (Librado and Rozas, 2009), and these were then employed in downstream analyses. Indels were coded as presence/absence characters except for gaps caused by mononucleotide repeats, which were excluded from analyses. To evaluate genetic relationships among haplotypes, a haplotype network was constructed using NETWORK 4.6.11 with default parameters (Bandelt et al., 1999). A phylogeny of plastid haplotypes was estimated using MrBayes 3.1.2 (Yang and Rannala, 1997). Thirteen closely related taxa were selected for outgroups based on Kilian et al. (2017; Appendix 1); please note that P. souliei is listed under the older synonym Melanoseris souliei). To identify the optimal nucleotide substitution model, MrModeltest 2.3 was used (Nylander, 2004) with default parameters under the Akaike Information Criterion (AIC). For MrBayes, a random starting tree was used, and one cold and three heated chains were run simultaneously in two independent runs for 10,000,000 generations, with trees sampled every 2,000 generations. After discarding the first 25% of trees as burn-in, chain convergence was examined using TRACER 1.5.1 ( In addition, to estimate genetic differentiation among P. souliei haplotypes during Pliocene to Pleistocene, molecular dating was performed under a Bayesian approach as implemented in BEAST 1.6 (Drummond and Rambaut, 2007). Analyses were run using the uncorrelated log-normal (UCLN) relaxed-clock model, and models of sequence evolution were the same as for the MrBayes analyses (GTR+I+G). The root node was constrained using a normal prior distribution [Mean: 6.9 Ma, Stdev: 0.5] by using secondary calibration based on the results of Kilian et al. (2017).

Genetic and Spatial Analyses

The distribution of each plastid haplotype was plotted on a digital elevation map (DEM) derived from Shuttle Radar Topography Mission (SRTM) data ( Population diversity parameters, including number of haplotypes, haplotype diversity and nucleotide diversity were calculated from the haplotype network and phylogenetic analyses using DnaSP 5.0 (Librado and Rozas, 2009). Average within-population diversity (HS) and total haplotype diversity (HT) were calculated in PERMUT ( to evaluate the level of genetic variation; GST and NST were used to estimate differentiation between populations. When NST is compared to GST using U-statistics, a significantly higher NST than GST usually is interpreted as indicative of phylogeographic structure. The differences between GST and NST were assessed using a permutation test with 1,000 random permutations. Samova 1.0 with default parameters (Dupanloup et al., 2002) was used to define the number of groups of populations (K) that are geographically homogenous and maximally differentiated from each other. To evaluate population genetic structure, an analysis of molecular variance (AMOVA) was used to examine the genetic variation within and between populations using Arlequin 3.1 (Excoffier et al., 2005). To investigate the correlation between genetic and geographic distance, a Mantel test was performed with 1,000 permutations in GenAlEx 6.3 (Peakall and Smouse, 2006).

Demographic History and Species Distribution Modeling

To evaluate whether populations experienced recent range expansion(s), Arlequin 3.1 (Excoffier et al., 2005) was used to perform pairwise mismatch distribution and neutrality analyses for the major population groups identified. Tajima's D and Fu's Fs (Tajima, 1989; Fu, 1997) were calculated to detect past demographic history using DnaSP 5.0 (Librado and Rozas, 2009), with 1000 simulated samples and a coalescent algorithm. Populations that have experienced expansion are expected to have a unimodal shape in the mismatch distribution. The sum of squared deviations (SSD) between observed and expected mismatches was also computed. Harpending's raggedness index (r) and its P-values were estimated to test for the significance of the population expansion model (Harpending, 1994).

Through fieldwork and using specimen records from three herbaria (PE, KUN, and HNWP; abbreviations follow Index Herbariorum:, a total of 159 species occurrence points were compiled for species distribution modeling (Figure S1, and these data were submitted to PANGAEA, accession number 10.1594/PANGAEA.887276). Nineteen bioclimatic variables were obtained from WorldClim 1.4, with a resolution of 1 km2 (Collins et al., 2004). The following three ecologically relevant bioclimatic variables were retained after eliminating highly correlated variables (i.e., r ≥ 0.75) (Hijmans et al., 2005): precipitation of the wettest quarter, mean temperature of the warmest quarter and precipitation seasonality (coefficient of variation). The following three layers from the Model for Interdisciplinary Research on Climate (MIROC) (Hasumi and Emori, 2004) were employed for distribution modeling: the current layer, the Mid-Holocene layer (c. 6 ka), and the Last Glacial Maximum (c. 22 ka) layer. All models were run in MaxEnt 3.2 with 10 replicates (Phillips et al., 2006), and we partitioned the locality data into training and testing data sets (75 and 25%, respectively) to evaluate the quality of the model.


Sequence Characteristics, Phylogenetic Relationships, and Haplotype Distribution

The length of the concatenated two-locus alignment was 1507 bp, including 77 parsimony-informative sites within P. souliei. A total of 86 substitutions and eight indels (not including poly-A and poly-T repeats) were detected (Appendix 2). MrModeltest identified GTR+I+G as the optimal nucleotide substitution model. A total of 37 haplotypes (H1–H37) were identified across all 417 individuals (accession number MH023959-MH023995), with 15 populations fixed for a single haplotype (Figure 2, Table 1). The most frequent and widespread haplotype (H1) was found in 141 individuals (33.9%) and 15 populations (41.7%), and the second most frequent haplotype (H3) was found in 91 individuals (21.9%) and 14 populations (38.9%). The MrBayes phylogeny of haplotypes revealed that all P. souliei haplotypes formed a well-supported clade (PP = 1.00). This is congruent with the result of the network and the BEAST-based chronogram. Four haplogroups were identified based on the phylogenetic and network analyses (Figures 3, 4). In haplogroup 1, the most frequent haplotype H1 was shared by approximately half of all populations (41.7%), and haplotype H3 was widely distributed, but most common in more northern regions (Figure 2). Two private haplotypes (H5 and H22) were unique to populations CN and LW, respectively. Haplotypes in haplogroup 2 (Figure 2) were restricted to six populations (GE, ZD, MX, JC, ZM, and YJ) in the NHDM, and five of them possessed more than two haplotypes (Figure 2). Haplotypes in haplogroup 3 were restricted to the SHDM (Figure 2), with all populations in this region containing three or more haplotypes. Haplotypes in haplogroup 4 were restricted to five populations (ZX, ZL, TT, DD, and SJ) in the QTP (Figure 2). These same four haplogroups were recovered as major clades in the haplotype phylogeny (Figure 4). The molecular dating analysis based on this tree recovered a Pliocene to Pleistocene diversification of these four P. souliei haplogroups (Figure 4).


Figure 2. The distribution of 37 cpDNA haplotypes and the four population groups identified within 36 populations of Parasyncalathium souliei. Four population groups were defined according to the haplogroups. The few populations (i.e., KW, AL, GE, MX, DX, HS, DD, TT, and ZX) that contained haplotypes from multiple haplogroups were assigned to population groups based on the proportion of haplotypes from each group.


Table 1. Details of sampling locations, number of individuals, geographic coordinates, elevation and haplotypes of Parasyncalathium souliei.


Figure 3. Network of relationships among the 37 cpDNA haplotypes found in Parasyncalathium souliei. Missing haplotypes are represented by white circles. The sizes of circles are approximately proportional to sample size. Branch lengths are indicated by the numbers along branches.


Figure 4. BEAST-derived chronogram for cpDNA haplotypes of Parasyncalathium souliei, indicating the four main haplogroups. For deeper nodes, the 95% highest posterior density interval (purple bars), the time of divergence (above branch) and posterior probabilities > 0.5 (below branch) are provided.

Because there was generally strong geographic localization of these haplogroups, four population groups were defined according to the haplogroups (Figure 2). The few populations that contained haplotypes from multiple haplogroups were assigned to population groups based on the proportion of haplotypes from each group. For example, populations KW and AL, GE and MX, DX and HS, and DD were assigned to population groups 1, 2, 3 and 4, respectively, based on the majority of individuals possessing haplotypes from those haplogroups. Population TT was assigned to population group 4 based on the large number of individuals with haplotypes from haplogroup 4, and population ZX was assigned to population group 4 based on the presence of 3 haplotypes from haplogroup 4.

Genetic Analyses

Total plastid genetic diversity in P. soulei was high (HT = 0.838), while average within-population diversity was low (HS = 0.290) (Table 2). Overall haplotype diversity (Hd) was 0.828 and nucleotide diversity (π) was 0.00685. Genetic diversity was relatively low within most populations (Appendix 3), and the highest genetic diversity was found in population GE (Hd = 0.769, π = 0.00987). The pairwise between-population differentiation (FST) was also highly variable (Appendix 4). Haplotype diversity of the four population groups (Hd) ranged from 0.604 to 0.888, while nucleotide diversity (π) ranged 0.0009–0.0140 (Table 2). Genetic differentiation was remarkably high (GST = 0.653 and NST = 0.797, P < 0.05), which is usually interpreted as indicative of phylogeographic structure. AMOVA revealed strong population variation among populations, with the highest FST value obtained when samples were divided into four population groups (FST = 0.8102, Table 3). The AMOVA revealed that 60.85% of the genetic variation was explained by differences among the four populations groups, suggesting a significant regional population substructure, and that 20.17% was partitioned among populations within population groups. Only 18.98% was partitioned within populations, consistent with the relatively low within-population diversity (HS = 0.290) (Table 3). SAMOVA failed to detect meaningful spatial population structure, as the FCT values fluctuated with increasing values of K rather than showing a distinct highest value (Appendix 5). The Mantel test (r = 0.009; P = 0.13) did not find a statistically significant pattern of isolation by distance for P. souliei.


Table 2. Results of genetic and mismatch analyses in Parasyncalathium souliei, including number of individuals (n) and kinds of haplotypes (Kh), haplotype diversity (Hd), nucleotide diversity (π), Fu's Fs (Fs), Tajima's D (DT),τ = 2ut, (where t is the expansion time and μ is the mutation rate per generation), pre-expansion population size (θ0), post-expansion population Size (θ1), sum of squared deviations (SSD, P-value), and Harpending's raggedness index (RAG, P-value).


Table 3. Structure of genetic variation within and among four population groups of Parasyncalathium souliei (FCT, differentiation among groups within the species; FSC, differentiation among populations within groups; FST, differentiation among populations within the species).

Estimates of Historical Demography and Species Distribution Modeling

Hierarchical mismatch analyses indicated that distributions of pairwise differences for population groups 2, 3, 4, and the entire species were ragged or multimodal (Figure 5), which rejected the hypothesis of demographic population expansion of these population groups. The SSD and RAG index (P < 0.05) also rejected the historical population expansion of population group 1 despite its unimodal pattern in the mismatch analysis (Table 3,Figure 5).


Figure 5. Mismatch distribution analysis plots across all populations and for the four population groups of P. souliei.

Species distribution modeling yielded an over-prediction of the geographic distribution of P. souliei to the north and west of its known range, albeit with weak confidence in these areas (Figure S2). Overall, the predicted ranges of P. souliei in the Mid-Holocene and Last Glacial Maximum were very similar to that of the predicted current distribution; the predicted range was slightly smaller during the Mid-Holocene and was shifted slightly northward during the Last Glacial Maximum but still with the high climatic suitability of HDM (Figure S2).



Compared to the genetic diversity of 170 plant species compiled by Petit et al. (2005), the plastid haplotype diversity of P. souliei (HT = 0.838) is higher than average (HT = 0.670). High plastid diversity has also been detected in other QTP-HDM species such as Eriophyton wallichii Benth. (Lamiaceae): HT = 0.983; Thalictrum squamiferum Lecoy. (Ranunculaceae): HT = 0.973; Paraquilegia microphylla (Royle) J. R. Drumm and Hutch. (Ranunculaceae): HT = 0.984; and Chionocharis hookeri (C. B. Clarke) I. M. Johnst. (Boraginaceae): HT = 0.935 (Luo et al., 2016). Chloroplast genetic differentiation among populations is relatively high (NST = 0.797; GST = 0.653, P < 0.05), which may result from the topographic complexity of the QTP-HDM. Considering the lack of spatial population genetic structure discernable from Samova suggesting such high inter-population differentiation may have been driven by genetic isolation and drift arising from dispersal limitation and topographic effects during Pleistocene full-glacial periods (Xing and Ree, 2017).

Our analyses imply that P. souliei diversified in the Pliocene to Pleistocene and that the current geographic range of the species may result from a complex pattern of long-term population stability. Pliocene to Pleistocene diversification has also been inferred in other QTP-HDM taxa such as Solms-laubachia Muschl. (Brassicaceae) (Yue et al., 2009), Baimashania Al-Shehbaz (Brassicaceae) (Koch et al., 2012), and Meconopsis Vig. (Papaveraceae) (Wen et al., 2014), suggesting that this period was an important period for species/lineage differentiation in the QTP-HDM. These dates coincide with later estimates for the age of the uplift of the QTP (Li and Fang, 1999) and the formation of the Hengduan Mountains (Shi et al., 1998; Akciz et al., 2008). The highest levels of within-population variation in P. souliei occur in the HDM, and particularly the SHDM, and this combined with the many unique haplotypes found there suggest that P. souliei populations may have been relatively stable through time in the HDM. The HDM is a region of highly heterogeneous microclimates and it has been suggested to be a refugium for many plant species during Pleistocene glacial periods, such as Taxus wallichiana Zucc. (Taxaceae) (Liu et al., 2013) and Quercus aquifolioides Rehder and E. H. Wilson (Fagaceae) (Du et al., 2017). Our species distribution modeling analyses (Figure S2) also suggest that climatic conditions favorable for the persistence of P. souliei existed across all modeled periods. In contrast with previous studies evaluating possible glacial effects on phylogeographic patterns in the QTP-HDM region (e.g., Zhang et al., 2005; Yang et al., 2008), our distribution modeling analyses indicate relative stasis of P. souliei populations during the later Pleistocene, rather than clear contraction (Figure S2). This relatively unusual demographic scenario may be explained by the wide ecological tolerance of P. souliei combined with the extreme topoclimatic complexity of the QTP-HDM region, which may have provided suitable areas for the in situ persistence of P. souliei even during glaciations. This contrasts to the general “contraction–expansion” scenario usually described in other plants (e.g., Yan et al., 2012). Away from the HDM, unique haplotypes are also found in other parts of the range of P. souliei, including the QTP, although within-population diversity is generally lower in these areas. It is important to note that the deepest phylogenetic split among haplotypes is between haplogroup 4 (group EHM, with haplotypes in the western part of the species range) and all other haplogroups (which are centered in the eastern part of the range). This implies that populations could have survived in other regions of the range away from the HDM as well, suggesting vicariance of at least some populations. However, the widespread nature of the relatively young haplotypes 1 and 3 (Figure 2) in QTP populations may reflect more recently established populations in this area.

Incorporating Haplotype Information in Phylogeographic Regionalization

Biogeographical boundaries delineate the basic macrounits of diversity in biogeography, conservation, and macroecology (Daru et al., 2016, 2017). Species composition in different biogeographic areas reflects the historical processes, such as orogeny and glaciations, that have helped to shape the present-day distribution of biodiversity (Mackey et al., 2007; Kreft and Jetz, 2010). Traditionally, biogeographical regions have been defined via qualitative and/or quantitative approaches that use geographical distribution and landscape features to delimit hierarchical floristic regions (Takhtajan, 1986; Wu and Wu, 1996; González-Orozco et al., 2014; Zhang et al., 2016). Recently, phylogenetic information has begun to be incorporated into regionalization schemes by considering phylogenetic dispersal, extinction, speciation, and niche conservatism to define assembly of species into distinct biogeographic units (Daru et al., 2017). Within the QTP-HDM hotspot region, two dominant methods have been used to delimit biogeographic regions: (1) based on the presence and absence data at the family, genus and species level (Wu and Wu, 1996), and (2) based on patterns of similarity in taxonomic composition (Zhang et al., 2016). We propose that the history of individual species, and in particular those that are widely distributed within biodiversity hotspots, may provide key information to help delineate phylogeoregions, such as in the QTP-HDM (Appendix 6). Such phylogeoregions have the advantage of incorporating information at relatively shallow (i.e., recent) time scales, and may or may not correspond to existing proposed biogeographic regions. For example, Amarilla et al. (2015) identified three clades of haplotypes within the grass Munroa argentina Griseb. that correspond to previously identified biogeographic regions in the South American Transition Zone (SATZ).

As with the Munroa example, our phylogeographic results in Parasyncalathium souliei are largely consistent with the biogeographic regions proposed by Wu and Wu (1996). Three of the four haplogroups detected in P. souliei (groups 2–4; Figure 2) are geographically localized within previously proposed subregions: haplogroup 2 was distributed through the NHDM; haplogroup 3 was located through the SHDM; and haplogroup 4 was mainly distributed in the eastern Himalayan region (hereafter called EHM), including the Three river-gorges and Southeast Xizang subregions (Figure 2). In contrast, haplotypes from haplogroup 1 are more widely distributed, ranging from the Paleotropical floristic kingdom (population CN in Figure 2) to the HDM (Wu and Wu, 1996), although most of the individual haplotypes in haplogroup 1 are highly localized (Figure 2). The wide distribution of sister haplotypes 1 and 3 is difficult to accommodate in hypotheses of phylogeographic regionalization (Figure 2). Most of the populations that are fixed for haplotypes in haplogroup 1 are at higher elevations and/or are further north, and it is possible that these populations have been established more recently by one or a few ancestral population(s) dominated by haplotypes in this haplogroup. However, it is difficult to pinpoint the cause for this pattern given available evidence.

Based on our results and a survey of the literature (e.g., Zhang et al., 2016), we propose that the distribution of private haplotypes in P. souliei may reflect a putative set of phylogeoregions that roughly follow the biogeographic regions of Wu and Wu (1996), with some important differences. Overall, we suggest that the QTP-HDM region is best divided into three phylogeoregions (I, II, and III; Figure 1B) rather than the more complicated subdivisions proposed by Wu and Wu (1996). More specifically, there are three main differences between the earlier qualitative division of Wu and Wu (1996; Figure 1A), which is based on the geographical distribution of endemic species, genera, and families, and our proposed phylogeographic divisions (Figure 1B). First, we find no support (i.e., private haplotypes; Figure 2) in P. souliei for two of the subkingdoms (Figure 1A; III: East Asia floristic kingdom; IV: Paleotropical floristic kingdom) of Wu and Wu (1996). Likewise, other phylogeographic studies that have involved this region (region IVG 24 in Figure 1A) have not found private haplotypes that distinguish these two regions (e.g., Li et al., 2011; Jia et al., 2012; Liu et al., 2013). It is therefore more reasonable to assign IVG 24 of Wu and Wu (1996) to the East Asia floristic kingdom (Figure 1).

Secondly, our results suggest that the Jinsha River may mark an important phylogeographic boundary between what we propose is a new Himalayan phylogeoregion (region III; Figure 1B) and the North and South HDM phylogeoregions (regions I and II; Figure 1B). The areas east and west of the Jinsha River are profoundly different in terms of ecology and species diversity and these differences are associated with the uplift of the QTP and the north/south barrier to the monsoon climate created by the Gaoligong and Nushan Mountains (Liu et al., 2013). This differentiation between the EHM and HDM regions is also supported by paleoclimatic evidence (Kou et al., 2006). Similarly, Zhang et al. (2016) have also assigned the Three river-gorges subregion of Wu and Wu (1996) to our proposed Eastern Himalayan phylogeoregion based on quantificationally hierarchical clustering of species endemism. Likewise, a recent phylogeographic study of Marmoritis complanata (Dunn) A. L. Budantzev (Lamiaceae) detected what the authors called the “Ward Line-Mekong-Salween Divide” (Luo et al., 2017) as an important floristic boundary, which is consistent with our placement of the old Three river-gorges area in the EHM phylogeoregion.

Thirdly, we suggest that the marked division in the geographic distribution of haplotypes between haplogroups 2 and 3 in P. souliei at approximately 29° N latitude (Figure 2) may reflect another important phylogeographic boundary in the QTP-HDM, specifically the division between the North and South HDM. Hence, we tentatively suggest that the boundary of the putative NHDM and SHDM phylogeoregions be placed at this latitude (red horizontal line in Figure 1B). Significant differences in climate and soil exist north and south of this line: the NHDM is characterized of by a drier, subfrigid plateau climate dominated by alpine meadows, whereas the SHDM is characterized by a wetter, subtropical montane climate dominated by broadleaf evergreen forest (Yu et al., 1989; Li et al., 2010). The two putative phylogeoregions have also been shown to be distinctly different in soil microbiota, especially in the number of radiobacteria (Shao et al., 2005).


Compared to qualitative methods of biogeographic regionalization based on the distribution of endemic plant species and quantitative methods that employ hierarchical clustering and ordination of plant communities, such as those of Wu and Wu (1996) and Zhang et al. (2016), our approach helps to identify potential barriers to dispersal and survival within individual taxa. However, we emphasize that our study cannot fully resolve the phylogeographic patterns in the QTP-HDM because it is based on a single species with a short life span. While such species may provide valuable insight into phylogeographic patterns because they have elevated mutation rates, a more complete picture of phylogeographic regions will require studies of many more taxa including a variety of growth forms in many different plant families. Such future studies will help bring about a fuller view of the biogeographic complexities of the QTP-HDM, and hence further improve our classification of biogeographic regions.

Author Contributions

HS, JZ, and HW conceived and designed this study; NL, TD, and YS analyzed the data, and MM, NL, TD, and YS wrote the manuscript; XH, WS, and DL aided in field collections. All authors read and approved the final manuscript.


This study was supported by grants-in-aid from the Strategic Priority Research Program of Chinese Academy of Sciences (XDA20050203 to HS), the Major Program of the National Natural Science Foundation of China (31590823 to HS), the National Natural Science Foundation of China (NSFC, 31000101, 31370004 and 31570213 to JZ), the National Key R&D Program of China (2017YFC0505200 to HS), the International Program of Chinese Academy of Science (151853KYSB2016021 to HS).

Conflict of Interest Statement

The 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.


We are indebted to Prof. J. Luo, Drs. Y. H. Zhang, J. G. Chen, X. H. Li, B. Song, B. Xu, Z. Zhou, and Mr. Q. B. Gong for help with sampling, and to C. R. Ouyang for lab work.

Supplementary Material

The Supplementary Material for this article can be found online at:

Figure S1. The 159 species occurrence points (black dots in figure) were for species distribution modeling.

Figure S2. Climate-based species distribution models of P. souliei for three geological time periods: present day, Mid-Holocene (c. 6 ka), and Last Glacial Maximum (c. 22 ka).

Appendix 1. List of sampling localities, herbarium voucher specimens, and GenBank accession numbers. Herbarium acronyms is as follows: KUN = Herbarium of Kunming Institute of Botany, Chinese Academy of Sciences.

Appendix 2. Chloroplast DNA sequence polymorphism detected at mutation sites at haplotypes in P. souliei.

Appendix 3. Genetic diversity within populations and correlations test. Genetic diversity within populations. (Hd indicates haplotype diversity; π indicates nucleotide diversity).

Appendix 4. Pairwise FST (blow diagonal) and P-values (above diagonal) among populations of P. souliei.

Appendix 5. Table and Figure of relationship between the fixation index (FCT; y-axis) and number of clusters (K; x-axis) in the SAMOVA of P. souliei.

Appendix 6. Several cases of biogeographic regionalization based on one species from Qinghai-Tibet Plateau (QTP) - Hengduan Mountains (HDM) region.


Akciz, S., Burchfiel, B. C., Crowley, J. L., Yin, J., and Chen, L. (2008). Geometry, kinematics, and regional significance of the chong shan shear zone, Eastern Himalayan Syntaxis, Yunnan, China. Geosphere 4, 292–314. doi: 10.1130/GES00111.1

CrossRef Full Text | Google Scholar

Amarilla, L. D., Anton, A. M., Chiapella, J. O., Manifesto, M. M., Angulo, D. F., and Sosa, V. (2015). Munroa argentina, a grass of the South American Transition Zone, survived the Andean uplift, aridification and glaciations of the Quaternary. PLoS ONE 10:e0128559. doi: 10.1371/journal.pone.0128559

PubMed Abstract | CrossRef Full Text | Google Scholar

Bandelt, H. J., Forster, P., and Hamburg, A. R. (1999). Median-Joining networks for inferring intraspecific phylogenies. Mol. Biol. Evol. 16, 37–48. doi: 10.1093/oxfordjournals.molbev.a026036

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y. S., Deng, T., Zhou, Z., and Sun, H. (2018). Is the East Asian flora ancient or not?National Science Review nwx156. doi: 10.1093/nsr/nwx156

CrossRef Full Text

Collins, W. D., Blackmon, M., Bitz, C. G. B., and Bretherton, C. S. (2004). The community climate system model version 3 (CCSM3). J. Climate 19, 2122–2143. doi: 10.1175/JCLI3761.1

CrossRef Full Text | Google Scholar

Daru, B. H., Elliott, T. L., Park, D. S., and Davies, T. J. (2017). Understanding the processes underpinning patterns of phylogenetic regionalization. Trends Ecol. Evol. 32, 845–860. doi: 10.1016/j.tree.2017.08.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Daru, B. H., van der Bank, M., Maurin, O., Yessoufou, K., Schaefer, H., Slingsby, J. A., et al. (2016). A novel phylogenetic regionalization of phytogeographical zones of southern Africa reveals their hidden evolutionary affinities. J. Biogeogr. 43, 155–166. doi: 10.1111/jbi.12619

CrossRef Full Text | Google Scholar

Drummond, A. J., and Rambaut, A. (2007). BEAST: Bayesian Evolutionary Analysis by Sampling Trees. BMC Evol. Biol. 7, 214. doi: 10.1186/1471-2148-7-214

PubMed Abstract | CrossRef Full Text | Google Scholar

Du, F. K., Hou, M., Wang, W., Mao, K., and Hampe, A. (2017). Phylogeography of Quercus aquifolioides provides novel insights into the Neogene history of a major global hotspot of plant diversity in south-west China. J. Biogeogr. 44, 294–307. doi: 10.1111/jbi.12836

CrossRef Full Text | Google Scholar

Dupanloup, S., Schneider, T., and Excoffier, L. (2002). A simulated annealing approach to define the genetic structure of populations. Mol. Ecol. 11, 2571–2581. doi: 10.1046/j.1365-294X.2002.01650.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Excoffier, L., Laval, G., and Schneider, S. (2005). Arlequin (version 3.0): an integrated software package for population genetics data analysis. Evol. Bioinform. 1, 47–50. doi: 10.1177/117693430500100003

PubMed Abstract | CrossRef Full Text

Favre, A., Päckert, M., Pauls, S. U., Jähnig, S. C., Uhl, D., Michalak, I., et al. (2015). The role of the uplift of the Qinghai-Tibetan Plateau for the evolution of Tibetan biotas. Biol. Rev. 90, 236–253. doi: 10.1111/brv.12107

PubMed Abstract | CrossRef Full Text | Google Scholar

Fu, X. Y. (1997). Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics 147, 915–925.

PubMed Abstract | Google Scholar

González-Orozco, C. E., Thornhill, A. H., Knerr, N., Laffan, S., Miller, J. T., and Richardson, D. (2014). Biogeographical regions and phytogeography of the eucalypts. Divers. Dist. 20, 46–58. doi: 10.1111/ddi.12129

CrossRef Full Text | Google Scholar

Hall, T. A. (1999). BioEdit: a user-friendly biological sequence alignment editor and analysis program for windows. Nucleic Acids Symp. Ser. 41, 95–98.

Google Scholar

Harpending, H. C. (1994). Signature of ancient population growth in a low-resolution mitochondrial DNA mismatch distribution. Hum. Biol. 66, 591–600.

PubMed Abstract | Google Scholar

Hasumi, H., and Emori, S. (2004). K-1 Coupled GCM (MIROC) Description. Tokyo: Center for Climate System Research, University of Tokyo.

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

Jia, D. R., Abbott, R. J., Liu, T. L., Mao, K. S., Bartish, I. V., and Liu, J. Q. (2012). Out of the Qinghai–Tibet Plateau: evidence for the origin and dispersal of Eurasian temperate plants from a phylogeographic study of Hippophae rhamnoides (Elaeagnaceae). New Phytologist 194, 1123–1133. doi: 10.1111/j.1469-8137.2012.04115.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Katoh, K., Kei, K. M., Kuma, K. I., and Miyata, T. (2002). MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Res. 30, 3059-3066doi: 10.1093/nar/gkf436.

PubMed Abstract | CrossRef Full Text | Google Scholar

Kilian, N., Sennikov, A., Wang, Z. H., Gemeinholzer, B., and Zhang, J. W. (2017). Sub-Paratethyan origin and middle to Late Miocene principal diversification of the Lactucinae (Compositae: Cichorieae) inferred from molecular phylogenetics, divergence-dating and biogeographic analysis. Taxon 66, 675–703. doi: 10.12705/663.9

CrossRef Full Text | Google Scholar

Koch, M. A., Karl, R., German, D. A., and Al-Shehbaz, I. A. (2012). Systematics, taxonomy and biogeography of three new Asian genera of Brassicaceae tribe Arabideae: an ancient distribution circle around the Asian high mountains. Taxon 51, 955–969.

Google Scholar

Kou, X. Y., Ferguson, D. K., Xu, J. X., Wang, Y. F., and Li, C. S. (2006). The reconstruction of paleovegetation and paleoclimate in the Late Pliocene of West Yunnan, China. Clim. Change 77, 431–448. doi: 10.1007/s10584-005-9039-5

CrossRef Full Text | Google Scholar

Kreft, H., and Jetz, W. (2010). A framework for delineating biogeographical regions based on species distributions. J. Biogeogr. 37, 2029–2053. doi: 10.1111/j.1365-2699.2010.02375.x

CrossRef Full Text | Google Scholar

Li, J. J., and Fang, X. M. (1999). Uplift of the Tibetan Plateau and environmental changes. Chin. Sci. Bull. 44, 2117–2124. doi: 10.1007/BF03182692

CrossRef Full Text | Google Scholar

Li, Y., Zhai, S. N., Qiu, Y. X., Guo, Y. P., Ge, X. J., and Comes, H. P. (2011). Glacial survival east and west of the 'Mekong-Salween Divide' in the Himalaya-Hengduan Mountains region as revealed by AFLPs and cpDNA sequence variation in Sinopodophyllum hexandrum (Berberidaceae). Mol. Phylogenet. Evol. 59, 412–424. doi: 10.1016/j.ympev.2011.01.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Z. X., He, Y. Q., Xin, H. J., Wang, C. F., Jia, W. X., Zhang, W., et al. (2010). Spatio-temporal variations of temperature and precipitation in Mts. Hengduan region during 1960-2008. Acta Geographica Sinica 65, 564–579. doi: 10.1063/1.4973107

CrossRef Full Text | Google Scholar

Librado, P., and Rozas, J. (2009). DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics 25, 1451–1452. doi: 10.1093/bioinformatics/btp187

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, J., Möller, M., Provan, J., Gao, L. M., Poudel, R. C., and Li, D. Z. (2013). Geological and ecological factors drive cryptic speciation of yews in a biodiversity hotspot. New Phytologist 199, 1093–1108. doi: 10.1111/nph.12336

PubMed Abstract | CrossRef Full Text | Google Scholar

Luo, D., Xu, B., Li, Z., and Sun, H. (2017). The 'Ward Line - Mekong Salween Divide' is an imporant floristic boundary between the eastern Himalaya and Hengduan Mountains: evidence from the phylogeographical structure of subnival herbs Marmoritis companatum (Lamiacae). Bot. J. Linn. Soc. 4, 482–496. doi: 10.1093/botlinnean/box067

CrossRef Full Text | Google Scholar

Luo, D., Yue, J. P., Sun, W. G., Xu, B., Li, Z. M., Comes, H. P., et al. (2016). Evolutionary history of the subnival flora of the Himalaya-Hengduan Mountains: first insights from comparative phylogeography of four perennial herbs. J. Biogeogr. 43, 31–43. doi: 10.1111/jbi.12610

CrossRef Full Text | Google Scholar

Mackey, B. G., Berry, S. L., and Brown, T. (2007). Reconciling approaches to biogeographical regionalization: a systematic and generic framework examined with a case study of the Australian continent. J. Biogeogr. 35, 213–229. doi: 10.1111/j.1365-2699.2007.01822.x

CrossRef Full Text | Google Scholar

Myers, N., Mittermeier, R. A., Mittermeier, C. G., and Kent, J. (2000). Biodiversity hotpots for conservation priorities. Nature 403, 8453–8858. doi: 10.1038/35002501

PubMed Abstract | CrossRef Full Text | Google Scholar

Nylander, J. A. A. (2004). MrModeltest v2. program distributed by the author. Genet. Soc. Am. 123, 597–601. doi: 10.1109/icdcs.1994.302463

CrossRef Full Text

Peakall, R. O. D., and Smouse, P. E. (2006). GENALEX 6: genetic analysis in Excel. Population genetic software for teaching and research. Mol. Ecol. 6, 288–295. doi: 10.1111/j.1471-8286.2005.01155.x

CrossRef Full Text

Petit, R. J., Duminil, J., Fineschi, S., Hampe, A., Salvini, D., and Vendramin, G. G. (2005). Comparative organization of chloroplast, mitochondrial and nuclear diversity in plant populations. Mol. Ecol. 14, 689–701. doi: 10.1111/j.1365-294X.2004.02410.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Phillips, S. J., Anderson, R. P., and Schapire, R. E. (2006). Maximum entropy modeling of species geographic distributions. Ecol. Modell. 190, 231–259. doi: 10.1016/j.ecolmodel.2005.03.026

CrossRef Full Text | Google Scholar

Ren, G., Mateo, R. G., Liu, J., Suchan, T., Alvarez, N., Guisan, A., et al. (2017). Genetic consequences of Quaternary climatic oscillations in the Himalayas: Primula tibetica as a case study based on restriction site-associated DNA sequencing. New Phytologist 213, 1500–1512. doi: 10.1111/nph.14221

PubMed Abstract | CrossRef Full Text | Google Scholar

Robert, L., Simon, Y. W. H., Davies, T. J., Angela, T. M., Lonnie, A., Nathan, G. S., et al. (2013). Taller plants have lower rates of molecular evolution. Nat. Commun. 4:1879. doi: 10.1038/ncomms2836

CrossRef Full Text | Google Scholar

Shao, B. L., Gong, G. S., Zhang, S. R., Yu, X., Yang, D., Zhang, H., et al. (2005). Soil microbial quantity and its relations with ecological factors in northern alp region of Hengduan Mountains. Chin. J. Ecol. 25, 885–890.

Shi, Y. F., Li, J. J., and Li, B. Y. (1998). Uplift and Environmental Changes of Qinghai-Xizang (Tibetan) in the Late Cenozoic. Guangzhou: Guangdong Science and Technology Press.

Sun, H., Zhang, J. W., Deng, T., and Boufford, D. E. (2017). Origins and evolution of plant diversity in the Hengduan Mountains, China. Plant Divers. 39, 161–166. doi: 10.1016/j.pld.2017.09.004

CrossRef Full Text | Google Scholar

Tajima, F. (1989). The effect of change in population size on DNA polymorphism. Genet. Soc. Am. 111, 597–601.

Google Scholar

Takhtajan, A. (1986). Floristic Regions of the World. Berkeley, CA: University of California Press.

Google Scholar

Wen, J., Zhang, J. Q., Nie, Z. L., Zhong, Y., and Sun, H. (2014). Evolutionary diversifications of plants on the Qinghai-Tibetan Plateau. Front. Genet. 5, 12–16. doi: 10.3389/fgene.2014.00004

PubMed Abstract | CrossRef Full Text

Wu, Z. Y., and Wu, S. G. (1996). “A proposal for a new floristic kingdom (realm): the East Asiatic Kingdom, its delineation and characteristics,” in Floristic Characteristics and Diversity of East Asian Plants, Proceedings of the First International Symposium on Floristic Characteristics and Diversity of East Asian Plants, eds A. L. Zhang and S. G. Wu (Beijing: China Higher Education Press), 3–42.

Xing, Y. W., and Ree, H. R. (2017). Uplift-driven diversification in the Hengduan Mountains, a temperate biodiversity hotspot. Proc. Natl. Acad. Sci. U.S.A. 114, E3444–E3451. doi: 10.1073/pnas.1616063114

PubMed Abstract | CrossRef Full Text | Google Scholar

Yan, H. F., Zhang, C. Y., Wang, F. Y., Hu, C. M., Ge, X. J., and Hao, G. (2012). Population expanding with the phalanx model and lineages split by environmental heterogeneity: a case study of Primula obconica in subtropical China. PLoS ONE 7:e41315. doi: 10.1371/journal.pone.0041315

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, F. S., Li, Y. F., Ding, X., and Wang, X. Q. (2008). Extensive population expansion of Pedicularis longiflora (Orobanchacea) on the Qinghai-Tibetan Plateau and its correlation with the Quaternary climate change. Mol. Ecol. 17, 5135–5145. doi: 10.1111/j.1365-294X.2008.03976.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, Z., and Rannala, B. (1997). Bayesian phylogenetic inference using DNA sequences: a Markov chain monte carlo method. Mol. Biol. Evol. 14, 717–724. doi: 10.1093/oxfordjournals.molbev.a025811

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, Y. D., Liu, L. H., and Zhang, J. H. (1989). Vegetational regionalization of the Hengduan Mountains regions. Mountain Res. 7, 47–55.

Yue, J. P., Sun, H., Baum, D. A., Li, J. H., Al-Shehbaz, I. A., and Ree, R. (2009). Molecular phylogeny of Solms-laubachia (Brassicaceae) s.l., based on multiple nuclear and plastid DNA sequences, and its biogeographic implications. J. Syst. Evol. 47, 402–415. doi: 10.1111/j.1759-6831.2009.00041.x

CrossRef Full Text | Google Scholar

Zhang, D. C., Ye, J. X., and Sun, H. (2016). Quantitative approaches to identify floristic units and centres of species endemism in the Qinghai-Tibetan Plateau, south-western China. J. Biogeogr. 43, 2465–2476. doi: 10.1111/jbi.12819

CrossRef Full Text | Google Scholar

Zhang, J. W., Boufford, D. E., and Sun, H. (2011a). Parasyncalathium, J.W. Zhang, Boufford and H. Sun (Asteraceae, Cichorieae): a new genus endemic to the Himalaya-Hengduan Mountains.Taxon 60, 1678–1684.

Zhang, J. W., Nie, Z. L., Wen, J., and Sun, H. (2011b). Molecular phylogeny and biogeography of three closely related genera, Soroseris, Stebbinsia, and Syncalathium (Asteraceae, Cichorieae), endemic to the Tibetan Plateau, SW China. Taxon 60, 15–26.

Google Scholar

Zhang, Q., Chiang, T. Y., George, M., Liu, J. Q., and Abbott, R. J. (2005). Phylogeography of the Qinghai-Tibetan Plateau endemic Juniperus przewalskii (Cupressaceae) inferred from chloroplast DNA sequence variation. Mol. Ecol. 14, 3513–3524. doi: 10.1111/j.1365-294X.2005.02677.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: Qinghai-Tibet Plateau (QTP)-Hengduan Mountains (HDM), Parasyncalathium souliei, phylogeography, genetic distribution, phylogeoregion

Citation: Lin N, Deng T, Moore MJ, Sun Y, Huang X, Sun W, Luo D, Wang H, Zhang J and Sun H (2018) Phylogeography of Parasyncalathium souliei (Asteraceae) and Its Potential Application in Delimiting Phylogeoregions in the Qinghai-Tibet Plateau (QTP)-Hengduan Mountains (HDM) Hotspot. Front. Genet. 9:171. doi: 10.3389/fgene.2018.00171

Received: 05 February 2018; Accepted: 27 April 2018;
Published: 17 May 2018.

Edited by:

Kangshan Mao, Sichuan University, China

Reviewed by:

Igor V. Bartish, Institute of Botany (ASCR), Czechia
Gonzalo Nieto Feliner, Real Jardín Botánico (RJB), Spain

Copyright © 2018 Lin, Deng, Moore, Sun, Huang, Sun, Luo, Wang, Zhang and Sun. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Hengchang Wang,
Jianwen Zhang,
Hang Sun,

These authors have contributed equally to this work.