Resistance and Resilience of Hyrcanian Mixed Forests Under Natural and Anthropogenic Disturbances

Biological disturbances are integral to forest ecosystems and have pronounced effects on forest resistance, resilience, and diversity. The Hyrcanian mixed forest, in northern Iran, is at risk of declining resistance, resilience, and diversity due to ongoing pressure from land use change, harvesting, and biological disturbances. We analyzed the resistance and resilience of this area under two biological disturbances (i.e., oak charcoal fungus, Biscogniauxia mediterranea, and alder leaf beetle, Galerucella lineola) and in concert with proposed harvesting. We used a simulation modeling approach whereby we simulated 12 combinations of biological disturbances and harvesting scenarios using the LANDIS-II landscape change model. We estimated the correlation between forest resistance and resilience and tree species diversity to harvesting and biological disturbance. We analyzed the full species composition and age class for 30 and 100 years after disturbances in order to assess resistance as the change in species composition over time. We considered resilience as the ability to recover from a disturbance back to a similar initial state. Results indicate a positive effect of biological disturbances and harvesting on diversity. Our simulations resulted in a negative relationship between diversity-resistance and diversity-resilience within high diversity areas. Our simulation of the Hyrcanian forest reveals that harvesting and biological disturbances, as tested, fulfill the goal of maintaining forest diversity. However, increasing diversity does not always follow by increasing forest resistance and resilience. Scenarios with oak charcoal fungus, both with and without harvesting indicate the lowest decrease in resilient and resistant.

Biological disturbances are integral to forest ecosystems and have pronounced effects on forest resistance, resilience, and diversity. The Hyrcanian mixed forest, in northern Iran, is at risk of declining resistance, resilience, and diversity due to ongoing pressure from land use change, harvesting, and biological disturbances. We analyzed the resistance and resilience of this area under two biological disturbances (i.e., oak charcoal fungus, Biscogniauxia mediterranea, and alder leaf beetle, Galerucella lineola) and in concert with proposed harvesting. We used a simulation modeling approach whereby we simulated 12 combinations of biological disturbances and harvesting scenarios using the LANDIS-II landscape change model. We estimated the correlation between forest resistance and resilience and tree species diversity to harvesting and biological disturbance. We analyzed the full species composition and age class for 30 and 100 years after disturbances in order to assess resistance as the change in species composition over time. We considered resilience as the ability to recover from a disturbance back to a similar initial state. Results indicate a positive effect of biological disturbances and harvesting on diversity. Our simulations resulted in a negative relationship between diversity-resistance and diversity-resilience within high diversity areas. Our simulation of the Hyrcanian forest reveals that harvesting and biological disturbances, as tested, fulfill the goal of maintaining forest diversity. However, increasing diversity does not always follow by increasing forest resistance and resilience. Scenarios with oak charcoal fungus, both with and without harvesting indicate the lowest decrease in resilient and resistant.

INTRODUCTION
Forest landscapes are subject to a broad range of disturbances. Biological disturbances, those caused by insect and disease outbreaks, are among the principle agents of forest change and cause tree mortality at scales ranging from individual trees to entire regions (Wood and Van Sickle, 1994;Unger and Stewart, 1995;MacLean, 2008;Sturtevant et al., 2004;Trumbore et al., 2015).
Insects and pathogens play an important role in regulating forest ecosystem functions, and there is increasing evidence that these organisms are one of the primary determinants of forest structure and composition (Lewis and Lindgren, 2000).
Biological disturbances are one of the most prevalent disturbances in the Hyrcanian forests bordering the Caspian Sea. These forests stretch like a green belt over the northern slopes of the Alborz mountain range and along the southern borders of the Caspian Sea. Their range extends across three provinces of Iran (i.e., Gilan, Mazandaran, and Golestan) from Astara in the northwest of Iran to the vicinity of Gorgan (Golidaghi) in the northeast of Iran, and also include a small western portion in the country of Azerbaijan ( Figure 1A). This covers an area of 1.85 million ha, and contains ∼15% of the total Iranian forests and 1.1% of the country's area (Sagheb-Talebi et al., 2004). The Hyrcanian forests contain a unique assemblage of endemic tree species (Tohidifar et al., 2016), yet are threatened by land-use change and mortality from biological disturbances. Several pathogen and insects species have caused extensive tree mortality in Hyrcanian forests (Hajizadeh et al., 2016), including oak charcoal disease [Biscogniauxia mediterranea (De Not.) Kuntze], and brown alder leaf beetle [Galerucella lineola Fabricius] (Sadeghi et al., 2004;Karami et al., 2018). Tohidifar et al. (2016) noted that since 1950, the total area of Hyrcanian forests has declined from 2,750,000 ha to 1,850,000 ha (32.7%) due to a wide range of perturbation (anthropogenic and biological disturbances).
Harvesting also occurs within the Hyrcanian forests through land use change and logging. Homogenous harvesting (e.g., clear-cut) can result in greater landscape scale susceptibility to other disturbances, reducing both resistance and resilience (Drever et al., 2006). Harvesting alters the disturbance legacy of a region and has the potential to remove the ecological memory (such as soil build-up, seedbanks, or heterogeneity) reducing its resilience to subsequent disturbance (Johnstone et al., 2016). Other silvicultural treatments, however, may increase the resilience of forests by increasing heterogeneity, reducing the drivers of initiation or spread of a disturbance, or increasing the likelihood of regeneration of a given species (DeRose and Long, 2014). It is therefore crucial to understand whether the Hyrcanian forests will remain ecologically resilient under this disturbance regime (the long term pattern of a disturbance acting on particular location), or whether tree species composition of the area will be permanently altered in the coming century.
Resilience is a key component of forest ecosystems (Yan et al., 2011). While the concept of resilience is ubiquitous in ecology, researchers are still uncertain on the definitive characterization or demarcation of ecosystem resilience. As mentioned by Klein et al. (2003), there are a vast number of definitions of resilience with various characteristics. Engineering resilience is defined as the time required for a system to return to predisturbance conditions, following a disturbance event (Holling, 1973;Pimm, 1984). This is how we will define "resilience" throughout this paper.
Some authors have regarded resistance as an integral part of resilience (Walker et al., 2004), while Grimm and Wissel (1997) have considered it as a separate property, which measures the degree to which a system remains unchanged despite perturbation, likewise, we use this definition of resistance here. Newton and Cantarello (2015) noted that if the maintenance of a forest ecosystem in its current state is a priority, despite a changing disturbance regime, then resistance could provide an appropriate objective. Resistance can also be regarded as synonymous with the stability of a landscape. Higher diversity may lead to greater stability and resistance (Elton, 1958;Tilman and Downing, 1994;McCann, 2000;Tilman et al., 2006;Moore et al., 2011). It is also important to understand the resistance and resilience of the Hyrcanian forests, given the magnitude of the damage of alder beetle and oak fungus outbreaks, and the economic significance of the host trees (Sadeghi et al., 2004;Karami et al., 2018).
Ecological resilience is hard to quantify (Grimm and Calabrese, 2011), and consequently it is still controversial how to measure it directly (Webb, 2007). Holling (1973) and Gunderson (2000) defined ecological resilience as the magnitude of the maximum disturbance the ecosystem can absorb before it changes into another state. Others have estimated the time for the system to recover back to the stable state following disturbance (Grimm and Calabrese, 2011) or the maximum stress the ecosystem can absorb without changing (Mitchell et al., 2000). In contrast, Martin et al. (2011) used engineering resilience, the time required following disturbance to come back to its initial state. Such analyses can be carried out with simulation models (Pimm and Lawton, 1977;Duveneck and Scheller, 2016).
It is important to quantify the resistance of the Hyrcanian forest to disturbances as it is a refuge for many quaternary species and considered as a national treasure mainly due to its high species diversity (Shakeri et al., 2012;Sagheb-Talebi et al., 2014). These Hyrcanian forests are considered as one of the oldest and most threatened temperate forests globally. These forests are a mixture of uneven-aged deciduous broadleaf tree species (Sagheb-Talebi et al., 2014). The mixed Hyrcanian forests rise from 20 m to an altitude of 2,800 m, and encompass a variety of different forest types. Fagus orientalis, Quercus castaneifolia, Carpinus betulus, Parrotia persica, Acer velutinum, Alnus subcordata, and Alnus glutinosa are among the most important tree species and can form either pure stands (i.e., if dominant tree covers more than 80% of the area) or mixed stands (if none of the dominant species can cover more than 80% of the stand). The management objective of Iranian forest organization is to preserve the current composition of forest stands and avoid any shift in the composition of commercial forest types (e.g., pure Fagus orientalis, Quercus castaneifolia, and other mixed types). Therefore, we explored how biological disturbances and harvesting may alter the long-term resistance and resilience of Hyrcanian forests via forest change simulations. Our main objective was to determine whether Hyrcanian forests will be resilient and resistant following insect and fungal outbreaks under ongoing harvesting operations.

Study Location
Our study area covered a section of the Hyrcanian mixed forests along the southern coast of the Caspian sea and northern slopes of the Alborz mountains ( Figure 1B). We selected this study area due to its detailed records of forest management practices and frequent reports by local experts about the alder leaf beetle and oak charcoal outbreak. Our focal watershed has an area of 2,567 ha ( Figure 1C). Elevation ranges from 220 to 700 m. The slope varies between 0 and 30%. The climate is temperate and humid with an average annual precipitation of 1188.6 mm and mean annual temperature of 15.2 • C.

Forest Types Classification
The dominant tree species are common hornbeam (Carpinus betulus), caucasian alder (Alnus subcordata), chestnut-leaved oak (Quercus castaneifolia), date-plum (Diospyros lotus) and persian ironwood (Parrotia persica). Caucasian alder constitutes over 7% the land of the northern forests of Iran, and is considered as the fourth most commercially important tree of the Hyrcanian forests after the oriental beech (Fagus orientalis), common hornbeam and chestnut-leaved oak trees (Sadeghi et al., 2004). Caucasian alder is a native species to the Hyrcanian forest and has an important ecologic and economic role in the region (Sagheb-Talebi et al., 2004). To classify forests we followed the forest type identification used by the Iranian forest service [modified from Küchler and Zonneveld (2012)]. If a dominant tree species makes up more than 80% of the stand composition, the stand is classified as a pure type, otherwise if the second species makes up more than 20%, the classification will be composed of the first and second dominant tree names. In case the stand contains more than three codominant tree species, it will be called mixed broadleaf type. Mixed broadleaf is the second main forest type after Fagus orientalis type in the Hyrcanian forest. It is both the most diverse group and contains the most economically valuable tree species such as chestnut-leaved oak, Velvet maple (Acer velutinum), Cappadocian Maple (Acer cappadocicum) and Caucasian alder. Other major remaining forest types include Carpinus, Carpinus-Diospyros, and Carpinus-Parrotia and Alnus-Carpinus (Forests, Rangelands and Watershed Management Organization).

Dominant Biological Agents
Oak charcoal fungus can survive in aerial organs of chestnutleaved oak and acts as a pathogen when chestnut-leaved oak is under stress, particularly during drought (Capretti and Battisti, 2007). Oak Charcoal disease results in cambium necrosis which can initiate a decline in tree health leading to the death of the host trees. The first symptoms of charcoal disease are the discolorations and browning of the leaves, the drying of foliage, viscous liquid exudates and lengthwise bark cracks. Subsequently, the outer bark begins to slough off the area of infection, this can be identified by recognizing pieces of bark at the base of the tree (Ragazzi et al., 2012). Cankers are usually localized in the bark of stems or twigs of trees and can rapidly enlarge lengthwise in humid weather conditions. The pathogen results in cambium necrosis and kills portions of the sapwood; hence increasing potential for stem breakage. If cankers girdle the stem, the branches will die. These weakened trees are predisposed to secondary infection through branch stubs and other wounds (Inácio et al., 2011). The infection and mortality rate caused by charcoal fungus were 69 and 3% per year, respectively, in the Hyrcanian forests (these values related to the study that was conducted in the east part of Hyrcanian forests of Iran) (Karami et al., 2018). Oak charcoal fungus may cause irreversible damages to the wood quality and subsequently to the forest (Inácio et al., 2011;Henriques et al., 2014).
Brown alder leaf beetles (Galerucella lineola Fabricius) are herbaceous defoliators that use species of Salicaceae family for feeding and oviposition (Hambäck et al., 2013). Both larvae and adults feed on the plants emerging buds leading sometimes to quite severe and extensive damage (Escherich, 1923;Kovačevič, 1957;Hunter, 1992;Hoeglund et al., 1999;Sage et al., 1999). High activity of the alder leaf beetles is known on alder trees (Alnus subcordata) in Iran (Sadeghi et al., 2004). Separately Gharadjedaghi (1997) studied the occurrence and dynamics of the abundance of G. lineola in the vicinity of Bayreuth (northern Bavaria). According to the author there was a heavy outbreak of the chrysomelid on alder trees.
The larger understanding of the brown leaf alder beetle life cycle is as follows: after over-wintering, adults emerge in April and start their feeding period, followed by copulation and oviposition (egg-laying). The larval stages are generally found from May-July. Then the new generation of adult beetles appeared in July-August and feeds before hibernating (Kendall and Wiltshire, 1998). However, in Iran imagoes activate already at the end of March (Sadeghi et al., 2004). The damage is mostly related to adult wintering insects because these insects are veracious feeders and consume a large amount of plant material while growing. When adults feed it leads to irregular holes in the leaf surface destroying the leaves and reducing photosynthesis and weaking the tree, therefore it increase the mortality rate (Sadeghi et al., 2004). Caucasian alder trees are among the most dominant species in Hyrcanian forests, and the incidence of the beetle could considerably alter forest composition (Sadeghi et al., 2004).

Simulation Model Framework and Parameterization
We quantified the effects of biological disturbances and harvesting on resilience and resistance using LANDIS-II (v6.1) a forested landscape succession and disturbance model (Scheller et al., 2007). LANDIS-II is a spatially explicit and stochastic model that simulates how ecological processes (including growth, succession, dispersal, and disturbances) affect a forested landscape over time. LANDIS-II requires unique life history parameters for each tree species, such as shade tolerance and dispersal distance ( Table 1). Tree species age classes are mapped onto a grid of sites representing the landscape. The sites were initialized to represent the contemporary tree species and age classes present (Scheller et al., 2007). For each species present on a forest site, individual trees classified into age classes. Here, we refer to age classes by the upper bound of its range. For example, a tree cohort aged 11 to 20 would be referred to as age class 20 (Scheller et al., 2007). Our age classes were spaced by 10 years.
LANDIS-II requires an initial community layer (tree species that are present on each site categorized by age class) and ecoregions (areas of similar abiotic characteristics). We developed the initial community layer from our field work that identified the species of each tree and estimated age by measuring tree diameter (see Supplementary Material, section 1 for sampling design). To prepare the initial community layer, we field inventoried the Hyrcanian forests using a randomsystematic sampling with a circular plot. The inventory included 667 systematic random (SR) sampling plots of size 1,000 m 2 , designed within a 200 × 200 m grid (the cells resolution in the software was reduced; cell length: 5 × 5 m = 0.001 ha). Existing stands were used to map the initial communities onto the landscape, however, each cell simulates succession independently when the model run begins.
Landscape ecoregions and species life history were defined using forest and natural resource organization data and forest expert consultation, respectively. We divided our landscape into four ecoregions as identified by soil type.
We used Age-only Succession extension v4.1 of LANDIS II to simulate forest succession (Mladenoff and He, 1999). The extension simulates the establishment, aging and mortality of each cohort on each site. To calculate establishment, it uses algorithms that consider species-specific life history attributes (e.g., longevity, shade tolerance, and age) ( Table 1). The Ageonly Succession Extension also simulates tree mortality caused by senescence and age.

Description and Parameterization of Biological Disturbances Extension
We used the Biological Disturbance Agent extension (BDA) v3.0 (Sturtevant et al., 2004) which simulates tree mortality following major outbreaks of insects and disease at every 5-year time step.  Biological disturbances in LANDIS II are probabilistic at the site (i.e., cell) scale, where each site is assigned a probability value called the biological disturbance probability (BDP). After the BDP is calculated, it is compared to a uniform random number to determine whether the site is disturbed or not. Disturbance causes species-and age cohort-specific mortality in the cell. In the simplest case, BDP is equal to the Site Resource Dominance. This number is calculated on the tree species and age cohorts present and ranges from 0.0 (no host) to 1.0 (most preferred host). In total, five main elements control the probability of biological disturbance within the extension: (1) local resource dominance; (2) host value modifiers that reflect environmental conditions and recent disturbance history; (3) host dominance within a user-specified neighborhood; (4) the temporal outbreak pattern characteristic of the BDA (as parameterized by the user); and (5) BDA dispersal. Additionally, the BDA allows us to simulate multiple biological disturbances concurrently. In this extension, species parameters for each biological agent comprise the host susceptibility classes (1, 2, 3) which are ranked in terms of the probability of sustaining an outbreak (class 1 being the least susceptible and class 3 the most). Each host species has susceptibility class ages indicate the minimum age at which a species enters a respective susceptibility class. These classes determine the age cohorts are subject to mortality if a site is disturbed. Cohorts younger than the minimum age for youngest susceptibility class are immune to the BDA. This helps differentiate young cohorts that are not yet susceptible to certain biological disturbance (Sturtevant et al., 2015).

Description and Parameterization of Harvest Extension
We also utilized Base Harvest Disturbance Extension v3 (Gustafson et al., 2000) at every 5-year time step to implement harvest prescriptions. The harvest simulations were carried out within management areas. Each of the management area consist of a unique collection of stands. Each of the management areas and stands also requires its own input map . To simulate harvesting, we set a 5-year time step. Three management areas were specified based on Forests, Rangelands and Watershed Management Organization classification which is formed from the quality, quantity and type of tree species in various parts of the forest. The management areas were classified as productive forest (Area 1), degraded forest (Area 2) and reforestation (Area 3) (Figure 2). Within each management area, we designated stands based on the type of tree species and their dominance on the forest landscape resulting in 193 stands.
Each management area had one harvesting prescriptions. Hence, the three prescriptions were parameterized according to the policy of the Forests, Rangelands and Watershed Management Organization (FRWO). For the Hyrcanian mixed forests the management goal is to conserve forests as a diverse and mixed uneven-aged structure. According to forest expert's views, selecting an appropriate cutting method with respect to characteristics of each management areas would be a proper tool for implementing the policy. The productive forest management area prescription was based on the single-tree selection system wherein harvest is based on the longevity of each species. Therefore, we ranked these stands by weighting stand selection based the number of the trees available for each species per hectare; the most trees per hectare were ranked highest. Cohorts within the stand were harvested based on their longevity. This harvesting prescription increases tree regenerating, thus helping to restore the forest to its natural state (mixed and uneven-aged forest). The reforestation management area was parameterized according to the forest organization plan for the thinning cuts. Stands were ranked to target shade intolerant and pioneer tree species planted in the reforestation areas (e.g., Caucasian alder and Velvet maple). Cohorts were removed based on the typical thinning regimes (ages from 15 to 30). This prescription allowed reproduction of other tree species to turn these areas into mixed and diverse forests. The forest organization's plan for the degraded management area was a clear-cutting prescription; all trees were removed from stands that were randomly selected. We simulated planting of shade intolerant tree species including Caucasian alder, Velvet maple and chestnut-leaved oak. In this harvest prescription, planting shade intolerant and pioneer tree species resulted in a more diverse forest. We fixed harvesting rate at 5% of each management area at each time step.

Experimental Design
Our experimental design contained 12 scenarios by which we examined whether each stand was resistant or resilient (or neither) following biological disturbances and harvesting operation ( Table 2). We tested oak charcoal fungus by itself, alder leaf beetle by itself, and both together. Because there is limited data on the ultimate mortality caused by LAB, we designed two scenarios for alder leaf beetle disturbance regime; one where the ultimate morality is slight (equilibrium state; only older cohorts would be killed), and the other where the alder leaf beetle results in a high level of mortality (non-equilibrium state; all age-class cohorts could be killed). We ran each of the two biological disturbance scenarios separately and combined, both with the presence and absence of harvest activity to evaluate the effect of harvesting on landscape resilience and resistance. All scenarios were simulated for 100 years. Our simulations had a 10 m 2 (0.001-hectare) resolution and our landscape was 2,567 ha. Each scenario was replicated three times. We analyzed the patterns with a classification scheme based on forest types as defined by FRWO of Iran.

Measuring Resilience and Resistance
To measure the resistance and resilience of our study area, we used a method similar to Duveneck and Scheller (2016).
We assessed simulated forest resistance by measuring the dissimilarity of initial species composition and age class under effects of biological disturbances outbreak over the course of 100 years. Following their method, we used the Hyrcanian forest manager point of views to define a proper resistance criterion for the Hyrcanian forest. If the forest type classification was the same and the maximum age did not decrease by more than 5%, the stand was considered resistant, otherwise it was not. To measure resilience, we evaluated the recovery of species composition from the initial point over a 30-year period. Each cell was checked against its future self every 10 years that followed. For this, we looked at each cell deemed not to be resistant then we looked at whether it recovers to its initial state over the next 30 years (or at all). Those that recovered to the original state any time in the next 30 years are classified resilient.

Species Composition
Species composition changed between the years 0 and 100 dependent on disturbance and harvest scenario (Figures 3, 4).
In the oak charcoal fungus scenario without harvest (OCF-NH), pure Carpinus type remained dominant and relatively unchanged by the year 80, and then it decreased rapidly and was replaced primarily by Alnus-Carpinus type in the year 100 ( Figure 3D). Mixed broadleaf type was the second most dominant at the end of simulation. Whereas, with forest management (OCF-H), Carpinus declined gradually by the year 50 and then plunged sharply. Carpinus was replaced by mixed broadleaf and Alnus-Carpinus ( Figure 3C). The ALB-NH (alder leaf beetle/no harvest, Figure 3H) and BD-NH (both disturbance/no harvest, Figure 3F) scenarios showed that in the absence of harvest management, the number of pure Carpinus type stands experienced a relatively steep decline from the start of the period and was nearly eliminated by the year 100. There was a corresponding rise in the number of Carpinus-Diospyros stands which exceeded that of pure Carpinus stands around the year 50. Carpinus-Diospyros stands peaked around the year 60, climbing to between 300 ha in ALB-NH and BD-NH scenarios. Following this, Carpinus-Diospyros declined and was supplanted with Carpinus-Parrotia and Alnus-Carpinus stands. Carpinus-Parrotia became dominant in ALB-NH and BD-NH scenarios (Figures 3H,F). By contrast, when the alder leaf beetle was parameterized in the equilibrium state under the same condition as ALBE-NH and BDE-NH scenarios, there was not a sharp and constant reduction in the amount of pure Carpinus stands (Figures 4A,C). However, these stands do decline by the year 80 and then pure Carpinus sites undergo a steep decline. Additionally, there was negligible growth in the number of other stand types especially in Carpinus-Diospyros and Carpinus-Parrotia stands unlike the non-equilibrium state. Harvesting accelerated the decline of pure Carpinus stand and slowed the increase of Carpinus-Diospyros in all scenarios (Figures 3A,C,E,G). Harvesting increased other species dominance, especially mixed broadleaf species. Harvesting generated similar behavior in the alder leaf beetle equilibrium/non-equilibrium scenarios (Figures 4B,D).

Forest Resistance
We analyzed each scenario to assess how the number of resistant cells changed through time ( Table 3). After 100 years, we forecast that the number of resistant sites will decline around 46% for the no disturbance (ND) scenario as baseline, suggesting that a majority of the forest will not be resistant under any of the scenarios examined. For the scenarios without harvesting, the greatest decline in resistant sites after 30 years was under BD-NH scenario. Compared to a no disturbances (ND) scenario, BD-NH (oak charcoal fungus plus alder leaf beetle without harvest) had 30% fewer resistant sites. The alder leaf beetle in non-equilibrium without harvest operation (ALB-NH) was about 19% less resistance than the ND scenario. The scenarios of ALBE-NH and BDE-NH are in the third and fourth place, respectively, with about 7 and 3% less resistant area ( Table 3). The least reduction in the number of resistant sites was OCF-NH scenario. In the year 100, the order of the decreasing effect of the scenarios on the number of resistant sites remained similar to year 30 (Table 3).
These patterns changed when we included harvesting ( Table 4). After 30 years, the greatest difference in the area of resistant sites was between the ALB-H and the scenario of no disturbances (ND), a 37% decline in resistant sites. This was followed by BD-H with a 37% decrease in resistant sites compared to the control scenario ( Table 4). Other important declines due to harvesting came from BDE-H and ALBE-H scenarios with 16 and 16% drop in resistant sites, respectively. The smallest reduction was under the OCF-H scenario with a 13% decline in resistant sites. After 100 years, this trend changed so that in a descending order from the largest decline to the lowest is as follows: scenario ALB-H with a 70% decrease, BD-H scenario with a 69% decrease, BDE-H scenario with a 67% decrease, ALBE-H and OCF-H scenarios with 65 and 57% reduction in the amount of resistant sites ( Table 4).

Forest Resilience
We forecast the number of resilient sites (measured as recovery from initial point after 30 years) and its changes under different biological disturbance. Throughout the simulated period the forest was resilient (>50% cells recovered) under the no disturbance (ND) and OCF-NH (oak charcoal fungus/NH) scenarios ( Table 5).
The simulated forest landscape was not resilient (<50% cells recovered) under any of the scenarios in the presence of harvest ( Table 6).

Effect of Biological Disturbances on Species Dissimilarity
The species dissimilarity remained unchanged under outbreaks of oak charcoal fungus (OCF; Figure 3D) from the control scenario of no disturbances (ND; Figure 3B). This suggests that forest demography along may lead to these successional shifts over time regardless of disturbance.
Unlike our hypothesis (decreasing trend of dissimilarity under biological disturbances in the future), the results also show that the most and second most dissimilarity will happen under outbreak of alder leaf beetle (ALB-NH; Figure 3H) and both of biological agent together (BD-NH; alder leaf beetle plus oak charcoal fungus, Figure 3F), respectively, in comparison of no disturbances scenario simulation (ND; Figure 3B). As the area of pure Carpinus forest type will decline, the area of Carpinus species combined with other tree species will increase simultaneously (Carpinus-Diospyros and Carpinus-Parrotia).  Table 2 for scenarios abbreviations. OCF/H-oak charcoal fungus with harvest operation, OCF/NH-oak charcoal fungus without harvest operation, ALB/H-alder leaf beetle with harvest operation, ALB/NH-alder leaf beetle without harvest operation.  Table 2 for scenarios abbreviations.

Biological Disturbances
The simulated forest landscape became not-resilient (<50% resilient) after 30 years except in the OCF-NH (oak charcoal fungus/no harvest) disturbance regime ( Table 5). Over the same time period the forest was resistant following all biological outbreak scenarios, despite the decreasing trend in the number of resistant stands, after 100 years the cumulative decrease in resistant stands resulted in a non-resistant landscape under all scenarios (Table 3). Pathogenic organisms and herbivorous insects are integral components of forest ecosystems, altering the diversity of those ecosystems in many ways (Goheen and Hansen, 1993;Hansen et al., 1997). It is the subject of debate whether or not biological disturbances increase diversity; it has been Not resistant (NR), Resistant (R), 95% CI in parentheses See Table 2 for scenarios abbreviations.  proposed that endemic pathogens increase the diversity of forests (van der Kamp, 1991;Goheen and Hansen, 1993;Thom and Seidl, 2016), On the contrary, exotic pathogens may destroy diversity (Hansen, 1999). Our results using an endemic biological disturbance suggest that they generate higher diversity through increasing tree species dissimilarity by decreasing the dominate forest type (pure Carpinus), thus leading to the increase the amount of other forest types (Carpinus-Diospyros and Carpinus-Parrotia).
The biological disturbances on the landscape increased the diversity in the forest via increasing the dissimilarity of tree species. The disturbances with the greatest mortality (BD-NH and ALB-NH scenarios; Figures 3F,H, respectively) drastically increased tree species dissimilarity by causing a significant decrease in the single dominant forest type and therefore generated higher total forest diversity. During severe outbreaks of the alder leaf beetle, the dominate forest type (Carpinus) experienced a noticeable reduction, leading to an increase of the less dominant (Carpinus-Diospyros and Carpinus-Parrotia) forest types (Figures 3F,H).
The less severe disturbance scenarios (BDE-NH and ALBE-NH scenarios) showed a similar trend but to a lesser affect (Figures 4A,C).
The response of the forest to the alder leaf beetle may depend on the intensity of the outbreaks. Because this is not known, we ran two alder leaf beetle disturbance regimes. In the first only older cohorts would be killed (equilibrium), and in the second all age-class cohorts could be killed (non-equilibrium). In the equilibrium state (alder leaf beetle equilibrium/No harvest; ALBE-NH; Figure 4A), forests showed less dissimilarity than the same as the scenarios in non-equilibrium state (ALB-NH, Figure 3H), but the reproduction of Caucasian alder was greater in the non-equilibrium state and the reproduction of young and middle-aged cohort of Caucasian alder tree allowed the species to survive (Figure 4A).
The difference between these two disturbance regimes is accounted for in the uncertainty in the alder leaf beetle mortality pattern. During the strongest outbreak (class 3), the scenarios with more severe mortality (ALB-NH, BD-NH) considered that Caucasian alder trees in all age classes (healthy and susceptible trees) would be killed ( Table 3). In the equilibrium model (ALBE-NH, BDE-NH), only age classes greater than 90 would be killed (only susceptible and weakened trees) (Supplementary Table 3). This in turn leads to a period between the years of 40 to 100 in which common hornbeam mixed forest becomes the dominant species. This same pattern can be seen in the harvest only scenario with the window for alder domination being shorter. This difference in susceptible ages allows changes the establishment of new species, which then contribute to the shift in forest type.
The selection pressure, and the timing of the decline in Carpinus trees, however, did not affect the final forest composition as Carpinus and Carpinus-Diospyros declines and the other forest types establish more readily (Figures 3, 4). Given that alder leaf beetle (ALB-NH) scenario was parameterized by local forest expert assumptions, entomologist knowledge and observational accounts, future scientific research are necessary to determine correlations of age, host preference, and susceptibility to improve these scenario parameterizations. This difference in heuristic understanding seems to determine a large amount of the forest composition difference.
The oak charcoal fungus did not have a substantial effect on forest type, likely because oak is only a component of the mixed broadleaf forest type. It additionally had a minimum impact on the resilience or resistance on the landscape, likely because it affected only a minor component.

Harvesting
Although harvesting created more diversity by increasing tree species dissimilarity via raising mixed broad leaf forest type (as the most diverse forest type) in all scenarios (Figures 3A,C,E,G,  4B,D), it did not increase the number of resistant and resilient sites over time and rather decreased them (Tables 4, 6). Harvesting did, however, increase the mixed broadleaf forest type, the most diverse group, and other composition types (Figures 3A,C,E,G, 4B,D). The impact of harvesting and total diversity was location and forest structure specific. Peltzer et al. (2000) found a positive correlation of harvesting and diversity while others emphasize negative association of diversity and harvesting (Brown and Gurevitch, 2004;Jafari et al., 2013). Given the intermediate disturbance hypothesis, mild disturbance would increase species richness while species richness is markedly decreased by the intense and low disturbances (Mackey and Currie, 2001;Mishra et al., 2004).
We had anticipated that harvest would increase the forest resilience (measured as recovery from initial point) by increasing the total landscape diversity over time, and the forest would therefore be more resilient under harvesting scenario. However, OCF-NH is the only scenario in which the forest was resilient (Table 5), the forest became not-resilient under the same scenario in the presence of harvesting (i.e., OCF-H) ( Table 6). This is likely because the overall trajectory of the landscape moving from pure Carpinus to a more diverse landscape type. The harvesting accelerates this transition. Perhaps over longer time frame this increase in diversity may increase the resistance and resilience.
This type of harvesting simulated has been performed for many years by FRWO to keep Hyrcanian forest as a mixed forest and promote diverse forest types. Our research supports the organization's goal at this scale. When biological disturbances and harvesting occur simultaneously (BD-H; Figure 3E and ALB-H; Figure 3G), the diversity peaked. Although biological disturbances and harvesting maintained diversity, harvest accelerated the shift in forest composition, lowering both the resistance and the resilience of the landscape over the next 100 years. The biological disturbance regime (scenario differences), however, had far more influence on the resilience numbers than the harvest differences. This suggests that while harvesting may reduce forest resilience in the short term, this difference is small in comparison to changes by other disturbances. Given that resistance scores over the next century are low under all likely scenarios, the relative difference in resilience through the addition of harvesting seems slight.

Diversity and Resistance/Resilience
The relationship between diversity and ecosystem stability has long been debated (Pimm and Pimm, 1991). This correlation was considered as positive until the early 1970s, assuming that more diverse ecosystems result in more energy and nutrient flows (Bengtsson et al., 2000). The key importance of diversity may be the capacity of a variety of species to continue ecosystem functions even if some of the present species should disappear (Folke et al., 1996). However, we found that higher diversity was more likely to display lower resistance and resilience (Figures 3F,H and Tables 3, 5; ALB-NH, BD-NH). Similarly, May (2019) argued that diverse communities were less stable than simple ones. On the contrary, other studies (Tilman and Downing, 1994;Naeem, 1998;Peterson et al., 1998;Bengtsson et al., 2000) noted that more diverse ecosystems are more resistant and resilient to disturbances. (Duveneck and Scheller, 2016) reported a positive relationship between diversity and resistance within low diversity areas, whereas we found a negative relationship between diversity, resistance and resilience within high diversity areas (Figures 3F,H and Tables 3, 5; ALB-NH, BD-NH). In many cases, tree diversity is also expected to provide resistance to natural disturbances including insect herbivores so that mixed forests are more resistant than monocultures (Folke et al., 1996;Jactel and Brockerhoff, 2007;Griess et al., 2012;Castagneyrol et al., 2014;Jactel et al., 2017). Our result contradict these trends to biological disturbances in Hyrcanian mixed forest.

STUDY LIMITATIONS
It is important to consider that our simulations looked only at dominant forest type and not total forest biomass. Given the amount of mortality experienced under the alder leaf beetle disturbances, total biomass would have declined over the next century. There are also other drivers shaping this landscape that were not considered in the study including cattle browsing of young trees, fire, and windthrow In addition, future research is necessary to determine correlations of age, host preference, and susceptibility of alder leaf beetle parameterization to improve the accuracy of the disturbance forecasts. There are also a vast number of biological disturbances and each of them can play an important role in resistance and resilience of the forest and its diversity. Our harvest management regimes were based on what has been conducted in the Hyrcanian forest. Harvesting is currently banned in the Hyrcanian forest for the next 10 years, and new decisions may be made in the future that would change these forecasts. Some of the factors that we did not include in our management strategies (e.g., planting extinct species) represent future opportunities for research. Our result suggests that the resistance and resilience of the study area will decrease as the tree diversity will increase by increasing dissimilarity of tree species forest type (dissimilarity inside the 2,500 ha study area). From this perspective, it is expected that the forest would be less resistant and resilient at the finer spatial scales due to more tree species dissimilarity (that is, in the finer scales, alder leaf beetle hosts (Alnus Sp.) may be as single dominate forest type. Therefore, it is possible that with their widespread mortality, there will be enough space for other species to regenerate and make the forest more diverse. But, at a scale larger than our study area (i.e., for the whole range of Hyrcanian forests), forests may be less vulnerable against the leaf beetle due to lower densities of the host tree species (alder trees are found in low-land areas), meaning that the defoliator insect would not have a significant effect on the tree species composition. Hence, at larger scales, tree species composition and forest resistance and resilience might not be substantially affected by alder leaf beetles. Our findings on resistance and resilience are based on age and dominant forest type classification, two common though coarse understandings of forest composition. Accessing the resilience of forests is inherently a value judgment as to which ecosystem structures, functions or services are measured by different metrics (Millar et al., 2007). Including additional metrics such as the resilience or resistance of biomass or carbon captured were beyond the scope of this study, given the available data. Future studies of the region may allow for a more wholistic understanding of the forest's resilience with regards to the ecosystem services this forest would provide a broader scope for managers to understand these trade-offs.

CONCLUSION
Our simulations of the Hyrcanian forest reveals that management activity and disturbances, as tested, fulfill the goal of maintaining forest diversity, however, current diversity is not correlated with increased resilience and resistance of species composition. Our approach provides an effective method to quantify forest resilience, resistance, and diversity across the landscape. Our simulations suggest that harvesting and biological disturbances enhance forest diversity in the region, although this diversity did not contribute to the area of resistant and resilient forest.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
MV, RS, SM, ZS, and ZR conceived the idea and designed the methods. MV collected the data. MV, RS, and ZR analyzed the data. MV, RS, ZR, ZS, and MF led the writing of manuscript. All authors contributed to the article and approved the submitted version.