Forests have been known an important terrestrial carbon (C) sink [ 1 - 6 ]. Generally, C sequestration in forests can be expressed by using CO2 flux measurements. Recent development and applications of eddy flux measurement technology reflect these metrics [ 7 - 11 ]. However, CO2 flux is highly impacted by changing environmental factors, including topography, climate, hydrogeology, soil, vegetation and various disturbances [ 12 - 14 ]. Therefore, there are large uncertainties in the flux estimation using this technology due to large differences in environmental conditions and inadequate equipment to cover large enough areas to do the assessments for large regions.
Although understanding C sequestration in forests is fundamental for REDD, it is impossible for us to conduct intensive inventories everywhere to know C stocks and long-term C dynamics in forests in large regions, such as nationwide or continent-wide, because of limitations in personnel, equipment, funds and complex environmental conditions. However, we can do this by means of computer models, which are developed from expert knowledge, long-term experiences, and observations. Recent applications of biogeochemical C models for assessing forests response to land use change and natural and human disturbances reflect their merits [ 2 , 15 - 21 ]. Accordingly, application of these effective tools is fundamental to evaluate the responses of forests to climate change and the role of forests in mitigating global warming, and to effectively assess C stocks and long-term dynamics for forest management and restoration.
The main objective of this study was to use the spatially-explicit and process-based model Forest- DNDC (FDNDC) to simulated spatiotemporal C dynamics, involving threefold: 1) to evaluate the model performance using observations from sample plots in Nez Perce-Clearwater National Forest in Idaho of USA, 2) to assess spatiotemporal C stocks for the forest, and 3) to estimate the responses of C sequestration to disturbances and climate variability. To achieve the first objective, FDNDC was evaluated using the data of FIA (Forest Inventory and Analysis) sample plot data to simulate C dynamics for Nez Perce-Clearwater National Forest. To achieve the other two objectives, the evaluated model was used to assess spatiotemporal C stocks in the forest and the responses of C sequestration to climate variability and disturbances using spatial climate data for a 33-year period from 1981 to 2013 and spatial disturbance information for fires, insects and harvest occurred in the period 1991-2011.
2. METHODS AND DATA
2.1. Site Description
The study site, Nez Perce-Clearwater National Forest (NPCNF), is a temperate forest, located in northern Idaho of USA (45.27˚N - 47.13˚N, 114.3˚W - 117.0˚W) (Figure 1). The forest area is over 16,000 km2 at present. Geographical gradient alters greatly; the slope ranges from 1 to 52% with a mean of 10%; the elevation varies between about 390 and 2710 m above mean sea level, with a mean of over 1400 m.
Climate in this forest area alters substantially with changing gradient; precipitation increases with an increment in the altitude, but air temperature is converse. Mean annual precipitation was 550 - 1730 mm over the forest in the 33-year period from 1981 to 2013 based on the Daymet data (http://daymet.ornl.gov/), with a regional mean of 1100 mm and an annual precipitation range of 400 - 2440 mm over space and time in the same period. Spatial difference in long-term mean temperature was large in the last 33-year period, ranged from 0.0˚C to 9.6˚C across this forest with a regional mean of about 5.1˚C in the same time period, and the annual mean over the time and space ranged from −4˚C to 12˚C. Because of changing in the altitude, the snow is different place-to-place, ranged from 5 to 760 kg m−2 yr−1, with an average of 300 kg m−2 yr−1 in the last 33 years.
Figure 1. Nez perce clearwater national forest in Idaho (purple area), USA.
There are 55 soil series in the area of NPCNF due to complex geographical topography and climate, they are distributed on different slopes. Soil type ranges from loamy sand to sandy clay loam, but the main soil type is loam; the clay content is between 2.5% and 50.0%, the mean is about 11.3%. The acidity of the soils is different, pH 4.8 - 8.2, with a mean pH of 6.2, based on the data obtained from NRC downloaded from www.nrc.org. The rock content in the soils varies largely, ranging from 0% to 97%, with a mean of about 30%, substantially increasing or decreasing with an increase or decrease in slopes (R2 = 0.11, n = 1056, P = 0.01), and the soil organic matter is between 0.25% and 50%, with an average of 2.15%, on the basis of the NRC data.
There are over 20 tree species in this forest, but the main species are coniferous, including fir (Pinaceae) family [Douglas fir (Pseudotsuga menziesii), grand fir (Abies grandis), and subalpine fir (Abies lasiocarpa)], spruce (Picea engelmannii), lodgepole pine (Pinus contorta) and ponderosa pine (Pinus ponderosa), western redcedar (Thuja plicata). Due to disturbances (see Data below), stand age ranged from seedling to 818 years old in 2010 based on the data of FIA plots in the NPCNF [ 22 ].
Climate data for a 33-year period from 1981 to 2013 were obtained from Daymet database (http://daymet.ornl.gov/). Spatial climatic data points used for this study were over 1000, about 4.2 km resolution on average. A poly-point map was created using ArcGIS for the spatial climate data. Soil data for NPCNF, including map and physiochemical soil properties, were obtained from NRC database, downloaded from www.nrc.org. The data of vegetation species/types and inventories used to calibrate and validate the model FDNDC for modeling C dynamics in this forest were obtained from the FIA database.
Disturbance map for NPCNF were derived from the Landsat TM time series using the Vegetation Change Tracker [ 23 ] (Figure 2). The disturbance types were fire, harvest/thinning and insect, identified using Monitoring Trends in Burn Severity (MTBS, http://www.mtbs.gov/), Forest Service Activity Tracking System (FACTS, http://www.fs.fed.us/nrm/index.shtml), and Aerial Detection Survey (ADS, http://www.fs.fed.us/). Disturbance magnitudes in NPCNF were initially identified using post-disturbance change in modeled live aboveground carbon (for harvests) or percent canopy cover (for fires, insects, abiotic), then calibrated with FIA sample plots (www.fia.fs.fed.us/tools-data). Disturbance magnitudes were grouped in four bins (0% - 25%, 25% - 50%, 50% - 75%, 75% - 100%) in canopy cover [ 24 ]. The 30-m
Figure 2. Disturbances occurred in 1991-2011 in nez perce clearwater national forest, Idaho, USA.
resolution gridded disturbance magnitude datasets were nearest-neighbor interpolated to 90-m resolution.
To simulate spatial C stocks in NPCNF, the spatial datasets, including soil properties, vegetation and climate data, were alimented to 90 m resolution, as same as the resolution used for disturbances. The polygonal vegetation and soil and poly-point climate and disturbance maps were sequentially joined (overlapped) to create spatial units (polygons) using ArcGIS 10 to assess spatial C stocks in this forest.
FDNDC is a spatially explicit biogeochemical model used to simulate forest growth and C and N dynamics in forest ecosystems and emissions of trace gases such as CO2, CH4 and N2O based on the balance of water, light, and nutrients in forest ecosystems [ 2 , 25 , 26 ]. Vegetation is divided into three layers, i.e., over-story (dominant canopy), understory (saplings) and ground-growth (sedges and moss/lichen). Each of the vegetation layers is simulated based on competition for energy and nutrients and plant physiological process responses to variations in environmental conditions. The soil profile is divided into multiple layers, 1 - 5 cm in thickness; soil conditions and the dynamics of C and N in each soil layer are modeled hourly. This model has been widely tested and used for estimating greenhouse gas (GHG) emissions from forested wetland and upland ecosystems and assessing C sequestration in forests in a wide range of climatic regions, from boreal to tropical [ 11 , 26 - 32 ]. The model structure, algorithms and parameterization can be found in the model manual (http://www.dndc.sr.unh.edu) and publications [ 25 - 28 ].
2.4. Model Performance Evaluation
FDNDC was calibrated and validated using 203 and 381 plots, respectively, to simulate C dynamics in the NPCNF. The information of these 584 sample plots was added to the joined map mentioned above, all information for vegetation, soil, and disturbances were verified, and then the model was independently parameterized using these plot-based information for evaluating the model performance. The key parameters are presented in Table 1.
Table 1. Key vegetation and soil parameters for Forest-DNDC*.
*: Values for vegetation, soil and climate are spatially based, data were obtained from FIA, NRC and Daymet (see Site Description and Data); #: Usually, the coefficient of stem density is the ratio of forest to the bare (non-forest) areas in each simulating unit, and from 0 - 1.
To prevent bias from individual method, model performance was evaluated employing four quantitative methods [ 32 , 33 ]: the coefficient of determination (R2, squared correlation coefficient), model performance efficiency (E) [ 34 ], percent bias (PBIAS), and the RRS [the ratio of the root mean squared error (RMSE) to SD (standard deviation)] [ 35 ].
The key variable used to evaluate model performance is E (−∞, 1), and it is calculated as
where , and are observed values, observed mean and simulated results, respectively. The evaluation variables, PBIAS and RRS are computed, respectively, as
where SD is the observation standard deviation; RMSE is the root mean squared error, the equation is
where n is the number of samples, or the pairs of the observed and simulated values; and the R2 is
Model performance rating ranges are: 0.25 ≤ E < 0.5 (satisfactory), 0.5 ≤ E < 0.75 (good), and E ≥ 0.75 (excellent or very good); a PBIAS between −25 and 25 and a RRS less than 0.7 are considered satisfactory [ 32 , 35 , 36 ].
FDNDC was run for a 225-year period starting from 1788 when the oldest one of the 584 FIA plots regenerated to evaluate the model for assessing C stocks in NPCNF.
Because there were no available climate data for the period before 1981 and disturbance data for the period before 1991, the 33-year climate data and 21-year disturbance information (see Section 2.2) were repeated for the simulated historical period based on the assumption that the climatic pattern and the type and magnitude of disturbances in history were similar with those for observation periods in this forest.
2.5. Model Parameterization for Assessing Spatial Carbon Stocks
The model was configured to simulate spatial C dynamics in NPCNF. Accordingly, the joined map mentioned above was used to generate spatial simulation units using the tools of FDNDC (see details in the Forest-DNDC User’s Guide, http://www.dndc.sr.unh.edu). The model was parameterized for those spatial simulation units using spatial vegetation, soil, climatic and disturbance information. The model was run to simulate the spatial C dynamics in NPCNF and C stocks for the period 1990-2012.
3.1. Model Evaluation
The result of the model calibration using 203 sample plots with a stand age range of 4 - 215 years old in 2012 and all main tree species in NPCNF showed that the simulated aboveground biomass using FDNDC (Figure 3(a)) was in a good agreement with the observation (R2 = 0.91, P = 0.001, n = 203) with a
Figure 3. (a) Biomass (Mg C ha−1) simulated vs. observed for the model calibration (red dash line, 1:1); (b) Biomass (Mg C ha−1) simulated vs. observed for the model validation (red dash line, 1:1).
slope (0.97) of the regression model between simulation and observation that was approximate 1.0 and an intercept (11.25 Mg C ha−1) that was about 13.8% of the observed mean (79.77 Mg C ha−1). The simulated mean (83.65) was about 4.9% higher than the observed.
The result of the simulation for model validation against observations from 381 sample plots is plotted in Figure 3(b), indicating that FDNDC performed properly well with a proper slope (0.94) and a small intercept (5.87 Mg C ha−1, about 7.4% of observation average) of the regression model between observation and simulation. The simulated mean (80.34 Mg C ha−1) was only about 1.9% higher than the observed (78.83 Mg C ha−1). The four quantitative evaluation variables (Table 2) concluded consistently that FDNDC performed excellently to estimate C sequestration in this temperate forest in Idaho in USA with high performance efficiency; based on the model performance rating, FDNDC performance for assessing C stocks in stands over the temperate forest was within the range of “very good” (E ≥ 0.75; RRS ≤ 0.7; −25 ≤ PBIAS ≤ 25) [ 35 , 36 ].
3.2. Spatiotemporal Carbon Sequestration
The C sequestration in the temperate forest varies considerably in space, aboveground biomass (Figure 4(a)) and net primary production (NPP) (Figure 4(b)) in this forest ranged from 1.1 to 368.1 Mg C ha−1 (Figure 4(a)) with a mean of 85.1 Mg C ha−1 and 33.9 to 649.5 g C m−2 with a mean of 238.7 g C m−2, respectively, in 2010. Mean aboveground biomass increased annually by about 0.6 Mg C ha−1 in this forest in the period 1990-2012, excluding the harvested wood products and dead woods due to disturbances of fires and insects.
The spatial divergence is mainly associated with the stand age over the forest, also related to climate
Table 2. Results from model performance*.
*: R2, coefficient of determination; E, model performance efficiency; PBIAS, percent bias; RRS, the ratio of the root mean squared error (RMSE) to SD (standard deviation); a and b, the slope and intercept of the regression model between observations and simulations.
Figure 4. (a) Spatial aboveground biomass (Mg C ha−1) in NPCNF (~16,000 km2) in 2010; (b) Spatial NPP (g C m−2) in NPCNF (~16,000 km2) in 2010.
and soils, but the impacts of climate and soils are smaller than the stand age. The biomass from spatial simulation for 855 FIA plots were extracted to compare the simulation with observations. The results showed that the simulated mean biomass (81.4 Mg C ha−1) was approximately to the observation (79.7 Mg C ha−1). The simulated belowground biomass, including live and dead roots, was 23.3 Mg C ha−1, about 11% higher than the value estimated (21.0 Mg C ha−1) by FIA.
C sequestration varies temporally, and there are differences among species in NPCNF. The results from the simulation showed that net primary productivity (NPP) reached the peak at different ages (Figure 5(a)). Ponderosa might reach the NPP peak at about age ≥10 years old, but Douglas fir and Grand fir
Figure 5. (a) Annual net primary productivity (NPP) vs stand age in NPCNF, Idaho, USA; (b) Aboveground biomass (Mg C ha−1) vs stand age in NPCNF, Idaho, USA; Subalpine species includes subalpine fir and subalpine spruce.
were at about age ≥18 years old. The NPP correlation to stand age was significant (P < 2.4E−73), described by
where Age is the stand age; and are coefficients; . However, the value of m is different specie-to-species in this forest, m = 2 for grand fir, subalpine fir, spruce, subalpine spruce, western larch and ponderosa, the m = 3 for western cedar and lodgepole pine, and m = 6 for Douglas fir, except for those coefficients that are different. This difference indicates that there are some differences in C sequestration rate among the species in the forest.
The interannual changes in NPP were caused principally by climate. The fluctuation exhibited that the sensitivity of different species to climate variability in this forest was different, Grand fir was more sensitive to climate in this forest than Douglas fir and ponderosa pine, indicating that climate change can impact NPP in this forest.
The relationship between the aboveground wood biomass and stand age (Figure 5(b)) is slightly different from NPP in this forest, a quartic polynomial (e.g., m = 4 in Equation (7)) may better describe the relationship for all main species in this forest (P < 9.3E−275), as follows
where AGB is accumulative aboveground wood mass (Mg C ha−1); others are as same as those used in Equation (6). However, the changes in AGB with stand age were different species-to-species although their growth curves were similar and nonlinear. Subalpine species (subalpine fir and subalpine spruce) have lower growth rates than others, but their AGB can continuously increase with an increment in stand age for a longer period.
3.3. Impact of Disturbances on Carbon Stocks
Disturbances affect substantially C stocks in NPCNF. Main disturbances occurred in this forest area were insects, fire and harvest. Total disturbed area was about 1062.14 km2 within the 21 years period from 1991 to 2011, about 859.88, 116.42 and 85.84 km2 disturbed by insects, harvest and fires, respectively. These factors caused a substantial decrease in NBP (net biome production) (Figure 6(a)), especially harvest
Figure 6. (a) Impact of disturbances on NBP, and C removal from the forest in NPCNF due to disturbances occurred in 1991-2011; (b) Interannual changes in carbon in mineral soils (ΔSOC, Mg C; solid black and grey dash lines) and forest floor (Δlitter, Mg C; dark and bright red lines).
and canopy fires that were main factors reducing C stocks in the forest. The total harvested C was 763.54 Gg C in the 21 years with a mean of 36.36 Gg C yr−1. However, fires removed about 1742.38 Gg C in the 21 years with a mean of 82.97 Gg C yr−1, indicated that fires plaid an important role in removing C from the forest, reduced forest floor. No matter what type the fires are, canopy or ground fires, they can cause a decrease in forest floor.
Changes in forest floor in the forest (Figure 6(b)) reflected disturbances that impacted substantially the C stocks in this forest; a substantial increase in forest floor was caused by harvest due to harvest residues, but a substantial decrease in forest floor reflected that a severe fire occurred. Total C loss to atmosphere from this forest due to fires in the 21 years was over 950 Gg C, including C loss from forest floor and canopy. Accordingly, there were substantial interannual changes in forest floor in this forest due to disturbance events.
There were about 4100 Gg C of deadwoods produced in the period 1991-2011 due mainly to fire and insects, excluded those deadwoods decomposed slowly in natural conditions. However, insects played very limited role in directly removing C from the forest, their main role was to produce a large amount of deadwoods. Accordingly, accumulative mean deadwood was about 34.3 Mg C ha−1 in 2003 based on FIA plot data.
4.1. Model Evaluation
The model performance evaluated using observations from 584 sample plots in NPCNF shows that the FDNDC can perform well to assess C stocks in this forest, with high model performance efficiency (E ≥ 0.83, R2 ≥ 0.85), proper slope (a ≥ 0.94) and intercept of the regression model between the observation and simulation. The mean error between the observation and simulation was 9.7%. The four evaluation variables (E, R2, PBIAS and RRS in Table 2) concluded consistently that this model was applicable for assessing C stocks in NPCNF with a proper model performance efficiency.
4.2. Differences between Observation and Simulation
This model over-estimated live aboveground biomass for this forest by about 9.7%. This error may be mainly because of unknown historical disturbances before 1991 that we assumed that the disturbances occurred before 1991 were as same as those in the 21-year recorded period 1991-2011, including disturbance types and magnitudes. The assumption can cause an overestimation of live aboveground biomass and underestimation of deadwoods in this forest. However, the total aboveground wood mass for simulation (110.7 Mg C ha−1), including live and dead woods, was only about 2.9% less than the observation (114.0 Mg C ha−1) for 855 FIA plots in 2003, indicated that the simulated mean biomass growth rate might be slightly lower than the actual value for this forest.
4.3. Impact of Disturbances on Carbon Stocks in Forests
Wild fires are not uncommon, and their impact on C stocks in forests may be different fire-to-fire. The results simulated for the fires in NPCNF showed that those fires had produced a large amount of dead woods and C loss to atmosphere. The wild fire occurred in NPCNF in 2000 burned the forest of 31.7 km2 with a mean canopy loss of about 35.2% based on the estimate of the canopy loss to fires. Based on the simulation using fire disturbance data (Figure 2), the fires left dead woods of about 94.4 Gg C and caused a loss of over 40 Gg C to air due mainly to a high fuel loading, including dead woods and forest floor in this forest, so that the forest floor in the burnt area was lower in 2001 than pre-fire in 2000; and total C loss to atmosphere from burning was over 50 Gg C in 2000 if the fine branches and leaves of live trees burnt were included.
The impact of insects on C sequestration may be slightly different from other disturbances. Generally, insects may mainly cause tree death or leaf loss. Accordingly, a large scale of insects may lead to a great amount of live trees died. The simulated result for NPCNF showed that the dead woods produced only in 1993-1994 were over 481.9 Gg C due to the outbreak of insects over 65.7 km2. Because the decomposition of dead woods, especially standing dead woods, is slow, there are a large amount of dead woods in this forest; total dead woods, produced by all factors, simulated and observed for NPCNF were 29.3 and 34.33 Mg C ha−1 in 2003, respectively, indicating that deadwoods are an important C pool in this forest.
Harvest is a common anthropogenic disturbance, substantially impacts C stocks in forests. Harvest removed forest product of about 763.5 Gg C from NPCNF in 1991-2011. However, harvest is different from other disturbances. If those harvested woods were used in longer life commercial products, C storage should increase substantially. However, harvest may increase wild fire severity due to fuel loading added by harvest based on the findings of Stone et al. [ 37 ], and harvest can leave a large amount of dead roots, those dead roots and residues left by harvesting can cause a subsequent increase in soil CO2 flux due to decomposition.
4.4. Effect of Climate Variability on Carbon Sequestration
The impact of climate variability on C stocks in NPCNF was estimated using the Daymet climate data for a 33-year period 1981-2013. The results show that climate variability can substantially impact C stocks in the forest in NPCNF (Figure 7). Annual mean aboveground biomass in NPCNF can increase substantially with an increase in temperature (Figure 7(a)), at a rate of about 248 kg C ha−1 yr−1 per ˚C on average.
Figure 7. (a) Impact of air temperature on annual mean aboveground biomass (kg C ha−1 yr−1) in Nez perce clearwater national forest within 100 years; (b) Impact of precipitation on annual mean aboveground biomass (kg C ha−1 yr−1) in Nez perce clearwater national forest within 100 years.
The impact of temperature on C sequestration in this forest is completely different from that in tropical semidry forests where the C sequestration can decrease with an increase in temperature [ 32 ]. This is due mainly to lower temperature in this forest (about 5.1˚C averaged from 1981 to 2013) than that in the tropical dry forests, such as the forests in the Yucatan Peninsula of Mexico (the mean of 26.5˚C for the period 1981-2012) [ 32 ]. Temperature rise can increase the length of growing seasons in this forest such that the biomass can increase with an increment in mean air temperature in this forest.
Biomass can substantially decrease small with an increase in precipitation in this forest (Figure 7(b); P < 0.001), at a rate of about 98.5 kg C ha−1 yr−1 per 100 mm. However, there are some divergences in the impact of precipitation on C stocks in this forest among the vegetation types. The biomass in the forestlands dominated by Ponderosa pine, subalpine fir and subalpine spruce does not decrease substantially with an increase in precipitation (R = 0.107, n = 204, P > 0.1 for subalpine fir and spruce; R = 0.145, n = 51, P > 0.2 for ponderosa pine). Compared to the tropical dry forests where the ratio of precipitation to potential evapotranspiration (RPPE) is less than 1.0, such as the forests in the Yucatan Peninsula in Mexico where the ratio is only about 0.8 [ 38 ], the biomass in this forest does not increase with an increment in precipitation due to the RPPE ≥ 1.0 in this forest.
The results from model evaluation indicate that the model Forest-DNDC can be applied for assessing C sequestration in NPCNF in Idaho of USA with a high model performance efficiency (E > 0.8) and a small error between observation and simulation. The slightly overestimated live biomass and small underestimated deadwood for this forest might be related to unavailable disturbance information for the time period before 1991 that was assumed as same as the observed period from 1991 to 2011.
There are some differences in C sequestration rates among species in the forest. The relationship between NPP and stand age can alter from quadratic to sextic polynomial. Although AGB in this forest can increase at a rate of quartic polynomial with an increment in stand age, there are substantial divergences in the growth rate among the species. The subalpine species have lower growth rates, but they can continuously increase for a longer time than other species.
Temperature rise can increase biomass in this forest. This is mainly related to low air temperature (long-term mean of 5.1˚C) in this forest. Accordingly, an increase in temperature may increase the length of growing seasons such that forest production can increase with an increase in temperature in this forest. However, precipitation increase cannot produce more biomass in this forest because this forest has adequate precipitation (long-term mean of 1100 mm).