Broad-scale distribution of epiphytic hair lichens correlates more with climate and nitrogen deposition than with forest structure

Hair lichens are strongly influenced by forest structure at local scales, but their broad-scale distributions are less understood. We compared the occurrence and length of Alectoria sarmentosa (Ach.) Ach., Bryoria spp., and Usnea spp. in the lower canopy of > 5000 Picea abies (L.) Karst. trees within the National Forest Inventory across all productive forest in Sweden. We used logistic regression to analyse how climate, nitrogen deposition, and forest variables influence lichen occurrence. Distributions overlapped, but the distribution of Bryoria was more northern and that of Usnea was more southern, with Alectoria's distribution being intermediate. Lichen length increased towards northern regions, indicating better conditions for biomass accumulation. Logistic regression models had the highest pseudo R2 value for Bryoria, followed by Alectoria. Temperature and nitrogen deposition had higher explanatory power than precipitation and forest variables. Multiple logistic regressions suggest that lichen genera respond differently to increases in several variables. Warming decreased the odds for Bryoria occurrence at all temperatures. Corresponding odds for Alectoria and Usnea decreased in warmer climates, but in colder climates, they increased. Nitrogen addition decreased the odds for Alectoria and Usnea occurrence under high deposition, but under low deposition, the odds increased. Our analyses suggest major shifts in the broad-scale distribution of hair lichens with changes in climate, nitrogen deposition, and forest management.


Introduction
Filamentous "hair" lichens in the genera Alectoria, Bryoria, and Usnea often dominate forest canopies throughout the boreal zone, as well as some temperate forests.Globally, Alectoria and Bryoria mainly occur in cool and cold climates, whereas Usnea occurs worldwide (Brodo and Hawksworth 1977;Thell and Moberg 2011).Hair lichens have important functions in forests.They participate in nutrient and water cycling, provide habitat and food for animals, and constitute a significant part of the winter diet for caribou and reindeer (subspecies of Rangifer tarandus (Linnaeus, 1758); Hauck 2011;Stanton et al. 2014;Esseen and Coxson 2015).Hair lichens are useful indicators of forest ecosystem integrity and have strongly declined in areas with atmospheric pollution (Kuusinen et al. 1990;Bruteig 1993) and intensive forestry (Esseen et al. 1996).Hair lichens such as Alectoria sarmentosa (Ach.)Ach., Bryoria nadvornikiana (Gyeln.)Brodo & D. Hawksw., and Usnea longissima Ach. are now red-listed in Fennoscandia (Kålås et al. 2010;Rassi et al. 2010;ArtDatabanken 2015).
Lichen abundance results from a balance between positive and negative factors that differ among species.Positive factors include availability of specific substrata (Ellis 2012), a suitable combination of water, light, and temperature, and a balanced availability of nutrients (Palmqvist et al. 2008).Negative factors include environmental stressors (e.g., pollution, suboptimal microclimate) and disturbances (e.g., fire, wind, herbivory, forestry; Hauck 2011).For example, the essential nutrient nitrogen enhances lichen growth in low to moderate doses, whereas high doses of nitrogen are detrimental (Geiser et al. 2010;Johansson et al. 2012).Environmental variability and dispersal limitation operating at tree, stand, and landscape scales also influence epiphytic lichens (Ellis 2012).Key factors at a tree scale include tree species, canopy height, branch size, tree age, and nutrient availability (Esseen et al. 1996;Ellis 2012).Important stand factors include horizontal and vertical distributions of the canopy, stand age, and microclimate (Coxson and Coyle 2003;Sillett and Antoine 2004).
Anthropogenic airborne sulphur and nitrogen strongly negatively affected forest ecosystems in Europe and North America (Bobbink et al. 2010;Hauck 2011).The SO 2 emissions were highly detrimental to epiphytic lichens in western Europe but have been reduced since the 1970s (Vestreng et al. 2007).Therefore, many lichens have recovered locally, but emissions of nitrogen are still high or increasing with a negative impact on lichens inhabiting oligotrophic environments such as Bryoria and Usnea (van Herk et al. 2003;Hultengren et al. 2004).Hair lichens are particularly sensitive to pollution and climate change, as their large surface area to mass ratios filter moisture and elements from the air (e.g., Knops et al. 1996;Stanton et al. 2014).In Sweden, nitrogen deposition shows a steep south-north gradient and exceeds 10 kg N•ha −1 •year −1 in southern regions (Pihl Karlsson et al. 2011).This is above or close to the critical load (threshold between harmless and harmful nitrogen deposition) for many lichens (Bobbink and Hettelingh 2011;Pardo et al. 2011;Johansson et al. 2012).
The regional distribution of hair lichens has received considerable interest (e.g., Brodo and Hawksworth 1977;Thell and Moberg 2011), but there are still knowledge gaps.Ahlner (1948) mapped the distribution of several hair lichens in Fennoscandia before the start of large-scale clearcutting in the 1950s.Hair lichen maps based on large-scale surveys are available for Finland (Kuusinen et al. 1990;Poikolainen et al. 1998) and Norway (Bruteig 1993).However, few have applied statistically rigorous methods to analyse how environmental factors affect the regional distribution of hair lichens (e.g., Bruteig 1993;van Herk et al. 2003;Berryman and McCune 2006;Shrestha et al. 2012).Recently, Boudreault et al. (2015) analysed forest characteristics influencing hair lichen distribution in ecosystems dominated by Picea mariana (Mill.)Britton, Sterns & Poggenb.across a large west-east gradient in Quebec, Canada.These studies help to explain mechanisms behind broad-scale distribution of hair lichens.However, we need more detailed data and better models to understand how climate, nitrogen deposition, and forestry interact and influence the distribution of different hair lichens.
Our study analyses factors correlating with the large-scale distribution of the hair lichen genera Alectoria, Bryoria, and Usnea in the lower canopy of Picea abies (L.) Karst. in Sweden, from temperate to boreal and subalpine forests.We base our analyses on data from the National Forest Inventory (NFI) that provides a large probability sample.Our aims are to (i) compare the distribution of the three genera in NFI plots across five regions differing in climate and in human impact; (ii) compare the thallus length (indicating growth conditions) among these regions; and (iii) use logistic regressions for quantification of links between the occurrence of hair lichens and macroclimate, nitrogen deposition, and forest structure.

Study area
The study area is in Sweden (55°N-69°N) with a length of 1500 km and a width up to 400 km (Fig. 1).The productive forests (site productivity ≥ 1 m 3 •ha −1 •year −1 ) cover 23 million ha; in addition, there are 5 million ha of unproductive forest and 2 million ha of other wooded land (Anon 2014).There are strong south-north gradients in climate, forest composition, land use, and element deposition.The climate ranges from humid warm temperate climate in the south to a humid snow climate with a cold summer in most of the country, with polar tundra in northwestern mountains.The temperate (nemoral) zone forms a narrow belt in the south and southwest.It not only has coniferous forest (P.abies), but also comprises substantial amounts of broad-leaved trees, particularly Betula spp.and Fagus sylvatica L., as well as Acer spp., Fraxinus excelsior L., Quercus spp., and Tilia cordata Mill.Most of southern Sweden is in the hemiboreal zone, which is a transition zone between the temperate and the boreal zones, where temperate deciduous trees together with P. abies dominate on nutrientrich soils, whereas Pinus sylvestris L. dominates nutrient-poor soils.The boreal zone, dominated by P. abies, P. sylvestris, and Betula spp., covers most of Sweden.Industrial forestry (mainly even-aged forest management) is the dominant land use, whereas agriculture occurs mainly in the southern and central regions.The mean volume on productive forest land in Sweden is 135 m 3 •ha −1 and is dominated by P. abies (42%), P. sylvestris (39%), and Betula spp.(12%; Anon 2014).The onset of forest exploitation started in the late 1700s in the southernmost parts, with a northward expanding timber frontier during the 1800s and 1900s (Östlund et al. 1997).Since the 1950s, forest management has largely been based on even-aged forestry, i.e., clear-cut harvesting followed by artificial planting of conifers (Östlund et al. 1997).Rotation periods range from 50 years in the south to 120 years in the north (Fries et al. 2015), resulting in a highly fragmented landscape, particularly in southern and central regions (Esseen et al. 2016).Only along the Scandinavian mountain range is the current forestry less intense.

Study species
Alectoria is here represented by only A. sarmentosa, a large, pendent species widely distributed across Eurasia and North America (Ahlner 1948;Brodo and Hawksworth 1977;McMullin et al. 2016).Alectoria sarmentosa is strongly associated with old-growth P. abies forests in Fennoscandia (Esseen et al. 1996).Bryoria and Usnea are species-rich and taxonomically difficult genera containing both pendent and shrubby species (Thell and Moberg 2011)

The NFI
The Swedish NFI is a multipurpose and multiscale monitoring program designed to provide information about forests for developing national-and regional-level policies and international reporting (Fridman et al. 2014).It covers issues such as wood resources for the forest industry, biodiversity, and emissions and removals of greenhouse gases.The NFI includes all forests in Sweden except subalpine birch forests in the Scandinavian mountain range.The design includes stratification into regions (Fig. 1), with different sampling intensities.The survey consists of clusters of sample plots (tracts; square-formed), where the length of tract side varies from 300 to 1200 m among regions.The plots (4-8 per tract in this study) are located around the tract perimeter; they are circular with a radius of 10 m (area, 314 m 2 ).More than 200 variables are recorded, including soil characteristics and presence and cover of different plant species.Characteristics of the trees such as species, diameter, and height are core data.Hair lichens are surveyed on permanent plots (ϳ500-700 per year; remeasured every 10 years).Here, we use data collected between 1993 and 2002 (a full inventory cycle; one measurement per plot) from plots in productive forest land (see above).

Lichen inventory
We recorded hair lichens on one randomly selected live P. abies with a diameter at breast height (DBH, 1.3 m) ≥ 150 mm per sample plot.Presence and maximum length of A. sarmentosa (Alectoria in the following), Bryoria spp., and Usnea spp.were recorded up to a height of 5 m on live and dead branches, as well as on the stem.The length of the longest thallus (maximum length) of each type of lichen was measured to the nearest centimetre.Length represents the actual length of bent and entangled thalli.Maximum length not only indicates growth conditions for lichens, but also correlates with epiphyte biomass (McCune 1990).

Explanatory variables
We focused on ecologically important variables (known to affect lichen growth and survival) and excluded location variables such as latitude and elevation, even though some had high explanatory power.Our goal was to explore possible mechanistic relationships and not to construct models for predicting lichen occurrence at particular locations.Among about 40 candidate variables, we selected 11 variables representing different ecological aspects that correlated with lichen occurrence: six forest variables (DBH, crown limit, stand height, stand basal area, stand age, and site quality), four macroclimate variables (temperature, continentality, precipitation, and rain index), and deposition of nitrogen.The DBH reflects substratum availability (branch size is correlated with DBH).Crown limit represents the height of the lowest live branch of the sample tree and indicates vertical extent of the lower canopy.Stand height is the mean height weighed by basal area and was only included to show differences between regions.Stand basal area (m 2 •ha −1 ), assessed by a relascope, indicates substratum availability at the sample plot level and also reflects canopy cover and thus light availability.Stand age (time for lichen development) is weighed by basal area for stands taller than 7 m.Site quality index for P. abies, an index of potential forest production capacity (m 3 •ha −1 •year −1 ), was estimated from forest and site characteristics.
We obtained climate data with monthly averages for the last full reference period 1961-1990 in a 4 km × 4 km grid from the Swedish Meteorological and Hydrological Institute (SMHI).We extracted mean annual temperature and mean total annual precipitation for each NFI plot using ArcGis version 10.3.An index of continentality was calculated as the difference in mean temperature between July and January.Continentality correlates with less rain and higher diurnal temperature amplitude during summer, implying more frequent dew formation (Gauslaa 2014).Much precipitation in northern areas falls as snow in seasons not supporting high lichen growth rates.Instead of using total annual precipitation with low explanatory power, we calculated a "rain index", which was the total precipitation in months with mean temperature ≥ 0 °C during which lichens can grow.Finally, we extracted data on annual deposition of inorganic nitrogen (dry plus wet depositions).We obtained nitrogen data (NO x + NH x ) in a 20 km × 20 km grid from SMHI using the Match model (available from http://www.smhi.se/klimatdata/miljo/atmosfarskemi).We calculated mean annual nitrogen deposition for 1998-2002 in each grid cell and then extracted data for each NFI plot.

Data analysis
The total data consisted of 5586 plots (but we excluded 105 plots; see below), each with lichen data from one P. abies tree.We calculated percent lichen occurrence, mean length, and standard error (SE) for each region.A general linear model was used to test whether thallus length differed by region.Length was log transformed to obtain approximately normal distributions, and the Tukey-Kramer posthoc test was used to evaluate differences among regions.We calculated means and SE for explanatory variables by regions, as well as their intercorrelations across all regions.
We used logistic regression (Hosmer et al. 2013) to predict occurrence of the three genera using climate, nitrogen deposition, and forest variables.To simplify model construction, we excluded plots classified as unstocked (basal area < 3 m 3 •ha −1 , n = 16), thicket (mean height < 1.3 m, n = 13), and young stands (mean height ≥ 1.3 m and < 10 cm DBH of dominant trees, n = 75).These plots mainly had residual trees from the previous stand (before logging), and the observed lichen data are not strongly linked to current forest structure.We also excluded one plot with no data for basal area, thus resulting in 5481 plots.To account for the possibility of nonlinear relationships, fractional polynomials of first and second degree were applied (Sauerbrei and Royston 1999).We fitted separate and multiple fractional polynomial logistic regression models (Appendix A) with the library mfp in R version 3.1.0(R Core Team 2014).Continentality was removed from the multiple model for Usnea to avoid multicollinarity.We used odds ratios to help interpret the derived logistic regression models.The odds ratio is widely used as a measure of association, as it approximates how much more likely or unlikely (in terms of odds) it is for the outcome (lichen) to be present among those subjects (P.abies) with a one-unit increment in an explanatory variable, i.e., x + 1 versus x (Hosmer et al. 2013).That is, the odds ratio is the relative change in the odds of occurrence when increasing the explanatory variable while holding other explanatory variables fixed.Several analogues to the linear regression R 2 have been proposed for logistic regression, but McFadden's pseudo R 2 (recommended by Menard (2000)) is used throughout this study.It is computed as follows: 1 Ϫ maximized log likelihood for the model containing only the intercept maximized log likelihood for the fitted model which, like linear regression R 2 , is on a [0,1] scale.

Distribution and length
Bryoria was most common (44.6% of the trees), followed by Usnea (37.5%) and Alectoria (16.7%).Alectoria was exclusive on 1.0%, Bryoria was exclusive on 11.3% and Usnea was exclusive on 9.7%, while all genera co-occurred on 7.7% of the trees.The distributions overlapped broadly, but Bryoria had a mostly northern distribution and Usnea's distribution was mostly southern, with Alectoria being intermediate (Fig. 2).Bryoria and Alectoria in particular had low frequency along the southwest coast, and all hair lichens were absent from the very south.Alectoria had its highest frequency in northern and central regions and was extremely rare in southern regions (Fig. 3A).Bryoria gradually increased towards northern regions, whereas Usnea peaked in region 3.

Explanatory variables
All explanatory variables except forest age and continentality decreased from south to north (Table 1).Site quality and nitrogen deposition showed the steepest gradients, DBH showed the weakest gradient, followed by precipitation.Southern regions had fewer live branches in lower canopy than northern regions: 61% of P. abies had a crown limit < 5 m in region 5 compared with 92% in region 1.Intercorrelations between all variables were significant due to the large sample size (Supplementary material, Table S2 1 ).Strong correlations occurred between ecologically important factors such as temperature, continentality, site quality, nitrogen deposition, and rain index, e.g., between site quality and temperature, as well as between nitrogen deposition and rain index (r = 0.88 in both cases).
Figure 4 shows hair lichen occurrence in relation to the combination of mean annual temperature and precipitation.Bryoria had the widest temperature amplitude but was rare in warm climates.The occurrence of Alectoria and Usnea decreased in colder climates.Bryoria occurred over a wide precipitation range, including the driest areas.In contrast, both Alectoria and Usnea steeply declined with decreasing precipitation below 550 mm.Usnea was the most frequent genus in warm and wet climates.

Logistic regression models for separate variables
Most models were nonlinear and all variables had highly significant slope coefficients (Supplementary material, Table S3 1 ) with highest pseudo R 2 value for Bryoria, followed by Alectoria and Usnea (Table 2).The sequence of variables was similar for all genera 1 Supplementary data are available with the article through the journal Web site at http://nrcresearchpress.com/doi/suppl/10.1139/cjfr-2016-0113.   when ranked by R 2 values.Temperature, nitrogen, continentality, and site quality had higher R 2 value than the rain index.Stand age had the highest R 2 value of the forest variables, whereas basal area, crown limit, and DBH had low explanatory power (in terms of pseudo R 2 values).Bryoria showed a steep decrease from cold to warm climates, with almost 100% estimated probability of occurrence at -2 °C to a 12% estimated probability of occurrence at +6 °C (Fig. 5).In contrast, probabilities for Alectoria and Usnea peaked at approximately +0.6 and +3 °C, respectively.The pattern for continentality was similar to that seen for temperature.For the rain index, Bryoria showed a steep decrease from dry to wet climates, and Usnea was more frequent in wetter climates than Alectoria.For nitrogen deposition, the probability of occurrence peaked at 3.5, 4, and at 5 kg N•ha −1 •year −1 for Bryoria, Alectoria, and Usnea, respectively.The probabilities decreased steeply at higher nitrogen levels, particularly in Alectoria.Patterns for site quality were similar to those of nitrogen, with very low probability in the most productive forests.For stand age, Bryoria and Usnea showed steep increases up to about 180 and 80 years, respectively.Usnea peaked at approximately 115 years, whereas Bryoria continued to increase with age.Alectoria increased slower with age up to >200 years, peaking at a higher probability than Usnea.

Multiple logistic regression models
The best multiple model, in terms of R 2 values, was found for Bryoria, followed by Alectoria and Usnea, with all variables highly significant (P ≤ 0.017; Table 3).However, only the Usnea model showed a large increase in R 2 value compared with the best singlevariable model.The Bryoria model was least complex with five variables.Climate and forest variables were included in all models, whereas nitrogen was only included for Alectoria and Usnea.
The odds ratios (Fig. 6) show how odds for lichen occurrence change when increasing one variable in the models.The odds increase significantly when the lower confidence bound is higher than 1.0 and decrease when the upper confidence bound is lower than 1.0.The odds ratio curves for increasing temperature with 1 °C differed among genera.Bryoria decreased at all temperatures.By contrast, warming increased odds for both Alectoria and Usnea occurrences in colder climates but decreased the odds in warmer climates.The thresholds were -0.1 and +1.8 °C for Alectoria and Usnea, respectively.Increasing continentality with 1 °C decreased the odds that Alectoria would occur.Adding 20 mm of rain decreased the odds for Bryoria occurring.Adding 1 kg of nitrogen at low nitrogen levels increased the odds up to approximately five times for Alectoria occurring, but odds decreased at >3.9 kg•ha −1 •year −1 .Usnea showed an even higher increase in odds at low nitrogen levels, whereas odds decreased at levels >5.7 kg•ha −1 •year −1 .Similarly, increasing site quality with 1 m 3 •ha −1 •year −1 increased odds for Usnea occurrence below 4.8 m 3 •ha −1 •year −1 but decreased odds above this threshold.An increase in stand age of 10 years increases the odds of Alectoria being present.Bryoria and Usnea had more than a 2.3-fold increase in odds for occurrence in young stands (15-25 years) but no significant change in odds at approximately 100-130 years.Odds decreased in stands older than 135 years.All genera decreased when adding 1 m 2 •ha −1 to basal area.By contrast, increasing DBH by 1 cm increased the odds for occurrence in both Alectoria and Usnea.Odds ratio curves for crown limit were similar in all genera.

Discussion
Compared with former studies on regional distribution of hair lichens, our study is based on the largest sample so far.Monitor-Fig.5.Estimated probability of occurrence of Alectoria, Bryoria, and Usnea in the lower canopy of P. abies from logistic regression models of the six explanatory variables with the highest R 2 value.Transformations and P values are shown in Supplementary material, Table S3 1 .Figure provided in colour online.ing programs such as NFIs have large potential to analyse factors affecting broad-scale distributions of forest organisms.Although based on simple measures of presence and absence and a crude abundance estimate (thallus length), the large number of sample plots enables highly valuable information on factors potentially shaping the distribution of hair lichen genera.

Factors affecting distribution
The distribution maps (Fig. 2) represent "probability maps" of hair lichen occurrence on randomly selected P. abies in NFI plots.Probability maps differ from standard distribution maps because they do not show the entire distribution range.We recorded data for lichens only on large P. abies; however, hair lichens occur on a broad range of tree species.For example, Bryoria is common on P. sylvestris, particularly in northern regions.Moreover, we only studied the canopy below 5 m height.Nevertheless, our study confirms general distribution patterns of hair lichens in Fennoscandia from previous large-scale surveys (Kuusinen et al. 1990;Bruteig 1993;Poikolainen et al. 1998) but adds considerably more detailed information.Bryoria was most northern, Alectoria was intermediate and rare in the hemiboreal zone of Ahti et al. (1968).Usnea was most southern with a marked decrease from the middle to the northern boreal zone, as in East Fennoscandia (Halonen et al. 1999).By comparing our Alectoria map (Fig. 2A) with the detailed map of Ahlner (1948, p. 25; based on observations in 1790-1947), we note a steep decline in southern Sweden.
Variability in environmental conditions affecting lichen colonization, growth, reproduction, and mortality shape the regional distribution of hair lichens.The spatial distribution of water sources (rain, humid air, and dewfall) also plays a fundamental role for lichen distribution (Gauslaa 2014).Confounded natural and anthropogenic factors complicate the separation and identification of mechanisms controlling the distribution of hair lichens in the study area.However, although forestry is intense throughout Sweden, air pollution is low in northern regions (Pihl Karlsson et al. 2011).This suggests that climate and forest structure are main drivers for hair lichens in northern regions.The more open, low forests and drier continental climate likely cause the high frequency of Bryoria in the two northernmost regions (cf.Gauslaa 2014).Among lichens, Bryoria growing in the lower canopy has one of the lowest capacities to store internal water (Esseen et al. 2015) but rapidly takes up water from humid air and dewfall.Dewfall is an important water source for Bryoria in open, continental forests with large diurnal temperature amplitude (Gauslaa 2014).Bryoria also has dark absorbing pigments (melanins) that efficiently screen and protect the underlying algae from excess light (Färber et al. 2014).These pigments allow Bryoria to grow in lower canopies of open forests, as well as in exposed upper canopies.Their high absorbance causes snowmelt, supplying Bryoria with water and allowing photosynthesis during winter (Coxson and Coyle 2003).By contrast, Alectoria and Usnea have the bright reflecting pigment usnic acid and are more susceptible to high light in the desiccated state than Bryoria (Färber et al. 2014).This mechanism may contribute to the lower frequency of Alectoria and Usnea in northern regions that generally have a drier climate and smaller areas of shady and moist P. abies dominated forests.Alectoria is associated with humid forests in Fennoscandia (Ahlner 1948) and North America (McCune and Geiser 2009;Boudreault et al. 2015).This also applies to many pendulous Usnea species.
The northward decrease of Usnea was correlated with lower temperatures (Figs. 4 and 5), suggesting that Usnea is sensitive to cold climates.However, low temperatures per se may not limit Usnea, as boreal lichens are tolerant to freezing (Hauck 2011).Instead, we hypothesize that the low amount of rain limits Usnea in northern regions (Figs. 4 and 5), where more than one-third of the precipitation falls as snow.Usnea dasypoga, the most common Usnea in northern Sweden, has a lower internal water storage capacity than Alectoria (Esseen et al. 2015).Compared with Alectoria, U. dasypoga needs more frequent hydration to maintain high growth, which may explain why Usnea becomes more restricted to humid sites in northern regions such as riparian forests and gullies.In Norway, Usnea is abundant in rainy coastal areas with high rainfall, whereas Bryoria is scarce (Bruteig 1993), as this genus is sensitive to excess hydration (Goward 1998).
The factors affecting hair lichens in southern regions are less known.Air pollution still reduces hair lichens in southern and southwestern Sweden, despite the strong reduction in SO 2 pollution since the 1970s.Only minor recovery in lichen communities has occurred, likely because nitrogen deposition is still high (van Herk et al. 2003;Hultengren et al. 2004;Grandin 2011).Nitrogen deposition in southern Sweden (9.8 and 14.1 kg N•ha −1 •year −1 in regions 4 and 5, respectively) exceeds critical loads for acidophytic lichens.Our odds ratio curves from multiple logistic models (Fig. 6) suggest that the critical load is 3.9 kg N•ha −1 •year −1 for Alectoria and 5.7 kg N•ha −1 •year −1 for Usnea.Nitrogen was not included in the multiple model for Bryoria, but the separate model (Fig. 5) indicates a somewhat lower critical load than for Alectoria and Usnea.Our critical loads are in the same range as for other nitrogen-sensitive species (Geiser et al. 2010;Bobbink and Hettelingh 2011;Pardo et al. 2011).In a whole-tree experiment in an oligotrophic boreal forest, Johansson et al. (2010Johansson et al. ( , 2012) ) found that Alectoria and Bryoria increased at 6 kg N•ha −1 •year −1 and remained stable or decreased at 12.5 kg•ha −1 •year −1 after 4 years.These levels are higher than our critical loads, but nitrogen loads that are neutral or positive in a short-term perspective could be negative over longer periods.McCune and Geiser (2009) classified most hair lichens as oligotrophs in nutrient-poor sites, but some were mesotrophs.They reported a peak detection frequency of 2.4 kg N•ha −1 •year −1 for A. sarmentosa, 2.0 kg N•ha −1 •year −1 (mean) for 12 Bryoria species, and 2.3 kg N•ha −1 •year −1 (mean) for 12 Usnea species in the Pacific Northwest, USA.However, U. subfloridana was eutrophic, with a peak detection at 6.1 kg N•ha −1 •year −1 (McCune and Geiser 2009).In fact, U. subfloridana probably contributed to our higher critical load for Usnea.It is one of the most widespread Usnea species in southern regions (Halonen et al. 1999;Thell and Moberg 2011), likely favoured by nitrogen on conifers with acidic bark.The lower importance of variables in our models follows partly from sampling confined to large trees and by excluding early-successional stands.However, the higher crown limit in southern regions (mean, 4.0-4.7 m; Table 1) suggests that the scarcity of live branches in the lower canopy contributed to the low occurrence of hair lichens in these regions.Moreover, the high basal area in southern regions suggests that low light limits hair lichens, particularly Bryoria in productive forests such as dense P. abies plantations.Bryoria has a higher light compensation point than Alectoria (Coxson and Coyle 2003).Overall, hair lichens in southern regions are negatively impacted by changes in land use, including clearcutting, planting of conifers, and denser forests.Bryoria and Usnea are also more common on deciduous trees than on P. abies along the western coast (S.Hultengren, personal communication).
Our results reinforce the fundamental importance of high forest age for the accumulation of hair lichens (Esseen et al. 1996;Horstkotte et al. 2011;Boudreault et al. 2015).Timber harvesting has negative impact on hair lichens throughout Sweden because rotation cycles are too short (50-120 years) for high biomass to accumulate (Dettki and Esseen 2003).Alectoria depended most on high forest age and prefers slow-growing forests considerably older than the current rotation cycles.Less efficient dispersal in Alectoria than in Bryoria (Esseen et al. 1996) is a likely cause.

Factors affecting length
Hair lichens have considerable growth potential and can grow several centimetres per year (Stevenson and Coxson 2003;Jansson et al. 2009).Thereby, the longer thalli in northern regions for all genera suggests better growth conditions in these regions.Mean maximum length of Alectoria across Sweden (19 cm) was half the length in northern, old-growth P. abies forests (41 cm, Esseen 2006).Alectoria can be 0.5-1 m long in optimal sites (Esseen and Renhorn 1998), thus the larger dimensions in Alectoria than in Bryoria are not surprising.The small size of Usnea (8.4 cm) across Sweden and of Alectoria and Bryoria in southern Sweden (regions 4 and 5) is remarkable.Short Bryoria and Usnea thalli on P. sylvestris in southern Norway result from air pollution (Bruteig 1993).Usnea dasypoga could be 30-50 cm long or even longer.Thereby, conditions for Usnea are hardly optimal in large parts of Sweden.However, we do not know how many observations are shrubby species such as U. subfloridana.Our length data cannot be converted to lichen biomass in the lower canopy.However, Esseen (2006) showed that a 2-fold increase in maximum length corresponded to a 5-fold biomass increase in Alectoria.Hence, regional differences in biomass are likely significantly larger than indicated by our length data.

Model performance
Many factors influence model performance, including accuracy of lichen observations, forest data, gridded climate, and nitrogen data, as well as variable selection and transformations used in model construction.We have no data on observer accuracy, but lichen identification was likely less problematic in this study focusing on genera and not individual species.However, small thalli are difficult to detect and possibly overlooked in the inventory.
All models were highly significant because of large sample sizes, but the Bryoria model had markedly higher R 2 value than the Usnea model despite small difference in occurrence.We hypothesize that the lower model performance for Usnea results from greater interspecific variability in habitat requirements than within Bryoria.Usnea species are highly variable in morphological and anatomical traits and include short, shrubby species such as U. hirta in drier microclimates, to long, pendulous species such as U. dasypoga and U. longissima in humid forests (Thell and Moberg 2011).The variability in habitat preferences is partly explained by greater variability in water storage capacity in Usnea than in Bryoria (Gauslaa 2014;Esseen et al. 2015;P.-A. Esseen, unpublished data).Future studies should address factors affecting regional distribution of individual Bryoria and Usnea species (Shrestha et al. 2012) and include higher canopy positions.
Our models could probably be improved by including light availability and landscape context variables (not available in the NFI) such as proximity to old forests functioning as a propagule source (Dettki et al. 2000).Distance to forest edge could also be of interest, as hair lichens are sensitive to edge microclimates (Esseen and Renhorn 1998).Sweden's forests, particularly the southern and central regions, are highly fragmented by forestry and agriculture with a large extent of abrupt forest edges (Esseen et al. 2016) that may strongly influence hair lichens.

Implications
If our models truly reflect species requirements, they suggest major range shifts in the distribution of hair lichens following changes in climate, nitrogen deposition, and forest structure.Climate scenarios for Sweden predict a 2-6 °C increase in mean annual temperature and a 10%-40% increase in precipitation by the end of the century, with the largest changes in northern regions (Sjökvist et al. 2015).In fact, the changes in the estimated odds of occurrence were surprisingly large following relatively moderate increases in some of the variables such as a temperature increase of 1 °C (Fig. 6).Our models suggest that a warmer climate would always be negative for Bryoria.By contrast, warming might be positive for Alectoria and particularly for Usnea in colder climates but negative in warmer climates.Such effects also depend on whether predicted altered precipitation patterns during warmer winters would affect lichens negatively (cf.Bjerke 2011).Results also suggest that a wetter climate facilitates Usnea (Figs. 4  and 5), whereas Bryoria would decrease (Fig. 6).If nitrogen deposition increased, the occurrence of Alectoria and Usnea, probably also of Bryoria, would decline in most of Sweden but not necessarily in northern regions with low deposition.Extending the rotation periods would increase hair lichen occurrence (Figs. 5 and 6), whereas shorter rotation periods would have a very negative effect (Dettki and Esseen 2003), particularly for Alectoria.If several factors change simultaneously, we could expect synergetic responses that could dramatically disrupt ecological functions of hair lichens in forests.This may have cascading negative effects on other organisms such as reindeer, canopy-living invertebrates, and passerine birds (Pettersson et al. 1995) and adversely impact water and nutrient cycling (e.g., Knops et al. 1996;Stanton et al. 2014).

Conclusions
Hair lichens are key components in forest canopies but face multiple environmental hazards threating their function in forest ecosystems.Our models may help to interpret and (or) predict current and future changes in regional distribution of hair lichens.Such knowledge is critical for choosing adequate measures to mitigate impacts of climate change, element deposition, and forestry on hair lichens.The correlations between hair lichen occurrence and explanatory variables presented here suggest mechanistic hypotheses for future experimental testing.We conclude that broadscale distribution of hair lichens in the lower canopy of P. abies is stronger coupled to macroclimate and nitrogen deposition than to forest structure in Sweden with strong south-north environmental gradients.The three hair lichen genera show broad similarities in response to environmental variables but differ in important aspects such as the shape of relationships with temperature, rain, nitrogen deposition, and stand age.The multiple logistic regression models suggest major range shifts in the regional distribution of the three hair lichen genera following future changes in climate, nitrogen deposition, and forest management.Our study also shows that odds ratios is a useful tool for interpreting complex logistic regression models of lichen occurrence.More studies are needed from other parts of Eurasia and North America to test our hypotheses at global scales.

Fig. 2 .
Fig. 2. Distribution and maximum length of (A) Alectoria, (B) Bryoria, and (C) Usnea in the lower canopy of P. abies in Sweden.Circle area is proportional to mean lichen length of all trees with lichen occurrence in each tract (cluster of plots).

Fig. 3 .
Fig. 3. (A) Occurrence and (B) maximum length (mean ± 1 standard error) of Alectoria, Bryoria, and Usnea in the lower canopy of P. abies in five regions.See Fig. 1 for depiction of the regions.

Fig. 4 .
Fig. 4. Occurrence of Alectoria, Bryoria, and Usnea in the lower canopy of P. abies in relation to mean annual temperature and mean annual precipitation (1961-1990).Black dots show plots with lichen occurrence, and blue dots show plots with P. abies present but without lichens.Figure provided in colour online.

Fig. 6 .
Fig. 6.Plots of odds ratio functions (solid lines, log scale; reflecting relative change in probability of occurrence) with 95% confidence bands (dashed lines) for an increase in one unit of variables in multiple logistic regression models for Alectoria, Bryoria, and Usnea.The grey horizontal lines indicate an odds ratio of 1.0.Note the different scales for the odds ratio.Figure provided in colour online.

Table 1 .
Means (±1 standard error) of explanatory variables in the five regions and range across all regions.
a For sample trees.b Stand height and precipitation were not included in the logistic regression models.

Table 2 .
Pseudo R 2 values for separate variables (after transformation; see Supplementary material, TableS31 ) in logistic regression models predicting the occurrence of Alectoria, Bryoria, and Usnea in the lower canopy of P. abies (N = 5481).

Table 3 .
Summary of multiple logistic regression models predicting occurrence of Alectoria, Bryoria, and Usnea in the lower canopy of Picea abies (N = 5481).See Table1for explanations of variable abbreviations. Note: