The structure and function of oak Quercus spp. savanna ecosystems in the North American Midwest were originally maintained by an active disturbance regime (often fire). Subsequent reductions in the frequency of disturbance after European settlement have facilitated rapid, regional conversion of these ecosystems to more closed-canopy forest. Hence, regional-scale management strategies are now needed to restore critical spatial gradients of light, temperature, soil moisture, and soil organic matter for recovery and sustenance of the unique mosaic of understory grass and forb species assemblages that define oak savannas. Tree species composition, distribution, mortality, basal area, and canopy cover are important forest structural parameters that are intrinsically linked to oak savanna restoration ecology. In this benchmark study, we seek to determine whether Landsat-based monitoring protocols can be developed as a tool to guide and monitor regional-scale restoration and management efforts. Using the Sherburne National Wildlife Refuge in central Minnesota as a test case, ground-based forest-structure data were collected and used with multitemporal Landsat sensor data and iterative exclusion partial least-squares regression to calibrate six predictive overstory structure models. Model calibrations produced moderate- to high-accuracy results with respective adjusted coefficient of determination and root mean-squared error values as follows: 0.859, 9.3% (canopy cover); 0.855, 2.95 m2 ha−1 (total basal area); 0.741, 11.6% (red oaks relative basal area); 0.781, 11.9% (bur oak relative basal area); 0.861, 3.20 m2 ha−1 (living oak basal area); and 0.833, 9.1% (dead oak relative basal area). We used the resulting structure models for the Sherburne test site to demonstrate how these data could be applied to help managers prioritize areas within management zones for restorative treatments. Although our Sherburne oak savanna test ecosystem is small (12,424 ha) compared with the size of a full Landsat scene (3.4 million ha), resulting structure models can be extended to the whole Landsat scene, which demonstrates how such modeling protocols can be used for repeated (e.g., annual to decadal), regional-scale analysis and assessment to improve management, planning, and implementation of oak savanna restoration efforts elsewhere.
The specific structure and function of oak Quercus spp. savanna ecosystems in the North American Midwest were maintained by frequent disturbances (often fire) in pre-European settlement times (Curtis 1959). Subsequent reductions in the frequency of disturbance after European settlement have facilitated rapid conversion of these ecosystems to more closed-canopy forest (Curtis 1959; Vogl 1974; White 1983; Faber-Langendoen and Tester 1993). Combined with land-use conversion and fragmentation, once abundant oak savannas are now rare (<0.02% of presettlement area) in this region (Nuzzo 1986), which increasingly threatens the survival of many savanna-specific ground-layer species (Packard 1988, 1993; Leach and Givnish 1996, 1999; Peterson et al. 2007). In the absence of a robust, natural disturbance regime among remaining oak savannas, prescribed burning as well as mechanical and chemical tree thinning are frequently recommended (often in combination) to control both overstory tree-cover composition and to limit woody plant encroachment including shade-tolerant and fire-intolerant tree species regeneration (Bowles and McBride 1998; Boyles and Aubrey 2006). Regional, cross-ownership management strategies of this nature are needed to restore critical spatial gradients of light, temperature, soil moisture, and soil organic matter (Ko and Reich 1993; Breashears et al. 1997) for recovery and sustenance of the unique mosaic of understory grass and forb species assemblages (Bowles and McBride 1998; Grundel et al. 1998; Bader 2001) and faunal associations that depend on them (Boyles and Aubrey 2006; Peterson et al. 2007; Olechnowski et al. 2009).
Tree basal area (BA, computed as the sum of cross-sectional area of trees measured at 1.37 m above ground per unit area) and canopy cover (CC) are common forest structural and ecological attributes that together affect ground-layer light environment in all forest ecosystems, including oak savannas (Apfelbaum and Haney 1991; Peterson and Reich 2001; Peterson et al. 2007). Basal area is typically used in forest management to describe occupancy of trees and to estimate wood volume (Ginrich 1967). However, BA is also an important ecological attribute used to estimate forest productivity (Burrows et al. 2003), understand ecosystem structure and function (Whittaker et al. 1974; Marshall and Waring 1986; Pacala and Deutschman 1995), as a surrogate for tracking forest carbon (Brown et al. 1989; Phillips et al. 1998), and as a metric to understand the ecology and habitat requirements of forest fauna (Blais 1958; Niemi and Hanowski 1984; Pastor et al. 1998; Sharov et al. 1999). Similarly, CC is used to distinguish different plant communities, assess forest-floor microclimate and light penetration (Jennings et al. 1999; Peterson et al. 2007), and as an indicator used for the conservation of endangered species (Packard 1993; Grundel et al. 1998).
Each of these forest parameters are routinely measured in situ, but are capable of being mapped regionally with the aid of satellite sensor data, such as those routinely acquired by Landsat. For local-scale monitoring and mapping efforts, it may be possible to reliably map some forest-structure parameters, such as CC (Liknes et al. 2010) and mortality (Allain et al. 2011), using National Agriculture Imagery Program (NAIP) data (http://datagateway.nrcs.usda.gov/). The NAIP effort began as a pilot program in 2001 and has been tasked with acquiring 1–2-m spatial-resolution imagery on a near-annual basis since 2003 for the continental U.S. NAIP's natural color imagery (three bands) is always collected during the peak agricultural growing season, with an option to add a near-infrared band if funding is available in a particular year. The true-color, county-wide NAIP mosaics are made freely available for download, but the complimentary near-infrared imagery, if available, is sold at a fee that scales with data volume requested.
Although NAIP is an excellent data source for many forest applications, regional discrimination of oak species is more complex as compared with CC or mortality mapping. The optimum time for spectral discrimination between hardwood tree species, including oaks, is during autumn leaf senescence when leaf pigmentation differences are at a maximum, which makes discrimination easy and unequivocal (Schwaller and Tkach 1985; Wolter et al. 1995). During peak summer months, red oaks and bur oaks are effectively indistinguishable in the electromagnetic wavelengths covered by NAIP. Compounding this, NAIP data are a composite of imagery collected by different contractors (varying sensor formats and acquisition methods), on different days (potential atmospheric differences), and at different times of day (shadow fraction differences). Collectively, this produces spectral artifacts in the NAIP mosaic products that require unique handling on a mosaic-to-mosaic basis (Allain et al. 2011), which effectively precludes automated, regional forest composition mapping that has become a hallmark of Landsat (Schriever and Congalton 1995; Wolter et al. 1995; Vogelmann et al. 2001; Reese et al. 2002; Huang et al. 2010). Although spatially coarse by comparison, Landsat data are frequently described as “synoptic”—capable of seeing large areas at once (one scene is 185 km × 185 km). Moreover, Landsat data are freely available and have spectral, radiometric, and temporal resolution (16-d repeat cycle) characteristic specifically designed for consistent, regional mapping and monitoring efforts. Hence, we explore this potential for regional mapping and monitoring of oak savannas in this research.
Spatially explicit estimates of BA and CC are increasingly sought for landscape analyses to inform regional-scale management decisions (He et al. 1998; Mäkelä and Pekkarinen 2004; McRoberts et al. 2007; Wolter et al. 2008). Satellite imagery can meet this need by providing spatially explicit, multispectral, reflectance measurements that can be related to tree species composition, tree size, stem density, and canopy cover (Franklin and McDermid 1993; Xu et al. 2003; Wolter et al. 2008; 2009; Wolter and Townsend 2011). Here we combine field measurements of oak BA and CC with satellite image data from the Sherburne National Wildlife Refuge (SNWR) in central Minnesota (Figure 1). The SNWR was chosen as a proof-of-concept test ecosystem to demonstrate the capability for repeated (annual to decadal), regional-scale mapping of oak savanna structure. We used iterative exclusion partial least-squares (xPLS) regression (Wolter et al. 2012; Singh et al. 2013) to evaluate 94 Landsat Thematic Mapper (Landsat-5) predictor variables to answer three main questions. First, can Landsat imagery be used to accurately predict critical aspects of oak savanna structure: forest composition, BA, CC, mortality? Second, given a set of hypothetical criteria, how much of the oak savanna within our proof of concept site may benefit from restorative treatments? Third, are Landsat data used in this fashion appropriate for guiding regional-scale management to prioritize and track prescriptive restoration treatments?
Oak savanna study site
The 12,424 ha (30,700 ac) SNWR is located on the gently rolling to flat Anoka Sandplain in Central Minnesota (Figure 1; N45.49° Lat., W93.73° Long.). Oak savanna habitat at SNWR is a fire-dependent, dynamic community characterized by scattered trees, or groves of trees, mostly consisting of oak species with absolute ranges in canopy cover from <1 to 100% and BA between <1–35 m2 ha−1; but these values are more typically between 33 and 78% for CC and 8–23 m2 ha−1 BA (Figure 2). Canopy cover can vary at small scales (stand-level) where areas of both scattered trees and areas with groves of trees are present within a single stand. Red oak Quercus rubra or northern pin oak Quercus ellipsoidalis, which commonly hybridize in this area (Swain 1972; Tester 1989), are typically the dominant tree species interspersed with bur oak Quercus macrocarpa (Law et al. 1994; Buchanan 1996). The oldest living trees on the refuge date back to ∼1800 and consist mainly of bur oak (Kipfmuller and Hepola 2007). Understory vegetation consists of both native grasses (25–100% cover) and forbs (5–50% cover; MNDNR 2005). Northern pin oak is occasionally present as a secondary tree species in the overstory or in the shrub layer. The shrub layer is typically composed of American hazelnut Corylus americana or oak sprouts and is patchy to continuous (30–100% cover). The ground layer is composed mostly of native forbs and woody plants characteristic of oak forests with some prairie and savanna forbs and graminoids as well (Wovcha et al. 1995). The spatial extent of prairie openings within this oak savanna ecosystem is <30% (USFWS 2005).
Savanna ecosystems were sampled in 2009 (60 plots) by the authors and other SNWR personnel (Figure 1). Field-plot data consisted of five variable-radius (Grosenbaugh 1952) subplots: one located at the intersection and one at each of the four end points of two crossing, orthogonal, 60-m transect lines placed near the center of large (≥4.4 ha), homogenous stands (Figure 3). We determined plot locations randomly using Geographic Information System (GIS) software and we used a digital map of relic oak savanna from the original, presettlement vegetation map of Minnesota compiled in 1930 (Marschner 1974) to stratify. These stratified random locations all fell in forested zones (BA >2.4 m2 ha−1), which resulted in the following numbers of plots within three BA ranges: 21 low-BA plots (range: 2.4–15.6 m2 ha−1), 16 medium-BA plots (range: 16–23 m2 ha−1), and 23 high-BA plots (range: 24–34.4 m2 ha−1). No two ground plots were closer than 60 m apart. Sufficient stand size and homogeneity assured that stand edge effects were minimized during analysis and that image georectification (image-to-earth) and coregistration (image-to-image) errors were inconsequential. We navigated to each plot center to within 3 m (2dRMS) using a Wide Area Augmentation System (satellite-based, real-time, differential correction) –enabled Garmin eTrex GPS receiver (Garmin International, Inc., Olathe, KS). We estimated overstory tree CC for the plot using 121 densitometer measurements (CC% = canopy hits/121) made in 1-m increments along the two plot axes (Figure 3). A densitometer is a hand-held instrument that allows one to view a discrete, vertical point in the canopy above a ground position, where observations are recorded as either canopy or sky. Living and dead overstory-tree-species BA, dominated by oak subgenre Leucobalanus (white oak group—mostly bur oak at SNWR) and Erythrobalanus (red oak group—northern red and pin oak), was measured at each of the five subplots using a metric basal area factor- (BAF-) 2 prism. Hereafter, members of oak subgenus Erythrobalanus are referred to simply as red oaks and Leucobalanus as bur oak. Averaged BA metrics from the five subplots became the respective BA estimates for the full plot. We performed conversion of bur oak BA, red oaks BA, and dead oak BA to relative BA (RBA) by dividing the respective partial BA measures by total plot BA to convert m2 ha−1 to a percent of total BA (bur RBA, red RBA, and dead RBA). We then used total BA, living BA, bur RBA, red RBA, dead RBA, and CC data as dependent variables to develop regression models to estimate and map these parameters using Landsat sensor data (Table S1).
Landsat-5 Thematic Mapper (hereafter referred to as Landsat) data have 30-m nominal pixel resolution in six optical bands (blue 450–520 nm [B1], green 520–600 nm [B2], red 630–690 nm [B3], near-infrared 760–900 nm [B4], and two shortwave infrared [SWIR] bands of 1550–1750 nm [B5] and 2080–2350 nm [B7]). We acquired Landsat images (one winter, two spring, one summer, and two autumn) for this study based on temporal proximity to field data collection (summer 2009). One older image (7 October 2003) was acquired to provide contrast between spectral signatures of living oak foliage from 2003 and signatures of leafless oaks that had recently died (prior to 2009). This basic strategy using satellite sensor data capturing leaf-on versus leaf-off (before and after) conditions is routine among studies involving aspects of forest change or mortality (Vogelmann and Rock 1989; Wolter et al. 1995; Collins and Woodcock 1996; Wolter and Townsend 2011). We included winter imagery because the spectral manifestation of tree shadows on snow ground cover was found to empirically relate to overall forest BA in this region (Wolter et al. 2012). Moreover, winter imagery has been used to accurately discriminate forest from nonforest (Wolter et al. 2008). Hence, we used a winter Landsat image (10 March 2008) to calculate the Tasseled Cap brightness index (Crist and Cicone 1984) for distinguishing forest from nonforest prior to structure model development according to Wolter et al. (2008). We selected October imagery to capture peak or near-peak autumn foliage coloration prior to leaf drop as an aid to distinguish oak species (Wolter et al. 1995). Senescent bur oak leaf coloration (typically pale brown) differs substantially from members of the red oaks group (e.g., Quercus velutina, Q. rubra, Q. ellipsoidalis, etc.), which are typically more reddish during autumn senescence because of greater expression of anthocyanin pigments (Boyer et al. 1988; Miller et al. 1991; Lee and Gould 2002). We included April imagery to capture potential leaf-flush timing difference between oak subgenres, while summer imagery provided stable growing-season conditions as spectral contrast to the autumn and spring imagery. Henceforth, the seven Landsat images are referred to using letter codes (Table 1): J (24 January 2009; Figure S1), P (14 April 2009; Figure S2), U (1 June 2009; Figure S3), A (17 August 2008; Figure S4), OF (4 October 2008; Figure S5), OS (7 October 2003; Figure S6), and OT (20 October 2008; Figure S7).
We downloaded all Landsat images from the U.S. Geological Survey's Earth Resources Observation and Science Center web site (source: http://glovis.usgs.gov/) in Universal Transverse Mercator zone 15 coordinates. These 30-m Landsat images are precision-orthorectified (corrected for image horizontal displacement due elevation variation) and georectified by Earth Resources Observation and Science Center using Global Land Cover Facility's GeoCover data (source: www.landcover.org). Landsat images used in this study exhibited excellent pixel coregistration with each other. Each image was converted to top of atmosphere reflectance using published sensor-calibration coefficients (Thome et al. 2004). Because topography is known to influence Earth surface reflectance (Olyphant 1984; Dozier 1989; Culvenor and Coops 1999), all images were topographically adjusted with a lambertian correction methodology (Riaño et al. 2003) using a 30-m digital elevation model (source: http://ned.usgs.gov/) and scene-specific sun position information (Table 1).
In addition to each image's six reflective bands, established indices (simple ratio [SR], Jordan 1969; soil adjusted vegetation index [SAVI], Huete 1988; moisture stress index [MSI], Rock et al. 1986; shortwave infrared visible ratio [SVR], Wolter et al. 2008; Autumn Index [AI], Wolter and Townsend 2011; and Tasseled Cap transformations (brightness, greenness, and wetness; Crist and Cicone 1984) were calculated and used as predictive variables of forest structure in regression analysis (Table 2). These indices are known to exhibit sensitivity to tree species composition, forest BA, tree canopy cover, and tree mortality (Crist and Cicone 1984; Rock et al. 1986; Xu et al. 2003; Wolter et al. 2009; Wolter and Townsend 2011).
Partial least-squares (PLS) regression is frequently used in forest ecology studies to calibrate models of biochemical and biophysical forest parameters using full-spectrum (i.e., hyperspectral) imagery (Coops et al. 2003; Townsend et al. 2003; Martin et al. 2008), multispectral imagery (Wolter et al. 2008, 2009, 2012), and multisensor data fusion (Wolter and Townsend 2011). Partial least-squares regression has become routine in such studies because the algorithm does not assume zero error in the predictor data—often falsely assumed in remote sensing-based studies (Curran and Hay 1986)—and addresses the problem of collinearity (dependence) among multiple independent and dependent variables (Helland 1988). For a detailed discussion of the theory and application of PLS regression, see Geladi and Kowalski (1986); Wold et al. (2009); and Wolter et al. (2008).
Partial least-squares regression was recently adapted to run in an iterative fashion to eliminate weak predictors and to find the most parsimonious model (i.e., fewest predictor variables) that accounts for the largest amount of variation among dependent variables (Wolter et al. 2008, 2009; Wolter and Townsend 2011). This variable-reduction strategy is now known as iterative exclusion PLS or simply xPLS (Wolter et al. 2012; Singh et al. 2013). Partial least-squares regression differs from stepwise and Akaike's Information Criteria or Bayesian Information Criteria methods in that it uses cross-validation to reduce the importance of (but not eliminate) intra- and intercorrelated variables. What distinguishes xPLS from PLS is that xPLS governs the execution of thousands of PLS regressions in an iterative fashion such that predictors showing no or low sensitivity to the dependent variable(s) are completely eliminated; iterations continue until model improvement ceases (reductions in model root mean-squared error; RMSE). Wolter et al. (2008, 2009) and Wolter and Townsend (2011) used xPLS regression with field-plot data and multispectral, multisensor satellite-image predictor variables to calibrate models for estimating forest composition and structure.
In this study, we use xPLS regression to evaluate 94 initial multispectral image predictor variables (Table 2) to calibrate models for estimating CC (Figure S8), total BA (living + dead, Figure S9), living BA (Figure S10), and the relative BA (RBA) of bur oak (Figure S11), red oaks (Figure S12), and dead oak (Figure S13), where RBA = specific BA/total BA. Basal area measurements include predominantly oak species (bur, northern red, and northern pin) with an assortment of other forest species, such as Populus tremuloides, Acer saccharum, Fraxinus pennsylvanica, Ulmus americana, Ostrya virginiana, Prunus serotina, with low, overall, relative occurrence at SNWR (Kipfmuller and Hepola 2007).
Spatial autocorrelation is known to weaken some statistical tests by lowering the effective degrees of freedom (Koenig 1999). If present, spatial autocorrelation among regression residuals produces spatially dependent trends that may violate basic assumptions of parametric statistical tests. However, model coefficients of determination (R2) or F-statistics are not rendered invalid (Kerr et al. 2001). We performed isotropic spatial autocorrelation analyses on the regression residuals for total forest BA to determine whether critical levels of spatial autocorrelation were present.
Leave-one-out cross-validation (Gong 1986), using the prediction residual sum-of-squares prediction residual sum-of-squares procedure (Myers 1986), was used to validate all model calibrations using two-thirds of the field-plot data (n = 40). This method of model validation is analogous to applying the equation to an independent data set as the prediction residual sum-of-squares residual is obtained using observations that are excluded during equation derivation (Sun et al. 2003; Woods et al. 2008). However, independent model validations were also performed using plot data (n = 20) that were withheld from the xPLS model-calibration process.
Forest-structure model calibration and validation
The isotropic range of spatial autocorrelation among the regression residuals of total forest BA was 30.0 m. Hence, spatial dependence among regression residuals was not a factor in these analyses, because the minimum ground-plot separation distance was 60 m. Leave-one-out cross-validation (henceforth referred to as simply cross-validation; Table 3) of forest-structure models calibrated using ground-plot data and multitemporal Landsat imagery resulted in adjusted coefficients of determination, R2, that ranged between 0.741 for red oaks RBA (RMSE = 11.6%) and 0.861 for living oak BA (RMSE = 3.198 m2 ha−1). Cross-validation results for bur oak RBA and red oaks RBA models were similar, with bur oak RBA slightly better (adj. R2 = 0.781, RMSE = 11.9%) than red oaks RBA. Living BA, CC, and total BA were the three best models according to cross-validation (adj. R2 = 0.861, 0.859, and 0.855; RMSE = 3.198 m2 ha−1, 9.3%, and 2.948 m2 ha−1, respectively). While the living oak BA model had the highest adjusted coefficient of determination overall, the RMSE (3.198 m2 ha−1) for this model was slightly greater than that reported for total BA (Figure 4). Dead oak RBA calibration and cross-validation resulted in estimation errors of <10% (Figure 4), while adjusted R2 (0.833) was intermediate among the other model calibrations (Table 3).
Independent validation of the forest-structure models using data from 20 plots withheld from the xPLS calibration procedure show slight deviations from leave-one-out cross-validation results (Table 3). The three RBA models (bur oak, red oaks, and dead oak) showed similar error levels between the two validations used, but the adjusted R2 values under the independent validation for the red oaks RBA model was 9.1% greater and the dead oak RBA model 7.2% lower than that reported under cross-validation. However, the adjusted R2 values under each validation method were similar for the bur oak RBA model (Table 3). Independent model validation showed lower adjusted R2 values (−7.0 to −14.2%) than cross-validation for living BA, total BA, and CC accompanied by slight to modest increases in error levels: CC, +4.9%; living BA, +0.69 m2 ha−1; and total BA, +1.39 m2 ha−1.
Model reduction and important image-based predictor variables
Of the original set of 94 multispectral, image-based predictor variables (Table 2), 23 were retained for final model calibrations. This represents a model reduction of 75.5% over the initial full model. The six dependent structure variables (ground-plot data) were all placed together in a matrix and regressed simultaneously against the independent, image-based, predictor variable matrix (Wolter et al. 2008; Wold et al. 2009); therefore, all of the 23 retained image variables are used in each of the six structure models. However, individual predictor variable coefficients vary in sign and magnitude among the different models, which gives an indication of the relative importance of the retained predictors among the six respective structure model calibrations (Figure 5). Image predictor variables are reported below with a + or − sign preceding the variable name to indicate direction of coefficient loading.
For CC, 1 June 2009 SVR (−SVR_U) and 4 October 2008 Tasseled Cap wetness (+OF9) were the two strongest image predictors, followed by 14 April 2009 Tasseled Cap greenness (−P8) and 1 June 2009 Tasseled Cap brightness (−U7). In terms of sign and magnitude, these top two variables were also the leading variables for estimation of total BA, living BA, and bur oak RBA (Figure 5). The top two predictor variables for red oak RBA (northern red and pin oak), were 20 October 2008 Tasseled Cap greenness (+OT8) and 14 April 2009 Tasseled Cap greenness (+P8), followed by 20 October 2008 shortwave infrared band 6 (−OT6) and 1 June 2009 near infrared (−U4). The strongest loadings for calibrating the dead oak RBA model were +SVR_U and −OF9, which were similar to that observed for other models, but with opposite signs (Figure 5). The next two strongest predictor variables were 1 June 2009 Tasseled Cap brightness (+U7) and 7 October 2003 moisture stress index (−MSI_OS). Of the stronger variables from the set of 23 for calibrating the dead oak RBA model, three were from the 7 October 2003 image (−MSI_OS, +OS9, and –SVR_OS). These 2003 image variables were stronger for calibrating dead oak RBA than for any other structure model (Figure 5).
Current oak savanna structure at Sherburne NWR
In all, 2,898 trees were tallied within our 60 field plots, with red oaks consisting of 33.6%, bur oak 39.7%, other species 4.0%, and dead oaks 22.6%. Most of the dead trees tallied were red oaks (68.2%), followed by bur oak (22.9%) and other species (8.9%). Within the 12,424 ha of land at SNWR, ∼45.3% supports some level of tree cover. Across 10 CC intervals (each representing a range in CC of 10%) between 1 and 100 (referred to hereafter as CC classes), the 41–50% CC class encompasses the greatest proportion of oak savanna at SNWR, with a mean total BA = 11.7 m2 ha−1 (Figure 6). Greater proportions of oak savanna at SNWR occur in the next three higher CC percent classes (51–60 CC class, 14.3%; 61–70 CC class, 13.4%; and 71–80 CC class, 12.8%) than that falling in the lower CC classes (11–20 CC class, 4.2%; 21–30 CC class, 9.0%; and 31–40 CC class, 13.4%). Bur oak RBA exceeds that of the red oaks RBA (northern red and pin) within CC levels from ∼21 to 57% (total BA: 5.7–15.5 m2 ha−1); otherwise, red oaks RBA is greater than that of bur oak (Figure 6).
Based on visual interpretation of modeled results, CC and total BA distribution patterns at SNWR appear to be similar (Figure 7), while distributions of bur oak and red oaks RBA appear to differ substantially (Figure 7). Distribution of bur oak with higher RBA levels is concentrated in the north-central and northwestern regions of SNWR, while red oaks at similar RBA levels are concentrated in the west-central and southeastern areas. The distribution of high red oaks RBA appears to correspond with high percent CC levels, while bur RBA distribution is less distinct with reference to percent CC. The RBA of dead oak is generally low across SNWR, with concentrations of higher RBA in the southeastern areas where forest total BA is generally low (Figure 7).
Landsat and oak savanna structure mapping
Landsat sensor data have been used in numerous studies to quantify forest structural parameters (Franklin 1986; Horler and Ahern 1986; Woodcock et al. 1994; Franco-Lopez et al. 2001; Xu et al. 2003; McRoberts et al. 2007; McRoberts 2008; Wolter et al. 2009), but use of iterative PLS regression (xPLS) as a data-fusion and model-calibration optimization technique for Landsat data analysis is relatively new (Wolter et al. 2009, 2012; Wolter and Townsend 2011; Singh et al. 2013). Wolter et al. (2012) used xPLS regression with multispectral, multitemporal winter Landsat data to model one element of oak savanna structure (total BA) with strong results (R2 = 0.898, RMSE = 2.73 m2 ha−1). However, our use of xPLS regression with multiseason Landsat sensor data for more comprehensive oak savanna structure mapping, including oak mortality, is novel. Moreover, simultaneous calibration of the six forest-structure models (CC, total BA, living BA, bur oak RBA, red oaks RBA, and dead oak RBA) using 23 of the original 94 image variables (75.5% model reduction) is notable because the process is fully automatic and repeatable. However, success of these analyses also hinges on the dates of Landsat imagery and indices used as inputs to the xPLS algorithm, which should always be coupled with knowledge of the biology and phenology of forest species under investigation (Boyer et al. 1988; Miller et al. 1991; Wolter et al. 1995; Wolter et al. 2012).
Automated selection of image variables by xPLS regression during model calibration highlights the strengths of several image variables for quantifying oak savanna structure at SNWR—specifically, the importance of SVR (Wolter et al. 2008) and Tasseled Cap wetness (Crist and Cicone 1984) for nearly all elements of structure modeled (Figure 5). Landsat wetness is a linear contrast of SWIR with the near-infrared and the visible portion of the electromagnetic spectrum and is known to provide highly accurate information on forest canopy structure due, in part, to relative insensitivity to topographical effects (Cohen and Spies 1992; Cohen et al. 1995; Hansen et al. 2001). Similarly, SVR is the ratio of SWIR to visible, while MSI (Rock et al. 1986) is a ratio of SWIR to near-infrared (Table 2); both are known to be strong predictors of forest basal area and tree species relative to basal area in subboreal forests (Wolter et al. 2008; Wolter and Townsend 2011). The SWIR component of these indices is known to vary inversely with forest biomass due to increases in shadowing with increasing heterogeneity in forest structure (Berner et al. 2012). Percent CC increases linearly with total BA (Jennings et al. 1999); therefore, it is not surprising that SWIR-based variables were of similar magnitude and sign for estimating these structural parameters (Figure 5).
Predictor variable loading (magnitude and sign) for bur oak RBA and red oaks RBA may indicate phenology differences between stands of the two oak subgenres at SNWR. For instance, of all the image variables retained during xPLS model calibration, each of the seven variables from 14 April 2009 (most for any one date) show opposite loading signs between the two oak subgenres (Figure 5). For living versus dead oak, the sign of respective loadings on almost all the retained image predictors were opposed, which is logical. The opposed April image (P) predictor loading signs for bur and the red oaks, however, is puzzling because both subgenres appeared, upon visual inspection of April imagery, as if they had not yet flushed leaves. Alternatively, it is possible that ground-layer green-up had begun prior to image acquisition. Hence, the Landsat sensor may be detecting potential differences in understory vegetation structure, composition, and phenology associated with the two oak subgenres (Leach and Givnish 1999; Schulte et al. 2011), which would be optimal in the absence of fully flushed oak leaves. However, we have neither ground-layer vegetation data nor specific phenology characteristics required to test this hypothesis.
Oak savanna management
Similar to oak savanna trends reported elsewhere in Minnesota (Margoles 2009), the apparent increase in the proportion of red oaks over bur oak at SNWR has resulted in a disproportionate spatial extent of canopy cover classes >61% (Figure 6), which agrees with tree composition (n = 55 point samples) reported by Kipfmuller and Hepola (2007). Similar to other oak savanna ecosystems (Leach and Givnish 1999), the distributions of bur oak at SNWR appears, based on visual inspection, to correspond with soils having higher silt content, such as Becker loam and Isanti loamy fine sand. Soils at SNWR are predominantly sandy, but the distribution of red oaks does not appear to be distinctly associated with any one soil type or series (see USFWS 2005). Patterns of red oak abundance at SNWR (especially northern pin) are likely aligned more closely to fire-disturbance history (White 1983; Crow 1988) and, potentially, ground-water availability (P. Drobney, U.S. Fish and Wildlife Service, Neal Smith National Wildlife Refuge, pers. comm.) than soil type alone. Whether bur oak regeneration has declined in favor of fire-intolerant northern red and northern pin oaks at SNWR over past decades due to changes in fire disturbance, competition for light from invasive understory species, herbivory from white-tailed deer Odocoileus virginianus, climate trends, or a combination of these factors, as observed elsewhere (Kipfmuller and Hepola 2007; Ziegler et al. 2008; Schulte et al. 2011), is unclear. What is clear is that red oaks are increasing in proportion and generally outcompeting bur oaks at SNWR.
Regardless of causal factors, there are no precise analogs of oak savanna remaining in this region from which to model or guide restorative efforts (Bowles and McBride 1998); although, a few highly diverse savannas do exist (e.g., Timberhill in southern Iowa; Wilhelm and Rericha 2007) that provide good benchmark examples of species–light relationships and ecosystem functioning. In the absence of extant savanna benchmarks, restoration targets are often identified by old aerial photos, General Land Office maps, surveyor's comments from original land surveys, herbaria, dendrochronology, evidence of plant communities in archeological collections and soil sedimentation, and other historical records (USFWS 2012).
In order to achieve regional restoration objectives based on historical references, appropriate management tools should be implemented broadly and, ideally, across ownership boundaries to mimic historical disturbances that influenced and maintained savanna ecosystems. For example, because many savannas are highly degraded, with fire-intolerant species dominating the overstory and understory, prescribed burning is used as a tool to mimic historical fires. Prescribed burning accomplishes several things: it promotes the growth of native wildflowers and native grasses; retards competition from exotic grasses and fire-intolerant trees and shrubs; promotes higher water-table levels; and converts heavily wooded or closed-canopy forest ecosystems to a more open-canopy condition that more closely approximates oak savannas in this region (USFWS 2005). The severity of degradation may warrant additional restoration techniques such as tree thinning operations. At our SNWR case-study ecosystem, tree-thinning operations are often implemented prior to prescribed burning to expedite restoration, focusing efforts on closed-canopy oak savanna forest with high BA values and favoring bur oak trees for retention. Removal of undesirable trees from the overstory creates open canopy conditions facilitating light penetration and changes species composition. Depending on specific management objectives or definitions (Bray 1955; Curtis 1959; Botts et al. 1994), these Landsat-based structure maps can be developed relatively quickly, and repeatedly (annual to decadal intervals), for large regions and serve as a framework from which to help target or monitor regional- or local-scale restorative treatments.
In our local-scale test ecosystem for example, one of the management objectives in the Comprehensive Conservation Plan for SNWR is conversion of existing woodlands to oak savanna (USFWS 2005). Basal area (BA) and CC attributes among oak savannas in this region typically ranged from 1 to 11.5 m2 ha−1 (5–50 ft2 acre−1) and 10–70%, respectively. If an upper BA limit of 11.5 m2 ha−1 is used as the threshold above which restorative treatments are applied, then 66.5% (i.e., 3,698 ha) of the oak savanna at SNWR qualifies for treatment. Rather, if a 50% upper CC threshold is used, then 57.7% (3,205 ha) qualifies for treatment. Using the intersection of these two ground areas (11.5 m2 ha−1 BA and 50% CC) qualifies 55.9% (3,106 ha) of the oak savanna at our SNWR test study for restorative treatments. Furthermore, these Landsat-based structure maps can provide spatial reference to assist managers in focusing restoration efforts, such as tree-thinning operations, to areas with high BA and CC values. Managers at SNWR may also wish to apply their resources to those locations where red oaks make up >50% of the total forest BA. In this case, 28% (1,556 ha) of forest at SNWR qualifies for savanna restoration treatment (Figure 7). Use of these Landsat-derived forest-structure data and other ancillary GIS data can effectively help managers prioritize oak savanna ecosystems for restorative treatments, at both the local scale (such as our SNWR study site), but also for regional-scope management. With this said, it is important to note that while overstory CC and BA are certainly salient factors affecting understory light environment, specific determination of whether a site is truly “degraded,” is also a function of the quality and composition of understory vegetation present, regardless of overstory factors (P. Drobney, pers. comm.). Hence, these Landsat data should be used as a template to guide local- and regional-scale restoration efforts.
Within our SNWR test-study ecosystem, managers (E. Berkley and G. Dehmer, U.S. Fish and Wildlife Service, pers. obs.) estimate that approximately 50% of the mortality mapped using 2003 and 2009 Landsat image data (Figure 7) was due to oak wilt disease caused by the fungus Ceratocystis fagacearum, which affects red oaks more readily than white oaks (Brinkman 1952; Boyce 1959; Jacobi and MacDonald 1980). Once oak wilt spreads to an area, such as SNWR, the amount of damage and rate of spread is a function of stem density, tree age, and site differences (Menges and Kuntz 1985). Although practical methods for retarding the spread of oak wilt outside of urban settings do not exist, the Landsat-based forest-structure maps produced in this study could be used to facilitate a better understanding of the variability in rates of disease spread (Juzwik 2000; Downing et al. 2009).
Tree-species composition, distribution, mortality, basal area (BA), and canopy cover (CC) are important forest-structure parameters that are intrinsically linked to oak savanna restoration ecology in the upper Midwest region of the United States. We have demonstrated that predictive forest-structure models, calibrated using xPLS regression, simple ground measurements, and readily accessible Landsat image data are capable of producing moderate to highly accurate, spatially explicit, estimates for these parameters. Although we chose a relatively small oak savanna ecosystem for this proof of concept, our intent is for repeated, regional application. With regard to these Landsat-derived oak structure estimates (in the context of savanna restoration), one can only speculate whether a site is truly “degraded,” because much also depends on specific understory vegetation composition. However, these Landsat-based structure estimates can be used both locally and regionally in a variety of ways to 1) visualize and quantify forest composition trends, 2) prioritize potential areas for oak savanna restoration efforts (e.g., identification of overgrown oak savanna for thinning), and 3) monitor BA and CC of oak savanna through space and time as restoration treatments are implemented. Forest mortality models predicting the relative BA of living versus dead oaks may also be useful for determining rates of mortality (whether from prescribed fire or disease). This Landsat-based method of forest-structure estimation is robust and may be extended as a large-scale assessment tool to guide and improve management, planning, and implementation of oak savanna restoration efforts at multiple scales.
Please note: The Journal of Fish and Wildlife Management is not responsible for the content or functionality of any supplemental material. Queries should be directed to the corresponding author for the article.
Table S1. Sherburne National Wildlife Refuge 2009 field-plot data—average values calculated from the five subplots per plot. Also, the associated Universal Transverse Mercator zone 15 (WGS84) plot center coordinates (meters) used to extract Landsat band and index values.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S1 (14 KB XLSX).
Figure S1. Subset of the 24 January 2009 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S2 (13.8 MB TIF).
Figure S2. Subset of the 14 April 2009 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S3 (13.8 MB TIF).
Figure S3. Subset of the 1 June 2009 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S4 (13.8 MB TIF).
Figure S4. Subset of the 17 August 2008 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S5 (13.8 MB TIF).
Figure S5. Subset of the 4 October 2008 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S6 (13.8 MB TIF).
Figure S6. Subset of the 7 October 2003 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S7 (13.8 MB TIF).
Figure S7. Subset of the 20 October 2008 Landsat-5 Thematic Mapper imagery that was used to calibrate biophysical models of oak forest structure for the Sherburne National Wildlife Refuge test site in central Minnesota. These Landsat data are provided in Universal Transverse Mercator zone 15 projection and WGS84 datum at 30-m spatial resolution.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S8 (13.8 MB TIF).
Figure S8. Modeled percent canopy cover estimates for the Sherburne National Wildlife Refuge test site in central Minnesota. Percent canopy cover estimates were generated at 30-m spatial resolution in Universal Transverse Mercator zone 15, WGS84 coordinates.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S9 (4097 KB TIF).
Figure S9. Modeled estimates of total oak basal area for the Sherburne National Wildlife Refuge test site in central Minnesota. Basal area estimates were generated at 30-m spatial resolution in Universal Transverse Mercator zone 15, WGS84 coordinates and represent the trunk cross-section area (square meters measured at 1.37 m above ground) of all forest trees per unit ground area (m2 ha–1).
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S10 (4097 KB TIF).
Figure S10. Modeled estimates of living oak basal area for the Sherburne National Wildlife Refuge test site in central Minnesota. Basal area estimates were generated at 30-m spatial resolution in Universal Transverse Mercator zone 15, WGS84 coordinates and represent the trunk cross-section area (square meters measured at 1.37 m above ground) of living oak trees (Q. rubra, Q. ellipsoidalis, and Q. macrocarpa) per unit ground area (m2 ha−1).
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S11 (4097 KB TIF).
Figure S11. Modeled estimates of relative basal area of bur oak for the Sherburne National Wildlife Refuge test site in central Minnesota. Relative basal area estimates were generated at 30-m spatial resolution in Universal Transverse Mercator zone 15, WGS84 coordinates and represent the ratio of trunk cross-section area (square meters measured at 1.37 m above ground) of bur oak Q. macrocarpa per unit area to the total cross-section area of all tree trunks (live + dead) in the same units (m2 ha−1).
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S12 (4097 KB TIF).
Figure S12. Modeled estimates of relative basal area of red oaks for the Sherburne National Wildlife Refuge test site in central Minnesota. Relative basal area estimates were generated at 30-m spatial resolution in Universal Transverse Mercator zone 15, WGS84 coordinates and represent the ratio of cross-section area (square meters measured at 1.37 m above ground) of red oaks Q. rubra and Q. ellipsoidalis per unit area to the total cross-section area of all tree trunks (live + dead) in the same units (m2 ha−1).
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S13 (4097 KB TIF).
Figure S13. Modeled estimates of relative basal area of dead oaks for the Sherburne National Wildlife Refuge test site in central Minnesota. Relative basal area estimates were generated at 30-m spatial resolution in Universal Transverse Mercator zone 15, WGS84 coordinates and represent the ratio of cross-section area (square meters; measured at 1.37 m above ground) of dead oak trees per unit area to the total cross-section area of all stems (live + dead) in the same units (m2 ha−1).
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S14 (4097 KB TIF).
Reference S1. Buchanan H. 1996. Oak savanna restoration at Myre–Big Island State Park. Student On-Line Journal 1(2):1–4.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S15; also available at http://conservancy.umn.edu/bitstream/58421/1/1.2.Buchanan.pdf (269 KB PDF).
Reference S2. Law JR, Johnson PS, Houf G. 1994. A crown cover chart for oak savannas. Columbia, MO: USDA Forest Service, North Central Forest Experiment Station. Technical brief TB-NC-2.
Reference S3. Kipfmuller KF, Hepola T. 2007. Fire history and age structure analysis in the Sherburne National Wildlife Refuge: establishing reference conditions in a remnant oak savanna woodland. Fort Collins, CO: USDA Forest Service, Rocky Mountain Research Station. USDA Forest Service Proceedings RMRS-P-46CD, 507–514.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S17; also available at http://www.fs.fed.us/rm/pubs/rmrs_p046/rmrs_p046_507_514.pdf (890 KB PDF).
Reference S4. [USFWS] U.S. Fish and Wildlife Service. 2005. Comprehensive conservation plan: Sherburne National Wildlife Refuge. Washington, D.C.: U.S. Department of Interior.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S18; also available at http://www.fws.gov/Midwest/Planning/sherburne/ccp/FinalCCP.pdf (2563 KB PDF).
Reference S5. Margoles SS. 2009. Perspectives on oak savanna restoration in Minnesota: a dendroecological approach. Doctoral dissertation. Minneapolis: University of Minnesota.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S19; also available at http://conservancy.umn.edu/bitstream/60108/1/Margoles_Sarah_December%202009.pdf (3641 KB PDF).
Reference S6. Wilhelm G, Rericha L. 2007. Timberhill savanna assessment of landscape management. Elmhurst, IL: Conservation Research Institute.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S20; also available at http://wwx.inhs.illinois.edu/files/1913/4021/3314/Timberhill_Final.pdf (5553 KB PDF).
Reference S7. [USFW] U.S. Fish and Wildlife Service. 2012. Working definition of “oak savanna” for restoration efforts at Sherburne NWR. Unpublished addendum to USFWS (2005) Comprehensive Conservation Plan.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S21 (144 KB PDF).
Reference S8. Downing MC, Thomas V, Juzwik J, Appel DN, Reich RM, Camilli K. 2009. Using classification tree analysis to predict oak wilt distribution in Minnesota and Texas. Pages 67–83 in Billings RF, Appel DN, editors. Proceedings of the national oak wilt symposium. Newtown Square, PA: USDA Forest Service, Northern Research Station. Texas Forest Service Publication 166:61–69.
Found at DOI: http://dx.doi.org/10.3996/022013-JFWM-010.S22; also available at http://texasoakwilt.org/Professionals/NOWS/conference_assets/NOWS_Proceedings.pdf#page=6 (2785 KB PDF).
This research was supported by funding from the U.S. Fish and Wildlife Service. The authors would also like to acknowledge the three anonymous Journal reviewers as well as the Subject Editor whose comments and suggestions greatly improve the quality of this manuscript.
Any use of trade, product, or firm names is for descriptive purposes only and does not imply endorsement by the U.S. Government.
Wolter PT, Berkley EA, Peckham SD, Singh A. 2014. Satellite-based management tool for oak savanna ecosystem restoration. Journal of Fish and Wildlife Management 5(2):252–269; e1944-687X. doi: 10.3996/022013-JFWM-010
The findings and conclusions in this article are those of the author(s) and do not necessarily represent the views of the U.S. Fish and Wildlife Service.