Gap models across micro- to mega-scales of time and space: examples of Tansley’s ecosystem concept

Abstract

Background

Gap models are individual-based models for forests. They simulate dynamic multispecies assemblages over multiple tree-generations and predict forest responses to altered environmental conditions. Their development emphases designation of the significant biological and ecological processes at appropriate time/space scales. Conceptually, they are with consistent with A.G. Tansley’s original definition of “the ecosystem”.

Results

An example microscale application inspects feedbacks among terrestrial vegetation change, air-quality changes from the vegetation’s release of volatile organic compounds (VOC), and climate change effects on ecosystem production of VOC’s. Gap models can allocate canopy photosynthate to the individual trees whose leaves form the vertical leaf-area profiles. VOC release depends strongly on leaf physiology by species of these trees. Leaf-level VOC emissions increase with climate-warming. Species composition change lowers the abundance of VOC-emitting taxa. In interactions among ecosystem functions and biosphere/atmosphere exchanges, community composition responses can outweigh physiological responses. This contradicts previous studies that emphasize the warming-induced impacts on leaf function.

As a mesoscale example, the changes in climate (warming) on forests including pest-insect dynamics demonstrates changes on the both the tree and the insect populations. This is but one of many cases that involve using a gap model to simulate changes in spatial units typical of sampling plots and scaling these to landscape and regional levels. As this is the typical application scale for gap models, other examples are identified. The insect/climate-change can be scaled to regional consequences by simulating survey plots across a continental or subcontinental zone. Forest inventories at these scales are often conducted using independent survey plots distributed across a region. Model construction that mimics this sample design avoids the difficulties in modelling spatial interactions, but we also discuss simulation at these scales with contagion effects.

Conclusions

At the global-scale, successful simulations to date have used functional types of plants, rather than tree species. In a final application, the fine-scale predictions of a gap model are compared with data from micrometeorological eddy-covariance towers and then scaled-up to produce maps of global patterns of evapotranspiration, net primary production, gross primary production and respiration. New active-remote-sensing instruments provide opportunities to test these global predictions.

Introduction and background

In this paper, we will provide examples of different models, all of which have are unified in their use of modeling forest dynamics, but operate over different time and space domains. These models simulate the physical structure of forests across their respective domains. Over the time of development of these models, there has been a parallel development of a remote-sensing capability to observe change are associated at the micro-, meso-, macro- and mega-scales shown in Fig. 1. In this paper, we present examples of individual-based forest models, notably “gap models” to utilize these new data, to test models and to generate forest-ecosystem predictions and theories.

Fig. 1

(a) Environmental disturbance regimes (including climate) across space and time scales with (b) biotic responses of forests. Disturbance events here include wildfire, wind damage, clearcut, flood, earthquake, etc. Blue shaded area is region of high predictability (sensu Wiens 1989), matching appropriate processes to disturbances at similar scales. Although landscapes may be viewed as scale independent, they are often viewed as a scale larger than stands in practice and encompass scales at which ecosystem processes connect forest mosaics (Lertzman and Fall 1998). (c) Classification of vegetation with scale. (d) Representative existing forest datasets by scale. Dashed lines divide conceptual scales from local to mega, including mesoscales (in red), which bridge the responses of vegetation from gaps and forest succession to species migration. Reprinted from Druckenbrod et al. (2019) with panels a-c adapted from original in Delcourt et al. (1983). See also Prentice (1992) for a similar analysis as in panel d

New technologies in remote sensing (RS) are providing rich challenges and opportunities to increase the understanding of forest ecosystems. These technologies can provide new observations of structural and functional traits to examine patterns and processes of ecosystems at different spatial and dimensional resolutions. The consequences of the interactions between pattern and processes, which is the yang and yin of terrestrial ecology, are at a relatively advanced level in forest ecosystem science, but there is still much more to learn, and results from research employing these new technologies can speed the learning processes. When Tansley originally coined the neologism, “ecosystem” in 1935, he made the point of that, “These ecosystems, as we may call them, are of the most various kinds and sizes.” (Tansley 1935, page 299). Diagrams of scale using panels with time-and-space scales as axes of (1) disturbances or drivers of ecosystem change at particular time-and-space scales, (2) processes that respond to these drivers at equivalent scales and (3) patterns in ecosystems that arise at these scales (e.g. Delcourt et al. 1983; Druckenbrod et al. 2019) represent the “kinds-and-sizes” part of Tansley’s definition. The union of drivers, processes, and patterns of responses define an ecosystem sensu Tansley — an ecosystem is a system as defined by its inputs, processes and outputs, all at commensurate space- and time-scales.

For example, consider a forested watershed in a given location as an ecosystem. Define precipitation, humidity and temperature as the inputs, transport and evapotranspiration as processes, and the variation in streamflow as the pattern. A second ecosystem might be a patch of forest inside the watershed where inputs of photosynthetically active radiation (PAR) drive photosynthesis to produce daily productivity patterns. While different from one another, both of these example forest ecosystems could be collocated and observed simultaneously. Sampling to characterize these two forest ecosystems would involve measuring different variables resolved across different time- and space-scales. Often, but not always, response times of the processes included in ecosystem formulations decrease with increasing with spatial scale (Druckenbrod et al. 2019).

The multiple-scale aspect of Tansley’s ecosystem implies a matching of the space-and-time domains of ecosystem drivers, processes and responses in spatial pattern (Wiens 1989), as in the examples just mentioned. Druckenbrod et al. (2019) recent review of these ecosystem concepts adds a fourth panel representing the time- and space-scales of available data sets (Fig. 1, also see Prentice 1992). Air- and satellite-borne remote-sensing represent different time-space domains from other data collections. These data often are collected across vast spatial scales with approximately daily-to-weekly sample-return intervals. The cumulative period of data collection for these instruments is at multiple decades in some cases (Fig. 1).

These remote-sensing (RS)-based data sets are an upscaled representation of processes normally studied at relatively small areas but not observed over areas (global-, regional- and landscape-scales) at a resolution for which there is no experimental comparison. Success in upscaling in these cases is difficult to test using traditional statistical procedures. Hence, correlation-related procedures, various pattern recognition techniques and, notably, ecosystem models arise as tools for the analysis of these data. Examples of earlier large-scale RS observations in addressing ecological questions of processes are determining how much PAR was being absorbed by a pixel on the terrestrial surface (Tucker 1979; Tucker et al. 2005), the extent of regional wildfires (Justice et al. 1996) and how climate change might be altering the global phenology of vegetation (Nemani et al. 2003). Shugart et al. (2015) saw the fusion of ecological modeling and remote sensing as a necessary synthesis needed to improve our modeling prediction of the ecological responses to global change for forests over regional and continental scales. Because the planet is warming (IPCC 2014), the need for such a capability could not be greater.

The United States National Aeronautics and Space Administration (NASA), the European Space Agency (ESA), and other national space agencies have developed and are beginning to launch a diverse array of new RS instruments. Many of these have already been tested from ground and airborne platforms. They are capable of distinguishing the vertical, horizontal, and 3-D structure of forests with either LiDAR (Light Detecting and Ranging) instruments (Lefsky et al. 2002), or RaDAR (Radio Detecting and Ranging) instruments (Shugart et al. 2010; Hall et al. 2011; Le Toan et al. 2011). Hyperspectral imaging spectroscopy can quantify leaf-and-species-level chemical and functional traits (Asner et al. 2012). At global and ecosystem scales, SIF (Solar-Induced Chlorophyll Florescence) has been demonstrated to be linearly correlated with GPP at the seasonal scale, and thus can potentially serve as an optical proxy for GPP (Frankenberg et al. 2011; Joiner et al. 2014; Yang et al. 2015; Coppo et al. 2017). A wide range of microwave and optical sensors are currently providing global observations of soil and canopy moisture, exchanges of water through evapotranspiration, capturing impacts of droughts on ecosystem function, tree mortality and carbon cycling (Saatchi et al. 2013; Zhou et al. 2014; Eswar et al. 2018; Fan et al. 2019).

The relationship between the RS observations and ecosystem models is also undergoing a paradigm shift. The different national aerospace agencies are currently in a major research phase of calibration and validation of the new instruments, which should provide a continuing string of independent tests for a priori forest model predictions, often made over large areas (Shugart et al. 2018). Many of these involve quantification of the physical structure of forests along with other attributes more usually measured at microscales (Fig. 1).

Methods

Our intent is to identify recent developments in forest gap models and to identify a fortuitous synergism with the capabilities of developing remote sensing technologies to evaluate these models. Since we are using models to compare change prediction at different scales, we divide the results section that follows to each of the ranges of time and space scales that we study. Forest gap models are a class of individual-based forest models (IBMs) that simulate the establishment, growth, and mortality of individual trees on independent plots, or forest patches, i.e. ‘gaps’, about the area of influence of a dominant canopy tree (Shugart 1984; Shugart et al. 2018). They are usually applied to forest of mixed tree sizes and species (more “natural” forests). The first of these models, JABOWA (Botkin et al. 1972) and FORET (Shugart and West 1977), were developed for use in the eastern United States. Follow-ons from these models are numerous and have been developed for use in the forests of China and Russia (Yan and Shugart 2005; Shuman et al. 2017), the western US (Bugmann 2001; Foster et al. 2017), the tropics (Huth and Ditzer 2000; Fischer et al. 2016), Europe (Bugmann and Solomon 2000), and boreal North America (Bonan 1989; Foster et al. 2019). Gap models simulate vegetation-soil interactions (Pastor and Post 1985; Bonan 1989; Foster et al. 2019), wildfire, windthrow and insect outbreak impacts on vegetation (Schumacher et al. 2006; Shuman et al. 2017; Foster et al. 2018), and volatile organic carbon (VOC) emissions from forests (Wang et al. 2017a).

Gap models compute individual tree growth, mortality, and regeneration through a combination of deterministic processes such as species-specific optimal-diameter-increment growth over time and individual-tree growth response to environmental conditions, and stochastic processes such as stress-related mortality, regeneration success, and disturbances. Each simulated plot represents a single, independent forest gap undergoing successional and gap-dynamics processes through time. At a single location or “site”, several hundred of such plots are generally run with similar starting conditions and site-wide parameters. Through the combination of deterministic and stochastic processes, individual plots differ from differences in mortality and regeneration events. Thus, output from a single simulated plot represents a potential outcome arising from the incorporation of these processes and interactions. Simulations typically are produced as Monte Carlo simulations. The average of an array of simulated plots represents the mean expectation of the characteristics of a forested landscape of indeterminate size, with the plots representing a dynamic mosaic of forest gaps, each with its own dynamical history in any given year (Shugart and Seagle 1985). The landscape-scale output from a gap model is similar to a random sampling of an actual landscape using forest inventory plots. Monte Carlo simulations on average produce properties of forest landscapes emerging over time as forest succession, cyclical dynamics, and forest response to shifting climate and disturbance regimes (Shugart and Woodward 2011; Foster et al. 2015; Shuman et al. 2015; Shugart et al. 2018).

Typically, individual-tree growth is simulated annually. Other processes such as soil moisture and decomposition dynamics are simulated at monthly- or daily-time scales. Individual trees differ in their tolerance to the ongoing environmental conditions on each plot based on their size, species, and current growth rate (Shugart 1984). Trees shade one another and compete for resources, and impact the soil conditions on the plot through changes in litter inputs and nutrient requirements (Pastor and Post 1985; Yan and Shugart 2005). Trees may die from prolonged low growth or by disturbances. Generally, disturbances such as fire or windthrow occur at the plot-level and do not spread to other plots within the same site. Regeneration of new trees is dependent on species-specific seed- and seedling banks, modified according to each species’ abundance on the plot, regeneration strategy, and environmental tolerances (Yan and Shugart 2005).

Through explicit simulation of individual trees interacting with one another and their environment, gap models reproduce forest dynamics, compositional change, biomass, and structure at a resolution comparable to forest inventory data across a wide range of ecosystem types. Tabulations of dozens of examples, mostly used as model performance testing are in a sequence of reviews with progressive updates (Shugart 1984, 1998; Shugart et al. 1992; Shugart and Woodward 2011; Larocque et al. 2016). Through simulation at sites spanning large regions or continents, gap models can provide large-scale estimates of forest characteristics and response to environmental change (Shuman et al. 2015, 2017).

Results

Micro-scale (10 m² to 10⁶ m²) models and observations

The scaling up of production and emissions of Volatile Organic Compounds (VOCs) from leaf to ecosystem level needs to confront a challenge of high interspecific variability in the emissions of these compounds. In contrast to the primary metabolic processes of photosynthesis and transpiration, which show shallow phylogenetic conservation, secondary metabolism of VOCs is relatively deeply conserved (Harley et al. 1999; Monson et al. 2013). In other words, while photosynthesis, respiration, and transpiration show variability across species within a system in the order of tens of percent, the variation in VOCs production capability across species in a plant community is often orders of magnitude (Lerdau and Slobodkin 2002). This contrast in heterogeneity between primary and secondary metabolisms is globally true across biomes; both emitter and non-emitters of VOCs co-exist in an ecosystem, and among the emitters the emission capacity varies significantly (Loreto and Fineschi 2014). For forests, from temperate forests in the eastern United States to tropical systems about one-third of tree species produce isoprene; even low diversity ecosystems, such as boreal forests, contain a mixture of emitting and non-emitting species (Lerdau 2007). This inter-specific heterogeneity in VOCs production intrinsically requires both species-level accounting and vertically explicit accounting for variation within forest canopy, as light, operating through both direct effects and indirect effects on leaf temperature, is a crucial factor influencing VOCs emissions. By contrast, considering a model complexity versus efficiency tradeoff, existing VOCs models primarily adopt the scheme of plant functional types that is extensively used in biosphere models (Sellers et al. 1986). Such models, notably MEGAN (Guenther et al. 1995, 2012), have undoubtedly made significant contributions to understanding spatial-temporal patterns and magnitudes of VOCs emissions. However, a representation with plant functional types circumvents ecological complexity (in terms of both functional and structural diversity) in ecosystems. Remaining huge uncertainties in magnitude of VOCs emissions and a lack of prognostic capability arguably suggest that they are far from robust and should be improved.

Current state-of-the-art gap models, e.g. UVAFME, are better positioned to tackle these problems of complexity (Wang et al. 2016, 2017a, 2018). In contrast to models based on plant-functional types, a gap model-based VOCs emission simulator is explicit in considering inter-specific emission variation and variation associated with vertical light change within canopy and tree crown (Fig. 2). Given the established forest gap modelling framework of being functionally and structurally explicit, the development is not challenging, especially for forest ecosystems in the eastern US where species-specific VOCs emission factors are readily available (Fig. 3a and b). With this individual-based gap model, for the first time Tansley’s ecosystem concept is enriched with a new dimension of secondary metabolism in terms of VOCs emissions. Thus, the interplay of ecosystem processes driven by primary and secondary metabolisms over time can be scaled and examined, facilitating a more complete elucidation of ecosystem dynamics and functioning.

Fig. 2

Schematic of light transmission through forest canopy and tree crown in gap model-based forest VOCs emissions model—UVAFME-VOC. Kdt, Kt, and Kb represent the extinction coefficient for diffuse, total beam, and beam, respectively. Qsc and Qabs denote the scattered and absorbed beam radiation, respectively. Reproduced from Wang et al. (2017a)

Fig. 3

Micro-scale simulations of forest dynamics in composition, productivity, and hydrocarbon emissions in Southeastern US using UVAFME-VOC. (a) Contemporary forest dynamics in the Southeastern US. (b) Forests around the 1900s with American Chestnut (pre-blight). The width of each color band represents the biomass of different species (see the color code in Wang et al. 2017a). (c) Unsuppressed forest biomass and (d) enhanced hydrocarbon emissions driven by ozone pressure, which is implemented by imposing species-specific sensitivity on the forest illustrated in panel a. (e) Schematic of diversity-mediated threshold of forest hydrocarbon production under warming, where ecological suppression offsets physiological enhancement. Prior to and beyond the threshold, distinctive feedback mechanisms between climate and atmospheric chemistry mediated by forests may occur. Figures reproduced from Wang et al. (2016, 2017a, 2018)

The relative roles of environmental pressure versus forest community in driving ecosystem functioning in terms of both forest productivity and VOC emissions have been investigated (Fig. 3) with this new micro-scale gap model. These studies clearly indicate the considerable roles of community-level processes in mediating ecosystem responses to environmental pressures (Wang et al. 2016, 2018). For example, ozone, a secondary air pollutant that is harmful to both human health and plant activity, is assumed responsible for global-vegetation-productivity decline and land-carbon-sink reduction (Sitch et al. 2007; Wang et al. 2017b). However, long-term simulations suggest that, over time, ozone does not necessarily dampen forest productivity and carbon stock (Wang et al. 2016). This system-level result emerges from a community-level process of ozone-resistant species replacing sensitive ones (Fig. 3c). Moreover, we found that ozone pressure enhances isoprene emissions by favoring isoprene-emitting species which are less ozone sensitive (Fig. 3d), a result of a potential plant metabolic tradeoff of resource consumption versus stress tolerance. Besides ozone, increasing temperature is another pressing driver faced by forest ecosystems. Contrary to the prevailing opinion of increasing VOCs emissions under climate warming, increasing temperature does not necessarily continuously enhance forest isoprene emissions because of forest compositional changes in relative abundance of emitters versus non-emitters (Fig. 3e). These results strongly support the role of forest community in mediating the forest ecosystem responses to global change agents, pointing to the deficiency of PFT-based models in scaling up physiology directly and to the advantage of forest gap model in incorporating community processes (Wang et al. 2016, 2018).

These initial micro-scale explorations and findings with a forest gap model regarding secondary metabolism-mediated processes warrant more investigations into their performance at differing locations and into their larger scale implications for biosphere-atmosphere interactions (Lerdau 2007; Wang et al. 2019). Achieving such scale-ups of these local scale emissions to larger scales can be inspired by the following discussions on applications to addressing patterns of forest composition and structure at scales from meso- to macro-scale.

Meso-scale (10⁶ m² to 10¹⁰ m²) models and observations

At the meso-scale, gap-models scale from dynamics at the individual plant- and stand-scale to those of whole landscapes. At the stand-scale, these models simulate exogenous and endogenous forest dynamics through explicit tracking of individuals throughout their life cycle, from initial regeneration following release by the death of a canopy dominant, through their growth and response to local-scale weather, site, and environmental conditions, and to their mortality because of low growth or disturbances (Shugart et al. 2018). These individuals additionally interact with one another via shading and competition for other resources as well as through impacts on their abiotic environment (e.g. soil depth, soil moisture, litter quality, etc.). Gap models thus capture how individual trees and forest stands respond to and interact with their changing environment. They are particularly useful in forest ecosystems where the disturbance regime and successional dynamics within the region lead to a heterogeneous landscape of mixed-age, mixed-species stands. Individual trees respond differently to their environment and disturbance events based on their size, age, species, and life history. Gap models, which track all of these variables on an individual-scale, potentially can more closely match actual forest response to change over a wider domain. In particular, gap models can aid in simulating multi-scale interactions between vegetation, disturbances, and climate.

An example landscape-scale gap model application

Foster et al. (2018) utilized the individual-based gap model, UVAFME, to predict interactions between vegetation, bark beetle infestation, climate change, and other disturbances in the subalpine zone of the US Rocky Mountains. Spruce beetles (Dendroctonus rufipennis (Kirby)) are an aggressive bark beetle species that infests spruce (Picea spp.) species throughout the western United States, Alaska, and Canada (Jenkins et al. 2014). Because spruce beetles preferentially attack older and larger spruce trees, the post-outbreak species composition and stand structure tend to be a mix between small (< 10 cm diameter) spruce trees as well as variable-sized non-host tree species (Veblen et al. 1991). This post-outbreak composition and structure is in contrast to the effects of stand-replacing fire in the western US, whereby even-aged, often even-species stands tend to arise.

Spruce beetles infest trees through “mass attacks”, whereby pheromones released by attacking beetles draw more and more beetles to a host tree. This strategy allows for successful infestations of otherwise healthy trees (Raffa et al. 2008). The success of a mass-attack is therefore predicated on the size, age, and condition of the host tree, as well as the local population size of beetles. Beetle population size is additionally impacted by the availability of suitable hosts and climate impacts on survivorship and population growth rate (Berg et al. 2006; Hansen et al. 2011, 2016; Hart et al. 2015). Through explicit simulation of individual trees and individual stands, UVAFME captured the multi-scale factors that influenced infestation rate under current and future climate scenarios (Foster et al. 2018). The infestation probability of a spruce tree was determined by the tree-level characteristics such as size, stress-level, and proximity to other infested trees, stand-level characteristics such as basal area of spruce and down woody debris, and site-wide climate characteristics. With this methodology, the model accurately produced the shifts in species composition and stand structure following an outbreak seen in field studies documenting such events (Veblen et al. 1991; Derderian et al. 2016). This ability was in large part due to the fine-scale nature of UVAFME. Without the representation of vegetation dynamics at their inherent scale – individual trees – the detailed response of forest stands to outbreaks cannot be simulated without strong assumptions about stand structure.

Through a Monte Carlo-style aggregation of several hundred plots, UVAFME additionally represented landscape-level forest properties and dynamics over time in response to spruce beetle outbreaks. At this scale the model produced, without prescription, emergent properties of these interactions, predicting the rising and falling of infestation rates over time (Foster et al. 2018) (Fig. 4). This periodicity is comparable to periodicities in spruce beetle outbreaks found in field studies (Veblen et al. 1994; Zhang et al. 1999; Berg et al. 2006; Hart et al. 2014), and arises as a result of the nature of the vegetation-spruce beetle system. During an outbreak, most if not all of the larger spruce trees are killed, leaving only trees too small to sustain high levels of beetle populations (DeRose and Long 2012). Beetle populations then decline, allowing the surviving trees to grow and eventually become suitable for infestation. Some outbreak-inciting event (e.g. drought, windthrow, etc.) then occurs and starts the cycle anew (O’Connor et al. 2015). This emergent property of tree- and stand-level infestations was not prescribed within UVAFME, yet the model was able to produce it via simulation of hundreds of plots as a mosaic of forest stands across a landscape, each undergoing separate fine-scale vegetation-disturbance dynamics but experiencing similar climate conditions.

Fig. 4

(a) UVAFME-simulated species-specific biomass (t C·ha^− 1) at the USDA Forest Service’s Glacier Lakes Ecosystem Site site in southern Wyoming. Spruce beetle infestations of Engelmann spruce (Picea engelmannii) begin at simulation year 400. Cyclical dynamics can be seen over time in the spruce biomass as a result of cycles of infestation-related mortality. (b) Spruce-beetle killed Engelmann spruce biomass (t C·ha^− 1·yr^− 1) over time for the same simulation averaged across all 200 plots. The red line corresponds to a 10-year running average. Redrawn from Foster et al. (2018)

This type of multi-scale modeling is crucial for capturing vegetation drivers and biotic-abiotic interactions under current and future climate scenarios. In ecosystems where such feedbacks are tied to species- and tree-size specific responses to disturbances, gap models can simulate the potential non-linear and cascading effects of shifting climate and disturbance regimes (Shugart and Woodward 2011; Seidl et al. 2017). Simulations at the stand (500 m²) and landscape scale (10⁶ m²) can then be applied across whole regions or continents to simulate forest change across very large scales (Shuman et al. 2017; Foster et al. 2019).

Scaling-up to regional applications

A recent application with UVAFME in boreal Alaska found that simulating the tight linkage at the tree- and stand-level between vegetation demography, soil characteristics, wildfire, and climate was necessary to represent forest dynamics, structure, and composition across the entire region (Foster et al. 2019). Vegetation interactions along with soil characteristics and the fire regime dominate the Alaskan boreal zone (Viereck et al. 1983; Chapin III et al. 2006b; Johnstone et al. 2010a). In these forests, different stable states of vegetation type arise because of species-specific vegetation-soil-fire interactions. Black spruce (Picea mariana) is able to grow and reproduce on deep, poor nutrient quality, moist soils with shallow permafrost layers (Viereck et al. 1983; Burns and Honkala 1990). The slow decay rate of black spruce litter (Flanagan and Van Cleve 1983; Vance and Chapin III 2001) leads to the buildup of a thick organic-moss layer, a shallow active layer, low nutrient contents in the soil, and the dominance of black spruce over other tree species that do not tolerate such conditions. In contrast, mixed white spruce (Picea glauca) and deciduous stands generally occur on warmer slopes without permafrost. These tree species are less tolerant of deep soils, permafrost, and low nutrients, and the faster decay rate of deciduous litter allows for these species’ more favorable conditions to persist (Johnstone et al. 2010a). Such self-perpetuating ecotypes additionally interact with the fire regime, as black spruce stands are more flammable than deciduous stands, and black spruce depends on fire for rapid post-fire reproduction via serotinous cones (Greene and Johnson 1999).

Foster et al. (2019) updated UVAFME to include daily freeze/thaw and seasonal active layer depth dynamics, which in turn impacted soil moisture dynamics, litter decay, and individual tree growth and reproduction. Litter decayed according to site/soil conditions as well as litter characteristics, which differed based on litter type (e.g. leaves, branches, boles, etc.) and genus (for leaves). Thus, the litter decay rate and litter influx were tied to the species composition, forest structure, and successional cycle of each simulated plot. The decay rate and litter influx in turn impacted the depth and nutrient content of the humus and litter layers, which fed back to soil moisture and permafrost dynamics as well as tree growth and regeneration (Foster et al. 2019). The updated model additionally tied fire intensity to litter content and characteristics, with thicker, drier soils burning at higher intensities. Fire intensity along with forest structure and species composition impacted fire mortality and post-fire regrowth on each stand (Shuman et al. 2017; Foster et al. 2019).

With these updates to UVAFME each individual tree’s growth, mortality, and decay influenced not only the surrounding trees on the plot but also the environmental conditions that those trees experience. Prior to these updates, the model could not accurately reproduce forest characteristics and dynamics within the study region (Foster et al. 2019). Even intermediate testing between model updates (e.g. following permafrost/soil moisture updates but before decomposition updates) did not produce accurate results when compared to inventory data or expected forest successional dynamics. It was only until all biotic-abiotic feedbacks were included in UVAFME that the model could simulate the forest successional dynamics and resulting interactions. Such a failure to simulate forest dynamics accurately within boreal Alaska without the tree-level links between vegetation, soils, wildfire, and climate is unsurprising, given the importance of these interactions in structuring the mosaic of forest types and ages within the region (Chapin III et al. 2006a, 2006b; Johnstone et al. 2010a).

The fine-scale interactions between trees and their environment in this study scaled up to influence region-wide changes in biomass, structure, and composition. Foster et al. (2019) found that species- and tree size-specific interactions drove changes in vegetation and soil conditions under warming temperatures (Fig. 5). Over the course of the climate change simulation, the most important growth stressors shifted from mainly low temperature and shade stress to drought and nutrient stress. Additionally, at the stand-scale, changes in vegetation in response to climate impacted the soil and fire regime – with historically cold, black spruce stands underlain by permafrost shifting to dry deciduous stands with thin soils and a deep active layer. These results indicate that along with drought stress, tree-tree competition for resources and vegetation-soil feedbacks will become increasingly important drivers of vegetation change. When applied across the Tanana River Basin (~ 115,000 km²), UVAFME was also able to predict the differential response of forests across the region to climate change depending on site characteristics and pre-climate change species composition and soil conditions.

Fig. 5

(a) UVAFME-simulated species-specific biomass for a site in interior Alaska under RCP 8.5 climate change scenario (+ 7 °C, + 146 mm precipitation); (b) simulated deciduous fraction (%) for the same simulation; (c) simulated organic layer depth (cm; brown) and active layer depth (cm; gray) for the same simulation; (d) simulated plant-available nitrogen (kg N·ha^− 1) for the same simulation. Red dashed lines indicate the period of climate change. Increasing temperatures result in a shift in the species composition and biomass, impacting the soil regime, which further influenced the changing vegetation. Redrawn from Foster et al. (2019)

Other ecosystem models investigating these biotic-abiotic interactions have been applied within the North American boreal region (Euskirchen et al. 2009; Johnstone et al. 2011; Genet et al. 2013; Fisher et al. 2014; Trugman et al. 2016; Mekonnen et al. 2019). However, most of these applications have represented vegetation at much broader scales in term of composition (groups of plant functional types rather than species) and structure (“big-leaf” canopies rather than individual trees), or have used a state-transition modeling system rather than a process-based simulation of forest dynamics. Given the importance of fine-scale interactions within the boreal region, as well as other forest ecosystem types (Keane et al. 2001; Chapin III et al. 2006a; Araujo and Luoto 2007; Purves and Pacala 2008; Shugart et al. 2018), this lack of fine-scale vegetation representation has implications for predicting both the current and future state of forests, worldwide.

Without representation of individual species and their response and interaction with their environment, much of the detail and potential for non-linear, interacting effects is lost. Fisher et al. (2014) compared simulations of annual carbon flux over Alaska across 40 broad-scale terrestrial biosphere models and found high variability (both in magnitude and in direction) across the model simulations, citing uncertainty in NPP, plant-functional type, and soil conditions as some of the more important factors driving this overall uncertainty (Fisher et al. 2018). Many of these broad-scale models represent boreal vegetation grouped into needle-leaved evergreen (including both black and white spruce), needle-leaved deciduous (larch), and broad-leaved deciduous functional types (Fisher et al. 2018; Mekonnen et al. 2019). The grouping of black and white spruce is problematic given the species’ differential tolerances and resource requirements (Burns and Honkala 1990; Chapin III et al. 2006a), as well as their differential impacts on local-scale soil conditions and the fire regime (Johnstone et al. 2010a). Under climate change scenarios, the compounding and non-linear vegetation, soils, and wildfire responses may alter the existing biotic-abiotic feedbacks (Johnstone et al. 2010b), resulting in new species mixtures and interactions which can only be predicted if species-specific effects are considered at a fine scale.

The representation of fine-scale forest structure is important for considering tree competition, biophysical impacts on climate, as well as how trees of different sizes respond to disturbances and their environment. Tree-tree competition is important for determining forest response to climate change (Purves and Pacala 2008) and implies the use of individual-based models. Tree size impacts disturbance mortality, response to environmental stressors, and the local-scale environment that a tree experiences throughout its life (Shugart 1998; Keane et al. 2001; McDowell and Allen 2015; Hood et al. 2018) (Fig. 4). Models that do not include simulation of a mixed-age, mixed-size forest lose this detail and under scenarios of climate change may not be able to fully capture how forests will respond and interact with their changing environment.

Macro-scale (10¹⁰ m² to 10¹² m²) and mega-scale (> 10¹² m²) models and observations

Many of the applications of gap models, have simulated landscape-scale units. Even the earliest gap models included site variables (soil depth, water capacity of the soil, site quality, etc.) that could vary over the simulated landscape. As was mentioned in the section above, the most straightforward simulations of landscapes are to exercise the model to predict the fates of survey plots the size of a simulated gap-model plot with these “plots” at a spacing that reduces contagion effects to a minimum. These contagion effects can be important and several extensions of gap models have simulated and tested against the effects of contagion. One approach to this problem was manifested as the ZELIG code, in which a gap-model is used as a computational “window” to change each tree in a dynamically changing gridded map (Urban et al. 1991). Windowing is used in other individual-based models in other fields. For example, the simulation of the development of galaxy geometry arising from billions of gravitational interacting stars is solved computing the effect of nearby stars on each single star. This effectively converts the computation of gravitational effect from a function of the square of the number of stars (when all stars interact with one another) to a function of number of stars.

Macro-scale (regional-scale) applications

Forest structure and composition at the meso- and macro-scale impact biophysical feedbacks to climate through changes in albedo, surface roughness, and latent heat flux (Bonan et al. 1992; Liu et al. 2006). Thus, changes in structure and composition have the capacity to impact the trajectory of climate regionally and even globally. Historically, gap models were excluded from use in global climate models due to limits in computer processing power (Shugart et al. 2018). Gap models have a rich and early history of simulating the dynamic changes in landscapes and regions during the Quaternary (Solomon et al. 1980, 1981; Solomon 1986). However, modern technology now allows these fine-scale models to be applied continentally (Sykes and Prentice 1996; Shuman et al. 2017) and provides the impetus for continental and global-scale applications that couple vegetation and climate.

Weishampel et al. (1992) included side-shading effects on canopy geometry in natural Douglas-fir (Pseudotsuga menziesii) of different ages. They tested their prediction against semivariance patterns in canopy heterogeneity that were independently collected using high-resolution photogrammetry. ZELIG has subsequently been used to predict expected forest patterns over time and spatial resolution (van Tongeren and Prentice 1986; Huston and Smith 1987; Smith and Urban 1988; Huston et al. 1988) and in several applications in Boreal forests (Larocque et al. 2006, 2011). Recent simulations for the Russian boreal forest with the SIBBORK model (Brazhnik and Shugart 2015) have simulated contagious effects including orographic shading — a south-facing slope in a deep valley is compositionally and structurally different from a south-facing slope without shading without the north-facing opposite valley side (Brazhnik and Shugart 2016). The SIBBORK model has also been applied to one of the principal contagion effects in boreal forests, wildfire (Brazhnik et al. 2017). Generally, the inclusion of spatial effects in individual-based models has been limited by data for parameter estimation and not by modeling limitations.

Mega-scale (global-scale) applications

Canopy process models simulate the flux of heat, CO₂ and water between plant canopies and their environment over time scales of a few seconds to a day. They represent the canopy as a single- or multi-layer unit with a fixed structure (i.e. leaf area). Photosynthesis and transpiration are simulated by estimating microclimatic variation and stomatal conductance for the canopy (or canopy layers). The inside of a leaf and its external environment are coupled by the stomata, microscopic pores in the leaf surface that can change under different environmental conditions and are controlled by the plant. If these stomatal openings are relatively large, resistance to molecular diffusion is low — H₂O diffuses out from the moist spaces inside the leaf and CO₂ diffuses into the same internal spaces to compensate for the CO₂ taken up by the plant through photosynthesis. As the stomatal aperture closes, the resistance to these diffusion-based transfers increases. Thus, the balance of CO₂ and H₂O are interwoven. Because the loss of water from the leaf is evaporative, the outward flux of water also removes heat from the leaf. In canopy process models, formulae relating the CO₂, H₂O and energy fluxes for leaves form bases for simulating the CO₂, H₂O and energy fluxes of plant canopies over large areas.

For example, Woodward (1987) developed a simple model of the energy and hydrological balance of a plant canopy using the Penman-Monteith equation (Penman 1948; Monteith 1981) to determine canopy transpiration. The water transpired by plant canopy of a given leaf area, along with evaporated water, is subtracted from the water held in the soil. If over the course of a year the soil dries too much or the soil is unable to recharge its water, then the leaf area is assumed too high and a lower value for leaf area used until the maximum sustainable leaf area is found. Provided with spatial data for the environmental conditions needed by the model (solar radiation, precipitation, temperature, etc.), the model simulates the expected leaf area at regional and global scales.

Woodward’s model illustrates several of the general features found in homogeneous landscape models:

• “An Appeal to Optimality” as seen in the procedure used to determine the expected leaf area at a location — that the vegetation will optimize its leaf area to use the available water in a region. This optimum balances the positive tendency for vegetation to add leaf area when water is available with the constraint that, if vegetation “mines” the soil water beyond the resupply rate, then plant death and leaf shedding will reduce leaf area. More recently, it has become apparent that this approach to optimization incorrectly assumes that optimization at the vegetation-assemblage scale is evolutionarily stable; in fact, canopy optimization is only stable when it occurs at the individual tree scale, which leads to an assemblage-scale canopy-cover that is lower than possible when optimization occurs at assemblage-scale (Anten and Hirose 2001). Another constraint is that vegetation cannot gain additional net leaf area once the lowest layer of leaves is sufficiently shaded.

• “Use of Limiting Factors,” particularly for light and water as in Woodward example. The heat balance of the vegetation links to the water evaporated from the vegetation by the Penman-Monteith equation.

• “An Expectation of Generality” allowing the models to be applied to vegetation worldwide. This is often based on physically-based heat-flux equations solved at equilibrium (“Equilibrium Seeking Behaviour”).

Inasmuch as gap models have included metabolism (Friend et al. 1993), they have tended to focus on photosynthesis, respiration, and transpiration because of the central role these processes play in the mass and energy balance of an ecosystem.

Simulation of global carbon fluxes from natural systems has traditionally involved models that scale up ecological processes from study sites to large, assumed homogeneous units of about ~ 0.1° to ~ 0.5° latitude × longitude blocks. Many of these aggregated models such as TEM (McGuire et al. 1992), CENTURY (Parton et al. 1993), CASA (Potter et al. 1993), IBIS (Foley, 1994), SiB2 (Sellers et al. 1996), LPJ (Sitch et al. 2003) and ORCHIDEE (Krinner et al. 2005) have been used to project ecosystem carbon fluxes. There are similar models that use stand-based dynamics or even individual-based models to include population dynamic processes such as plant establishment, mortality and resource competition, notably: Hybrid (Friend et al. 1997), LoTEC (King et al. 1997), LPJ-GUESS (Smith et al. 2008), TRIPLEX (Peng et al. 2002) and INTCARB (Song and Woodcock 2003).

A direct application of a gap model to simulate global forests is the use of the FORCCHN model to evaluate the global forest carbon fluxes (Ma et al. 2017). The FORCCHN model has several significant aspects. It can produce daily estimates of gross primary production (GPP), ecosystem respiration (ER) and net ecosystem production (NEP), using canopy modeling approaches we discussed in the section on micro-scale models (above). This allows the physiological aspect of the model to be tested against 37 forest eddy-covariance sites, which were drawn from the daily data in the LaThuile FluxNet free use data set (http://www.fluxdata.org), AmeriFlux (http://ameriflux.ornl.gov), CarboEuroFlux (http://www.carboeurope.org), ChinaFlux (http://www.chinaflux.org) and FFPRI FluxNet (http://www2.ffpri.affrc.go.jp/labs/flux/index.html). Figure 6 shows four examples using monthly data from four of these study sites against the FORCCHN model. Across all 37 sites, the daily correlation coefficients averaged 0.72 for GPP, 0.70 for ER and 0.53 for NEP.

Fig. 6

Monthly variation of observed and simulated (a) GPP, (b) ER and (c) NEP at four forest eddy correlation sites: (from top to bottom) Changbaishan site in China (CN-Cha), Loobos site in Netherlands (NL-Loo), Collelongo-Selva Piana site in Italy (IT-Col) and Tumbarumba site in Australia (AU-Tum). The grey histogram denotes observation and black dot denotes simulation results

After inspecting the ability of the FORCCHN model to reproduce GPP, NEP and ER for the 37 eddy correlation sites, the model was then applied to the problem of estimation of global carbon fluxes from forests. To do so one must develop a data set to drive the FORCCHN model:

• 1) The climatological forcing was from Princeton University over the 1982~2011 period at a grid resolution of 0.5° × 0.5° (http://hydrology.princeton.edu, Sheffield et al. 2006). Derived through a combination of reanalysis data and observations, these variables include the daily maximum and minimum air temperature (°C), precipitation (mm), relative humidity (%), wind speed (m·s^− 1), atmospheric pressure (hPa) and total solar radiation (W·m^− 2).

• 2) Soil parameters are soil organic matter (carbon and nitrogen pool in units of kg C·m^− 2 and kg N·m^− 2, respectively), soil physical parameters and litter pool decomposition parameters. The soil physical parameters are strongly dependent on the geographical position and include the soil field capacity (mm), wilting point (mm), bulk density (kg·m^− 3), sand content (%), silt content (%) and clay content (%). The Global Gridded Surfaces of Selected Soil Characteristics (Global Soil Data Task Group 2000) coupled with Harmonized World Soil Database (Nachtergaele et al. 2012) provide resources for the soil organic matter and physical parameters. The litter pool decomposition parameters are calculated according to Kirschbaum and Paul (2002). Moreover, to insure that the allocation proportion of organic carbon in ten soil pools are in equilibrium in FORCCHN, the model is run for 300 years at each grid point and then the new allocation proportions are used as model input data to simulate C fluxes for the past 30 years.

• 3) Global forest types, used to select which tree functional types are used, are derived from International Geosphere Biosphere Program-Data and Information Service (IGBP-DIS) DISCover land cover classification system, with a spatial resolution of 0.5° × 0.5° (Loveland et al. 2009). The 8-day 5-km LAI of Global LAnd Surface Satellite (GLASS) in 1982 (Liang et al. 2013) is also used to drive the model. Product developers ascribe quality control flags based on LAI to screen and reject poor quality data. The 8-day LAI are composited into the yearly maximum and minimum values. Note that satellite-derived LAI datasets are resampled to the geographic projection and spatial resolution of the global climatological forcing.

The resultant global simulations for forest GPP and ER (Fig. 7) are consistent with Model Tree Ensemble-based GPP estimates (Jung et al. 2011) — except that FORCCHN-derived GPP is about ~ 300 g C·m^− 2·yr^− 1 smaller in most tropical rain forest and ~ 900 g C·m^− 2·yr^− 1 larger in parts of south-central Africa. Over the simulated interval (1982–2011) both GPP and ER significantly increased (P < 0.01, see Ma et al. 2017) across all forest types. The forest type with the greatest CO₂ uptake by photosynthesis and CO₂ release by ecosystem respiration were evergreen broadleaf forests with multi-year averaged values for GPP of 2631 ± 233 g C·m^− 2·yr^− 1 (mean ± 1 standard deviation) and ER of 2513 ± 216 g C·m^− 2·yr^− 1. Deciduous broadleaf forest (GPP of 1428 ± 183 g C·m^− 2·yr^− 1; ER of 1346 ± 184 g C·m^− 2·yr^− 1) and mixed forest (GPP of 961 ± 84 g C·m^− 2·yr^− 1; ER of 917 ± 84 g C·m^− 2·yr^− 1) were the next two most important types of forests.

Fig. 7

Spatial distribution of mean (a) GPP, (b) ER and (c) NEP (g C·m^− 2·yr^− 1) for global forest ecosystems during 1982 to 2011 as simulated by the FORCCHN model

The role of forests as sources or sinks in global carbon budgets is a consequence of NEP. NEP is the difference between two relatively large numbers, GPP-ER. For global forest ecosystems FORCCHN gave an annual total GPP and ER of 58.83 ± 5.61 and 55.77 ± 5.18 Pg C·yr^− 1 for global forest ecosystems during 1982–2011, such value is within the range reported by other GPP models. Global forest ecosystems as simulated by FORCCHN contribute a substantial C sink for the same period, with total NEP being 3.06 ± 0.67 Pg C·yr^− 1. This is also comparable to the results from studies using observation-base estimations of 52.61–67.54 Pg C·yr^− 1 (Beer et al. 2010) and from satellite-based observations of 37.59–59.77 Pg C·yr^− 1 (Cai et al. 2014). These are initial global results from a globally distributed forest gap model. Ma et al. (2017) mention several improvements that would be invaluable in the model and in the available data. Nonetheless, the agreement between the FORCCHN results and predictions from other global models using different assumptions are robust convergence of outcomes.

Conclusions

In this paper, we have focused on individual-based gap models of forests across the time and space domains of their applications. These models draw strongly from information on the silvics, allometry and environmental responses of individual tree species. These descriptive parameters are usually estimated from descriptions of the tree species and not fitted to data. For example, a parameter for tree maximum height for a species might have some variability in its estimate, but it could not be estimated to be 500 m, even if this produced better statistical fits to a data set. Other significant parts of the models include, for example, standard biophysical models for evapotranspiration of biogeochemical models for nutrients released by decomposition. Here again, the parameters are not assigned arbitrary values for best statistical fit. What these features mean is that these models are more likely to be tested for agreement with data, rather than fit to data.

A new generation of satellite observation systems have the capability to provide new observations to test existing model predictions. This class of models can efficiently interact in a hypothesis–testing mode in at least three ways:

• 1) One can use the models for predicting the expected patterns of variability among important ecological variables as simulated by the models. These predictions include a priori estimation of expected composition of forests (Asner et al. 2012), remote sensing observations implying physical structure of forests and their relationships to biomass (Köhler and Huth 2010; Le Toan et al. 2011; Saatchi et al. 2011; Lobo and Dalling 2014), productivity (Yang et al. 2015) and VOCs (Fu et al. 2019).

• 2) One can use remote sensing and gap models to predict how ecosystem processes produce patterns across scales from micro to global. While previous field studies typically needed to specify a spatial scale of interest prior to sampling forests, the fine-scale and global coverage of remote sensing products, along with the available computing power for model simulations, enables the simulation of fine-scale processes across large spatial extents, encompassing all of the intermediate scales along that spectrum. While ecologists have long been interested in the relevance of scale on ecological processes, these approaches now provide a clear means to explicitly evaluate the impact of scale on patterns, such as suggested by hierarchy theory (O’Neill 1989).

• 3) Our microscale example shows how ecological processes may also affect multiple scales. While the canopy positions of tree species at a gap scale affects the competitive growth rates of individual trees, it also contributes to the regional production of ozone at larger scales, which in turn impacts fine-scale competitive interactions by favoring ozone-resistant species. Such scale-spanning interactions are also inherent in studies involving forests and climate change as we discuss in our global example. While the uptake of carbon dioxide through photosynthesis is ultimately dependent on leaf-level physiology, the global scale consequences of that uptake influence climate depending on whether forests serve as a carbon source or sink. The global summation of that physiological process returns to affect the competitive interactions between tree species at fine scales, as species’ competitive ability will change along with a changing climate. Remote sensing and model data will be needed to test predictions on how these processes impact forest ecosystems across scales.