Next Article in Journal
Reading between the Lines: RNA-seq Data Mining Reveals the Alternative Message of the Rice Leaf Transcriptome in Response to Heat Stress
Next Article in Special Issue
Different Roles of Introgression on the Demographic Change in Two Snakebark Maples, Acer caudatifolium and A. morrisonense, with Contrasted Postglacial Expansion Routes
Previous Article in Journal
Pressurized Hot Water Extraction of Okra Seeds Reveals Antioxidant, Antidiabetic and Vasoprotective Activities
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Environmental Heterogeneity Leads to Spatial Differences in Genetic Diversity and Demographic Structure of Acer caudatifolium

Department of Life Science, National Taiwan Normal University, Taipei 116059, Taiwan
*
Author to whom correspondence should be addressed.
These authors contributed equally to this work.
Plants 2021, 10(8), 1646; https://doi.org/10.3390/plants10081646
Submission received: 22 June 2021 / Revised: 2 August 2021 / Accepted: 9 August 2021 / Published: 10 August 2021
(This article belongs to the Special Issue Plant Molecular Evolution and Population Ecology)

Abstract

:
Under climate fluctuation, species dispersal may be disturbed by terrain and local climate, resulting in uneven spatial-genetic structure. In addition, organisms at different latitudes may be differentially susceptible to climate change. Here, we tracked the seed dispersal of Acer caudatifolium using chloroplast DNA to explore the relationships of terrain and local climate heterogeneity with range shifts and demography in Taiwan. Our results showed that the extant populations have shifted upward and northward to the mountains since the Last Glacial Maximum. The distributional upshift of A. caudatifolium is in contrast to the downward expansion of its closest relative in Taiwan, A. morrisonense. The northern populations of A. caudatifolium have acquired multiple-source chlorotypes and harbor high genetic diversity. However, effective gene flow between the north and south is interrupted by topography, geographic distance, north-south differences in October rainfall, and other climate heterogeneities, blocking southward genetic rescue. In addition, winter monsoon-driven rainfall may cause regional differences in the phenological schedule, resulting in adaptive effects on the timing of range shift and the genetic draft of chlorotype distribution. Terrain, distance, and local climate also differentiate the northernmost populations from the others, supporting the previous taxonomic treatment of Acer kawakamii var. taitonmontanum as an independent variety.

1. Introduction

Both environmental conditions and dispersibility may be limiting factors for changes in the distribution ranges of species. The measurement of ecological niche is one way to link the suitable environmental conditions for species in which their fitness can be maintained or increased. These theoretical niche ranges are the fundamental (Grinnellian) niches [1] and could be predicted by the ecological niche modeling (ENM) under the assumption of phylogenetic niche conservatism (PNC). However, this approach does not consider life-history traits and limitations of dispersibility [2]. For sessile organisms, dispersibility could be relatively important for range change. As representative sessile organisms, plants can overcome spatial constraints through seed dispersal with variable strategies [3,4]. Consequently, seed dispersibility is one of the determinants of the ability of plants to cope with climate change.
Gene flow and demographic dynamics may be subject to spatial and environmental influences [5,6]. Given the rugged terrain distribution and steep environmental differences in mountain regions, the seed dispersal and pollen flow of mountain plants and their temporal range dynamics may jointly reflect topography and climate [7]. According to the center-periphery hypothesis [8,9], range change (e.g., spatial expansion) may be more adaptively constrained at range margins [10]. The core-to-edge decline of an adaptive genetic variation supports seed-based dispersal of mountainous maples in Taiwan rather than pollen flow [7]. However, whether the seed dispersal of maples in the same area but at lower elevations is also spatially and climate-related under climate change remains unknown.
In this study, an island maple, Acer caudatifolium Koidzumi, was selected as the research object. Acer caudatifolium is widely distributed at low to high elevations, is endemic to the island of Taiwan [11], and is a relative of temperate maples in continental Asia [12,13,14]. Taiwan is a mountainous continental island situated off the southeastern Asian Continent. Due to high seed dispersibility, the population genetic differentiation of A. caudatifolium is expected to be low. However, according to the mountain barrier hypothesis, mountains may hinder long-distance seed dispersal of low-elevation plants [7], resulting in population differentiation among mountain islands. In addition, the varied mountainous environment could affect population structure, interrupting gene flow. Given its massif distribution from low to high elevation and high dispersibility, A. caudatifolium is a suitable object for exploring the issue of spatial-environmental variation in demographic and range change.
This study asked the following two questions: (1) What are the spatial patterns of genetic diversity and population differentiation of A. caudatifolium? (2) Are the population genetic structure and demographic dynamics of A. caudatifolium associated with geographic or environmental factors? To answer these questions, we reconstructed the current and paleo distributions of A. caudatifolium using maternally inherited chloroplast DNA (cpDNA) to track seed dispersal and demographic dynamics, which were then correlated with geographic and environmental factors. The results illustrate how spatial and environmental heterogeneities determine the population genetic structure of mountain trees under climate change.

2. Materials and Methods

2.1. Ecological Niche Modeling

To estimate the potential distribution of Acer caudatifolium, we performed ecological niche modeling (ENM). We used the current distribution and projected with the same climatic variables to the middle Holocene (Holocene thermal optimal, HTO, ~6 kya) and the Last Glacial Maximum (LGM, ~21 kya) to construct the putative paleodistribution. Since the available paleoclimatic variables in the open database were limited and needed to be consistent with current variables while constructing paleo-ENM, we extracted current (1960–1990) climatic factors (19 bioclims, 24 monthly temperatures (Tmin and Tmax), and 12 monthly precipitations) in spatial resolution of 2.5 arc-minutes (approximately 4.5 km × 4.5 km at the equator) from WorldClim (https://www.worldclim.org/; accessed on 26 July 2021). The corresponding paleoclimatic variables in the HTO and LGM were also downloaded from the WorldClim database.
Variables with variance inflation factor (VIF) >6 were removed to avoid multicollinearity using the R package usdm [15]. Only six climatic variables (bio3, bio7, prec2, prec6, prec9, and tmax2) were retained (Table S1). ENM was conducted under the maximum entropy model in MaxEnt GUI v.3.4.1 [16]. The occurrence data referred to the Global Biodiversity Information Facility (GBIF) database and our sampling records. Unreliable records (e.g., in urban regions or coasts) were discarded. To avoid extreme weighting, only one of multiple records within a grid was retained. Twenty percent of the sampling data were treated as the testing dataset. We run the ENM with 10 bootstrapped replicates by default settings (convergence threshold of 1 × 10−5, 10,000 background points, a maximum of 500 iterations, and a prevalence of 0.5). To evaluate the performance of the predicted niche model, we considered the area under the receiver operating characteristic curve (AUC) value from MaxEnt and also estimated the partial ROC value using NicheToolBox (http://shiny.conabio.gob.mx:3838/nichetoolb2/; accessed on 26 July 2021) with 1000 bootstrapped replicates and an E-value of 0.05.

2.2. Sampling and Chloroplast DNA Sequences

The population sampling included a distribution range of A. caudatifolium from low (381 a.s.l.) to high elevation (2987 a.s.l.), an individual has collected at least 10m away from each other without replication. For every individual, the leave was dried in silica gel and stored at 4 °C for DNA extraction. Total genomic DNA was extracted following the modified cetyltrimethylammonium bromide (CTAB) method [17] and stored at −20 °C with 1 × TE buffer. Two cpDNA fragments, trnH-psbA spacer and rpl16 intron, were sequenced and amplified: the former with primers 5′-GTTATGCATGAACGTAATGCTC-3′ and 5′ -CGCGCATGGTGGATTCACAATCC-3′ and the latter with primers 5′-GCTATGCTTAGTGTGTGACTCGTTG-3′ and 5′-CTTCCTCTATGTTGTTTACG-3′. The sequencing was conducted by ABIPRISMH® 3730XL DNA Sequencer (Perkin-Elmer, Foster City, CA, USA) with ExoSAP-IT (Thermo Fisher Scientific Inc., Waltham, MA, USA) and the ABI BigDye 3.1 Terminator Cycle Sequencing Kit (Applied Biosystem, Foster City, CA, USA). The BioEdit [18] was applied for sequence quality check.
Concatenated chloroplast DNA (cpDNA) sequences with a total length of 1579 bp (440 bp from trnH-psbA and 1139 bp from rpl16) were generated from 294 samples in 19 populations (Table 1). Indels ≥ 2 nucleotides were recorded as a single mutation event, and 11 variant regions (sites) were ultimately generated. The obtained sequences are deposited in NCBI GenBank (accession numbers: trnH-psbA: MZ275974–MZ276267 rpl16: MZ275679–MZ275972).

2.3. Genetic Diversity and Haplotype Network

Indices of genetic diversity (FST, θs, and θπ) and the site-frequency spectrum (Tajima’s D and Fu’s Fs) were estimated by Arlequin v3.5.1.3 [19]. Tajima’s D and Fu’s Fs estimated the deviance between segregating sites and nucleotide diversity to infer the demographic change against the null constant size model. The increasing rare alleles after population expansion will inflate θs relative to θπ and resulted in negative Tajima’s D and Fu’s Fs. To show the relationships between cpDNA haplotypes (chlorotypes), a pairwise-distance haplotype network was built with the minimum spanning tree (MST) algorithm in Arlequin v3.5.1.3 [19].

2.4. Mismatch Analysis

Mismatch analysis was performed in every polymorphic population under both demographic and spatial expansion models in Arlequin v3.5.1.3 [19]. The former model assumes population size variation (θ0 and θ1) at time τ with no gene flow (m = 0), while the latter assumes gene flow (m) among demes at time τ with constant population size (θ). The mismatch analysis took advantage of pairwise differences between haplotypes within populations and estimates the skewed distribution of mismatch sites to infer demography change against the null expansion model. The sum of the square deviation (SSD) and the raggedness index (Rag) was used to test the deviation from the expectation of 1000 simulations. Since the method of evaluating “expansion” in the spatial expansion model assumes m > 0, the “expansion” is regarded as the “range change”. In addition, because expansion times may be overestimated by orders of magnitude [20], we did not calculate the exact time of expansion t with the formula t = τ/2μk (where μ and k are the mutation rate and sequence length) but discuss the relative time (τ) only.

2.5. Genetic Barriers

Genetic barriers were reconstructed based on geographic networks and genetic distances under Monmonier’s algorithm in the R package adegenet [21]. The method utilized Monmonier’s algorithm to identify between-group differences (i.e., barrier) based on the geographic networks and genetic distances. The K Nearest Neighbors (NN), Delaunay Triangulation (DT), and Gabriel Graph (GG) models were applied to build the population geographic connections. In the NN model, we used one-third of the data points (sampling sites) as the criterion to set the neighbors (i.e., K = 6). The genetic distance among populations was calculated by Nei’s distance. The optimize.monmonier function was used to compute the boundaries with 20 independent runs. The threshold of local difference was set to 0 to seek all potential barriers.

2.6. Genetic Differentiation across Geography and Environment

The maximum likelihood population effects mixed-effects model (MLPE) in the R package lme4 [22] was used to test whether population differentiation (gen in model) was affected by geographic distance (geo in model) (i.e., gen~geo + (1|pop), isolation by distance, IBD), adaptability to different altitudes (alt in model) (i.e., gen~alt + (1|pop), isolation by altitudinal difference, IBAlt) or different environments (env in model) (i.e., gen~env + (1|pop), isolation by environment, IBE). The combined effects among IBAlt, IBD, and IBE were also considered in model selection (Table 2). Besides the environmental variables used for ENM, we added the monthly average temperature (Tavg), solar radiation (srad), windspeed (wind), vaporization (vapr) from WorldCLim and mean global annual aridity index (GAI), actual evapotranspiration (AET), potential evapotranspiration (PET) from Global Aridity and PET Database. Among these 106 variables (Table S2), we also removing multicollinear factors with VIF values > 6 [19] and left only five factors: precipitation in October [prec10], solar radiation in June and July [srad6 and srad7], mean annual actual evapotranspiration [AET], and the global annual aridity index [GAI] (Table S1).
Genetic distance was calculated as FST/(1 − FST) (Table S3) [23]. The Euclidean and Canberra distances were used to calculate the geographic and altitudinal distances and the environmental differences (Tables S4–S6), respectively. The Akaike and Bayesian information criteria (AIC and BIC) values were ranked, and the model with the smallest AIC (or BIC) was selected as optimal. If two or three models had very similar AIC (or BIC) values, the likelihood ratio test was conducted to test whether the more complicated model rejected the simpler one. Subsequently, the Mantel test was performed to test whether the predictor (i.e., geographic distance) of the optimal model (i.e., IBD, see Results) significantly explained the genetic divergence. One thousand time permutations were conducted to test the significance of the correlation between the genetic and geographic distances. We also performed the partial Mantel test to exclude environmental interference when testing IBD.

2.7. Factors Affecting Demographic Dynamics

Population genetic patterns can also be influenced by demographic dynamics, which can be spatially and environmentally related. We, therefore, tested whether the demographic change was associated with geographic and environmental factors. The generalized linear model (GLM) was conducted in the R package stats [24], and Tajima’s D and the demographic or spatial expansion time (τ) were selected as the responses. The latitude, longitude, and altitude were taken as the geographic predictors, and the five climatic factors (prec10, srad6, srad7, AET, and GAI) were used as the environmental predictors. The optimal model was chosen using backward selection by the stepAIC function in the R package MASS [25].

3. Results

3.1. Genetic Diversity and Population Structure

The cpDNA estimation revealed that northern Taiwan is a hotspot of genetic diversity, whereas the southern populations are almost monomorphic in chlorotype (Hap6). The northern population comprises Hap4, Hap6, and their derived chlorotypes (Figure 1). Mismatch analysis showed that most of the polymorphic populations could not reject both demographic and spatial expansions according to Rag and SSD, exception the spatial expansion of the YMS and MC populations (SSD = 0.026 and 0.025, p = 0.001 and 0.012, respectively, Table 3 and Table S7). The extent of population expansion varied, as did the expansion times (τ = 0.643–8.350 and 0.135–7.924 in the demographic and spatial expansion models, respectively). These estimates do not include monomorphic populations, which does not imply constant demography of these populations but rather an inability to perform mismatch analysis.

3.2. Genetic Barriers

Inference of population genetic barriers by the Monmonier algorithm suggested that most of the geographic barriers are in northern Taiwan (Figure 1C–E). The NN model allowed the most corridors (connections) among populations, suggesting the fewest barriers; the DT and GG models allowed relatively fewer corridors and indicated more and similar barrier patterns. The main barrier was the same in the three models and surrounded the northeast population RF in Taiwan, although this population was still connected with the YMS population in the NN model. Both the DT and GG models isolated YMS from the other populations. RF and YMS directly face the winter northeast monsoon, which brings cold and humid air and precipitation in winter. The GG model suggested that the MF population in central Taiwan was isolated from the south. Most populations south of this barrier are monomorphic, with Hap6 as the dominant chlorotype, whereas most of the northern populations are polymorphic, with Hap4 as the dominant chlorotype (Figure 1A).

3.3. Geographic Distance as the Source of Population Differentiation

Model selection in MLPE indicated that IBD had the lowest AIC and BIC, followed by IBD + IBE, which had a very similar AIC value (ΔAIC = 0.4, Table 2). IBD could not be rejected by IBD + IBE in LRT (2ΔL = 1.549, df =1, p = 0.213, Table S8), suggesting that IBD is the optimal (best) model explaining population genetic differentiation. Although the influences of the combined effects of geographic with environmental differences (IBD + IBE), altitudinal differences (IBD + IBAlt), and both differences (IBD + IBE + IBAlt) did not exceed geographic distance only (IBD), their AIC and logLik were close to the best IBD model (Table 2). Neither pure IBE nor IBAlt can reject IBD. However, we weighted all factors equally to explore their influence on genetic differentiation, but these factors may have different weights in reality. It must be noted that the operation of weighting factors may lead to uncertainty in the judgment of similar models.
In analyses of Mantel and partial Mantel tests, the positive correlation between the genetic and geographic distances was only marginally significant (p = 0.096 and 0.059, respectively, Figure S1). Although far geographic distance does not necessarily indicate high differentiation, populations with higher genetic distances must have longer geographic distances, meaning that a long geographic distance is a necessary rather than sufficient condition for high genetic differentiation.

3.4. Demographic Dynamics Are Spatially and Environmentally Related

Model selection in GLM suggested that the null (empty) model could not be rejected by any alternative model in demographic expansion time. However, the optimal models for predicting the extent of population dynamics (Tajima’s D) and spatial expansion time comprised both spatial and environmental factors, although certain factors were excluded under backward selection (Table S9). Under these optimal models, latitude, prec10, srad6, srad7, and AET significantly predicted Tajima’s D (p = 0.003, 0.07, 0.021, 0.008, and 0.013, respectively), whereas only prec10 significantly predicted spatial expansion time (p = 0.010) (Table 4). Together with the Mantel test, these results indicate that although gene flow may occur between remote populations, spatial expansion is still constrained by the environment.

3.5. Upward and Northward Expansion of the Distribution Range

ENM predicted a potential distribution very similar to the exact distribution, with the largest percent contribution from Tmax2 (maximum temperature of February, 70.8%) (Figure 2A,B). Both AUC (mean AUC = 0.88) and partial ROC (mean = 0.913, p < 0.001) indicates the predicted models of Acer caudatifolium are reliable. The greatest contributing factor, Tmax2, was highly correlated with temperature-related Bioclims (Bio1, 6, 9, and 11), Tavg, Tmax, Tmin, vapr, AET, and PET (r > 0.9, p < 0.05) in Pearson’s correlation (Table S10). According to this map, the Yushan Mountain Range (populations TTC and ALS) is the most suitable habitat for the current distribution. When projected to the paleoclimates, the paleodistribution was restricted to the west side of the southern Central Mountain Range in the HTO, and to the southernmost Central Mountain Range in the LGM (Figure 2C,D). These predictions suggested that A. caudatifolium may be upward and northward to the higher-altitudinal mountainous regions with range expansion since the LGM.

4. Discussion

4.1. Paleodistribution and Climate Change Affect the Geodistance-Related Genetic Structure

The present analysis demonstrated upward and northward expansion of the distribution of A. caudatifolium into higher mountain ranges from the narrow habitats in the southwestern area of Taiwan. The ENM speculates that the south was more suitable for A. caudatifolium growth than north-central Taiwan, although the northeast may also have potentially suitable habitats (Figure 2). Restricted range in LGM may be responsible for the spatial difference in genetic diversity between northern and southern Taiwan. The north-south separation of ancestral distribution is also reflected in the genetic differentiation related to the geographical distance between the north and south populations (IBD model). The high genetic variation and differentiation in the north and monomorphic characteristics in the south indicate a correlation of latitude with demography. In addition, the population structure and spatial-temporal expansion are affected by multiple climatic factors, revealing an environmental impact of spatial-genetic structure on island mountain trees, i.e., environmentally driven genetic draft [26,27]. Climate-related population structure has also been observed in Taiwan for skullcap flowers [28,29], Rhododendron [30,31], and cow-tail fir [32]. However, the gene flow of A. caudatifolium is still related to geographic distance, albeit marginally. However, geographic distance alone may not be sufficient to explain the population differentiation caused by genetic drift in the short term. Other factors, such as the spatial and temporal differences of local climate, may also limit population distribution and gene flow. In other words, climate change determines species distribution, and geographic constraints and environmental heterogeneity affect genetic distribution.

4.2. Climate Change Facilitates the Contact of Two Closely Related Maples with Divergent Grinnellian Niches

According to PNC, the distribution range of the most suitable habitat varies temporally with climate change. ENM in this study indicated that A. caudatifolium upward shifted to mountain ranges after the LGM (Figure 2), in contrast to the downward expansion of A. morrisonense from high mountains (Figure 5 of [7]). Among all maples in Taiwan, these two species have the closest phylogenetic relationship but are not sister species [12,13,14]. It can be inferred that phylogenetic niche divergence between these two species caused their respective ancestors to occupy different territories when entering Taiwan. Subsequent climate warming might facilitate their distributional contact, increasing competitive interaction, although such competition is likely to be modest [33]. However, further study is required to establish the direction of the causal relationship, that is, whether niche divergence promotes distributional contact leading to competitive pressure [33] or the competition caused by distributional contact leads to niche divergence [34].
Like A. morrisonense, the chlorotypes of A. caudatifolium are differentiated from north to south (Figure 2b of [7]). However, the north-south differentiation of A. caudatifolium resulted from the collection of different ancient populations (Figure 2), whereas divergent alpine refugia were responsible for the differentiation of A. morrisonense (Figure 5 of [7]). During the LGM, the west side of the southern Central Mountain Range of Taiwan were most suitable for A. caudatifolium, and northeastern Taiwan was another appropriate range separated from the southwest by the Central Mountain Range. The ancient northeastern populations harbored diverse chlorotypes, whereas the ancient southern populations had low genetic diversity. Under climate warming, these temperate-origin maples shifted to higher altitudes accompanied by range expansion. Therefore, the extant northern populations have chlorotypes from both the ancient southwestern and northeastern populations. By contrast, the southern populations only retain the original southwestern chlorotype due to geographic barriers.

4.3. Environmentally Biased Dispersal Constrains the Northernmost Populations

The distinct chlorotype composition suggests the northernmost populations as a unique evolutionary significant unit from the southern populations. The formation of the northeast-southwest geographic barrier may be attributable to environmentally based dispersal bias. The disruption of gene flow by terrain according to the mountain-barrier effect is more significant in the northern populations, particularly for blocking gene flow between the southern and northernmost populations on and around the Datun volcano group (i.e., populations YMS and RF). Environmentally determined migratory characters, e.g., seed dispersibility, may constrain the extent of range shift. Although the seeds of A. caudatifolium have wings that aid spreading, maple seed dispersal is unexpectedly distance-limited and aggregated [35]. The warming-driven downslope range shift is also more limited than the upshift [36], resulting in the limited passage of mountainous maples through the lower basins and lowlands and interrupting gene flow between the Datun volcano group and the southern populations.
Northeastern Taiwan is affected by the cold and humid winter northeast monsoon from October to April of the next year. The local climate difference (e.g., prec10 and Tmax2) between the rainy northeast and arid southwest in Taiwan has affected forest phenology [37,38]. Xylem embolism of maples under water deficiency may cause leaf senescence [39], facilitating leaf color change and defoliation. For deciduous trees, the duration of green leaves affects nutrient accumulation and reproductive yield in the following year [40,41]. The leaves of A. caudatifolium typically turn yellow in October, and differences in late autumn rainfall between northeastern and southwestern Taiwan may result in differences in leaf duration and next-year seed yield (i.e., fitness). The significant correlation between prec10 and the timing of range change implies that autumn rainfall was one of the key factors affecting adaptive migration under paleoclimate change. In addition, since selection may eliminate maladaptive migrants, such a phenological difference may affect the success of colonization, thereby interrupting north-south gene flow. The greatest contributing factor of ENM, Tmax2, are highly correlated with other temperature-related factors (Bio1, 6, 9, 11, Tavg1~5, 11, 12, Tmax1~5, 10~12, Tmin1~3, 12, vapr1~5, 10~12, AET, and PET, Table S10). These climate factors are mostly related to the winter temperature, highlighting the impact of the temperature in the winter northeast monsoon period on the current distribution of A. caudatifolium and on the maintenance of genetic differences between southern and northern populations.

4.4. Local Climate Heterogeneity Underlies Genetic Draft of Chlorotype Distribution

Our results show that local climate heterogeneity combined with the influence of geography (i.e., terrain and distance) reduces the ability of the northernmost chlorotypes to enter the south, resulting in genetic differentiation from the south. Prec10, srad6, srad7, and AET were related to Tajima’s D, indicating that autumn rainfall, summer sunshine, and evapotranspiration jointly determine the extent of population renewal and demographic structure. Precipitation, solar radiation, and evapotranspiration may affect photosynthesis, nutrient accumulation, growth, and yield simultaneously. These environmental factors related to survival and reproduction are undoubtedly crucial to population genetic structure via local adaptation and affect the chlorotype distribution via the genetic draft, particularly for populations in central to northern Taiwan.
Seasonality is more pronounced in northern Taiwan than in southern Taiwan. The higher genetic diversity of the northern populations promotes their adaptability and enables local adaptation to long-term seasonal and climate changes. By contrast, the southern margin may be more susceptible to climate change [42], as reflected in the reduced distributions after the warming in the LGM. Warming can even reduce the competitiveness of organisms [43], leading to more severe genetic losses in the south. Therefore, in the face of global warming, the potential of the lack of genetic variation in the south to weaken resilience deserves further attention.

Supplementary Materials

The following are available online at https://www.mdpi.com/article/10.3390/plants10081646/s1, Table S1: Variance inflation factors (VIFs) of the retained environmental variables; Table S2: Environmental factors of the 19 populations of Acer caudatifolium used in this study (csv file); Table S3: Population Pairwise genetic distances FST/(1-FST) matrix (csv file); Table S4: Population pairwise geographic distances matrix (csv file); Table S5: Population pairwise altitudinal distances matrix (csv file); Table S6: Population pairwise environmental differences matrix (csv file); Table S7: Details of the estimated indices of mismatch analysis; Table S8: Likelihood ratio test comparing the best two models in MLPE; Table S9: The optimal model used for generalized linear model (GLM) testing of the factors affecting demography; Table S10: Pearson’s correlation of environmental factors (csv file); Figure S1: Correlations between geographic distance and genetic differentiation according to the Mantel test and partial Mantel test conditioned on environmental differences.

Author Contributions

Conceptualization, P.-C.L.; data curation, M.-X.L. and H.-P.L.; formal analysis, M.-X.L., H.-P.L., M.-W.C., J.-T.C., and P.-C.L.; Funding acquisition, P.-C.L.; writing—original draft, P.-C.L.; writing—review and editing, M.-X.L., M.-W.C., and J.-T.C. All authors have read and agreed to the published version of the manuscript.

Funding

We are grateful to the staff of the National Center for Genome Medicine of the National Core Facility Program for Biotechnology, Ministry of Science and Technology (MOST) of Taiwan for technical and bioinformatics support. This research was financially supported by MOST (109-2621-B-003-003-MY3 and 109-2628-B-003-001) and was also subsidized by National Taiwan Normal University (NTNU).

Data Availability Statement

The cpDNA sequences are available on NCBI GenBank (accession numbers: trnH-psbA: MZ275974–MZ276267 rpl16: MZ275679–MZ275972).

Acknowledgments

We thank Bing-Hong Huang, Huan-Yi Hsiung, and Chien-Ti Chao for their assistance in field sampling and Dawn Schmidt for English proofreading and editing.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Grinnell, J. The niche-relationships of the California Thrasher. Auk 1917, 34, 427–433. [Google Scholar] [CrossRef]
  2. Fordham, D.A.; Resit Akçakaya, H.; Araújo, M.B.; Elith, J.; Keith, D.A.; Pearson, R.; Auld, T.D.; Mellin, C.; Morgan, J.W.; Regan, T.J.; et al. Plant extinction risk under climate change: Are forecast range shifts alone a good indicator of species vulnerability to global warming? Glob. Chang. Biol. 2012, 18, 1357–1371. [Google Scholar] [CrossRef]
  3. Nathan, R.; Muller-Landau, H.C. Spatial patterns of seed dispersal, their determinants and consequences for recruitment. Trends Ecol. Evol. 2000, 15, 278–285. [Google Scholar] [CrossRef]
  4. Yang, Y.-Z.; Zhang, R.; Gao, R.-H.; Chai, M.-W.; Luo, M.-X.; Huang, B.-H.; Liao, P.-C. Heterocarpy diversifies diaspore propagation of the desert shrub Ammopiptanthus mongolicus. Plant Species Biol. 2021, 36, 198–207. [Google Scholar] [CrossRef]
  5. Gamelon, M.; Grøtan, V.; Nilsson, A.L.K.; Engen, S.; Hurrell, J.W.; Jerstad, K.; Phillips, A.S.; Røstad, O.W.; Slagsvold, T.; Walseng, B.; et al. Interactions between demography and environmental effects are important determinants of population dynamics. Sci. Adv. 2017, 3, e1602298. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  6. Evans, L.M.; Allan, G.J.; DiFazio, S.P.; Slavov, G.T.; Wilder, J.A.; Floate, K.D.; Rood, S.B.; Whitham, T.G. Geographical barriers and climate influence demographic history in narrowleaf cottonwoods. Heredity 2015, 114, 387–396. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  7. Chang, J.-T.; Luo, M.-X.; Lu, H.-P.; Tseng, Y.-T.; Liao, P.-C. Population genetics under the Massenerhebung effect: The influence of topography on the demography of Acer morrisonense. J. Biogeogr. 2021, 48, 1773–1787. [Google Scholar] [CrossRef]
  8. Pironon, S.; Papuga, G.; Villellas, J.; Angert, A.L.; García, M.B.; Thompson, J.D. Geographic variation in genetic and demographic performance: New insights from an old biogeographical paradigm. Biol. Rev. 2017, 92, 1877–1909. [Google Scholar] [CrossRef] [PubMed]
  9. Sexton, J.P.; McIntyre, P.J.; Angert, A.L.; Rice, K.J. Evolution and Ecology of Species Range Limits. Annu. Rev. Ecol. Evol. Syst. 2009, 40, 415–436. [Google Scholar] [CrossRef] [Green Version]
  10. Angert, A.L.; Bontrager, M.G.; Ågren, J. What Do We Really Know About Adaptation at Range Edges? Annu. Rev. Ecol. Evol. Syst. 2020, 51, 341–361. [Google Scholar] [CrossRef]
  11. Li, H.-L.; Lo, H.-C. Aceraceae. In Flora of Taiwan, 2nd ed.; Taiwan, Editorial Committee of the Flora of Taiwan, Ed.; National Science Council of the Republic of China: Taipei, Taiwan, 1993; Volume 3, pp. 589–598. [Google Scholar]
  12. Zhang, Z.; Li, C.; Li, J. Conflicting phylogenies of Section Macrantha (Acer, Aceroideae, Sapindaceae) based on chloroplast and nuclear DNA. Syst. Bot. 2010, 35, 801–810. [Google Scholar] [CrossRef]
  13. Gao, J.; Liao, P.-C.; Huang, B.-H.; Yu, T.; Zhang, Y.-Y.; Li, J.-Q. Historical biogeography of Acer L. (Sapindaceae): Genetic evidence for Out-of-Asia hypothesis with multiple dispersals to North America and Europe. Sci. Rep. 2020, 10, 21178. [Google Scholar] [CrossRef] [PubMed]
  14. Areces-Berazain, F.; Hinsinger, D.D.; Strijk, J.S. Genome-wide supermatrix analyses of maples (Acer, Sapindaceae) reveal recurring inter-continental migration, mass extinction, and rapid lineage divergence. Genomics 2021, 113, 681–692. [Google Scholar] [CrossRef]
  15. Naimi, B.; Hamm, N.A.; Groen, T.A.; Skidmore, A.K.; Toxopeus, A.G. Where is positional uncertainty a problem for species distribution modelling? Ecography 2014, 37, 191–203. [Google Scholar] [CrossRef]
  16. Phillips, S.J.; Dudík, M. Modeling of species distributions with Maxent: New extensions and a comprehensive evaluation. Ecography 2008, 31, 161–175. [Google Scholar] [CrossRef]
  17. Doyle, J.J.; Doyle, J.L. A rapid DNA isolation procedure for small quantities of fresh leaf tissue. Phytochem. Bull. 1987, 19, 11–15. [Google Scholar]
  18. Hall, T. BioEdit: A user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symp. Ser. 1999, 41, 95–98. [Google Scholar]
  19. Excoffier, L.; Lischer, H.E. Arlequin suite ver 3.5: A new series of programs to perform population genetics analyses under Linux and Windows. Mol. Ecol. Resour. 2010, 10, 564–567. [Google Scholar] [CrossRef]
  20. Grant, W.S. Problems and cautions with sequence mismatch analysis and Bayesian skyline plots to infer historical demography. J. Hered. 2015, 106, 333–346. [Google Scholar] [CrossRef] [Green Version]
  21. Jombart, T.J.B. adegenet: A R package for the multivariate analysis of genetic markers. Bioinformatics 2008, 24, 1403–1405. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  22. Bates, D.; Mächler, M.; Bolker, B.; Walker, S. Fitting linear mixed-effects models using lme4. J. Stat. Softw. 2015, 67, 1–48. [Google Scholar] [CrossRef]
  23. Rousset, F. Genetic differentiation and estimation of gene flow from F-statistics under isolation by distance. Genetics 1997, 145, 1219–1228. [Google Scholar] [CrossRef] [PubMed]
  24. R Core Team and Contributors Worldwide. The R Stats Package ver. 3.5.0; R Foundation for Statistical Computing: Vienna, Austria, 2017. [Google Scholar]
  25. Ripley, B.; Venables, B.; Bates, D.M.; Hornik, K.; Gebhardt, A.; Firth, D.; Ripley, M.B. Package ‘mass’. Cran R 2013, 538, 113–120. [Google Scholar]
  26. Neher, R.A. Genetic draft, selective interference, and population genetics of rapid adaptation. Annu. Rev. Ecol. Evol. Syst. 2013, 44, 195–215. [Google Scholar] [CrossRef] [Green Version]
  27. Orsini, L.; Vanoverbeke, J.; Swillen, I.; Mergeay, J.; De Meester, L. Drivers of population genetic differentiation in the wild: Isolation by dispersal limitation, isolation by adaptation and isolation by colonization. Mol. Ecol. 2013, 22, 5983–5999. [Google Scholar] [CrossRef] [PubMed]
  28. Huang, B.-H.; Huang, C.-W.; Huang, C.-L.; Liao, P.-C. Continuation of the genetic divergence of ecological speciation by spatial environmental heterogeneity in island endemic plants. Sci. Rep. 2017, 7, 5465. [Google Scholar] [CrossRef] [Green Version]
  29. Hsiung, H.-Y.; Huang, B.-H.; Chang, J.-T.; Huang, Y.-M.; Huang, C.-W.; Liao, P.-C. Local climate heterogeneity shapes population genetic structure of two undifferentiated insular Scutellaria species. Front. Plant Sci. 2017, 8, 159. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  30. Huang, C.-L.; Chen, J.-H.; Chang, C.-T.; Chung, J.-D.; Liao, P.-C.; Wang, J.-C.; Hwang, S.-Y. Disentangling the effects of isolation-by-distance and isolation-by-environment on genetic differentiation among Rhododendron lineages in the subgenus Tsutsusi. Tree Genet. Genomes 2016, 12, 53. [Google Scholar] [CrossRef]
  31. Huang, C.-L.; Chen, J.-H.; Tsang, M.-H.; Chung, J.-D.; Chang, C.-T.; Hwang, S.-Y. Influences of environmental and spatial factors on genetic and epigenetic variations in Rhododendron oldhamii (Ericaceae). Tree Genet. Genomes 2014, 11, 823. [Google Scholar] [CrossRef]
  32. Shih, K.M.; Chang, C.T.; Chung, J.D.; Chiang, Y.C.; Hwang, S.Y. Adaptive genetic divergence despite significant isolation-by-distance in populations of Taiwan Cow-Tail Fir (Keteleeria davidiana var. formosana). Front. Plant Sci. 2018, 9, 92. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  33. Krosby, M.; Wilsey, C.B.; McGuire, J.L.; Duggan, J.M.; Nogeire, T.M.; Heinrichs, J.A.; Tewksbury, J.J.; Lawler, J.J. Climate-induced range overlap among closely related species. Nat. Clim. Chang. 2015, 5, 883–886. [Google Scholar] [CrossRef]
  34. Freeman, B.G. Competitive interactions upon secondary contact drive elevational divergence in tropical birds. Am. Nat. 2015, 186, 470–479. [Google Scholar] [CrossRef] [PubMed]
  35. Gómez-Aparicio, L.; Gómez, J.M.; Zamora, R. Spatiotemporal patterns of seed dispersal in a wind-dispersed Mediterranean tree (Acer opalus subsp. granatense): Implications for regeneration. Ecography 2007, 30, 13–22. [Google Scholar] [CrossRef]
  36. Mamantov, M.A.; Gibson-Reinemer, D.K.; Linck, E.B.; Sheldon, K.S. Climate-driven range shifts of montane species vary with elevation. Glob. Ecol. Biogeogr. 2021, 30, 784–794. [Google Scholar] [CrossRef]
  37. Chang-Yang, C.-H.; Lu, C.-L.; Sun, I.F.; Hsieh, C.-F. Flowering and fruiting patterns in a subtropical rain forest, Taiwan. Biotropica 2013, 45, 165–174. [Google Scholar] [CrossRef]
  38. Kuo, C.-C.; Bain, A.; Chiu, Y.-T.; Ho, Y.-C.; Chen, W.-H.; Chou, L.-S.; Tzeng, H.-Y. Topographic effect on the phenology of Ficus pedunculosa var. mearnsii (Mearns fig) in its northern boundary distribution, Taiwan. Sci. Rep. 2017, 7, 14699. [Google Scholar] [CrossRef] [Green Version]
  39. Sperry, J.S.; Donnelly, J.R.; Tyree, M.T. Seasonal occurrence of xylem embolism in sugar maple (Acer saccharum). Am. J. Bot. 1988, 75, 1212–1218. [Google Scholar] [CrossRef]
  40. Chabot, B.F.; Hicks, D.J. The Ecology of Leaf Life Spans. Annu. Rev. Ecol. Syst. 1982, 13, 229–259. [Google Scholar] [CrossRef]
  41. Montserrat-Martí, G.; Camarero, J.J.; Palacio, S.; Pérez-Rontomé, C.; Milla, R.; Albuixech, J.; Maestro, M. Summer-drought constrains the phenology and growth of two coexisting Mediterranean oaks with contrasting leaf habit: Implications for their persistence and reproduction. Trees 2009, 23, 787–799. [Google Scholar] [CrossRef] [Green Version]
  42. Franco, A.M.A.; Hill, J.K.; Kitschke, C.; Collingham, Y.C.; Roy, D.B.; Fox, R.; Huntley, B.; Thomas, C.D. Impacts of climate warming and habitat loss on extinctions at species’ low-latitude range boundaries. Glob. Chang. Biol. 2006, 12, 1545–1553. [Google Scholar] [CrossRef]
  43. Reich, P.B.; Sendall, K.M.; Rice, K.; Rich, R.L.; Stefanski, A.; Hobbie, S.E.; Montgomery, R.A. Geographic range predicts photosynthetic and growth response to warming in co-occurring tree species. Nat. Clim. Chang. 2015, 5, 148–152. [Google Scholar] [CrossRef]
Figure 1. Spatial-genetic distribution of Acer caudatifolium. (A) CpDNA haplotype distribution; (B) haplotype network; (CE) genetic barriers (blue lines) predicted by the Monmonier algorithm under (C) nearest neighbor (K = 6), (D) Delaunay triangulation, and (E) Gabriel graph. The black lines indicate the gene-flow paths allowed by the models.
Figure 1. Spatial-genetic distribution of Acer caudatifolium. (A) CpDNA haplotype distribution; (B) haplotype network; (CE) genetic barriers (blue lines) predicted by the Monmonier algorithm under (C) nearest neighbor (K = 6), (D) Delaunay triangulation, and (E) Gabriel graph. The black lines indicate the gene-flow paths allowed by the models.
Plants 10 01646 g001
Figure 2. Ecological niche modeling (ENM) of the potential distribution of Acer caudatifolium. (A) Current spatial distribution; (B) percentage contribution of the predicted climatic factors used for current ENM; (C) potential distribution projected in the middle Holocene (Holocene Thermal Optimal, HTO); (D) potential distribution projected in the Last Glacial Maximum (LGM).
Figure 2. Ecological niche modeling (ENM) of the potential distribution of Acer caudatifolium. (A) Current spatial distribution; (B) percentage contribution of the predicted climatic factors used for current ENM; (C) potential distribution projected in the middle Holocene (Holocene Thermal Optimal, HTO); (D) potential distribution projected in the Last Glacial Maximum (LGM).
Plants 10 01646 g002
Table 1. Sampling site information and summary statistics of genetic diversity.
Table 1. Sampling site information and summary statistics of genetic diversity.
PopLongitude (E)Latitude (N)Altitude (m)NPolymavgFSTΘs *SD(θs) *Θπ *SD(θπ) *Tajima’s DPFu’s FsP
YMS121.53925.181762–10212920.7240.3230.2390.3650.3380.2730.6880.3260.527
RF121.78725.065381–484910.4320.2330.2330.1410.207−1.0880.198−0.2630.190
LLS121.40124.6871166–1317700.4620N.A.0N.A.0N.A.0N.A.
JS121.2724.6671197–14681780.2231.4990.7211.4020.910−0.2290.461−0.2950.454
MC121.48624.6281106–12182110.4360.1760.1760.3260.3191.5050.9501.4740.733
SKR121.14324.561546–19752000.3000N.A.0N.A.0N.A.0N.A.
TPS121.52124.5231566–17641410.2380.1990.1990.2780.2970.8420.8500.9440.572
SY121.31224.3381848–20001080.5091.7910.9241.3370.919−1.0940.1760.7130.637
DXS120.97724.2361680–20262310.1700.1720.1720.1900.2290.1860.7690.6120.438
TRK121.40724.1922336–2987520.3380.6080.4800.5070.505−0.9730.189−0.8290.106
MF121.17624.0932087–22832680.3291.3280.6110.4350.380−2.1150.0020.6100.586
DD121.16923.7872196–2396400.2510N.A.0N.A.0N.A.0N.A.
RL120.92223.7081362–16421580.3671.5580.7611.4230.930−0.3180.4363.0810.935
TTC120.91323.531604–23881700.2890N.A.0N.A.0N.A.0N.A.
ALS120.85523.4821930–24041500.8860N.A.0N.A.0N.A.0N.A.
LD120.99523.2452033–23092600.3390N.A.0N.A.0N.A.0N.A.
TJ120.74123.061467–1490200.3480N.A.0N.A.0N.A.0N.A.
JBS120.75722.7261306–20402100.2640N.A.0N.A.0N.A.0N.A.
JSY120.7522.412571310.3520.2040.2040.0970.163−1.1490.169−0.5370.128
N, sample size; Polym, number of polymorphic sites (including indels). * These values have been multiplied by 1000.
Table 2. Summary results of model selection in MLPE. The models are listed in order of AIC value; the IBD model has the best performance.
Table 2. Summary results of model selection in MLPE. The models are listed in order of AIC value; the IBD model has the best performance.
ModelFormulanparAICBIClogLikDeviance
IBDgen~geo + (1 | pop)41103.31115.8−547.6301095.3
IBD + IBEgen~geo + env + (1 | pop)51103.71119.4−546.8601093.7
IBD + IBAltgen~geo + alt + (1 | pop)51105.11120.8−547.5401095.1
IBD + IBE + IBAltgen~geo + env + alt + (1 | pop)61105.71124.5−546.8501093.7
IBEgen~env + (1 | pop) 41110.71123.3−551.3501102.7
IBAltgen~alt + (1 | pop)41111.51124.0−551.7401103.5
IBE + IBAltgen~env + alt + (1 | pop)51112.61128.3−551.3001102.6
npar, number of parameters; gen, genetic distance; geo, geographic distance; env, environmental difference; alt, altitudinal difference.
Table 3. Summary results of mismatch analysis under the demographic and spatial expansion models. Standard deviation (SD) and 95% confidence intervals (95% CI) are displayed in Table S7.
Table 3. Summary results of mismatch analysis under the demographic and spatial expansion models. Standard deviation (SD) and 95% confidence intervals (95% CI) are displayed in Table S7.
Demographic ExpansionSpatial Expansion
Popτθ0θ1SSDRagτθMSSDRag
YMS0.74409990.0260.1990.7420.0039990.0260.199
RF2.9300.9003.6000.3070.3580.2600.0089990.00050.358
JS8.35001.6420.0430.1127.0991.3810.310.0350.112
MC0.7890.0029990.0250.2650.7850.0069990.0250.265
TPS0.64309990.0120.2080.6450.0019990.0120.208
SY1.0310.0049990.0370.1297.8971.690.1590.0630.129
DXS2.9820.9003.6000.2380.2500.3870.0039990.0020.250
TRK1.03709990.0650.3501.0350.0039990.0650.350
MF300.2470.0070.4397.6460.2080.0790.0050.439
RL0.7250.0109990.0580.1667.9241.0090.3260.0570.166
JSY2.9650.4500.4500.0280.5030.1350.1102.7050.00020.503
LLS, SKR, DD, TTC, ALS, LD, TJ, and JBS are genetically monomorphic, and thus mismatch analysis could not be conducted.
Table 4. Summary results of the generalized linear model for testing significant factors in the extent of demographic dynamics (Tajima’s D) and spatial expansion time (τ).
Table 4. Summary results of the generalized linear model for testing significant factors in the extent of demographic dynamics (Tajima’s D) and spatial expansion time (τ).
Tajima’s DSpatial Expansion Time (τ)
EstimateSEtPEstimateSEtP
Lat2.8560.7413.8530.003 *−6.6263.613−1.8340.090
Long−2.6811.697−1.5800.14317.7108.3582.1190.054
AltNANANANA−0.0060.003−1.9550.072
prec100.0130.0043.2790.007 *−0.0510.017−3.0250.010 *
srad60.0050.0022.6910.021 *NANANANA
srad7−0.0040.001−3.2200.008 *NANANANA
AET−0.0150.005−2.9720.013 *−0.0350.028−1.2500.233
GAI−0.00010.00004−1.7750.104NANANANA
* p < 0.05.
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Luo, M.-X.; Lu, H.-P.; Chai, M.-W.; Chang, J.-T.; Liao, P.-C. Environmental Heterogeneity Leads to Spatial Differences in Genetic Diversity and Demographic Structure of Acer caudatifolium. Plants 2021, 10, 1646. https://doi.org/10.3390/plants10081646

AMA Style

Luo M-X, Lu H-P, Chai M-W, Chang J-T, Liao P-C. Environmental Heterogeneity Leads to Spatial Differences in Genetic Diversity and Demographic Structure of Acer caudatifolium. Plants. 2021; 10(8):1646. https://doi.org/10.3390/plants10081646

Chicago/Turabian Style

Luo, Min-Xin, Hsin-Pei Lu, Min-Wei Chai, Jui-Tse Chang, and Pei-Chun Liao. 2021. "Environmental Heterogeneity Leads to Spatial Differences in Genetic Diversity and Demographic Structure of Acer caudatifolium" Plants 10, no. 8: 1646. https://doi.org/10.3390/plants10081646

APA Style

Luo, M. -X., Lu, H. -P., Chai, M. -W., Chang, J. -T., & Liao, P. -C. (2021). Environmental Heterogeneity Leads to Spatial Differences in Genetic Diversity and Demographic Structure of Acer caudatifolium. Plants, 10(8), 1646. https://doi.org/10.3390/plants10081646

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop