Imputing forest carbon stock estimates from inventory plots to a nationally continuous coverage.

5.6 MB·English

by  Wilson, Barry Tyler

Wilsonetal.CarbonBalanceandManagement2013,8:1
http://www.cbmjournal.com/content/8/1/1
RESEARCH Open Access
Imputing forest carbon stock estimates from
inventory plots to a nationally continuous
coverage
Barry Tyler Wilson*, Christopher W Woodall and Douglas M Griffith
Abstract
TheU.S.hasbeenprovidingnational-scaleestimatesofforestcarbon(C)stocksandstockchangetomeetUnited
NationsFrameworkConventiononClimateChange(UNFCCC)reportingrequirementsforyears.Althoughthese
currentlyareprovidedasnationalestimatesbypoolandyeartomeetgreenhousegasmonitoringrequirements,
thereisgrowingneedtodisaggregatetheseestimatestofinerscalestoenablestrategicforestmanagementand
monitoringactivitiesfocusedonvariousecosystemservicessuchasCstorageenhancement.Throughapplicationof
anearest-neighborimputationapproach,spatiallyextantestimatesofforestCdensityweredevelopedforthe
conterminousU.S.usingtheU.S.’sannualforestinventory.Resultssuggestthatanexistingforestinventoryplot
imputationapproachcanbereadilymodifiedtoproviderastermapsofCdensityacrossarangeofpools(e.g.,live
treetosoilorganiccarbon)andspatialscales(e.g.,sub-countytobiome).Comparisonsamongimputedmaps
indicatestrongregionaldifferencesacrossCpools.TheCdensityofpoolscloselyrelatedtodetritalinput(e.g.,dead
wood)isoftenhighestinforestssufferingfromrecentmortalityeventssuchasthoseinthenorthernRocky
Mountains(e.g.,beetleinfestations).Incontrast,livetreecarbondensityisoftenhighestonthehighestqualityforest
sitessuchasthosefoundinthePacificNorthwest.Validationresultssuggeststrongagreementbetweenthe
estimatesproducedfromtheforestinventoryplotsandthosefromtheimputedmaps,particularlywhentheCpool
iscloselyassociatedwiththeimputationmodel(e.g.,abovegroundlivebiomassandlivetreebasalarea),with
weakeragreementfordetritalpools(e.g.,standingdeadtrees).Forestinventoryimputedplotmapsprovidean
efficientandflexibleapproachtomonitoringdiverseCpoolsatnational(e.g.,UNFCCC)andregionalscales(e.g.,
ReducingEmissionsfromDeforestationandForestDegradationprojects)whileallowingtimelyincorporationof
empiricaldata(e.g.,annualforestinventory).
Keywords:Forest,Carbondensity,Imputation,UnitedStates,Forestinventory,Rastermaps
Background at sub-national scales including states (e.g., California) and
Forest ecosystems represent the largest terrestrial carbon ownerships (e.g, National Forest System climate change
(C) sink on earth [1,2], such that the United Nations scorecard). Forest C stocks in the U.S. are estimated using
Framework Convention on Climate Change [3] has recog- data from the national forest inventory conducted by the
nized their management as an effective strategy for offset- USDA Forest Service, Forest Inventory and Analysis (FIA)
ting greenhouse gas (GHG) emissions [4,5]. As part of the program [7]. Broad forest ecosystem components (e.g.,
Convention,theU.S.hasbeensubmittingnationalreports, aboveground live biomass) have been delineated to
the National Greenhouse Gas Inventory (NGHGI), detail- generalize C stocks to meet international reporting agree-
ing emissions and removals of GHGs [3] on an annual mentspursuanttorefiningunderstandingofglobalcarbon
basis for many years [6]. In addition to international cycling [2,3]. Carbon estimates for the ecosystem compo-
reportingrequirements,GHGbudgetsarebeingdeveloped nents of forest floor (inclusive of litter, fine woody debris,
and humic soil horizons), down dead wood, belowground
*Correspondence:[email protected] (BG) biomass, and soil organic matter are calculated by
USDAForestService,NorthernResearchStation,ForestInventoryand FIA using models based on geographic area, forest type,
Analysis,St.Paul,MN55108,USA
©2013Wilsonetal.;licenseeBioMedCentralLtd.ThisisanOpenAccessarticledistributedunderthetermsoftheCreative
CommonsAttributionLicense(http://creativecommons.org/licenses/by/2.0),whichpermitsunrestricteduse,distribution,and
reproductioninanymedium,providedtheoriginalworkisproperlycited.
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page2of15
http://www.cbmjournal.com/content/8/1/1
and, in some cases, stand age [6,8]. Estimates of above- accuracyof C maps may often be of equal importance as
ground (AG) standing live and dead tree C stocks are the need for spatial accuracy. Woodall et al. [12] found
based on biomass estimates obtained from inventory tree that actual standing dead tree C stocks were often sig-
data[6,9].AlthoughforestCstockestimates,suchasthose nificantly different than those modeled for the same in-
from FIA, are readily available at national and regional ventory plots. Despite the measurement/model error
scales [6,7], there is increasing interest in disaggregating associated with annual forest inventory programs, the
theselarge-scalenumericalestimatesintomapsofcontinu- temporally dynamic nature of forest ecosystems (e.g.,
ous estimates to enable strategic forest management and wildfires and wind events) necessitates the incorporation
monitoring activities geared toward offsetting GHG emis- of annual data into map products employed by scientists
sions[10]andadvancingCdynamicsresearch. andstakeholdersalike.
Secondary to the need for spatially continuous forest Wilson et al. [13] developed a methodology (hereafter
C maps, numerous constituents (e.g., managers, policy referred to as Phenological Gradient Nearest Neighbor,
makers, and scientists, forest analysts) require an effi- or PGNN, for convenience) for producing maps of tree
cient methodology for incorporating annual monitoring species occurrence and relative abundance over large
information into C maps. Sophisticated approaches to areas by utilizing information collected on FIA field
mapping forest C stocks may provide robust estimates plots in conjunction with 250 m pixel resolution raster
of stocks [11], but lack the flexibility to rapidly incorpor- data in a k-nearest neighbor (kNN) imputation frame-
ate annual monitoring information. As numerous forest work. The PGNN approach builds upon the Gradient
C pools may change on annual time steps, especially in Nearest Neighbor (GNN) work of Ohmann and Gregory
response to stochastic disturbance events, temporal [14], who integrated nearest-neighbor imputation of FIA
Figure1Totalforestecosystemcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.Includesabove-and
belowgroundlivetrees,downeddeadwood,forestfloor,soilorganiccarbon,standingdeadtrees,understoryabove-andbelowgroundpools.
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page3of15
http://www.cbmjournal.com/content/8/1/1
plots with ecological ordination via canonical corres- to a spatially continuous raster grid in order to produce
pondence analysis (CCA). PGNN is best described as a mapped estimates of the conterminous U.S.’s forest C
hybrid of the kNN and GNN approaches, since it also densitywiththesespecificobjectives:1)toproduceandin-
makes use of CCA but utilizes k nearest neighbors dur- terpret maps of forest carbon density by individual pools
ing imputation rather than only a single neighbor. An- andcombinationsthereof(totalforestecosystemCdensity,
other distinguishing characteristic is that it utilizes live tree AG, live tree BG, live understory AG and BG,
vegetation phenology information derived from multi- standing dead tree AG, downed dead wood, forest floor,
temporal satellite imagery, as well as climate, topo- soil organic carbon, and the pool that has the highest
graphic, and ecoregion data compiled at a 250 m pixel proportion of total forest ecosystem C density); 2) to con-
resolution. One of the most attractive features of this duct validation of the C mapping approach by comparing
approach is the efficiency with which a plot identifica- map-based and field plot-based estimates using a variety
tion map can be produced at the national scale. In of metrics; and 3) to suggest future research directions
other words, every pixel is assigned a forest inventory and applications.
plot label as well as the attributes of the labeled plot’s
nearest neighbors, as defined by the CCA model. In the Results
case of forest C accounting, every pixel could be TheimputedrastermapsoftotalforestecosystemCstocks
assigned C stock estimates in a rapid fashion on an an- (i.e., sum of all pools) suggest a rather disparate distribu-
nual time-step. tionoflargetotalCstocksacrosstheU.S.(Figure1).While
Given the need for C maps at the national scale and the most forested areas of the U.S. have moderate C stock
possibleapplicationofPGNN,thegoalofthisstudywasto density (< 100 Mg/ha) (e.g., lower elevations of the Rocky
apply PGNN for imputing national forest inventory plots Mountains, Central, and Plains states), there are other
Figure2Livetreeabovegroundcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page4of15
http://www.cbmjournal.com/content/8/1/1
areasthathaveCdensitiesinexcessof200Mg/hasuchas density exceeding 2 Mg/ha. The highest downed dead C
the Pacific Northwest and the upper Great Lakes. As total stock densities (> 12 Mg/ha) are almost exclusively
Cstocksarecomprisedofdiverseforestecosystemcompo- found in the Pacific Northwest and West Coast/Sierra
nents, examining the distribution of C stock density by Nevada (Figure 6). The detrital components of forest
component can refine understanding of C density dynam- floor and SOC have spatial distributions fundamentally
icsatanationalscale. different from woody biomass C stock distributions
Live AG and BG (Figures 2 and 3) C stock density is (Figures 7 and 8). The highest C stock densities for
highest in the Pacific Northwest, northwest California, forest floor (> 15 Mg/ha) are found in the Pacific
northern Rockies, and Appalachian Mountains. The Northwest, California, Rocky Mountains, upper penin-
highest live AG C stock density often exceeds 80 Mg/ha. sula of MI, and New England. The highest C stock
As BG C stocks are modeled as a function of AG C densities for SOC (> 80 Mg/ha) are found in the
stocks, their spatial distributions are closely aligned. Live upper Lake States, Pacific Northwest, northern New
understory AG and BG C density follow spatial patterns England, and coastal areas of the Southeast. In order
in allocation similar to live tree distributions, albeit at a to better appraise areas subject to varying C stock dy-
much lower density (< 1 Mg/ha) (Figure 4). Standing namics, each pixel was assigned to one of three cat-
dead tree C stock densities are highest (> 8 Mg/ha) in egories, indicating where it had the largest proportion
the Olympic Mountains, Cascade Range, and North and of its total C stocks apportioned: 1) live biomass (live
Central Rocky Mountains (Figure 5). In comparison to tree and understory AG and BG), 2) SOC, and 3) dead
the western U.S., eastern standing dead tree C stock wood and forest floor (Figure 9). Live biomass is the
density is minimal with only the Adirondacks and iso- dominantCstockalongtheWestCoastandAppalachian
lated areas of the Appalachian Mountains having a stock Mountains. In contrast, SOC is the dominant C stock
Figure3Livetreebelowgroundcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page5of15
http://www.cbmjournal.com/content/8/1/1
Figure4Liveunderstoryaboveandbelowgroundcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
along Southeastern coastal areas, New England, Great 50km(216,500ha).Thedifferencesappeartobedistrib-
Lakes, and Great Plain’s forests. The dead wood and for- uted in a rather spatially unbiased manner across the
estfloorisdominantinareasoftheRockyMountains. conterminous U.S., although there is a tendency for the
Validation metrics suggest good agreement between model to overestimate forest C at the edges of the
map-based and field plot-based estimates of C density forested extent. The same phenomenon was present and
across pools and spatial scales (Table 1). The strongest discussed in the earlier study of species relative abun-
agreement according to all three validation metrics is at dance [13]. This is most likely an effect of the spatial
the coarsest spatial scale (200 km) with slight reductions mismatch between pixels and plots, as well as a mis-
in agreement statistics down to the finer spatial scale match between the “forest” stratum used during imput-
(50 km). At 25 km, most pools demonstrate a more sub- ation and the FIA definition of forest land, that includes
stantial drop in agreement, albeit most statistics still in- a component based on land use not readily detected
dicate strong agreement (e.g., agreement coefficient usingremote sensingdata.
above 0.90). The one exception is the standing dead tree
C pool which had good agreement at the spatial scale of Discussion
100-200 km, but demonstrated a marked decline in This study demonstrated that a spatially explicit imput-
agreement down to thefinest spatial scale of25 km (e.g., ation approach may be applied to a standard forest in-
agreement coefficient=0.67). The distribution curves ventory to efficiently produce continuous maps of forest
and CI maps (Figures 10, 11, 12, 13, 14, 15, 16 and 17) C stock estimates in a timely manner. Across spatial
reinforce the fit statistics: there is fairly robust agree- scales ranging from 25 to 200 km, imputed C stock esti-
ment between the map-based and field plot-based esti- mates closely matched those derived from the forest in-
mates of C density by pool at the finer spatial scale of ventory data that serve as the basis for the U.S.’s
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page6of15
http://www.cbmjournal.com/content/8/1/1
Figure5Standingdeadtreeabovegroundcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
NGHGI. Thus, the opportunity exists to down-scale a highlighted the need to provide temporally continuous
NGHGI to finer scales (e.g., sub-county scale). As vari- monitoring of forest C stocks. Whether the PGNN ap-
ous studies [5,10] have recently highlighted the role that proach or alternative is used, spatially explicit imput-
forest management may play in mitigating climate ation approaches enable the use of annual forest
change, the accurate carbon assessment of “project- inventories to inform real-time and real-world C man-
scale”forest managementactivitiesisparamounttofore- agement situations.
casting potential mitigation benefits, if any. A project Although an imputed forest C map may provide a
scale C inventory that is consistent with regional and reasonable down-scaling of a NGHGI, the situation
NGHGI’smayallowforrobustverification,adesirableat- remains that most NGHGI’s have high variability at
tribute of forest C projects [15]. In addition, consistency fine scales such that the statistical power to detect
in monitoring forest C across spatial scales may be inte- stand-level change is limited [20]. Therefore, while
graltofullyunderstandingthedynamicsofforestCfrom an imputed C estimate for an individual pixel may
stands to ecosystems [16]. Beyond the spatial scale, such be consistent with a NGHGI, the level of uncertainty
imputed maps may be rapidly updated annually to in- associated with that one pixel will be very high. It
corporate recent disturbance information. It is expected may not matter how consistent a map is with a
that large scale disturbance events may increasingly NGHGI, if the NGHGI itself contains tremendous
affect forest C stocks in the face of climate change [17]. uncertainty. Another limitation to imputing NGHGI
Disturbance events may occur on annual time steps data to maps is that dedicated analytical staff is
such that forests become net emitters of C [18] over a needed to produce these outputs on an annual time
relatively short period of time. Recent wildfires and in- step. Indeed, if the value in such an approach is its
sect outbreaks [19] in western North America have sensitivity to recent disturbance events, then likewise
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page7of15
http://www.cbmjournal.com/content/8/1/1
Figure6Downeddeadwoodcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
forest analysts will need to develop and apply imput- the map in this study (Figure 7) highlights the reduced
ation models within the same time frame. Finally, an uncertaintythatcouldberealizedbybothadoptingem-
important component of the PGNN approach used in pirical measurements of C pools within a NGHGI and
this study is that the live tree attribute of basal area adapting imputation models to fit those unique ecosys-
was a central dependent variable [13]. The imputed temcomponents.
maps of standing dead tree C density had the poorest Can imputed forest C maps inform our knowledge
agreement with the empirical inventory present in regarding the dynamics and attributes of C across the
the NGHGI. It has already been demonstrated in the U.S.? Perhaps the dynamics of forest C is best illustrated
U.S.’s NGHGI that detrital C models based on live by Figure 9. It is obvious that the predominant forest C
tree attributes may substantially over/under-estimate pool varies by ecosystem across the U.S.. While at high
actual C stocks at finer scales (e.g., plot-level) [12]. latitudes or in coastal/wetland areas it may be the SOC
Thismaylikewiseoccurwithspatialimputationmodels and forest floor pools that require the most attention
and it is suggested that future research explore alterna- when itcomesto managementandmonitoring,itcanbe
tive imputation models for the non-living C stocks in the AG biomass (whether dead or alive) that should be a
forest ecosystems. There was reasonable agreement focal point in most other areas. Most telling was the
between the imputed estimates of SOC and the forest dominance of detrital forest C pools in most areas of the
floor C densities compared to the NGHGI because Intermountain West. The value of imputed forest C
the NGHGI currently uses models based partially on maps may be beyondmonitoring C monitoring efforts at
live tree attributes to estimate these stocks [6]. How- scales ranging from national to sub-county, rather it
ever, a comparison between the first empirical inven- may be in identifying emerging research areas and eco-
tory of forest floor C stocks across the U.S. [21] and logical “hotspots” [21] such as areas where forest C
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page8of15
http://www.cbmjournal.com/content/8/1/1
Figure7Forestfloorcarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
stocks may be responding to climate change events or Methods
how detrital andsoilorganicCstocksmayberelated. Data
The FIA program is the primary source for information
Conclusions about the extent, condition, status, and trends of forest
Down-scaling forest NGHGI’s to finer scales (i.e. project resources intheUnitedStates[22].FIA appliesanation-
level) is needed to provide project verification that is re- ally consistent sampling protocol using a systematic de-
gionally consistent while at the same time refining the sign covering all ownerships across the U.S., at a base
science of forest C monitoring. A map-based imputation national sample intensity of one plot per 2,428 ha. Land
approach, such as the PGNN technique applied in this area is stratified using aerial photography or classified
study, affords an efficient and timely method for produ- satellite imagery to increase the precision of estimates.
cing spatially continuous maps of diverse forest C pools Remotely sensed data may also be used to determine if
thatare consistentwithaNGHGI.Theuncertainty asso- plot locations are forested and should be measured in
ciated with each imputed pixel is dependent on not only the field. FIA defines forested land as areas that have at
the imputation model, but also on the models that least 10 percent tree canopy cover, are at least 0.4 ha in
underlie the NGHGI. We suggest that further research size, and are at least 36.6 m wide [23]. FIA inventory
in both modeling areas be undertaken to refine forest C plotsconsist offour,7.32-m fixed-radius subplotsspaced
maps in the future. Until such refinements occur, the 36.6 m apart in a triangular arrangement with one sub-
maps produced in this study can only provide rough plot in the center [8,24]. All trees (live and standing
guidance for forest C projects while highlighting the re- dead) with a diameter at breast height (dbh) of at least
gional differences invarious C poolsand their associated 12.7 cm, are inventoried on forested subplots. Within
dynamicsthatinturnmayguidefutureresearch. each subplot, a 2.07 m microplot offset 3.66 m from
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page9of15
http://www.cbmjournal.com/content/8/1/1
Figure8Soilorganiccarbondensityimputedfromforestinventoryplots,conterminousU.S.,2000-2009.
subplot center is established where only live trees with a others have chosen to intensify their sampling intensity
dbhbetween 2.5and12.7cmareinventoried. by 2-3 times the base intensity. FIA field data, with ap-
The field data for this study were taken entirely from proximate plot locations, are freely available for down-
the FIA database [8] using the most recent annual col- load fromthe program’swebsite [25].
lection (what FIA refers to as an “evaluation”) of forest These field data were used in conjunction with data
inventory plots available at the time of the study for the extracted at each plot location from a 250 m pixel reso-
conterminous 48 states (N.B., annual plots from western lution raster stack. This predictor dataset included vege-
Oklahoma, New Mexico, and Wyoming were not avail- tation phenology information derived from a time series
ableatthetimeofthisstudy).Thedata collectionperiod of vegetation indices based on MODIS satellite imagery
for the annual inventories conducted in most of the (2002-2008), mean monthly climate characteristics from
states used in this study was initiated since 2000 and ran the Daymet climatological model of interpolated climate
through 2009, with a full cycle of plots being collected station observations (1980-1997) [26], topographic metrics
over a5-year period in the East (roughly2005-2009) and from the Elevation Derivatives for National Applications
a 10-year period in the West (roughly 2000-2009). The digitalelevationmodel,andOmernik’sLevelIIIecoregions
collection of annual plots used for each state in this (orecologicalzones)[27].Foramorecompletedescription
study contained all plots associated with that “evalu- ofthesedatasetsandhowtheywereusedinthestudy,see
ation” from the FIA database. This includes non- Wilsonetal.[13].
forested plots and, by definition, only the data collected
from the most recent observation of each plot. Sample Plot-levelcarbonestimates
intensities vary by state, since some have not yet com- Plot-level estimates of forest C stocks are a combination
pleted data collection for a full cycle of plots, while of empirically measured tree/site attributes combined
Wilsonetal.CarbonBalanceandManagement2013,8:1 Page10of15
http://www.cbmjournal.com/content/8/1/1
Figure9Majorforestcarbonpoolswiththepluralityoftotalforestcarbonstockforeachpixelimputedfromforestinventoryplots,
conterminousU.S.,2000-2009.Majorpoolsare:1)livingbiomass(aboveground,belowground,andunderstory),2)deadwoodandforestfloor
(includingstandingdead,downdead,andforestfloor),and3)soilorganiccarbon.
Table1ValidationmetricsbyscaleandforestCstock
Metric Scale Total Above-ground Below-ground Understory Standingdead Downdead ForestFloor Soilorganic
Agreement 200km 0.9911 0.9942 0.9944 0.9909 0.9399 0.9932 0.9947 0.9828
Coefficienta
100km 0.9943 0.9926 0.9925 0.9866 0.9379 0.9924 0.9943 0.9933
50km 0.9866 0.9785 0.9785 0.9697 0.8635 0.9844 0.9850 0.9850
25km 0.9610 0.9295 0.9300 0.9220 0.6707 0.9534 0.9605 0.9639
KSStatisticb 200km 0.0364 0.0424 0.0485 0.0364 0.0485 0.0364 0.0364 0.0424
100km 0.0274 0.0316 0.0316 0.0247 0.0521 0.0288 0.0233 0.0288
50km 0.0866 0.0882 0.0882 0.0863 0.1135 0.0866 0.0863 0.0863
25km 0.1246 0.1260 0.1260 0.1244 0.1591 0.1244 0.1244 0.1244
RMASlopec 200km 1.0073 1.0039 1.0050 0.9941 1.1515 1.0122 1.0003 1.0089
100km 1.0012 1.0025 1.0027 0.9995 1.0466 1.0050 1.0059 0.9987
50km 1.0102 1.0171 1.0174 1.0140 1.1405 1.0151 1.0114 1.0073
25km 1.0248 1.0405 1.0411 1.0403 1.2158 1.0319 1.0271 1.0206
aAgreementCoefficient(largervaluesindicatebetteragreement,min=0,max=1).
bKolmogorov–SmirnovStatistic(smallervaluesindicatebetteragreement,min=0,max=unbounded).
cReducedMajorAxisSlope(valuescloserto1indicatebetteragreement,min,max=unbounded).