Evaluating Southern Appalachian Forest Dynamics without Eastern Hemlock: Consequences of Herbivory by the Hemlock Woolly Adelgid

Eastern hemlock (Tsuga canadensis Carriére) and the Carolina hemlock (Tsuga caroliniana Engelmann) are ecologically important tree species in eastern North America forests that are currently threatened by the hemlock woolly adelgid (HWA, Adelges tsugae Annand, Hemiptera: Adelgidae). HWA has spread rapidly from its original introduction site into new areas. Once present, HWA kills its hosts over a period of 4 to 10 years leading to a phenomenon that is known scientifically and colloquially as hemlock decline. To date, quarantine, chemical management, and biocontrol efforts have failed to curb the spread of the HWA. As such, forest management efforts are now being redirected towards developing an understanding of the effects of hemlock removal on vegetation dynamics, changes in forest composition, and changes in ecosystem function. In this study, we parameterize a spatially explicit landscape simulation model LAN-DIS II for a specific forested region of the southern Appalachians. Parameterization involves defining the life-history attributes of 37 tree species occupying 11 ecological zones and is based on knowledge of: current vegetation composition data, recent historic management and fire regimes, and life-history traits of each species. The parameterized model is used to explore a simple scenario of catastrophic hemlock mortality likely to occur as a result of HWA herbivory. Our results emphasize that hemlock is an important foundation species. When hemlock is removed from the system, forest composition changes considerably with a greater presence of shade intolerant pine and oak species. Additionally, hemlock removal leads to a period of transient, relatively unstable vegetation dynamics as the forest communities restructure.


Introduction
The eastern hemlock (Tsuga canadensis Carriére) and the Carolina hemlock (Tsuga caroliniana Engelmann) are foundation species in eastern US forests (Ellison et al., 2005).Hemlocks provide habitat for terrestrial wildlife, but also influence surrounding vegetation, soil processes, and aquatic environments.Currently, eastern hemlock forests are under threat from the hemlock woolly adelgid (Adelges tsugae Annand; Hemiptera: Adelgidae) (HWA), an exotic invasive sap sucking insect that feeds on both eastern and Carolina hemlock.Since its introduction in 1953 in Richmond, Virginia, HWA populations have expanded and currently occupy approximately half of the range of native eastern hemlock forests (USDA, 2010).When local HWA infestations occur, trees typically lose vigor over a 4 to 10-year period and then die.Among scientists, professional foresters, and increasingly, the general public, this pattern of expanding HWA populations and gradual tree mortality is often referred to as hemlock decline.The potential impact of HWA on eastern North American forests has been compared to that of chestnut blight (Cryphonectria parasitica) on the American chestnut (Castanea dentata) and Dutch elm disease (Ophiostoma ulmi) on the American elm (Ulmus americana) (Small et al., 2005).
HWA infests hemlocks of all ages and size classes.The insect feeds on sap at the base of needles leading to reduced needle growth, reduced tree vigor, and ultimately host mortality (McClure, 1991).Through the course of an infestation, rapid population growth leads to an excess of highly dispersible offspring that forms new infestations.The susceptibility and mortality rates of hemlocks infested by HWA may be driven by site characteristics and/or needle chemistry but are most likely related to HWA population density (Orwig et al., 2002).Consequently, hemlock decline is best understood hierarchically.First, HWA populations expand into new geographical areas.Second, as local population densities increase, infestations occur gradually through space and time.Third, trees in each infestation lose vigor and die at different rates.Observationally, this results in a gradual decline in hemlock trees and forests, but hemlock mortality is virtually guaranteed once HWA becomes endemic to an area.Over time, it is probable that both eastern hemlock and Carolina hemlock will become extinct in southern Appalachian forests.In turn, hemlock decline is likely to lead to significant changes in the ecology of the surrounding area (Ellison et al., 2005;McClure, 1991;Ross et al., 2003;Snyder et al., 2002) including changes to the decomposition processes in the forest floor of hemlock stands (Cobb et al., 2006;Nuckolls et al., 2009), and changes in the composition and structure of vegetation in its immediate neighborhood (Orwig & Foster, 1998).
Given that large scale hemlock mortality is inevitable, two important forest management research questions become evident: 1) what successional pathways will occur after hemlock mortality?and 2) what type of forest landscape will emerge as a result of hemlock mortality?Researchers are beginning to address these questions at different spatial and temporal scales (e.g., Heard & Valente, 2009;Nuckolls et al., 2009).However, given the dynamic nature of forested landscapes and the large spatio-temporal scale of forest processes, it is essential to provide quantifiable links between the ecology of individual tree species, the dynamics of vegetation through time, and the resulting state of the forest.The goal of this study was to use a spatially explicit forest landscape model (LANDIS II) to simulate and understand the long term impact of HWA on forest composition and succession in a representative southern Appalachian forest landscape (Grandfather Ranger District, Pisgah National Forest, NC).Our specific approach to answering these questions involves two tasks: 1) To parameterize LANDIS II to simulate competition between 36 of the most important tree species in the region and to understand the current or recent historic pattern of vegetation dynamics.
2) To compare and contrast the results from 1) with the vegetation dynamics of simulations where hemlocks are removed from the system.
These comparisons allow evaluation of changes in the vegetation dynamics to follow hemlock removal by HWA and provide insights to into the successional pathways of forest development and the emergent forest landscape after the HWA impact.

Study Area
Our study area was the Grandfather Ranger District (latitude 35˚55'N, longitude 81˚56'W, Figure 1), one of the four ranger districts of Pisgah National Forest in the southern Appalachian mountains.This area covers 777 km 2 and is one of the most diverse forest landscapes in the eastern United States.Simon et al. (2005) identified eleven major ecological zones, ranging from mesophytic assemblages in low-elevation moist coves to oak and pine stands on dry ridges (Figure 1).Three ecological zones (Acidic Cove, Xeric Pine-Oak Heath and Oak Heath, and Mesic Oak-Hickory) together occupy more than 92% of the study site.Acidic Coves occupy low-elevation sites (300 -500 meters) near perennial streams, Xeric Pine-Oak Heath and Oak Heath cover upper slopes and ridges at low elevations (500 -900 meters), and Mesic Oak-Hickory forests inhabit steep slopes at high elevations (900 -1500 meters).Eastern hemlock is one of the most abundant species in these three major ecological zones.The spread of HWA in this forest landscape is expected to cause large scale mortality of eastern hemlock and to generate novel successional dynamics in all ecological zones in subsequent decades.

LANDIS II
LANDIS II (He & Mladenoff, 1999;Mladenoff & He, 1999;Scheller et al., 2007) is a raster-based model designed to simulate the spatial and temporal changes in large landscapes.The model has been used to simulate vegetation changes in a variety of forest landscapes under different disturbance regimes.Specific applications to the southern Appalachians include studies by Cairns et al. (2008aCairns et al. ( , 2008b)); Lafon et al. (2007); Waldron et al. (2007); Xi et al. (2009).LANDIS II consists of three interacting components: maps of the underlying biophysical landscape, vegetation dynamics, and disturbance.The vegetation succession component simulates the birth, death, dispersal, and growth of each species within and between each site (cell) in the landscape.The disturbance module consists of components to simulate wind, fire, insect, and harvest effects.Both the succession and disturbance modules are driven by maps of land types, fire disturbance regimes, vegetation communities, and the initial vegetation structure that describes the spatial heterogeneity of the processes represented by the simulation.For our study area, the maps constitute raster files with a cell size of 90 m × 90 m and 731 × 746 rows and columns respectively.LANDIS II simulates vegetation dynamics and tracks the abundance (and age cohorts) of each tree species in each cell using a 10-year time step.

Model Parameterization
The first step in the LANDIS II parameterization process is to define a landscape based on 11 ecological zones defined for our study area by Simon et al. (2005).Since fire disturbance is one of fundamental driving forces of existing southern Appalachian vegetation patterns, we then parameterized LANDIS II's fire module so that simulated burn areas matched (+/−10%) the historic and current values for each ecological zone based on research by Flatley et al. (2011) and Lafon and Grissino-Mayer (2007).
Additionally, each tree species requires a parameter called an "establishment coefficient".This parameter represents a species' ability to establish under the specific environmental conditions (climate, soil, vegetation composition, and topography) in each ecological zone.In contrast to the basic life history and fire related parameters, the species establishment coefficients are a more abstract parameter that represent the likelihood of a species establishing in a particular ecological zone, based on multiple interacting environmental factors such as climate, soil, topography, and early stage competition.This abstraction makes establishment coefficients more difficult to parameterize directly from the literature.Instead, we estimated the establishment coefficients through an iterative process of model calibration/validation described in the following section.

Estimation of Establishment Coefficients and Model Calibration
An essential objective of this study is to parameterize LAN-DIS II such that the output of the model is representative of the current and recent vegetation composition and dynamics of our study region.This objective allows us to develop baseline scenarios against which changes to the system such as those driven by HWA herbivory can be compared.
A fundamental issue for any complex simulation model is that small differences in parameter values may interact to generate much larger differences in model output.Given the integrative, large scale nature of LANDIS II, we have taken the approach that any parameter estimates from the literature, expert knowledge, or any other data source are subject to error; but that some are likely more accurate than others.In the case of LANDIS II, many of its parameters have a direct, ecological meaning such that they can be measured or estimated using literature or expert knowledge.Table 1 illustrates ecologically meaningful parameters used in LANDIS II that were estimated from the literature and/or expert knowledge.However, LAN-DIS II also requires a set of parameters called "establishment coefficients" that are essential to the validity of the model, but are difficult to estimate directly from the literature.The objective of the calibration procedure is to infer establishment coefficients based on a priori information about known parameters (including their likely error); known recent historical disturbances; and data outlining current forest composition.The goal is to derive a fully parameterized LANDIS II model that can be used to simulate southern Appalachian vegetation dynamics.
To calibrate process LANDIS II, the model was run iteratively, making fine-scale adjustments to the establishment coefficients for each species.To maintain the integrity of our a priori parameter values, we only adjusted life-history parameters for species that were poorly defined by information in the scientific literature.Carolina Vegetation Survey (CVS) data (266 plots located in the study region) were used to quantify current species abundances.We considered the model calibrated when differences between the observed and the simulated abundances of 70% of the species in this landscape were within +/−15% in the three major ecological zones.All calibrations (and simulations) were performed using a randomized initial species composition map to ensure that our simulated vegetation dynamics consistently tended towards real world species assemblages, rather than just maintaining an initial forest composition.
Although previous LANDIS studies used a calibration procedure, only a few studies have calibrated vegetation communities against actual vegetation data (e.g., Shifley et al., 1998).As a result, relative to other applications of LANDIS II, the parameterization approach used in this study is among the most thorough and systematic reported.In large part, this process is necessitated by the diversity of southern Appalachian forests and the use of a defined and real study area.The net result of this process is a realistic, data driven baseline scenario for current and recent historical southern Appalachian vegetation dy-namics against which the effects of removing hemlock from the system can be compared.

Simulating Scenarios and Processing Model Outputs
To run a simulation the following input data are required: a map of the ecological zones (Simon et al., 2005); a map of randomized initial vegetation abundance; life history parameters for each species; the establishment coefficients of each species in each ecological zone; and fire parameters for each ecological zone.Though species in all eleven ecological zones are simulated and results were analyzed in our study, we focus on the three extensive ecological zones: Acidic Cove, Xeric Pine-Oak Heath and Oak Heath, and Mesic Oak-Hickory.Detailed results of these three zones are presented along with concise results from other ecological zones.The rational for focusing on the three ecological zones was that 1) together they represented approximately 92% of the forestland in the study site and hemlock was one of the most abundant species present, 2) the results from the other ecological zones were qualitatively similar (Figure 2), and 3) emphasis was placed on general vegetation dynamics impact occurring at a broad landscape scale.
Once parameterized and calibrated, scenarios were simulated with and without hemlock.To model the second scenario, hemlock establishment coefficients for each ecological zone were set to zero, thereby simulating instantaneous hemlock mortality.This basic scenario was chosen for simplicity and because detailed hemlock decline data is currently lacking.Simulations were run for 1200 years, at a time-step of 10 years resulting in a run time of approximately 12 minutes per simulation.The output of LANDIS II is a set of raster maps, each showing the maximum age of each species cohort in each cell of the simulated area at each time step (10 years).The first 400 years of the simulations were discarded as these represent a transient period during which initial pattern of random species abundance moves towards a consistent repeating pattern of nontransient vegetation dynamics.To extract temporal and nontemporal species abundance for the relevant ecological zones, simulation output was summarized using LANDISVIEW (Birt et al., 2009) and Microsoft Excel.

Results
The average difference in the abundance of each species within the landscapes with hemlock and without hemlock is illustrated in Figure 2. Here, species abundance is defined as the proportion of cells in each ecological zone occupied by each species (LANDIS II does not provide a means for calculating the abundance of each species in each cell).Figure 2 is derived by averaging the abundance of each species in each ecological zone through time (between years 400 to 1200, i.e., once the simulation has reached a steady state); and then averaging the percent increase in abundance between hemlock and nonhemlock forest.Figure 2 also shows the average abundance of hemlock found in baseline simulations.The removal of hemlock from the simulations results in an increase in the abundance of all species in all ecological zones roughly in propor- tion to the abundance of hemlock in the baseline scenario.
In Figure 2, ecological zones marked with an asterisk represent 95% of the study area.Figure 3 provides a more detailed view of species abundance within these three ecological zones.In simulations without hemlock, the abundance of most species increase between 25% to 74%.Together, Figures 2 and 3 show that the greatest change in species abundance occurs in Mesic zones (Acidic Cove and Mesic Oak-Hickory zones) with least change in Xeric zones (Xeric Pine-Oak Heath and Oak Heath).In the Acidic Cove and Xeric Pine-Oak Heath and Oak Heath ecological zones, hemlocks are mostly replaced by 4 species (white pine, chestnut oak, red oak, and yellow poplar).In contrast, the Mesic Oak-Hickory ecological zone shows relatively uniform increases in abundances across 13 species.Figure 3 also illustrates differential response to hemlock removal for two important understory shrubs, great rhododendron and mountain laurel.In all three ecological zones, great rhododendron shows very little change in abundance between simulations with and without hemlocks.In contrast, in all ecological zones, mountain laurel increases by 74%, 22% and 91% in the Acidic Cove, Xeric Pine-Oak Heath and Oak Heath, and Mesic Oak-Hickory ecological zones respectively.
Figures 4-6, illustrate the temporal pattern of vegetation dynamics in simulations with and without hemlock in the three major ecological zones.In simulations with hemlock, the relative changes of each species follow a smoother trajectory and reach a stable state more quickly than in the non-hemlock scenario.In scenarios without hemlock, the abundance of each species fluctuates over greater amplitudes and for longer time periods before reaching a stable state.This phenomenon is related to the structure and parameterization of the model.By removing hemlock from the system, the relative competitive abilities of the remaining species change and radically alter the vegetation dynamics.The extent to which hemlocks exert a stabilizing effect on the forest landscape is illustrated in Figures 4-6.In addition to this general trend, the figures also illustrate that, for any species, fluctuations in abundance occur over relatively long time-periods (between 200 and 400 years).

Discussion and Conclusions
Removal of hemlock leads to an increase in the abundance of the remaining 36 species in all ecological zones.This general increase in species abundance is approximately proportional to the baseline abundance of hemlock in each ecological zone and   leads to a substantial change in the forest landscape.The simulation results highlight two important dimensions to this disturbance.The first involves the long-term, stable dynamics of the landscape and the types of species that are expected to become dominant in the absence of hemlock.As Figure 3 shows, the removal of hemlock from the system generally leads to an increase of shade intolerant species such as pines and oaks.The second dimension of this change involves the temporal fluctuations in forest composition that may occur as the forest responds to this disturbance (Figures 4-6).Catastrophic removal of hemlock from the system leads to a restructuring of forest composition that occurs over 200+ years.In forests without hemlock, short term fluxes in species abundance are much greater.This relative instability illustrates the importance of models such as LANDIS II for forest planning, especially for forests with diverse communities of interacting species.Without the benefit of integrative models such as LANDIS II, it is tempting to interpret forest patterns in a single dimension (for example by looking at current landscape composition and configuration).An important contribution of models such as LANDIS II is to illustrate the link between individual species ecology, and the short term and/or long term vegetation dynamics that may arise under different forest management and/or disturbance scenarios.
In the broadest sense, these results show that hemlock is a dominant, foundation species that significantly affects forest composition and structure.Through a combination of longevity and large size, it inhibits the establishment of shade-intolerant species such as oak and pine.The model results also suggest that this dominance is exacerbated by suppressed fire regimes.In comparison to pine and oak species, hemlocks are relatively fire intolerant.Fires create gaps in forest canopies, and light conditions suitable for the growth and establishment of shadeintolerant pine and oak species.However, our results show that nearly all species (not just shade intolerant pine and oaks) increase as a result of hemlock removal.It is likely that this general increase in species abundance, drives increased competetion, and contributes to the long term instability of vegetation dynamics in simulations without hemlock.
The results of this analysis largely concur with other studies.For example, Nuckolls et al. (2009) proposed that, over short temporal scales, hemlocks will either be replaced by a mix of advanced regeneration and early successional tree species such as red maple, black birch, and yellow poplar, or great rhododendron.Although, in our study, mountain laurel increases in abundance, great rhododendron does not.Heard and Valente (2009) studied a catastrophic hemlock decline in the mid-Holocene using fossilized pollen, and found increases in maple, birch, beech, and oak.HWA hemlock decline is already being compared to the impact of chestnut blight on American chestnuts between 1900 and 1940.In line with our results for hemlock, Day and Monk (1974) studied the changes in species dynamics following chestnut mortality in a southern Appalachian watershed and found that in most cases, chestnuts were replaced by multiple species, rather than a single co-dominant species.They also suggest that these large scale disturbances may also have allowed "pest" species such as great rhododendron and mountain laurel to establish and spread.
In our future work, we intend to use the model to explore scenarios under more realistic patterns of hemlock mortality, and to investigate management scenarios with the potential to mitigate hemlock decline.Given the destabilizing effect of hemlock removal in this study (which uses the simplest assumption of catastrophic, instantaneous hemlock mortality), it is important to investigate how different spatio-temporal patterns of hemlock decline may affect short and long-term vegetation dynamics.Our current model structure provides the opportunity to explore such scenarios, and ongoing, independent research to measure and characterize the current pattern of Hemlock decline will provide the inputs necessary to address this issue (Clark et al., 2012).Similarly, the diversity of species, and the complex social and ecological function of forests in this area suggest that the model can be used to address forest management options with potential to mitigate the effects of Hemlock decline.For example, the HWA resistant Chinese hemlock (Tsuga chinensis (Franch.)E. Pritz.)(Tredici and Kitajima, 2004) has been proposed as a direct replacement for native hemlocks.Similarly, modification of the current fire suppression regime has been proposed as a method to control the expansion of mountain laurel and possibly increase the dominance of fire-tolerant species such as oaks and pines.Both mitigation efforts could be assessed by building upon the parameterization of LANDIS II outlined in this paper.Furthermore, they can be tested by using a specific culturally, socially, and ecologically important landscape (i.e., the Grandfather Ranger District).

Figure 1 .
Figure 1.Location of the study area, Grandfather Ranger District, Pisgah National Forest, North Carolina.

Figure 2 .
Figure 2. Average difference in species abundance between scenarios with hemlock and without hemlock for each ecological zone (except Spruce-Fir zone).Gray bars show the original (baseline) abundance of hemlock.Black bars show the average change in abundance of all species without hemlock.

Figure 3 .
Figure 3.Comparison of species abundance between forest landscapes with and without hemlock.A. in Acidic Cove; B. in Xeric Pine-Oak Heath and Oak Heath; C. Mesic Oak-Hickory.Note that only species that have greater than 20% abundance in any ecological zone are included.Gray bars show the original (baseline) abundance of hemlock.Black bars show the average change in abundance of all species without hemlock.

Figure 5 .
Figure 5. Species abundance in Xeric Pine-Oak Heath and Oak Heath.

Table 1 .
Species life history characteristics.
b Non-native invasive species; c Shrubs.