Simulated effects of forest management alternatives on landscape structure and habitat suitability in the Midwestern United States

Understanding the cumulative effects and resource trade-offs associated with forest management requires the ability to predict, analyze, and communicate information about how forest landscapes (1000s to>100,000 ha in extent) respond to silviculture and other disturbances. We applied a spatially explicit landscape simulation model, LANDIS, and compared the outcomes of seven forest management alternatives including intensive and extensive even-aged and uneven-aged management, singly and in combination, as well as no harvest. We also simulated concomitant effects of wildfire and windthrow. We compared outcomes in terms of spatial patterns of forest vegetation by age/size class, edge density, core area, volume of coarse wood debris, timber harvest, standing crop, and tree species composition over a 200-year simulation horizon. We also used habitat suitability models to assess habitat quality for four species with diverse habitat requirements: ovenbird (Seiurus aurocapilla), prairie warbler (Dendroica discolor), hooded warbler (Wilsonia citrina), and gray squirrel (Sciurus carolinensis). Management alternatives with similar levels of disturbance had similar landscape composition but different landscape patterns. The no-harvest scenario resulted in a tree size class distribution that was similar to scenarios that harvested 5% of the landscape per decade; this suggests that gap phase replacement of senescent trees in combination with wind and fire disturbance may produce a disturbance regime similar to that associated with a 200-year timber rotation. Greater harvest levels (10% per decade) resulted in more uniform structure of small or large patches, for unevenor even-aged management, respectively, than lesser levels of harvest (5% or no harvest); apparently reducing the effects of natural disturbances. Consequently, the even-aged management at the 10% level had the greatest core area and least amount of edge. Habitat suitability was greater, on average, for species dependent on characteristics of mature forests (ovenbird, gray squirrel) than those dependent on disturbance (prairie warbler, hooded warbler) and habitat suitability for disturbance dependent species was more sensitive to the management alternatives. The approach was data-rich and provided opportunities to contrast the large-scale, long-term consequences for management practices from many different perspectives. Crown Copyright # 2006 Published by Elsevier B.V. All rights reserved.


Introduction
Management of forest landscapes can benefit from the ability to predict and assess the long-term, large-scale consequences of natural and anthropogenic disturbances (or their absence) on forest structure, species composition and the related spatial patterns of forest vegetation. Such information is key to understanding how management alternatives are likely to affect wildlife habitat, timber, recreation opportunities, species diversity, landscape diversity, and a host of other products, amenities, and ecological services that forests provide. All these factors are affected by the current and future condition of forest vegetation at site, stand, and landscape scales. On public lands there is an additional need to effectively communicate the expected outcome of various management alternatives.
Simultaneous consideration of all these factors mandates working at the landscape scale -typically thousands to tens of thousands of hectares in spatial extent. In most situations, landscape-scale field experiments are impractical, but spatially www.elsevier.com/locate/foreco Forest Ecology and Management 229 (2006) [361][362][363][364][365][366][367][368][369][370][371][372][373][374][375][376][377] explicit computer simulation models can effectively provide such a landscape perspective (Mladenoff and Baker, 1999). Parameterization of such models for specific ecological conditions can be difficult, and data requirements for model implementation are often demanding. Landscape simulation models, however, are often the best tools available to predict future forest conditions and provide perspective on long-term, large-scale outcomes of management decisions. Maps of projected forest conditions can illustrate general landscape patterns through time and provide data needed to assess future impacts on wildlife habitat, aesthetics, large-scale biodiversity, and a host of other factors that depend on the spatial arrangement of landscape features. Maps of projected conditions are also useful for illustrating and discussing management alternatives (e.g., Gustafson et al., 2000;Zollner et al., 2005).
Selection of forest management methods and the level of harvest or rotation period are major factors affecting forest landscapes. For example, the application of even-aged versus uneven-aged management has implications for tree species composition, stand structure, landscape structure and wildlife populations in Midwestern, oak-dominated forests (Thompson et al., 1995;Johnson et al., 2002;Dey, 2002). Forest management planning on public lands in the United States often is contentious because of the important consequences of management decisions on forest landscapes and implications for the array of benefits that society expects from public forests. Although harvest (or its absence) is the greatest contemporary disturbance process affecting U.S. forests, other ever-present factors can also greatly affect landscape change and wildlife (e.g., weather, fire, invasive insects, disease, etc.) (Brawn et al., 2001).
In this paper we use a spatially explicit landscape simulation model, LANDIS (Mladenoff et al., 1996;He et al., , 2005Mladenoff and He, 1999), to simulate seven forest management alternatives for a 71,142-ha forested landscape in the Missouri Ozarks. We apply landscape-scale habitat suitability models (Larson et al., 2003(Larson et al., , 2004 to assess the effects of these landscape changes on wildlife habitat. This region is one of few locations in the Midwest where detailed information on ecological land types, forest type and size class, and wind and fire disturbance patterns exist for a large landscape. Our objectives are to (1) demonstrate the utility of this approach to forest planning and management, (2) to draw general conclusions about long-term and large-scale effects of forest management alternatives in oak-dominated forests, and (3) quantify the impact of alternative forest management practices on wildlife habitat quality. The spatially explicit nature of the model allowed us to contrast management alternatives over time in terms of forest size structure, patch size, length of edge habitat, spatial juxtaposition, timber harvest, residual timber, down wood, wildlife habitat suitability, and other metrics.

Study site
The study area is a 71,142 ha portion of the Mark Twain National Forest in the Missouri Ozarks (Fig. 1). We chose this region because it provided an extensive mapped landscape including hypsography, ecological land types (Miller, 1981), stand boundaries (i.e., contiguous operational management units 1-20 ha in size), management area boundaries (i.e., thematic management zones thousands of hectares in size and often spatially discontinuous), and an inventory of initial vegetation conditions (age and forest cover type) ( Table 1, Fig. 1). The area also has a well-documented fire history (Westin, 1992;Guyette, 1995;Guyette et al., 2002) and local information on wind disturbance (Rebertus and Meier, 2001).

Simulating vegetation change
We applied the LANDIS software (version 3.6) (Mladenoff et al., 1996;He et al., , 2001He et al., , 2005Mladenoff and He, 1999) to simulate forest vegetation response to disturbance by timber harvest, wind, and fire. In LANDIS, a landscape is organized as a mapped grid of cells (or sites), with vegetation information stored as attributes for each cell. Within each cell LANDIS represents the forest vegetation as a matrix with the presence or absence of tree species (or species groups) recorded by 10-year age classes. LANDIS simulates four spatial processes (fire, windthrow, harvesting, and seed dispersal) that affect the projected species composition and age structure of individual cells and, in aggregate, of the landscape as a whole. LANDIS and its various modules are described elsewhere in greater detail Gustafson et al., 2000). We followed the general approach of Shifley et al. (1997 who previously applied LANDIS to simulate forest landscape change on a relatively small (3216 ha) Ozark landscape embedded within the landscape used in this study. Overstory vegetation on any given 30 m by 30 m pixel on the landscape was represented by the presence or absence of trees in four species groups in 10-year age classes. We used the following four species groups that in combination comprise nearly 80% of the basal area of mature forests in the region: the white oak group (Quercus alba L., Q. stellata Wangenh., Q. muehlenbergii Engelm.), the black oak group (Q. velutina Lam., Q. coccinea Muenchh., Q. rubra L.), the shortleaf pine group (Pinus echinata Mill. and Juniperus virginiana L.,), and the maple group (Acer rubrum L. and A. saccharum Marsh). Hickories (Carya spp.) comprise the majority of the remaining basal area; they occur ubiquitously across the landscape at low frequency, and they were not modeled explicitly.
We initially populated each pixel in the landscape with one of the four species groups based on a random draw from observed species probability distributions by age class (seedling or sapling, age 0-29 years; pole, age 30-59 years; and sawtimber, age 60 years), forest cover type (shortleaf pine, oak-pine, oak-hickory, black-scarlet oak, oak-gumcypress, elm-ash-cottonwood, maple-beech), and ecological land type (south and west slopes, north and east slopes, ridge tops or upland flats, upland waterways, floodplains or low terraces, side slopes on limestone, or glades). We derived those species probability distributions from two other sources of detailed field inventory data collected in close proximity to our study area: the Missouri Ozark Forest Ecosystem Project (Shifley and Brookshire, 2000) and forest inventory and analysis data collected by the USDA Forest Service (Miles, 2005;Miles et al., 2001). Based on these data sources, for any given forest cover type, stand age class, and ecological land type on the initial landscape map we were able to estimate the relative frequency of trees in the white oak, black oak, shortleaf pine, and maple species groups and place them on the map in the proper proportion.
We examined seven forest management alternatives ( Table 2) that encompass the range of timber harvesting practices likely to be considered for public forest management. These included no harvesting (no harvest), extensive and intensive even-aged management (EAM 5%, EAM 10%, respectively), extensive and intensive uneven-aged management (UAM 5%, UAM 10%, respectively), and an extensive and intensive mix of even-and uneven-aged management (mixed 5%, mixed 10%) ( Table 2). The percentages refer to the proportion of the management area regenerated by timber harvest each decade. Moreover, through recognition of established management areas, the simulations take specific account of areas of the Mark Twain National Forest that are permanently reserved from any type of timber harvest. For the mixed, even-and uneven-aged management alternatives we varied harvest techniques by ecological land type in accordance with local practices (i.e., even-aged management on ridges, south slopes, west slopes and upland drainages; uneven-aged management on other land types) ( Table 2). In reality, the landscape will be managed using a mixture of practices. The practices that we examined bracket the likely range of outcomes and this approach efficiently contrasts differences among the management practices. Timber harvest in LANDIS is simulated using algorithms described by Gustafson et al. (2000).
We set the mean fire-free interval for the landscape to approximately 415 years. Thus, on average a given point on the landscape would burn once every 415 years. We based this level of fire disturbance on the reported frequency of wildfires on state and federal lands including reported flame heights (USDA Forest Service wildfire database, Westin, 1992) coupled with a published model predicting tree mortality based on tree diameter and flame height (Loomis, 1973). This fire-free interval corresponds to the approximate frequency of fires with flame heights greater than 1.2 m. These fires are likely to kill trees 18 cm dbh in the black oak or maple species groups and trees 10 cm dbh in the white oak and shortleaf pine groups. Even in sawtimber age class such fires are likely to create openings in the forest canopy large enough for new trees to regenerate or for advance reproduction to grow into the forest canopy. We did not simulate fires of lesser intensity (flame heights less than 1.2 m) because we assumed they would have relatively little impact on forest structure for the majority of forest which is predominantly in the pole and sawtimber age classes. This is a compromise that was necessitated by imperfect knowledge of fire effects and the related complexity of modeling surface fires that are not stand replacing. This scenario assumes continuation of the current practice of active fire suppression.
We set the mean return interval for wind disturbance at 800 years based on data from Rebertus and Meier (2001). This corresponds roughly to the interval between blowdowns creating openings greater than 0.05 ha in size (i.e., greater than half of the 0.09 ha pixel size in the study). Windthrow of individual trees (or small groups) occurs frequently, but we could not effectively model events that are smaller than the 30 m by 30 m (0.09 ha) pixel resolution used to depict the landscape.
Processes such as disturbance and regeneration are stochastic (probabilistic) in LANDIS, and they are simulated by random draws from probability density functions that define the range and frequency of possible outcomes for a particular disturbance or regeneration event. Repeated simulation runs based on a different sequence of random draws (i.e., based on different random number ''seed'') will result in a different simulation outcome. Consequently, we ran five simulations of each management alternative and evaluated how they differed due to inherent stochasticity of the modeling process.

Landscape analysis
For each management alternative (Table 2), we simulated landscape characteristics for 200 years and retained output maps for tree species and age classes. We combined the 10-year age classes into four size classes: seedling (age 0, 10 years), sapling (age 20, 30 years), pole (age 40, 50 years), and sawtimber (age 60 years). We calculated total area for each forest size class and expressed it as a proportion of the total forest area. We used the program FRAGSTATS (McGarigal and Marks, 1995) to calculate patch size, total core area (the area of a given patch that is greater than 60 m from the patch edge), and edge density (m/ha) based on maps of the four forest size classes. We also recomputed edge for two size classes representing open or young forest (seedling plus sapling size classes combined) versus closed forest (pole plus sawtimber classes combined). We also computed some of these measures for the seedling size class versus all older size classes combined so we could focus on patterns of canopy gaps created by recent disturbance or timber harvest. We summarized area disturbed by harvest, fire or wind from LANDIS output files. We estimated timber harvest volume and residual timber volume using a local, age-based volume table described in . We estimated volume of down wood with the age-based formula from Spetich et al. (1999). We used maps, graphs, and simple summaries to compare changes in stand and landscape characteristics over the 200-year simulations among the seven management alternatives. Except where otherwise specified, reported results are the mean of five simulations for each alternative.

Habitat suitability
We use landscape-scale habitat suitability models to assess the impacts of management alternatives on four wildlife species: ovenbird (Seiurus aurocapilla), prairie warbler (Dendroica discolor), hooded warbler (Wilsonia citrina), and gray squirrel (Sciurus carolinensis). They represent a late-successional, edge-sensitive species; an early-successional species; a gap-dependent species; and a mast-dependent species, respectively. The habitat suitability models are described by Larson et al. (2003) and utilize raster GIS files such as those generated by LANDIS. Habitat suitability index (HSI) values represent an index of habitat quality that is assumed to vary linearly from 0 (non-habitat) to 1 (the best habitat). No direct relationship can be made between HSI values and density or viability without other supporting data, and we do not attempt this in this manuscript (but see Larson et al., 2004). Habitat suitability models compute a HSI value for each pixel, a number between zero and one that indicates relative habitat quality. We present means and medians of HSI values for the entire landscape as well as selected maps of individual pixel values.
To demonstrate the impacts of forest management on habitat quality we report habitat conditions for simulation years 20 and 200 for four management alternatives: no harvest, EAM 10%, UAM 10%, and mixed 10% (Table 2). We selected management alternatives that represented the greatest departure from the noharvest alternative because they are most likely to demonstrate

Forest size class distribution
The initial landscape was predominantly populated by forest in the pole and sawlog size classes. In the early decades of simulation the proportion of area in the sawlog size class increased as sites initially in the pole size class matured and moved to the sawlog size class. Over several decades as the disturbance regimes (harvest, wind, fire) were consistently implemented, the proportion of area by size class equilibrated. This occurred after about 70 years of simulation for management regimes that harvested 10% of the area each decade and after about 120 years of simulation for management regimes that harvested 5% of the area per decade (Fig. 2). For a given percent harvest per decade (5 or 10%) the proportion of the landscape in various size classes (seedling, sapling, pole or sawtimber) over time was similar, regardless of the management practice (even-aged, uneven-aged or mixed). Despite these similarities in the total area by age class for the landscape (Fig. 2), there were obvious differences in the spatial arrangement of size classes on the landscape (Fig. 3). The evenaged management regimes produced even-aged patches averaging about 7 ha in size; the uneven-aged regimes produced a landscape of intermixed age classes with individual age cohorts generally smaller than 0.3 ha in size. These differences are clearly visible in the mapped results of tree size classes (Fig. 3). Differences in the spatial arrangement of vegetation size classes are reflected in the greater edge and lesser core area associated with uneven-aged treatments, as illustrated in the following sections.
The no-harvest management alternative produced an age class distribution with proportions nearly identical to those where 5% of the area was harvested each decade (Fig. 2). Under those low-impact scenarios roughly one-third of the total area was split among the seedling, sapling and pole size classes. The no-harvest alternative maintained younger size classes on the landscape by a combination of wind disturbance, fire disturbance, and gap-scale replacement of senescent trees with new trees in young age classes.

Edge
The greatest length of edge resulted from the uneven-aged management scenario that harvested 10% of the area each decade (UAM 10%, Fig. 4a). At the end of 200 years of simulation, that alternative produced three times as much edge habitat per hectare as the even-aged management alternative with the same harvest intensity (EAM 10%).
After 200 years of simulation, the least edge (257 m/ha) was associated with the most intensive even-aged management practice (EAM 10%). The no-harvest management alternative produced the least edge after the first 20 years of simulation ( Fig. 4a) because the majority of the landscape was in the sawlog size class (Fig. 2). But over time under the no-harvest alternative, the amount of edge increased due to increases in dispersed patches of forest regeneration resulting from wind disturbance, fire disturbance, and gap-scale replacement of senescent trees. After 200 years of simulation the quantity of edge per hectare for the no-harvest scenario was similar to the UAM 5% and mixed 10% scenario. When we compared the length of edge for young forest (seedling and sapling size class) versus mature forest (pole and sawlog size classes) the total length of edge was less but the relative values among treatments were similar to the results based on all four size classes (Fig. 4b).

Core area
Core area generally decreased as the amount of edge increased. Timber harvest practices had a large influence on the core area statistics. After 200 years of simulation, the two evenaged scenarios produced the largest quantity of core area (Fig. 5). In contrast, the no-harvest scenario had relatively little core area after 200 years of simulation (Fig. 5) because aging tree cohorts died (i.e., on individual pixels) and were replaced by younger cohorts. This created large expanses of unevenaged forest that, due to the interspersion of multiple age classes, did not meet our definition of a core area (i.e., total area in a patch of a single forest size class that is >60 m from an edge). The situation was similar for the uneven-aged management scenarios. Through group selection harvesting, the unevenaged management scenarios generated many small patches of trees each decade and reduced the occurrence of large areas in a single age class.

Tree species composition
We initially populated the landscape with one species per pixel which was a relatively simple way to establish the dominant tree cover in the correct proportion for a wide range of initial stand conditions across the landscape. Over the course of the simulation new species and age cohorts were regenerated on the landscape. Following simulated regeneration, many pixels supported multiple species, which is realistic for a 0.09 pixel size and typical of comparison sites such as the Missouri Ozark Forest Ecosystem Project (Shifley and Brookshire, 2000). Under most scenarios, trees in the white oak group expanded to occur on >90% of the pixels (Fig. 6). The proportion of sites with trees in the black oak group and pine groups also increased over time. The increase for these species was greater on sites with 10% harvesting per decade than on those with 5% harvest or no harvest. The proportion of sites with the shade-tolerant maple group increased slightly under the no-harvest scenario, remained nearly constant under the EAM 5% and UAM 5% scenarios, and declined under the EAM 10% and UAM 10% scenarios. By year 200 of every scenario the white oak group, the black oak group and the pine group occurred (represented by at least one tree of each species group) on most of the 790,462 pixels (0.09 ha in size) tracked by the model.

Timber harvest and down wood
The total timber harvest per decade and the residual volume of standing timber varied with the intensity of harvest (Fig. 7.) The scenarios that harvested 10% of the area per decade removed roughly 400,000 m 3 of growing stock per decade and left 4.2 million m 3 of standing volume. For the 5% harvest scenarios the harvest volume dropped to roughly 200,000 m 3 per decade with a residual standing volume of about 4.5 million m 3 . Standing wood volume for the no-harvest scenario peaked at 5.1 million m 3 in decade 6 and stabilized at a level of roughly 4.7 million m 3 for the remaining 140 years of simulation.
The harvested volume decreased in decades 8 through 10 for the uneven-aged 10% harvest regime and the mixed 10% harvest regime (a mixture of even-and uneven-aged harvest practices). A similar pattern of declining harvest volume emerged at decade 18 for the 5% uneven-aged and mixed harvest scenarios.
The model of down wood volume per acre that we applied to the age classes on the landscape predicted high levels of down wood for stands in the years immediately following a standinitiating event (e.g., harvest, fire, or blowdown). That volume then decreased through stand age 80; after age 100 this trend reversed and the volume of down wood increased rapidly for age classes older than 100 years. The estimated volume of down wood (10 cm in diameter) was greatest under the no-harvest scenario (with the oldest mean age across the landscape) and least under the 10% harvest scenarios (youngest mean age across the landscape) (Fig. 7).

Fire and wind disturbance
In the simulations, young forest and forest experiencing simulated blowdown events had a relatively high probability of fire-related tree mortality. In all situations, younger trees had a greater probability of fire-caused mortality, and trees killed by simulated fire could resprout. Older forests experienced fewer simulated fire events and had less damage when fires occurred. The scenarios with 10% harvest per decade had mean fire-free intervals of roughly 400 years (Table 3). For scenarios with 5% harvest per decade the mean fire-free interval increased to approximately 500 years. For the no-harvest scenario the firefree interval exceeded 700 years. Total burned area varied Fig. 5. Core area by size class in simulation year 20 (white bars) and 200 (shaded bars) under seven forest management scenarios in southern Missouri. Core area is the total area by size class (seedling, sapling, pole, or sawlog) that is at least 60 m from an edge boundary between any two adjacent size classes. among scenarios (Table 3), but it also varied for repeated runs within a single scenario. For example, among five repeated runs of the no-harvest scenario, the total burned area varied from 212 to 395 ha and the corresponding mean fire-free interval varied from 703 to 816 years. The size of individual simulated fire events ranged from 1 to more than 700 ha; individual fire sizes followed a negative exponential frequency distribution with many small fires and few large fires.
The mean wind-damage-free interval generally ranged from approximately 1200-2200 years. Mean overstory blowdown for each decade of simulation was 440 ha (Table 3). Like fire disturbance, wind disturbance was modeled as a stochastic process and the location and total area of wind disturbance events varied among treatments and among repeated runs for a single treatment (Table 3). In general, older forests were subject to greater wind damage and less fire damage than younger forests.

Variation among multiple simulation runs
The LANDIS model simulates seedling success, fire, sprouting, wind damage, and a variety of other processes as stochastic events with the specific outcome in each case determined by drawing from a probability distribution. The five repeated simulation runs for each management alternative were initiated using a different sequence of random numbers and had Fig. 6. Percent of pixels where a tree species group is present (regardless of age) under seven simulated forest management scenarios in southern Missouri. We populated the initial landscape (year 0) with one species per 0.09 m pixel. Over time the number of species on most sites increased to two or more through simulated regeneration, succession, and response to disturbance. Hence, as the simulations progressed species in the black oak, white oak, and pine groups were each present across 60 to >90% of the landscape. different simulation outcomes due to chance. Generally the differences among multiple runs for a single treatment were most apparent in the spatial location of disturbance and regeneration events and in the spatial arrangement of species. For the variables that were summarized and expressed as landscape means, the coefficient of variation for repeated runs was generally less than 5% (Table 4).

Habitat suitability
The mean HSI is a measure of the average habitat quality across the landscape while the median HSI value gives additional insight into the distribution of HSI values (Table 5). Maps of HSI values provide additional insight into factors affecting suitability. The distribution of high-suitability pixels for ovenbirds and gray squirrels was greatly affected by the distribution of older forest (Fig. 8). Prairie warblers are earlysuccessional species and were most affected by the distribution of young forest that met their minimum requirements for habitat patch size (Fig. 8). Hooded warblers required juxtaposition of old and young forest so high-suitability pixels occurred around the edges of disturbance patches (Fig. 8).

Anticipated versus emergent results
Many of the general trends that were simulated (e.g., that greater harvest intensity produces a younger mean landscape age as well as greater susceptibility to fire damage or that uneven-aged management results in more edge per unit area) are intuitive. That is comforting because the results qualitatively validate the design and calibration of the model with respect to the ecological processes simulated. Virtually none of the modeled outcomes for a given management alternative (e.g., the relative change in area by age class or size class, fire susceptibility, or edge habitat) could be estimated quantitatively without the use of a simulation model. Model results quantify the relative change associated with alternative disturbance (harvest) regimes. The relationships serve to provide (1) provisional guidelines for management decisions and (2) hypotheses of effects that will be subject to field testing via long-term experiments such as the Missouri Ozark Forest Ecosystem Project (Shifley and Kabrick, 2002).
Some results that emerged from the simulation would have been difficult to anticipate without the simulation methodology. One such outcome was the similarity of the size class distributions for the no-harvest scenario to the EAM 5%, UAM 5%, and mixed 5% scenarios. After 120 years under the  The mean fire-free interval is the number of years it would take to burn an area equivalent to then entire landscape or, equivalently, the mean number of years between repeated fires at one location. The mean wind-free interval is analogous to mean fire-free interval.
no-harvest scenario the seedling, sapling, pole and sawlog size classes equilibrated at about 7, 10, 11, and 72%, respectively (Fig. 2). These proportions were similar to those for the 5% harvest regimes. This suggests that gap phase replacement of senescent trees in combination with wind and fire disturbance may produce a disturbance regime similar to that associated with harvesting 5% of the landscape per decade -roughly equivalent to a 200-year timber rotation. This is not unrealistic given that relatively few trees in the region's remaining remnant old-growth forests exceed 200 years in age (Parker, 1989). Our simulation of low intensity uneven-aged management is a process similar to what we anticipate for gap phase replacement. Harvest openings are small (0.09-0.27 ha or 1-3 pixel in size), scattered, and the oldest forests are affected first.
The simulation results also illustrate the effect of small but omnipresent disturbances of wind and fire. Over the short term, windthrow and fire disturbances affect a relatively insignificant portion of the landscape. Over the course of our 200-year simulation, however, fire and wind disturbances affected approximately half the landscape and had a notable impact on the forest size class distribution. Wind and fire disturbances add diversity to the landscape age structure, and they also reduce the area suitable for timber harvest. The net effect is a reduction in the harvested area or in the volume of harvested material due to these disturbances. Long-term management plans rarely anticipate the changes in structure and composition that inevitably result from wind and fire disturbances over time. Simulation modeling helps put this issue into context.
Effects of wind and fire disturbance on the age and species composition of affected pixels carry forward indefinitely from one decade to the next if there is no harvesting. Simulated harvesting of any type resets pixel age to zero and establishes new regeneration on harvested sites. When previously fire-or wind-damaged pixels are harvested, those pixels become indistinguishable from any other site harvested during the same decade. Thus, the historical legacy of the wind or fire disturbance is lost from those harvested pixels. Management Table 4 Standard error of the mean and coefficient of variation (in parentheses, expressed as a percent) at simulation year 200 of landscape characteristics calculated from five repeated simulation runs of the same forest management scenario alternatives with more intensive timber harvest tend to reduce the cumulative effects of wind and fire disturbance over time because they frequently regenerate forest patches. This was an unanticipated emergent result.

Harvest effects
A second unexpected outcome was the drop in harvest area and volume from decades 8-12 in the UAM 10% scenario and, to a lesser extent, for the mixed 10% scenario (Fig. 7). A similar drop in harvest volume began at decade 18 in the UAM 5% and mixed 5% scenario. This outcome resulted from the group selection harvest algorithm that simulated harvest of trees in groups ranging from 0.09 to 0.27 ha (1-3 pixel) in size. The harvest algorithm (Gustafson et al., 2000) actually tracks the location of group harvests over time and does not perform a second harvest in an earlier group opening until all locations in a stand have been harvested once or the rotation age is reached (i.e., after about 100 years for the UAM 10% scenario or the uneven-aged portion of the mixed 10% scenario, or after about 200 years for the UAM 5% or uneven-aged portion of the mixed 5% scenario). Near the end of the rotation period it became increasingly difficult for the harvest algorithm to identify previously unharvested locations of sufficient size to contain the new group openings. The issue was exacerbated because at each decade of the simulations the harvest algorithm continued until it either exactly met or slightly exceeded the desired harvest area (i.e., 10 or 5% of the stand area per decade). The accumulation of small overruns in harvest area further reduced the availability of unharvested sites at the end of the rotation period. This outcome is a direct result of how we chose to apply the LANDIS harvest algorithm, but it corresponds to a very real issue associated with the group selection silvicultural method. As areas undergo long periods of uneven-aged management, tracking locations of group openings in the field becomes difficult or impossible. Placement of new groups without overlapping prior openings becomes increasingly difficult over time and failure to do so effectively shortens the rotation length and/or reduces the harvest volume.

Landscape legacy
The simulation results show the enduring influence of the initial landscape conditions on the future landscape. On our landscape very little area was initially in the seedling size class, and over time this led to decreases in the area in the sapling size class and later the pole size class as the seedling cohorts matured and moved to a successively older age classes. Under regimes with a low disturbance rate (5% harvest or no harvest) this effect of the initial conditions lasted more than a century. For regimes with a greater disturbance frequency (10% harvest) this initial size class effect was smaller and shorter in duration. But in all cases the proportion of area in the sapling and pole size classes decreased over the first few decades and the proportion in seedling and sawlog size classes increased. Simulation modeling indicated how this legacy is perpetuated through time for alternative management practices.

Context specific data analysis and reanalysis
With the information tracked in LANDIS and recorded in the raster data structure it is possible to define and analyze patches, edges, and core areas in a variety of ways, and this is often necessary depending on the landscape characteristics of interest. For example, we measured core area as the area of patches greater than 60 m (2 pixel) from an edge, and a disturbance patch as small as 1 pixel (0.09 ha) created a surrounding edge. As a result, small but widely dispersed disturbance events had a large effect on core area and edge statistics, while intensive, even-aged management (with fewer, larger harvest disturbances) resulted in less edge and more core area than uneven-aged or no-harvest management. We considered this approach appropriate for assessing disturbance impacts on songbirds that respond positively or negatively to small disturbance patches (i.e., hooded warbler and ovenbird, respectively). Alternatively, we could have used an algorithm that required a larger disturbance patch (e.g., >1 ha) to delineate forest edge habitat, and the related landscape statistics would change. For example, a heterogeneous, uneven-aged mixture of many small age cohorts across a large area can function as a homogenous (uneven-aged) forest habitat for some wildlife species (e.g., white-tailed deer [Odocoileus virginianus]). Thus, in some contexts a tally of the edge habitat associated with thousands of small regeneration openings in a managed, uneven-aged, mature forest matrix may be irrelevant, and definitions of core areas and openings must be redefined accordingly. The LANDIS output provides the opportunity to redefine categories and reanalyze results via post-processing. For example, when we revised our analyses to compare the edge and core area for the combined seedling and sapling size classes (representing forest openings) versus the combined pole or sawlog size classes (representing closed forest cover), the length of edge under the UAM 10% scenario at year 200 decreased from 760 to 242 m/ha and the core area increased from 1129 to 4011 ha (Figs. 4 and 5).
The data-rich, spatially explicit modeling approach provides the opportunity to estimate values for a variety of attributes that can be linked to forest age and species composition. The volume of timber and down wood (Fig. 6) can be readily derived from information about forest age and analyzed spatially and temporally if desired. Other supplemental models can be quite complex and can integrate other sources of information. We were able to apply GIS-based habitat suitability models to simulate changes in wildlife habitat over time from the outputs produced by LANDIS and other GIS layers. Sullivan (2001) previously linked patterns of weather variability with LANDIS output to simulate the temporal and spatial distribution of hard mast production -a characteristic important to wildlife species and to the process of oak regeneration. Similarly Fan et al. (2003Fan et al. ( , 2004 were able to link a model of cavity tree abundance to LANDIS output. The ability to map and view simulation results (including values derived via post-processing) provides a useful way to spatially evaluate and communicate management implications that often get lost in tabular summaries.

Variation among repeated simulation runs
From a scientist's perspective there is a desire to replicate this simulation experiment on other landscapes and search for general trends that persist for many different landscapes. From a manager's perspective, however, a population of simulated outcomes based on multiple landscapes may be less desirable; managers are frequently most concerned with the outcomes for a specific landscape (including its unique set of initial conditions). In practice it may be reasonable to (a) look at the differences among repeated simulation runs on one landscape for characteristics of interest (e.g., for values reported in Figs. 2-7 and Table 3); (b) evaluate in practical terms which of those differences are likely to be relevant to management objectives (e.g., to the response of specific wildlife populations or humans); and (c) then determine if the variation (uncertainty) among multiple runs is large relative to the magnitude of differences presumed to be of practical importance. For many variables the variation among multiple runs for the same treatment was small (Tables 3 and 4). Treatment effects (e.g., Figs. 2-7) that differ by only a few percent are often not considered different from a practical standpoint, nor could such small differences be shown to be statistically different based on the variation reported in Table 4. In our simulations the harvested sites varied by location among repeated runs, but not in the disturbance patterns they created at the landscape scale. Simulated wind and fire events, which tend to be more variable in the patterns they create, disturbed only a small part of the landscape relative to timber harvest and age-dependent tree mortality. Moreover, given the size of our landscape (71,132 ha) differences among repeated runs in spatial patterns at specific locations are masked in the landscape means.

Habitat suitability
HSI estimates generally confirmed differences among the management alternatives that we would expect based on our knowledge of silvicultural systems and habitat relationships of the species we investigated (Thompson et al., 1995). For example, ovenbirds utilize a wide range of mid-to latesuccessional deciduous forest, but are edge sensitive. After 200 years of simulation mean HSI values for ovenbird varied by only a factor of 1.1 (0.54-0.60, Table 5); they were greatest in the landscape with the most late-successional forest (no harvest), and lowest in the landscapes with tree harvest and the most edge (uneven-aged management). Gray squirrels are dependent on mast and mature trees for cavities. They had the greatest mean HSI values in the landscape with no harvest followed by uneven-aged management, mixed management, and even-aged management. We suspect this pattern of HSI values is the net result of older trees being retained under noharvest and uneven-aged management regimes relative to even-aged management regimes and the greater estimated capacity for mast production associated with the older trees. Nevertheless, HSI values for gray squirrels after 200 years varied by only a factor of 1.2 (0.34-0.42, Table 5) among alternatives. Apparently all scenarios sustained enough attributes of older forests for both ovenbirds and gray squirrels so habitat suitability did not vary greatly.
Prairie warblers and hooded warblers are dependent on disturbance to create large or small patches, respectively, of early-successional habitat. Habitat suitability, on average, was much lower for these species than the mature forest species discussed above; most of the landscape provided no habitat (Table 5, Fig. 8). The prairie warbler and hooded warbler were also much more sensitive to the management alternatives with HSI values varying by a factor of 4.3-5.4, respectively, after 200 years of management (Table 5). The pattern of HSI values for prairie warblers (greatest in even-aged management, lowest in no-harvest) reflects their dependence on early-successional forest, preference for large patches, and avoidance of edge. In contrast, hooded warblers had the greatest mean HSI under the uneven-aged management scenario because they utilize small gaps within mature forest (Table 5, Fig. 8).
Interpretation of HSI values is currently limited to comparing differences in the index and how it varies between 0 (non habitat) and 1 (the best habitat); efforts are currently underway by the authors and other to relate HSI values directly to density or viability. There are multiple ways to summarize and report habitat quality at the landscape scale. In addition to means and medians, the distribution of HSI values (low to high) can be summarized over time or mapped to summarize them over space.

Modeling issues
A spatially explicit landscape framework for simulating forest composition and age structure is particularly powerful for integrating information. Software, such as LANDIS, automates many of the tedious record keeping and mapping functions needed to analyze results. Models can accommodate large landscapes (e.g., >75,000 ha) that encompass public agency planning units. The base data layers of forest vegetation by age class through time provide a mechanism for integrating other information and overlaying models and analyses for other resources. Examples include models for mast production (Sullivan, 2001), cavity trees (Fan et al., 2003(Fan et al., , 2004, and habitat suitability for a variety of wildlife species (Larson et al., 2003).
Landscape simulation models are not without limitations. For application in new regions the onus of model calibration falls on the user. The information required for LANDIS model calibration includes successional dynamics, wind and fire disturbance rates, and harvest regimes. These elements must be derived from external sources of information that are often limited in availability. Calibration of LANDIS to realistically simulate successional dynamics and response to disturbance for multiple species is a time-consuming process that requires protracted cycles of parameter setting, test runs, evaluation, and revision (e.g., Shifley et al., 1997Franklin et al., 2001). Once model calibration is complete and initial landscape conditions are mapped, comparing alternatives is relatively straightforward, but it requires tracking and analyzing large volumes of data describing future landscape conditions. Model application is limited to landscapes where information on land types and initial vegetation conditions are available or can be estimated. Currently such landscapes are relatively few in number. For some applications the requisite data layers can be developed via remote sensing (e.g., He et al., 1998He et al., , 2002Shao et al., 1996). Detailed maps with stand or management unit boundaries are often limited to public or corporate forest lands or areas where large investments have been made in forest mapping and inventory. Although we used a relatively simple (minimalist) algorithm to populate the initial tree species composition (one tree species per pixel), the simulations over time realistically increased the number of species per pixel in a manner consistent with observations at the adjacent Missouri Ozark Forest Ecosystem Project (Shifley and Brookshire, 2000). This dynamic based on minimal estimates of initial forest cover was encouraging because establishing the initial forest cover layer is generally difficult and data intensive under the best of circumstances.
Despite the difficulties typically associated the building initial data layers and calibrating a model, after a landscape model is initialized and calibrated it often opens new opportunities for collaboration, synthesis and development. Recent additions to the LANDIS fire modeling capability and ability to model biological disturbances (e.g., invasive insects, disease, oak decline, etc.) provide avenues build on past work and explore other dimensions of landscape change (Sturtevant et al., 2004;Yang, 2005;Yang et al., 2004). Incorporating effects of land use change is another obvious, albeit complex, way to build upon forest landscape simulation capabilities.

Conclusions
Spatially explicit landscape simulation models such as LANDIS are useful tools for exploring the potential effects of timber harvest, wind, and fire on future landscape conditions. Although they cannot predict the time and location of individual disturbance events, such models can describe the patterns that disturbances (natural or anthropogenic) are likely to create on a forest landscape. This modeling approach is well suited for examining emergent, spatially explicit properties relevant to evaluation of wildlife habitat, biodiversity, and the inevitable tradeoffs among multiple forest commodities, amenities, and services.
Landscape models are too coarse for site-specific planning, but they provide important, quantifiable indicators of the largescale, long-term consequences of management practices and natural disturbances. For example, management plans rarely account explicitly for natural disturbances that are likely to occur over time, but those disturbances are inescapable and over time they will affect a large portion of the total landscape.
The amount of data that can be produced by a LANDIS simulation run is both an asset and a liability. The output is voluminous -sometimes to the point of creating data collection, storage, and management issues. We used a 30 m pixel size, but pixel size can be scaled differently for other applications. Pixel sizes from 10 to 1000 m have been used elsewhere, and there are tradeoffs. For example, a 10 m pixel size corresponds roughly to the crown size of a mature tree in our study area. That resolution adds a degree of realism to the pattern of canopy openings when modeling single-tree selection harvests. However, the processing time and data storage requirements increase exponentially with decreasing pixel size. For large landscapes and long time horizons these are significant issues, especially for some of the HSI processing algorithms that perform data-intensive ''moving-window'' summaries for every point on the landscape. Despite these data processing issues, the spatially explicit approach provides the ability to overlay models for numerous forest attributes and to integrate other sources of information. In addition to its value for analysis, the ability to map simulated forest change over space and time is also important in displaying and communicating the outcomes of alternative management practices.
The value of this approach to forest management planning lies in the ability to analyze characteristics that depend on the spatial arrangement of forest vegetation. The ability to map and visualize forest characteristics through time is a great asset for communication and discussion. For many species of wildlife, assessment of habitat characteristics is highly dependent on knowledge of the spatial arrangement of forest types and age classes. Analysis of aesthetic considerations is dependent on spatial data as is analysis of landscape diversity. Many traditional aspects of timber management related to harvest quantity and harvest scheduling can be addressed without spatial data, but they often are much easier to visualize and communicate when summarized in map form. Other timber issues including transportation networks, riparian buffers, or adjacency constraints related to the pattern of harvest treatments are inherently spatial and are avenues ripe for exploration.
Several important and not entirely intuitive results emerged from the comparison of forest management alternatives. Alternatives with similar levels of disturbance had similar landscape composition but different landscape patterns. The no-harvest scenario resulted in a tree size class distribution that was similar to scenarios that harvested 5% of the landscape per decade; this suggests that gap phase replacement of senescent trees in combination with wind and fire disturbance may produce a disturbance regime similar to that associated with a 200-year timber rotation. Natural disturbance and mortality under the no-harvest scenario produced a landscape pattern that was most similar to the UAM 5% scenario. The greater harvest levels (10% per decade) regenerated more of the landscape each decade and obscured much of the impact of natural disturbances. The most intense management scenario, EAM 10%, created large uniform blocks of forest and resulted in the greatest core area and least edge, as we defined those habitat characteristics. Habitat suitability for late-successional wildlife species varied the least among alternatives, likely because all alternatives sustained enough attributes of mature forest in the landscape to accommodate these species. Habitat suitability for early-successional species was much more sensitive to the management alternatives, likely because there was less earlysuccessional forest than older forest, and these species were sensitive to disturbance patch size.