Evolutionary Pattern of High Variation Traits in Subfamily Zygophylloideae (Zygophyllaceae)

Patterns in traits and trait combinations reflect how organisms cope with their environment. Owing to different degrees of variability, the performance of traits varies during adaption to the changing environment. In this study, we focused on a taxon dominant in arid regions – the subfamily Zygophylloideae. We analyzed the evolutionary patterns of functional traits to clarify the impact of trait attributes on niche shifts. The results of phylogenetic signal analysis of traits revealed that quantitative traits, such as plant height, were not evolutionarily conserved. Phylogenetic regression pointed out that there are synergistic changes in environmental factors and in some traits within a phylogenetic context. These traits can meet the requirements of different environments more easily, possibly owing to their high variability. As a result, species in the subfamily Zygophylloideae showed clustering in some phenotypic spaces. Thus, the adaptive evolution of traits reduced niche restrictions of related environmental factors on species. The evolutionary analysis of functional traits in subfamily Zygophylloideae proved that high variability of traits allows fine-tuning according to the related environmental factors during the evolutionary process and promotes niche shifts.


Introduction
Distribution and abundance of species are governed by many environmental factors that together constitute the niche of a given species (Grinnell, 1917). Different niche factors restrict species in different aspects and select for independent trait dimensions (Laughlin, 2014). Thus, patterns among traits and trait combinations reflect how organisms cope with their environments (Winemiller et al., 2015). Many organismal traits have well-known functions, and their associated environmental characteristics are predictable (Segar et al., 2013). For example, many convergent phenotypes span diverse taxonomic groups, as stated by McGhee (2011). In most cases, an organism inhabits a suitable environment where a population can maintain a positive net rate of increase (Buckley et al., 2014). However, a series of environmental events, such as habitat and climate changes, invasion of exotic species, etc., can result in a population existing outside its niche envelope (Holt, 2009;Matsuzaki et al., 2016;Tingley et al., 2009). In such a situation, to avoid extinction, the population must migrate to other suitable habitats or its traits must evolve to adapt to the new environment (Brunet & Larson-Rabin, 2012;Chevin et al., 2010). Moreover, adaptive evolution is required to be faster than changes in the environment (Merilä & Hendry, 2014;Ning et al., 2010). Although the evolutionary direction of traits is subject to environmental factors, the rate of evolution is determined by its own properties.
Zygophylloideae is the largest subfamily of Zygophyllaceae, and it consists of about 180 species of shrubs, subshrubs, and herbs (Sheahan & Chase, 2000). The widespread genus Zygophyllum was found to be polyphyletic with other genera, such as Augea, Fagonia, and Tetraena, which were found embedded in it in phylogenetic analyses (Sheahan & Chase, 2000). A substantial revision of this subfamily was made by Beier et al. 2003, who suggested six monophyletic genera (Augea, Fagonia, redefined genera Tetraena and Zygophyllum, new genera Melocarpum and Roepera) that make up the Zygophylloideae. Species of Zygophylloideae are found in all major warm deserts of the world, except in South America, and exhibit different strategies for drought adaptation, hence representing an ideal subject to analyze the adaptability of morphological traits to different environments (Lauterbach et al., 2016).
Owing to different degrees of conservatism, the performance of different traits should be diverse during adaptation to environmental changes. Traits with higher variability are likely to undergo adaptive evolution in response to the environment. In this study, we focus on the impact of trait variability on niche shifts. Evolutionary patterns of functional traits of members of subfamily Zygophylloideae were analyzed using phylogenetic methods to verify the following two hypotheses: (i) selective pressure led to the emergence of a convergent phenotype that spanned diverse taxonomic groups; (ii) there are some phylogenetically correlated synergistic changes between climatic factors and biological traits. In addition, we also briefly discuss how associated limiting factors effectively shifted and restricted the ecological niche by influencing traits.

Material and Methods
In this study, a total of 75 individuals representing 75 species and six major clades of Zygophylloideae were used to establish a phylogenetic tree. Three individuals of another subfamily of Zygophyllaceae (Guaiacum angustifolium, Larrea tridentata, and Tribulus terrestris) were designated as outgroups. DNA sequences of three regions, including nuclear ribosomal internal transcribed spacer region (ITS), and plastid DNA (cpDNA) trnL and trnL-trnF, were combined and used for further analysis. Voucher information and GenBank accession numbers are provided in Table S1.
Bayesian inference (BI) and maximum likelihood (ML) methods were used to construct phylogenetic trees. The best fitting substitution model (general time reversible; GTR+Γ) was determined by the Bayesian information criterion (BIC) using jModelTest 2.1.10 (Darriba et al., 2012). ML analysis was implemented in MEGA version 10.0.5 (Kumar et al., 2018), with 1,000 bootstrap replicates. BI analysis was performed using MrBayes 3.2.5 (Ronquist et al., 2012). Markov chain Monte Carlo (MCMC) runs were performed with the following settings: sampling one tree every 1,000 generations for 10,000,000 generations; 1/4 of the samples were discarded as burn-in. Nodal robustness was estimated using posterior probability (PP). Tracer version 1.6 (Rambaut et al., 2014) was used to assess the stationarity of the runs.
We obtained morphological traits from various floral records, such as Flora of China (Wu et al., 2008), Flora of Australia (Wilson, 2013), collection of online floras network (http://www.efloras.org/index.aspx), and some additional literature (Beier, 2005;Sheahan, 2007). In total, one binary and five continuous characters were taken into consideration (Table S3). We downloaded the GPS latitude/longitude data for each species from the GBIF database (http://www.gbif.org/). The 30 arc-seconds resolution climate data was downloaded from WorldClim (http://www.worldclim.org/), and bioclimatic variables of species habitat were extracted using ArcGIS version 10.0.
We used D statistics of the phylo.d function in the caper package to calculate the phylogenetic signal of binary characters (Orme et al., 2012), and K statistics of the multiPhylosignal function in the picante package to calculate the phylogenetic signal of continuous characters (Kembel et al., 2010). The prcomp function was used to implement a principal components analysis of bioclimatic variables from all sample localities. Phylogenetic regression between principal components and biological traits was investigated using functions from the ape and nlme packages (Pinheiro et al., 2013). The continuous characters were ranked from large to small according to the values, and species in the top 20% and bottom 20% were selected as the focal species. They were then analyzed for convergent evolution using the Wheatsheaf index in the windex package (Arbuckle & Minter, 2015). All these analyses were performed on the phylogenetic framework obtained from the BI analysis.

Results
Phylogenetic relationships of Zygophylloideae obtained from BI and ML analyses of the three-marker dataset are shown in Figure S1 and Figure S2. The effective sample sizes (ESSs) for all traces were greater than 200. The results of BI and ML analyses were relatively consistent with previous research, with six monophyletic groups being identified within Zygophylloideae (Augea, Fagonia, Melocarpum, Roepera, Tetraena, and Zygophyllum). However, relationships between the six major clades were different from those obtained in previous studies (Beier et al., 2003;Wu et al., 2018).
The D statistics of the phylogenetic signal for the binary character "life form" (herbaceous or woody) was 0.16. This result indicated that life form evolution can be better fitted using the Brownian motion model compared to randomization of data, i.e., there was a certain degree of phylogenetic signal in the data. The opposite was observed for continuous data, where K statistics values lower than one in all cases indicated a higher level of homoplasy than expected from the Brownian model, i.e., there was little phylogenetic signal in the data (Table 1).
Phylogenetic generalized least squares (GLS) is a linear model in which the covariance (correlation) structure between species is permitted to match with that expected under a Brownian motion process of evolution on the tree. There were some relationships between the principal components of bioclimatic variables and plant traits, and slopes of some of these relationships were significant with p < 0.05 (Table 2). This implies that plant traits change in response to changes in environmental factors. Although evolutionary studies are a type of observation, the results of this study support the hypothesis that environmental selection pressure plays a role in the evolution of plant traits.
The Wheatsheaf index provides quantitative measures of the strength of convergence. This method penalizes phenotypic distances across species by the phylogenetic distance before investigating similarity so that close phenotypic similarity is weighted higher for distantly related species, then calculates a ratio of the mean (penalized) distances between all species to the mean (penalized) distances between focal species. The results of this analysis are shown in Table 3. The bounds of 95% confidence interval for some Wheatsheaf indices were much greater than 1, indicating that the focal species showed a tighter clustering to each other. The Wheatsheaf index detected stronger convergence in the analysis of the following   phenotypic spaces: lower plant height and larger sepal and petal. In other words, there were selective forces that resulted in the similarity of these traits between species.

Discussion
DNA-based species delimitation may be compromised by limited sampling effort and species rarity, including "singleton" representatives of species, which hampers estimates of intra-versus interspecies evolutionary processes (Ahrens et al., 2016). Zygophylloideae includes about 180 species, and this study sampled only 75 species (<50%). The small sample size resulted in a difference in the position of Augea capensis between the ML and Bayesian analyses, as shown in Figure S2 and Figure  S3. Nevertheless, ML and Bayesian analyses based on nucleotide sequences resulted in relatively consistent phylogenies, and bootstrap support values indicated that the resulting phylogeny was comparatively robust.
Since Darwin established that natural selection is an essential evolutionary mechanism, adaptive evolution of species has been the core issue of evolutionary biology, and has attracted the attention of evolutionary biologists (Berenbaum et al., 1996;Kocher, 2004). Past adaptations and differentiation of species determine their current distribution, richness, and functional traits. However, the traits (or genes) that vary in response to natural selection to allow adaptation remain to be determined. Therefore, consideration of the precise evolutionary process has become an essential part of ecological adaptation research (Feder et al., 2012;Jablonski, 2008). In the case of key innovations, the processes seem obvious. For example, nectar spurs of the genus Aquilegia, detoxification mechanisms in the butterflies of Pieridae, and hypocone of mammals open specific routes of diversification (Hunter & Jernvall, 1995;Wheat et al., 2007;Whittall & Hodges, 2007). These innovations are usually conservative, in the sense of being gained few times but being lost several times (Agrawal et al., 2009;Ree, 2005). Allometric changes follow different trajectories with easier transitions; however, they are potentially equally important for natural selection. For example, the positive relationship between mimetic fidelity and body size suggests that reduced predation pressure on less profitable prey species limits the selection for mimetic perfection (Penney et al., 2012).
In this study, phylogenetic signals of quantitative traits indicated that these traits cannot be very well explained by stochastic evolution according to the Brownian motion model. Phylogenetic GLS pointed out that environmental factors are associated with plant traits. For example, plant height and PC2 showed a significant phylogenetic negative correlation. Owing to selection pressure, species in subfamily Zygophylloideae exhibited clustering in some phenotypic spaces. For example, some species showed a tendency of convergent evolution towards smaller height.
The above results suggest that these traits have higher variability, allowing for faster adaptation to different environments. Thus, changes in these traits allow adaptive fine-tuning over short time scales and reduce niche restrictions of related environmental factors on species (Walsh et al., 2016).
Organisms are composed of different interrelated traits, which respond to various environmental factors. Adaptive evolution of traits maintains the ability of organisms to acquire outside resources in a changing environment and sustain themselves (Reich et al., 2003). In a stable environment, conservatism can reduce the generation of unfavorable traits and save resources. However, when the environment changes, conservative traits are less likely to adapt to the new environment, which can result in localized extinction of the population (Jump & Penuelas, 2005;Vincenzi, 2014). In this study, we demonstrated that high variability of traits allows fine-tuning according to the related environmental factors in the evolutionary process and promotes niche shifts. This research helps to understand the law of adaptive evolution and the adaptive changes necessary to allow survival of Zygophylloideae members during climate change.

Supporting Material
The following supporting material is available for this article: • Figure S1. Phylogenetic relationships within Zygophylloideae obtained from BI analysis based on the three-marker dataset. • Figure S2. Phylogenetic relationships within Zygophylloideae obtained from ML analysis based on the three-marker dataset. • Table S1. Voucher information and GenBank accession numbers of sequences used in Bayesian inference. • Table S2. Summary of principal component analysis. • Table S3. Morphological trait information and GBIF ID for each species used in this study. • Appendix S1. The R code used in this study.