Non-Destructive Lichen Biomass Estimation in Northwestern Alaska: A Comparison of Methods

Terrestrial lichen biomass is an important indicator of forage availability for caribou in northern regions, and can indicate vegetation shifts due to climate change, air pollution or changes in vascular plant community structure. Techniques for estimating lichen biomass have traditionally required destructive harvesting that is painstaking and impractical, so we developed models to estimate biomass from relatively simple cover and height measurements. We measured cover and height of forage lichens (including single-taxon and multi-taxa “community” samples, n = 144) at 73 sites on the Seward Peninsula of northwestern Alaska, and harvested lichen biomass from the same plots. We assessed biomass-to-volume relationships using zero-intercept regressions, and compared differences among two non-destructive cover estimation methods (ocular vs. point count), among four landcover types in two ecoregions, and among single-taxon vs. multi-taxa samples. Additionally, we explored the feasibility of using lichen height (instead of volume) as a predictor of stand-level biomass. Although lichen taxa exhibited unique biomass and bulk density responses that varied significantly by growth form, we found that single-taxon sampling consistently under-estimated true biomass and was constrained by the need for taxonomic experts. We also found that the point count method provided little to no improvement over ocular methods, despite increased effort. Estimated biomass of lichen-dominated communities (mean lichen cover: 84.9±1.4%) using multi-taxa, ocular methods differed only nominally among landcover types within ecoregions (range: 822 to 1418 g m−2). Height alone was a poor predictor of lichen biomass and should always be weighted by cover abundance. We conclude that the multi-taxa (whole-community) approach, when paired with ocular estimates, is the most reasonable and practical method for estimating lichen biomass at landscape scales in northwest Alaska.


Introduction
Lichen biomass is an important indicator of grazing impact and the availability of winter forage for caribou, reindeer, muskox, and other animals in northern regions [1], [2]. In northwestern Alaska, these animals rely on forage lichens including many species of Cladonia (''reindeer lichen'', previously Cladina), Alectoria, Bryocaulon, Bryoria, and Cetraria [3]. Lichens provide the main winter sustenance for the Western Arctic Caribou Herd, which is one of the largest caribou herds in North America and is important in the subsistence economy of native Alaskans [1]. Aside from utility as wildlife forage, terrestrial lichen biomass can also be used as a vegetation monitoring metric to assess the impact of disturbances such as fire [3]- [5], climate change [6], [7] and air pollution [8], [9]. However, estimating lichen biomass by destructive sampling is very time consuming and does not allow for assessment of the same area over time. Researchers studying epiphytic lichens have previously approached this problem in places such as Norway [10], China [11], British Columbia [12] and the U.S. Pacific Northwest [13] by developing regression equations, yet very few have estimated biomass for grounddwelling, terrestrial lichens.
Biomass estimation requires accurate volumetric measurements for lichen height and area cover. The method by which cover is measured represents a trade-off between the potential bias of estimates (if estimated subjectively) and the potential to underrepresent rare or patchily distributed taxa (if measured quantitatively) [14]. Point-count estimates have the potential to be less biased given a high enough density of points [15], while ocular cover estimates can assess larger areas rapidly and can integrate patchiness. Therefore, ocular and point count methods must be comparatively assessed for accuracy and for their influence on biomass estimates.
Climate and vegetation type are important predictors of lichen biomass [16] and should be accounted for in landscape-level assessments. The landscapes of northwestern Alaska span several climatic ecoregions and a wide variety of vegetation and landcover types including wetlands, spruce woodlands, and tundra dominated by lichens, low shrubs or graminoids [17], [18]. Within these landcover types, there are strong gradients of substrate pH and Table 1. Habitat type codes [17]: L = Dwarf Shrub-Lichen Dominated; M = Mesic/Dry Herbaceous; P = Open Low Shrub Dwarf Birch/Ericaceous; S = Sparse Vegetation. ''Cladonia'' includes all species except C. stellaris. Ecoregion codes [18]: vascular vegetation physiognomy that influence lichen community composition [19]. Because species differ in growth forms and levels of dominance, it could be expected that variability in lichen biomass and bulk density would be associated with distinctive landcover types. Lichens in northwestern Alaska occur in very diverse, highly mixed assemblages typified by high species richness [3], [19]. Yet, most estimates of terrestrial lichen biomass completed elsewhere have been limited to only a single [5] or a very restricted number of species [2], [20], mainly focusing on the critical forage lichen genus Cladonia. Some workers [20] have suggested that in diverse lichen communities, the most rapid and effective way to estimate biomass would be to integrate estimates over entire communities rather than parsing out unique biomass relationships among individual taxa. Whole-community surveys [21] are also the most realistic way to monitor biomass across landscapes that have heterogeneous or patchy distributions of a large number of species, or for landscapes that have large species turnover (beta-diversity) among sites.
This study explored the feasibility of using height and cover measurements to estimate terrestrial lichen biomass in northwestern Alaska. Our inferences focused mainly on forage lichens (including Alectoria, Bryocaulon, Cetraria and Cladonia spp.) considered both as individual taxa and as the dominant members of aggregate, multi-taxa lichen communities. Here, we extended the range of species to include multi-taxa ''community'' samples that integrated all observed species, and we presented regression slope estimates to use as conversion factors applicable to northwestern Alaska. Our primary objectives in this study were: 1) to estimate the slope of zero-intercept linear regressions for biomass vs. volume; 2) to examine model fit between single-taxon vs. multi-taxa sampling techniques; 3) to examine model fit between ocular vs. point count cover estimation techniques; 4) to determine whether observed and predicted biomass in multi-taxa communities differed significantly among landcover types within ecoregions; and 5) to examine whether height or cover by themselves could each be good predictors of lichen biomass. We also place our results in the context of previous studies and make recommendations for best survey practices.

Ethics Statement
Field sampling was conducted on public lands administered by the US Bureau of Land Management by permission of that agency. We complied with all national and international rules regarding ethics; the research did not involve measurements on humans or animals. Plant material collected for this study was sampled on a very limited scale and therefore had negligible effects on landscapes. We declare no commercial interests or conflicts of interest.

Field Sampling
Field data are available as Dataset S1. We measured lichen height, cover and biomass (both for single taxa and for all taxa aggregated) among two different ecoregions [22] and four landcover types on the Seward Peninsula, northwestern Alaska (Table 1). Although both ecoregions converge on the Seward Peninsula, they differ sharply in rainfall and climate [23]. The first ecoregion, the Northern Seward Peninsula/Kotzebue Sound Lowlands (hereafter, ''Kotzebue Lowlands''), has a Chukchi Seainfluenced and more continental climate with mean annual precipitation of 250-500 mm. The second ecoregion, the Southern Seward Peninsula/Bering Sea Coast (hereafter ''Bering Sea Coast'') is moist with mean annual precipitation ranging from 400-800 mm. Within each ecoregion we sampled four landcover types [17] which formed the basis of long-term lichen monitoring in the Bering Land Bridge National Preserve [24]; these were: Dwarf-Shrub and Lichen Dominated (L), Mesic/Dry Herbaceous (M), Open Low Shrub -Dwarf Birch/Ericaceous (P), and Sparse Vegetation (S) [17]. While there was a newer landcover classification available at the time of sampling, we chose an older data set compatible with NPS's Arctic Network's long-term lichen monitoring studies [24].
In the study area, we assigned plot locations non-randomly. For single-taxon sampling, we targeted sites with continuous coverage of a dominant forage lichen species. For multi-taxa sampling, we targeted sites with continuous coverage of lichens in general, avoiding non-lichen vegetation. Because multi-taxa communities were often dominated by high proportions of forage lichens (Alectoria spp., Bryocaulon divergens, Cetraria spp., or branchedfruticose Cladonia), we recorded the dominant lichen taxon within  (Fig. 1A), multi-taxa samples (Fig. 1B), and single-taxon samples (Fig. 1C). Lines indicate fitted model slopes (see Table 2 for estimates), while grey polygons indicate 95% confidence intervals. doi:10.1371/journal.pone.0103739.g001 Table 2. Estimated parameters and goodness-of-fit statistics for fitted models of forage lichen biomass in northwestern Alaska. each community to facilitate later comparison with single-taxa samples, and also recorded the dominant vascular vegetation type for landcover classification according to [17].

Method
On the ground at each sampling location, we placed a quadrat frame made of PVC pipe (internal area 25625 cm) divided into twenty-five cells (565 cm). In the center of each cell, we measured the height of lichens (if present) using a narrow, ruled metal rod (3 mm diameter) lowered to the base of the lichen until firm resistance at the duff/soil layer was encountered, but not pushing into layers below. We calculated lichen height as the average of the twenty-five points within each plot (including ''no-hit'' zero values) because we wanted to facilitate adaptation to other techniques that might include many lichen absences (e.g., other point count or line intercept surveys). For single-taxon samples, we recorded only the height of the target lichen, while for multi-taxa samples we also recorded the identity of the lichen species and visually estimated its cover to the nearest one percent. Based on ground surface similarity of growth form (i.e., without dislodging the lichen) and bulk density, species found in single-taxon plots were aggregated to the genus level if several were present in the same sample. For example, Alectoria ochroleuca was grouped with A. nigricans, and Cetraria laevigata was grouped with C. islandica. Morphologically similar Cladonia species (C. arbuscula and C. mitis; C. rangiferina and C. stygia) were also grouped, although we measured C. stellaris separately because it had apparently different bulk density and a morphology that was readily recognizable in the field.
After recording height and cover, lichens were misted with water, then manually harvested by cutting vertically along the edge of each quadrat and gently separating intact mats from the ground. We removed non-lichen plant parts, non-target species, obviously dead or decaying lichen parts, and any dirt, gravel or extraneous debris, then placed samples in dry paper bags for transport to the lab. Each sample was completely oven-dried at 80uC until reaching a stable, unchanging mass (minimum 8 hours, no more than 12 hours). Note that this is considerably less drying time than the 76 hours used in a similar study [20]. We weighed each sample with a digital scale (precision 6 0.01 g) within just a few minutes of removal from the drying oven to prevent moisture uptake that might bias dry mass estimates. This also differs from [20], in which lichens were allowed to cool to room temperature before weighing.

Calculations
Our goal was to assess the relationship between biomass and lichen volume (calculated as the product of cover6height). We used two methods to calculate volume; the first was based on the ocular (visual) method of estimating cover described above, and the second was based on the point count estimation of cover. For the point count method, we estimated percentage cover by dividing the number of lichen presences by the total number of points (twenty-five) in each plot. For both methods, we multiplied percentage cover by the area of the plot to convert all measurements to cubic centimeters.
We regressed lichen biomass on volume (or height or cover) using ordinary least squares regression with forced zero intercepts implemented in R version 3.0.1 (R Core Development Team 2013) and compared model slope parameters. Although weighted regressions have been used elsewhere [20] to estimate forage lichen biomass, examination of model residuals (not reported here) suggested this was not necessary given our data. Because goodnessof-fit R 2 statistics are artificially inflated and are essentially meaningless for zero-intercept models, we instead report confidence intervals for each slope, and compare models using three metrics: percent difference in slope, Bayesian Information Table 3. Observed and predicted lichen biomass (g m 22 ) for each sample type and landcover type in northwestern Alaska.  Criterion (BIC) [25] and log-likelihood ratio (LLR) tests [26]. The latter two are likelihood-based methods that do not attempt to assess whether two model parameters differ more than expected by chance (as would a probability-based F-test, for example), but instead they provide a relative measure of how much better or worse each model fits, given the data at hand. For our application, this is more useful than F-tests because the sign and magnitude of difference between two methods provides more information than simply whether or not any such difference occurs.
To determine whether mean observed and predicted biomass in multi-taxa ''communities'' differed significantly among landcover types within ecoregions, we used nested analysis-of-variance Ftests (landcover type nested within ecoregions), accompanied by Tukey's Honest Significant Difference test for pairwise comparisons; this procedure used only the multi-taxa data to simulate realistic conditions.

Biomass vs. Volume Slope
From the regression of biomass on volume, estimates of the slope ranged from 0.0143 to 0.0203 g cm 23 for lichen communities, depending on the method used (Table 2 and Fig. 1A). Individual taxa ranged from 0.0056 g cm 23 (Alectoria) to 0.0221 g cm 23 (Cladonia stellaris). All slopes differed significantly from zero (p,0.05 from t-test of linear slope coefficients). For both observed mass and bulk density, mean values differed among all taxon groups (F-tests p,0.0001) except the pairwise comparisons among Bryocaulon -Cetraria and Cladonia-C. stellaris (for all others, p,0.05 from Tukey test).

Single-taxon vs. Multi-taxa Sampling
The single-taxon method consistently underestimated biomass relative to the multi-taxa method (Fig. 1B vs. 1C). Model fit was always better for the multi-taxa method (Table 2), which in all cases also had less variation in biomass than the single-taxon method (Table 3).

Ocular vs. Point Count
The point count method had lower biomass estimates relative to ocular methods (Fig. 1, filled vs. open symbols). This is directly related to the point count method estimating higher volume ( Fig. 2A) and lower bulk density (Fig. 2B) relative to the ocular method. Model fit was usually but not always better with the ocular method compared to the point count method (Table 2). Among individual taxa, differences in slope between methods (Fig. 3) were more pronounced for samples dominated by taxa with diffuse, filamentous growth forms like Alectoria, rather than for denser taxa like Cladonia that tended to fill entire sample frames.

Ecoregions and Landcover Types
When slope estimates were used to predict standing biomass for each of the four landcover types, the point count method had generally greater deviation from observed values than the ocular method (Table 3), though this was not true in all cases; magnitudes of deviation ranged from 0.01-19.5%. The biomass values reported here for the four landcover classes studied represents the high end of the lichen biomass distribution for each cover type because sites were selected based on high lichen cover to allow for better biomass-to-volume estimates. Of the seven combinations of landcover types nested within ecoregions (Table 1), there were 21 possible pairwise comparisons. Of the 21 comparisons, Tukey tests revealed very few significant pairwise differences for observed mass (3 significant comparisons), mass predicted using the ocular method (3 comparisons), mass predicted using the point count method (2 comparisons), and mass predicted using only height (0 comparisons). Each of these nominally ''significant'' pairwise comparisons had p-values that were close to the p = 0.05 significance threshold, so we do not report them here.

Biomass vs. Height or Cover
Deviation of estimated values from observed values (Table 3) was similar but of slightly greater magnitude when using height as the sole predictor of biomass in multi-taxa communities than when using volume; magnitudes of deviation ranged from 0.5-18.1% for different landcover types. When using cover as the sole predictor of biomass, the deviation of estimated from observed values was far greater than when using volume (about 3-5 times greater: Table 3); this was true of the ocular method (absolute deviation: 3.1-28.4%) as well as the point count method (absolute deviation: 8.0-31.7%).

Biomass as a Function of Volume
When lichen biomass is regressed on volume, the resulting linear slope is equivalent to a mass-to-volume ratio or a measure of the lichen mat's bulk density. Our findings indicate that taxa have different bulk densities as a result of their different growth forms, a trend that is reflected in multi-taxa communities dominated by each growth form. For individual taxa, Cladonia was most dense, while Bryocaulon, Cetraria, and Alectoria were successively less dense. There were congruent patterns in multi-taxa communities dominated by each of those genera, with Cladonia-dominated communities being densest.
The study most comparable with our methods [20] used volumetric estimation of biomass for four species (three Cladonia spp. and Cetraria islandica) in Sweden using 50650 cm plots. Their mean Cladonia bulk densities agree closely with our estimates (within 12% of each other for ocular methods), but were slightly lower, perhaps because they manually pulled lichen mats off the ground (rather than cut them cleanly as we did), which could leave lichen material attached to the ground. By contrast, Cetraria spp. bulk densities from our study were nearly 50% less than that of C. islandica in Sweden [20]. We assume that unmeasured factors (e.g., nutrient availability, moisture, growing season, fire regime, successional status, herbivory patterns) allow Cetraria in Sweden to occur in much denser colonies than on the Seward Peninsula of Alaska.
In Finland, Kumpula et al. [2], reported a much lower bulk density than either our results or those from Sweden [20], possibly resulting from differing definitions of what constitutes a ''living'' portion of a lichen thallus. Dunford et al. [5] did not report a bulk density because they relied solely on lichen cover (not height or volume) for the single species Cladonia mitis. Though their scope of inference extended to 25 sites, Dunford et al. [5] also had a drastically smaller sample size than ours (n = 8 vs.144).
To the extent that we had outlier observations, these were likely due to the unintentional inclusion of non-target vegetation in some plots, especially where lichen thalli were layered over top of unseen moss or shrub tissues beneath. While our volumetric calculations probably included portions of non-target vegetation, we note that this is consistent with all situations in which non-destructive estimates of cover and height are used for applying conversion factors to larger study areas.

Single-taxon vs. Multi-taxa Sampling
Individual genera had generally lower bulk density than the multi-taxa samples. This is probably because lichens on the Seward Peninsula grew in mixed assemblages with intertwined and overlapping layers composed of different growth forms and bulk densities. Thus what superficially appeared to be a continuous cover of one lichen from the surface was in fact usually a tangled mat of several different lichen taxa below. This may be why we found it difficult to locate monotypic mats for every target species, even using a plot size of just 25625 cm.
Estimating biomass-volume ratios for every species found in northwestern Alaska would be a prohibitive task. Given that multitaxa surveys require only general taxonomic training and that they capture variation in biomass at levels equivalent to or better than individualistic single-taxon assessments, we conclude that a multitaxa focus is preferable to single-taxon surveys when lichens across entire landscapes must be faithfully represented.

Ocular vs. Point Count
The use of point-count estimations of cover gave consistently higher cover (and volume) estimates than ocular methods, leading to lower bulk density estimates; this is consistent with other findings [20]. Despite its initial appeal for reducing observer bias, we noted that the point count method often seemed to overestimate cover by ''catching'' or touching disproportionately small pieces of lichen thallus in areas of the plot that lacked appreciable coverage, especially for taxa that were filamentous (e.g., Alectoria) or patchily distributed. Though the ocular method seemed initially less objective, it had the benefit that it allowed observers to visually integrate patchy or diffuse spatial coverage of Regressions were fitted using either single-taxon sampling (top two rows) or multi-taxa samples (bottom two rows). Lines indicate fitted model slopes (see Table 2 for estimates), while grey polygons indicate 95% confidence intervals. In the bottom panel, the distribution of estimated bulk density (all methods, grey boxes) for 144 lichen samples from northwestern Alaska is shown, where dark bars in boxes are median values, boxes represent the interquartile range of data values, and whiskers are the maxima/minima within each group. Mean bulk density differs among all species groups (F-test p,0.0001) except the pairwise comparisons among Bryoria -Cetraria and Cladonia -C. stellaris (for all others, p,0.05 or less from Tukey HSD test). doi:10.1371/journal.pone.0103739.g003 lichens in a more realistic manner. Point count methods overestimate cover of vascular plants, while visual estimates are more accurate [27], and counting enough points to be useful would be logistically limiting in remote Arctic settings. We note that point count methods may still have utility for plots larger than our 25625 cm size because ocular cover becomes less precise in larger plots [15]. Point count methods may also be preferable for long-term monitoring where repeatability is paramount [28] and where results may differ among observers [14]. Future work is needed to resolve variation among observers and identify sources of error.

Ecoregions and Landcover Types
The two ecoregions we examined on the Seward Peninsula appeared to have similar bulk densities (g cm 23 ) of lichens. Bulk densities also differed far less than expected among the four landcover types. We expected that differences in substrate availability, substrate chemistry, moisture and vascular vegetation would be expressed in lichen community biomass because growth forms differ according to those conditions. For example, filamentous groups such as Bryocaulon and Alectoria dominate the Sparse Vegetation (''S'') cover type on the Seward Peninsula, whereas the low-elevation L, P, and M landcover types tend to be dominated by denser Cladonia spp. We suspect that ecoregion and landcover would have been more influential had we sampled across a broader latitudinal gradient (see [29] for one example of an Arctic gradient spanning 1800-km north-south). Note that while the biomass values in Table 3 are a starting point for landscape modeling, they portray only one possible scenario (i.e., one with lichen cover averaging ,80% due to our preferential plot assignment). Future landscape estimates must account for patchy or sparse lichen cover, and therefore require a truly random sample of sites from within each landcover type.

Biomass vs. Height or Cover
In some instances, a direct measure of biomass as a function of average height of the lichen can be useful, avoiding the need to estimate cover. For example, line intercept methods can provide copious amounts of height data in a short time frame and may be more applicable than closed-frame methods in the open, unforested landcover types that are commonly found in northwestern Alaska. However, we found a large amount of variation in the estimates that relied on height as the sole predictor, and we advise that height-based estimates must be always weighted by some measure of cover to prevent overestimation. Given that caveat, transect-based line intercept methods could provide biomass estimates if the total area of the study site is known.
If height alone appeared to be an inadequate predictor of lichen biomass, then cover alone was even worse. This makes logical sense -two different lichen patches could cover an equal surface area yet differ by an order of magnitude in height; observers would have no way of knowing their true volume. Because we found that deviations from true values were 3 to 5 times greater when using cover alone as compared to the other methods, we do not recommend using cover as the sole predictor of lichen biomass.

Extending Spatial Coverage
Remote sensing is one way that estimates can be scaled up to landscape-level inferences. Field surveys can be coupled with remote sensing imagery to construct biomass estimates via several algorithms [4]. Other workers [30] have similarly used Landsat imagery to derive leaf area index and biomass for multiple vegetation groups including fruticose, foliose and crustose lichens (in aggregate). Spectral decomposition of remote imagery is a further refinement that can distinguish among taxonomic or functional groups of lichens [31]. Though location-specific correction factors must always be derived and validated on the ground, remote sensing is a promising technology that should allow rapid estimation of biomass for very large, remote or difficult-to-access locations in northern ecosystems.

Management Applications
The biomass functions presented here promise to have utility in several different areas of applied ecology in arctic Alaska. Studies of ungulate biology should be able to use these equations to help scale plot-level vegetation data (where both lichen cover and height measurements are available) to remote sensing imagery, as in [32]. Long-term studies of ungulate grazing exclosures are underway on Alaska's Seward Peninsula [33], [34], and conversion of lichen cover/height data to forage biomass will be an essential component of these studies. The National Park Service has other long-term lichen/vegetation studies in which both lichen cover and height are measured [33]; these may similarly be coupled with remote sensing modeling using our conversion factors for applications including changes in vegetation community structure with climate change effects [35], vegetation response to changing fire regimes [36], and carbon accounting [37]. The National Park Service's Arctic Network has sufficient density of ground-based measurements to model lichen biomass at a landscape scale [24]. Repeat measurements would allow for temporal and spatial depictions of lichen biomass on the landscape. Improved burn severity geospatial layers and historic fire perimeters would increase the resolution of these layers by allowing for the modeling of successional status, and hence lichen height and cover. Lichen winter range for the Western Arctic Caribou Herd is projected to suffer some decline with climatedriven shrub increase [1], [35], more frequent wildfire [36], and increasing inputs of nitrogen and sulfur from regional development [38]. Continued monitoring of lichen biomass will be critical for detecting and addressing ongoing changes.

Conclusions
To summarize, we found that multi-taxa (''bulk'') sampling, when coupled with ocular estimation, was the preferred method for biomass estimation because it yielded the most accurate estimates and was logistically most efficient in the field. This method required no more than general taxonomic knowledge in the field and was tractable for statistical purposes. The zerointercept regression equation (biomass = 0.0203 g cm 23 ) which we derived by this method is readily applied to lichen communities in northwestern and arctic Alaska, though we caution that this does not necessarily apply to other geographical regions where locationspecific equations will be required. Height alone was not the preferred predictor of biomass, and cover alone was a poor predictor; the best estimates should ideally include both height and cover. Because bulk densities varied among the forage lichen species we measured, we suggest that in cases where there is only a single species of interest, investigators may benefit from using separate regression equations for target species. Lichen biomass estimation has a wealth of applications that will help managers estimate wildlife forage, understand successional trends, detect climate and air quality signals, and account for landscape-level carbon. Therefore, continued monitoring will be vital for understanding how ongoing changes in lichen biomass and distribution affect other elements of Arctic ecosystems.