Using Spectral Indices Derived from Remote Sensing Imagery to Represent Arthropod Biodiversity Gradients in a European Sphagnum Peat Bog

: Monitoring of peatlands is an important conservation issue. We investigated communities of soil mites (Acari: Oribatida, Mesostigmata) inhabiting a relatively undisturbed European boreal mire characterized by a mosaic of oligotrophic and meso-eutrophic areas. We assess the potential of using remote sensing approach as a mapping and predictive tool for monitoring productivity and arthropod biodiversity in a peat bog. In georeferenced plots, Acari biodiversity, water table level, water pH and plot productivity class on the oligotrophic-eutrophic gradient were recorded. Data from the Landsat 8 OLI sensor were used to calculate several spectral indices known to represent productivity and surface moisture gradients in terrestrial ecosystems. We then explored the relationship between spectral indices, environmental gradients and biodiversity of mites. We found that several spectral indices were signiﬁcantly and consistently correlated with local environmental variables and biodiversity of soil mites. The Excess Green Index performed best as a predictor of plot trophic class on the oligotrophic-eutrophic gradient and showed signiﬁcant relationship with Oribatida diversity in 2016. However, following hot summer in 2019, there was no signiﬁcant relationship between abundance and species richness of Oribatida and remotely sensed data; there was a weak correlation between abundance of Mesostigmata and spectral indices which represent surface moisture gradient (e.g., Normalised Difference Moisture Index). We discuss advantages and challenges of using spectral indices derived from remote sensing imagery to map biodiversity gradients in a peatland.


Introduction
Ombrotrophic (rain-fed) peat bogs are characterized by high water table, acidic nutrient-poor conditions and dominance of Sphagnum mosses and occur in the boreal and temperate zone on most continents [1]. Worldwide, peatland environments are important as providers of ecosystem services and as long-term carbon storage reservoirs [2][3][4].
The plant communities of peat bogs provide information on ecosystem processes such as primary production and carbon sequestration and are sensitive indicators of environmental change [5][6][7][8][9][10][11]. Based on plant communities, bog environments are usually classified along the oligotrophy-eutrophy 'productivity' gradient [12]. Although not strictly a nutrient gradient, as plant communities reflect predominantly the gradient of pH values, the oligotrophy-eutrophy gradient frequently also reflects the availability of nutrients [13,14]. If a bog receives an additional nutrient input, either globally in the form of nutrient pollution from atmospheric sources, or locally from mineral-rich ground water seepage or from surface streams, the resulting areas of nutrient enrichment are characterized by higher pH and characteristic changes in vegetation [15]. The second major gradient in bogs is the ground water level, driven by the micro-topography of hummocks and hollows, which gives bogs their characteristic patterning and which is linked to hydrology and carbon sequestration [1,16]. Both productivity and microtopography exert significant controls on patterns of plant and animal diversity in peat bogs.
Bog environments are a habitat for many species; among invertebrates which are abundant in peat bogs are mites (Arthropoda: Acari). Mite taxa such as Oribatida are well studied, abundant and diverse, and are known to respond to a wide range of environmental and anthropogenic stressors in bogs and elsewhere [17][18][19][20][21][22][23][24][25]. At a site level, the diversity patterns of mites are influenced by the same gradients that are recognized as richness drivers for broader groups of bog organisms, such as pH, nutrient availability and ground water level [25,26].
Peat bogs worldwide have been affected by nutrient pollution, peat extraction, drainage and other types of development, leading to degradation and loss of biological diversity [2,3,27]. Consequently, monitoring of environmental conditions and productivity of peatlands is important. Due to access difficulties and expense of ground sampling in remote peat bogs, coupled with the need for spatially-explicit landscape-scale information, various remote sensing methods have been used for monitoring of peatlands [28][29][30][31][32]. For example, remote sensing has been used to obtain information on bog hydrology and water table depth, functional types and phenology of plants, and to map vegetation classes and fertility gradients [9,28,29,31,33,34]. Other relevant values, such as the Leaf Area Index or the Normalized Difference Vegetation Index, both important parameters linked to terrestrial ecosystem productivity, are routinely mapped using spectral indices derived from remote sensing imagery [10,[35][36][37][38][39].
Because satellite data can be remotely obtained and are free for satellites such as Landsat 8 or Sentinel-2, they are useful to enhance our understanding of response of bog arthropod communities to productivity gradients in less accessible areas. In this study, we investigated communities of free-living mites inhabiting a relatively undisturbed European boreal mire "Shichengskoe". Shichengskoe mire system is characterized by a mosaic of oligotrophic and meso-eutrophic areas, and our previous data for this system show that the pH gradient and the nutrient availability gradient in this bog are linked [25]. The trophic class of sampled plots (classified on a oligotrophy-eutrophy gradient based on plant community features) was one of the best predictors of abundance and species richness of non-aquatic oribatid mites in this bog [25]. For aquatic Oribatida, water table depth was of significant importance [25]. Mesostigmatid mites of the bog have not been studied previously. Here, we aimed to assess the potential of using remote sensing data as a mapping and predictive tool for monitoring arthropod biodiversity of a peat bog. As the connection between environmental data and mite biodiversity can been demonstrated directly, we aimed to confirm the link between spectral data and environmental data and to check if the link between spectral data and biodiversity data is consistent with the ground survey results.

Study Site
The mire system "Boloto Shichengskoe" in Vologda region in the north-western Russia (59 • 56 30.4 N, 41 • 16 57.1 E, 120 m a.s.l.) is a large (15,900 ha) wetland system of predominantly lacustrine origin (Figure 1a). The central part of the mire is occupied by a shallow oligotrophic lake Shichengskoe (1060 ha). Extensive area of the mire is the oligotrophic peat bog, dominated by Sphagnum mosses. Within the oligotrophic bog, there are several ground water seeps, associated with forested (Picea-Pinus-shrublets) islands. Areas near the seeps and the south-eastern part of the mire system are meso-and eutrophic [40]. The climate in the region is humid continental (Dfb in Köppen climate classification) with long moderately cold winter (mean temperature of January −12 • C) and short warm summer (mean temperature of July 16-17 • C). Annual precipitation is 500-650 mm; snow cover lasts 165-170 days of the year [40].

Sampling on the Ground
Two ground datasets were used in this study. The first dataset (training set) was collected in July 2016 in Sphagnum-dominated communities in the ombrotrophic western part of the mire system ( Figure 1b). Samples were collected from 48 plots dominated by Sphagnum species; each plot was 1 m 2 in size. The minimum distance between plots was 5-7 m, but usually over 30 m. In each sampling plot, the following was recorded: GPS coordinates, Sphagnum moss identity, water table level, water pH and Sphagnum nutrient content (C, N, P, K) (see details in [25]). Based on plant community features, each plot was assigned to one of the five qualitative productivity classes ("trophic classes") on the oligotrophic-eutrophic gradient (oligotrophic, oligo-mesotrophic, mesotrophic, meso-eutrophic, eutrophic) [41,42]. Such trophic classes do not directly reflect nutrient status, but instead are based on plant indicator species of the poor-rich (acidity-alkalinity) gradient-however, in Shichengskoe mire these trophic classes are correlated with pH and Sphagnum nutrient content measurements [25]. The distribution of trophic classes among plots was not homogeneous, as meso-eutrophic and eutrophic patches are rarer on the ground. There were 21 oligotrophic plots, 18 oligo-mesotrophic and mesotrophic plots, and 9 meso-eutrophic and eutrophic plots. Sphagnum moss for mite extraction was collected as 10 × 10 cm samples to the depth of living moss plants (including capitula and the length of stems) (one sample per plot). Mites from moss samples were extracted in modified Berlese funnels until samples were fully dry (at least for five days). Adult Oribatida were identified to a species level using published keys and original species descriptions and classified into two functional groups following [22,43]-'aquatic' species, living on submerged vegetation in freshwater habitats, and 'terrestrial' (all other species). Oribatida juveniles were excluded from abundance and richness counts.
The second dataset (validation set) was collected in August 2019 using the same methods, except Sphagnum nutrient content was not measured. Sampling was again limited to 1 m 2 plots in Sphagnum-dominated communities but covered a larger extent of the mire (Figure 1b). In each sampling plot GPS coordinates, Sphagnum moss identity, trophic class on the oligotrophic-eutrophic gradient, ground water level and pH were recorded. Mites were sampled and extracted using the same methods as above. Adult Oribatida and Mesostigmata were identified to a species level. Aquatic and terrestrial Oribatida were counted separately, with juveniles excluded. Mesostigmata are all terrestrial; their juveniles were excluded from abundance and richness counts. Some of the records in this dataset were excluded from the analysis as plots were located too close to each other-at the end, the dataset comprised 42 records (see Data S1 in Supplementary Materials).

Remote Sensing Data
The U.S. Geological Survey Landsat 8 OLI/TIRS scenes for path 177, row 18 covering Vologda region were used with the 2016 (scene LC81770182016151LGN01, 30 May 2016) and 2019 (scene LC08_L2SP_177018_20190608_20200828_02_T1, 6 August 2019) data sets. The satellite imagery search criteria were (a) acquisition date close to the date of sampling and the period of maximum vegetation activity and (b) minimum cloud cover (5% threshold). The 24-Jul-2016 scene, also available, was rejected due to the larger % cloud cover and the presence of haze.
Remote sensing images require radiometric and geometric correction. The path 177 row 18 scenes for 30 May 2016 and 6 August 2019 are "level T1", indicating that they have been geometrically corrected. As a rule, systematically corrected T1 Landsat 8 data have geodetic accuracy of ≤12 m [44].
The remote sensing data were processed using the ArcMap 10.6, Esri Inc., Redlands, CA, USA. Landsat 8 scenes were cropped to the area of interest and co-registered. The 16-bit digital numbers (DN band ) recorded in the Landsat 8 OLI spectral radiance bands 2-7 were converted to the top of the atmosphere band reflectances (r band ) using solar zenith angle and band-specific scaling coefficients provided in the Landsat 8 OLI metadata file [44]. Spectral bands were significantly correlated to local environmental data (Table 1). However, such correlations are difficult to interpret biologically, so we focused on known spectral indices which represent productivity gradients and surface moisture gradients in terrestrial ecosystems. The following spectral indices representing productivity gradient and ground water level gradient were calculated for the area of interest and sampled using ground plot coordinates: (i) Normalized Difference Vegetation Index (NDVI), computed using band 4 (red) and band 5 (NIR) reflectances as (r5 − r4)/(r5 + r4). NDVI is an indicator of the amount of vegetation; it approaches 1.0 if a pixel contains vegetation; 0 if a pixel contains soil; and −1.0 if a pixel contains water. NDVI is commonly used in remote sensing as a proxy for productivity [36]; (ii) Excess Green (ExG), calculated as 2r3 − r4 − r2; (iii) Excess Green minus Excess Red (ExG−ExR), calculated as 1.4r4 − r3 using normalized reflectances of band 2 (blue), band 3 (green) and band 4 (red) [45]. Both ExG and ExG−ExR have been found useful as proxies of gross primary productivity [45,46]; (iv) Normalized Difference Moisture Index (NDMI), computed using band 5 (NIR) and band 6 (SWIR1) reflectances as (r5 − r6)/(r5 + r6) [47]. NDMI is used to determine vegetation water content; it is sensitive to changes in liquid water content and in spongy mesophyll of vegetation canopies [47,48], otherwise known as NDWI (normalized difference water index); (v) Moisture Stress Index (SWIR1/NIR), computed as r6/r5; this index is negatively correlated with surface water content and has been suggested as a broad-band index of surface moisture (reflective of water table position) in peatlands [32,49]. Moisture Stress Index is used for canopy stress analysis, productivity prediction and biophysical modeling [50].

Data Analysis
The values of spectral indices were calculated and extracted from the satellite imagery for the GPS coordinates of plots using the ArcMap 10.6 software package. The values of spectral indices were compared with local-scale environmental and biodiversity parameters for that year using Pearson correlation in R 4.0.5. Variables at a plot level in 2016 were: water table level; water pH; Sphagnum nutrient content (C, N, K, P); trophic class; Oribatida abundance and species richness (see [25] for detailed analysis of 2016 data). Variables at a plot level in 2019 were: water table level; water pH; trophic class; Oribatida abundance and species richness; Mesostigmata abundance and species richness.
The 2016 dataset was used as a training set-using the 2016 data, the spectral indices with the significant correlation to environmental and biodiversity parameters were selected to develop the regression model. Random forest regression (randomForest package in R) [51] was employed for variable selection and to predict the trophic class of ground plots. The output provides the total % variance explained, and the importance score for each explanatory variable. The significance of random forest model was tested using permutation procedure (rfUtilities package in R 4.0.5; [52]).

Results
Previous results showed that plot trophic class and water table depth were the two best predictors of Oribatida biodiversity in Shichengskoe mire [25]. Higher abundance and species richness of oribatid mites in 2016 were correlated with higher productivity and lower water table depth (Table 2). In 2019, we observed no relationship between plot trophic class and abundance, or species richness of mites (Table 2) and the total abundance of mites collected in 2019 sampling was low (Appendix A, Tables A1 and A2), possibly due to the heat wave anomaly-June and July 2019 were the hottest months on record for the region. The relationship between water table level and mite diversity, on the other hand, remained consistent between 2016 and 2019, with higher abundance and richness of terrestrial Oribatida and Mesostigmata associated with lower water table (Table 2). Table 2. Pearson's correlations (r, + p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001) for local environmental data and mite diversity in Shichengskoe mire; na-data not collected in 2016. Landsat 8 OLI spectral indices representing productivity gradient and ground water level gradient were consistently correlated to local environmental gradients in Shichengskoe mire in both 2016 and 2019 data sets (Tables 3 and 4). ExG and NDVI performed well as productivity indices, showing significant correlative relationship with trophic classification of sampling plots. The NDMI and SWIR1/NIR showed significant correlation with water table depth in both data sets. ExG−ExR was correlated to both trophic class and ground water level in 2016, but only to trophic class in 2019.  Spectral indices representing productivity gradient were significantly correlated with oribatid mites biodiversity in 2016, consistent with the ground data; the strongest relationship was between Oribatida biodiversity and the ExG index (Table 3, Figure 2a,b). However, spectral indices representing ground water level gradient showed no correlation with Oribatida diversity in 2016, even though the ground data (Table 2) suggest that this gradient was important. In the 2019 dataset, a significant correlation was seen between abundance and richness of Mesostigmata and spectral indices representing water table gradient, which is consistent with the ground data (Tables 2 and 4). Random forest regression model predicting trophic class using spectral indices was significant at p < 0.001, with ExG index again showing as the best predictor of trophic class (Table 5). The random forest model based on 2016 satellite and ground data predicted the trophic class of plots visited in 2019 with 69.0% accuracy; this of course is relying on the assumption that the trophic class did not change between 2016 and 2019.

Discussion
The use of data derived from remote sensing imagery for assessment and prediction of biodiversity is receiving increasing attention [53][54][55][56]. Unlike vegetation, soil communities cannot be directly linked to spectral indices, but close links between landscape forms, vegetation and soil biota allow us to link information from remote sensors to soil biodiversity patterns. There are only a few such studies-but, for example, satellite-derived spectral information has been used to predict patterns of biodiversity in soil microbiome communities, soil mesofauna (springtails) and earthworms [57][58][59].
Our results demonstrate both the potential and the limitations of using freely available satellite data for generating information on environmental gradients, productivity and biodiversity in a peatland. Among spectral productivity indices, the ExG index was the best predictor of trophic class and could be used as a broad-scale variable representing productivity gradient in the Shichengskoe mire, with a significant relationship with Oribatida diversity. Both ExG and ExG−ExR have been used successfully elsewhere as a proxy for gross primary productivity [45,46]. On the other hand, NDVI is known to perform less well in peatland environments due to atypical near-infrared reflectance of Sphagnum mosses [60,61].
The application of spectral indices derived from remote sensing imagery to map ecological gradients in peatlands presents several challenges. The first of these is the scale-for medium-scale sensors such as Landsat 8 OLI or Sentinel-2, it results in aggregation and averaging of patterns and features which are smaller than the pixel size of the sensor [28]. The use of fine-scale sensors such as IKONOS improves the accuracy of peatland land cover mapping [31] but is more costly. The second challenge is the difficulty of deriving the details of bog microtopography from the satellite sensor data. Peat bog microtopography (hummocks and hollows) are linked to hydrology and biodiversity [1,24,26]. In the Shichengskoe mire, two major environmental drivers which explained the abundance, species richness and community composition of Oribatida were trophic class, linked to acidity-alkalinity (pH) gradient and nutrients (N-P-K) availability, and water table level, linked to microtopography [25]. For aquatic Oribatida, the water table was the single most important variable [25]. Our 2019 data (this study) suggests again that microtopography and moisture are significant drivers for diversity of both Oribatida and Mesostigmata mites. A satellite sensor, by definition, records a "flat' pixel reflectance value and loses structural information, which may result in a poor discrimination of microtopography and its effects. Soil moisture indices (such as NDMI) may have limitations in peat bog environment, especially during dry conditions [32,62]. Crichton et al. [28] used Landsat ETM+ to map the spatial distribution of six ecohydrological classes on a lowland ombrotrophic peatland in the UK; their accuracy with Landsat ETM+ bands alone was 74%; including the texture brightness layer derived from band 5 increased accuracy of prediction.
Other option includes supplementing satellite imagery with airborne Light Detection and Ranging (LiDAR) data [31,63,64] or using the Synthetic Aperture Radar (SAR) backscatter in Sentinel-1 satellite [65], which adds information on surface texture and water table depth to increase accuracy.

Conclusions
To conclude, Landsat 8 OLI spectral indices representing productivity gradient and ground water level gradient were consistently correlated to local environmental gradients in Shichengskoe mire. Spectral indices were significantly correlated with mite biodiversity parameters, and the patterns were consistent with the ground data; the strongest relationship was between Oribatida biodiversity and the ExG index in 2016. Changing weather conditions can override local environmental gradients, as we have seen with lower abundance and diversity of Oribatida after the abnormally hot summer of 2019.
There is an opportunity to use freely available medium-scale remotely sensed data for soil biodiversity monitoring in space and time, once the links between biodiversity of specific taxa and land features which affect land reflectance values (e.g., soil type, soil moisture, vegetation) are established. Ground data should be used to validate the accuracy of predictions. Surface topography (e.g., digital elevation model) and texture (roughness, microtopography) data, if possible, should be used to supplement spectral data.