Phylogenetic diversity of ferns reveals different patterns of niche conservatism and habitat filtering between epiphytic and terrestrial assemblages

: Much attention has been directed to understanding species richness patterns, but adding an evolutionary perspective allows us to also consider the historical processes determining current diversity patterns. We analyzed phylogenetic patterns of fern species assemblages in 868 plots along a wide range of elevational (0-4000 m) and latitudinal (0°-23°N) gradients in the Neotropics to allow a deeper understanding of evolutionary processes underlying current patterns of diversity and community assembly. Overall, we found that phylogenetic mean pairwise distance (sMPD) and mean nearest taxon distance (sMNTD) decreased with increasing latitude and elevation, but that these geographical factors per se were weak explanatory variables. Incorporating environmental variables strongly enhanced the power of the predictive model, indicating that fern assemblages are phylogenetically more diverse under wet and warm to cool conditions at low latitudes and elevations. Further, whereas epiphytic fern assemblages were strongly influenced by climatic factors, this was not the case for terrestrial ones, suggesting that edaphic conditions and vegetation structure may have a stronger influence on the evolution and diversification of terrestrial ferns. We conclude that fern assemblages are strongly influenced by phylogenetic niche conservatism and environmental filtering. This has also been found for angiosperms, but the direction of the environment-phylogenetic relationship is often opposed in the two groups, suggesting that the older age of many fern lineages includes historical signals that are not evident in the more recent angiosperm radiation. This paper is part of an Elevational Gradients and Mountain Biodiversity Special Issue. Abstract Much attention has been directed to understanding species richness patterns, but adding an evolutionary perspective allows us to also consider the historical processes determining current diversity patterns. We analyzed phylogenetic patterns of fern species assemblages in 868 plots along a wide range of elevational (0-4000 m) and latitudinal (0°-23°N) gradients in the Neotropics to allow a deeper understanding of evolutionary processes underlying current patterns of diversity and community assembly. Overall, we found that phylogenetic mean pairwise distance (sMPD) and mean nearest taxon distance (sMNTD) decreased with increasing latitude and elevation, but that these geographical factors per se were weak explanatory variables. Incorporating environmental variables strongly enhanced the power of the predictive model, indicating that fern assemblages are phylogenetically more diverse under wet and warm to cool conditions at low latitudes and elevations. Further, whereas epiphytic fern assemblages were strongly influenced by climatic factors, this was not the case for terrestrial ones, suggesting that edaphic conditions and vegetation structure may have a stronger influence on the evolution and diversification of terrestrial ferns. We conclude that fern assemblages are strongly influenced by phylogenetic niche conservatism and environmental filtering. This has also been found for angiosperms, but the direction of the environment-phylogenetic relationship is often opposed in the two groups, suggesting that the older age of many fern lineages includes historical signals that are not evident in the more recent angiosperm radiation.


Introduction
Biodiversity is distributed very unevenly around the globe and explaining geographical patterns of species richness remains one of biogeography's major challenges. An example of this unevenness is the latitudinal gradient of increasing species richness from the polar regions towards the tropics (Darwin 1859, Pianka 1966, von Humboldt [1828 2004, Mittelbach et al. 2007, Brown 2014. Despite over 200 years of research, the factors driving the observed patterns are still imperfectly understood (e.g., Hillebrand 2004, Qian and Ricklefs 2011, Oliveira and Scheffers 2019. Recently considered and partly interlinked explanations for such large-scale richness patterns appear to involve mainly climatic variables (e.g., Hawkins et al. 2003, Currie et al. 2004, Weigand et al. 2019, Abotsi et al. 2020, and historical and evolutionary processes (e.g., Ricklefs 2004, Wiens and Donoghue 2004, Fine 2015. Elevational gradients have increasingly gained attention for assessing the relative contribution of the plethora of proposed drivers of species richness patterns. For a long time, there was an unquestioned belief that elevational gradients simply mirror latitudinal gradients because both span a transition from warm to cold climatic conditions. However, Rahbek (1995) showed that this perception was the result of overemphasizing a few studies but that across nearly all taxonomic groups the majority of studies find unimodal richness patterns with maximum richness at some intermediate point of the gradient, observation that has been confirmed by many studies since then (e.g., Grytnes and Beaman 2006, Grau et al. 2007, Pandey et al. 2020. Elevational gradients further differ from the latitudinal gradient in that they are more concise (usually tens of km) and have many more replications (Rahbek 2005, McCain 2009), which allows for disentangling the effects of different ecological factors that vary in different ways along different elevational gradients (Lomolino 2001). Both latitudinal and elevational gradients thus represent large-scale biogeographical gradients along which biodiversity patterns change in clear and often predictable ways. Because of the different spatial extents and specific ecological transitions of the two types of gradients, a combination of latitudinal and elevational gradients offers invaluable opportunities to disentangle the relative roles of different drivers of species richness patterns (Qian et al. 2020).
There is increasing evidence that local diversity is strongly connected to regional diversity, raising the question of the historical component of local diversity patterns (e.g., Wiens and Donoghue 2004, Ricklefs 2005, Roy and Goldberg 2007, Weigand et al. 2019. Biotic assemblages along ecological gradients may thus show phylogenetic patterns reflecting historical differences between regions as well as inherited traits of clades contributing to the assemblages (Hawkins et al. 2005, 2006, Ricklefs 2005. For example, the tropics may harbor higher diversity across many taxa than temperate regions because the taxa have originated there ("out of the tropics", Jablonski et al. 2006) and had a greater extent in the past (Behrensmeyer et al. 1992, Kissling et al. 2012, so that most extant clades are originally tropical, leading to greater time and space availability for speciation. Moreover, a higher amount of energy in the tropics leads to shorter generation times and thus higher diversification rates ("speciation-extinction", Cardillo 1999, Mittelbach et al. 2007). Turning towards outer-tropical regions, adaptations are necessary to disperse and persist in cold and climatically seasonal regions, and these have evolved only in some taxa ("habitat filtering", Harper 1977, "niche conservatism", Wiens and Donoghue 2004), that represent clusters of closely related species that tend to be ecologically more similar and should therefore show higher competition (Darwin 1859, Cavender-Bares et al. 2009).
As species assemblages may show non-random patterns of species co-occurrences (Webb et al. 2002), the degree of relatedness between co-occurring species in terms of phylogenetic similarity may shed light on prevailing ecological processes such as environmental filtering (Weiher andKeddy 1995, de Bello et al. 2006) and niche differentiation (Stubbs and Bastow Wilson 2004). The processes behind these assumptions state that either the probability of a species to persist in an assemblage may increase if it is distinct from other species in the assemblage ("overdispersion": more efficient exploitation of the available resources to reduce competitive interactions, see Diáz andCabido 1997, Hooper et al. 2005), or that the environment sets limits on the occurrence of species unsuited to these ambient conditions and admit only species that are viable under the given environmental conditions (e.g., climate, disturbance regimes).
Ultimately, neither latitude nor elevation directly influence organisms. Rather, the actual drivers of the diversity patterns are ecological and historical factors that are associated with latitude and elevation (Körner 2007). For this reason, to truly understand how diversity patterns have developed, it is useful not only to analyze patterns against elevation and latitude, but also in relation to the covarying factors. The challenge here is identifying the crucial factors among the dozens of potential ecological drivers, whose influence is also likely to interact among each other and to differ in importance along different parts of the geographical gradients (Perrigo et al. 2020). Historical factors such as the timing of mountain uplift (Hoorn et al. 2013) are even more difficult to quantify. For this reason, a combination of geographical (latitude, elevation), ecological, and (if available) historical factors should be included to understand the origin of diversity patterns.
In order to elucidate large-scale richness patterns, we investigated an evolutionarily well-delimited plant group (ferns) within a region of high diversity along both latitudinal and elevational gradients. Ferns are the second-most diverse lineage of vascular plants on Earth, with a rather complex evolutionary history and an estimated age of 431 ma, which also includes major recent lineages with an approximate age of 40-80 ma Sundue 2016, Lehtonen et al. 2017 ferns are largely independent of biotic interactions, thus reducing the number of covarying factor as a potential explanation of the biodiversity patterns. For these reasons, they represent a suitable model group for the investigation of phylogenetic aspects of biogeographical questions. Ferns are also useful for understanding differences in evolutionary pathways between different life forms, in this case, terrestrial (ground-living) and epiphytic lifestyles. The epiphytic lifestyle differs from the terrestrial one in that due to a lack of direct contact to the soil, plants are more strongly limited by water and nutrient availability, while at the same time being more exposed to climatic extremes (Zotz 2016). Accordingly, numerous adaptations are necessary to conquer the epiphytic habitat, limiting the number of lineages that have become epiphytic. The ancestral condition of ferns is terrestrial, but epiphytism has evolved in six major lineages, in most cases alongside the rise of angiosperms to ecological dominance during the Cenozoic (Schuettpelz and Pryer 2009). Today, epiphytes make up about 27% of total fern diversity (Zotz et al. 2021), although locally they can contribute over half of all species (Kessler 2001). Because epiphytic ferns belong to a lower number of phylogenetic lineages, and because of the more stressful environmental conditions of their habitat, one might expect that they show different patterns of phylogenetic metrics than terrestrial ferns.
By evaluating the relationship between climate and phylogenetic metrics on latitudinal and elevational gradients, our aim was to provide insight into the mechanisms driving the assembly of species from regional pools into local assemblages (Qian et al. 2020, in press). The Neotropics offer a unique opportunity for the analysis of gradients of biodiversity because they contain one of the greatest centers of fern species diversity (Tryon and Tryon 1982), numerous mountain ranges, and because the mountain ranges are largely arranged in a north to south fashion.
We gathered a dataset of 868 study plots and 922 fern species using a uniform sampling protocol, spanning the latitudinal gradient from the equator to 23°N, between sea level and mountain tops up to 4000 m. Thus, the latitudinal gradient covered by our study runs from the inner tropics through the transition zone of subtropical regions, which allowed us to integrate seasonal and dry climates, which are known to limit the diversity of ferns as humidity-dependent organisms (Page 2002, Brodribb andMcAdam 2011); therefore we expected a response with respect to species diversity and the phylogenetic metrics.
In previous studies, we have shown that along our studied elevational gradients species richness (number of species per plot) of ferns follows the characteristic hump-shaped pattern , Salazar et al. 2015, Hernández-Rojas et al. 2020, and that overall richness of the gradients decreases from the equator northwards (Hernández-Rojas et al. 2018). In light of the phylogenetic hypotheses outlined above, for the present study we expected that the values of phylogenetic metrics of fern assemblages should decrease with increasing environmental stress. Within tropical humid mountain forests, where temperatures are moderate and humidity is consistently high, ferns do not experience physical but biotic interactions (Karger et al. 2014), whereas towards low and high elevations (heat and frost) and towards high latitudes (increasing frost, drought, and seasonality), physical limitations for fern growth and reproduction emerge, not only changing richness but also community composition and hence phylogenetic richness.
Specifically, we expected the next: H1: Phylogenetic metrics (MPD and MNTD) show a decrease (i) latitudinally from south to north (i.e., tropics to subtropics) and (ii) elevationally towards extreme ends of the elevational gradient; H2: These trends are closely related to the physical environment, especially amount and seasonality of precipitation; H3: These trends differ between terrestrial and epiphytic fern assemblages, with the latter showing stronger climate-related patterns.

Study Area
We studied 11 elevational gradients spanning Ecuador (0°) to northern Mexico (23°N; Fig. 1, Table 1). Elevational gradients extended from almost sea level to timberline (where present) and were selected to cover natural zonal forest with limited human disturbance patterns. With the exception of the transect in Honduras, the gradients have been described before, and further details may be found in previous publications: Costa Rica: Kluge and Kessler (2005)

Vegetation sampling
On each gradient, we sampled the fern assemblages at regular elevational intervals of 100-300 m (every 500 m in Ecuador and Perote), depending on the accessibility of zonal natural forest. At each elevation, up to 8 plots of 20 m x 20 m (400 m 2 ) were sampled with a consistent, standardized methodology (Kessler andBach 1999, Karger et al. 2014), avoiding special topographic and thus microclimatic changes like deep ravines and ridges. This is a small enough size to keep environmental factors and forest structure more or less homogeneous within the plots and is the minimum area required for representative pteridophyte surveys in humid tropical forests (Kessler and Bach 1999). Within each plot, all terrestrial and epiphytic ferns were recorded. The following climbing, Since we mostly did not climb the host trees, crown and high trunk epiphytes were recorded using binoculars and collecting poles, and by searching recently fallen trees and branches within the plot or adjacent locations (Gradstein et al. 2003, Sarmento-Cabral et al. 2015.

Phylogenetic reconstruction and molecular dating
To account for evolutionary history in our analyses, we used a supermatrix approach to generate a Phylogenetic diversity of ferns along gradients of latitude and elevation densely sampled phylogeny of ferns. We first used the sequence data of Testo and Sundue (2016) and then added sequence data for an additional 1068 taxa from Genbank and sequences newly generated for this study. The combined matrix included sequence data from seven chloroplast regions: atpΑ, atpΒ, matK, rbcL, rpl16, rps4, and trnL-trnF. Each region was aligned using the MAFFT (Katoh et al. 2002) plugin in Geneious v. 9.0.1 (Biomatters, Ltd.), and the resulting alignments were visually inspected and edited manually. Best-fitting models of nucleotide substitution were obtained for each region using jModelTest2 (Darriba et al. 2012) and the markers concatenated into a single matrix of 24572 base pairs. Several iterations of preliminary trees were generated using RAxML 8.2.10 (Stamatakis 2006) on the Cipres Science Gateway portal (Miller et al. 2010) to identify rogue and misidentified taxa. Following the removal of these problematic accessions, our matrix included 5075 taxa, including 27 outgroup-representing bryophytes, lycophytes, and seed plants, as well as 5048 ferns, representing approximately 50% of global fern diversity (PPG 1 2016; species number per family in phylogeny: see Supplementary Material Table S1). We generated a time-calibrated phylogeny using a multi-step approach. First, we generated a Maximum Likelihood (ML) tree using RAxML on the Cipres Science Gateway Portal (Miller et al. 2010) and retained the best tree as a constraint tree for dating analyses. We then used treePL (Smith and O'Meara 2012) to time-calibrate this best tree, using the fossil constraints employed by Testo and Sundue (2016) and a penalized-likelihood smoothing value of 0.001. Both our final ML tree and the time-calibrated tree are available (Supplementary Material Table S1a and S1b).
We measured the simple evolutionary relatedness within ecological communities or species assemblages using the most compressive available time-calibrated phylogeny that includes the sequence data from six chloroplast regions of 5057 tips of fern species (more than half of the extant ferns, Noben et al., unpubl. data). We excluded lycophytes because they are a sister group to the euphyllophytes (ferns and seed plants) and thus do not belong to the ferns (PPG I 2016). Only taxa determined to the species level, and some determined to genus (mainly for Ecuador, if they were the only species of the genus recorded), were included. Of the total of 922 species, 632 (68.5%) were represented in the phylogenetic tree (Supplementary  Material Table S2). For the 290 species (31.5%) not represented in the phylogeny, we searched for the putative closest relative in the phylogenetic tree (based on morphology) using it in the place of the focal species and by avoiding using the 632 species that had already been selected.
We analyzed three different data sets separately: all species, terrestrials, and epiphytes, because the two major fern life forms represent different ecological groups and thus react differently to environmental drivers in terms of richness and composition. We used presence-only data, because a few plots had no abundance information, to avoid personal biases in abundance estimation and because preliminary analyses showed that the phylogenetic patterns were similar using both presence-absence and abundance data. We used two different measures for calculating the phylogenetic structure of the assemblages: mean pairwise distance (MPD), which is an estimate of the average phylogenetic relatedness on basis of the branches between all possible pairs of taxa in a local community (plot assemblage; Kembel 2010), and MNTD or mean nearest taxon distance, which is the mean distance separating each species in the community from its closest relative. MPD is generally thought to be more sensitive to tree-wide patterns of phylogenetic clustering and evenness, whereas MNTD is more sensitive to patterns of evenness and clustering closer to the tips of the phylogeny (Kembel 2010, Cadotte and Davis 2016). Since these metrics of phylogenetic structure are not statistically independent of species richness, we used the standardized effect sizes (SES) of these metrics (sMPD and sMNTD), which compare empirical values against randomized null communities  to standardize for unequal richness across samples (Faith 1992; for constructing an appropriate Null model see next section). Positive SES values indicate phylogenetic overdispersion, or a greater phylogenetic distance among co-occurring species than expected by chance and thus high phylogenetic values. Negative SES values indicate phylogenetic clustering, or smaller phylogenetic distances among co-occurring species than expected by chance ).
All phylogenetic metrics and null model estimation were obtained by using the R package picante ).

Null model selection
Selecting an appropriate null model is crucial and requires knowing the model assumptions and choosing the appropriate species pool (Gotelli and Graves 1996, Gotelli 2000, Webb et al. 2011. Because all transects belong to a coherent biogeographical region, we analyzed all transects together by taking the whole species list as the basic species pool. We tested two procedures to gain 9999 randomly constructed assemblages using a fully constrained model that randomizes the community data matrix keeping richness and occurrence frequencies constant (constrained model: "trial swap algorithm", Miklós and Podani 2004, which is a modification of the "independent swap algorithm", Gotelli andGraves 1996, Gotelli 2000, but performs much faster using the same constraints). This model assumes that the ability of a species to colonize a plot is proportional to its frequency in the community matrix (Kembel and Hubbell 2006). We also used a model without any constraints in the randomization procedure, randomly shuffling the tips labels across the tips of the phylogeny (unconstrained model, "taxa labels"). Both models showed highly congruent results in the phylogenetic pattern, differing only in the magnitude (Supplementary Material Fig. S1). Hereafter, we therefore present only the results for the null assemblages constructed with the trialswap algorithm.

Explanatory variables
As environmental variables, we included energyrelated variables (temperature) and humidity-related variables since the study group is closely dependent on climatic factors related to humidity, because their sexual reproduction is linked to the presence of water (Page 2002) and because of their poor stomatal control (Kessler 2001, Brodribb andMcAdam 2011). Besides mean temperature and precipitation and their temporal variability and extreme values, we additionally included cloud cover as a variable because clouds reduce solar radiation and provide extra 'occult' precipitation (Bruijnzeel and Veeneklaas 1998). We extracted the following climatic variables per plot from the global climate data base set CHELSA ): Annual mean temperature and precipitation (Bio1, Bio12), as well as temperature and precipitation seasonality (Bio4, Bio15), the seasonal extremes maximum temperature of the warmest and coldest month, and the precipitation of the wettest and driest month (Bio5, Bio6, Bio13 and Bio14). Annual cloud cover and its seasonality were extracted from the 'EarthEnv'-dataset (Wilson and Jetz 2016). We checked for correlations between these variables to control for multicollinearity in models and excluded variables from the dataset that were highly correlated to already included variables (Pearson's R>75). As a result, we used the following variables to construct the full model: mean annual temperature TEMP, temperature seasonality TEMPs, annual precipitation PREC, and mean annual cloud cover CLOUD. Also, we ran models using latitude, elevation, and elevation to the square (Elevation + Elevation 2 ) because of the mostly humped nature of fern elevational gradients ).

Statistical procedures
Because spatially proximate plots often share many species and to account for differences in plot numbers between transects, for the analyses we combined the phylogenetic values of all plots in steps of 500 m and used these means for further analyses. Since our individual data points (plots) were not randomly distributed between the latitudinal and elevational limits, we used linear mixed effects models (LMMs) to control for the non-independence among data points in assessing changes in the phylogenetic variables with elevation, latitude, and in relation to climatic variables (fixed effects), because these models allow for spatial autocorrelation between neighbours (Crawley 2013). We used Likelihood Ratio Tests (LRT) or 'deviance tests' to compare between a null model without the term of interest and the model including this term to determine if one is a better fit to the data than the other (Luke 2017, Winter 2019. For the models including climatic variables or many fixed effects, we performed a full suite of likelihood ratio tests for all fixed effects in a model and constructed the correspondent comparison model providing p values for all fixed effects in a model (Singmann et al. 2020, Winter 2019). All variables used in the models were standardized (mean = 0 and SD = 1). To account for the non-independence of our data, the analyses were performed using the mean per elevational band every 500 m in elevation using the random structure "1|TRANSEC" indicating that the intercept is different for each transect, and "1" stands for the intercept here. In other words, there is going to be multiple responses per transect and these responses will depend on each transect baseline level. To decide whether such a simplified model was an enhancement to the previous model, we calculated the cAIC (conditional Akaike Information Criterion, Saefken and Ruegamer 2018), with a lower cAIC indicating a better model. Additionally, for each model, we calculated the amount of variation explained by the fixed effect only (marginal R 2 ) and the random effect (conditional R 2 ), which characterizes the variance described by fixed and random effects together (Winter 2019). All analyses were performed with the statistical platform R (R Core Team 2020) using the packages 'lme4' (Bates et al. 2015), 'afex' (Singmann et al. 2020), 'MuMIn' (Barton 2019), 'cAIC4' , 'MuMIn' (Nakagawa andSchielzeth 2013, Barton 2019), and 'ggtree' (Yu et al. 2017, Yu et al. 2018, Yu 2020 to construct the phylogentic trees showed here.

Results
Along the latitudinal gradient, species richness decreased with increasing distance from the tropics, with concentrated losses in particular phylogenetically old clades (Fig. 1, large branches; Supplementary Material Fig. S2). Typical tropical fern groups like tree ferns (Cyatheaceae) and filmy ferns (Hymenophyllaceae), and also families like Gleicheniaceae, Lomariopsidaceae, Nephrolepidaceae, Marratiaceae, and Schizaceae faded out, whereas phylogenetically young families like Dryopteridaceae and Polypodiaceae kept their dominance increasing latitude (Fig. 1, short branches).
The core tropical gradients in Ecuador and Costa Rica showed overdispersion of sMPD almost throughout their elevational ranges, with decreasing values towards high elevation sites around tree line (clustering, Fig. 2). Towards higher latitudes, overdispersion was increasingly limited to mid-elevations and clustering dominating in the northernmost transects. In contrast, values of sMNTD rarely revealed overdispersion and showed unspecific patterns. Separation for major life forms showed that this overall pattern was strongly driven by epiphytic species, whereas in terrestrials the pattern was less pronounced.
Model relations of sMPD and sMNTD values to spatial and climatic variables showed that when considering all species, latitude and elevation per se did not have a strong effect on the phylogenetic variables (Table 2). In terms of explained variance (marginal R 2 ), models of sMPD only including spatial variables (latitude, elevation) had R 2 m values of 0.13-0.17, whereas when we included climatic variables, values were much higher (R 2 m = 0.43) and showed positive relationships. sMNTD showed no patterns, neither with latitude or elevation, nor with climatic factors. Terrestrial assemblages only showed weak latitudinal and elevational patterns of sMPD (R 2 m = 0.09-0.13), but none in relation to climatic factors, whereas sMNTD showed no significant patterns at all. Finally, the patterns of epiphytic assemblages largely mirrored the overall patterns for sMPD; whereas for sMNTD, they showed additional, albeit weak, relationships with latitude and elevation (R 2 m = 0.13-0.20), as well as with climate (R 2 m = 0.10). Table 2. Results of the likelihood ratio tests by a phylogenetic variable with p-values for all fixed effects in the linear mixed models (coefficients), using the transects as the random factor for ferns along 11 elevational transects (0-23° N) in America. We calculated three different models, one with elevation only, one with latitude only, and one with climatic factors only, but only the models with delta = 0 are shown here (function dredge). TEMP and TEMPs: annual mean temperature and its seasonality; PREC: annual precipitation; CLOUD: mean annual cloud cover. *** P<0.001, **P<0.01, *P<0.05.

Discussion
Incorporating phylogenetic information in ecology is critical because it provides the most complete perspective of evolutionary processes in a context where the priority should be minimizing the loss of biodiversity with unique evolutionary qualities (Vane-Wright et al. 1991). With these first phylogenetic analyses of a large dataset from a latitudinally wide range of elevational transects in the Neotropics, we are able to detect relatively low, but clear and significant phylogenetic trends of the taxonomic compositions of ferns and their environmental correlates. These relationships do not only shed light on the environmental limits to the distributions and niches of this group of vascular plants. More importantly, they unravel the evolutionary processes and abilities of major life forms and taxonomic groups to originate and adapt to increasingly stressful environments.
An underlying consideration of our study is that neither of the two geographical gradient types, latitude and elevation, explain anything per se. Rather, it is their climatic and spatial structure which are relevant. While this is well-known in principle (Körner 2007, Qian et al. 2020), many studies still treat latitude and elevation as direct explanatory factors. In the following, we instead focus our interpretations more on the underlying ecological factors, also because these showed stronger relationships with the phylogenetic patterns.
The first clear pattern was that sMPD showed clear positive relationships with mean annual precipitation, temperature, and cloud cover, meaning that assemblages in hot, wet, and cloudy habitats tended to have higher relative phylogenetic complexity than those under cold, dry, and sunny conditions. Interestingly, this pattern was only seen for sMPD and not for sMNTD, indicating that it was driven by the presence or absence of deep phylogenetic lineages, rather than by more recent radiations within distinct lineages. Indeed, examinations of the distribution of individual families showed that this pattern is driven by early divergent families such as Cyatheaceae, Dicksoniaceae, Dennstaedtiaceae, Hymenophyllaceae, Lomariopsidaceae, Marratiaceae, Nephrolepidaceae, Plagiogyriaceae, Saccolomataceae, and Schizaceae ( Fig. 1), which all showed a more or less sudden decline in diversity when reaching subtropical climates in Mexico (~15°N) and dry and cold habitats in general. In contrast, evolutionary younger families such as Dryopteridaceae and Polypodiaceae gained dominance in colder and dryer habitats (Fig.1).
Overall, this pattern suggests a strong effect of phylogenetic niche conservatism (Jablonski et al. 2006, Jansson et al. 2013, in which early divergent fern lineages evolved under humid tropical conditions, which were globally prevalent over much of the time period of the original diversification of major fern lineages (Page 2002, Wiens and Donoghue 2004, Becker et al. 2016. The extant lineages from this original fern radiation clearly still favor such conditions, and their number decline towards colder and dryer habitats and thus has a strong effect on phylogenetic measures emphasizing deep phylogenetic splits such as MPD (Kembel 2010, Webb et al. 2011. The underlying climatic gradients in both directions tell the wellknown story of diversity gradients in ferns, which are poorly adapted to withstand extreme temperatures and drought because of passive stomatal control (Barrington 1993, Brodribb et al. 2009, Brodribb and McAdam 2011. This pattern of decreasing phylogenetic diversity with increasing latitude and elevation and towards dry and cold habitats has been previously documented in ferns at local  and regional scales (Qian et al. in press), as well as in other groups of organisms (e.g., hummingbirds: Graham et al. 2009;trees: Qian and Ricklefs 2011;microbes: Wang et al. 2012;herbs: Masante et al. 2019), and it is one of the major explanations for high tropical species diversity in these groups. However, there appear to be differences between study groups and geographical regions that may deserve closer examination. For instance, among forest trees in eastern North America, Qian et al. (2020) found that phylogenetic patterns were clearer along elevational than latitudinal gradients, which is not the case in our study, and that tip-weighted phylogenetic measures (MNTD) showed clearer patterns than base-weighted measures (MPD), which is the opposite of what we found. It is too early to derive generalizations from such a contrast, but unsurprisingly it would appear that ferns and angiosperm trees show different patterns of phylogenetic community structure reflecting their different evolutionary histories, with ferns including many more old lineages than angiosperms. The same conclusion was reached at a regional scale for China by Qian et al. (in press).
Consistent with our hypothesis number three, we found that terrestrial and epiphytic fern assemblages showed very different phylogenetic patterns. Whereas epiphytic assemblages largely mirrored the overall pattern; this was not the case for the terrestrial assemblages. We believe that these different patterns may at least partly be driven by the fact that epiphytic ferns only occur in six major lineages: the evolutionarily old Hymenophyllaceae and Psilotaceae, and the younger Pteridaceae (in the vittarioid ferns), Aspleniaceae, Dryopteridaceae, and Polypodiaceae (Schneider et al. 2004, Schuettpelz andPryer 2009), although single epiphytic species occur in many other families (Zotz et al. 2021). Of these, the Psilotaceae and vittarioids are almost entirely tropical, the Hymenophyllaceae is mainly so, and the Aspleniaceae, Dryopteridaceae, and Polypodiaceae are well represented in temperate regions, although mainly with non-epiphytic species. Thus, it is not surprising that a very clear pattern emerges: the phylogenetically early divergent lineages are tropical, and the phylogenetically younger lineages are widespread. As a result, the loss of many families, especially the Hymenophyllaceae (which is very abundant in wet tropical forests), towards arid and cold habitats, results in a strong phylogenetic signal that presumably drives the pattern of sMPD. Indeed, in the most arid and cold environments, only the family Polypodiaceae remains Frontiers of Biogeography 2021, 13.3, e50023 © the authors, CC-BY 4.0 license 10 as epiphytes (Krömer et al. 2005, Sylvester et al. 2014, Sundue et al. 2015. This young (in terms of ferns) family (90-30 my old, Schneider et al. 2004, Schuettpelz and Pryer 2009, Testo and Sundue 2016, has evolved numerous adaptations to both dry and cold conditions, including water storage, poikilohydry, and frost tolerance (Barrington 1993, Watkins et al. 2007, Sundue et al. 2015. Among terrestrial ferns, in contrast, the pattern is much less clear. Essentially, all fern families have terrestrial representatives, and many fern families are exclusively terrestrial. Thus, although some tropical fern families clearly avoid dry and cold habitats (e.g., Cyatheaceae, Lomariopsidaceae, Lygodiaceae, Marattiaceae, Nephrolepidaceae, Oleandraceae, Schizeaceae, Tectariaceae), there are enough other phylogenetically old families that do occur under such conditions (e.g., Anemiaceae, Equisetaceae, Ophioglossaceae), thus a wide range of phylogenetic lineages is represented in the majority of fern assemblages. As a result of this, we did not find a clear phylogenetic structure as seen for the epiphytic ferns.
An even more striking finding was that whereas among epiphytes sMPD showed a much stronger signal relative to climatic factors (R 2 m = 0.43) than relative to latitude (R 2 m = 0.13) or elevation (R 2 m = 0.17), the reverse was true for terrestrial assemblages, with low values for elevation (R 2 m = 0.13) and even lower values for climatic factors (temperature only; R 2 m = 0.09). This suggests that the community composition of epiphytic ferns is much more strongly determined by climatic factors than that of terrestrial ferns. This makes ecological sense since epiphytic plants are directly exposed to climatic factors and strongly depend on water availability and drought avoidance to be able to survive in the canopy habitat (Zotz 2016). Terrestrial ferns, in contrast, appear to be much less dependent on regional climatic factors, but rather are influenced by such factors as soil conditions or vegetation structure (Krömer et al. 2013). For instance, the distribution of many fern species is well-known to be closely linked to soil conditions (Tuomisto et al. 2002(Tuomisto et al. , 2014 and microtopography, which influences microclimate and light availability (Jones et al. 2011). Because these factors are unavailable at the scale of our study plots, we were unable to include them in our study, but we suspect that if such data were available, we would find that the latitudinal and elevational patterns of the phylogenetic structure of terrestrial ferns, as seen by us, are in fact driven by ecological factors that are still undisclosed.
We found that, overall, fern assemblages in the Neotropics show a pattern of phylogenetic niche conservatism, with phylogenetically more diverse assemblages under wet and warm to cool conditions at low latitudes and elevations. This suggests that environmental filtering (Donogue 2008, habitat filtering sensu Harper 1977) is the main process for structuring fern assemblages at the scale of our study. This is in accordance with findings that the global diversification of ferns was strongly influenced by extinction events driven by unfavorable climatic conditions in the past (Lehtonen et al. 2017), thus limiting many fern lineages to the remaining suitable habitats. However, whereas epiphytic fern assemblages appear to be strongly influenced by regional climatic factors, this is not the case for terrestrial ones, suggesting that edaphic conditions and vegetation structure may have a stronger influence on the evolution and diversification of terrestrial ferns. Thus, although the overall pattern of niche conservatism and environmental filtering is in accordance with that observed in many groups of organisms, the striking divergence of the actual factors determining this among epiphytic and terrestrial ferns opens exciting research opportunities for unravelling generalities of how community assembly processes can be driven in similar (or different) directions by different sets of environmental factors. This is particularly relevant in comparison with angiosperms, which appear to show opposite environment-phylogenetic metric relationships (Qian et al. 2020). This may well reflect the fact that many extant fern lineages are much older than the majority of angiosperm lineages, and thus include different historical signals (Qian et al. in press).

Data accessibility
The data on which this study was based will be available after publication in Dryad.

Supplementary Material
The following materials are available as part of the online article from https://escholarship.org/uc/fb Figure S1. Comparison of the standardized effect sizes of the mean pairwise distance (MPD) using two null models along 11 elevational transects in America. Figure S2.Detailed version of the phylogenetic trees of all fern species recorded along 11 elevational transects in America. Figure S3.Patterns of standardized effect sizes of the mean nearest taxon distance (sMNTD) for all fern species, terrestrial and epiphytic ones separately against latitude, elevation and climatic variables along 11 elevational gradients in America. Table S1a. Maximum Likelihood (ML) tree using RAxML on the Cipres Science Gateway Portal for achieving the best tree as a constraint tree for dating analyses. Table S1b. Time-calibration of the best tree (Table S1a) with treePL. Table S2. Species list including the closest relatives for the absent species in our phylogeny.