Spectral Response Analysis: An Indirect and Non-Destructive Methodology for the Chlorophyll Quantification of Biocrusts

Chlorophyll a concentration (Chla) is a well-proven proxy of biocrust development, photosynthetic organisms’ status, and recovery monitoring after environmental disturbances. However, laboratory methods for the analysis of chlorophyll require destructive sampling and are expensive and time consuming. Indirect estimation of chlorophyll a by means of soil surface reflectance analysis has been demonstrated to be an accurate, cheap, and quick alternative for chlorophyll retrieval information, especially in plants. However, its application to biocrusts has yet to be harnessed. In this study we evaluated the potential of soil surface reflectance measurements for non-destructive Chla quantification over a range of biocrust types and soils. Our results revealed that from the different spectral transformation methods and techniques, the first derivative of the reflectance and the continuum removal were the most accurate for Chla retrieval. Normalized difference values in the red-edge region and common broadband indexes (e.g., normalized difference vegetation index (NDVI)) were also sensitive to changes in Chla. However, such approaches should be carefully adapted to each specific biocrust type. On the other hand, the combination of spectral measurements with non-linear random forest (RF) models provided very good fits (R2 > 0.94) with a mean root mean square error (RMSE) of about 6.5 µg/g soil, and alleviated the need for a specific calibration for each crust type, opening a wide range of opportunities to advance our knowledge of biocrust responses to ongoing global change and degradation processes from anthropogenic disturbance.


Introduction
Water scarcity restricts vegetation growth and productivity in most desert regions around the world. Under these hostile conditions, other life forms, such as biological soil crusts (biocrusts) have been described to play a key role in ecosystem primary productivity [1]. Biocrusts are inconspicuous communities, mainly dominated by photosynthesizing cyanobacteria, algae, lichens, or bryophytes that can survive during long drought periods in a dormant state and rapidly become active after small water pulses [2]. These traits allow biocrusts to fix a large amount of atmospheric carbon and nitrogen [3,4], which becomes incorporated within the upper layer of the soil to be used by heterotrophic fungi, bacteria, archaea [5], and other soil surface inhabitants [6,7]. As a result, they improve soil fertility [8] and other soil properties, such as surface stability and soil water retention capacity [9][10][11][12], preventing degradation processes and maintaining the ecosystem capacities to provide multiple services to society [13].
Biocrusts' capacity to fix atmospheric carbon and nitrogen, as well as their role in all other key soil surface properties and processes, is mainly driven by the biomass of photosynthetically active organisms within the biocrust community [14] that have been demonstrated to be one of the best indexes of biocrust development [15,16]. Moreover, it is known that chlorophyll can be used to assess biocrust states in response to seasonal dynamics [17,18], anthropogenic disturbances [19,20], ongoing climate changes [21,22], and restoration activities [23,24]. For all these reasons, large efforts have been devoted during the past decades to obtain reliable measurements of biocrust chlorophyll concentrations [25][26][27]. In a similar way to wet chemical methods used for the extraction of chlorophyll a from plant tissue, most of these techniques require extraction in a solvent followed by spectrophotometric determinations of the chlorophyll-related absorbance peaks in the supernatant, and then conversion from absorbance to concentrations using standard published equations [25][26][27]. Though accurate, these laboratory methods require destructive field sampling and are expensive and time consuming. This thereby limits the applicability for punctual measurements that may not represent the inherent spatial and temporal variability of biocrust communities within the landscape [28,29]. For this reason, although biocrusts have been described as one of the main photosynthetic communities in dryland regions around the world [30], there exist only a few rough estimations of their contribution to the total photosynthetic biomass pool at regional [31] or global scales [32][33][34] and these estimates rarely account for seasonal dynamics in biocrust biomass.
Chlorophyll a, as well as other pigments found in phototrophic organisms, absorbs incoming solar radiation at specific wavelengths. These specific spectral absorptions can be detected by means of field spectroscopy techniques that have been proven as an alternative for the non-destructive estimation of leaf and canopy chlorophyll concentrations [35][36][37]. Moreover, the combination of field spectroscopy techniques with airborne and satellite images, which provide spatially explicit information at high temporal resolutions, represents an ideal solution for the derivation of chlorophyll spatial distributions within the landscape (mapping and scaled-up, site-specific measurements for a whole ecosystem or region), and for large-scale monitoring of chlorophyll dynamics at a reasonable cost [38]. Indeed, during the last few decades, there have been numerous studies using remote sensing methods for vegetation and phytoplankton chlorophyll analysis at different spatial and temporal scales [39][40][41]. Vastly different methodologies and indexes have been developed for these analyses using both multi-and hyper-spectral information [42][43][44][45][46][47]. Contrastingly, though there are numerous studies demonstrating that biocrusts have spectral traits related to the presence of photosynthetic pigments [2,[48][49][50][51][52], remote sensing applications for biocrust chlorophyll a quantification are mostly limited to a few local analyses using the normalized difference vegetation index (NDVI; e.g., Karnieli et al., 2001 [53]) or surface reflectance at 680 nm [2]. This occurs for two main reasons. First, monitoring drylands and biocrusts from remote sensing platforms has been quite challenging because the spectral signal of these sites is often the results of the interaction between different components with different spectral properties and life cycles such as soil, biocrusts, and vascular plants [2,52,53]. Second, when dry, the biocrust spectral signal is very subtle [52], which increases the difficulty in detecting chlorophyll-related absorption peaks, and therefore, the biocrust's chlorophyll quantification. Therefore, any application of spectroscopy for the retrieval of Chla in the biocrust should first examine the spectral response of different biocrust communities colonizing soils with different biocrust coverages and soil conditions, such as differing soil texture, color, or soil quality, before this technology can be directly transferred to produce reliable Chla estimations in the field.
The goals of this study were: (i) to analyze the effect of Chla on the spectral response of different natural and artificially induced biocrust communities in different soils, and with different biocrust cover; and (ii) to explore whether this information, the spectral transformations, and the most widely used spectral indexes proposed in the literature are appropriate for non-destructive chlorophyll a quantification in biocrusts. Finally, we discuss the applicability of the different methods for non-destructive monitoring applicable to long-term studies aimed at analyzing the effects of ongoing climate change on the biocrust communities, as well as their usefulness to monitor biocrust recovery after disturbances under natural conditions and on assisted restoration activities, depending on the available information and requirements.

Sample Collection
Selected with contrasting differences in soil composition and color (Supplementary Table S1), a total set of 93 samples of natural biocrust communities were collected from three semiarid experimental areas within the province of Almeria (Spain) (Figure 1)  At the same experimental areas, soil substrates were collected at a 2 cm depth, air-dried, and sieved to 2 mm in order to remove the small herbs and rock fragments. Three native filamentous nitrogen-fixing cyanobacteria species, Nostoc commune (CANT2 UAM 817), Scytonema hyalinum (CAU6 UAM 820), and Tolypothrix distorta (CANT7 UAM 825) that had been previously isolated from undisturbed biocrusts [57] and cultured ex situ were inoculated separately and in a mixture (1:1:1) onto Petri dishes containing soils from the three different study sites. Furthermore, non-inoculated samples from each soil type were set to serve as a control (hereafter "bare soil"). A total set of 119 inoculated soils were placed in a Fitotron Plant Growth Chamber (University of Almería, Almería, Spain) at a constant temperature (25 ± 1 °C) and a light intensity of 60 µmol photons m −2 s −1 with a 16 h photoperiod for 90 days before the spectral and chlorophyll a measurements were conducted (see details in Román et al., 2018 [58]).
Combining all these samples, we obtained a final dataset of 212 samples (93 natural and 119 artificially induced) that represented a wide range of biocrust coverage and developmental stages in Las Amoladeras is a semiarid steppe dominated by the alpha grass species Macrochloa tenacissima (L.) Kunth (= Stipa tenacissima L.) and sandy loam soils with a 10YR 4/6 color (dry). Most of the open spaces between plants are covered by three main biocrust types that together represent more than 30% of the total area: (a) cyanobacteria-dominated biocrust, composed of well-developed communities of cyanobacteria along with some dark lichens; (b) lichen-dominated biocrusts, mostly Diploschistes diacapasis (Ach.) Lumbsch and Squamarina lentigera (Weber) Poelt; and (c) moss-dominated biocrusts that also contain a small proportion of cyanobacteria (≈15%) (see Chamizo et al., 2012 [8] for further details).
El Cautivo is a badlands landscape developed on Tortonian gypsiferous mudstones with silty loam textured soils [54] that have a 2.5Y 6/3 color (dry). The vegetation is mainly M. tenacissima that grows along with other dwarf shrubs and a dense layer of biocrust mainly dominated by two different types: (a) cyanobacteria-dominated biocrusts that also contained numerous pioneer lichens, and (b) lichen-dominated biocrusts mainly composed of the light-colored Diploschistes diacapasis and Squamarina lentigera [55]. Moss-dominated biocrusts are also representative in some topographical positions, such as in the steep northern aspect hillslopes where low incoming solar radiation levels maximize water availability. A detailed description of the area can be found in Cantón et al., 2003 [54].
The third study site (Gádor limestone quarry experimental area), described in Luna et al., 2018 [56], was chosen as an example of a highly degraded ecosystem where biocrusts are present in an early successional stage. The soil texture in the area is clay loam with a 2.5Y 8/3 color (dry). Two main biocrust types were identified at this site: (a) cyanobacteria-dominated biocrusts, mostly incipient cyanobacteria mixed with small, dark pioneer lichens including Collema sp.; and (b) moss-dominated biocrusts in spaces close to the shrub patches in undisturbed areas.
From each biocrust community described at the three different sites (cyanobacteria-dominated biocrust from Las Amoladeras, El Cautivo, and Gádor quarry; lichen-dominated biocrusts from Las Amoladeras and El Cautivo; and mosses-dominated biocrust from Las Amoladeras, El Cautivo, and Gádor quarry), a minimum of 12 replicates, representing different compositions and coverage values, were collected between 5-7 December 2017. Before sampling, the crust coverage and composition was visually estimated and samples were selected representing different compositions and coverage values. For sampling, a 90 mm (diameter) per 15 mm (depth) Petri dish was used. The Petri dish was pressed 15 mm deep into the soil, a spatula was pushed underneath, and then both were pulled out and the Petri dish was sealed with the cap. Then, samples were transported to the University of Almería, air dried, and stored at room temperature at low light intensities for less than one week until the surface spectra and chlorophyll a extraction measurements were conducted.
At the same experimental areas, soil substrates were collected at a 2 cm depth, air-dried, and sieved to 2 mm in order to remove the small herbs and rock fragments. Three native filamentous nitrogen-fixing cyanobacteria species, Nostoc commune (CANT2 UAM 817), Scytonema hyalinum (CAU6 UAM 820), and Tolypothrix distorta (CANT7 UAM 825) that had been previously isolated from undisturbed biocrusts [57] and cultured ex situ were inoculated separately and in a mixture (1:1:1) onto Petri dishes containing soils from the three different study sites. Furthermore, non-inoculated samples from each soil type were set to serve as a control (hereafter "bare soil"). A total set of 119 inoculated soils were placed in a Fitotron Plant Growth Chamber (University of Almería, Almería, Spain) at a constant temperature (25 ± 1 • C) and a light intensity of 60 µmol photons m −2 s −1 with a 16 h photoperiod for 90 days before the spectral and chlorophyll a measurements were conducted (see details in Román et al., 2018 [58]).
Combining all these samples, we obtained a final dataset of 212 samples (93 natural and 119 artificially induced) that represented a wide range of biocrust coverage and developmental stages in different soil types with contrasting spectral properties.

Spectral Measurements and Biocrust Coverage Estimation
Thirty minutes before the spectral measurements, all samples were irrigated with 30 mL of distilled water. Next, the surface reflectance of each sample was measured using an ASD FieldSpec ® Handheld portable spectroradiometer (ASD Inc., Boulder, CO, USA) with a 3.5 nm optical resolution from 325 nm to 1075 nm. All measurements were conducted in stable light conditions provided by two 14.5 V, 50 W ASD Pro Lamp (ASD Inc., Boulder, CO, USA) halogen lamps, with the optical fiber placed 16 cm above the sample surface, resulting in a field-of-view (FOV) area of 66 cm 2 corresponding to the entire surface of the Petri dish. A total of three replicate measurements, each consisting of the average of the three individual spectra, were collected for each sample and then averaged to produce the sample surface spectrum. A Spectralon (ASD Inc., Boulder, CO, USA) white reference sample was also measured before each set of three biocrust spectra.
In addition, after the spectral measurements, a zenithal picture of each sample was acquired using a digital Canon EOS 600D (Canon Inc., Tokyo, Japan) placed at nadir, 25 cm above the sample surface. Then, the biocrust coverage (%) on the soil substrates was calculated using the supervised maximum likelihood classification in the ENVI 4.7 software (ITT VIS, Boulder, CO, USA).

Chlorophyll a Extraction
After the spectral measurements, a 3-mm thick biocrust subsample was collected from each Petri dish for determining the biocrust chlorophyll a concentrations (Chla). Samples were air-dried at ambient temperature, crushed with a roller, and then ground with a mechanical agate mortar. Chlorophyll a extraction from soil samples was conducted by following the ethanol double-extraction method developed by Castle et al., 2011 [25]. The supernatants from both extractions were combined and the absorbance was measured at 649, 665, and 750 nm. As the inoculated soils were only comprised of cyanobacteria-inoculated strains, a specific equation developed by Ritchie, 2006 [59] for these organisms was applied in order to calculate Chla, as follows: where the absorbance coefficients have the dimensions µg·cm/mL·A, A is the absorbance value at the specific wavelength, V is the volume of the extract (mL), and L is the optical path length of the spectrophotometer cuvette (cm). For natural biocrust communities, comprised of a mix of different photoautotrophic organisms, Chla was calculated by applying the equation for chlorophytes and bryophytes reported by Ritchie, 2006 [59] (Equation (2)): where the absorbance coefficients have the dimensions µg·cm/mL·A, A is the absorbance value at the specific wavelength, V is the volume of the extract (mL), and L is the optical path length of the spectrophotometer cuvette (cm).

Spectral Data Analysis
Before the spectral analyses, all data were subjected to a pre-treatment that consisted of removing the noisy bands (350-400 nm and above 950 nm) and applying a cubic polynomial smoothing filter with a 17-nm window size [60]. Moreover, all spectra were resampled to a Sentinel-2 multispectral resolution, which included broad bands at the blue (490 nm), green (560 nm), red (665 nm) and near-infrared (NIR, 842 nm) regions of the spectra, as well as some narrower bands within the red edge region (740, 783, and 865 nm). This procedure was done by applying the spectral resampling functions included within the "hsdar" package in R version 3.3.3 (R-project; [61]).

Data Transformation
As the reflectance at a single determined wavelength (R) is very sensitive to variable irradiance, background effects, and the geometric arrangement of the sensor and surface [62], we explored the usefulness of well-known spectra transformations that have proven to be very useful for isolating the particular absorption features of different pigments, such as the continuum removal (CR; [63,64]) and the first derivative of the spectra ( ; [65][66][67]). The CR and were calculated for each biocrust type using the pre-treated spectra from 400 to 950 nm by 1 nm. The CR was calculated by applying the convex hull fit over the top of a spectrum using straight-line segments that connect the local spectra maxima so that values equal to 1 indicate the absence of absorption peaks (the spectrum matches the convex curve), while values lower than 1 indicate the presence of absorption peaks [63]. The was calculated using the Savitzky and Golay [60] method with third-order polynomials and a 17-nm window size, based on previous exploratory analysis. We also estimated the band depth (1-CR; BD) for the complete 212-sample dataset and from this we calculated several well-known spectral properties related to two main chlorophyll spectral absorptions at 460 nm and 680 nm including: (i) the area under the absorption peak (area), (ii) the maximum absorption value (max), (iii) the wavelength at which maximum absorption is observed (maxwl), (iv) the wavelength positions of upper and lower half-max values (up-wlhm and lo-wlhm), (v) the difference between wavelength positions of upper and lower half-max values (width wlhm), and (vi) the root mean square error (RMSE) between a Gaussian distribution function and the feature values between the lower half-max and the maximum position and between the upper half-max and the maximum position (gaus-low wlhm and gaus-up wlhm). All analyses were done using the "hsdar" package in R (R-project; [61]).

Band Ratios and Standard Spectral Indexes
Biocrust spectra at hyper-and Sentinel-2 multispectral resolutions were also used to calculate the normalized difference indexes (ND) for all possible band combinations between 450-900 nm, following Equation (3): where ND is the normalized difference index between bands 1 and 2, and RB1 and RB2 are the reflectance values at spectral bands B1 and B2, respectively. Finally, a set of 38 existing published Chla-related narrow-and broad-band vegetation indexes were also calculated for each sample (Supplementary Table S2). The selected indexes were chosen from previous studies that reviewed leaf chlorophyll quantification from spectral information [36,67,68]. Specific details about the formulations and reference sources are summarized in Supplementary  Table S2.

Statistical Analysis
The differences in Chla between the different biocrust types and sample sites were analyzed using the two-way analysis of variance (ANOVA) and Tukey's honestly significant difference (HSD) post hoc test. Moreover, a linear regression analysis was performed to explore the relationships between Chla concentration and R, CR, and values of different bands of the spectrum, as well as the relationships between Chla and all possible ND and standard vegetation indexes. This was repeated for the complete 212 sample dataset and for each crust type separately: (i) artificially induced cyanobacteria biocrusts ("inoculated") on soils from Las Amoladeras, El Cautivo and Gádor quarry; (ii) natural cyanobacteria biocrust communities ("cyanobacteria-dominated") from Las Amoladeras, El Cautivo, and Gádor quarry; (iii) lichen biocrusts ("lichen-dominated") from Las Amoladeras and El Cautivo; and (iv) mosses biocrust ("moss-dominated") from Las Amoladeras, El Cautivo, and Gádor quarry. For all regression models, the differences between predicted values and observed values was assessed using the root mean squared error (RMSE) and the coefficient of determination (R 2 ) of the chlorophyll quantifications were estimated.
Finally, we developed a random forest (RF) multivariate model for Chla estimation, which is applicable to a wide range of biocrust coverage and compositions using all spectral properties derived from the two main Chla absorption features corresponding to the first two peaks of the BD values (at 460 nm and 680 nm). To do this, we used the complete dataset in a 20-fold cross-validation strategy. Moreover, in order to test the influence of biocrust coverage on our Chla estimations, we repeated this model also considering the coverage of incipient cyanobacteria, well-developed cyanobacteria and some pioneer lichens such as Collema cristatum (L.) F.H. Wigg, Endocarpon pusillum Hedw, Fulgensia fulgida (Nyl.) Szat, and Fulgensia desertorum (Tomin), lichens, and mosses. The number of predictor variables performing the data partitioning at each node and the total number of trees to be grown in the model run were set to 20 and 100, based on a previous tuning experiment. Model fit was tested using the mean RMSE obtained during the cross-validation process, and importance of predictor variables was measured using the Gini decrease in node impurity measure, which is computed by permuting the predictor variables with the out-of-bag data in the RF validation approach (details in Liaw and Weiner, 2002 [69]). Finally, in order to get an estimation of linear fit between estimated versus measured Chla concentrations that were comparable with this obtained by linear regression, we applied the obtained RF model to the complete dataset.

Chla of Different Biocrust Types
From the different natural biocrust types considered in this study, the moss-dominated biocrusts showed the highest Chla amongst the three different study sites. The values ranged from 8 µg/g soil in the moss samples with a surface coverage below 28%, to 54 µg/g soil on fully covered samples from Las Amoladeras (  Figure 3 shows the mean reflectance, CR, and ρ of the different natural biocrust communities from the three different study sites, as well as the values of the spectral variables for the artificially induced cyanobacterial biocrusts on the soils from the three sites. Overall, the biocrust samples dominated by white lichens showed higher reflectance values than the bare soils, whereas the mossand natural cyanobacteria-dominated crusts showed the opposite effect ( Figure 3). The artificially  Figure 3 shows the mean reflectance, CR, and of the different natural biocrust communities from the three different study sites, as well as the values of the spectral variables for the artificially induced cyanobacterial biocrusts on the soils from the three sites. Overall, the biocrust samples dominated by white lichens showed higher reflectance values than the bare soils, whereas the moss-and natural cyanobacteria-dominated crusts showed the opposite effect ( Figure 3). The artificially inoculated cyanobacteria biocrusts also produced a slight decrease in the overall surface reflectance of the three different reference soils, but this effect was lower than the effect from the natural communities ( Figure 3). There were also some important differences in the spectral shape between the bare soils and biocrusts. For example, all biocrust types showed an important absorption peak at about 670 nm that was not pronounced in the reflectance spectrum but could be easily identified in the CR, , and BD ( Figure 3). Overall, samples with higher biocrust coverage or Chla had a more accentuated spectral absorption with a minimum in the CR curve at 670 nm and a spike in the (Figure 3), with some exceptions occurring on the lichen and natural cyanobacteria samples (Figure 4b). Subsequently, this spectral absorption was more pronounced on the moss-dominated biocrusts than on the lichen, cyanobacteria, and artificially inoculated cyanobacteria biocrusts (Figure 3), exhibiting a similar pattern to the Chla. Another important spectral absorption at about 500 nm was identified from the CR and of the moss, cyanobacteria, and bare soil samples that was masked in lichen-dominated biocrusts. When the spectral signatures of the different biocrusts were resampled to Sentinel-2 multispectral resolution, the spectral curves showed a similar shape, but differences were subtler than those observed at a hyperspectral resolution ( Figure 5).

Direct Chla Predictions from the Surface Reflectance, CR, and ρ
The relationships between the biocrust spectra and Chla are shown in Figure 6. Overall, hyperspectral information provided more accurate results in Chla predictions than multispectral data, as shown by slightly higher values of the coefficient of determination and lower RMSE (Figure

Direct Chla Predictions from the Surface Reflectance, CR, and ρ
The relationships between the biocrust spectra and Chla are shown in Figure 6. Overall, hyperspectral information provided more accurate results in Chla predictions than multispectral data, as shown by slightly higher values of the coefficient of determination and lower RMSE ( Figure   Figure 5. Sentinel-2 resampled reflectance for different biocrust types in Las Amoladeras, El Cautivo, and Gádor quarry study sites.

Direct Chla Predictions from the Surface Reflectance, CR, and
The relationships between the biocrust spectra and Chla are shown in Figure 6. Overall, hyperspectral information provided more accurate results in Chla predictions than multispectral data, as shown by slightly higher values of the coefficient of determination and lower RMSE ( Figure 6, Table 1). The soil surface reflectance within the Vis-NIR region decreased as Chla increased, with the best regression coefficient found at about 674 nm ( Figure 6, Table 1). However, the regression fits were very low at both multi-and hyperspectral resolutions (R 2 = 0.21 and 0.19, respectively). The CR in the red region and transformation showed a better relationship with Chla than surface reflectance, yielding a coefficient of determination of up to 0.68 and 0.51, respectively. Overall, the CR values decreased as Chla increased, with the best regression fits between CR and chlorophyll observed at 655 nm ( Figure 6, Table 1). On the other hand, showed the best regression fits with spectral bands at the inflection point of the red-edge (about 700 nm).   When the different biocrust communities were analyzed separately, we found that the correlation patterns between the surface reflectance, CR, and and Chla strongly varied amongst the biocrust communities (Figure 7a,b). Whereas the CR within the red region showed the best correlation with Chla on the artificially induced cyanobacteria biocrusts (646 nm, R 2 = 0.71), and at the inflection point of the red edge showed the best results on the moss and natural cyanobacteria dominated biocrust communities (723 nm and R 2 = 0.65, and 720 nm and R 2 = 0.27, respectively) ( Figure 7b, Table 1). Chla within the lichen-dominated biocrusts best fit to the values within the green region of the spectra (581 nm, R 2 = 0.41).

Chlorophyll a Predictions from the Normalized Band Ratios
The normalized band combination ratios were as sensitive to chlorophyll changes as CR and at both the hyperspectral and Sentinel-2 multispectral resolutions (Table 1). For example, when all biocrust samples were analyzed together, the best fits were obtained by combining the chlorophyll absorption band at 687 nm and the reflectance from 731-736 nm for the hyperspectral data, and band 4 and band 5-6 for the Sentinel-2 spectral resolution, yielding a coefficient of determination of 0.69 (Table 1; Supplementary Figures S5 and S10). However, similar to the results found when the Chla-R, CR, and regressions were analyzed, the best band selection and overall fit strongly varied between the crust types. For artificially inoculated cyanobacteria, the best results were obtained using the combination of the surface reflectance at 704 nm (band 5 of Sentinel-2) and the surface reflectance at 730 nm (band 6 of Sentinel-2), yielding a coefficient of determination of about 0.72 (Table 1; Supplementary Figures S1 and  S6). The moss chlorophyll showed a better correlation with band ratios based on spectral reflectance at 518 and 519 nm, yielding a coefficient of determination of 0.65. However, these subtle spectral traits were not correctly characterized at the Sentinel-2 spectral resolution, at which the best correlation was found using band 4 and bands 5-6 (R 2 = 0.62) ( Table 1; Supplementary Figure S9). Band ratios within the red or NIR regions were not sensitive to changes in Chla in the natural cyanobacteria-dominated crusts that showed the best relationship with band ratios combining the surface reflectance at 466 nm and 455 nm ( Table 1; Supplementary Figure S2). As previously shown for R, CR, and , the absolute values were again very low (R 2 = 0.26). On the lichen-dominated biocrust samples, no significant relationship was found between the normalized bands and Chla for most bands. In these crusts, only the region around 680 nm showed a poor relationship with Chla (R 2 = 0.19; Table 1; Supplementary Figure S3).

Sensitivity of the Standard Spectral Indexes to the Chlorophyll a Determination
Similar to the normalized band ratios, the standard spectral indexes also were sensitive to changes in the biocrust chlorophyll concentration. Of the 28 narrow band and 10 broad band spectral indexes (Supplementary Table S2), 6 were robustly sensitive to the biocrust Chla (R 2 > 0.6; Figure 8). The best retrievals for the general biocrust communities (all samples) were provided by the broad band greenness indexes using red and NIR information, such as the NDVI, simple ratio (SR), modified simple ratio (MSR), soil adjusted vegetation index (SAVI), and the optimized soil-adjusted vegetation index (OSAVI), yielding a coefficient of determination above 0.6 at both the hyper and Sentinel-2 spectral resolutions ( Figure 8, Table 1). Biocrust specific indexes crust index (CI) and biological soil crust index (BSCI), on the other hand, did not provide good results.  The broad band greenness indexes also appeared to be good predictors of Chla for the artificially induced cyanobacteria and the natural moss biocrust samples. However, when these biocrust types were analyzed separately, some narrow band pigment-related indexes, such as the modified chlorophyll absorption in reflectance index (MCARI[705,750]) (for mosses) and the ratio of first derivative maxima in red-edge region and green region (EGFR) (for artificially induced cyanobacteria) showed better relationships with the chlorophyll concentrations than the broad band indexes, yielding coefficients of determination up to 0.69 and 0.72, respectively. In the case of the lichen-dominated biocrusts, only the narrow band normalized band ratio between 700-722 nm (BND; [67]) showed a good fit, whereas none of the narrow or broad band ratios fit well with Chla on the natural cyanobacteria-dominated biocrusts. Figure 9 summarized the results of the random forest model performance. Overall, both RF models (the model that exclusively considered spectral properties derived from the two main absorption peaks and the RF model also including biocrust coverage) provided a very good fit with a mean RMSE about 6.5 µg/g soil. Moreover, when they are applied to the complete dataset, in order to get an accuracy estimation comparable with that obtained using linear regression, the prediction accuracy increased (coefficient of determination between predicted and observed Chla values = 0.95; Figure 9). The three most important predictor variables for Chla estimation in both models were maximum depth within the 680 nm absorption peak, gaus-low wlhm, and gaus-up wlhm, whereas the coverage of the different biocrust components played only a secondary role (Supplementary Figure S11).  Figure 9 summarized the results of the random forest model performance. Overall, both RF models (the model that exclusively considered spectral properties derived from the two main absorption peaks and the RF model also including biocrust coverage) provided a very good fit with a mean RMSE about 6.5 µg/g soil. Moreover, when they are applied to the complete dataset, in order to get an accuracy estimation comparable with that obtained using linear regression, the prediction accuracy increased (coefficient of determination between predicted and observed Chla values = 0.95; Figure 9). The three most important predictor variables for Chla estimation in both models were maximum depth within the 680 nm absorption peak, gaus-low wlhm, and gaus-up wlhm, whereas the coverage of the different biocrust components played only a secondary role (Supplementary Figure S11). , and the minimum and maximum values (whiskers) of the mean root-mean-square error (RMSE; µg/g soil) obtained by the 20-fold cross-validations of the model that exclusively considered spectral properties derived from the two main absorption peaks (red) and the RF model also including biocrust coverage (green). (B) Linear regression between measured the chlorophyll a concentration (µg/g soil) and estimations obtained when the RF model that exclusively considered spectral properties derived from the two main absorption peaks (red) and the RF model also including biocrust coverage (green) were applied to the complete dataset.

Discussion
Chla provides valuable information about the photosynthetic capacity of biocrusts [27] and have been widely used as a surrogate for biocrust development and growth of both natural and artificially induced biocrust communities [15,28,58,70,71]. However, until now, most analyses were based on punctual measurements and required destructive sampling that hindered the analysis of seasonal dynamics and growth rates in these spatially heterogeneous communities during long term monitoring experiments. The present study demonstrates that natural and artificially induced biocrust communities present certain common spectral traits related to the presence of chlorophyll that, although less accurate than traditional laboratory analyses, represent an alternative non-destructive and reliable method to quickly quantify Chla and biocrust development.
The spectral analysis of our experimental dataset corroborated that, according to several studies analyzing different biocrust communities around the world [2,48-52], biocrusts modify soil surface , and the minimum and maximum values (whiskers) of the mean root-mean-square error (RMSE; µg/g soil) obtained by the 20-fold cross-validations of the model that exclusively considered spectral properties derived from the two main absorption peaks (red) and the RF model also including biocrust coverage (green). (B) Linear regression between measured the chlorophyll a concentration (µg/g soil) and estimations obtained when the RF model that exclusively considered spectral properties derived from the two main absorption peaks (red) and the RF model also including biocrust coverage (green) were applied to the complete dataset.

Discussion
Chla provides valuable information about the photosynthetic capacity of biocrusts [27] and have been widely used as a surrogate for biocrust development and growth of both natural and artificially induced biocrust communities [15,28,58,70,71]. However, until now, most analyses were based on punctual measurements and required destructive sampling that hindered the analysis of seasonal dynamics and growth rates in these spatially heterogeneous communities during long term monitoring experiments. The present study demonstrates that natural and artificially induced biocrust communities present certain common spectral traits related to the presence of chlorophyll that, although less accurate than traditional laboratory analyses, represent an alternative non-destructive and reliable method to quickly quantify Chla and biocrust development.
The spectral analysis of our experimental dataset corroborated that, according to several studies analyzing different biocrust communities around the world [2,[48][49][50][51][52], biocrusts modify soil surface reflectance with a main absorption detected in the reflectance red peak at about 680 nm, and along the red edge ( Figure 3) that corresponds to the well-known chlorophyll a spectral absorption [72][73][74][75]. However, contrary to previous studies [2,76], we found that the use of single band reflectance values at this wavelength were not able to accurately estimate Chla of biocrust communities, as shown by the lower coefficient of determination and higher RMSE values obtained compared to those obtained with the , CR, ND, and standard spectral indexes (Table 1). This may be due to the fact that our dataset was composed of different biocrust communities dominated by a variety of phototrophic organisms with different spectral traits [2,50] that, as previously reported by Escribano et al., 2017 [77], interact with the soil reflectance in multiple ways, modifying the final reflectance values independently of the biocrust chlorophyll concentration. For example, whereas dark cyanobacteria-and green moss-dominated biocrusts reduced the overall soil surface reflectance when compared to the bare soil spectra, white lichens that occur in the experimental zones caused the opposite effect, especially in the dark-colored soils from Las Amoladeras (Figure 3). A chlorophyll-a-related absorption peak at about 680 nm was enhanced by the biocrust spectra normalization using the CR and spectral transformations. In this study, all biocrust types showed a deep absorption peak in the CR spectra at about 670 nm and a subsequent spike at 700 nm in the spectra (Figure 3). CR and are commonly used to isolate particular absorption features because they are less sensitive to differences in the overall reflectance and color between different samples [36,64,66] and they have been found to be a good predictor of Chla in plants and leaves [64][65][66]78,79]. However, previous studies on natural biocrust communities that described this chlorophyll-related spectral feature were focused on identification and mapping (see Weber and Hill, 2016 [2] and references within) and thereby did not explore the potential of the CR and to quantify biocrust Chla. We found that absorption values within the red and the red-edge spectral region extracted from the CR and normalization, respectively, improved chlorophyll a estimations when compared with reflectance values, yielding a coefficient of determination above 0.6 ( Figure 6). Moreover, the depth of the absorption peaks increased from the very incipient artificially induced cyanobacteria-dominated biocrust to the well-developed, moss-dominated biocrusts as Chla increased (Figure 3), highlighting the potential of spectral transformations to be used as an alternative biocrust development index.
Normalized difference band ratios including the red and NIR spectral bands, as well as the standardized broad-and narrow-band spectral indexes using the same spectral information (e.g., MSR, OSAVI, NDVI, SAVI), were also strongly correlated with the biocrust Chla, giving a similar prediction accuracy as the CR and spectral transformations (Table 1). Thus, results obtained two decades ago by Karnieli et al., 2001 [53] who found a good relationship between the NDVI value of biocrust communities representing different levels of development and chlorophyll, are now confirmed for a dataset composed of 212 samples from different ecosystems that comprised a wide range of biocrust coverage and compositions with chlorophyll a values that varied from 0-54.63 µg/g soil on different soils with contrasting soil textures and physicochemical properties, and thus, spectral properties. This approach may reduce the obtained accuracy in comparison with local applications that only considered one soil or biocrust type, but it increases the transferability of the results to a wider range of applications that consider different soil types and biocrust compositions (including monitoring, analysis and comparison of natural biocrust dynamics, and artificially induced biocrust communities).
The good relationships found between the hyperspectral data and Chla are promising and suggest that the best indexes found in this study (Table 1) can be used for reliable estimations of biocrust Chla. Moreover, in a similar way as Chla, spectral absorption at 670 nm and the value of the different spectral indexes values increased as biocrust development did. This positive relationship allowed us to propose biocrust reflectance as a quick, cheap, and non-destructive sampling method to monitor biocrust development over long-time spans, with several potential applications for field and laboratory studies. Furthermore, the good results obtained with the use of hyperspectral information denote a high potential of space-based hyperspectral imagery, as the already launched DESIS or PRISMA hyperspectral missions [80] or several upcoming missions, such as EnMAP and HySPIRI, with hundreds of spectral channels covering the RED and NIR regions with a high potential for biocrust monitoring and biomass retrieval over a wide variety of locations and scenarios. Of note are also the good results obtained by the application of broad band greenness indexes, yielding models that were as reliable and accurate at predicting Chla as those obtained at hyperspectral resolution ( Figure 8). The good fits obtained via the application of broad band greenness indexes at the Sentinel-2 multispectral resolution, combined with the high temporal and spatial resolution that characterize this satellite, open a wide range of opportunities for biocrust monitoring at regional scales that would allow for their use in monitoring programs aimed at understanding biocrust responses to ongoing global change and degradation processes. However, before imaging spectroscopy can be applied in the quantification of chlorophyll a in biocrusts, more effort is needed in order to deal with additional difficulties related with atmospheric effects, shadows effects, effect of soil water content on spectral response, or interactions with other surface components such as vegetation that absorb solar radiation at similar wavelengths to biocrusts. Moreover, a specific calibration would be necessary as marked differences are found when the biocrust behaviors are analyzed separately (Figure 7a). Chla for the mosses and artificially induced cyanobacteria showed good correlations with most indexes that fit well for the whole dataset, whereas worse results were obtained for lichen-dominated biocrusts and natural cyanobacteria communities. This can be explained by several factors. To cope with high UV levels, natural cyanobacteria communities produce a wide variety of sunscreen pigments, such as scytonemins and carotenoids [81,82], that also absorb solar radiation and interact with the chlorophyll spectral absorption profile [83]. Other species, such as Microcoleus vaginatus or Oscillatoria sp., reside below the soil surface during droughts and glide upwards only during long favorable periods when soil pore water is available [84,85], a trait that reduces detection accuracy. Moreover, in dryland soils with numerous vesicular pores where solar radiation can penetrate the soil profile, recent studies revealed that high levels of activity beneath the soil surface would not be observable using remote sensing technologies [86]. In the case of lichen-dominated biocrusts, the photobiont is usually sheltered in tight fungal structures that may be difficult to detect in chlorophyll-a-related spectra. Moreover, the lichen species analyzed in this study, such as Diploschistes sp., have a rugged surface that can create an additional habitat between the thalli and the soil for other organisms, such as cyanobacteria, that cannot be observed using optical techniques. In these cases, when the chlorophyll variability is below the detection limit, accurate laboratory analyses are necessary, and other approaches such as the use of coverage values, surface albedo, or roughness could be an alternative.
As an efficient machine learning algorithm, our RF model successfully sorted out all difficulties related with spectral dissimilarities between different biocrust types and provided the most accurate estimation of all methodologies tested in this study without the need for a specific calibration for each biocrust type (Figure 9). This can be explained by the capacity of RF to model non-linear relationships [87], but also because the RF model implies the combination of Chla absorption values at a specific wavelength with other properties related with the shape and specific position of Chla absorptions that are specific to the different biocrust components (Figure 3). This represents a step toward the usability of spatially-explicit spectral information with coarser spatial resolution, such as that obtained from future satellite hyperspectral images, as they showed good transferability during the cross-validation process (Figure 9a). However, as in the case of the application of specific absorption features or vegetation indexes, more effort is needed in order to deal with the additional difficulties presented when airborne or satellite images are analyzed.

Conclusions
Our results demonstrate that biocrust reflectance represents a valuable alternative for non-destructive and reliable retrieval of general biocrust chlorophyll a concentration and biocrust development that could be incorporated into monitoring programs. From the different methods analyzed, spectral transformation, such as CR and , as well as normalized band ratios and standard hyperspectral and broad band indexes, provided accurate results. However, such approaches should be carefully adapted to each specific case as these methods were not able to correctly predict subtle differences in chlorophyll a concentration in lichen-dominated biocrusts and natural cyanobacteria communities, and therefore, specific calibrations prior to broader application are recommended. These marked differences in the biocrust spectral behaviors can be solved by the application of RF models based on different spectral properties derived from Chla or by a combination of absorption spectral properties and coverage values.
Supplementary Materials: The following are available online at http://www.mdpi.com/2072-4292/11/11/1350/s1, Figure S1: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at hyperspectral resolution, for cyanobacteria artificially inoculated. Only the significant values (P < 0.05) are represented. Figure S2: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at hyperspectral resolution, for natural cyanobacteria-dominated subsamples. Only the significant values (P < 0.05) are represented. Figure S3: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at hyperspectral resolution, for natural lichen-dominated subsamples. Only the significant values (P < 0.05) are represented. Figure S4: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at hyperspectral resolution, for moss-dominated subsamples. Only the significant values (P < 0.05) are represented. Figure S5: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at hyperspectral resolution, for the entire dataset. Only the significant values (P < 0.05) are represented. Figure S6: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at Sentinel-2 spectral resolution, for artificially inoculated cyanobacteria subsamples. Only the significant values are represented (P < 0.05). Figure S7: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at Sentinel-2 spectral resolution, for natural cyanobacteria-dominated subsamples. Only the significant values (P < 0.05) are represented. Figure S8: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at Sentinel-2 spectral resolution, for lichen-dominated subsamples. Only the significant values (P < 0.05) are represented. Figure S9: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at Sentinel-2 spectral resolution, for moss-dominated subsamples. Only the significant values (P < 0.05) are represented. Figure S10: 2-D correlation plot illustrating the coefficient of determination (R 2 ) of the normalised difference indices for all possible band combinations between 450-900 nm at Sentinel-2 spectral resolution for the entire dataset. Only significant values (P < 0.05) are represented. Figure S11: GINI importance of each variable in Random Forest model. Two options were tested: (a) with cover (green bars) and, (b) without cover (pink bars). IC: Incipient cyanobacteria, C + PL: mix of cyanobacteria and pioneer lichens. Table S1: Soil texture, pH, electrical conductivity, total organic carbon (TOC) and total nitrogen (TN) of the of the three soils employed in this study: Las Amoladeras, El Cautivo and Gádor quarry (from Román et al., 2018). Table S2: Summary of the different spectral indices used in this study. R: reflectance; : the first derivative of reflectance; RBLUE: reflectance in the blue region, RGREEN: reflectance in the green region, RNIR: reflectance in near-infrared region.