Low-level Adelges tsugae infestation detection in New England through partition modeling of Landsat data

Remote Sensing of Environment 190 (2017) 13–25

Contents lists available at ScienceDirect

Remote Sensing of Environment journal homepage: www.elsevier.com/locate/rse

Low-level Adelges tsugae infestation detection in New England through partition modeling of Landsat data Justin P. Williams a,⁎, Ryan P. Hanavan b, Barrett N. Rock c, Subhash C. Minocha d, Ernst Linder e a

Department of Natural Resources and the Environment, James Hall, University of New Hampshire, Durham, NH 03824, USA U.S. Forest Service, Forest Health Protection, 271 Mast Rd., Durham, NH 03824, USA Earth Systems Research Center, Morse Hall, University of New Hampshire, Durham, NH 03824, USA d Department of Biological Sciences, Rudman Hall, University of New Hampshire, Durham, NH 03824, USA e Department of Mathematics and Statistics, Kingsbury Hall, University of New Hampshire, Durham, NH 03284, USA b c

a r t i c l e

i n f o

Article history: Received 9 August 2016 Received in revised form 23 November 2016 Accepted 12 December 2016 Available online xxxx Keywords: Landsat Adelges tsugae Invasive species NDVI Decision tree classification

a b s t r a c t The hemlock woolly adelgid (HWA) (Hemiptera: Adelges tsugae Annand) is an invasive insect causing damage to Eastern (Tsuga canadensis (L.) Carr.) and Carolina (Tsuga caroliniana Engelm.) hemlock trees in the eastern United States. Maine and New Hampshire are currently the northernmost front of HWA's range. Developing methods to locate newly infested stands is paramount in the effort to monitor HWA range expansion; presently, the most reliable method of detection requires extensive on-the-ground manual surveying. Field surveys for invasive pests like HWA consume multiple resources, limiting the amount of area that can be surveyed, and the results often misrepresent the true extent or patchiness of an invasion. Satellite based remote sensing, using vegetation indices, enables us to detect insect driven changes in forest health at a landscape scale. Our objectives were to classify HWA infested hemlock stands along the infestation front and demonstrate how the classification product could be used to improve HWA survey planning. Our workflow consisted of 1) modeling hemlock habitat suitability using a Maximum Entropy algorithm; 2) developing decision tree rules to classify likely infested stands from a Landsat time series, using the habitat suitability model to mask out non-hemlock areas; and 3) field check the final classification product. The hemlock habitat suitability model attained an overall accuracy of 68.2%. Partitioning of leaf-on multi-year (1995–2013) Landsat 5 and 8 data resulted in seven probability of infestation classes with a training R2 = 0.782. Agreement between the classification and previously reported HWA infestations was 75.0% in conifer forests, 33.3% in mixed forests and 50.0% in deciduous forests. Agreement between the classification and test-survey locations was 78.6%; verified new infestations were detected as far as 19 km away from previously reported infestations. The methods presented outline how Landsat could be used to detect low-level HWA infestations. Classification products such as this ultimately could be used by federal and state agencies to target specific areas for efficient survey, suppression, and eradication efforts. Published by Elsevier Inc.

1. Introduction Adelges tsugae Annand (Hemiptera: Adelgidae), the hemlock woolly adelgid (HWA), is an invasive and damaging insect in the eastern United States. Both eastern (Tsuga canadensis (L.) Carr.) and Carolina (Tsuga caroliniana Eng.) hemlock are at risk and approximately half the combined natural range of these two species is currently infested. The HWA's biannual and asexual reproductive biology, in addition to its 1–2 mm size and windborne dispersal vector (Evans and Gregoire, 2007), has resulted in rapid range expansion since its first detection in ⁎ Corresponding author at: U.S. Forest Service Southern Research Station, 775 Stone Blvd., Box 9861, Mississippi State, Mississippi 39762, USA. E-mail addresses: [email protected] (J.P. Williams), [email protected] (R.P. Hanavan), [email protected] (B.N. Rock), [email protected] (S.C. Minocha), [email protected] (E. Linder).

http://dx.doi.org/10.1016/j.rse.2016.12.005 0034-4257/Published by Elsevier Inc.

Richmond, VA in 1951 (Souto et al., 1995). Currently, Maine and New Hampshire are the northernmost front of the HWA's range. Developing methods to remotely detect infested stands prior to significant damage is paramount in effectively monitoring and managing the spread of this pest. Hemlock woolly adelgid initiates host tree and stand level decline by feeding on starch reserves and nutrient solution in the xylem ray parenchyma (Young et al., 1995). Damage consists of premature needle loss, branch dieback and subsequent mortality; this damage can occur over a 15 + year period (McClure, 1991), the timing and severity of which has been shown to be negatively correlated with latitude (Orwig et al., 2012) on account of winter temperatures in the north which reduce HWA populations (Trotter and Shields, 2009). Infestations and resultant damage are often patchy throughout the landscape (Mayer et al., 2002; Orwig et al., 2002). Along the infestation front, ground based surveys of high risk areas are the most reliable early detection effort.

14

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

Since the first detection of HWA in Portsmouth, NH, in the year 2000, efforts to monitor HWA have consisted of cooperative state and federal field surveys. The states of New Hampshire and Maine publicize the resulting infestation maps at the town level; federal HWA infestation maps are produced at the state county level. In New Hampshire and Maine, state run surveys using the methods of Costa and Onken (2006) consist of surveying two branches for each of 100 trees at five or more high risk hemlock stands within a town. High risk stands are located near water bodies, travel corridors, housing developments and infestations in neighboring towns (Struble et al., 2011). Although this is a reliable method that has detected multiple new infestations, it may also be inadequate because HWA may exist at locations other than the selected high risk survey sites, thereby resulting in false negative infestation statuses. In addition, maps produced through these surveys only indicate whether a town or county is infested and therefore they do not truly characterize the patchy distribution of infestations across the landscape. Remotely detecting areas that have a high probability of being infested could therefore decrease the time and cost of surveys by increasing the chances of locating previously unknown infested areas while providing a landscape-scale perspective of the patchiness and extent of HWA infestations. Several studies have investigated the use of Landsat Thematic Mapper (TM) reflectance data, transformations, and vegetation indices (VIs) to detect hemlock decline caused by HWA. Royle and Lathrop (2002) used the normalized difference vegetation index (NDVI) and change detection techniques to identify four classes of hemlock defoliation in northern New Jersey with 82% overall accuracy. Bonneau et al. (1999) used the modified soil-adjusted vegetation index (MSAVI2) and clustering techniques to identify four classes of hemlock health in Connecticut with 82% accuracy. In both studies the damage caused by HWA ranged from light defoliation (b 25%) in newly infested areas to widespread mortality. Jones et al. (2015) used anniversary leaf-off Landsat TM data, transformations, and image difference values, along with ancillary soil and topographic layers, in MaxEnt species distribution modeling software to predict five classes of HWA damage at the Delaware Water Gap National Recreation Area (AUC = 0.986). Our study aimed to improve early infestation symptom detection by focusing on the northern front of HWA spread where defoliation and mortality of hemlock has been minimal. In our previous work, Williams et al. (2016) reported significantly less reflectance for HWA infested current year foliage in the visible (400–680 nm) and red edge (680–750 nm) spectral regions, compared to non-infested hemlock, during June and July. This was thought to be the result of greater chlorophyll concentrations absorbing greater light energy. The data suggested that wavelengths or indices correlated with leaf chlorophyll content and turgidity would likely be useful in detecting early stage HWA infestations. More importantly, the data suggested that the timing of remote sensing data collection was critical because the greatest differences in leaf reflectance occurred in the weeks immediately following HWA sisten attachment at the base of new foliage (late June through July in the northeast). Hypotheses for the perceived increases in chlorophyll and moisture content included a compensatory response of the host to HWA (shifting resource allocation), or a manipulation of the host by HWA (Gómez et al., 2012). During the planning stages of this study change detection methods such as image differencing and discriminant function analysis were investigated as possible means of separating HWA infested from noninfested forests. These methods were abandoned when it became clear that the changes happening to hemlock forests in this area were occurring too slowly in time and space to be separated from other sources of change (see Vogelmann et al. (2016) for discussion on monitoring gradual change). Subsequent review of publications regarding decision tree analysis of inter-seasonal Landsat data for forest inventory purposes (Doucette et al., 2009; Koch et al., 2005; Lefsky et al., 2001; Schriever and Congalton, 1995) gave rise to the idea of using partition (decision tree) modeling of multi-year Landsat data to detect HWA

infested hemlock stands, using the ‘time of year’ and spectral guidelines outlined in Williams et al. (2016). Classification learning via partition modeling is useful because it can statistically identify non-parametric patterns in data to develop hard classification rules (Doucette et al., 2009). Although the algorithms used can vary by software or program proprietor, in general, decision trees predict class membership by recursively partitioning a data set into more homogenous subsets (Defries and Chan, 2000; Jensen, 2005). The hypotheses, rules and conditions derived from the analysis can then be entered as classification rules in an expert classification system. The goal of this study was to develop methods for using Landsat data to detect HWA infested hemlock stands on a landscape scale prior to significant defoliation or mortality. Our objectives were to 1. Classify HWA infested hemlock stands along the infestation front, and 2. Demonstrate how the final product could be used to coordinate ground surveys for likely infested stands. 2. Methods 2.1. Study region The region of interest for this study encompassed approximately 8800 km2 of southern New Hampshire and Maine (Fig. 1). Elevation ranged from 0 m at the coast to 889 m above mean sea level. Forest types of this region include the Northern Mesic Hardwood and Conifer type which consists of hardwoods such as oak (Quercus sp.), hickory (Carya sp.) maple (Acer sp.), birch (Betula sp.), and beech (Fagus grandifolia), and softwoods such as pine (Pinus sp.), and hemlock (Tsuga canadensis); and the higher elevation Red Spruce and Fir forest type consisting of spruce (Picea sp.) and fir (Abies balsamea) (“U.S. National Vegetation Classification”, 2014; Eyre, 1980). Within this region of interest four publicly owned properties were utilized so that classification accuracy could be assessed at both large (property level) and small (regional) scales (Fig. 1). The properties used in this study were the Massabesic Experimental Forest in Alfred, Maine (N43.43927, W70.67931) which has been infested with HWA since 2012; the Rachel Carson Wildlife Refuge on Cutts Island, Kittery, Maine (N43.09558, W70.67142) which has been infested with HWA since 2003; the RussellAbbott State Forest in Wilton, New Hampshire (N42.79048, W71.76155) which has been infested with HWA since 2011; and the Northwood Meadows State Park in Northwood, New Hampshire (N43.20466, W71.19827) where HWA presence had not yet been detected as of 2015. 2.2. Reference data collection and treatment Reference data points used in the partition modeling and in accuracy assessments were collected during the summer of 2013 (see Fig. 2). For the regional reference dataset, 150 random points were established (minimum distance = 1 km) within forested areas as classified by the 2006 National Land Cover Dataset (NLCD) (Fry et al., 2011) using ArcMap 10.1 (ESRI, 2012) create random points tool (Fig. 3). Field crews located these points using handheld GPS units; 40 of those plots were inaccessible. At each point (N = 110) a 10-factor basal area prism (BAF-10) was used to estimate the percent basal area occupied by hemlock. Presence of HWA was determined by visually checking the underside of the terminal 1 m of all branches within reach for each hemlock in the plot. If no branches were within reach, understory and plot perimeter hemlock were checked instead. In all, 13 of the regional reference data points were found to be infested with HWA (Fig. 3). Establishment of the property scale reference points at the Russell-Abbott State Forest (N = 29), Rachel Carson Wildlife Refuge (N = 29), Massabesic Experimental Forest (N = 29) and Northwood

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

15

Fig. 1. A map indicating the region of interest and the locations of the study plots.

Meadows State Park (N = 30) and the data collected therein followed the same protocols. A minimum distance between points was determined by the area of the property and ranged from 50 m at Rachel Carson Wildlife Reserve to 150 m at Massabesic Experimental Forest.

cultivated and wetland classes was used in accuracy and agreement analyses (see Table 1) (Anderson et al., 1976).

In addition to these reference data points, a subset of point data of known HWA infestations from state surveys and public reporting were used in classification validation procedures. A total of 331 confirmed points of infestation were within the study region. Five points in highly developed areas were excluded due to the lack of tree cover. Landsat's 30 m spatial resolution often results in pixels having mixed percentages of different land cover types (Fisher, 1997); without a high percentage of hemlock tree cover within a pixel the likelihood that the pixel would be identified as infested decreases. Twelve points that were within 30 m of a waterbody were excluded to account for a ± 1 pixel missregistration that may have occurred during image processing and masking procedures. Surveyed areas that contained a high density of HWA detection points were thinned so that there was at least 300 m distance between points, resulting in 131 points being excluded. If there were a choice between two or more points, points in conifer or mixed forests were preferred over residential or deciduous forest points (Table 1). The point thinning exercise resulted in N = 183 valid infestation reference points. Using ArcMap, the 2006 NLCD (Fry et al., 2011) Level II class value was extracted for each HWA detection point (Table 1). A modified classification system (Levels I and II) that aggregated the multiple developed (20s), cultivated (80s), and wetland (90s) class into single developed,

The distribution of suitable hemlock habitat within the study region was modeled using a Maximum Entropy algorithm in MaxEnt software version 3.3.3k (Phillips et al., 2006). Thirteen environmental variables (Table 2) and 648 hemlock presence points photo-interpreted (Hershey and Befort, 1995) from Google Earth leaf-off imagery (2006– 2013) were used in the model. Environmental variables were transformed into 30 m2 raster grid files in New Hampshire State Plane (NAD83) coordinate system. The final model was an average of three individual model iterations; iterations were randomly seeded bootstrapped models with a 50/50 training and validation split of the presence data and used 10,000 random background points. Linear, quadratic and product feature types were selected; the regularization multiplier was set at 0.25 and model output was set to logistic data format. Model statistical outputs included area under the receiving operator curve (AUC), environmental variable percent contribution and permutation importance, and environmental variable response curves. To create a binary map of suitable or unsuitable hemlock habitat a logistic threshold of 0.33 was calculated from the model test data using the maximum difference of cumulative frequency of probabilities method (Browning et al., 2005; Clark et al., 2012; Fei et al., 2007). Pixel values greater than the threshold value were deemed suitable hemlock habitat. The final model was assessed for accuracy at small and large scales using the regional reference data (N = 110), locations of known HWA

2.3. Hemlock habitat suitability

16

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

infestations (N = 183), and property level (N ≤ 30 points) data sets (Fig. 2). A 30 m buffer around each point was used to account for spatial errors associated with using handheld GPS units and the environmental variable datasets. Error matrices were constructed with overall, producer's, and user's accuracies calculated (Congalton et al., 1983).

2.4. Remote sensing

Fig. 2. This flowchart summarizes where and how the multiple datasets were applied in the modeling and validation processes.

Landsat TM 5 and Landsat 8 Operational Land Imager (OLI) data were acquired though Earth Explorer and Glovis USGS websites (Table 3). Cloud free (≤10% cloud cover) July imagery from Landsat 5 TM was acquired from 1995 (pre-HWA invasion) to 2011. Previous work had indicated that this was an optimal time of year for spectral discrimination between infested and non-infested hemlock (Williams et al., 2016). In 2013 Landsat 8 OLI was successfully launched to orbit and was fully operational by spring; however, images from June and July were obscured by clouds. The best image available was September 23, 2013 which had minor cumulus cloud cover (b10%). Although September was not within the optimal time period for using satellite data, spectral differences between infested and non-infested foliage were still apparent (Williams et al., 2016). Furthermore, using Landsat data from 2013 was important because the reference data were collected in 2013, meaning that the reference data and Landsat data would be directly related. All image processing was performed in ERDAS Imagine 2013 (Hexagon Geospatial, Norcross, GA). For each image date six Landsat spectral bands within the visible, near-infrared (NIR), and middle-infrared (MidIR) reflectance regions in digital number format were layer stacked into a single image file. Image stacks were projected to New Hampshire State Plane (NAD83) coordinate system. ERDAS's Georeferencing

Fig. 3. General location and distribution of random reference data points within the study region. Reference points marked with a star indicate presence of HWA.

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25 Table 1 Description of the 2006 NLCD Level II classes assigned to each HWA detection point. Definitions are from http://www.mrlc.gov/nlcd01_leg.php. Only classes found with the region of interest are listed. Value Class

Description

11

Open water

21

Developed Open space

Areas of open water, generally with b25% cover of vegetation or soil. Areas with a mixture of some constructed materials, but mostly vegetation in the form of lawn grasses. Impervious surfaces account for b20% of total cover. These areas most commonly include large-lot single-family housing units, parks, golf courses, and vegetation planted in developed settings for recreation, erosion control, or aesthetic purposes. Areas with a mixture of constructed materials and vegetation. Impervious surfaces account for 20% to 49% percent of total cover. These areas most commonly include single-family housing units. Areas with a mixture of constructed materials and vegetation. Impervious surfaces account for 50% to 79% of the total cover. These areas most commonly include single-family housing units. Highly developed areas where people reside or work in high numbers. Examples include apartment complexes, row houses and commercial/industrial. Impervious surfaces account for 80% to 100% of the total cover. Areas of bedrock, desert pavement, scarps, talus, slides, volcanic material, glacial debris, sand dunes, strip mines, gravel pits and other accumulations of earthen material. Generally, vegetation accounts for b15% of total cover. Areas dominated by trees generally N5 m tall, and N 20% of total vegetation cover. N75% of the tree species shed foliage simultaneously in response to seasonal change. Areas dominated by trees generally N5 m tall, and N 20% of total vegetation cover. N75% of the tree species maintain their leaves all year. Canopy is never without green foliage. Areas dominated by trees generally N5 m tall, and N 20% of total vegetation cover. Neither deciduous nor evergreen species are N75% of total tree cover. Areas dominated by shrubs; b5 m tall with shrub canopy typically N20% of total vegetation. This class includes true shrubs, young trees in an early successional stage or trees stunted from environmental conditions. Areas dominated by gramanoid or herbaceous vegetation, generally N80% of total vegetation. These areas are not subject to intensive management such as tilling, but can be utilized for grazing. Areas of grasses, legumes, or grass-legume mixtures planted for livestock grazing or the production of seed or hay crops, typically on a perennial cycle. Pasture/hay vegetation accounts for N20% of total vegetation. Areas used for the production of annual crops, such as corn, soybeans, vegetables, tobacco, and cotton, and also perennial woody crops such as orchards and vineyards. Crop vegetation accounts for N20% of total vegetation. This class also includes all land being actively tilled. Areas where forest or shrubland vegetation accounts for N20% of vegetative cover and the soil or substrate is periodically saturated with or covered with water. Areas where perennial herbaceous vegetation accounts for N80% of vegetative cover and the soil or substrate is periodically saturated with or covered with water.

22

Developed Low intensity

23

Developed Medium intensity

24

Developed High intensity

31

Barren land (rock/sand/clay)

41

Deciduous forest

42

Evergreen forest

43

Mixed forest

52

Shrub/scrub

71

Grassland, herbaceous

81

Pasture, hay

82

Cultivated crops

90

Woody wetlands

95

Emergent herbaceous wetlands

wizard was used to co-register all image years to the 1995 reference image using affine or polynomial model settings (Table 3). Coregistered images were converted from digital number format to surface reflectance using the COST method (Chavez, 1996), then stretched to 8 bit radiometric resolution (0–255). Although Landsat 5 had a native 8 bit radiometric resolution, Landsat 8 products are delivered at 16 bit

17

Table 2 Source data for the predictor variables used in the MaxEnt modeling platform. Predictor variable

Source

Elevation Aspect Slope Average annual precipitation Average minimum temperature Average maximum temperature Soil drainage Soil taxonomic class Parent material Soil pH Soil water Topographic wetness index Land cover type

LANDFIRE LANDFIRE LANDFIRE USDA NRCS Geospatial Data Gateway USDA NRCS Geospatial Data Gateway USDA NRCS Geospatial Data Gateway USDA NRCS Web Soil Survey USDA NRCS Web Soil Survey USDA NRCS Web Soil Survey USDA NRCS Web Soil Survey USDA NRCS Web Soil Survey USDA FS Remote Sensing and Application Center USDA NRCS Geospatial Data Gateway

radiometric resolution. The Landsat 8 data was rescaled to match the Landsat 5 data to limit possible modeling confusion and error. Principal components, VIs (NDVI, IR/R, MSI, MSAVI2), and Tasseled Cap (except 2013) were calculated for each year of imagery and layer stacked with the original six spectral bands (Table 4). Reflectance values for all spectral bands, VIs, and transformations were extracted for each of the 110 points of the regional reference data set for partition analysis in JMP (JMP Pro Version 10.0. SAS Institute Inc., Cary, NC, 1989–2007) (Fig. 2). Prior to partition analysis the regional reference data set was divided into training and validation sets. Of the 110 points included in the reference data set fifteen were excluded due to cloud interference in the imagery (including one HWA infested point), leaving 95 points remaining in the data set. Twenty percent (n = 19) of the remaining points were randomly selected to be reserved for validation procedures; however, three of these points were confirmed infested sites. Due to the limited number of confirmed infested sites in the reference data set (n = 12) it was important that all of them be used in the training portion of the analysis, therefore the three infested points that were chosen for validation were put back into the training set. This reduced the validation set to 16 non-infested points. The classification hypotheses and rules developed through the partition analyses were re-constructed in ERDAS Imagine Knowledge Engineer. The hemlock habitat suitability map was then used to mask out areas that potentially do not support hemlock. Finally, a supervised classification of the September 2013 Landsat image was used to mask out cloud, cloud shadow and highly developed areas (Fig. 2). 2.5. Accuracy and survey assessments Regional detection capabilities of the classification model were assessed using location data of 126 known HWA infestations that agreed with the hemlock habitat model (Fig. 2). The agreement assessment was separated by the modified 2006 NLCD (Fry et al., 2011) Level I and II land cover classes (see Table 1) to provide a more comprehensive examination of how well infestations could be detected under certain

Table 3 Landsat image acquisition dates and results of geometric corrections. Landsat sensor and image date

Geometric correction model

RMSE

Landsat 5 Landsat 5 Landsat 5 Landsat 5 Landsat 8

Reference image Affine 3rd order polynomial Affine Affine

0.000 0.271 0.498 0.489 0.579

07/04/1995 07/15/1999 07/10/2009 07/16/2011 09/23/2013

18

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

Table 4 Descriptions and formulae for vegetation health indices and Landsat data transformations. Transformation or vegetation index

Description

Tasseled cap

(Crist and Cicone, 1984) Brightness, greenness, and wetness bands calculated using linear combinations of Landsat TM bands multiplied by coefficients. Transformation of original highly correlated bands into orthogonal (uncorrelated) bands; the first 3 components account for N98% of the variance from the original data set. (Rouse et al., 1974) Normalized difference vegetation index: sensitive to the amount of biomass in a pixel, signal saturates in high-biomass conditions. (Band 4 − Band 3)/(Band 4 + Band 3) (Birth and McVey, 1968) Sensitive to the amount of biomass in a pixel, but signal does not saturate in high-biomass conditions. (Band 4/Band 3) (Rock et al., 1986) Moisture stress index: scaled from 0 to 1 with higher values indicating varying degrees of vegetative moisture stress. (Band 5/Band 4) (Qi et al., 1994) Modified soil adjusted vegetation index 2: pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi 2

Principal components NDVI

IR/R MSI MSAVI2

2NIRþ1−

ð2NIRþ1Þ −8ðNIR−REDÞ 2

Biomass index adjusted for background noise from soils.

land cover conditions. For this dataset 90 m buffers around the points were used to account for spatial errors associated with handheld GPS units and Landsat data preparation and analysis. Property level accuracy assessments (N ≤ 30 points at each of four properties) were conducted with no regard to the hemlock habitat model or land cover type, and 30 m buffers were used to account for spatial error (Fig. 2). Traditional error matrices were constructed with overall, producers, and user's accuracies calculated. The predictive capability of the classification model and its usefulness to assist surveyors was assessed by sending field crews to public lands at fourteen locations scattered throughout the region (Fig. 2). Each of the 14 locations had undocumented HWA infestation statuses, were classified as infested, and were determined to have hemlock through photo interpretation. Five of the fourteen sites were located in towns already known to be infested; the remaining nine locations were in towns that did not have a previously reported HWA infestation. Surveys were completed using modified methods from Costa and Onken (2006). At each survey location four 1 ha blocks were established from a center point. In each block 25 trees (at least 2 branches per tree) were sampled; tree selection was randomized by having the surveyor walk 20 paces in a semi-random cardinal direction. Once an infested tree was located the survey stopped. If no HWA was found after sampling 100 trees, one could state with 75% reliability that b2% of the stand was infested. Note that not finding HWA during the survey did not imply that the area was adelgid free but that the percentage of infested trees was less than the detection threshold of 2%.

Accuracy of the habitat suitability model was evaluated at regional and property scales. Overall accuracy from the regional 110 point reference data set was 68.2% (Table 6). A second regional evaluation examined overall habitat model agreement with the list of known HWA infestations within the study region; 78.3% of known HWA infested points within areas classified as forest by the NLCD were modeled as suitable hemlock habitat (Table 7). Agreement between the habitat model and known HWA infestation points was substantially reduced in developed and other areas (Table 7). Overall accuracy of the hemlock habitat model at the property scale ranged from 33.3 to 80% (Tables 8-11).

3.2. Remote sensing Partition analysis of the Landsat data set resulted in seven probability of infestation classes (Fig. 5). Entropy R2 of the training data was 0.78 (RMSE = 0.165) with a misclassification rate of 0.038. Validation of the model was one sided since all 16 validation points were non-infested; all 16 points were correctly classified. The model relied on vegetation indices (NDVI and MSI) and individual band values (Bands 1, 2, 3, 7) from four of the 5 years (1995, 2009, 2011, and 2013) of imagery (Fig. 5). After masking procedures only pixels classified into the highest probability of infestation class (0.886) were included in the final map (Fig. 6A; Supp. Figs. 1 and 2). Of the area modeled as suitable hemlock habitat in the study region approximately 5% was classified with a high probability of being infested with HWA.

3. Results 3.1. Hemlock habitat suitability The final model of hemlock habitat suitability resulted from the average of three model iterations. Area underneath the receiver operating curves (AUC) for the training data was 0.831 (± 0.005) and 0.788 (±0.019) for the test data. Response curves produced by MaxEnt indicated how each predictor variable influenced the model and showed the best potential environments for hemlock success. The preferred hemlock habitat consisted of moderate (~20%) north facing slopes on well drained, acidic (~4.5 pH), loamy till inceptisol soils derived from granite, gneiss and schist parent materials. Land cover type, soil taxonomic class, soil drainage class and soil parent material were the most important variables in building the model (Table 5). Implementing a logistic threshold of 0.33 resulted in a binary hemlock non-hemlock map (Fig. 4). The model predicted that nearly one third of the total study region (300,000+ ha) was suitable hemlock habitat.

Table 5 Percent contribution and permutation importance of the predictor variables in development of the hemlock habitat suitability model. Variable

Percent contribution

Permutation importance

Land cover Soil taxonomic class Soil drainage class Soil parent material Slope Topographic wetness index Annual precipitation Elevation Soil pH Aspect Maximum annual temperature Minimum annual temperature Soil water availability

60.6 13.6 13.3 8.5 1.1 0.8 0.6 0.3 0.3 0.3 0.3 0.2 0

59.5 11.6 9.9 9.9 2.1 0.4 1.3 1.7 1.8 0.5 0.3 0.6 0.3

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

19

Fig. 4. Suitable hemlock habitat as predicted by the MaxEnt software through thirteen environmental variables (Table 5).

3.3. Accuracy and survey assessments Classification agreement was evaluated using point data of 126 known HWA infested areas that were within the modeled hemlock habitat (Fig. 6B). Overall 57.1% of current known infestations were detected through the classification; 61.1% of all known infestations in forested areas were detected (Table 12). In coniferous forests 75.0% of known infestations were detected; classification agreement in forested areas declined as the hemlock component of forests decreased (Table 12). In mixed forests the classification detected 33.3% of known infestations, and 50.0% of infestations in deciduous forests were detected (Table 12). Classification agreement with known infestations in developed and other areas were 44.8% and 57.1% respectively (Table 12).

Table 6 Regional accuracy assessment results of the hemlock habitat suitability model. Producer's accuracy (error of omission) indicates the likelihood of a pixel being correctly classified by the map producer; user's accuracy (error of commission) indicates the likelihood that a pixel classified on the map correctly represents that category on the ground (Congalton, 1991).

Reference points

Model

Hemlock Non-hemlock Producers

Hemlock

Non-hemlock

Users

41

19

68.3%

16

34

68.0%

71.9%

64.2%

75/110

Overall

68.2%

Infestation detection capabilities were also evaluated at each of the four properties within the study region. Plot surveys at Northwood Meadows State Park (Table 15) and at Massabesic Experimental Forest (Table 14) found no HWA infestations and overall accuracy was 86.6% and 72.4% respectively. The Rachel Carson Wildlife Refuge (Table 13) and Russell-Abbott State Forest (Table 16) were both highly infested properties; overall accuracy was 70.0% and 62.1% respectively. Predictive strength of the classification was tested by surveying fourteen locations where HWA had not been previously identified and reported. All five locations in towns with previously reported infestations were found to be infested (Fig. 7). Of the nine locations in towns not previously reported to have HWA infestations, 6 were found to be infested (Fig. 7). Overall 78.6% of the sites surveyed that

Table 7 Agreement of the hemlock habitat suitability model with point data of known HWA infestations, categorized by 2006 NLCD Level II class value. Land cover

Agree

Total

Percent

Forested Conifer Mixed Deciduous Developed Other Pasture Shrub Wetland Overall

90 52 18 20 29 7 2 1 4 126

115 56 28 31 55 13 4 2 7 183

78.3% 92.9% 64.3% 64.5% 52.7% 53.8% 50.0% 50.0% 57.1% 68.9%

Data in bold are major land cover catagories and overall agreement values. Data not in bold are subcatagories of the major catagories.

20

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

Table 8 Localized accuracy assessment results of the hemlock habitat suitability model at the Rachel Carson Wildlife Refuge in Kittery, Maine.

Table 11 Localized accuracy assessment results of the hemlock habitat suitability model at RussellAbbott State Forest in Wilton, New Hampshire.

Reference

Reference points Hemlock Model

Non-hemlock

Hemlock

4

20

Non-hemlock

0

6

Producers

100.0%

Users 16.7 %

Non-hemlock 100.0% 23.1 % Overall

10/ 30

Reference points Hemlock

Non-hemlock

Users

Hemlock

18

2

90.0%

Non-hemlock

5

4

44.4%

Producers

78.3%

Producers

Hemlock

Non-hemlock

Users

14

5

73.7%

7

3

66.6%

30.0% 37.5%

17/29

Overall

58.6%

33.3 %

Table 9 Localized accuracy assessment results of the hemlock habitat suitability model at the Massabesic Experimental Forest in Alfred, Maine.

Model

Model

Hemlock

points

66.6%

22/29

Overall

75.9%

were identified as infested through the classification were found to be infested with HWA. 4. Discussion and conclusions The MaxEnt model relied most on land cover and soils data in its prediction of suitable hemlock habitat. Land cover type was likely the most important variable because the presence points located through photointerpretation were often located in areas of dense hemlock which were classified as evergreen forest by the NLCD; in all 57% of the hemlock presence points were located in areas classified as evergreen forest in the 2006 NLCD (28% mixed forest; 12% deciduous forest; 1% each urban, wetland, and shrub/scrub). This variable, although helpful in building the predictive model by restraining it to coniferous forest areas, also introduced substantial error into the model by identifying most coniferous and mixed forest habitats as suitable hemlock habitat. Soil taxonomic class, drainage class and parent material type were the next most important variables. Interpretation of the modeled preferred hemlock habitat (moderate (~ 20%) north facing slopes on well drained acidic (~ 4.5 pH) loamy till inceptisol soils derived from granite, gneiss and schist parent materials) was similar to other published studies

(Clark et al., 2012; Godman and Lancaster, 1990). It is important to note that the interpretation of preferred hemlock habitat was confined to the study region. In contrast to Dunckel et al. (2015), where hemlock occurrence and percent basal area were modeled for the entire state of Maine, minimum and maximum average temperature and elevation were not important in predicting suitable hemlock habitat; this was likely due to the size of our region of interest where hemlock occurred over the entire range of values for those variables. The overall distribution of hemlock habitat mapped through the MaxEnt model was similar to existing host species datasets except for some areas along the seacoast. The National Insect and Disease Risk Map (NIDRM), produced by the U.S. Forest Service's Forest Health and Technology Enterprise Team (FHTET), and the Landscape Fire and Resource Management Tools (LANDFIRE), produced through the U.S. Department of Agriculture and the U.S. Department of the Interior, provide freely available host species raster data at 250 m and 30 m spatial resolution, respectively. Besides the difference in spatial resolution between the NIDRM and MaxEnt hemlock products, the most apparent difference between the MaxEnt model and the other host layers was that the MaxEnt model indicated greater suitable hemlock habitat in southeast New Hampshire and along the southern Maine seacoast. Using the NIDRM or LANDFIRE host layers in the masking procedures

Table 10 Localized accuracy assessment results of the hemlock habitat suitability model at the Northwood Meadows State Park in Northwood, New Hampshire.

Reference points

Model

Hemlock Non-hemlock Producers

Hemlock

Non-hemlock

Users

23

5

82.1%

1

1

95.8%

50.0% 16.7%

24/30

Overall

80.0%

Fig. 5. Partitioning of multi-year Landsat data resulted in this decision tree. The numbers within the rectangles are the band or vegetation index values upon which the data was split. The numbers within the diamonds are the resulting probability of HWA infestation under those reflectance based conditions. Only pixels with the greatest probability of being infested (0.886) were included in the final model.

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

21

Fig. 6. A. Final classification map indicating areas that are likely infested with HWA in red. B. Landsat 8 OLI false color composite (Band 7, Band 6, Band 4; RGB) with the HWA classification and known points of infestation overlaid. These figures highlight the patchiness of infestations across the landscape, the accuracy of the classification, and how the classification may be used to coordinate survey, suppression and eradication efforts.

would have excluded these important areas, many of which are infested with HWA, and mischaracterized the true range of infestation in the final classification map. Building a species model specific to this study region allowed us to retain these important habitat

areas, but unfortunately at the expense of over predicting suitable hemlock habitat. Although the MaxEnt model achieved a test AUC of 0.788 (± 0.019), a traditional accuracy assessment of the hemlock

22

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

Table 12 Regional agreement between point data of known HWA infestations and model predicted HWA infestations, categorized by 2006 NLCD Level II class value. Land cover

Agree

Total

Percent

Forested Conifer Mixed Deciduous Developed Other Pasture Shrub Wetland Overall

55 39 6 10 13 4 2 1 1 72

90 52 18 20 29 7 2 1 4 126

61.1% 75.0% 33.3% 50.0% 44.8% 57.1% 100.0% 100.0% 25.0% 57.1%

Table 14 Property-level agreement between plot survey data and model predicted HWA infestations at the Massabesic Experimental Forest in Alfred, Maine.

Reference points HWA

No-HWA

0

8

No-HWA

0

21

Producers

0.0%

HWA

Users 0.0%

Model 100.0%

72.4%

21/29

Overall

72.4%

Data in bold are major land cover catagories and overall agreement values. Data not in bold are subcatagories of the major catagories.

habitat model was sub-optimal regionally (b 70%) and inconsistent at the property scales (33.3–80%). For example the Rachel Carson Wildlife Refuge had the lowest overall accuracy of 33.3%, mostly due to white pine (Pinus strobus L.) existing in areas that were modeled as suitable hemlock habitat. User's accuracy of the hemlock class at Rachel Carson was 17% (100% for non-hemlock), meaning that of the areas identified as suitable hemlock habitat only 17% of those areas actually contained hemlock. In comparison at the other three properties the user's accuracies for the hemlock class were 73.7–90%, and 30–50% for the non-hemlock class. These inconsistencies indicated that MaxEnt's over-reliance on the NLCD (60.6% contribution) likely fostered the over-prediction of suitable hemlock habitat in some areas. However, it is also important to note that MaxEnt models the theoretical distribution of species based on occurrence data, but that the realized distribution of species may be quite different due to factors such as forest management and competition (Phillips et al., 2004). At the regional scale, user's and producer's accuracies of hemlock and nonhemlock were similar (64% and 71% respectively) and the model was generally a good representation of the distribution of suitable hemlock habitats. The purpose of the logistic threshold (0.33) was to maximize accuracy and to limit errors of omission and commission. Improvements in both large and small scale accuracy could be attained by incorporating remote sensing data, rather than remote sensing products, into the model. Output of the Landsat partition model consisted of seven probability of infestation classes. Given that recent HWA infestation in this region had relatively little impact in terms of defoliation and mortality, it seemed logical that only pixels with the highest probability of being infested should be retained in the final classification. The highest probability class of 0.866 was defined by a high 2013 NDVI (≥0.901), a high 1995 blue band (≥2), and a low 2013 MSI (≤0.379) value. This agreed with our previous research (Williams et al., 2016) in which laboratory based reflectance measurements indicated that infested hemlock in

this study region had greater NDVI and lower MSI values than noninfested hemlock. Williams et al. (2016) noted that increased leaf chlorophyll and moisture contents in the new flush of foliage on infested hemlock trees may have resulted from either a compensatory resource re-allocation by the host, or a manipulation of the host tissue by HWA (Gómez et al., 2012). On a separate but related note, Mayer et al. (2002) reported an initial increase in new growth on lightly infested trees, and during periods of hemlock recovery after a crash in HWA populations, in New Jersey. If hemlocks in this region respond to HWA infestation in a similar way, this could also explain the higher NDVI and MSI values as being important indicators of developing infestations. The high blue band values were likely helpful to the model in distinguishing between coniferous and deciduous pixels; coniferous forest pixels typically exhibit greater blue band reflectance. As discussed in Vogelmann et al. (2016), anniversary date (before and after) change detection may not be sufficient when attempting to identify gradual changes in remotely sensed data. Vogelmann et al. (2016) defined gradual changes as those occurring over several years or more, involving gradual spectral changes over that time. Hemlock woolly adelgid infestations in this region may be categorized as a source of gradual change; cold weather limitations on populations (Trotter and Shields, 2009) have restrained widespread major impacts and subdued infestation driven changes in hemlock canopy reflectance. The major advantage to using multiple years of Landsat data that covered the temporal range of HWA establishment within the study region, as opposed to comparing pixel values between points on a single image or a change in pixel values between two dates, is that a single observation point became a repeated observation point (Goodwin et al., 2008) upon which HWA establishment could be characterized. The fact that two indices from Landsat 8 (September 2013) defined the highest probability class validated the importance of using imagery that was concurrent with the reference data collection, but should not diminish the importance of the data from the previous year's in helping the model to isolate the infestation signature. Additionally, Landsat 8′s advances in land imaging such as an

Table 13 Property-level agreement between plot survey data and model predicted HWA infestations at the Rachel Carson Wildlife Refuge in Kittery, Maine.

Table 15 Property-level agreement between plot survey data and model predicted HWA infestations at the Northwood Meadows State Park in Northwood, New Hampshire.

Reference points

HWA

Reference points

HWA

No-HWA

Users

2

2

50.0%

Model

HWA HWA

0

No-HWA 4

Users 0.0%

Model No-HWA Producers

7 22.2%

19

73.1%

No-HWA

0

90.5%

21/30

Producers

0.0%

Overall

70.0%

26 86.7% Overall

100.0% 26/ 30 86.7%

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25 Table 16 Property-level agreement between plot survey data and model predicted HWA infestations at the Russell-Abbott State Forest in Wilton, New Hampshire.

Reference points HWA

No-HWA

HWA

6

3

No-HWA

8

12

Users 66.7 %

Model Producers

42.9%

60.0%

80.0%

18/29

Overall

62.1%

improved signal-to-noise ratio, narrower bandwidths, and increased radiometric resolution (Loveland and Irons, 2016), cannot be ignored as another possible reason for those indices defining the highest probability class. Regionally, agreement between the classification and known HWA infestation point data was greatest in coniferous forest areas (75.0%) and substantially lower in mixed (33.3%) and deciduous (50.0%) forests. This trend was expected since the masks applied to the classification should have eliminated most nonconiferous forested pixels and because of the inherent reduction of the hemlock spectral signal. Known HWA infestation points in deciduous areas would only have been included in the assessment if they were located within 30 m of modeled hemlock habitat and were within 90 m of modeled HWA infestations. This is likely a contributing factor in why greater classification accuracy was attained in the deciduous forest category than in the mixed forest category. Overall classification accuracy between the properties ranged from 62.1 to 86.7% (72.8% average). Although the overall accuracy values

23

indicated moderate to high classification accuracy at large scales, the producer's and user's accuracies indicated that the overall accuracy values were biased by the No-HWA class. Among the properties the average user's and producer's accuracies for the No-HWA class was 83.3% and 82.4%, respectively. In contrast, the average user's and producer's accuracies in the HWA class were 31.1% and 16.3%, respectively. Two factors may account for this. First, a very small percentage of the region of interest was classified as infested (b5%) leaving a much greater chance for areas not infested to be classified correctly and thereby biasing the traditional error matrix. This becomes most apparent when using a random point sampling design; in hindsight using a stratified point sampling design may have produced a less biased accuracy assessment. Second, using the hemlock habitat model in the classification masking procedure likely increased errors at the property scale. As stated above the hemlock habitat map had a wide range of accuracies at the property scale, that error would therefore be incorporated through the masking procedure into the final classification. For example the Rachel Carson Wildlife Refuge and the Russell-Abbott State Forest both had the lowest overall accuracies in the hemlock habitat model assessment, the same trend was observed in the assessment of the HWA classification. Like the hemlock habitat model, the HWA classification likely depicts a moderately accurate estimate of HWA infestations in the region of interest, but its accuracy can vary at larger scales. With that limitation in mind, the true test of the classification model was whether it was (and would be) useful to land managers and surveyors. In order for the HWA classification product to be useful to land managers and surveyors it must not only be accurate, but also improve survey efficiency and reduce costs. Improving survey efficiency would require an increase in HWA detections relative to the number of areas surveyed. In Struble et al. (2011) it was reported that over the course of two years 125 towns (378 sites) were surveyed in Maine, and 115 towns (575 sites) were surveyed

Fig. 7. A map of the study region indicating towns currently known to be infested with HWA.

24

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25

in New Hampshire (between 3 and 5 survey sites per town). The surveys resulted in 23 new towns being designated as HWA positive in each state. Hypothetically, survey planners using a remotely sensed product to direct surveys to sites classified with a high probability of being infested could reduce survey time and costs upwards of 50% (surveying 1–3 sites per town) depending on the survey objectives; however, this would assume a highly accurate classification. In reality, classification maps contain errors of omission and commission. In an effort to reduce errors only pixels with the highest probability of infestation were retained in the final classification product, and the masking procedures limited this further to areas that could support hemlock. Since the hemlock habitat model over-predicted areas of suitable hemlock habitat (i.e. omission errors of ‘non-hemlock’ were greater than ‘hemlock’) it is likely that commission error of the HWA class would be high (i.e. high number sites with false positives). However, the effect of the commission errors could be reduced during the planning of surveys by cross referencing likely-infested areas on the classification map with aerial photographs to confirm hemlock presence. Using the final classification map to identify potentially infested areas for surveying was successful. Surveyed areas were limited to public access areas classified as infested and where a strong hemlock component was confirmed through photointerpretation of Google Earth imagery. It was apparent that strategically targeting hemlock forests identified with a high probability of being infested would likely increase the chances of finding HWA infestations because the greatest accuracy between the classification model and known infestations occurred in coniferous forest types. These targeted surveys resulted in six new towns being designated as infested, two in Maine and four in New Hampshire; one infestation was found nearly 19 km away from any previously reported infestation. Results showed that this product would be useful during survey planning stages by indicating areas that have a high probability of being infested and that targeting those areas would likely increase the success of early HWA detection as well as survey efficiency. In regard to the survey technique, there is likely a lag time of several years from when a stand becomes infested to the time when the HWA population would be large enough to be detectable in the understory. At each test survey location where HWA was found, woolly masses were observed on the previous 2 years of growth. Taking into account that infestations likely start in the canopy from windborne or other dispersion (phoresy) it is probable that the majority of the newly detected infestations were well established prior to the survey. At test locations identified as likely infested, but where no HWA was found through the survey, infestations may have been present in the canopy but HWA populations were not large enough to be found in a ground level survey. To this point, efforts to manage and monitor the spread of HWA solely through ground based surveying will likely always lag behind several years unless remote sensing based early detection technologies and methods are developed and utilized. The methods presented in this paper need refinement, but the results indicate that they are a significant step forward in using satellite based remote sensing to detect low level HWA populations on a landscape scale. Ways to improve the methods include identifying better predictor variables for the hemlock habitat model, increasing reference data sampling intensity and using a stratified sampling design, and surveying (or sampling) the hemlock canopy as well as ground level hemlock branches. In addition, measuring the range of infestation intensities that could be detected through these methods would help to guide what survey techniques work best. These and other improvements would aid in developing tools that could be used by federal or state agencies to monitor the regional expansion of HWA infestations and identify lowlevel or outlier infestations for suppression and eradication efforts.

Acknowledgements This research was a portion of a thesis submitted to the Graduate School at the University of New Hampshire as part of the degree requirements for a Master's of Science in Natural Resources. Funding for this project was provided through the United States Department of Agriculture Forest Service (Forest Health Protection), New Hampshire Space Grant Consortium, University of New Hampshire Farrington Fund Fellowship and College Woods Scholarship, and the Bearcamp Valley Garden Club Scholarship. Marc DiGirolomo, Angela Hammond, Michael Simmons, Michael Bohne, Robert Cooke, Edward Jordan, Kevin Dodds, and Garret Dubois of the USFS provided field support. We would also like to thank David Orwig and the anonymous reviewers for their comments on this manuscript. Appendix A. Supplementary data Supplementary data to this article can be found online at http://dx. doi.org/10.1016/j.rse.2016.12.005. References Anderson, J.R., Hardy, E.E., Roach, J.T., Witmer, R.E., Peck, D.L., 1976. A land use and land cover classification system for use with remote sensor data. U.S. Geological Survey Circular 671 (964), 41. Birth, G.S., McVey, G., 1968. Measuring the color of growing turf with a relfectance spectroradiometer. Agron. J. 60, 640–643. Bonneau, L.R., Shields, K.S., Civco, D.L., 1999. Using satellite images to classify and analyze the health of hemlock forests infested by the hemlock woolly adelgid. Biol. Invasions 1, 255–267. Browning, D.M., Beaupré, S.J., Duncan, L., 2005. Using partitioned Mahalanobis D 2 (K) to formulate a GIS-based model of timber rattlesnake hibernacula. J. Wildl. Manag. 69 (1), 33–44. Chavez, P.S.J., 1996. Image-based atmospheric corrections — revisited and improved. Photogramm. Eng. Remote. Sens. 62 (9), 1025–1036. Clark, J.T., Fei, S., Liang, L., Rieske, L.K., 2012. Mapping eastern hemlock: comparing classification techniques to evaluate susceptibility of a fragmented and valued resource to an exotic invader, the hemlock woolly adelgid. For. Ecol. Manag. 266:216–222. http:// dx.doi.org/10.1016/j.foreco.2011.11.030. Congalton, R.G., 1991. A review of assessing the accuracy of classifications of remotely sensed data. Remote Sens. Environ. 37, 35–46. Congalton, R.G., Oderwald, R.G., Mead, R.A., 1983. Assessing Landsat classification accuracy using discrete multivariate analysis statistical techniques. Photogramm. Eng. Remote. Sens. 49 (12):1671–1678 Retrieved from. http://agris.fao.org/agris-search/ search.do?recordID=US19840105628. Costa, S., Onken, B., 2006. Standardized Sampling for Detection and Monitoring of Hemlock Woolly Adelgid in Eastern Hemlock Forests. US Department of Agriculture, Forest Service, Washington, D.C. Crist, E.P., Cicone, R.C., 1984. A physically-based transformation of thematic mapper data—the TM tasseled cap. IEEE Trans. Geosci. Remote Sens. GE-22 (3):256–263. http://dx.doi.org/10.1109/TGRS.1984.350619. Defries, R.S., Chan, J.C., 2000. Multiple criteria for evaluating machine learning algorithms for land cover classification from satellite data. Remote Sens. Environ. 74, 503–515. Doucette, J.S., Stiteler, W.M., Quackenbush, L.J., Walton, J.T., 2009. A rules-based approach for predicting the eastern hemlock component of forests in the northeastern United States. Can. J. For. Res. 39 (8):1453–1464. http://dx.doi.org/10.1139/X09-060. Dunckel, K., Weiskittel, A., Fiske, G., Sader, S.A., Latty, E., Arnett, A., 2015. Linking remote sensing and various site factors for predicting the spatial distribution of eastern hemlock occurrence and relative basal area in Maine, USA. For. Ecol. Manag. 358:180–191. http://dx.doi.org/10.1016/j.foreco.2015.09.012. ESRI, 2012. ArcMap 10.1. Redlands, CA. Evans, A.M., Gregoire, T.G., 2007. A geographically variable model of hemlock woolly adelgid spread. Biol. Invasions 9 (4):369–382. http://dx.doi.org/10.1007/s10530006-9039-z. Eyre, F.H., 1980. Forest Cover Types of the United States and Canada. Society of American Foresters. Fei, S., Schibig, J., Vance, M., 2007. Spatial habitat modeling of American chestnut at Mammoth Cave National Park. For. Ecol. Manag. 252 (1–3):201–207. http://dx.doi.org/10. 1016/j.foreco.2007.06.036. Fisher, P., 1997. The pixel: a snare and a delusion. Int. J. Remote Sens. 18 (3):679–685. http://dx.doi.org/10.1080/014311697219015. Fry, J., Xian, G., Jin, S., Dewitz, J., Homer, C., Yang, L., Wickham, J., 2011. Completion of the 2006 national land cover database for the conterminous United States. Photogramm. Eng. Remote. Sens. 77 (9), 858–864. Godman, R.M., Lancaster, K., 1990. Tsuga canadensis (L.) Carr. Eastern Hemlock. Silvics of North America: Volume 1, Conifers. Agriculture Handbook 654 Vol. 1. U.S. Department of Agriculture Forest Service, Washington, D.C., pp. 604–612. Gómez, S., Orians, C.M., Preisser, E.L., 2012. Exotic herbivores on a shared native host: tissue quality after individual, simultaneous, and sequential attack. Oecologia 169 (4): 1015–1024. http://dx.doi.org/10.1007/s00442-012-2267-2.

J.P. Williams et al. / Remote Sensing of Environment 190 (2017) 13–25 Goodwin, N.R., Coops, N.C., Wulder, M.A., Gillanders, S., Schroeder, T.A., Nelson, T., 2008. Estimation of insect infestation dynamics using a temporal sequence of Landsat data. Remote Sens. Environ. 112 (9):3680–3689. http://dx.doi.org/10.1016/j.rse. 2008.05.005. Hershey, R.R., Befort, W., 1995. Aerial Photo Guide to New England Forest Cover Types (USDA Forest Service Gen. Tech. Reoprt NE-195). Randor, PA. Jensen, J.R., 2005. Introductory Digital Image Processing. third ed. Pearson Prentice Hall, Upper Saddle River, NJ. Jones, C., Song, C., Moody, A., 2015. Where's woolly? An integrative use of remote sensing to improve predictions of the spatial distribution of an invasive forest pest the Hemlock Woolly Adelgid. For. Ecol. Manag. 358:222–229. http://dx.doi.org/10.1016/j. foreco.2015.09.013. Koch, F.H., Cheshire, H.M., Devine, H.A., 2005. Mapping hemlocks via tree-based classification of satellite imagery and environmental data. In: Onken, B., Reardon, R. (Eds.), Third Symposium on Hemlock Woolly Adelgid in the United States. USDA Forest Health Technology Enterprise Team, Asheville,North Carolina, pp. 104–114. Lefsky, M.A., Cohen, W.B., Spies, T.A., 2001. An evaluation of alternate remote sensing products for forest inventory, monitoring, and mapping of Douglas-fir forests in Western Oregon. Can. J. For. Res. 31:78–87. http://dx.doi.org/10.1139/x00-142. Loveland, T.R., Irons, J.R., 2016. Landsat 8: the plans, the reality, and the legacy. Remote Sens. Environ. 185:1–6. http://dx.doi.org/10.1016/j.rse.2016.07.033. Mayer, M., Chianese, R., Scudder, T., White, J., Vongpaseuth, K., Ward, R., 2002. Thirteen years of monitoring the Hemlock Woolly Adelgid in New Jersey forests. Hemlock Woolly Adelgid in the Eastern United States Symposium. Agricultural Experiment Station and Rutgers University, New Brunswick, NJ, pp. 50–60. McClure, M.S., 1991. Density-dependent feedback and population cycles in Adelges tsugae (Homoptera: Adelgidae) on Tsuga canadensis. Environ. Entomol. 20 (1), 258–264. Orwig, D.A., Foster, D.R., Mausel, D.L., Forest, H., 2002. Landscape patterns of hemlock decline in New England due to the introduced hemlock woolly adelgid. Sciences-New York 1475–1487. Orwig, D.A., Thompson, J.R., Povak, N.A., Manner, M., Niebyl, D., Foster, D.R., 2012. A foundation tree at the precipice: Tsuga canadensis health after the arrival of Adelges tsugae in central New England. Ecosphere 3 (1), 1–16. Phillips, S.J., Dudik, M., Schapire, R.E., 2004. A maximum entropy approach to species distribution modeling. Proceedings of the Twenty-First International Conference on Machine Learning, pp. 655–662. Phillips, S.J., Anderson, R.P., Schapire, R.E., 2006. Maximum entropy modeling of species geographic distributions. Ecol. Model. 190 (3–4):231–259. http://dx.doi.org/10. 1016/j.ecolmodel.2005.03.026.

25

Qi, J., Chehbouni, A., Huete, A.R., Kerr, Y.H., Sorooshian, S., 1994. A modified soil adjusted vegetation index. Remote Sens. Environ. 48 (2):119–126. http://dx.doi.org/10.1016/ 0034-4257(94)90134-1. Rock, B.N., Vogelmann, J.E., Williams, D.L., Vogelmann, A.F., Hoshizaki, T., 1986. Remote detection of forest damage. Bioscience 36 (7):439–445 Retrieved from. http:// www.jstor.org/stable/10.2307/1310339. Rouse, J.W., Haas, R.H., Schell, J.A., Deering, D.W., 1974. Monitoring vegetation systems in the great plains with ERTS. Proceedings of the 3rd Earth Resouce Technology Satellite Symposium, pp. 48–62. Royle, D.D., Lathrop, R.G., 2002. Discriminating Tsuga canadensis Hemlock forest defoliation using remotely sensed change detection. J. Nematol. 34 (3), 213–221. Schriever, J.R., Congalton, R.G., 1995. Evaluating seasonal variability as an aid to covertype mapping from Landsat thematic mapper data in the Northeast. Photogramm. Eng. Remote. Sens. 61 (3), 321–327. Souto, D., Luther, T., Chianese, R., 1995. Past and current status of Hemlock Wolly Adelgid in Eastern and Carolina Hemlock Stands. In: Salom, S.M., Tignor, T.C., Reardon, R.C. (Eds.), Proceedings of the First Hemlock Woolly Adelgid Review. U.S. Department of Agriculture Forest Service, Forest Health Technology Enterprise Team, Charlottesville, Virginia, pp. 9–15. Struble, D., Kanoti, A., Searles, W., Lombard, K., Weimer, J., Burns, B., Souto, D. (2011). Slow-the-Spread Management of Hemlock Woolly Adelgid in Northern New England: A Coordinated Regional Response. Retrieved from http://na.fs.fed.us/fhp/ hwa/reports/final-tristate-redesign report.pdf Trotter, R.T., Shields, K.S., 2009. Variation in winter survival of the invasive hemlock woolly adelgid (Hemiptera: Adelgidae) across the eastern United States. Environ. Entomol. 38 (3):577–587 Retrieved from. http://www.ncbi.nlm.nih.gov/pubmed/19508766. U.S. National Vegetation Classification, 2014. Retrieved from. http://usnvc.org. Vogelmann, J.E., Gallant, A.L., Shi, H., Zhu, Z., 2016. Perspectives on monitoring gradual change across the continuity of Landsat sensors using time-series data. Remote Sens. Environ. 185:258–270. http://dx.doi.org/10.1016/j.rse.2016.02.060. Williams, J.P., Hanavan, R.P., Rock, B.N., Minocha, S.C., Linder, E., 2016. Influence of Hemlock Woolly Adelgid infestation on the physiological and reflectance characteristics of Eastern Hemlock. Can. J. For. Res. 46:410–426. http://dx.doi.org/10.1139/cjfr-20150328. Young, R.F., Shields, K.S., Berlyn, G.P., 1995. Hemlock Woolly Adelgid (Homoptera: Adelgidae): stylet bundle insertion and feeding sites. Ann. Entomol. Soc. Am.