The current capacity of northern high-latitude forests to sequester carbon has been suggested to be undermined by the potential increase in fire and insect outbreaks. Here, we use a process-based terrestrial ecosystem model for the province of British Columbia (BC), Canada to investigate the response of the provinceâs ecosystems to continually increasing atmospheric CO2 concentration and changing climate up to 2100, as well as the recent large mountain pine beetle (MPB) outbreak that started in 1999, in a combined framework. Model simulations suggest that the recent MPB outbreak results in BCâs forests accumulating 490 teragrams less carbon over the 1999-2100 period. Over this same period the changing climate and increasing atmospheric CO2 concentration, however, yield enhanced carbon accumulation that amounts to around 4000-8000 Tg C, depending on the future climate change scenario, which can potentially more than overcome the carbon loss associated with any future large insect outbreaks.
Changes in climate due to rising concentrations of greenhouse gases in the atmosphere affect terrestrial ecosystem processes in a variety of ways. Current evidence suggests that high-latitude forests are sequestering carbon (Ciais et al., 2013; Gourdji et al., 2012). Increasing atmospheric CO2 concentration, which increases photosynthesis through the CO2 fertilization effect, and the associated changes in climate which is gradually getting warmer and wetter, on average, in high-latitude regions (Hartmann et al., 2013; Mekis and Vincent, 2011) are considered to be the primary reasons for this high-latitude carbon sink. Climate change, however, can potentially also increase occurrences of both insect outbreaks and fire and these disturbances have been suggested to undermine the ability of northern forests to take up and store atmospheric carbon. Kurz et al. (2008b), for example, project the carbon balance of the managed Canadian forests to year 2022 by taking into account stochastic future disturbances and using the Carbon Budget Model of the Canadian Forest Sector (CBM-CFS3). They project that the managed Canadian forests will remain a source of carbon to the atmosphere until 2022. The CBM-CFS3 model, however, does not explicitly take into account any effect of increased growth of atmospheric CO2 and the associated climate change on forest growth rates, other than what is implicitly included in the empirical growth and yield data that are used to drive the model.
Here, we use a process-based terrestrial ecosystem model to explicitly simulate the effects of the recent large mountain pine beetle (MPB) outbreak that started in 1999 and the historical and future changes in climate and atmospheric CO2 concentration on the terrestrial carbon balance of the province of British Columbia (BC), Canada, in a combined framework. Our objective is to put the carbon loss associated with the recent MPB outbreak in the context of the carbon gain associated with the increase in atmospheric CO2 and climate change over historical and future periods. The inherent temporal and spatial scales associated with insect disturbances are, of course, different from those of climate change associated with increasing concentrations of greenhouse gases. Major insect disturbances are periodic, often cyclical, events that typically affect a fraction of the landscape over 2-15 year time periods. The effects of climate change and increasing CO2, in contrast, are more wide spread, gradual and expected to last over decades to centuries. Nevertheless, a modelling framework that simulates the natural response of terrestrial ecosystems to an insect disturbance and increasing atmospheric CO2 together with associated changing climate, in a combined framework, provides the means to compare the effect of both in a consistent manner.
2. Model, methods and experimental set up
Our study is based on simulations made with the Canadian Terrestrial Ecosystem Model (CTEM, version 1.1) which is a process-based ecosystem model, coupled to the Canadian Land Surface Scheme (CLASS, version 3.5) (Peng et al., 2014). The coupled CLASS-CTEM modelling framework has been used for the province of BC, Canada over the historical 1900-2010 period driven with observation-based climate and atmospheric CO2 concentration data at a spatial resolution of around 40 km (Peng et al., 2014) (see Figures S1 and S2 in supporting information), but without the MPB disturbance. In this framework, the model simulates vegetation growth and calculates time-varying carbon storage in three living vegetation pools (leaves, stems and roots) and two dead carbon pools (litter and soil organic matter) for four plant functional types (PFTs) that are present in the province of BC (coastal and interior needleleaf evergreen trees, broadleaf cold deciduous trees and C3 grasses). These PFTs occupy specified fractions in each 40 km grid cell, which do not change over time, and are based on a 25-m resolution land cover dataset (corresponding to year 2000) from the Canadian Forest Services Earth Observation for Sustainable Development of Forests (EOSD) (Wulder et al., 2003). The primary terrestrial ecosystem processes that are modelled in CTEM for this study are photosynthesis, autotrophic and heterotrophic respiration, allocation, phenology, turnover of leaves, stems and roots, and fire.
Here, we extend the historical simulation from 2006 onwards to 2100 for three representative CO2 concentration pathways (RCP 2.6, 4.5 and 8.5) which correspond to the three future climate change scenarios using climate data from the Max Planck Institute for Meteorology (MPI) Earth System Model (MPI-ESM-LR), and also include the effect of the MPB disturbance that started in 1999. The year 2100 atmospheric CO2 concentration in the RCP 2.6, 4.5 and 8.5 scenarios is around 420, 540 and 940 ppm, respectively (see supporting information). We used climate data from the Max Planck Institute for Meteorology (MPI) Earth System Model (MPI-ESM-LR) (Giorgetta et al., 2013) because the seasonality of climate data from this model matched well with observation-based estimates compared to other models for BC, Canada and as a result provided a smooth transitioning from the historical observation-based climate data to model-based climate data in 2006. The climate data from the MPI-ESM-LR were, however, adjusted for biases in annual mean values as well as for inter-annual variability. The historical climate data were obtained from the CRUNCEP data set for the period 1901-2005, which are based on the National Centre for Environmental Prediction (NCEP) reanalysis (Kanamitsu et al., 2002) with monthly means adjusted to match the Climate Research Unit (CRU) observations. The CRUNCEP data (surface temperature, pressure, precipitation, wind, specific humidity, shortwave and longwave radiation fluxes) are available at a resolution of 0.5 degrees and at a six hourly time interval. Data were extracted for BC and spatially interpolated to the ~40 km Ã 40 km grid used in this study (see Figure S2). The province area of 1,005,388 km2 used in the model is about 6% greater than the actual area of 944,700 km2, since some grid cells at the provinceâs boundary lie partially outside its borders. The six-hourly data were disaggregated to half-hourly time resolution. The resulting annual time series of province-wide average climate variables, together with the atmospheric CO2 concentration, used to drive the CLASS-CTEM model are shown in Figure S1 for the historical and future periods.
The MPB disturbance is implemented on the basis of spatially-distributed data for severity index and cumulative pine volume killed available at 450 m resolution from the BC Ministry of Forests, Lands and Natural Resource Operations (MFLNRO) for the period 1999 to 2011 (http://www.for.gov.bc.ca/ftp/HRE/external/!publish/web/BCMPB/Year9/). In the absence of spatially-distributed infested fraction (I, fraction) in each 450 m grid cell we relate severity index to fraction of killed trees observed during aerial overview surveys using Table 2 of Westfall and Ebata (2011). The cumulative pine volume killed data are available as fraction of trees killed in a grid cell and used to estimate annual fraction of trees killed (K, fraction/year). The rate of change of infested area is given by , and when I and K are known the annual recovery fraction (R, fraction/year) can be estimated. Finally, K and R are regridded from the 450 m resolution to the 40 km resolution at which the CLASS-CTEM model is implemented.
The fraction of interior needleleaf evergreen trees killed in a given year is gradually implemented during a year with the killed fraction increased from May 1 to Sep 30. The trees of the killed fraction stop photosynthesizing and their needles, stem and root contribute to the litter pool at a rate faster than the normal turnover rate for these components. As a result the ecosystem loses carbon. Conversely, as recovery occurs the recovered fraction begins to photosynthesize and slowly the ecosystem returns to being a carbon sink. The dead needles are assumed to have a leaf life span of 0.5 years. The stem component in the model for MPB affected trees has a half-life of around eight years. The half-lives for leaf and stem components determine the rate at which these components turnover and contribute to the litter pool. The half-life of eight years used for the stem component of the MPB affected trees is consistent with the 10-year half-life used in another modelling study (Edburg et al., 2011) and the observation-based estimates (Lewis and Hartley, 2006) which suggest in unthinned stands, 50% of trees had fallen within nine years and 90% within 14 years after their death. In any case, Figure S3 shows that while the simulated peak impact of the MPB disturbance is sensitive to the chosen half-life of the stem component of the MPB affected trees, the cumulative impact is not. The needle life span for the healthy interior needleleaf evergreen trees is 5 years and the default half-life for their stem component is 48 years. The recovered fraction each year, which is also increased gradually from May 1 to Sep 30, is assigned to a new age cohort that starts growing vegetation biomass as an interior needleleaf evergreen PFT. The simulated response to the MPB disturbance, which is the result of reduction in photosynthesis and increase in heterotrophic respiration following the disturbance, is illustrated for a single grid cell in Figure S4.
The modelâs response to increasing atmospheric CO2 and changing climate has been evaluated against observation-based stemwood growth rate in coastal British Columbia (Peng et al., 2014). For the period 1959-1998, CTEM simulated a rate of increase of stemwood growth of 2.7% and the observed inventory-based rate of increase of stemwood growth (Hember et al., 2012) was found to be 3.0% in response to extrinsic factors of climate change and CO2 fertilization (which stimulates plant growth through increased photosynthesis rates).
Figure 1a shows that the time series of province-wide infested area compares reasonably well with the observation-based infested area. The newly infested area has been decreasing in recent years as the MPB is invading new areas but at a diminishing rate. Figures 1b through 1e display the spatial distribution of cumulative fraction of trees killed (K) due to the recent MPB outbreak for years 1999, 2003, 2007 and 2011 at 40 km resolution at which the model is implemented and show how the MPB spread through the interior of the province of BC over time.
Figure 2a quantifies the effect of the recent 1999-2011 MPB disturbance using results from the historical (1900-2005) and future (2006-2100) simulations for the three climate change scenarios (RCP 2.6, 4.5 and 8.5) that are performed with and without the MPB disturbance. The results show the modelled net land-atmosphere CO2 flux for the RCP 4.5 scenario. The simulated net land-atmosphere CO2 flux represents the model response to changing climate and increasing CO2 and, when included, the MPB disturbance. Positive values of net atmosphere-land CO2 flux indicate that land is gaining carbon. The simulated net land-atmosphere CO2 flux is lower in the simulation with the MPB disturbance after 1999 because of lower photosynthesis and higher heterotrophic respiratory losses. The strength of the sink over the 1990s and early 2000s is around 44 g C m-2 year-1 and may be compared to an inversion-based estimate for 2003 of 38Â±66 g C m-2 year-1 (Deng et al., 2007), indicating that model response to warming climate and increasing atmospheric CO2 is reasonable. The decrease in net atmosphere-land CO2 flux after year 2000, even in the absence of the MPB disturbance, is due to climate variability. Other aspects of the simulated terrestrial carbon cycle, including gross primary productivity, vegetation biomass and stem wood growth rate have been assessed by Peng et al. (2014) and compare reasonably with observation-based estimates for the historical period. The difference in the net land-atmosphere CO2 flux between simulations with and without the MPB disturbance yields an estimate of the effect of the MPB disturbance (Figure 2b, the green line corresponds to the RCP 4.5 scenario). Figure 2b also shows the estimated effect of the MPB disturbance when the future simulations (with and without MPB disturbance) are driven with climate and atmospheric CO2 concentration corresponding to the RCP 2.6 and 8.5 scenarios. The estimate of the effect of the MPB disturbance is broadly insensitive to the future scenario chosen. These results indicate that at its peak in 2006 the MPB disturbance reduced the net atmosphere-land CO2 flux averaged over the total province area by around 24 g C m-2 year-1, or equivalently by 36 g C m-2 year-1 over the treed area. The effect of the MPB since then has diminished as the rate of MPB-related tree mortality is declining and forests are recovering. This recovery is simulated to continue in the future. Regardless of the changing climate and increasing CO2, forests recover naturally from disturbances. Changing climate and increasing CO2 imply that the rate of this recovery is higher than the rate if the forests were to recover in the absence of changing climate and increasing CO2. The cumulative loss over the period 1999-2020 is around 326 gC m-2 which yields a cumulative reduced carbon uptake of 328 Tg C when multiplied by the provincial area of 1,005,388 km2 used in the model. This estimate is in broad agreement with the estimate of 270 Tg C calculated by Kurz et al. (2008a) using the Carbon Budget Model of the Canadian Forest Sector, CBM-CFS3, for the same time period. The cumulative reduced carbon uptake due to the MPB outbreak over the period 1999-2100 is 492 Tg C.
In Figure 2b, the model simulates that by around 2060 the net effect of the MPB disturbance, averaged over the whole province, will turn positive. That is, by around this time the provinceâs forests will be taking up more carbon compared to the simulation without disturbance. This is consistent with the observed response of terrestrial ecosystems following a disturbance where the period of net loss is followed by a period of first increasing and then decreasing carbon gain as the system recovers (e.g. see Figure 1 of Kurz et al., 2013).
The effects of climate change and increasing CO2 concentration are summarized in Figure 3a which shows the simulated net atmosphere-land CO2 flux for the historical and future periods for the three climate change scenarios. These simulations include the effect of the MPB disturbance. The net province-wide averaged atmosphere-land CO2 flux is simulated to return to its 1990s value in the 2030s, depending on the future scenario, and the provinceâs ecosystems are simulated to continue to take up carbon in response to increasing CO2 and a changing climate that gets warmer and wetter (future temperature and precipitation increase in all scenarios, see Figure S1 in supporting information). Differences in the scenarios emerge after 2040 and cumulative uptake is highest in the RCP 8.5 scenario and lowest in the RCP 2.6 scenario. The sink generated due to changing climate and increasing atmospheric CO2, as reflected by positive values of the net atmosphere-land CO2 flux in Figure 3a, is the result of a larger increase in net primary productivity than in heterotrophic respiration (see Figure S5 in supplementary information). The net result is an increase in total land carbon, consisting of carbon in the modelâs vegetation, litter and soil components, as shown in Figure 3b for simulations with and without the MPB disturbance. The change in total land carbon which is equivalent to cumulative atmosphere-land CO2 flux, shown on the right hand side y-axis of Figure 3b, ranges from an increase of about 900 Tg C for the RCP 2.6 scenario to 1060 Tg C for the RCP 8.5 scenario, for the period 1999-2020, in simulations without the MPB disturbance (light coloured lines) indicating that the reduced carbon uptake by land due to the MPB disturbance of around 328 Tg C is already surpassed by 2020. The cumulative carbon uptake in simulations with the MPB disturbance (dark coloured lines) is consequently lower by 328 Tg C (~570 Tg C in RCP 2.6 scenatio and ~730 Tg C in RCP 8.5 scenario), for the period 1999-2020, because of the reduced carbon uptake associated with the MPB disturbance.
For the period 1999-2100, the increase in total land carbon is around 4000 Tg C for the RCP 2.6 scenario and about 8000 Tg C for the RCP 8.5 scenario. Litter and soil carbon increase by about 2000 Tg C in all scenarios so the remaining 2000 Tg C (RCP 2.6) to 6000 Tg C (RCP 8.5) increase is due to increase in vegetation biomass (see Figure 4). This province-wide increase in vegetation biomass by 2100 is equivalent to about 22% and 58% increase in diameter at breast height of trees in the RCP 2.6 and 8.5 scenarios, respectively, if vegetation height, stem density and wood density do not change substantially. At 2100, the dark coloured lines for each scenario that correspond to simulations with the MPB disturbance are about 490 Tg C lower than the values corresponding to the simulations without the MPB disturbance (light coloured lines), consistent with the 492 Tg C cumulative effect of the MPB disturbance over the period 1999-2100, and due to the fact that net primary productivity (NPP) is continually increasing (see Figure S5, panel a).
The model also simulates area burned and CO2 emissions associated with forest fires based on a fire parameterization of intermediate complexity (Arora and Boer, 2005). The area burned in the province of BC is, however, small (~0.08% area burned annually over the period 1970-2010). The simulated area burned in the model for the same time period is higher (~0.28% burned annually) which results in forest fire CO2 emissions of ~0.8 Tg C/year. The simulated annual area burned and emissions, averaged across all scenarios, for the future 2070-2100 period increase to 0.40% and 1.9 Tg C/year, respectively, averaged across all scenarios. However, the cumulative increase in fire emissions only contributes to a source of 82 Tg C (over the 1999-2100 period) – much less than the cumulative impact of the MPB disturbance and the response to climate change and increasing atmospheric CO2.
4. Discussion and conclusions
Limitations remain in our modelling framework. First, the model version used in our study does not dynamically simulate the fractional coverage of its PFTs. Studies that have used bioclimatic envelopes for projecting future distribution of tree species in BC (Hamann and Wang, 2006; Wang et al., 2012) suggest that climate envelopes for relatively productive species that currently exist in coastal and mild-climate interior regions will expand over much of BC at the expense of less productive sub-boreal, subalpine, and alpine ecosystems. The implications for the resulting carbon balance of these projected changes, however, are unclear since tree species do not migrate as quickly as the projected shifts in climate envelopes. Second, the terrestrial ecosystem model used does not account for the age distribution of forest stands. The modelled response to climate change and increasing atmospheric CO2 is based on that of an average-aged tree in the landscape, without an explicit representation of self-thinning that would increase mortality as biomass increases. The average forest age in BC is increasing and so the age-dependent reduction in tree growth counteracts the environmentally-driven growth enhancement (Hember et al., 2012). Third, although CTEM includes a parameterization of down-regulation of photosynthesis (Arora et al., 2009), due to nutrient limitations, as atmospheric CO2 concentration increases (based on results from plants grown in ambient and elevated CO2 environments), it does not include an explicit coupling of terrestrial carbon and nitrogen cycles. Consideration of age-class structure and an explicit coupling of terrestrial carbon and nitrogen cycles are both expected to reduce the simulated response to future climatic change and increasing CO2. Finally, our modelling framework does not include the effect of wood harvest and other smaller insect disturbances that have occurred over the historical period. Recent harvest related transfers from the forest to harvested wood products are estimated to be around 18 g C m-2 year-1 when averaged over the whole province (based on updated estimates from Stinson et al. (2011) from Natural Resources Canada, National Forest Carbon Monitoring, Accounting and Reporting System). Explicit modelling of harvest related carbon transfers will change the absolute values of the model response to the MPB disturbance and changing climate and increasing CO2. However, we do not expect inclusion of harvesting to substantially change the large response of the model to climate change and increasing atmospheric CO2 relative to the MPB disturbance.
The response of the CLASS-CTEM modelling framework to changing climate and increasing CO2, expressed in terms of stem wood growth rate, compares well with observation-based estimates in coastal BC (Peng et al., 2014) and so does the magnitude of the simulated sink of around 44 g C m-2 year-1 during the 1990s and early 2000s. Our estimate of the cumulative effect of the MPB disturbance for the period 2000-2020 of 325 Tg C (based on Figure 2b) is also in broad agreement with an existing estimate from Kurz et al. (2008a). The response of the CLASS-CTEM model to climate change and increasing CO2 and the recent MPB insect disturbance is thus consistent with observation-based and other estimates over the historical period. The extended simulations to 2100 then suggest that the enhanced carbon uptake by the forests of BC in response to increasingly warmer and wetter climate and gradually increasing CO2 can potentially easily overcome the carbon losses due to the recent as well as any future MPB outbreaks given their characteristic return interval of around 40 years in the province (Alfaro et al., 2010). The extent to which factors not included in our simulations will modify these increases is subject of ongoing research and therefore the results obtained from application of the CLASS-CTEM model for future climatic and atmospheric CO2 conditions must be interpreted in that context.