Assessing climate change impacts, benefits of mitigation, and uncertainties on major global forest regions under multiple socioeconomic and emissions scenarios

Our study uses an ensemble of climate change projections simulated by the MIT Integrated Global System Model-Community Atmosphere Model (IGSM-CAM) modeling framework (Monier et al 2013a). The climate ensemble is composed of different emissions scenarios (unconstrained versus stabilized radiative forcing), different global climate system response (climate sensitivity), and different realizations of natural variability (initial conditions) (table 1). The ensemble was prepared for the US Environmental Protection Agency's Climate Change Impacts and Risk Analysis (CIRA) project (Waldhoff et al 2015), which examines the benefits of global mitigation actions to the United States. In the core analysis of our study, we consider two emissions scenarios: a reference scenario (REF) that represents unconstrained emissions similar to the Representative Concentration Pathway RCP8.5 scenario (Riahi et al 2011), and a greenhouse gas (GHG) mitigation scenario (POL3.7) that stabilizes radiative forcing at 3.7 W m⁻² by 2100. The POL3.7 scenario was designed to fall between the RCP2.6 (van Vuuren et al 2011) and RCP4.5 (Thomson et al 2011) scenarios and is consistent with a 2 °C global mean warming from pre-industrial by 2100. We focus our analysis on the simulations with a climate sensitivity of 3 °C, and for each emissions scenario we use a five-member ensemble with different initial conditions, thus limiting the total number of climate simulations to ten and keeping our core analysis to a manageable number of MC2 simulations. We analyze the mean over the different initial conditions in order to obtain robust estimates of the anthropogenic signal and filter out the noise from natural variability, and identify the changes due to GHG mitigation. This approach allows us to account for the significant uncertainty in natural variability, highlighted in a number of studies (Hawkins and Sutton 2009, Deser et al 2012, Monier et al 2013b, 2015, 2016).

We further expand upon our uncertainty analysis by analyzing the range of climate impacts on the world's forests associated with different global climate system responses, analyzing simulations with a climate sensitivity of 2.0 °C and 4.5 °C, obtained via radiative cloud adjustment (see Sokolov and Monier 2012). We also analyze a slightly less stringent GHG scenario (POL4.5), similar to a RCP4.5.

Although the climate ensemble used in this study is derived using a single climate model, it accounts for the uncertainty in the emissions scenarios, the global climate response, and natural variability, which account for a substantial share of the uncertainty in future climate projections (Monier et al 2015). More details on the emissions scenarios can be found in Paltsev et al (2015), details on the climate projections for the US can be found in Monier et al (2015), and the implication for future changes in extreme events is given in Monier and Gao (2015).

Dynamic global vegetation models (DGVM) simulate terrestrial biosphere's response to climate by modeling vegetation biogeography, vegetation dynamics, biogeochemistry, and biophysics (Fisher et al 2014). DGVMs are the best tools for representing vegetation dynamics at global scales (Quillet et al 2010), and have been used by many to study vegetation dynamics at global scales (e.g. Cramer et al 2001, Friedlingstein et al 2006, Sitch et al 2008). MC1 DGVM (Bachelet et al 2001) has been applied at many scales, including continental and global scales (e.g. Bachelet et al 2015, Beach et al 2015, Drapek et al 2015, Gonzalez et al 2010). MC2 is MC1 DGVM re-written in C++ to improve computing speed and code organization. The design of MC1 and MC2 is comparable in complexity with other DGVMs (Fisher et al 2014, Quillet et al 2010). MC2 design is detailed elsewhere (Bachelet et al 2001, Conklin et al 2016), thus we highlight only the essential features and limitations of MC2 here.

MC2 represents land surface as a grid. It reads as input elevation, soil, and climate data and runs on a monthly time step. In each grid cell, the terrestrial ecosystem is represented as a web of above- and below-ground carbon pools. Plant growth, carbon and water fluxes are calculated monthly, using CENTURY Soil Organic Matter Model (Parton 1996) as a submodel. In each grid cell, a representative tree and grass compete for light and water. Monthly temperature, precipitation, and vapor pressure data drive calculations of plant productivity, decomposition, and soil water balance. Net primary productivity (NPP) is calculated directly as a function of temperature and available soil water. MC2 identifies the representative tree annually using a set of biogeography rules, recognizing a total of 35 plant functional types (table S2). Simulations require extensive computing resources, and are run on a high performance parallel computing platform. MC2 simulates CO₂ effects on NPP and potential evapotranspiration (PET) as simple multipliers, which vary linearly from 350 and 700 ppm of CO₂ concentrations.

As noted in the Introduction, we recognize that forests will be shaped by a complex interaction among multiple driving forces, including an array of disturbance regimes, including land cover change, logging, fire, and insect and pathogen outbreaks. Simulating all types of disturbances is ideal, but it remains a goal yet to be achieved by the earth system modeling community. Although climate change may significantly alter disease and insect outbreaks in forests (Dale et al 2001), specific mechanisms of how forest, insects, and fire interact are poorly understood (e.g. Andrus et al 2016, Harvey et al 2014). More importantly, insect and pathogen outbreaks are highly variable depending on organisms involved and are difficult to model on a global scale. We are not aware of any DGVM implementation that has successfully simulated interactions among forest, insects, and fire on a global scale. In this study, we focus on evaluating the role of wildfire as a major disturbance regime. MC2 is able to simulate grazing effects on grass, but it was disabled for this study. This is a limitation common to all DGVMs (Fisher et al 2014, Quillet et al 2010).

MC2 simulates fire occurrence as a function of the current vegetation type and fuel conditions. MC2 calculates fire consumption of vegetation and the associated ecosystem carbon pools based on the current weather. Conklin et al (2016) provide a detailed description of the MC2 fire module.

2.2.1. Model protocol and calibration

2.2.2. Model validation

For model validation, we ran MC2 for the full 0.5° global grid from 1901 to 2012, and compared the output with benchmark datasets resampled to the 0.5° grid. This represents two-fold cross-validation (Jopp et al 2011, Kleijnen 2008), appropriate when computational costs are heavy (Schwartz 2008). Since MC2 was modified to stochastically simulate the occurrence of fires, we ran 12 replicates for the 1901–2012 period, and calculated the mean and mode statistics of output variables. The comparisons with benchmark data were tabulated for the 16 major forest regions of the world (Sohngen et al 2001, Sohngen 2014) (figure 1). We compared MC2's estimates of NPP for 2000–2011 with MODIS NPP MODIS Terrestrial Gross and Net Primary Production Global Data Set, version MOD17 (Zhao et al 2005). We compared MC2's projection of locations of vegetation biomes with ISLSCP II Potential Natural Vegetation Cover (Ramankutty and Foley 2010), as well as the GLC2000 global land cover dataset (Bartholomé and Belward 2005). The land cover types in each dataset were translated to the biome types MC2 simulates (desert, shrubland, grassland, woodland, and forest) and Cohen's kappa was calculated in comparison to MC2 output. ISLSCP II did not distinguish between woodland and forest, so for that comparison the two MC2 biomes were combined. Finally, we compared MC2 estimates of burned area with GFED4 data (Giglio et al 2013). All comparisons, except the comparison with ISLSCP II, were made after agriculture and developed areas shown in the GLC2000 land cover data were masked out. For comparing NPP and burned area, MC2 outputs (g m⁻² y⁻¹ and %, respectively) were multiplied by the area of each grid cell. The validation results are described in the Results section below.

Zoom In Zoom Out Reset image size

To run simulations with future climate realizations, we compared 30 yr averages of total live vegetation carbon stock (C_veg) from the 12 replicates of the full grid, 1901–2012 MC2 simulations, and selected the replicate with C_veg the most similar to the ensemble average value of C_veg. We used the 1979 state of the selected simulation as the starting state for the future simulations. We first ran IGSM-CAM to produce global climate realizations, then used the climate realizations to drive MC2. The IGSM-CAM climate realizations were downscaled from their native resolution to 0.5° degree resolution using the delta method (Fowler et al 2007). For each climate realization, we ran 10 replicates of MC2, and, as for the validation analysis, agriculture and developed areas were masked out for analysis.

For a recent historical period (1983–2012) the global proportions of biomes projected by MC2 were comparable to the proportions derived from GLC2000 and ISLSCP II datasets (figure 2). MC2 projected 52% of land grid cells to be forest and 6% to be woodland, while the proportions derived from GLC2000 was 34% for forest and 5% for woodland. The proportion of woodland and forest combined projected by MC2 was 58%, while the combined proportion derived from ISLSCP II was 43%. The proportions for grassland and desert varied widely among all three data sources, while the proportions for shrubland were within 3 percentage points. Cohen's kappa for the globe between MC2 and ISLSCP II was 0.46, and with GLC2000 it was 0.43. Global average annual NPP simulated by MC2 for 2000–2011 was only 0.99 Pg (2%) over the MODIS MOD17 estimate, while the area burned simulated by MC2 for 1996–2011 was 109 Mha (31%) below the value reported by GFED4.

Zoom In Zoom Out Reset image size

We compared MC2 simulated global net biome production (NBP) and NPP with the values generated by 10 DGVMs included in Piao et al (2013) model inter-comparison study. The average NBP simulated by MC2 for the 1980–2009 period was 1.7 PgC y⁻¹, in the middle of the range of values generated by the ten DGVMs, and falls within the residual land sink value range reported by Friedlingstein et al (2010). The average NPP simulated by MC2 for the 1980–2009 period was 56.2 PgC y⁻¹. MC2 does not calculate gross primary productivity (GPP) published in Piao et al (2013), but assuming NPP is 50% of GPP (Waring et al 1998), MC2's GPP value falls within the range of values generated by the 10 DGVMs.

A region-by-region comparison of MC2 output for biome biogeography, NPP, and burned area with benchmark datasets is tabulated in table 2. The levels of agreement of biome biogeography with benchmark datasets were similar to the level of agreement for the globe, with kappa values for the majority of the regions ranging between 0.39 and 0.53. Kappa was 1.0 in Korea-Taiwan, because both MC2 and ISLSCP II projected 100% forest. Eastern Europe and Australia-New Zealand regions had particularly low agreement, with kappa values as low as 0.09 and 0.12. The average annual NPP for 2000–2011 simulated by MC2 were within 25% of the NPP reported by the MODIS MOD17 dataset for 10 of the 16 world forest regions. In North Africa and Australia-New Zealand, both of which contain much arid and semi-arid biomes, NPP simulated by MC2 far exceeded the values reported by MODIS MOD17 dataset. Area burned simulated by MC2 was within 33% of GFED4 values for five of the 16 regions: Western Europe, Russia, Korea-Taiwan, South Asia, and Australia-New Zealand. For another seven regions, the values deviated less than 100% from GFED4 values. For the remaining four regions—USA, Central America, India and China—the burned area simulated by MC2 deviated over 100% from GFED4. The worst of these was India, where burned area was 12.4 Mha (803%) over the value obtained from GFED4.

For each climate change scenario (i.e. REF and POL3.7), we calculated the probability of forest gain as the percentage of simulation replicates that indicated a biome conversion from non-forest to forest from a recent historical period (1980–2009) to the end of the century (2070–2099). For this analysis we aggregated all the forest vegetation types simulated by MC2 into a single 'forest biome' type. We calculated the probability of forest loss in a similar way, by calculating the percentage of simulation replicates that indicated conversion from forest to non-forest. A map of the combined percentages for each climate change scenario is shown in figure 1. The biogeography of forest biomes projected by MC2 was stable across most of the globe. Where there were shifts, MC2's projections were highly consistent across the multiple realizations and replicates within a scenario. That is, the areas showing forest gain and forest loss generally had high and low percentage values, with only a limited number of grid cells showing intermediate levels of agreement.

For both climate change scenarios, MC2 simulated poleward migration of forest biomes, where, in the leading-edge of the migration, lower-productivity biomes (e.g. grassland and woodlands) convert to forests under climate change; while at the trailing edge, forests convert to another biome (e.g. shrubland, grassland, or woodland) due to lower productivity or frequent fires. Large expanses of boreal forests in Canada and Russia shifted northward under the REF scenario, and to lesser degrees under the POL3.7 scenario. In the southern hemisphere, forest expanded southward in Southern Africa. Poleward migration of forests was not distinct in Western South America, where there the forest contracted along elevation gradients. In Australia, the tropical forests in the north contracted northward as they lost productivity and became woodlands. Simultaneously, MC2 simulated increased growth of trees in the woodlands in western Australia, converting those areas to forests.

Although MC2 simulated forest biomes to be geographically stable across much of the globe, the total live forest carbon stock increased dramatically and consistently under both climate change scenarios, gaining 447 Pg C (59%) and 410 Pg C (54%) under the REF and POL3.7 scenarios, respectively. MC2 simulated the vast majority of the total live forest carbon gain to occur in the southern hemisphere: Western South America, South America, and South Asia (figure 3(a)). Although the other regions—with the exception of Russia—also gained total live forest carbon, the amount gained were less than 20 Pg. Russia lost as much as 17 Pg (23%) of the total live forest carbon under the REF scenario, due to conversion of forest biomes to non-forest biomes. For Asia and Europe only small increases were projected, while Russia was projected to see a significant decline.

Zoom In Zoom Out Reset image size

To analyze the global impact of GHG mitigation on the world's forests, we show the range of changes among the 5-m ember initial condition ensemble for each emissions scenario for important metrics of global forest conditions: changes in forest carbon, NPP, forest carbon consumed by fire, forest area and burned area (figure 4). The analysis reveals large increases in forest carbon under both scenarios, along with increases in NPP, increases in wildfire, burned area and forest carbon consumed by wildfire, along with decreases in forest areas. The magnitude of the climate change effects are reduced by GHG mitigation under the POL3.7 scenario compared to the REF scenario, both the positive effects (increases in carbon stocks and NPP) and the negative (increase in wildfire and decreases in forest areas).

Zoom In Zoom Out Reset image size

Our analysis also estimates the uncertainty associated with natural variability, and thus provides a basic signal-to-noise analysis to test the robustness of the impact of climate mitigation. The role of natural variability on changes in forest carbon and NPP is small (figure 4(b)), but it is substantial for changes in forest carbon consumed by wildfire (figure 4(c)) and burned area (figure 4(e))—to the point where the ranges of the two emissions scenarios overlap—and to a lesser degree for changes in forest area (figure 4(d)). These results highlight the importance of relying on an ensemble of climate simulations with perturbed initial conditions to quantify the noise associated with natural variability and identify the robust impacts of climate policy. However, while identifying the aggregated impact of climate mitigation provides useful information for decision-making, it does not capture potentially heterogeneous responses at the regional scale. For this reason, we analyze the regional drivers of change next.

The exposure of the regions of the globe to climate change, as represented by change in mean annual temperature and precipitation from 1980–2009 to 2070–2099, are muted under the POL3.7 scenario but vary range widely under the REF scenario (figure 5(a)). For the majority of the 16 global forest regions, the REF scenario projects a temperature rise >3.5 °C and precipitation increase >12%. All of the major temperate and boreal forest regions (USA, Canada, and Russia) are exposed to a >5 °C warming under the REF scenario. For a small set of regions—Australia-New Zealand, Western Europe and Eastern Europe—the REF scenario projects the same magnitude of warming but only a small increase in precipitation (2%–4%). For Central America, the REF scenario actually projects a small reduction (−2%) in precipitation.

Zoom In Zoom Out Reset image size

Forest responses to those climate change exposures are projected to vary regionally, without simple correlations to the intensity of those exposures. As noted above, the greatest increase in forest carbon stock are projected to take place in the southern hemisphere, even though the greatest exposure to climate change is projected for North Africa, USA, Canada and Russia (figure 5). Russia is, however, projected to undergo the largest contraction of forest area (210 Mha) under the REF scenario (figure 3(b)), with 14 Mha lost to fire (figure 3(c)). Significant contraction of forest area is also projected for China (83 Mha, figure 3(b)), but area burned by fire is projected to decrease for China (figure 3(c)). See figure S3 for a complete set of regional change projections.

Benefits of mitigation are also not evenly distributed across the global regions. For example, for Canada the POL3.7 scenario results in a net gain of 3 Pg of forest carbon stock when compared to REF scenario (figure 3(a)), while for USA it results in a net loss of 2 Pg. Nearly half of the regions are projected to benefit from mitigation through increases in forest area (figure 3(b)) and reduction in forest area burned by fire (figure 3(c)). The greatest cost of mitigation—that is, a negative impact on forests—is projected to be borne by the southern hemisphere regions (W. South America, South America, South Africa and South Asia) where the greatest carbon gains are projected under both POL3.7 and REF scenarios.

The drivers of forest changes also vary by region. Different regions are projected to experience changes in forest carbon stock of similar magnitude but associated with differing mechanisms: 1) expansion or contraction of forests, with further loss of acreage to wildfire; and 2) changes in vegetation productivity. Plotted as a two-dimensional grid, these mechanisms have different levels of importance for the world's forest regions (figure 5(b)). The large increases in forest carbon stock projected for the southern hemisphere regions are driven primarily by increases in NPP, with little changes projected to the forest area or area burned. In contrast, the large decreases in forest carbon stock projected for Russia, and, to a lesser extent, for China, are both driven primarily by forest contraction, with only small changes projected in forest productivity. For USA, the increase in forest carbon stock is driven by a combination of forest expansion and increase in productivity. For Canada, forest contraction is balanced by an increase in productivity.

To frame the uncertainty in our estimate of climate impacts on the world's forests, we examine the range of impacts using three different ensembles: the range over 5 initial conditions (for REF and climate sensitivity 3.0 °C), the range over 3 climate sensitivities, namely 2.0 °C, 3.0 °C and 4.5 °C (for REF and initial condition member #1) and the range over the three emissions scenarios (for climate sensitivity 3.0 °C and initial condition member #1) (figure 6).

Zoom In Zoom Out Reset image size

This analysis identifies two major findings. First, increased levels of climate change along the emissions scenario dimension are associated with a larger increase in total forest carbon, but along the climate sensitivity dimension it is associated with a smaller increase in total forest carbon. The major difference between these two dimensions of climate change is the role of CO₂ fertilization. Under the REF scenario, CO₂ concentrations (827 ppm by 2100) are substantially higher than under the POL3.7 (462 ppm by 2100), and therefore so is the CO₂ fertilization effect. Meanwhile, the simulations with different climate sensitivities have the same CO₂ concentrations, which allows distinguishing the role of climate change versus the role of increases in CO₂ concentrations and CO₂ fertilization. This analysis shows that, at the global mean level, forest productivity under the REF scenario benefits substantially from the CO₂ fertilization effect and that higher warming alone does not necessarily increase global forest carbon. Higher levels of climate change, under fixed CO₂ concentrations, have a negative impact on global forest carbon, likely caused by more wildfires, and climatic effects like droughts.

Second, the analysis shows there are substantial uncertainties associated with our estimates of the benefits of GHG mitigation on the world's forests, highlighted by the large range of outcomes between different levels of global climate system response (i.e. climate sensitivity) and different representations of natural variability (i.e. initial conditions). The role of natural variability is even larger at the regional level, as shown in figure 3 (and figure S3), which shows that the range of outcomes between REF and POL3.7 can overlap when the range over different initial conditions is taken into account. This finding further highlights the need to account for natural variability when trying to obtain robust estimates of the impact of climate mitigation on forests, at both the global and regional scale.

Confidence in model projections can only be founded on an objective evaluation of model skill, as its ability to reasonably simulate past conditions is a necessary, though not sufficient, requirement for simulating future conditions. Comparing MC2 output with empirically obtained datasets requires some caution, because MC2 simulates potential natural vegetation without simulating the effects of various anthropogenic effects on the landscape. We evaluated our calibration of MC2 DGVM by analyzing output variables from each of MC2's three main internal modules: NPP for the biogeochemistry module, burned area for the fire module, and biome biogeography for the biogeography module (table 2). MC2's global NPP output compared closely with MODIS MOD17 NPP (Zhao et al 2005), and the global estimates are within the wide range of values reported in literature (Field et al 1998, Kicklighter et al 1999). For 10 of the 16 world forest regions considered, MC2's NPP values were comparable to MODIS values. With the exception of Russia, the regions where MC2's NPP values compared poorly were regions with smaller timber production. The problematic areas, however, highlight the many challenges remaining in DGVM design (Quillet et al 2010). Although there is broad agreement among the models, large uncertainties remain across models (Friedlingstein et al 2006, Piao et al 2013, Sitch et al 2008).

Our reformulation of the fire algorithm and its calibration appears to underestimate the prevalence of fire globally, although the geography of fire is comparable to previous versions (Gonzalez et al 2010). In the key temperate forests of Canada and USA, MC2 appears to overestimate the prevalence of fire (table 2), which may lead to an underestimation of forest C stock. In key tropical forests of Western South America, South America, South Africa and South Asia, MC2 both over- and underestimates fire prevalence reported by GFED4, and the mixture of errors may balance each other.

The agreement between MC2 biome biogeography and the two benchmark datasets may be called 'fair to good' (Banerjee et al 1999) for the globe and a majority of the regions. For comparison, kappa between ISLSCP II and GLC2000 was 0.51. Some of the disagreements arise from the disparate systems used to classify land surface, from aggregation to the 0.5° grid, and from translation to MC2 biomes.

MC2 projects northward migration of boreal forests in Canada, Russia and Alaska (USA), with massive forest losses at the trailing edges. Natural range shifts, however, may be disrupted by the rapidity of climate change and land use changes (Davis and Shaw 2001, Soja et al 2007). Meanwhile, productivity in tropical forests is already increasing (Lewis et al 2009, Pan et al 2011) and MC2 projects large increases in NPP under the REF scenario, with only a little change in fire activity. This result contrasts with decline of tropical forests simulated by some studies (Brienen et al 2015, Cramer et al 2001). Also, deforestation may play a key role in the future for tropical forests (Cramer et al 2004), a process not simulated in our study. This may increase the expansion of the short-rotation plantation wood market globally (Sohngen et al 2001), as temperate forests enjoys a relatively small gain in productivity.

MC2 simulates a shifting balance in global forest condition under the two climate change scenarios. Warmer temperatures, along with higher CO₂ concentrations and fertilization effects, drive forest C stock gains in all regions except Russia under REF and POL3.7 (figure 3(a)). The general trend simulated is consistent with many other terrestrial biome simulations (Friedlingstein et al 2006, Zscheischler et al 2014). However, because of distinct regional differences in climate change exposure (figure 5(a)), the drivers of change (figure 5(b)), and the resulting forest conditions (figures 3(a)−(c)), our models simulate divergent benefits and costs of mitigation for the global forest regions. For several regions the REF scenario is projected to bring significant increases in fire, while the POL3.7 scenario mitigates a significant fraction of those increases. Western Europe and Russia are projected to see significant increases in fire activity under the REF scenario. For Western Europe, the increase is likely driven by the singular increase in temperature without any increase in precipitation, consistent with CMIP3 and CMIP5 projections in the region (Christensen et al 2007, Collins et al 2013). Western Europe and Russia benefit from mitigation by reducing burned forest area by 3.1 and 12 Mha, respectively (figure 3(c)). For these two regions, the mitigation of fire ultimately contributes to the forest carbon stock gains seen under the POL3.7 mitigation scenario (5 and 12 Pg respectively for Western Europe and Russia, figure 3(a)). In contrast, for Western South America, the POL3.7 scenario mitigates burned forest area by 2.6 Mha (figure 3(c)), but the total forest carbon stocks are also reduced ultimately by 7 Pg under this scenario. The reduction in forest carbon stocks are also driven by lower forest productivity (figure 3(a)) and forest contraction (figure 5(b)). For this region, then, climate mitigation has both benefits and costs: mitigation reduces wildfires but also results in reduced forest carbon stocks. The US is under a similar dynamic, where mitigation reduces burned forest area by 0.7 Mha compared to the REF scenario, but the total forest carbon stock is also reduced by 2 Pg. In the US, higher fire suppression costs associated with increases in fire activity (Flannigan et al 2009, Mills et al 2015) may be particularly important in weighing the benefits and costs of mitigation.

Models are abstractions of the natural system, and their accuracy is limited by many types of uncertainties (Uusitalo et al 2015). MC2 simulations abstract the complex global land surfaces into discrete, coarse (0.5°) grid cells, where vegetation is represented by a limited set of plant functional types. Second generation DGVMs may resolve some of the uncertainties arising from coarse representation of vegetation at each grid cell (Scheiter et al 2013). Land use change and forest management practices can have large-scale effects on the forest carbon cycle (Houghton and Hackler 2000, Houghton et al 2000). While we excluded current developed and agricultural areas in our analysis, we also did not simulate the complex history of land use change and forest management practices that occurred on the natural lands. Nor did we simulate the effect of insect and pathogen outbreaks, which can have significant impact on forests, often through interaction with fire (Dale et al 2001). Our study used a single model (MC2) to simulate climate impact on forests. A multi-model ensemble approach could provide results with higher levels of confidence (Littell et al 2011).

All the limitations notwithstanding, running MC2 with a large ensemble of climate simulations allows us to confront the implications of our knowledge (Botkin 1977), and quantify uncertainties due to model formulation (Uusitalo et al 2015). Simulated global forest carbon stock responded in different directions when climate change was mitigated versus when climate sensitivity was decreased (figure 6). Two sources of uncertainty are at interplay here: the CO₂ fertilization effect and climate sensitivity. The CO₂ fertilization effect simulated by MC2 is considered to be moderate (Sheehan et al 2015), although CO₂ fertilization effect still remains poorly understood at the global scale (Körner 1993, Hickler et al 2008). Improving estimates of CO₂ fertilization effect for major vegetation types around the globe and improving estimates of climate sensitivity are needed to reduce uncertainties in projections of forest response. In addition, this study highlights the significant role of natural variability in future projections of vegetation productivity, fire activity, and biome biogeography, a finding consistent with Mills et al (2015). A reformulation of the fire occurrence and spread algorithm may further reduce uncertainties (Parisien and Moritz 2009, Mouillot and Field 2005, Thonicke et al 2001). Finally, the study's results are for realizations from a single climate model. The effects of mitigation policies on the forest carbon stock may be sensitive to climate model selection (Mills et al 2015).