The opportunity for disturbance interaction occurs when one event follows another. Kane et al. (2017) observed that the likelihood of such occurrence depends on both the frequency and size of each disturbance, where opportunities for disturbance interactions increase by chance alone as their respective frequency and size increase. The nature of the specific interaction may be determined by a range of factors such as the disturbance mechanisms at play, the precise sequence, extent, and arrangement of the overlap, and the time since the previous disturbance (Kane et al., 2017; Figure 3). Interactions may include one or more types defined by the constituent components of the disturbance event—its incidence/extent, its intensity, and (or) its effects (Kane et al., 2017). Any of these components may be positively or negatively affected, or otherwise unaffected by the previous event. The nature of the interaction can also take multiple forms—for example the effect of one disturbance on another, or the combined effect of two disturbances on an ecosystem property or function.

From a disturbance-interaction standpoint, the legacy from one disturbance becomes the new system state for the subsequent disturbance (Peterson, 2002). A “linked” disturbance interaction occurs when the legacy of one disturbance is a critical driver for another, affecting the likelihood, extent, or severity of the second disturbance (Buma, 2015). As elaborated by Kane et al. (2017), the specific components underlying the interaction are important to the nature of that interaction. For example, logging disturbance does not generally affect either the incidence (frequency and extent) or intensity of windstorms. In this regard, a wind event following logging disturbance might be considered a “random co-occurrence” (sensu Burton et al., 2020). However, the severity (i.e., effects) of a given wind event can be strongly linked to logging practices – for example, taller trees (generally older) tend to be more susceptible to wind damage than shorter (generally younger) trees, hard edges created by specific logging practices can increase the susceptibility of trees along that edge, and tree species with different rooting habits or wood densities can have differing responses to winds of a given intensity (Quine et al., 1995). Hence, patterns of wind effects may well be linked to the legacy of logging practices. Conversely, if salvage logging is a commonly applied practice, then the incidence and extent of salvage logging will be linked to the disturbance patterns created by wind events.

The occurrence, relative strength, direction, and duration of linked disturbance interactions are mediated by both internal ecosystem processes and external drivers. For example, insect disturbance leading to tree mortality can create pulses of heavy fuels that can set the stage for higher severity fire once the new fuels become combustible, but the period over which that legacy persists as a factor affecting fire depends on fuel decomposition rates (Fleming et al., 2002; but see Section “Empirical Approaches, Challenging Long-standing Dogma”). In some boreal systems, harvesting can transition highly flammable conifer forests to fire-resistant deciduous forest (Carleton and Maclellan, 1994; Cumming, 2001)—a legacy mediated by internal ecosystem recovery processes. Opportunities for disturbances to interact at a given location are therefore sensitive to the timing and sequence of when and where a given disturbance event overlaps a subsequent one, contingent on any relevant time-lag effects (Burton et al., 2020). For example, while tree-killing disturbances can influence soil stability underlying the likelihood of landslide disturbances, this linked disturbance interaction is constrained by topography (i.e., slope angle and position) (Carabella et al., 2019).

Certain disturbances are more sensitive to forest composition or age structure than others, and this sensitivity will determine the relative strength of the feedback between vegetation, the specific disturbance, and (by extension) linked interactions with other disturbances. For instance, most outbreaking insects are limited to a few host tree species, leading to relatively strong feedbacks between vegetation and insect disturbance impacts across all three interaction types (incidence, intensity, and effects; Figure 3) (Ohmart, 1989; Hennigar et al., 2008; Bentz, 2019). When landscape spatial pattern is shaped by the cumulative effects of one disturbance (e.g., logging; Sturtevant et al., 2014), the consequent spatial pattern can either enhance or limit subsequent disturbances in terms of their intensity, duration, and size (e.g., insect outbreak; Robert et al., 2018). Such disturbance “legacies” often accumulate across space and time, producing “landscape memory” (Peterson, 2002) with important consequences for future disturbance regimes and consequent landscape dynamics (Foster et al., 1998; Buma and Wessman, 2011).

“Looped disturbances” form a subset of linked disturbance interactions in the form of recurrent disturbances with feedbacks that create a metastability of the system (Burton et al., 2020). Classic looped disturbances are often associated as a component of disturbance regimes supporting a systems’ historic range of variability (HRV; Landres et al., 1999). For example, boreal conifers with fire-dependent life-history traits (e.g., serotiny) are both adapted and conducive to infrequent but extensive crown-fire events. Likewise, many temperate pine systems once dominant across much of the southern United States (Hanberry et al., 2020) are adapted (e.g., thick bark, high canopy) and conducive (e.g., grassy flammable understory) to frequent surface fire events (Figure 4). Looped disturbances of the same type can occur wherever life history traits of the primary vegetation are adapted to given type of disturbance so that it replicates itself – such as the build-up of fir regeneration in the understory that is released by the periodic defoliation by eastern spruce budworm (Choristoneura fumiferana Clemens) in the North American boreal forest (Baskerville, 1975). More complex looped disturbances across disturbance types are also possible when those types help reinforce each other in sequence (as illustrated by Burton et al., 2020).

“Compound” disturbance interactions occur when the combination of two or more disturbances have a multiplicative effect over the state of the system (Buma and Wessman, 2011), often via the recovery processes affecting ecosystem resilience (Buma, 2015; Coop et al., 2020). In the above example for boreal conifers, two fires in rapid succession can burn the tree regeneration of a serotinous species before they have sufficient time to produce new seeds—dramatically delaying subsequent system recovery (Figure 4). Indeed, starkly different vegetation communities resulted from different combinations and sequences of fire and wind in Minnesota as the life history strategies of different tree species were affected quite differently by each disturbance type and their compound disturbance effects (Johnstone et al., 2016). By contrast, regeneration failure due to more frequent fire is a feature of pine savannah (Figure 4) (Hanberry et al., 2020). Decreasing fire frequency enhances regeneration success, leading to increased stem density and canopy connectivity that make the systems more susceptible to catastrophic disturbances such as bark beetles, crown fire, or both in combination (Negrón and Fettig, 2014; Clarke et al., 2016). The consequence of the restriction of a fundamental structuring agent in this example (frequent surface fires) is therefore a bifurcation of the system to the extremes: either closed canopy forest or grassland (Figure 4). Analogous patterns have been observed in many regions (e.g., western North America; Coop et al., 2020). The potential for compound disturbance is likely increasing as climate change can directly and indirectly affect disturbance processes (Soranno et al., 2014; Becknell et al., 2015), human land use intensification is disrupting natural disturbance regimes on a global scale (Foley, 2005), and increasing mobility and transportation increases novel introductions of pests and disease (Liebhold et al., 2017). Such broad-scaled changes can have non-linear “cascading effects” resulting in disturbance regime shifts and transforming the system to alternative states (Buma and Wessman, 2011; Buma, 2015; Johnstone et al., 2016; Ramsfield et al., 2016; Burton et al., 2020).

The effects of the disturbances and their interactions depends on a larger context of system dynamics. This larger context includes (1) external drivers operating at regional scales [e.g., biome constraints (Burton et al., 2020), anthropogenic pressure (Turner, 2010)] (Figure 1), (2) states of the system (e.g., vegetation composition, age structure, soil conditions, etc.) that determine the behavior and effects of (3) disturbance mechanisms (e.g., defoliation, stem breakage, and combustion) that in turn determine the new state of the system (Peters et al., 2011; Figure 5). The new state is considered the legacy of the disturbance, upon which internal ecosystem processes (e.g., recolonization, growth, and succession) enable subsequent system recovery (Peterson, 2002, Johnstone et al., 2016). The type of ecosystem response, as defined by Kane et al. (2017; Figure 3), depends upon how the disturbance legacy relates to the internal processes of the system. If the essential structure (i.e., tree density, size, and composition) remains effectively intact, the ecosystem response is one of resistance. If the system requires a longer recovery period but remains within the natural variability of the broader ecosystem, the ecosystem response to the disturbance interaction is one of resilience. If the recovering system falls outside that natural variability the disturbances compounded to result in an alternate system state (Figure 5). External drivers such as climate and climate change interact with disturbances via three pathways (Seidl et al., 2017) – directly via its effect on a disturbance property such as incidence or intensity, indirectly via effects on internal processes that modify system properties that in turn affect the disturbance process, and through a secondary interaction with another disturbance process that in turn affects the focal disturbance (i.e., the focus of this paper).

To understand processes affecting multiple disturbances and their interactions, integration of processes and disturbances operating at more than one scale is needed (O’Neill et al., 1996; Wu and David, 2002; Peters et al., 2007; Soranno et al., 2014) (Figure 1). Living and non-living elements of a system tend to interact most strongly with other elements that have characteristic process rates and spatial scales similar to their own (Allen and Hoekstra, 2015), leading to more closely-linked feedbacks within a single domain of scale (Messier et al., 2016). For example, at tree and stand-levels, life history traits of tree species affect demographic (mortality and recruitment), ecosystem (nutrient cycling), and disturbance (e.g., flammability, plant defense, etc.) processes that can include neighborhood effects that may also be self-reinforcing (Frelich and Reich, 1999). Nonetheless, the statistical properties of disturbance regimes (i.e., Figure 2) are not defined by processes at any one scale, but rather they emerge from the interactions of processes across scales (Figures 1, 5). So while compartmentalizing system dynamics into scale-specific components is a useful way to disentangle system complexity (Allen and Hoekstra, 2015), cumulative changes at one scale can either amplify or attenuate system behaviors across scales—a phenomenon known as a “cross-scale interaction” (Peters et al., 2004, 2007; Figure 3). For example, when local fuels become interconnected across landscapes to the point that firestorms (fire events that create their own weather) are possible, then discontinuities that previously served as fuel breaks may become irrelevant (Peters et al., 2004).

Cross-scale interactions often underlie dramatic change as a system shifts from a strong negative effect to a strong positive one. More precisely, different processes operating at different scales can work synergistically to amplify disturbances across scales. For example, tree-killing bark beetles subject to strong, non-linear and density-dependent feedbacks (e.g., Raffa et al., 2008) are often linked to other disturbance types. Under low-density endemic situations, individual healthy trees easily defend against attacking beetles. Outbreaks are generally triggered by drivers that enable beetle populations to grow sufficiently to overwhelm the resources of healthy trees. In Europe, windstorm events provide a pulse of food resources in the form of weakened or fresh-killed trees (i.e., a disturbance legacy) that enables European spruce bark beetle (Ips typographus L.) populations to increase to the point they can overcome adjacent healthy tree defenses, such that a strong negative feedback (strong defense) changes to a strong positive feedback (ample food) (Kausrud et al., 2012). Such an effect can be further amplified if wind events overlap periods of drought, where trees are already stressed and defenses weakened (Kausrud et al., 2012). Climate warming in central Europe serves as a common driver increasing the frequency and intensity of windstorms, the frequency and intensity of drought events, and decreasing generation time of I. typographus to dramatically increase population growth potential (Senf and Seidl, 2018, Figure 6). Climate-driven synergies underlying explosive outbreaks of this insect have been further exacerbated by a long history of spruce plantation establishment in the region (Jansen et al., 2017), leading to relatively continuous food resources (Seidl et al., 2016b, Figure 6). While this example appears like a “perfect storm” event, it is a pattern replicated across many forest insect systems across the globe, particularly in bark beetle systems subject to similarly strong non-linear feedback processes (Burton et al., 2020; Kneeshaw et al., 2021).

Current conceptual advances in disturbance interaction investigation include formal disaggregation of disturbances into their constituent components, embedding the disturbance process into system dynamics affected by both external drivers and internal processes (Figure 5), and explicitly recognizing the interaction factors, types, and ecosystem responses (Figure 3). While early modeling investigations of disturbance interactions focused only simple binary interactions (Figure 7; Mladenoff et al., 1996), the desire to include more processes structuring landscapes has spurred development and application of more complicated and complex models including many interacting disturbances having direct and indirect effects (Figure 7) – a theme we will return to in Section “Modeling Approaches.” Furthermore, under the Anthropocene, socio-economic drivers are crucial to include as well in forest disturbance models (e.g., Figure 6).

Importantly, the same drivers affecting disturbance mechanisms may also influence the internal processes affecting ecosystem response (Figure 5). The sum total of the drivers’ effects on disturbances, internal processes, and their interactions can lead to changes in system properties and system behaviors that can be either gradual or abrupt, depending on the direction and form of the effect (Ratajczak et al., 2018). Empirical disturbance interaction studies (see Section “Empirical Approaches”) should therefore take care to define the nature and form of the interaction explicitly, while models (see Section “Modeling Approaches”) should effectively capture the form, direction, and magnitude of interaction processes, including accounting for processes that transcend spatiotemporal scales.

The study of disturbance interactions in forested ecosystems has unique challenges because it is nearly impossible to manipulate for the overlap of two (or more) disturbances in the field. Consequently, there are many serendipitous correlational studies in the disturbance interaction literature. For example, the extensive bark beetle epidemics affecting western North America in the last two decades, coincident with both widespread drought (Raffa et al., 2008), and highly active fire seasons (Hart et al., 2015) enabled a broad range of opportunistic disturbance-interaction studies (Kulakowski and Veblen, 2007; Kulakowski et al., 2012; Stevens-Rumann et al., 2015; Agne et al., 2016; Millar and Delany, 2019), including the widespread study of salvage logging and its effect on ecosystem response and recovery (Leverkus et al., 2018). Such studies are often descriptive in nature, documenting patterns rather than understanding the processes producing them. Nonetheless, as the Anthropocene increasingly enhances disturbances under climate change, such disturbance overlap may become increasingly frequent (Seidl et al., 2016b; Burton et al., 2020; Kim et al., 2020).

One common approach for investigating disturbance interactions is “retrogressive” (sensu Simard et al., 2011), where investigators reconstruct disturbances over a longer period of time, using a variety of data sources such as archived satellite imagery, air photos, dendrochronology, and disturbance records (e.g., Bebi et al., 2003; Howe and Baker, 2003; DeRose and Long, 2012; Sturtevant et al., 2014; Hart et al., 2015; Meigs et al., 2015; Robert et al., 2018, 2020). A challenge facing retrogressive approaches is that the quality of data defining underlying drivers, such as forest composition or age structure, generally degrades the further back in time one investigates (McKenzie and Perera, 2015; Robert et al., 2020) (Table 1). Conversely, empirical data defining a given disturbance legacy (e.g., fuel) are often used to drive models to infer future disturbance behavior (e.g., flame length) (DeRose and Long, 2009; Simard et al., 2011; Meigs et al., 2015; Seidl et al., 2016a). Such “prospective” studies, including those referenced above, typically apply models at the plot or stand scale as their method to infer consequent system behavior.

A number of recent manipulative experiments have focused on disturbance mechanisms and their interactions on a variety of ecosystem responses. For example, Royo et al. (2010) examined the interactive effects of canopy gaps, ungulate browsing pressure, and fire via girdling, fencing, and prescribed fire, respectively. The authors found substantial synergistic effects leading to greater herbaceous diversity under the combined disturbances relative to any one disturbance alone. Quentin et al. (2012) applied a combination of irrigation and mechanical defoliation to investigate the interactive effects of water supply and leaf consumption on eucalyptus (Eucalyptus globulus Labill) plantations in Tasmania. Sparks et al. (2018) investigated the interactive effects of drought stress and fire on western larch (Larix occidentalis Nutt.) seedlings and found that severely drought-stressed seedlings had higher survival following fire relative to moderately drought-stressed seedlings, in part due to the physiological response to drought. Cannon et al. (2014) simulated windthrow disturbance by the static winching of trees to investigate consequent fire behavior, burn intensities, and consequent effects of those intensities related to elevated fuels on vegetation composition. Such studies can support the development of more process-based models enabling the simulation of disturbance interactions based on first principles (Figures 3, 5). The challenges in their application relate to their inherent limitations, such as the extrapolation of disturbance or stressor responses across tree life stages (Máliš et al., 2016; Dayrell et al., 2018), or the scale or intensity of experimental vs. natural disturbances (Cannon et al., 2014) (Table 1). Indeed, it is often the spatial dynamics of disturbances and their interactions that simulation models are designed to investigate, precisely because empirical measurement of such processes are difficult to quantify explicitly in space. Nonetheless, experimental studies can inform even spatial processes, such as the spatial concentration of ungulate herbivory associated with burned locations, and the cascading effects on plant species diversity across space (Westlake et al., 2020).

Despite the value of process identification and quantification, a whole systems approach is often necessary to fully understand the nature and behavior of disturbance interactions. For example, Anderegg et al. (2015) offered a data-driven conceptual model to capture the processes underlying interactions between temperature, precipitation, tree stress, and the different physiological, chemical, and ecological constraints influencing herbivory by bark beetles and defoliators that can lead to widespread mortality under climate change. Likewise, Cobb and Metz (2017) adapted the heuristic “disease triangle”, traditionally applied to forest and plant pathology, to understand the landscape ecology of tree diseases. The disease triangle consists of the pathogen, its host, and the environment conducive to the infection of the host by the pathogen. From the perspective of disturbance interactions, past disturbance events from various agents might contribute to the “environment” point of the triangle that can either facilitate or impede infection rates. Such syntheses can serve as the foundation for process-based landscape simulation models.

Case study syntheses further enable the investigation of disturbance interactions from a more holistic systems perspective. Evidence for the relative strength of specific interactions are built from diverse data sources—often at widely different scales—and gaps in knowledge may be addressed through a combination of prospective modeling and speculation based on systems and scaling theory (Peters et al., 2007). For example, Allen (2007) investigated cross-scaled interactions among drought, bark beetles, fire, and erosion in northern New Mexico based on a conceptual model that included non-linear dynamics that could either amplify or ameliorate the propagation of the different disturbances at play. Synthesized data sources spanned seven distinct studies ranging in spatiotemporal scale from weekly dendrometer (i.e., physiology) measurements to regional maps of inter-annual forest dieback developed via remote sensing. Such syntheses may be uniquely poised to address the sharp non-linearities associated with cross-scaled disturbance interactions that can lead to fundamental disturbance regime shifts. However, broader inference from such case studies will remain limited until the greater body of underlying theory can be rigorously tested (Allen et al., 2014) (Table 1).

Analysis of ecological data in recent decades has shifted from the use of analytically tractable statistical methods constrained by parametric assumptions, such as ordinary least-squared regression and analysis of variance, to more advanced but computationally intensive methods such as mixed effects models, Bayesian inference, and machine learning techniques (Touchon and McCoy, 2016; Holmes and Huber, 2018). Such recent methods are better able to address the complexities of ecological datasets that are rapidly accumulating, including data collected via remote sensing and other automated sensors, as well as time-series data from long-term monitoring and experiments.

The combination of advanced statistical techniques and rich data resources enabled the modeling of underlying processes, legacy interactions, and driver roles related to disturbance interactions in time and space. For example, For example, Mezei et al. (2017) examined temporal patterns in disturbance rates of the European spruce beetle as a function of wind disturbance, prior beetle disturbance, and climatic variables affecting brood frequency. The authors used generalized additive models to allow detection of non-linear relationships, and applied multi-model inference (Burnham et al., 2011) to select from a set of candidate models reflecting plausible alternative hypotheses of underlying processes. Hierarchical and hidden process modeling better accommodate the variability and uncertainties inherent in ecological research (Gimenez et al., 2014), as well as non-linear responses due to the interactions among disturbances and climatic drivers (Brice et al., 2020), leading to more complete understanding of complex phenomena like disturbance interactions. As such, Itter et al. (2019) applied Bayesian hierarchical modeling to tree ring data in eastern and western regions of boreal Canada to evaluate the separate and combined impacts of tree defoliation and drought on tree growth patterns. They found a synergistic effect of defoliation plus drought on the “ecological memory” of the disturbances in terms of their temporal recovery period. Thus, advanced statistical modeling techniques can be applied to complex datasets to determine the relative importance of disturbance types and factors that affect the state of forest ecosystems.

Similarly, machine learning techniques are often more appropriate to modeling complex and non-linear relationships (Olden et al., 2008). For example, Candau et al. (2018) applied random forest methods to define spatial domains in central Canada where fire and eastern spruce budworm defoliation were likely to interact, based on compositional (i.e., host) and climatic gradients (i.e., moisture). Deep learning methodologies such as convolutional neural networks have been applied to extract useful patterns from complex data, such as imagery and LiDAR, although their application has focused on extracting patterns rather than understanding the processes underlying those patterns (Brodrick et al., 2019). Machine learning techniques are both diverse and evolving, and do have limitations such as propensity to over fit models and (particularly in the case of neural nets) lower transparency in the model-building process (Olden et al., 2008; Liu et al., 2018). Nonetheless, such methods may detect relationships or domains of scale that are poorly understood and worth more explicit attention.

Regardless of the empirical approach, a critical role of empirical studies is to challenge long-standing assumptions of disturbance interactions that often find their way into landscape models. For example, the persistent dogma of “insect disturbance increases fuel load and therefore fire risk”—has been confronted with conflicting empirical evidence. Analysis of spruce eastern budworm outbreaks suggested that they increased forest fire probability in the first 20 years post-outbreak (Fleming et al., 2002), while Meigs et al. (2015) found that fire probability was lower after a western spruce budworm (Choristoneura freemani Razowski) event. These two studies were similar in that they were retrospective studies that evaluated whether the area burned increased during periods of budworm outbreaks, using a range of lag periods (i.e., time since outbreak). By contrast, James et al. (2017) investigated a similar question for spruce budworm in a similar geographic area to that of Fleming et al. (2002) but included weather drivers as covariates in the analyses. They found that the magnitude of defoliation effects was lower than those associated with weather. Likewise, a literature review of bark beetle-fire interactions in western North America concluded the effects of mountain pine beetle (Dendroctonus ponderosae Hopkins) disturbance on fire occurrence and severity were ambivalent (Simard et al., 2008). A subsequent study applied prospective modeling to focus on a particular fire process (i.e., probability of active crown fire; Simard et al., 2011). The authors considered the changes in fuel characteristics with time since disturbance (i.e., dead needles in canopy, deposited needles, fallen trees, and understory development), weather characteristics affecting fire behavior (wind speed and moisture) and interactions between the disturbance legacy and weather drivers (e.g., micrometeorology below canopy). As in Simard et al. (2011), James et al. (2017) found that the potential amplifying effects of the fuel disturbance legacy from the outbreak were constrained by fire weather conditions. Hicke et al. (2012) found that categorizing studies according to a conceptual model of expected lags in different fuel components as they related to different stages of bark beetle outbreaks led to more consistency in results and consequent fire implications across studies.

Empirical study of processes and drivers governing disturbance interactions are prerequisites of plausible models that simulate their dynamics. The different empirical approaches to the investigation of such interactions each have their role to play in the development and application of such models. Empirical study should justify when and where the added complexity of disturbance interactions is warranted, and confront persistent dogma based on perception rather than data. While interest in this line of research has increased in recent decades, disturbance interactions remain complex and challenging to quantify. Recent conceptual advances (see Section “Conceptual Advances”) are helping to refine the questions and terminology to channel empirical research, while advances in statistical methods, coupled with increasingly detailed data sources, provide additional tools to probe the characteristically messy data for useful insights. As evidenced by a series of recent reviews, research in this area is active but also unevenly distributed across geography, biomes, and disturbance types. By extension, such deficiencies limit the scope of landscape modeling in this line of research.

A specific class of models known as forest landscape models (FLMs) dominated the literature of landscape disturbance interactions, representing 92% of studies from our systematic review (Table 3). FLMs were first developed in the late 20th century (Baker, 1989; Mladenoff and Baker, 1999; Gustafson, 2013) and have greatly diversified since then (He and Mladenoff, 1999; Schumacher et al., 2004; Lischke et al., 2006; Scheller and Mladenoff, 2007; Keane et al., 2011; Seidl et al., 2012; Wang et al., 2014, among others). Landscape modeling software (Fall and Fall, 2001; e.g., Bolte et al., 2007) further enabled the rapid development of customized FLMs. Despite this diversification, a common FLM structure has emerged integrating vegetation dynamics (succession, growth, mortality, senescence, etc.) with spatially-explicit disturbances affecting those dynamics. In some fashion, both vegetation and disturbance dynamics may be influenced from the “bottom up” by biophysical drivers such as soil characteristics and terrain, and from the “top down” by climatic drivers such as precipitation and temperature patterns, and in some cases atmospheric conditions (e.g., CO₂ concentrations, flux of photosynthetically active radiation). Model development during the last two decades improved the architecture of the models, in terms of modular design and scaling methodologies, the relative balance between phenomenological patterns and simulated processes, and the diversity of disturbance agents that can affect system dynamics. Active development and widespread application of FLMs precipitated several reviews in the last decade (Seidl et al., 2011; Keane et al., 2013, 2015; Shifley et al., 2017).

Disturbance processes simulated within FLMs and other landscape models share common elements of initiation, spread or contagion, intensity, and effects, where opportunities for interaction can occur within each of these processes (Figure 3 and Table 3). Direct interactions between disturbances are often approximated through some surrogate for underlying process. For example, time since fire may serve as a proxy the process of fuel build-up. Indeed, Mladenoff et al. (1996) extended this approximation to enable one of the first formal disturbance interactions within an FLM, where a recent wind event increases the intensity of a subsequent fire event. Such logic becomes increasingly complicated, however, as one simulates multiple overlapping disturbances or other disturbance types, such as harvesting or insects (Figure 7). To what degree does the surrogate (in this case, either time since fire or time since disturbance) approximate the process (fuel accumulation) when multiple disturbances are at play? There is the potential for compounding error as additional direct disturbance interactions are added (Figure 7). This issue relates to any of the interaction types that can occur within a system (Figure 3). The complexity of state and transition models, for example, multiply with the number of disturbance processes that are included, because transition probabilities must be specified for every combination of possible disturbance sequences (James et al., 2011a; Daniel et al., 2017).

The solution to the increasing intractability of “direct” interactions via surrogates is to model the underlying process, such as fuel dynamics, as a state variable. Here, the interaction becomes indirectly mediated by the internal processes underlying the dynamics of the state variable, that in some cases, may be more elegantly modeled than attempting to keep track of overlapping surrogate variables (Figure 7). Such methods have other benefits, such as allowing for external drivers such as climate to further mediate the interaction (Figure 5), and to model the lagged effects appropriate for a given process and location (Burton et al., 2020). Simulation of multiple interacting disturbances is likely related to the decline of the state and transition model approach within FLMs in favor of tracking individual tree species that can respond individualistically to a wider range of disturbance types and their various combinations. Only 10 percent of the studies in our review applied the state and transition model approach (Table 3). A more recent innovation in FLMs is the simulation of tree stress via carbon balance (e.g., Keane et al., 2011) or even more explicitly in the form of non-structural carbohydrate reserves (Seidl et al., 2012; DeBruijn et al., 2014). Such innovations enable trees to respond to the accumulation of multiple stressors that include both disturbances and environmental drivers.

A related model design choice relates to trade-offs in “static” versus “dynamic” disturbance regimes. A static disturbance regime is one parameterized according to its statistical properties (e.g., rotation, event size distribution, intensity distribution, etc.; Figure 2). The focus of such applications is typically the response of the system to a given disturbance regime, without feedback to the regime itself. A fully dynamic disturbance regime is one where the statistical properties of the regime are an emergent property of the system dynamics (Perera et al., 2015). In practice, any of the individual components of a given disturbance type (i.e., Figure 3) can be simulated as either static or dynamic. Most contemporary FLMs simulate disturbance effects dynamically, where mortality is contingent on the relative susceptibility of forest types, tree species, age classes, or some combination of these variables. Many introduce stochastic elements that may be related to process (e.g., likelihood of insect damage or fire ignitions related to vegetation characteristics). Some models offer different options where the user can specify the relative balance of static vs. dynamic components (e.g., Sturtevant et al., 2009; Keane et al., 2011; DeJager et al., 2017).

We developed a classification scheme to address the degree of process integration (i.e., balance of static vs. dynamic components) within a specific disturbance agent within a given landscape disturbance application (Table 4), and then applied that scheme to the different disturbance agents simulated across the studies identified by our systematic review (Table 3). While the classification scheme was both logical and readily applied given proper documentation – we found this supporting documentation was sometimes lacking within individual studies. Many referenced earlier model applications, or original model publications, requiring more extensive review to fully understand the nature of the disturbance interactions simulated. Other applications required assumptions on our part where the documentation was insufficient, including vague descriptions of the specific model options applied, or references to model code requiring literacy of the underlying programming language. Nonetheless, we were able to document ranges of process integration across disturbance types and model software platforms (Table 3), including dominant modes of process integration, and instances where more detailed process integration was possible via the coupling of models (see Section “Coupled Models”). We found that while the trend in disturbance modeling is toward more explicit modeling of disturbance processes and their emergent effects on system dynamics and landscape structure (Gustafson, 2013; Keane et al., 2015), static disturbance regimes are still commonly applied. The relative balance of static vs. dynamic components in application was also unevenly distributed across different disturbance types and model approaches – including non-FLM applications (Table 3).

Some contemporary FLMs are better described as model frameworks, where different trade-offs in process vs. surrogates, static vs. dynamic disturbance regimes, and the number and types of disturbance mechanisms simulated may be customized according to both the question addressed and the available supporting data. The most frequently applied model software from our review was LANDIS-II, an FLM model framework accounting for 59% of all studies (Table 3). A case study pair of applications of this model in north central Minnesota (United States) illustrates the framework approach (Lucash et al., 2017, 2018) (Figure 8). The authors applied a highly mechanistic succession option (Scheller et al., 2011a) where tree-species cohorts respond individualistically to different disturbance types and compete for above and belowground resources (light, water, nitrogen) via a set of life history traits and parameters. The first study focused on how the combination of climate change and forest management affected forest resilience following wind events (i.e., forest recovery rates), applying a factorial experiment that crossed climate scenarios with forest management scenarios (Lucash et al., 2017). The second study focused on the degree to which multiple disturbance regimes affected tree mortality rates (i.e., effects). Results of the latter supported antagonistic interactions, where total tree mortality caused by multiple disturbance regimes applied simultaneously was less than the sum of the mortality from each of the respective disturbance types applied independently. Closer examination of the simulation methods reveals a combination of process integration methods were applied to different disturbance types (Lucash et al., 2018). Wind disturbance was simulated as “mostly static” regimes with characteristic size, shape, and rotations of three classes of wind patterns (microburst, tornado, and derecho), each with species and age-specific impact patterns. Such static regimes are consistent with the lack of feedback between vegetation and wind event frequency or intensity, but cannot address spatial factors such as the effect of forest edges on tree susceptibility to windthrow (see Section “Disturbance Interaction Process and Terminology”). Fire regimes were simulated as mostly dynamic functions of the vegetation, fuel dynamics, and climate-related fire weather. Limitations of this approach in terms of emergent fire behavior are shared across virtually all FLMs: rules governing fire spread do not account for critical feedbacks between fire behavior and weather that fundamentally change those rules (such as long-distance dispersion of fire embers) (see Section “Cross-Scale Interactions”). Recent wind or insect disturbance could modify fuel attributes as direct interactions (Figure 7). The three different insect disturbance agents included intermediate blends of dynamic and static regime elements, where the extent of each disturbance was directly related to availability of vulnerable host, but the frequency and intensity of outbreaks was predefined. The exception was the patterns of mortality in oak (Quercus spp.) species by the two-lined chestnut borer (Agrilus bilineatus Weber) that occurred in instances where forest tent caterpillar (Malacosoma disstria Hübner) defoliation overlapped in time with drought events, reflecting a surrogate for tree stress defined by direct interactions with defoliation and climate, respectively (Figure 8).

We summarized the class of explicit disturbance interaction questions addressed across the different disturbance types (Table 5), using the conceptual advances summarized in Figures 3, 5. Frequencies of studies across questions and disturbance types were limited to those who either quantified the effect via their experimental design, or explicitly defined the nature of the simulated interaction within their methods. “Ecosystem response” was the one exception (see below). Linked disturbance interactions were specified according to which disturbance types were affected by another type. While looped disturbance types of the same disturbance (e.g., the effect of fire frequency on future fire likelihood) were not included in Table 5, instances where disturbance of a given type included multiple agents of that type (e.g., more than one insect agent) were included in the tabulation. Linked disturbances generally related to questions of disturbance resistance (Buma, 2015), but some studies aggregated responses across multiple disturbance types. In these cases, we tallied the study under “resistance (general)” (Table 5). Articles addressing resilience interaction questions were defined as studies that quantified the rates of return to a reference condition, shifts to alternative vegetation states, direct comparisons with HRV, or some combination of these responses. Climate-mediated disturbance questions refer to those studies where a given disturbance type was either modified to be consistent with future climate, or dynamically linked to climate. In the case of land management, climate-mediated disturbance involved some form of intentional climate adaption, such as assisted migration via planting. Cross-scaled disturbance frequencies refer to those studies that quantified cross-scaled interactions via their experimental design. Ecosystem response frequencies refer to those studies quantifying state variables other than those expressly linked to the disturbance interaction questions above (e.g., compositional trends). We note that most studies, whether self-identified as disturbance interaction studies or not, have an interest in capturing system dynamics under multiple interacting disturbance regimes. Frequencies of ecosystem responses in Table 5 therefore reflect the number of studies that included ecosystem responses regardless of whether the specific effect of a given disturbance type could be quantified in those response variables via the experimental design.

The most commonly investigated linked disturbance was fire (Table 5), with most investigations focused on how land management actions such as fuel treatments or logging affected fire severity or extent. Far fewer studies investigated fire as linked to insect disturbance. Insect disturbance was investigated as linked to a broad range of disturbance types in half the studies for which it was included. Wind was less often investigated as linked to other disturbances, with the largest exceptions from studies in central Europe. Studies including the major disturbance types (management, fire, wind, and insects) investigated questions related to system resilience in roughly a quarter of their respective cases. For less-studied disturbance types, questions addressed varied widely, though the sample size was low. Drought was primarily investigated as a climate-mediated disturbance (16 of 22 studies) and related questions of forest resilience (6 of 22 studies). Six of 15 studies including mammalian browse investigated questions of ecosystem resilience, often related to the loss of forest cover. Almost all of the studies including mass movement disturbances focused on linkages to other disturbance types that either removed forest cover or otherwise destabilized soil in steep terrain.

Half of the studies from our review included some form of climate mediation of disturbance regimes (Table 5). Of those, 61 studies evaluated the degree to which forest management scenarios could mitigate the less desirable effects of climate change, such as elevated mortality (i.e., resistance) or transition to alternative system states (i.e., resilience). For example, Cadieux et al. (2020) investigated how projected climate change interacted with logging and fire disturbance to affect bird species habitat in the boreal forests of Alberta, Canada. In this case, the focus was on the vegetation response rather than on disturbance processes. They therefore applied mostly static fire and drought regimes parameterized differently for each climate scenario, based on empirical studies (Boulanger et al., 2014; Chen et al., 2016, respectively). The combined effects of climate change in combination with the associated elevated fire and drought mortality resulted in transitions from conifer-dominated forests to either deciduous forest or open habitats, with important implications for bird habitat. How this change to alternative forest states might interact with the fire regime itself could not be quantified using this approach. Keane et al. (2019) explicitly focused on emergent fire dynamics related to climate change, and how such behavior might be mitigated via either fire suppression or fuel treatment rates. In this study, they used simulations under historic climate and without human intervention to define the multivariate HRV for three contrasting landscapes in the northern Rocky Mountains (United States). This approach, combined with their process-based modeling framework (Keane et al., 2011), enabled them to draw inferences with regards to tipping points in system behavior that indicated loss of system resilience relative to the benchmark HRV, and the degree to which human action might affect system resilience in the face of climate change.

The relative balance of static to dynamic process in disturbance agents define the degree to which reciprocal interactions among those agents may be investigated, and therefore must be specified carefully with regard to the questions asked and addressed. For example, Chen et al. (2011) investigated the effects of larch caterpillar (Dendrolimus superans Butler) defoliation of different intensities on fire regime. However, key assumptions governing the fire regime, such as whether fire size and rotation was an input or an output of the model, were not provided, making it impossible to evaluate their results regarding insect effects on response variables such as fire frequency. Sturtevant et al. (2012) investigated reciprocal interactions between spruce budworm and fire disturbance and demonstrated that budworm defoliation disturbance reduced live ladder fuels, and therefore decreased fire extent over long time scales (centuries), despite increasing fire spread and intensity at shorter (decadal) scales. Specifically, the authors examined how an intermediate static-dynamic insect disturbance regime affected a mostly dynamic fire regime. While they could directly address the question of budworm effects on fire regimes, they could not fully address the reciprocal question (i.e., how does fire influence the insect disturbance regime?). Temperli et al. (2015) coupled a spatially-explicit bark beetle model with the LandClim FLM (Schumacher et al., 2004) to investigate interactions among spruce bark beetles (Dendroctonus rufipennis Kirby), fire, wind, and climate. In this case, wind disturbance was simulated as a mostly static regime, fire was intermediate (i.e., ignitions and size distributions were static inputs, but spread and effects were dynamic), drought stress was linked to climate inputs, and beetles responded dynamically to vegetation, drought, host abundance, and windthrown host resources. Analogous to Sturtevant et al. (2012), the authors could address the effects of wind and fire on bark beetle disturbance behavior, but not the reverse. They found that climate change-driven drought events would promote beetle outbreaks, but the combination of the two mortality agents progressively reduced long-term susceptibility to bark beetle outbreaks. As a final example, Seidl and Rammer (2017) simulated emergent dynamics of bark beetles responding to wind disturbance and climate change in the Austrian Alps, concluding that climate change amplified the outbreaks initiated by wind events, in large part by shifting the optimal development rates of the beetles to higher elevations where its spruce host was more concentrated. Since the authors applied disturbance agents that each responded dynamically to both climate drivers and ecosystem state variables, they could make inferences about reciprocal interactions among disturbance agents, as well as the mediating role of climate.

Some of the constraints on the types of disturbance interaction questions that may be addressed via landscape modeling are imposed by the design constraints of the selected model software (e.g., Table 3). Nonetheless, model designs are ultimately constrained by supporting knowledge and data. In the insect examples above, there is greater consensus on the landscape factors controlling bark beetle dynamics than there are for defoliator dynamics (Kneeshaw et al., 2021). Indeed, bark beetles were included in all four of the disturbance interaction studies that explicitly addressed cross-scaled interactions (Table 5). Of those, one approximated non-linear responses of bark beetles to host abundance via percolation theory (Seidl et al., 2016a). The remaining three explicitly modeled emergence and spread of bark beetle agents as they interacted with environmental drivers such as temperature, food resources including both pulses provided by windthrow and host concentrations impacted in part by drought stress (Temperli et al., 2015; Dobor et al., 2020; Honkaniemi et al., 2020) (see Figure 6). While a few model designs may have allowed for cross-scale interactions in other studies (Table 3), they were not explicitly quantified so their importance cannot be evaluated.

Forest management was frequently investigated as a linked factor affecting other disturbance types, however, the reverse (effects of natural disturbance on land management) was far less common (Table 5). Eleven of the studies linking management to other disturbances included explicit salvage logging and/or sanitary logging, representing post-disturbance harvest of felled trees and preventative treatments to reduce insect impacts, respectively. Given the quantity of empirical studies investigating salvage logging, it is somewhat surprising this process is less investigated via landscape modeling.

FLMs most often model forest management practices as partially dynamic functions of landscape conditions (Table 3). Harvest locations are often selected via ranking algorithms such as oldest harvested first or biased toward species or forest types of interest (Gustafson et al., 2000). Harvest intensity and effects can likewise be defined to reflect a wide range of management prescriptions. However, the rates (i.e., extent) of any given treatment is generally predefined within scenarios. By contrast, forest planning models represent a different class of models that approach the question of how much to harvest dynamically based on the temporal changes and sometimes spatial distribution in the timber resource (Pukkala, 2013). Disturbance interaction questions from a forest optimization point of view often revolve around accounting for natural disturbances that can add uncertainty to the optimization process (Savage et al., 2011). There is a well-developed literature on spatial forest planning models that focus on the drivers of timber management such as growth and yield, access, road-building, and especially cost-benefit trade-offs that FLMs generally lack (Baskent and Keles, 2005; De Pellegrin Llorente et al., 2017). Nonetheless, such approaches have struggled with the stochastic and contiguous nature of natural disturbances, and therefore either simplistically account for such processes as aspatial rates of disturbance or ignore them altogether (De Pellegrin Llorente et al., 2017). Occasionally these divergent approaches have been blended to account for both human and ecological drivers of disturbance interactions. For example, forest optimization models have been used to define harvest rates for scenarios implemented in an FLM (e.g., Daniel et al., 2017). Raulier et al. (2014) used the SELES modeling platform to apply an iterative harvest optimization process – one that included spatially contagious processes such as road building and salvage logging where existing roads enabled rapid access. They used the model to simulate reciprocal interactions between fire and timber harvest to better quantify risk to timber supply under climate change in boreal Canada.

More recent methods for integration of socioeconomic drivers that contribute to both reciprocal interactions with forest landscape systems include demand and allocation models and also agent-based models (An, 2012; Egli et al., 2019). Developers recently incorporated agent-based logging and salvage logging within the iLAND model (Rammer and Seidl, 2015). Such development enabled simulation experiments examining reciprocal interactions among wind, European bark beetle, and forestry practices in the face of climate warming (e.g., Dobor et al., 2020) as illustrated in Figure 6, and account for the studies linking land management to wind and insect disturbance in Table 5. Agent-based forest management approaches are now being implemented in other FLM platforms (Spies et al., 2017; Sotnik et al., 2021). Land-use intensification and other anthropogenic disturbance activities have profound effects on landscape legacies and disturbance regime shifts (Watson et al., 2014). An emerging hotspot of research activity therefore lies in the investigation of reciprocal interactions between human activities and ecosystem processes—falling broadly under the scope of coupled human and natural systems (Wimberly et al., 2015). Anthropogenic and natural disturbance interactions might be modeled using integrated modeling platforms (Verburg and Overmars, 2009), where the combined effects of disturbance events and climate change can be explored via experimental designs that include either contrasting or a range of climatic and land-use intensification scenarios (Becknell et al., 2015). While six studies from our review used outputs from human land use change models to drive land cover change within FLMs (e.g., Duveneck and Thompson, 2019), none coupled land used change dynamically with forest harvest or natural disturbances to investigate their reciprocal interactions. One study applied an FLM to estimate forest productivity under alternative climate change scenarios and coupled those outputs with that of yield models from other economic sectors (agriculture and livestock) to drive an economic land allocation model to simulate land conversions under the twin drivers of climate change and human incentives (Briner et al., 2012). This example represents the closest approximation to the complete system portrayed in Figure 1 of any included in our review.

As illustrated by the coupling of FLMs with land cover change models in the previous section, some authors of articles in our review extended the capabilities of a given model by coupling them with one or more additional models. Different models may be linked via a “loose coupling,” where the output of one model is used as input for another (i.e., the most common form, see Section “Modeling Anthropogenic Forest Disturbance”), or “tight coupling”, where the dynamics of one model interact with the dynamics of another (Sturtevant et al., 2007). Loose couplings of models are widely used as a scaling technique, for example by deriving the input parameters of a coarse-scaled model from a finer-scaled process model (Urban et al., 1999; Xu et al., 2012; Boulanger et al., 2018; Huang et al., 2018) or down-scaling vegetation or disturbance behavior from coarse-scaled models such as DGVMs (Halofsky et al., 2014; Turner et al., 2015) (Table 3). FLM outputs are also commonly used as inputs for wildlife habitat or viability models that may be considered a dimension of ecosystem response to disturbance interactions (e.g., Longru et al., 2010; Scheller et al., 2011b; Nitschke et al., 2020).

Apart from human factors such land cover change or forest optimization, use of coupled models to address disturbance interactions centered on better representation of hydrologic and mass movement processes currently lacking within FLMs (Table 3). Research by Saksa et al. (2020a, b) coupled an ecohydrology model with a stand growth and yield model (for vegetation dynamics) and a fire behavior model to evaluate the integrated effects of fuel treatment patterns on fire risk and hydrologic outputs from small watersheds of the Sierra Nevada Mountains in the western United States. A planar surface model was used to define the flood disturbance regime for an FLM application to a Mississippi floodplain system in the central United States (DeJager et al., 2019). The combination of forest cover, disturbances, and impervious surface outputs from a coupled land cover change and FLM model system were used to project flood risk via a simple hydrologic model (Thompson et al., 2016). Scheidl et al. (2020) used the projections of an FLM to drive a physically-based spatial hydrologic model coupled with a slope-stability model to evaluate the integrated effects of climate change, forest management, and natural disturbance processes on shallow landslide frequency within the eastern Alps of Austria.

Hybrid empirical-modeling studies leverage the synthetic knowledge within process models to relate data representing measurable system components to processes underlying disturbance interactions that are either difficult or impossible to measure directly. For example, McGuire and Youberg (2019) empirically investigated the effects of repeat burns on both soil water infiltration rates and the occurrence of debris flow disturbances in dry montane forests of southeastern Arizona, United States. They then coupled remote sensing and terrain inputs with a high-resolution hydrologic model to demonstrate that differences in soil infiltration between once and twice burned soils were sufficient enough to affect threshold rain intensities necessary to generate debris flows. Similarly, Hancock et al. (2012) measured windthrown trees and consequent erosion patterns following a cyclone impacting northern Australia. They then used this data within the geomorphic landform model SIBERIA (Willgoose et al., 1991) to investigate ecohydrological controls (including windthrow patterns) on long-term erosion and catchment evolution. Økland et al. (2016) investigated a critical transition in bark beetle outbreak dynamics where beetle populations move from trees injured or killed by wind events to attack healthy trees (i.e., windfall-driven to patch-driven dynamics). They compared alternative models representing different hypothesized drivers within an individual-based model framework and found that the model driven by beetle aggregation was the most consistent with their empirical data on the spatial progression of an outbreak in Tatra National Park in Slovakia. Analogous hybrid modeling by Seidl and Rammer (2017) now serves as the foundation for the individual-based model of beetle disturbance contained within the iLAND FLM.

A primary challenge to conducting a systematic and consistent review of landscape models featuring disturbance interactions was the lack of standards by which to classify the implementation of a given disturbance agent in terms of its relative degree of process integration - or for that matter, documentation standards for FLMs more generally. There is an increasing call for standardization of methods in analogous classes of models, including agent-based model designs (Grimm et al., 2010) and species distribution models (Zurell et al., 2020). Table 6 represents a partial example of how documentation standards relevant to disturbance interactions can convey the nature of the interactions being simulated and investigated. Minimum standards should include the disturbance interaction type simulated (Figure 3), the response variable quantifying its effect, and the dynamic components of a given disturbance agent that dictate which interaction types are justifiably investigated. This documentation would facilitate a rapid classification – like that in Table 4 – so that the nature of the simulated interaction is transparent to the reader. Other essential documentation should include a clearly articulated study design and basic information such as the spatial and temporal grain and extent of simulated processes (Table 6). The latter is coming increasingly important to clarify as the temporal and spatial resolution (grain) has become more specific to individual processes as model designs become more refined (Table 6). Additional relevant documentation includes the level at which disturbance effects are applied (patch/cell, age class, species, cohort, tree, etc.), the spatial implementation of the disturbance agents (spread, dispersion, etc.), and how the various disturbance agents are interspersed in time. All such documentation should be specific to the model application at hand, rather than the potential capabilities of a given model software, as this was sometimes conflated in the articles we reviewed. Outputs from spatial models are complex enough without confusing artifacts of modeling assumptions with emergent behavior based on underlying process.

We do not assert that complex models are better than simple models. Rather, that one should account for the processes underlying the question being investigated to the extent that data and knowledge can support them. The contemporary argument for design of models based on first principles is that such models should make more robust predictions under the increasingly novel conditions of changing climate, atmospheric composition, and disturbance agents (Cuddington et al., 2013; Gustafson, 2013). We note also that many landscape modeling studies did not explicitly investigate disturbance interactions per se, but rather included multiple disturbance agents to simply account for all the relevant structuring processes within their forest system. Regardless, it is important to acknowledge and account for the cumulative effects of model error propagation as the number of interacting components increase and from one scale to another (Rastetter et al., 1992). To be more explicit about the relative importance of multiple processes and errors, Cressie et al. (2009) promoted the use of hierarchical statistical modeling that allows the quantification of the relative importance of uncertainty stemming from the data model (i.e., data error), the process model (i.e., model specification error), and the parameter model (i.e., parameter estimation error). There are, however, multiple additional sources of uncertainties—ranging from the data, models, processes, knowledge of the processes, and climate scenarios—that will affect the precision and accuracy of model projections.

To quantify the relative importance of these various types of uncertainty, one could compare the outcomes of several models (Cheaib et al., 2012). A recent comparison of four “state of the art” FLMs demonstrated wide divergence in model outcomes, in large part due to design choices affecting vegetation response to disturbance (Petter et al., 2020). This result suggests a certain degree of humility when reporting the results of landscape models – particularly where interactions by multiple disturbance agents may compound error associated with cumulative vegetative responses.

Scheller (2018) suggested that parameterization and “model confidence” are twin challenges facing contemporary landscape modeling studies. The increased burden of input parameter specification that generally comes with replacement of simple statistical relationships with more mechanistic detail may potentially undermine the perceived robustness of process-based models if those parameters are either poorly understood or quantified (McKenzie et al., 2019). Relevant knowledge acquisition necessary for parameter estimation will require multidisciplinary research crossing different domains in scale, integrating experimental research, monitoring of earth systems, advanced statistical analyses, and next generation models assembled from these insights to determine process-based responses to future climatic conditions. Knowledge gaps may further require more novel hybrid modeling frameworks (e.g., Talluto et al., 2016) that integrate process-based and statistical models to include both bottom-up processes and drivers and top-down drivers and account explicitly for feedbacks across scales (Figure 1). Such integrative models require knowledge about ecological responses based on targeted experiments in the field and carefully designed data analysis of existent data (Poisot et al., 2016).

Model confidence refers to the degree to which independent readers trust the utility of the model results – in this regard independent validation data remains a fundamental issue underlying confidence in FLMs (Shifley et al., 2017; Scheller, 2018). Validation issues may be particularly acute in the area of disturbance interactions given system sensitivity to the timing and order of disturbance mechanisms and the individualistic responses of the vegetation. While there were examples within our review that evaluated their results using empirical observations (Henne et al., 2013; Karam et al., 2013; Schwörer et al., 2014; Thom et al., 2018), these examples were rare. Validation of model components in combination with targeted uncertainty and sensitivity analyses may also improve model confidence (Cressie et al., 2009; Xu et al., 2012; Scheller, 2018; McKenzie et al., 2019). More specifically, model uncertainty refers to the degree to which model inputs are known as well as their associated measurement error, whereas sensitivity analyses focuses on the degree to which model inputs have influence on the model outputs (McKenzie et al., 2019). Guidance on effective uncertainty and sensitivity analyses can be found across disciplines for simulation modeling more generally (e.g., Saltelli et al., 2006; Wagener and Pianosi, 2019).

We found a broad spectrum of questions related to disturbance interactions in forested landscapes have been addressed by an equally diverse range of modeling approaches (Tables 3, 6). Notably, it is the more intractable of these questions that are also the most challenging to quantify empirically, and these require the support of simulation models to investigate. Integration of recent conceptual advances (i.e., linked disturbances, compound disturbances, cross-scaled interactions, and cascading effects) often require more detailed processes to simulate. Among the most important are non-deterministic tree-species competitive interactions and disturbance responses via life history traits, explicit modeling of stand structural characteristics, physiological traits of tree species that enable dynamic responses to climate and biophysical variables, and various integrations of these factors (Keane et al., 2011; Scheller et al., 2011a; Seidl et al., 2012; DeBruijn et al., 2014; Wang et al., 2014). Depending on the questions asked, modeling disturbances require ecological detail necessary to effectively capture their incidence, extent, intensity, and/or effects (see examples within Section “Questions and Approaches Across Disturbance Types”). We are approaching a point where some of the most process-based disturbance models, such as physical fire models, individual-based insect models, and agent-based human disturbances are possible either within FLMs or coupled model architectures (Rammer and Seidl, 2015; Bentley and Penman, 2017; Seidl and Rammer, 2017). Nonetheless, trade-offs between model elegance and complexity remain (McKenzie and Perera, 2015). For example, the parameter demands of highly complex models may exceed the data available to parameterize them (Scheller, 2018) (see also Section “Model Transparency, Uncertainty, and Confidence”). We also observe that the most mechanistic examples reviewed here were comparatively limited in the extent at which they were applied (e.g., thousands versus millions of ha, years versus centuries), illustrating practical trade-offs between ecological detail and scope remain despite advances in computing power (Table 6).

The study of disturbance interactions through landscape modeling reflects the asymmetry in systems and disturbance agents of the empirical studies (see Section “Empirical Approaches”). Studies were concentrated primarily and temperate, boreal, and alpine regions of North America and Europe (Table 2) and were almost absent from the southern hemisphere. Perhaps it is no coincidence that the prevalent disturbance regimes of subtropical and tropical regions (e.g., mass movement, hurricane/cyclones, mammalian grazers; Burton et al., 2020) were less commonly investigated in the context of disturbance interactions (Table 5). We found that human disturbance – by far the most commonly simulated disturbance type – was primarily simulated as a driver rather than an interactive component of the simulations. Cross-scale interactions remain an important concern that are particularly challenging to model and quantify. Each of these deficiencies suggests fruitful areas of future research. Finally, we echo the plea made by previous authors to improve model transparency – perhaps via documentation standards such as that developed for other classes of models.