Biophysical variables have both direct and indirect effects on the uptake and release of carbon dioxide (CO2) within tundra ecosystems. Arctic landscape shows high levels of spatial heterogeneity. High spatial-resolution remote sensing data has the ability to capture the fine-grain spectral response of various biophysical variables at the landscape scale. To accurately model CO2 flux patterns using remote sensing data we first need to model the relationships between biophysical variables and their spectral response. In this study we model percent vegetation cover (PVC), aboveground biomass (AGB), and soil moisture using high spatial-resolution (IKONOS 4 m) normalized difference vegetation index (NDVI) values. At two non-overlapping Arctic landscape sites statistically robust landscape-scale sampling procedures were used to characterize the biophysical variables. NDVI values were extracted from IKONOS data, and linear bivariate regression models were calibrated and validated using a k-fold cross-validation technique. PVC and percent soil moisture produced the strongest and most consistent results (r2≥ .84 and .73, respectively). Analysis of covariance tested the use of common models for each site. The models were not coincidental—combining data from various sites should be done with caution—but illustrated parallelism in that NDVI responds to each biophysical variable equally, regardless of site.
Introduction
Global climate is changing; warming at high latitudes is expected to exceed global mean warming by more than 40% (IPCC, 2007). Satellite remote sensing data indicate that Arctic tundra vegetation is changing at different rates, sometimes rapidly, throughout the Arctic (Myneni et al., 1997; Jia et al., 2003; Bhatt et al., 2010; Raynolds et al., 2012). Ecosystem structure, function, and processes in Arctic regions of Alaska have already been affected by climate warming (Oechel et al., 1993; Serreze et al., 2000; Sturm et al., 2001; Hinzman et al., 2005; ACIA, 2005). Warming within tundra ecosystems has the potential to: (1) promote plant growth and sequestration of carbon from the atmosphere; and (2) increase soil microbial respiration rates, releasing additional carbon to the atmosphere (Stieglitz et al., 2000). If ecosystem respiration (ER) exceeds gross ecosystem production (GEP), a positive feedback mechanism could develop, thus intensifying global climate change (Nobrega and Grogan, 2008). GEP and ER have been shown to vary in relation to ecosystem patterns. This variance is a result of many factors, including vegetation structure (cover and biomass), soil organic matter, soil moisture, N and/or P availability, macro- and microtopography, air temperature, and thaw depth (Shaver et al., 2007; Dagg and Lafleur, 2011). Therefore, there exists a need to understand the spatial distribution of ecosystems and their associated biophysical patterns.
Remote sensing provides spatially continuous data depicting vegetation and terrain patterns at a range of spatial, spectral, and temporal resolutions. These data provide a unique opportunity for investigating the biophysical properties of vegetation over space and time (Tieszen et al., 1997; Stow et al., 1998, 2000; Laidler et al., 2008). Remotely sensed data can: (i) provide baseline data to delineate vegetation community patterns (Spjelkavik, 1995; Hope et al., 1995; Rees et al., 1998; Walker et al., 1982, 2005; Stow et al., 1989; Mosbech and Hansen, 1994; Muller et al., 1999); (ii) examine community structure (e.g., estimate aboveground biomass) (Hope et al., 1993; Spjelkavik, 1995; Shippert et al., 1995; Walker et al., 1995; Epstein et al., 2012); and (iii) predict CO2 flux patterns at a variety of spatial scales (Stow et al., 1993b; Ostendorf and Reynolds, 1998; McMichael et al., 1999; Shaver et al., 2007). Stow et al. (1998) stated that for carbon storage and flux patterns to be predicted from remote sensing data, measurements of these variables are required to calibrate and validate the models. Before this can be attempted, the following requirements first need to be addressed: (i) unique electromagnetic signatures need to exist and correspond to variations in biophysical variables (i.e., vegetation patterns, structure, and ecological site factors); and (ii) one or more models are needed to transform remotely sensed data into derivative values pertaining to the biophysical variables (Stow et al., 1998).
Understanding biophysical patterns using remote sensing data involves modeling relationships between spectral vegetation indices (VIs) and biophysical variables, i.e., how they vary across landscapes, and how these variations are related to vegetation composition, biomass, and ecological site factors (Walker et al., 1995; Boelman et al., 2003). It has been demonstrated for a range of vegetation communities that simple transformations of band reflectance are closely correlated with plant biophysical properties and are generally less sensitive to external variables such as solar zenith angle than to individual image bands (Laidler and Treitz, 2003).
The normalized difference vegetation index (NDVI) (Rouse et al., 1974) has been the most commonly employed spectral vegetation index (VI) in biophysical analyses, including those conducted in Arctic regions (e.g., Shippert et al., 1995; Walker et al., 1995; Rees et al., 1998; Stow et al., 2004; La Puma et al., 2007; Laidler et al., 2008). The index is derived from the difference in reflectivity of the land cover in the near-infrared (NIR) band, where vegetation structure reflects strongly, and the red (R) band, where vegetation absorbs strongly as a function of chlorophyll concentration. It is calculated as follows:
NDVI has shown to have correlations to aboveground biomass (AGB) and percent vegetation cover (PVC) in Arctic environments (Hope et al., 1993; Shippert et al., 1995; Boelman et al., 2003; Walker et al., 2003; Laidler et al., 2008; Epstein et al., 2012; Raynolds et al., 2012).
The availability of remote sensing data at a range of spatial, spectral, and temporal resolutions has increased greatly in the past decade, though the challenge of matching the appropriate scale of data to the appropriate ecological scale remains. The majority of remote sensing studies modeling biomass and carbon flux in the Arctic have applied coarse spatial resolution data most suited to regional-scale ecosystem assessments (Shippert et al., 1995; Stow et al., 1998; Walker et al., 1998, 2003; Jia et al., 2003). In the Arctic at the landscape scale (within a common climatic zone), the soil hydrological regime driven by topography is the principle control on the distribution of vegetation communities and is highly spatially variable (Bliss and Matveyeva, 1992; Walker, 2000; Nobrega and Grogan 2008). Higher spatial resolution NDVI measurements are therefore essential to estimate the response of each tundra vegetation type to changing climatic conditions (Boelman et al., 2003). The use of high spatial-resolution multispectral data (e.g. IKONOS data [4 m spatial resolution]) has been minimally investigated in Arctic environments (Laidler et al., 2008; Fuchs et al., 2009). The improvement in satellite spatial resolving power provides enhanced capacity for the estimation of biophysical variables and the prediction of carbon flux values (Laidler et al., 2008).
Warming and its associated changes have the potential to change the carbon status of the Arctic tundra from a sink to a source, but the regional uncertainty of this varies with landscapes (ACIA, 2005; Hinzman et al., 2005). Spatially continuous data characterizing vegetation and terrain patterns at regional and landscape scales can be acquired with satellite-based optical remote sensing and various spectral derivatives (Boelman et al., 2003; Laidler et al., 2008; Fuchs et al., 2009). With a greater understanding of spectral response to variations in biophysical variables (i.e., vegetation patterns, structure, and ecological site factors), remote sensing data may be able to predict and monitor changes in CO2 flux patterns over time (Stow et al., 1998). High spatial resolution remote sensing may be the appropriate data source for estimating the biophysical variables at the landscape scale (Boelman et al., 2003). Before we can address the prediction of carbon flux patterns from remote sensing data at the these scales, we need to first answer the following question: Can various biophysical variables such as peak season aboveground biomass, fraction of vegetation cover, and soil moisture be accurately modeled with high-resolution spectral vegetation indices? The following hypotheses were tested at two non-overlapping Arctic landscapes along a latitudinal gradient: (1) robust landscape-scale sampling procedures that are paired with high spatial resolution imagery will enhance the capacity for the estimation of biophysical variables; (2) biophysical variables of species composition (i.e., vegetation functional groups), PVC, AGB, and soil moisture are strongly correlated with NDVI derived from high spatial resolution optical data and can be modeled using linear regression; and (3) regression equations derived at different study areas will be statistically similar in that one equation can be applied to multiple locations.
Methods
STUDY AREA
For this study two areas were selected: (1) the Cape Bounty (CB) Arctic Watershed Observatory on the south-central coast of Melville Island, Nunavut, Canada (74°55′N, 109°35′W); and (2) Lord Lindsay River watershed, west of Sanagak Lake (SL) (70°11′N, 93°44′W) on Boothia Peninsula, also in Nunavut, Canada (Fig. 1). These sites were selected to compare biophysical variables across a latitudinal climate gradient and allow for a more robust examination of biophysical modeling efforts. The northern study site at CB is approximately 150 km2 in size and is comprised of two adjacent watersheds that drain into two separate lakes, and then south into Viscount Melville Sound. The area is underlain by steeply dipping sandstone and siltstone of the Devonian Weatherall, Griper Bay, and Hecla Bay Formations and overlain by late Quaternary glacial and marine sediments (Hodgson et al., 1984). The climate is cold throughout the year with a melt season extending from June through to August (Lamoureux and Lafrenière, 2009). The mean daily July temperature in 2004 was 3.1 °C with infrequent, and typically low intensity, rainfall. Low stratus clouds and fog are common during the summer months.
The southern study site at Sanagak Lake (SL) is confined to the regions between the Lord Lindsay River to the south and an adjacent tributary to the north. The dimensions of the study area are approximately 15 km long and the width varies from 500 m in the southeast to over 5 km in the northwest. The area is comprised of sandy glacial outwash plains and plateaus with evidence of former oxbow lakes and channels. The area is underlain by dipping limestone formations with extensive outcroppings of granitic rocks to the north. In 2005 the mean daily July temperature was 6.0 °C, with only trace precipitation through the entire growing season.
At each study area, continuous permafrost with an active layer of 0.5–1 m dominates, and vegetation cover varies by drainage conditions along a mesotopographic gradient (Billings, 1973; Walker et al., 2002). According to Walker et al. (2005), CB is located in bioclimatic zone B and has a vegetation classification of G2—graminoid prostrate dwarf shrub, forb tundra, while SL, in bioclimatic zone C, is classified as P1—prostrate dwarf shrub, herb tundra.
FIELD DATA
Field sampling coincided with the peak growing seasons at each site, approximately the last three weeks of July (2004 at CB, and 2005 at SL). Sampling procedures were designed to satisfy the requirements for remote sensing and ecological data analyses through the collection of a variety of in situ data. Large contiguous areas of distinctly homogeneous vegetation cover were identified: (i) initially through an unsupervised spectral classification of Landsat TM data; and (ii) subsequently by a visual inspection of the study area. Sample vegetation plots were defined within these large homogeneous areas. An attempt was made to have at least two plots of each vegetation cover type measuring one hectare (100 m × 100 m). Although this study uses high spatial resolution imagery (4-m spatial resolution IKONOS data), 1-ha plots were sampled to allow for “up-scaling” to mid-resolution imagery (i.e., Landsat ETM + ) in future studies.
FIGURE 1.
Study areas: (A) Cape Bounty Arctic Watershed Observatory, Melville Island, Nunavut, Canada (CB); and (B) Sanagak Lake, Boothia Peninsula, Nunavut, Canada (SL). IKONOS images are displayed as color infrared composites: i.e., bands 4 (near infrared), 3 (red), and 2 (green).
Each plot was sampled to characterize the community in terms of species composition, PVC, AGB, and soil moisture at peak growing season (Lewis, 1998). In this study aboveground biomass was used to determine the sample size. As the primary variable that requires the most post field processing, it was used to determine sample size for all associated variables. Small sample size has been an issue for many studies in remote locations, and in many instances has limited their statistical reliability (Shippert et al., 1995; Laidler et al., 2008; Chen et al., 2009). This study took a statistically rigorous approach to determining sample sizes for each vegetation plot. Sample sizes were predetermined using pilot study results from 2001 (SL) and 2003 (CB) and were then adjusted in the field based on the mean and variability of field weight biomass samples. The sample size for each community is based on the sample size required to achieve a 95% confidence level (p = 0.05) within ±20% of the mean aboveground biomass for each plot (Wein and Rencz, 1976). For each sample plot the required number of samples (N req.) for a given level of confidence can be computed as follows:
where t* = 1.96 (n = α); d is the desired margin of error (±20% of the sample mean); and S is the standard error. There is a different t distribution for every sample size. Since we do not initially know the sample size we estimate it to be very large and select a t* value close to infinity (t = 1.96 [p = 0.05] where n = α). Upon calculation of the sample size we get a better estimate of the degrees of freedom and a new t* can be obtained and N req. is recalculated. This procedure is repeated until N req. no longer changes (Wein and Rencz, 1976). For the spatial location of the samples, a stratified random sampling technique without replacement was implemented with an equal number of N req. samples randomly located within each quadrant (25 m × 25 m) of the 1-ha plots. Randomly generated co-ordinates were uploaded into a Garmin GPSmap 76 so they could be located quickly within a plot. GPS precision was high (i.e., within 3 m) due to excellent line of sight and satellite signal reception. Sample sizes at CB ranged between n = 12 and n = 80, while at SL the range was from n = 5 to n = 72. The variation in sample size related directly to the variability of AGB within the plot; homogeneous plots required very few samples, while heterogeneous plots required a much higher sampling effort.
A 50 cm × 50 cm (0.25 m2) quadrat was used to delineate each sample location within the 1-ha plots. For each location a species list was generated and cover abundance was estimated using a modified Braun-Blanquet scale. Cover was estimated for biotic species and abiotic components, such as rock and till (Lewis, 1998). In some cases, species identification was limited to the genus level (e.g., Carex spp.), and moss species were classified into Sphagnum spp. and non-Sphagnum mosses (Chapin et al., 1996). Total aboveground biomass was clipped within the quadrat to the soil or brown moss layer (Norum and Miller, 1984). A 5-cm3 soil sample was collected at every third quadrat (with a minimum of n = 12 per plot) to calculate gravimetric soil moisture. Biomass and soil moisture samples were weighed wet in the field and returned to the lab for oven drying.
SATELLITE REMOTE SENSING DATA
IKONOS multispectral data (4-m spatial resolution) were collected for Cape Bounty on 22 July 2004, corresponding to the vegetation sampling for that site. A 23 July 2001 IKONOS image was collected for Sanagak Lake. Though the image does not correspond with the sampling year, it is representative of the peak growing season. As well, there was no identifiable land cover change or ecological disturbance during this period. IKONOS data (i.e., image channels blue: 0.45–0.52 mm; green: 0.51–0.60 mm; red: 0.63–0.70 mm; and near infrared: 0.76–0.85 mm) were calibrated to top-of-atmosphere reflectance following procedures outlined by NASA (2002) and Taylor (2005). Full atmospheric correction was not possible due to the lack of atmospheric condition data for both sites at the time of image acquisition. The images were geo-referenced to Universal Transverse Mercator (UTM) coordinates to correspond to the 1:50,000 Canadian National Topographic Survey (NTS) for each of the study areas and were corrected using ground control points (GCPs) collected concurrently with the field data collection. The overall Root Mean Squared Errors (RMSEs) for each of the corrected images were less than 3 m horizontally. NDVI images were derived for each site, and individual sample location values were extracted. Sample NDVI values were interpreted from the adjacent cells (n = 12) using bilinear interpolation to decrease any effects of misregistration errors that result from matching imagery with sample locations (Fuchs et al., 2009).
REGRESSION CALIBRATION
Linear bivariate regression of Arctic biophysical variables and vegetation spectral indices, specifically NDVI, from various satellite sensors is the most commonly applied and easily understood method of analysis (Stow et al., 1993a; Shippert et al., 1995; Laidler et al., 2008; Chen et al., 2009). Biophysical samples and NDVI values were aggregated to the plot quadrant level for analysis with each site (i.e., CB and SL) being analyzed independently. All variables were tested using the Shapiro—Wilk (n < 50) or Kolmogorov-Smirnov test (n > 50) to see if they meet the assumption of a normal distribution (SPSS, 2011). With the exception of NDVI, the biophysical variables of aboveground dry weight biomass (AGB) and soil moisture were transformed using a natural log function (In), while PVC required a square root transformation (Tucker and Sellers, 1986; Walker et al., 2003; Raynolds et al., 2012). Linear bivariate regression analyses were performed with the NDVI values as both dependent and independent variable. In the calculation of spatial models, as in Laidler et al. (2008) and Shippert et al. (1995), NDVI is used as the independent variable even though NDVI would be considered a function of the biophysical variable.
Although a withheld-data validation technique would have been ideal, the relatively small sample size of the aggregated data were better suited to a k-fold cross-validation technique (Snee, 1977). The data are portioned into k-sized folds with subsequent k iterations of calibration and validation calculated. At each iteration, k-1 folds are used to calibrate the model while the remaining fold is held out for validation. A tenfold cross-validation of all data sets was performed using Matlab R2012a.
Many statistics express the error of estimations, each requiring statistical criteria to define their validity. Ease of interpretation, clarity of presentation, reliability, and robustness are all characteristics of a good error statistic (Tayman et al., 1999). Measures of error can be either scale-dependent or independent (Hyndman and Koehler, 2006). The advantage of a scale-dependent estimate of error, such as the commonly used RMSE, is that its scale is the same as the data. RMSEs are representative of the size of an “average” error. Non-scale-dependent error estimates adjust for the population size by using a percentage error. Since percentage errors are not scale-dependent, they can be used to compare performance across data sets (Swanson et al., 2011). As in Chen et al. (2009), two relative statistics, the Relative Mean Absolute Error (RMAE) and the Median Absolute Percentage Error (MedAPE), were selected as the non-scale-dependent estimates of error. RMAE is computed as the Mean Absolute Error (MAE) divided by the corresponding average value of the dependent variable and is more sensitive to errors if the dependent variable has a large value. In contrast, the value of MedAPE generally gives more weight to percentage errors corresponding to low values for the dependent variable (Chen et al., 2009; Swanson et al., 2011). The coefficient of determination (r2) is the most common statistic given describing the fit of an empirical model and can be interpreted as the proportion of the total variation in the dependent variable explained by the independent variable.
TEST FOR COINCIDENCE (ANALYSIS OF COVARIANCE)
With regression equations calculated for each study site, we then tested for the hypothesis of coincidence, i.e., that the calculated relationships for NDVI versus the biophysical variable are exactly the same. This analysis is referred to as Analysis of Covariance (ANCOVA) (Wildt and Ahtola, 1978). If the hypothesis of coincidence is accepted, i.e., the empirical models exhibit equal slopes (H0: β1CB = β1SL) and equal intercepts (Ho: β0CB = βoSL), we can pool the data from each site to derive a single model (Y = β0 + β1X + E) to describe the relationship at each sites. Next, the hypothesis of equal slopes is tested. The null hypothesis is that the slopes of the regression lines are the same. If the slopes differ significantly, then coincidence is rejected. If the slopes are the same, we then test to determine if the intercepts are the same; if they are, then the models are coincident. If the intercepts are significantly different but the slopes are equal, we have parallelism. In parallel equations the effect of the explanatory variable on the response variable is the same for each site, but the “base-line” values are different (when X = 0). If parallelism is rejected, then the explanatory variable and the response variable are responding differently at each site.
Results
BIOPHYSICAL RESULTS
In the two years of field sampling, 21 independent 1-ha vegetation community sample plots were established, and 885 vegetation samples were collected and processed. Overall, a total of 44 species/ genus types were identified. All SL sample plots (12 plots and 398 samples) and all but one sample plot (i.e., A-CB) at CB (9 Plots and 487 samples) achieved their sampling goal of ±20% of the mean for aboveground dry biomass (p = 0.05); A-CB sampled to within ±33% (p = 0.05). Comparing the two sites (Fig. 2), CB had a mean aboveground dry biomass of 1019 g m-2, while SL was 43% lower with a mean of 580 g m-2. Figure 2 illustrates that CB, along with having a higher mean, also displayed a much larger variability and was skewed; 50% of the samples fell between the median value of 810 g m-2 and the upper quartile of 1542 g m-2. SL had a much narrower range of values (310 g m-2 as compared to 986 g m-2 between the upper and lower quartile) and was only slightly skewed.
Mean plot level gravimetric soil moisture values for CB and SL are presented in Figure 3. To illustrate the range of values and to allow for comparison between sites, the plots were arranged along a moisture gradient; i.e., driest to wettest for each of the study areas. At CB and SL there appeared to be three main groupings of plots in regards to soil moisture: <20% (dry), 20–55% (mesic), and >55% (wet). Hereafter, figures depicting percent cover and functional type will be arranged according to moisture. Further, and based on nomenclature described by Walker et al. (2005), each plot is characterized in terms of soil moisture and bioclimatic zone (Table 1).
Percent vegetation cover by functional group for CB and SL is presented in Figure 4. The dominant functional groups at CB were bryophytes with mean coverage of 45.2%, followed by graminoids at 26.8%. Non-vegetated surfaces (i.e., a combination of rock, till, bare soil, etc.) averaged 28.9% cover at CB. Each sample plot at CB has some degree of exposed surface, ranging from a high of 77.8% at C-CB to a low of 1.4% at J-CB, with the drier sites having a greater amount of non-vegetated area. The driest plots at SL exhibited high percentages of non-vegetated surface area (e.g., 78.5% at X-SL), but these declined quickly as soil moisture values increased with many plots having <1% of non-vegetated surface area and an overall site mean of 17.7%. Shrubs dominated SL (35.6%), particularly Dryas octopetala, which was found in 82.7% of the SL samples and in each plot. Bryophytes and graminoids represented 24.3% and 21.5% of the ground cover, respectively, at SL. Shrub cover was limited at CB (2.4%), consisting mostly of Salix spp. Based on Walker et al. (2005), the species composition and PVC of dry plots at CB and SL are most closely classified as P1—Prostrate dwarf-shrub, herb tundra; and P2—Prostrate/Hemiprostrate dwarf-shrub tundra, respectively. Non-vegetated surfaces were dominant at each study area followed by bryophytes and forbs at CB and D. octopetala at SL. The mesic plots at CB (G1—Rush/ grass, forb, cryptogam tundra) had high PVC of Nostoc commune, a nitrogen-fixing cyanobacteria that forms a black crust on the soil (Liengen and Olsen, 1997). N. commune was present at SL mesic plots (G2—Graminoid, prostrate dwarf-shrub, forb tundra) but represented a much lower PVC; instead these plots were dominated by Dryas spp., Salix arctica, and Cassiope spp. The only vegetation class that was shared by CB and SL is that of the wetlands variety, i.e., W1—sedge/grass moss wetland (wet sedge). For a more detailed description of the vegetation classifications the reader is referred to Atkinson and Treitz (2012). It can be seen that these sites had similar functional composition of primarily bryophytes and sedges, though again there was some shrub cover at SL that was not present at CB. Total aboveground biomass generally followed the previously discussed moisture gradient, though it was much better defined at SL. The driest sites (P1, P2) with the most nonvegetated surfaces had the lowest average aboveground biomass, i.e., 265.9 g m-2 for C-CB and 125.9 g m-2 for X-SL. CB and SL had similar high biomass values in the wet sedge group with 2225.5 g m-2 and 2129.6 g m-2, respectively.
FIGURE 2.
Boxplot of aboveground dry biomass samples for CB and SL, showing the lower, median, and upper quartiles, minimum and maximum sample values, and outliers.
FIGURE 3.
Average gravimetric soil moisture for each plot at CB and SL, arranged from driest to wet for each site.
TABLE 1
Research plots (n = vegetation sample size) along a moisture gradient of dry to wet with dominant cover and associated vegetation classification (Walker et al., 2005). CB = Cape Bounty Arctic Watershed Observatory, and SL = Lord Lindsay River watershed west of Sanagak Lake. See text for explanation of classification codes.
RELATIONSHIPS BETWEEN NDVI AND BIOPHYSICAL VARIABLES
Linear bivariate regressions were performed to examine relationships between satellite data (i.e., NDVI) and biophysical variables in this environment in order to determine the potential of satellite data for extrapolating these variables over large areas (Fig. 5). All reported regressions were significant (p < 0.01) (Table 2). The relationship between NDVI and PVC had the strongest relationship at CB and SL with r2values of 0.88 and 0.84, respectively. The error statistics for SL were the lowest with RMAE and MedAPE at 13.1% and 12.0%, respectively, and slightly higher at CB (i.e., 15.6% and 16.5%, respectively). Though soil moisture was not a causal variable, there was a strong relationship with NDVI at both sites with r2values of 0.72 and 0.73 at CB and SL, respectively. The error statistics for both sites were again relatively low, albeit slightly higher at CB (∼24% for RMAE and MedAPE) versus 17.4% and 13.0%, respectively, at SL. Both sites had similar RMSE values of 0.045 and 0.040 for CB and SL, respectively. The results for total dry aboveground biomass for each site varied more than previous relationships observed. SL showed a strong relationship between NDVI and biomass with r2values of 0.76 and error statistics less than 18%. CB demonstrated a moderate relationship (i.e., r2= 0.55) and higher error statistics at 30.4% and 28.7.8% for RMAE and MedAPE, respectively. CB had several plot quadrats with low NDVI values but high dry biomass values; thereby reducing the strength of the relationship.
FIGURE 4.
Percent vegetation cover (PVC) for each plot at CB and SL, indicating the percent per functional group and the average aboveground biomass (AGB).
FIGURE 5.
Relationships of normalized difference vegetation index (NDVI) to the biophysical variables of AGB, PVC, and soil moisture.
TABLE 2
Bivariate linear regression (Y = β0 + β1X + E) results for biophysical variables and normalized difference vegetation index (NDVI).
The ANCOVA results for these three regression equations confirmed that the slopes for each regression equation were statistically similar (p > 0.05) (Table 2). However, when the intercepts of the equations were examined, it was found that they were significantly different (p < 0.05). These equations were parallel in that NDVI varied equally with each independent variable, be it PVC, percent soil moisture, or AGB at each site. The analysis of the intercepts indicated that there are differences in the baseline values for each variable at each site.
FUNCTIONAL GROUPINGS AND NDVI
Species cover data were aggregated into functional groups; i.e., graminoids, bryophytes, forbs, shrubs, and lichens. Additionally, non-vegetated surfaces, such as mineral soil, till, and rock, were aggregated to create a non-vegetated surface cover group. The results of the linear bivariate regression of the functional groupings are presented in Figure 6 and Table 3. The regression results of shrubs as well as lichens to NDVI did not return statistically significant results (p> 0.05) and so are not displayed. The best results were observed for non-vegetated surfaces and NDVI; i.e., CB and SL each had r2≥ 0.80 and low error statistics. When an ANCOVA was performed to compare the two sites the results demonstrated that the slopes and the intercepts for relationships at both sites were statistically similar (i.e., p-values > 0.05). Since the relationships at each site were statistically similar, the data from each were combined to create one regression, applicable to both sites. The non-vegetated surface regression was strong with an r2 of 0.81 (RMSE = 0.035; RMAE and MedAPE scores ∼18%). Forbs were also aggregated after ANCOVA results determined that the relationships for each site were statistically similar. However, the regression results were poor (i.e., r2= 0.08). The slope for the combined forb group was negative, indicating that as NDVI values increased the percent cover of forbs decreased. Graminoids and bryophytes showed stronger relationships at CB (r2= 0.72 and 0.80, respectively) than at SL (r2= 0.62 and 0.44, respectively). For graminoids the ANCOVA results also demonstrated that the slopes were statistically similar (p > 0.05) at CB and SL though the intercepts were not, again indicating that the response of NDVI to changes in graminoid cover was similar at both locations.
Discussion
BIOPHYSICAL SITE COMPARISON
Moisture regime appeared to be the principal control on vegetation distribution and abundance at CB and SL, each of which appeared to have a natural moisture gradient, at peak season, that generated three distinct moisture categories across the landscape: dry, mesic, and wet (Bliss and Matveyeva, 1992; Walker, 2000; Nobrega and Grogan, 2008). Variations along this gradient may not have been simply limited to vegetation composition and structure but extended to ecosystem function and processes; e.g., carbon flux (CO2, CH4) and nitrogen availability/mobilization (Dagg and Lafleur, 2011). In this study, soil moisture was examined during the peak growing season. However, it should be noted that the manner in which the soil moisture regime changed within the vegetation communities throughout the growing season (i.e., phenological cycle) may also have had important impacts on ecosystem form and function.
The examination of a latitudinal gradient can serve as a proxy for predicted changes to the Arctic climate, particularly since increases in temperature are expected to enhance plant growth (Stieglitz et al., 2000; Epstein et al., 2012; Elmendorf et al., 2012). When examining the AGB data for CB and SL, we observed that CB, the more northern and cooler of the two sites, had ∼30% higher mean AGB than SL, albeit not uniformly distributed. There was high variability in the AGB samples from CB and a larger amount of non-vegetated surfaces. SL had more uniform vegetation cover, lower percentage of non-vegetated surfaces, and was dominated by prostrate shrub growth, particularly Dryas octopetala. At both sites microscale topographic features (e.g., frost cracks) created moisture gradients that influenced the pattern and distribution of vegetation (Laidler et al., 2008). Plots that had a large percentage of non-vegetated surfaces (those classified as P1 and P2) exhibited high within-plot variability. Plant growth was restricted to frost cracks and depressions where organic matter accumulation increased nutrient availability and, with reduced wind speeds, enhanced moisture retention (Oberbauer and Dawson, 1992). With vegetation limited to small-scale features the impact on the spectral signature of these plots will often be dominated by the exposed soil. Though the PVC values for these dry sites at SL and CB are comparable, CB did have higher AGB in the microtopographic features. The frost cracks at CB contained higher amounts of moss, while the SL sites were limited to shrubs and lichen. The amount of moss increased the AGB values for those sites.
FIGURE 6.
Relationships between NDVI and PVC of functional groups: graminoids, bryophytes, forbs, and non-vegetated surfaces.
TABLE 3
Bivariate linear regression results (Y = β0 + β1X + E) for PVC1 and NDVI (by functional group).
Within Arctic ecosystems, the primary nitrogen input originates from biological nitrogen fixation by cyanobacteria such as Nostoc commune (Alexander, 1974; Lennihan et al., 1994; Solheim et al., 1996). The mesic plots of SL and particularly CB contained the most PVC of Nostoc commune. Though nitrogen was not examined within this study, further understanding of the connections between soil moisture, nitrogen, and vegetation cover at these and other sites and the potential of remote sensing could prove interesting.
Wet sedge (W1) was common to both study areas. These plots contained the highest levels of biomass and the highest levels of soil moisture (often saturated soils). Wet sedge meadows have been shown to be comparable in their CO2 uptake but because of lower overall ecosystem respiration they tend to be sinks for carbon (Nobrega and Grogan, 2008). However, these communities had the highest soil moisture levels, often a function of large, often-perennial snow banks located upslope from these communities (Lamoureux et al., 2006). CB appeared to have a higher area of wet sedge due to the macrotopography of the adjacent watersheds. This may have contributed to the overall higher levels of biomass at CB. Changes to the hydrological regime, through quicker melt of these snow banks will directly impact these communities, perhaps altering their carbon sink status.
BIOPHYSICAL AND FUNCTIONAL RELATIONSHIPS WITH NDVI
The strong regression relationships observed between biophysical variables, functional groups, and NDVI, along with low error statistics, suggest that the increased sample size and sampling procedures in conjunction with high resolution imagery can produce enhanced estimations of biophysical variables, thereby supporting hypothesis (1). PVC produced the strongest and most consistent results for both CB and SL with high coefficients of determination and low error scores (< 16.5%). Percent soil moisture also showed strong results supporting the connection between soil moisture and vegetation within this environment. The relationship between AGB at CB was the only relationship that did not follow the pattern of PVC and soil moisture, having similar modeled results for both study areas. CB had a lower coefficient of determination (i.e., r2= 0.55) and the largest error statistics for AGB. In a 1:1 comparison of the estimated and actual NDVI values, when AGB was used to predict NDVI, CB samples that were classed as dry were over predicted, while the mesic and wet samples were under predicted (Fig. 7, part a); such a trend did not exist with the SL samples (Fig. 7, part b). When NDVI was used as the independent variable for the derivation of spatial models, the error resulted in an underestimation of biomass in the dry plots. This is potentially linked to the higher proportion of exposed soil in these plots lowering NDVI values, while there can still be considerable biomass within the microtopographic features. In this case a soil adjusted vegetation index may prove more useful (Laidler and Treitz, 2003).
Overall, and across both study areas, the bivariate linear regression relationships for the transformed PVC, soil moisture, and AGB values and the high spatial resolution NDVI were strong with very low error statistics, thereby supporting hypothesis (2). However, the ANCOVA analysis for these variables did not show strict coincidence, thereby refuting hypothesis (3). However, the variables did show parallelism in that NDVI was responding to changes in the variables equally at CB and SL. This indicates that if baseline values are obtained for a site, changes to NDVI over time can be used to monitor changes for each variable. When NDVI is used as the independent variable, and parallelism exists, the intercept represents the value of the biophysical variable when NDVI = 0. That variable may not necessarily be equal to zero. The differences in intercept at each site reflected differences in the amount of aboveground biomass, PVC, or soil moisture (when NDVI = 0), and the reflectance of exposed soil to the overall NDVI signature. The intercept values, though similar, were statistically different, thus one equation for both sites was not appropriate. The finding that the equations demonstrated parallelism has not been shown in other studies and so combining data from various sites should be done with caution.
FIGURE 7.
1:1 comparisons of actual and modeled NDVI values from aboveground dry biomass (g m-2): (A) CB—dry sites are overpredicting NDVI values resulting in underestimates of biomass for dry plots and vice versa; (B) SL—there is a stronger correspondence between modeled and actual NDVI.
PVC and percent soil moisture images modeled from IKONOS NDVI data for CB and SL are presented in Figure 8 along with the average modeled values for each study area in Table 4.(Shippert et al., 1995; Laidler et al., 2008). Modeling PVC and soil moisture for large areas created visualizations of the overall vegetation structure and pattern. However, it also depicted relationships between biophysical variables that would otherwise be difficult to discern. The findings of this study support Stow et al.'s (1993a) suggestion that high spatial resolution (i.e., <10 m) data should strengthen NDVI correlations to biophysical variables. These improved spatially explicit models can improve: (i) the identification of initial conditions for patch-scale models by inventorying landscape conditions and their relative proportions; (ii) the stratification of landscapes into relatively homogeneous response units for spatially distributed modeling of material and energy transport; (iii) the extrapolation of model simulations by mapping areas that are potentially sensitive to particular disturbances; and (iv) the assessment of landscape- and regional-scale model simulations by comparative spatial pattern analyses (Stow et al., 1993a; Laidler et al., 2008). In addition, the methodologies and results of this study could be paired with in situ measures of CO2 to better predict and monitor changes in CO2 flux patterns at landscape scales during the peak growing season. Further, these data (and methodology) can be integrated into regional-scale studies and coarser-resolution imagery to better understand the scaling of biophysical variables and CO2 flux patterns.
The relationships observed between NDVI and functional groups were mixed. Overall, hypothesis (2) is supported for the most prevalent functional groups: i.e., graminoids and bryophytes, whereas hypothesis (3) is refuted for all but non-vegetated surfaces. Lichens and shrubs did not return statistically significant results and so were not reported. At CB there was very low shrub cover, while at SL there was almost ubiquitous cover of Dryas spp. The lack of a gradient at either site contributed to the poor regression results. The lack of significance with lichen may be linked to the limited contribution that these groups made to the NDVI values (i.e., spectral response) (Shaver et al., 2007). Graminoids and bryophytes had the strongest relationships at CB where these two functional groups were the most abundant, along with non-vegetated surfaces. SL also had significant relationships with these functional groups. The PVC of graminoids and bryophytes at both sites increased with soil moisture. Graminoids showed parallelism for both sites, while bryophytes did not, indicating that the response of NDVI at CB and SL to changes in bryophyte coverage was statistically different. Forbs had a very weak but negative relationship to NDVI at both sites. Forbs are generally found in dry, barren tundra where there is little AGB. Though it exhibited a very weak relationship, ANCOVA did show coincidence for the NDVI and forb functions in that a single function could be applied at both study areas. Chen et al. (2009) developed field metric methods for assessing AGB for both graminoids and bryophytes; whereas the equations presented here have lower RMAE and MedAPE values, indicating the potential of high spatial resolution remote sensing data for modeling these variables. The strongest result was observed for nonvegetated cover. ANCOVA for this quasi-functional group showed the relationship was similar for both sites, indicating the models were transportable between the two study areas. Non-vegetated surfaces are different from PVC in that with a layered structure, PVC often has values greater than 100%, while non-vegetated surfaces cannot exceed 100%. Looking at the combined equation for this variable illustrates that NDVI responded equally at both sites to changes in vegetated cover. This equation has value in monitoring, as changes to NDVI in sites that have exposed soil may indicate increased plant growth and overall cover. Additionally, when modeling carbon flux, the amount of exposed soil is an important component, particularly when examining the contribution of soil respiration to ecosystem respiration and net ecosystem exchange (Nobrega and Grogan, 2008).
Conclusions
The Arctic is changing, and with it the composition, structure, function, and processes of terrestrial Arctic ecosystems (ACIA, 2005; Hinzman et al., 2005; Oechel et al., 1993; Serreze et al., 2000). Vegetation is both an integrator and indicator of climate and the properties of the ecosystem. In this study, NDVI data, derived from high spatial resolution IKONOS imagery, paired with appropriate and statistically based sampling techniques produced strong empirical models of biophysical variables at landscape scales, thus illustrating the value of remotely sensed data as a tool for evaluating changes within the tundra ecosystem. At CB and SL, similar but disparate locations along a latitudinal gradient, vegetation communities were strongly aligned along a peak-season moisture gradient. With each site having differing vegetation communities and bioclimatic subzones (Walker et al., 20005), strong linear bivariate regression equations for PVC, AGB, and gravimetric soil moisture were produced. Though most of the equations were not transferable between these specific study areas, they exhibited parallelism, indicating that NDVI responded to changes in the biophysical variables equally at these two locations (i.e., across 5° of latitude change), regardless of vegetation community type and bioclimatic subzone.
These models address two of the prerequisites established by Stow et al. (1998) for predicting CO2 flux and storage rates: (i) the existence of unique electromagnetic signatures that correspond to variations in biophysical variables; and (ii) models that transform remote sensing data into derivative values related to biophysical variables. What is needed further are carbon flux rates and storage amounts to calibrate and validate spatial models. The results of this study indicate a strong relationship between vegetation composition, structure, and function with peak-season soil moisture regime. The soil moisture regime is highly variable both spatially and temporally through the melt season. It has been established that soil moisture is related to available nitrogen and carbon flux (Nobrega and Grogan, 2008; Dagg and Lafleur, 2011). Hence, further research into the relationships between biophysical variables, remote sensing, and how they change within and between seasons should be further explored along an extended latitudinal gradient in order to better understand and predict the impacts of a warming climate at high latitudes.
FIGURE 8.
Spatial models of CB (upper) and SL (lower) for PVC (left) and soil moisture (right). Fine detail can be seen in each high-spatial resolution (4 m) image. Dense and multi-layered vegetation corresponds with wet drainage areas in each image.
TABLE 4
Average values and standard deviation, derived from the spatial models for the entire image areas of CB and SL.
Acknowledgments
The authors would like to gratefully acknowledge financial support from ArcticNet, International Polar Year (IPY), Northern Science Training Program (NSTP), NSERC and Queen's University. Logistical support provided by the Polar Continental Shelf Project (PCSP) was instrumental in supporting this research at two distant study areas. The authors would like to thank the Drs. Scott Lamoureux and Melissa Lafreniere for the support of this research at Sanagak Lake and the Cape Bounty Arctic Watershed Observatory. The principal author would like to thank field assistants Jake Wall, Freyja Forsyth, and Shanley Thompson for their sampling assistance in the field.
References Cited
lesund area, Svalbard.
Polar Record
, 31: 257–269. Google Scholar
