Contrasting altitudinal variation of alpine plant communities along the Swedish mountains

Abstract Changes in abiotic factors along altitudinal and latitudinal gradients cause powerful environmental gradients. The topography of alpine areas generates environmental gradients over short distances, and alpine areas are expected to experience greater temperature increase compared to the global average. In this study, we investigate alpha, beta, and gamma diversity, as well as community structure, of vascular plant communities along altitudinal gradients at three latitudes in the Swedish mountains. Species richness and evenness decreased with altitude, but the patterns within the altitudinal gradient varied between sites, including a sudden decrease at high altitude, a monotonic decrease, and a unimodal pattern. However, we did not observe a decline in beta diversity with altitude at all sites, and plant communities at all sites were spatially nested according to some other factors than altitude, such as the availability of water or microtopographic position. Moreover, the observed diversity patterns did not follow the latitudinal gradient. We observed a spatial modularity according to altitude, which was consistent across sites. Our results suggest strong influences of site‐specific factors on plant community composition and that such factors partly may override effects from altitudinal and latitudinal environmental variation. Spatial variation of the observed vascular plant communities appears to have been caused by a combination of processes at multiple spatial scales.

However, despite general declines in diversity and species richness with altitude, the patterns with which species communities varies along altitudinal gradients are not uniform (Rahbek, 1995). Although the importance of spatial scale for altitudinal variation in species communities has been recognized (Rahbek, 2005), we still have limited understanding of what scales different processes act and interact in shaping altitudinal gradients in species diversity.
Alpine vegetation, that is, vegetation found above the tree line (Körner, 2000), is generally structured into three vegetation zones which have a vertical extent of a couple of hundred elevation meters, although their extent depends on the level of maritime influences on local climate (Körner, 2003;Rydin, Snoeijs, & Diekmann, 1999).
However, the altitude at which alpine vegetation occurs also varies with latitude. The highest tree line is found close to the equator, at approximately 4,800 m.a.s.l., from which the alpine zone is found at successively lower altitudes toward the poles (Hoch & Körner, 2005;Körner, 2003). Understanding how the spatial variation in community structuring of alpine vegetation varies with latitude is of particular relevance in the face of the ongoing climate change, since both alpine and northern environments are predicted to experience the biggest ecological impacts from global warming (Bekryaev, Polyakov, & Alexeev, 2010;Beniston et al., 1997;IPCC , 2013;Pepin et al., 2015).
Biodiversity is frequently divided into alpha, beta, and gamma components (Whittaker, 1960), where alpha diversity describes the diversity within locations, gamma diversity describes the diversity within landscapes, and beta diversity describes the difference in diversity between locations within landscapes. Most research on biodiversity variation along altitudinal gradients has either focused on variation in alpha diversity (e.g., Rahbek, 1995) or on variation in beta diversity along the gradients (which has been used to quantify the gradient itself, Anderson et al., 2011). Along altitudinal gradients, species richness (a measure of taxonomic alpha diversity) most frequently displays a hump-shaped pattern (Rahbek, 1995). This pattern has generally been attributed to productivity-associated competitive exclusion at low altitudes and a productivity-driven filtering of species communities at high altitudes (Huston, 1999).
Less attention has been given to spatial variation in communities across altitude (Rahbek, 2005). Although beta diversity is highly influenced by the spatial scale over which it is measured, the dilution of communities at high altitudes is expected to give rise to a spatial homogenization and a subsequent altitude-related decline in beta diversity (Naud et al., 2019). Such patterns have been observed, but do not appear to be uniform across all spatial scales (Naud et al., 2019). Alpha and beta diversity, and subsequently gamma diversity, is driven by primary productivity and increasing altitude decreases productivity. Thus, we would expect that spatial diversity components and their relationship would change along altitudinal gradients.
While the spatial diversity components proposed by Whittaker (1960) describe several important characteristics of spatial community variation, they do not elucidate all spatial patterns of community structuring (Dalerum, Vries, Pirk, & Cameron, 2017). In particular, spatial nestedness and modularity are important structures that are not captured by alpha, beta, and gamma diversity (Leibold & Mikkelson, 2002). Although rarely used in evaluation of altitudinal effects on diversity, nested and modular patterns are highly informative for the processes shaping potential altitudinal gradients in diversity. In spatially nested patterns, the most species poor locations are inhabited by species that also exists in more species-rich ones, for example, generalists species (Almeida-Neto, Frensel, & Ulrich, 2012;Thébault, 2013;Wright, Patterson, Mikkelson, Cutler, & Atmar, 1997). Modular spatial community structures, on the other hand, result from species turnover and consist of local communities consisting of unique combinations of species, for example, by specialists species. Modularity is a well-established concept in nonspatial networks such as food webs and social networks (Newman, 2006), but are less investigated in biogeography partly because of previous lack of sufficient theoretical work and statistical methods (Hausdorf & Hennig, 2007;Thébault, 2013). Nested and modular structures are highly suitable to evaluate community transition along environmental gradients (Naud et al., 2019), where nested patterns suggest that communities at low diversity sites are maintained by generalist species. Modular patterns, on the other hand, suggest that each level of the gradients is inhabited by species specialized to its specific conditions. It should be noted that nested and modular structures are not mutually exclusive, for example, spatial nestedness can be found within spatial modules (Fortuna et al., 2010).
In this study, we investigate variation in alpha diversity, measured as species richness and species evenness, beta diversity, and spatial community structure of plant communities along altitudinal gradients at three sites at different latitudes in the Swedish mountains.
While Naud et al. (2019) recently highlighted variation in altitudinal gradients of alpine vegetation within a single alpine landscape in northern Sweden, we here evaluate if large-scale latitudinal contrasts in altitudinal gradients may override such landscape level variations.
We focus on strict quantification of patterns of community variation without evaluating potential hypotheses for what drives these potential patterns, since we regard it appropriate to first identify patterns before rigorous tests are developed to evaluate their causes (Dalerum, 2012). Based on empirical observations, we expect a unimodal or monotonically decreasing alpha diversity (species richness and species evenness) with increasing altitude (Bruun et al., 2006;Naud et al., 2019;Rahbek, 1995Rahbek, , 2005. We also expect a decline in beta diversity with increasing altitude, due to a spatial homogenization of species poor communities at high altitudes. Furthermore, we expect both of these two patterns to be weaker at higher latitudes, due to lower primary productivity at lower altitudes. Finally, we aim to test whether altitudinal variation in plant community composition is driven by species turnover and subsequently by altitudinal specialists or by continuous species loss where high-altitude communities consist of altitudinal generalists, as well as whether such patterns are influenced by latitude.

| Study area
The Swedish mountains form parts of the Scandes mountain range and have a north-to-south orientation with an extent of approximately 950 km. The main part of the Swedish mountains consists of Caledonian bedrock, formed ca 450 million years ago. The Caledonian bedrock mostly overthrust the older basal granitoid bedrock beneath, and forms a diverse mix of different layers of metamorphic and sedimentary bedrock. Botanically important calcareous bedrock occurs in patches and mainly in central Scandes (SGU, 2019;Sjörs, 1999). The upper tree line varies greatly with latitude and local climatic conditions along the mountain range, from ~600 m.a.s.l. in the north to 1,000-1,200 m.a.s.l. in the south (Odland, 2015). In Scandinavia, the alpine vegetation starts at the tree line and is divided into three different vegetation zones. The tree line is natural but defined by complex interactions between abiotic and biotic factors (Körner, 2003;Moen, Cairns, & Lafon, 2008). The low alpine zone ranges from the tree line up to where Vaccinium myrtillus ceases and is dominated by heath vegetation with dwarf shrubs and meadows with graminoids. The middle alpine zone is increasingly more dominated by grasses, sedges, and rushes, with fewer and lower postured dwarf shrubs and extends to where there is no continuous vascular plant cover left.
The upper border of the middle alpine zone is often diffuse and dependent on substrate and microclimatic conditions such as the accumulation and longevity of snow beds. The high alpine zone is dominated by rocky open ground, and the vegetation is sparse with single individuals or small groups of vascular plants, while bryophytes and lichens can form small dense communities (Barry, 2008;Nilsson, 1991;Rydin et al., 1999). Generally, the extent of the alpine vegetation zones is narrower in the west as these areas are more dominated by oceanic climate with colder summers; thus, the upper limit of the alpine zones sinks westward in the Scandinavian mountains (Rydin et al., 1999). Dominant vertebrate F I G U R E 1 Positions of the three study sites in northwestern Sweden (a). The sampling plots, 25 m contour lines and a hillshade background derived of data from Lantmäteriet GDS-höjddata 2 m, at the northern site (b), central site (c), and southern site (d). A generalized layout of the sampling design, 14 horizontal transect each consisting of six 1 m 2 plots (e) herbivores are willow grouse (Lagopus lagopus) and rock ptarmigan (Lagopus muta), mountain hare (Lepus timidus), several species of microtine rodents, and at lower elevations also moose (Alces alces).
All but the southernmost part of the Swedish mountains are subjected to substantial grazing by semidomestic reindeer (Rangifer tarandus) (Moen, 2008). The mountains are popular recreation areas with activities such as hiking during summer and skiing and snowmobiling during winter (Heberlein, Fredman, & Vuorio, 2002).
However, areas above the tree line are almost exclusively without roads, with limited human infrastructure and no resident human settlements.

| Study sites
We collected field data at three southern-facing mountain sites: one in the northern (Kebnetjokka), one in the central (Hemavan), and one in the southern (Helags) part of the Swedish mountains ( Figure 1 & Table 1). The sites abiotic conditions such as bedrock, solar radiation, slope, and aspect were compared using ArcMap and QGIS (Table 1; ESRI, 2015;QGIS Development Team, 2018). We focused the study on plant diversity along south-facing altitudinal gradients, partly to allow for a certain topographic consistency, but also since south-facing aspects have been shown to respond more strongly to climate variations under certain conditions (Danby & Hik, 2007). The sites were distributed from 62°53 to 67°52, and field sampling took place above the tree line at each site. The sites were dominated by alpine heath vegetation with dwarf shrubs and with increasing proportion of graminoids with higher altitude. The distances between the sites are 267 km between the northern and central and 360 km between the central and the southern site, respectively, covering a north-south gradient of over 550 km.

| Data collection
Data collection was conducted between June 27 and July 17 in 2016 starting at the southern site and successively sampling toward the north. At each study site, we used a stratified random sampling design in which 14 transects, each at a fixed and consistent altitude, were placed along an altitudinal gradient. Transects were spaced 50 elevation meters apart and spanned a total of 650 elevation meters.
The first transect was placed in the beginning of each mountain's southern slope, above the tree line (forest line) but below the tree species limit as defined by Körner (2003), and continued to the high alpine vegetation zone. Since the mountain slopes started at different altitudes, the sampled altitude range differed between sites (Table 1).
Each transect consisted of six square sampling plots of 1 m 2 , spaced between 30 and 150 m apart, resulting in a total of 84 plots per site.
To minimize the risk of biased results, we avoided areas in immediate proximity to waterways as well as recently snow-covered areas. We also restricted plot placements to locations with at least one vascular plant. The elevation at each plot was measured using a handheld GPS, and the plots were aligned along each mountain slope.
Within each 1 m 2 plot, all vascular plants were identified down to lowest possible taxonomic level using the flora by Krok, Almquist, Jonsell, and Jonsell (2013). To enable calculations of biodiversity metrics based on relative abundance, we counted the number of plant species in five smaller 20 × 20 cm 2 within each sample plot, one in the middle and one in each corner. A species could thus get a relative abundance value from 0 to 5 in each plot. All taxonomic TA B L E 1 Study site description with site name, location, coordinates, mountain altitude, altitude of closest treeline, sampled elevation range, bedrock conditions (SGU, 2019), and average monthly temperature and precipitation for the summer and winter months, respectively, for the period 1961-1990(SMHI, 2018 (2020). Vascular plant specimens that could not be identified to species level were identified to taxonomic level of genus or in rare cases family (Tables A1-A3).

| Data analysis
Species richness was defined as the number of species or taxa identified within each 1 m 2 sampling plot. Species evenness was calculated using the Shannon diversity index (Hill, 1973a;Shannon, 1948), based on the relative abundance value of identified species within each plot. We used a Euclidean distance to centroid based metric of beta diversity, which has been regarded as appropriate for evaluating turnover along environmental gradients (Anderson, 2006). We calculated this metric in two steps. First, we constructed a pair-wise distance matrix for each site using Bray-Curtis dissimilarities based on binary presence-absences of species in each plot. These were converted into principal coordinate space, in which we calculated the Euclidian distance from each plot to the centroid of all six plots along the same transect. Hence, this value gives a relative beta diversity measurement at each altitude along the slope.
We used generalized linear models with a Poisson error structure and a log link to evaluate the effect of altitude on plant species richness, and linear models to evaluate the effect of altitude on species evenness and beta diversity. We fitted separate models for each site and community metric, all which contained altitude as a continuous predictor.
Because we observed tendencies of nonlinear relationships between plant community metrics and altitude, we evaluated a series of models from a simple linear relationship to segmented regressions with up to three break points. Segmented regression has been suggested as a robust and practical alternative to nonlinear regression or more sophisticated nonlinear techniques such as generalized additive models (Muggeo, 2003). We used an algorithm was the number of break points were prespecified a priori, and the value of the breakpoints and the slope between them were estimated using a bootstrap algorithm (Muggeo, 2008;Wood, 2001). We used Akaike information criterion (AIC, Akaike, 1974) to evaluate which model provided the best fit to each site and plant community metric. Models within two AIC units were considered to have equal empirical support (Burnham & Anderson, 2002).
To provide a heuristic evaluation if any nested or modular structure followed the altitudinal gradient, we evaluated two separate methods of sorting the matrices. First, we calculated the indices using a sorting of rows and columns that generates optimal nestedness and modularity. For nestedness, these matrices were constructed by sorting according to row (number of species per plot) and column (number of occurrences per species) totals (Ulrich et al., 2009), and for modularity, by sorting according to plot and species scores obtained by reciprocal averaging (Hill, 1973b). Reciprocal averaging is an indirect ordination technique that finds ordination scores that optimize the correlation between a set of species and the locations in which they are sampled, and is particularly useful for evaluating environmental gradients (Gauch, Whittaker, & Wentworth, 1977).
We also sorted the columns in the plot/species matrices by altitude.
For each matrix, we compared observed values of nestedness and modularity against expected values from 1,000 null models. We used a null model designed for describing plant communities along gradients (Jonsson, 2001). This null model preserved the average column (i.e., species) frequencies, and hence, also on average the same matrix fills as the observed data. All analyses were made in the

| RE SULTS
In total, we identified 79 vascular plant taxa of which 63 (82.3%) could be identified to the taxonomic level of species (Tables A1-A3).
The northern site was most species-rich with a total of 64 identi- Sorting plots by species richness generated a more spatially nested structure of community composition than expected by chance (Table 3, Figure 5a-c). However, none of the sites showed a nested pattern if the plots were sorted by altitude (Figure 5d- in which case all sites instead had a less nested structure than expected by chance (Table 3). All sites had more modular spatial structures of plant community composition than expected by chance (Table 3), both when the plots were sorted by reciprocal averaging scores (Figure 6a-c) and by altitude (Figure 6d-f). Sorting plots by reciprocal averaging generated 11, 12, and 11 modules at the northern, central, and southern sites, and sorting plots by altitude generated 13 modules for the northern site and 11 modules for the central and southern sites (Table 4).

| D ISCUSS I ON
In line with our first expectation, we observed an expected overall decline in alpha diversity, measured both as species richness and F I G U R E 3 Relationships between vascular plant species evenness and altitude at three study sites along the Swedish mountains: (a) Kebnetjokka, a northern site, (b) Hemavan, a central site, and (c) Helags, a southern site. Each point reflects the Shannon index of species evenness for 1 m 2 plots, which were calculated using the presence and absence of plant species in 5 20 × 20 cm subplots. The lines show the linear model and segmented linear models which provided the best fit for the data, and dashed lines show nonsignificant linear relationship evenness, between the lowest and highest altitude at each site.
However, our second expectation, an equivalent altitudinal decline in beta diversity, was only supported by data from the southern site.
Instead, we observed an altitudinal increase or no significant altitudinal variation in beta diversity at the northern and central sites, respectively. These observations suggest that the overall altitudinal decline in species richness did not result in a spatial homogenization of community composition. The patterns of altitudinal variation in alpha diversity did not vary predictably according to latitude-related primary productivity (e.g., Currie, 1991). We also observed spatial community structures that were nested according to other environmental characteristics than altitude. Similar findings have been reported both for vascular plants (Jacquemyn, Honnay, & Pailler, 2007;Naud et al., 2019) and arthropods (Dalerum, Retief, Havemann, Chimimba, & Janse van Rensburg, 2019;Dalerum et al., 2017), and we suggest that spatial variation of species communities in a given location and area likely are formed by complex interactions of processes across several spatial and temporal scales (Willig, Kaufman, & Stevens, 2003). Such altitude-related species turnover have been found previously (e.g., Naud et al., 2019), and have mainly been attributed environmental regulation of local species communities, combined with evolutionary restrictions in the regional species pool (Laiolo, Pato, & Obeso, 2018). We suggest that site-specific variation in abiotic factors such as bedrock, local hydrology, and microtopography (Nilsson, 1991;Rydin et al., 1999), combined with local variation in biotic factors such as herbivory (Moen, 2008;Moen & Danell, 2003), interacts with the general altitudinal and latitudinal gradients in net available energy in shaping species communities across the observed landscapes. One such abiotic factor could be the presence of calcareous rocks in the northern area, leading to locally favorable botanical abiotic conditions (which also is supported by species findings in the area), disrupting both altitudinal and latitudinal patterns.
Our results were largely confirmatory with regards to alpha diversity, with observations of a unimodal altitudinal pattern (southern site), a monotonic altitudinal decline (central site), and a monotonic decline, but not along the whole gradient (northern site). Unimodal patterns have previously been observed for alpine plant communities in Fennoscandian mountains (Austrheim, 2002;Bruun et al., 2006;Naud et al., 2019;Rahbek, 1995), while a monotonic response to altitude has been observed in other alpine areas (Šťastná, Klimešová, & Doležal, 2012). However, we stress that we have observed the different patterns within a single mountain range. Unimodal patterns have been attributed to a productivity-driven increase in competition and subsequent competitive exclusion at low altitudes, and a corresponding productivity-driven environmental filtering at high altitudes (Callaway et al., 2002;Klanderud & Birks, 2003;Wilson & Nilsson, 2009). This explanatory model suggests that we should have found the strongest unimodal pattern at the southern site, where productivity is highest, and the weakest at the northernmost site.
Our observations only partly supported these predictions, which

| CON CLUS IONS
In conclusion, the patterns of species richness, species evenness, and beta diversity, as well as spatially nested community structures that did not follow altitudinal gradients, point to strong effects of local characteristics for site-specific community composition.
However, we found consistent spatial modularity along the investigated altitudinal gradients, which support that altitude do impose limits to distribution ranges of individual species, and imply that such limits are at least to some extent governed by abiotic regulation. We suggest that spatial variation of vascular communities along altitudinal and latitudinal gradients is likely caused by a combination of interactions of processes that occur across spatial scales. We propose that further studies focus at identifying interactions of biotic and abiotic processes, and at which scales these processes are ecologically relevant. We stress the importance of detailed information of plant diversity patterns and community structure, particularly in areas such as the Swedish mountains, since such information is paramount for our ability to investigate ecological consequences of climate-induced range shifts.

ACK N OWLED G M ENTS
We would like to thank Dr. Helle Skånes for assistance with data processing and the Bolin Centre for Climate Research for logistic support. Funding was provided by Formas, grant number FR-2017/0009 and by the Spanish Ministry of Economy and Competitiveness, grant number RYC-2013-14662. In addition, we would like to thank three anonymous reviewers for providing valuable comments on the manuscript.

CO N FLI C T O F I NTE R E S T
The authors declare no conflict of interest.

E TH I C A L A PPROVA L
The study is in compliance with ethical standards.

DATA AVA I L A B I L I T Y S TAT E M E N T
Data on species presences, Shannon evenness, and centroid distances used for estimation of beta diversity, as well as binary pres-