A consistent species richness-climate relationship for oaks across the Northern Hemisphere

Aim: Although the effects of climate on species richness are known, regional pro‐ cesses may lead to different species richness–climate relationships across continents resulting in species richness anomalies, especially for tropical groups. Phylogenetic niche conservatism may also influence species richness–climate relationships of dif‐ ferent lineages. Here, we tested whether regional effects also exist for temperate lineages using the genus Quercus . Northern Hemisphere. Time period: Present day. Major taxa studied: Quercus (Fagaceae). Methods: We used a dated phylogeny and distribution data for Quercus to evaluate its global species richness patterns and phylogenetic niche conservatism. To evaluate the consistency in species richness–climate relationships across continents of the genus Quercus as a whole and the temperate subgenus Quercus , we conducted analy‐ ses of covariance with continent as the categorical variable and climate variables as the covariate. We calibrated four widely used models using the global data or data from each continent separately and evaluated the predictive power of each model for different continents using the root mean squared error. Results: The relationships between species richness and climate were not signifi‐ cantly different among continents for both the genus Quercus as a whole and the subgenus Quercus . Unlike the models based on European data, those based on North American and eastern Asian data predicted both the global species richness and the richness in other continents. The species richness of a subtropical subgenus Cyclobalanopsis was most influenced by water availability, while that of a temperate subgenus Quercus was most influenced by environmental temperature. conclusions: In contrast to the subtropical subgenus , our results showed a consistent richness–climate relationship and absence of regional effects on species richness across continents for the temperate subgenus Quercus , that the species richness of temperate lineages, for example subgenus Quercus , may have reached equilibrium with the contemporary climate in the Northern Hemisphere.

Based on the strong correlation between species richness and environmental variables, ecologists have argued that contemporary species richness patterns are determined primarily by contemporary climate and habitat heterogeneity (Field, O'Brien, & Whittaker, 2005; O'Brien, Field, & Whittaker, 2000) and that there is a globally consistent species richness-climate relationship. For example, Currie and Paquin (1987) built a logistic model between species richness and actual annual evapotranspiration using the distributions of North American trees, and they demonstrated that this model could predict the tree species richness in Great Britain and Ireland. O'Brien and colleagues (Field et al., 2005;O'Brien, 1998O'Brien, ,2006 proposed the water-energy dynamic hypothesis to predict tree species richness in South Africa. Francis and Currie (2003) further evaluated this hypothesis and showed that it could explain over 87% of the variation in angiosperm family richness and is consistent globally and independent of evolutionary history and plant functional assemblages. Wang, Fang, Tang, and Lin (2011a) indicated that winter coldness was the primary driver of woody plant richness in eastern Asia because winter coldness exerted strong effects on species with tropical affinity and proposed a model including winter coldness and other variables representing water, temperature seasonality and habitat heterogeneity. This model predicted the tree species richness in North America with high accuracy. Although these models differ in the choice of variables, they all have high explanatory power in describing species richness variation and suggest a globally consistent correlation between species richness and climate.
Species richness anomalies occur when regions with similar environments show significantly different species richness (Latham & Ricklefs, 1993;Latham et al., 1994). A typical example is eastern Asia, which has a similar climate to that of eastern North America but has significantly higher species richness after accounting for environmental effects (Guo & Qian et al., 2017;Qian & Ricklefs, 2000;Ricklefs, Qian, & White, 2004). The mechanisms underlying the eastern Asian anomaly have been attributed to its distinct geological, climatic and evolutionary history (Qian et al., 2017;Ricklefs et al., 1999;Ricklefs & Renner, 1994;Xiang et al., 2004). For example, the uplift of the Qinghai-Tibetan Plateau and the resulting heterogeneous topography in south-west China and monsoon system in East Asia (Molnar, Boos, & Battisti, 2010) have provided ample opportunities for species in eastern Asia to diversify and survive the dramatic climatic fluctuations since the mid-Cenozoic (Wen, Zhang, Nie, Zhong, & Sun, 2014;Xing & Ree, 2017). Lineage speciation rates have been found to be higher in eastern Asia than in North America due to complex topography and long-term vegetation stability in the former (Hughes & Atchison, 2015;Qian & Ricklefs, 2000;Xiang et al., 2004). In contrast, the increased seasonality and aridity in North America during the Tertiary (i.e., the Eocene-Oligocene climate cooling) and the Quaternary may have led to the extinction of many tropical groups due to their phylogenetic niche conservatism constraining them to warm and wet climate (Latham et al., 1994;Prothero, 1994;Qian, 2002;Qian & Ricklefs, 2000). The large extent of ice sheets in North America during the Last Glacial Maximum (LGM) could also have lasting effects on the equilibrium between species richness and climate. The species richness anomaly in eastern Asia has been found primarily in old tropical groups (Latham et al., 1994;Qian et al., 2017), and whether it also exists for temperate groups remain controversial (but see Latham et al., 1994).
Here, we used the widely distributed temperate genus Quercus L. to test the consistency of the species richness-climate relationships and the species richness anomaly among continents after accounting for environmental effects. The genus Quercus is classified into two subgenera, the subtropical subgenus Cyclobalanopsis and temperate subgenus Quercus (Wu, Raven, & Hong, 1999 (Nixon, 1997;Wu et al., 1999). Phylogenetic studies suggest that the group Quercus, group Lobatae and group Protobalanus form a monophyletic clade, which is sister to the clade containing group Cerris, group Ilex and subgenus Cyclobalanopsis (Deng, Jiang, Hipp, Manos, & Hahn, 2018;Denk & Grimm, 2009Hipp et al., 2018;Manos & Stanford, 2001). Many Quercus species are important components of forest ecosystems in temperate regions in the Northern Hemisphere (Nixon, 1997;Wu et al., 1999). Previous studies on the species richness patterns of Quercus indicate that the species richness-climate relationship of this group is rooted in its evolutionary history due to niche conservatism (Xu, Wang, Rahbek, Lessard, & Fang, 2013) and that a water-energy interaction significantly influences Quercus species richness (Xu, Wang, Rahbek, Sanders, & Fang, 2016).

K E Y W O R D S
climate, eastern Asia, energy, historical effects, niche conservatism, North America, Quercus, regional effects, richness anomaly, water In the present study, we analysed the species richness patterns of the genus Quercus, especially of the temperate subgenus Quercus in different continents across the Northern Hemisphere. We assessed the effects of contemporary climate, climate change since the LGM and topography on richness patterns at continental and global scales as well as testing the consistency of the species richness-climate relationship across continents. We assessed the regional effects on species richness using continents as a variable to reflect the long isolation and independent geological history of each continent.

| Species distribution data
Quercus species distribution data were obtained from the global Quercus distribution database (Xu et al., 2013(Xu et al., , 2016, which includes over 20,000 distribution records for c. 400 Quercus species collected from published floras, checklists, online databases and herbarium specimens. The intraspecific taxa were merged to species level. Island species and hybrids [hybridization status according to Govaerts et al. (2011)] were excluded. The final dataset included distribution data for 376 species (Supporting Information Appendix S1). Among the 376 species, 91 belong to the subgenus Cyclobalanopsis and 285 belong to the subgenus Quercus (Denk & Grimm, 2009Menitsky, 2005;Nixon, 1997;Rangel, 2006;Rangel, Carlos, Zenteno, de Lourdes, & Enríquez, 2002;Wu et al., 1999). Of the 285 species in the subgenus Area of geographical units is not significantly correlated with Quercus species richness at continental or global scales (p > 0.05, Figure S1b in Supporting Information Appendix S2). Species richness was estimated as the total number of species in each geographical unit.

| Environmental data
To estimate the effects of contemporary climate, habitat heterogeneity and past climate change on Quercus species richness, we selected variables that have been widely used in previous studies and are easily interpretable in a biological context. Contemporary climate variables were grouped into three types representing environmental energy availability, water availability and climate seasonality (Table 1) MAT, MTCQ, MAP, altitude and LGM climate data were downloaded from the WorldClim database (www.worldclim.org) at a spatial resolution of 2.5 × 2.5 arc minutes. PET and AET with spatial resolutions of 30 arc seconds were obtained from the Consortium of International Agricultural Research Centers consortium for spatial information (http://www.cgiar-csi.org/). PET was calculated using the Hargreaves method, which requires less parameterization and reduces sensitivity to climatic input errors (Fisher, Whittaker, & Malhi, 2011;Trabucco, Zomer, Bossio, Straaten, & Verchot, 2008). WD was calculated as the difference between PET and AET.
PETmin was calculated as monthly minimum PET of a year. Rainfall was calculated as the sum of monthly precipitation values when mean monthly temperature was above 0°C (Francis & Currie, 2003).  Table S1 in Supporting Information Appendix S2. Correlation coefficients between variables are shown in Table S2 in Supporting Information Appendix S2.

| Phylogenetic niche conservatism of Quercus
To determine whether the species richness-climate relationship is mediated by climatic niche conservatism, as suggested Cyclobalanopsis in the consensus supertree were extracted from the North American tree, which covers more species and was dated using older fossils compared to the Eurasian tree.
Ancestral climatic niches were reconstructed for group Lobatae, group Quercus, group Cyclobalanopsis and the whole North American clade (group Lobatae, group Quercus and group Protabalnus) using the "ace" function in the R package ape (Paradis, Claude, & Strimmer, 2004). Blomberg's K (Blomberg, Garland, & Ives, 2003), as implemented by phyloSignal in the R package phylosignal (Keck, Rimet, Bouchez, & Franc, 2016), was employed to quantify the phylogenetic signal of each climatic variable because it is insensitive to sample size (Munkemuller et al., 2012). Blomberg's K was estimated using the subtrees for each Quercus group. Higher K-values indicate stronger phylogenetic niche conservatism (i.e., close relatives are more similar in their climatic niches than expected based on phylogenetic relationships). The phylogenetic signal of all climatic variables was compared among clades/groups. Because evolution mode could influence the conservatism of phenotypic traits, we fitted three evolution models to each climatic variable to investigate the mode of evolution of climatic niches (Blomberg et al., 2003) using the R package Geiger (Harmon, Weir, Brock, Glor, & Challenger, 2008). The three models were the (a) Brownian motion model (BM) (Felsenstein, 1973), (b) Ornstein-Uhlenbeck model (OU) (Butler & King, 2004) and (c) early-burst model (EB) (Harmon et al., 2010).
Climatic variables of each species were calculated as the mean values within their distribution ranges (Supporting Information Appendix S3). The mean values of mean annual temperature, mean annual precipitation, annual range of temperature and mean elevational range for each species were calculated and shown on the tree using the barplot.phylo4d function in the R package phylosignal (Keck et al., 2016).

| Richness-climate relationships among continents
To evaluate the consistency of species richness-climate relationships among North America, Europe and eastern Asia, we built both bivariate and multiple regression models for each continent separately and at global level. Bivariate regression models help detect the best single variable for explaining Quercus species richness on each continent. The variables that showed a significant effect on Quercus species richness in bivariate regressions at global level were further tested for their significant differences in slopes and intercepts among continents using analysis of variance (ANOVA). In particular, ANOVA was performed for each climatic variable, in which species richness was used as the response variable, continent (represented by eastern Asia, North America and Europe) as the categorical variable, and each climate variable as the covariate, respectively. We also included the interaction term of the climate variable and continent in the ANOVA models. The spatial autocorrelation in species richness could inflate type I error and, consequently, the significance levels of the statistical tests for regression and ANOVA models. Therefore, a modified t test was conducted in R following Dutilleul, Clifford, Richardson, and Hemon (1993) for significance testing of the models. We also conducted simultaneous autoregressive (SAR) models with spatially structured errors to account for spatial autocorrelation in our data (Kissling & Carl, 2008). A series of models with spatial weights at distances ranging from 400 to 2,500 km at 100-km intervals were conducted, and the one with the lowest spatial autocorrelation in the residuals and minimum Akaike information criterion (AIC) value was selected as the best model. Spatial autocorrelation in model residuals was estimated using Moran's I. Spatial statistics were performed with the "spdep" library in R (Bivand & Wong, 2018;R Core Team, 2018). SAR models showed consistent results (see Table S3 in Supporting Information Appendix S2) with the bivariate regression models.
We constructed separate multiple regression models based on our global Quercus dataset and the datasets for Europe, North America and eastern Asia. We selected four previously proposed models for the relationships between species richness and climate. The first three models were ( To test for consistency in the richness-climate relationship, we calibrated each of the six models at the global scale and separately for each continent and used them to predict the global species richness and that of the other continents. We then calculated the root mean squared error (RMSE) and adjusted r 2 for each prediction to assess the predictive power of the models. Lower RMSEs indicate higher predictive power and a higher level of consistency in the richness-climate relationship between the predicted continent and the continent used to build the model.

| RE SULTS
The six groups of oaks showed distinct distributions and species richness patterns in eastern Asia, Europe and North America ( Figure 1). The oak species richness in eastern Asia was dominated by the subtropical subgenus Cyclobalanopsis, but was dominated by temperate groups in North America and Europe. The total number of species in North America was higher than in eastern Asia. Ancestral climate reconstruction indicated that the ancestral MAT for the group Cyclobalanopsis was approximately 19°C, which was c. 3-5°C higher than that of the other groups ( Figure 2). Compared to the climatic niche of the subgenus Cyclobalanopsis, the climatic niches of the group Quercus and group Lobatae were significantly conserved with most K values larger than 1 (Table 2).
Simple regressions showed that water availability explained most of the variations in the species richness of the subtropical subgenus Cyclobalanopsis (Table 3). For the temperate subgenus Quercus, climate seasonality and past temperature change (TA) explained most variations in species richness in eastern Asia, while energy and habitat heterogeneity explained most variations in species richness in North America and Europe, respectively (Table 3). After accounting for the spatial autocorrelation in species richness, the slopes or intercepts of these simple regressions were not significantly different among continents for both the subgenus Quercus and the whole genus, which suggests that there is a consistent species richness-climate relationship among continents for oaks ( Figure 3).
Multiple-regression models showed that all models have a lower r 2 at the global scale than at the continental scale (Table 4). All four models showed a relatively high explanatory power in different continents (r 2 adj. from 0.571 to 0.854). Moreover, the r 2 of the four models was higher in eastern Asia (r 2 adj. ≥ 0.721) than in North America (r 2 adj. ≥ 0.571) and Europe (r 2 adj. ≥ 0.691). The Francis & Currie model (the second model) best predicted species richness of eastern Asia (r 2 adj. 0.822). The water-energy interaction model best predicted species richness in North America (r 2 adj. 0.799), and the water-energy dynamic model best predicted richness in Europe (r 2 adj. 0.854). However, for the subgenus Quercus, the r 2 of the four models (r 2 adj. ≥ 0.343) are lower than that of the whole genus. In eastern Asia, the r 2 of the four models for subgenus Quercus are higher than those for subgenus Cyclobalanopsis (Table S5 in Supporting Information Appendix S2).
Among the four models, the best global model was the Xu et al. (2016) model, which gave the highest explanatory power of species richness variation (whole genus, r 2 adj. = 0.710, p < 0.001, subgenus Quercus, r 2 adj. = 0.581, p < 0.001, Table 4 in maintext and Table S5 in Supporting Information Appendix S2) and could predict species richness in other continents with the lowest RMSEs compare to other three models ( Figure 4). However, when we used EU models

| Consistent oak richness-climate relationships across the Northern Hemisphere
Despite differences in the species composition of Quercus flora among continents, we did not find significant regional effects on the species richness of the entire genus and the temperate F I G U R E 2 Dated phylogenetic tree of the genus Quercus compiled from Hipp et al. (2000) and Deng et al. (2009). Species group membership, continental distribution, the mean annual temperature (MAT, °C), mean annual precipitation (MAP, mm), annual range of temperature (ART, °C) and range of altitude (RALT, m) are shown for each species. Range of altitude of a species was calculated as the mean of the altitudinal ranges of all geographical units in which the species occurs. The red figures besides the deep nodes represent the reconstructed ancestral temperature of these nodes [Colour figure can be viewed at wileyonlinelibrary.com] subgenus Quercus. This finding suggests that the eastern Asian species richness anomaly is not present in this group. Species anomalies between eastern Asia and eastern North America mainly exist for Tertiary relict lineages, which were extirpated from North America in the late Tertiary due to climate cooling and extensive glaciation but have survived in eastern Asia (Guo & Ricklefs, 2000;Qian et al., 2017;Qian & Ricklefs, 2000;Xiang et al., 2004). Temperate lineages, such as the subgenus Quercus, which diversified during this cooling period, have relatively higher species richness in North America than in eastern Asia. Other temperate genera (e.g., Carya) also have higher species richness in North America than in eastern Asia (Latham et al., 1994;Xiang et al., 2004). Our results, combined with previous findings, suggest that the species richness anomaly between North America and eastern Asia might not exist for temperate lineages.
Our results showed that Quercus species richness was a function of energy, water and their interaction (Table 4), which conforms with the tree species richness-climate models proposed in previous studies (Field et al., 2005;Francis & Currie, 2003;Hawkins et al., 2003;O'Brien, 1998;O'Brien et al., 2000;Whittaker, Nogues-Bravo, & Araujo, 2007;Xu et al., 2016). For example, in a series of classic studies, O'Brien and colleagues found that the species richness of woody plants in South Africa was a function of rainfall, PETmin and PETmin 2 due to higher evapotranspiration rates in warm regions (Field et al., 2005;O'Brien, 1998;O'Brien et al., 2000). However, we acknowledge that our study was limited to a typical temperate woody plant group in the Northern Hemisphere and therefore, further studies are needed to better understand the regional effects on species richness patterns of temperate clades. Although these studies support the importance of water, energy and their interaction on individual growth and distribution, the underlying physiological processes and the mechanistic link between individual growth and the general species richness pattern remain unknown.

| Phylogenetic niche conservatism and species richness-climate relationships
The current consistent Quercus richness-climate relationship across continents may be the result of the evolutionary history of this genus (Xu et al., 2013). According to the niche conservatism hypothesis, species tend to track their ancestral niches and have difficulties in adapting to new environments (Wiens et al., 2010). When climate deviates from its ancestral climate, the constraints on species richness induced by climate tend to be strengthened, and the species richness of a clade normally decreases. Moreover, groups with stronger niche conservatism have difficulties in adapting to new climates and may be more strongly limited by climate (Wiens & Donoghue, 2004). Recent studies on China's woody plants have found that winter temperature determines the latitudinal gradient of species richness by filtering out species with tropical affinities from the floras north of the tropics due to the tropical niche conservatism of these groups. The explanatory power of winter temperature on species with tropical affinity is higher than for those with temperate affinity (Wang et al., 2011b). Previous studies on Quercus have also found that the effect of environmental energy on the species richness of different groups depends on their ancestral niches; environmental energy limits species richness of tropical lineages at high latitudes but of temperate groups at both low and high latitudes (Xu et al., 2013). Here, we found stronger phylogenetic conservatism in the climatic niches represented by energy, water and seasonality variables for the temperate groups Quercus and Lobatae and the entire genus than for tropical group Cyclobalanopsis. Thus, niche conservatism may be one explanation for the globally consistent Quercus richness-climate relationship.
For a group with strong niche conservatism, adaptation to similar climate across its species tends to lead to a consistent species  Note. Subgenus Cyclobalanopsis was only distributed in Eastern Asia. Corrected p-values calculated using modified t tests (Dutilleul et al., 1993). Please see Table 1 for definitions of abbreviations. ***p-value < 0.01; **corrected p-value < 0.05; *p-value < 0.1.
subgenus Quercus, fossil records suggested its diversification in the Oligocene in North America when climate became cooler and dryer. Species from this subgenus might have dispersed from North America to Eurasia through the Bering and the potential Atlantic land bridges during glaciation periods (Denk, Grimsson, & Zetter, 2010. Its strong niche conservatism may have constrained the species of this group to similar climates across continents, leading to consistent richness-climate relationships in different continents. The high species richness of the groups Quercus and Lobatae is likely the result of the expansion of temperate forests and grasslands in North America due to continuous cooling and drying since the Oligocene (Axelrod, 1983). Dating results showed that the crown age of these two groups was c. 34 Ma, which corresponds to the Eocene-Oligocene climate cooling and drying. Fossil records also provide evidence that species from group Quercus and group Lobatae were widely distributed in high latitude regions of the Northern Hemisphere (Denk, Grimsson, & Zetter, 2010. When compared to Cyclobalanopsis, the species from these two groups show generally higher phylogenetic niche conservatism, which probably reflects their adaptations to seasonal climate, resulting in a strong species-climate relationship. Extinctions within a group would lead to low species richness of this group in some continents, and hence species diversity anomalies across continents. The earliest fossils with affinities to the group Cyclobalanopsis were found in North America in deposits from the Middle Eocene, when forests covered North America and the climate was warm and wet (Manchester, 1994). Subsequent cooling and drying during the Oligocene probably led to the regional extirpation of this group from North America (Axelrod, 1983). Currently, this group is confined to wet and warm subtropical and tropical eastern Asia where water availability is the main limiting factor for its species richness. It seems that adaptations to wet and warm climate may have give rise to the richness-climate relationship for subgenus Cyclobalanopsis. However, controversy exists between fossil and molecular evidence on the origination of the subgenus Cyclobalanopsis (e.g., Deng et al., 2018) and further studies are needed to solve the evolutionary history of this clade.
We did not detect significant phylogenetic signals for climatic niche variables in the subgenus Cyclobalanopsis. Currently, we have F I G U R E 3 Comparisons of the intercepts and slopes of species richness-climate regressions among continents for the genus Quercus and subgenus Quercus. Different letters on top of the bars indicate significant differences in slopes or intercepts between continents. After accounting for spatial autocorrelation, all slopes of richness-climate regressions were not significant among continents. In Europe, the intercepts of the regressions between richness and annual temperature range and range of mean annual precipitation are significantly lower than those in North America. MAT = mean annual temperature (°C); MTCQ = mean temperature of coldest quarter (°C); PET = annual potential evapotranspiration (mm); PETmin = minimum monthly potential evapotranspiration ( Note. All four models show a high explanatory power for all continents. Please see Table 1  North America, which has been regarded as one of the reasons for high species richness in eastern Asia (Latham et al., 1994;Qian & Ricklefs, 2000). In contrast, our study showed that the explanatory power of temperature anomalies since the LGM for species richness of the subgenus Quercus was higher in eastern Asia than in North America and Europe. One possible explanation of this pattern is that it was easier for species to track climate change through the southnorth mountain ranges in North America during the Quaternary climatic oscillations compared to the east-west mountain ranges in eastern Asia and Europe. Pollen fossil studies showed that Quercus and other dominant tree species have expanded their ranges following climate change since the LGM with a lag of <1,500 years (Jackson et al., 2000). This probably also explains the stronger species richness-climate relationship for the subgenus Quercus in North America than eastern Asia and suggests that species richness in North America is more likely to be influenced by recent climate change.

| Effects of habitat heterogeneity on oak species richness
The greater physiographical heterogeneity in eastern Asia compare to other parts of the world is often regarded as the main driver of the region's higher species richness because of the role of this heterogeneity in providing refuges and cradles for species evolution and offering greater habitat variability, enabling species coexistence (Andrews & O'Brien, 2000;Kerr & Packer, 1997;Tews et al., 2004). However, our results revealed that habitat heterogeneity represented by elevation and climatic range explains the species richness variation only in Europe but not in eastern Asia or North America. The relatively high isolation of the mountains in southern Europe from southern landmasses limits the southward distribution of species and increases the species pool, possibly contributing to high species richness in southern Europe and the significant effect of habitat heterogeneity on species richness in Europe (Svenning, Normand, & Skov, 2008).

| CON CLUS IONS
In contrast to the subtropical subgenus Cyclobalanopsis, our results support a consistent species richness-climate relationship for oaks in the temperate subgenus Quercus across the Northern Hemisphere, representing the effects of water, energy and their interaction. We did not find significant regional effects on species richness patterns of the temperate subgenus Quercus, which suggests that the species distribution and richness of temperate groups (e.g., subgenus Quercus) may have reached an equilibrium with the contemporary climate. The high level of niche conservatism and adaptation to cool and dry climate adopted by the temperate groups of Quercus may have played an important role in establishing its globally consistent richness-climate relationship across the Northern Hemisphere. Our results suggest that the distributions of groups with stronger niche conservatism are more strongly limited by climate, and hence are more likely to be affected by future climate change. F I G U R E 5 Latitudinal patterns of (a) the log-transformed species richness, (b) minimum monthly potential evapotranspiration (PETmin) and (c) rainfall. (d) The changes in species richness along the gradients of both minimum monthly potential evapotranspiration and rainfall. The size of dots is scaled by species richness [Colour figure can be viewed at wileyonlinelibrary.com]

AUTH O R CO NTR I B UTI O N S
Z.W., X.X. and C.R. designed the research; X.X. collected the data; X.X. and Z.W. performed the analyses; X.X. led the writing with contributions from all of the authors.

DATA ACCE SS I B I LIT Y
All data necessary to reproduce the analyses presented in this study are included in the Supporting Information (Appendices S1 to S4).