High-resolution Mapping of Forest Carbon Stocks in the Colombian Amazon

Introduction Conclusions References


Introduction
Tropical forests store roughly 475 billion tons of carbon (Pan et al., 2011), so retaining this carbon through conservation and increasing its stock through management activities that promote forest growth will play a major role in curbing a principal driver of climate change (Angelsen, 2008).Acknowledging this opportunity, the United Nations Framework Convention on Climate Change agreed to encourage reductions in greenhouse gas emissions from forests via the program for Reducing Emissions from Deforestation and Forest Degradation (REDD+) (UNFCCC, 2009).However, a technical barrier to REDD+ rests in monitoring carbon stocks and emissions.Although guidelines and discussion on the topic abound, few studies have delivered synoptic scale, high-resolution estimates of forest carbon stocks with spatially explicit uncertainty.

G. P. Asner: High-resolution mapping of forest carbon stocks in the Colombian Amazon
Sub-national or jurisdiction-scale mapping approaches, often in the multi-million hectare range, are fast advancing in the REDD+ development process, and they are a key steppingstone toward international REDD+ implementation (Herold and Skutsch, 2011).A case-in-point is the Governors' Climate and Forests Task Force (GCF; http://www.gctaskforce.org),which links major jurisdictions in Brazil, Indonesia, México, Nigeria, Peru and the United States to advance the role of forests in climate change mitigation.GCF could be the first international program to support compliance regulations with carbon emission offsets and sequestration, and this will require high-resolution mapping of forest carbon stocks and emissions to achieve its goals (Asner, 2011).Other regional programs are underway, mostly in voluntary or pre-compliance contexts (Parker et al., 2008) and most rely on general maps of forest cover and plot-based estimates of ACD extrapolated to large regional scales.These approaches result in high uncertainty of carbon stock estimates (Asner, 2009;Avitabile et al., 2011); hence, the cascading uncertainty in estimated carbon emissions are then very large (Pelletier et al., 2011).
Many approaches have been described for estimating aboveground carbon density (ACD; units of Mg C ha −1 ) at different scales in tropical forests.Plot-level inventories provide localized information on ACD, usually at a scale of 1 ha or less (e.g., Malhi et al., 2006).Developing large-scale, high-resolution ACD maps requires estimation of forest carbon stocks as accurately as with field plots, but over large gradients of climate, topography, hydrology, soils, and biological diversity (Goetz et al., 2009).The measurements must also resolve highly variable effects of land use on carbon stocks.Both natural gradients and land use impart a profound influence on ACD levels (Tian et al., 2000;Clark et al., 2002;Saatchi et al., 2007;Mascaro et al., 2011a), which are non-randomly distributed across the landscape (Loarie et al., 2009).As a result, ACD must be resolved spatially to support REDD+ emissions monitoring.
A major contributor to solving the carbon mapping challenge is airborne light detection and ranging (LiDAR), a technology that images forest canopies in three-dimensions using emitted laser light pulses (Lefsky et al., 2002b).LiDAR can be used to examine forest architecture in fine detail (Asner et al., 2008;García et al., 2010).When combined with field plots, LiDAR provides high-resolution, spatially contiguous estimates of ACD (e.g., Drake et al., 2002;Gonzalez et al., 2010), which can be used to map thousands of hectares of forest per day to quantify environmental controls over forest carbon storage (Asner et al., 2009a).Airborne LiDARbased estimates of tropical forest ACD are improving, with per-hectare errors recently becoming indistinguishable from those derived in field plots (Mascaro et al., 2011b).
To scale up to landscape level, LiDAR studies have often been closely tied to field calibration plots distributed throughout the mapping coverage.However, the logistical and cost burden of establishing an extensive plot network may limit the utility of LiDAR for carbon mapping, particularly in forests that remain very remote, either by distance or by difficult terrain.To address this problem, Asner et al. (2012b) recently developed a "universal" equation to estimate tropical ACD from airborne LiDAR.Based on data collected in Panamá, Hawaii, Perú and Madagascar, spanning a wide range of forest ages and floristic types, the universal LiDAR equation was developed to predict ACD with relatively high precision (r 2 = 0.80) and accuracy (RMSE = 27.6 Mg C ha −1 ).The equation provides estimates of ACD that are comparable in predictive power to locallycalibrated models, yet it relies only on limited basal area and wood density information for a given region, rather than traditional plot inventories.This approach has the potential to greatly reduce the time required to calibrate airborne LiDAR data, however it requires testing in new regions.
Despite the accuracy and extensive mapping capability of airborne LiDAR, it too reaches geographic limits due to cost and logistics, so methods are required to extend LiDARbased ACD samples to even larger scales.Two general approaches have emerged.One is high-resolution stratification of a region by vegetation type, topography and other environmental datasets, along with high-resolution mapping of vegetation losses and gains from deforestation, degradation and land abandonment (Asner, 2009;Helmer et al., 2009).Following stratification, the region is sampled to develop carbon density statistics for each stratum.LiDAR-assisted mapping of stratified regions can produce robust ACD statistics, which has proven useful in a wide range of forest types (Asner et al., 2010(Asner et al., , 2011(Asner et al., , 2012a)), but highly stratified maps can be difficult to obtain for many tropical regions (Herold and Skutsch, 2009;Pelletier et al., 2011).An alternative or decision-tree approach is to upscale field and LiDAR-based ACD estimates with spatially contiguous, regional correlates derived from satellite imagery (Baccini et al., 2008;Saatchi et al., 2011).Such a regression-based approach is elegant, and may allow for more rapid upscaling of field and LiDAR-based ACD measurements.However, regression approaches may also miss local or sub-regional controls over carbon stocks that can be resolved using stratification.
These two issues -LiDAR applicability with few field plots, and upscaling of LiDAR data to larger regions -remain critically important challenges to making high-resolution carbon stock and emissions monitoring possible.The Colombian Amazon is enigmatic of these challenges, where lowland to montane forests remain virtually unexplored in terms of carbon stocks and their environmental controls.The region not only has potential for carbon offsets and climate change mitigation work, but surveys indicate it is a major part of the western Amazon biodiversity hotspot.However, vegetation maps of the Colombian Amazon currently rely on coarse biological information (Forero, 1988;Armenteras et al., 2006), thereby lacking the definition needed for stratification and subsequent LiDAR sampling.We therefore developed a new top-down, high-resolution analytical approach for determining forest carbon stocks in this region.In doing so, we address three questions pertinent to carbon mapping efforts in remote, inaccessible tropical forests: (i) Using available satellite imagery and airborne LiDAR sampling, what are the principal determinants of aboveground carbon density detectable throughout the region?(ii) Despite limits to acquiring field inventory data on the ground, what are the estimated uncertainties associated with applying the universal LiDAR approach to the Colombian Amazon?(iii) What are the uncertainties associated with the stratification and regression approaches, and what are their advantages and disadvantages?
2 Methods and materials

Study area
The study region covers 16 561 695 ha (> 40 %) of the Colombian Amazon (Fig. 1), stretching from the Andean foothills in the west to the Brazilian border in the east, with unknown variation in natural forest carbon storage and very limited documentation of deforestation and forest degradation.The region combines a vast, forested sedimentary plain of Tertiary and Quaternary age, with large Paleozoic sandstone plateaus and remnant Precambrian surfaces derived from the Guiana Shield (Duivenvoorden and Duque, 2010).Andean-derived sediments within the plain are nutrient rich, while Guyana Shield soils are highly leached and nutrient poor (Quesada et al., 2012).Lowporosity basement rock within portions of the plain are associated with extensive swamps and inundation.Mean annual temperature is ∼ 25 • C, and mean annual precipitation ranges from ∼ 2000 mm yr −1 in the northwest to more than 3000 mm yr −1 in the southeast portions of the study region.Tree species diversity is thought to be among the highest in the Amazon (Duque et al., 2009).
The region is designated as a REDD+ pilot project area of the Colombian Institute for Hydrological, Meteorological, and Environmental Studies (IDEAM), stretching from the northwestern departments of Meta and Caquetá to the remote lowlands of Vaupés and Amazonas.The area is largely inaccessible due to a lack of roads and navigable rivers, and ongoing security issues prevent the widespread use of forest inventory plots.

Preliminary stratification
To guide LiDAR sampling, we performed pre-flight landcover stratification using a decision tree with input variables known to influence carbon stocks (Fig. 2a).With the CLASlite forest monitoring system (Asner et al., 2009b), we mapped forest cover for the year 2010 using 16 Landsat TM and ETM+ images at 30 m resolution.With 46 additional TM and ETM+ images from 1990, 2000, and 2005, we used CLASlite to map regrowth following historic deforestation and degradation with the techniques described in Asner et al. (2010).Finally, we segmented the region into elevation bands of 50 m, derived from the NASA Shuttle Radar Topography Mission (SRTM) digital elevation model (DEM) at 90 m spatial resolution (Fig. 2b).The SRTM DEM and other inputs are shown in the Supplement (Fig. S1).

LiDAR sampling
With the preliminary stratification map, we planned airborne LiDAR sampling to achieve at least 1 % coverage of each stratum, while targeting additional terrain variation within the forested classes that could contribute to regional variation in ACD.The LiDAR data were collected in January 2011 using the Carnegie Airborne Observatory (CAO) Alpha system (http://cao.ciw.edu).The CAO Alpha LiDAR (see Asner et al., 2007) was operated at 2000 m above ground level with 1.12-m spot spacing, 30-degree field of view, beam divergence customized to 0.56 mrad, and 50-kHz pulse repetition frequency, for which the aircraft maintained a ground speed ≤ 95 knots.With these flying parameters, CAO collected in continuous laser coverage without gaps between laser spots on the ground.In addition, all flights were planned with 100 % repeat coverage (50 % overlap of each swath to each adjacent swath) and therefore LiDAR pulse density averaged 2 points per 1.12-m spot.
LiDAR sampling totaled 465 622 ha (i.e., cloud-free, usable data), achieving coverage of 2.8 % of the region in www.biogeosciences.net/9/2683/2012/Biogeosciences, 9, 2683-2696, 2012 38 areas ranging in size from 9000 to 30 000 ha (Fig. 1).Due to security concerns within the study region, all flight operations were conducted at night -an option afforded by airborne LiDAR.Canopy three-dimensional structure, as detected by Li-DAR, was analyzed by binning discrete LiDAR returns into volumetric pixels (voxels) of 5 m spatial resolution and 1 m vertical resolution, yielding histograms representing the vertical distribution of vegetation in each 5 × 5 m spatial cell.These data were further reduced to mean canopy profile height, or MCH, which is the volumetric vertical center of the canopy (as opposed to top-of-canopy height).LiDAR MCH has been used in a large number of studies to estimate ACD with demonstrably high precision and accuracy (Lefsky et al., 2002a;Asner et al., 2010Asner et al., , 2011;;Mascaro et al., 2011a, b).

LiDAR-to-ACD conversion
This study is the first to apply a streamlined approach to convert LiDAR MCH measurements to ACD in tropical forests, which we summarize here.Asner et al. (2012b) calibrated a single model based on 482 field plots spread across four distinct tropical regions (Hawaii, Madagascar, Panamá, Perú).The architecture of the universal model follows basic tree allometry studies (Chave et al., 2005;Niklas, 2006), demonstrating that dry tree biomass (and thus carbon content, which is ∼ 48 % of biomass by weight; Martin and Thomas, 2011), can be estimated as a function of tree diameter, height, and wood density.Thus, the universal model relied on plot-level, rather than tree-level, inputs of field-measured diameters in the form of basal area (BA), wood density (a weighted average based on each tree's basal area), and LiDAR-derived MCH to estimate ACD in units of Mg C ha −1 .Upon calibration, the model was ACD = 2.04 MCH 0.436 BA 0.946 WD 0.912  (1) which explained 95 % of the variation in ACD across the four tropical regions.Asner et al. (2012b) further demonstrated that plot-level BA could be estimated from LiDAR according to a linear MCH-to-BA conversion specific to each region, hereafter termed the "stocking coefficient" (SC), and that regional WD could be substituted for plot-level WD.When these regional substitutions were made -that is, two numerical constants in place of exhaustive inventories of tree diameters, heights and wood densities -the model explained 81 % of variation in ACD.In this way, the model could be adjusted for a new study by simply inputting SC and WD for any tropical region.
Here, we used eleven 0.28-ha field plots to estimate the two regional input constants needed for the calibration, and subsequently validated this method using traditional inventory techniques (described below).We derived a regional SC of 1.52 (Fig. S2a), and thus we substituted 1.52 • MCH for the BA term in the LiDAR calibration equation (Eq.1).We determined that WD was not significantly related to MCH, and therefore we applied the basal-area weighted mean of 0.61 (Fig. S2b).For comparison, the SC and WD values for the southwestern Peruvian Amazon were 1.53 and 0.56, respectively (Asner et al., 2012b).With these inputs, Eq. ( 1) was simplified to a regionally calibrated LiDAR equation: (2)

Model validation
We validated the universal LiDAR calibration equation using traditional forest inventory techniques with allometric regression equations (Table S1).Although we used the same plots, this validation was critical to determine whether two simple regional constants (i.e., SC and WD) could substitute for traditional field inventory.For dead trees and palms, we utilized growthform-specific equations.Lianas were considered, but none were detected over the minimum size class (10 cm dbh).For all other trees, we utilized a general moist forest model of Chave et al. (2005).We corrected for local height variation by directly measuring the heights of the largest trees in all plots (> 50 cm dbh) and additional trees spanning a range of stem diameters using a laser hypsometer (Impulse-200, LaserTech Inc., USA).For the remaining trees, we produced a model relating height and diameter using maximum likelihood analysis (Fig. S3; R Development Core Team, 2011).Wood density values were assigned based on genus-(33 %) or family-level (42 %) identifications according to Chave et al. (2009), and are detailed in Table S2.These values were determined by averaging all listings at the genus-and familylevel within the database.In the absence of such an identification (25 %), a regional estimate of 0.58 was applied (ter Steege et al., 2006).Plot centers were determined using a global positioning system (GPS) with differential correction (Leica GS-50, Leica Geosystems Inc., Switzerland), which provided < 1 m positional uncertainty in most cases.

Regional upscaling based on stratification
To upscale LiDAR-derived ACD, we relied on the notion that nested sources of variation in carbon stocks exist throughout most tropical regions.These may occur, for example, due to localized gradients in riparian vegetation, mesoscale variation in canopy cover, and large-scale controls such as elevation, with overlapping human influences at all scales.Stratification allows for parsing a region based on these factors, and the finer the strata are delineated; the greater the possibility of detecting ACD differences among those strata with airborne LiDAR.First, however, to improve our coverage of the study region, and to minimize error resulting from atmospheric contamination of satellite imagery, we supplemented the original satellite data used for pre-flight stratification with an additional 85 Landsat TM and ETM+ images to derive maps of fractional photosynthetic vegetation (PV), non-photosynthetic vegetation (NPV) and bare soil cover using CLASlite (Supplement).These results were combined with our original historical forest change maps (deforestation, disturbance, regrowth) and variables from SRTM including elevation, slope, aspect and a terrain ruggedness index (TRI).The TRI was computed as the square root of the sum of the squared differences between each pixel and its surrounding 11 × 11 pixel kernel (modified from Riley et al., 1999), and thus it is an index of topographic variability at local scales.With these input satellite and airborne LiDAR data, we sought to maximize parsimony and minimize ACD variation within each stratum.In areas of airborne LiDAR coverage, we assessed relationships between environmental variables from satellite, both individually and in combination, and forest ACD from LiDAR using correlation and regression analyses, respectively.Based on these findings, the final stratification approach followed an extension of the decision tree used for preliminary stratification (Fig. 2c) in 136 final classes (Supplement).We then calculated ACD statistics from the LiDAR coverage of each class.Following Asner et al. ( 2010), we determined that the median value best represented the ACD distribution within each class, and then populated all classes in the map with its corresponding median ACD value.

Regional upscaling based on regression
We also upscaled the LiDAR-based ACD estimates using an alternative regression-based approach (Saatchi et al., 2011).However, to take advantage of the historical forest cover change information (deforestation, degradation, regrowth) afforded by the CLASlite analyses described earlier, we applied the regression approach only to the forested classes, and then embedded the non-forest results derived earlier into the regression-based map.Thus, the regression technique utilized here was a hybrid of regression and stratification, although forested pixels accounted for the vast majority of the study area.The regression employed the same potential environmental parameters considered for the stratified approach.A combination of the parameters that explained significant variation in LiDAR-based ACD results was ultimately used for regional mapping.

Uncertainty analyses
We estimated uncertainties in each step of the project, including measurement errors in LiDAR (type 1) and field allometric equations underlying the universal LiDAR model (type 2), as well as prediction errors in application of the universal model (type 3) and upscaling to the regional level (type 4).Errors of types 1-3 were quantified empirically by previous studies, and we summarize each of these in the results.For errors of type 4 (upscaling), we subset our LiDAR coverage to 75 % of its full extent and evaluated the performance of each upscaling technique on the remaining 25 %.To subset the data, we restricted each flight polygon to a width of 2600 m, simulating a more limited flight campaign.With the 75 % subset, we completed each step of the analysis using the same methodology.We then compared predicted ACD to LiDAR-derived ACD at the pixel level in the remaining 25 % of the original LiDAR extent.Using this 75 % training/25 % validation split, we calculated average pixellevel uncertainty as the root mean square error between predicted and observed ACD.Because spatial resolution needs vary across data applications (e.g., research versus carbon accounting), we examined pixel level uncertainty at both 30-m and 1-ha resolution.
We also considered the possibility of spatial dependence in the errors produced by both upscaling approaches.We used Moran's I, an index of spatial autocorrelation, to determine whether the stratification and regression-based upscaling produced clustering of errors, which would indicate controls over ACD not captured by the models (Moran, 1950).Moran's I ranges from −1 to 1, where positive values indicate clustering and negative values indicate values arranged in opposition (e.g., a checkerboard type pattern).Due to computational limitations, we assessed Moran's I on randomly selected continuous subsamples of the upscaling residuals, accounting for 0.3 % and 2 % of the study area at 30-m and 1-ha resolution, respectively (R Development Core Team, 2011).

Validation of the universal LiDAR approach
Traditional field inventory produced a LiDAR-to-ACD regression that was in close agreement with the universal Li-DAR equation (compare Eqs. 2 and 3): ACD = 2.1043 MCH 1.363 .
(3) Specifically, the universal and traditional approaches were nearly identical in terms of both slope (r 2 values were the same) and predictive power (RMSE differed by only 0.3 Mg C ha −1 ) (Fig. 3).Thus, the universal model (Eq.2) was deemed sufficient for application to the 38 LiDAR sampling polygons used for mapping ACD.Validation of the universal LiDAR approach demonstrates for the first time that plot-scale calibration of airborne Li-DAR can be accomplished from a very small number of field plots for a new region.Furthermore, it demonstrates that LiDAR data can be calibrated using simpler field measurements of basal area (for determination of the stocking coefficient; Fig. S2a), and genus-level taxonomic sampling (for determination of regional WD; Fig. S2b).These very limited inputs required to calibrate the universal LiDAR model (as in Eq. 2) can lead to a major reduction in cost and effort needed to effectively characterize the LiDAR-to-ACD relationship for a region.

LiDAR-scale carbon patterns
Airborne LiDAR-based estimates of ACD revealed pronounced variation within and among sampling areas totaling 465 622 ha (Fig. 4).Landscapes impacted by human activity exhibited abrupt changes in ACD at forest edges and in degraded zones, particularly in the northwest region along the base and foothills of the Andes.In forested areas of the remote lowlands, ACD fluctuations aligned spatially with a range of natural environmental controls, mainly terrain features associated with mesas and local highlands, depressions in swamps and inundated areas, and riparian zones.Mapping and quantifying this local and sub-regional ACD variation throughout the region was requisite to the analysis of environmental controls over carbon density.

Environmental controls over ACD
Using environmental variables derived from the SRTM and CLASlite datasets against the LiDAR-based ACD maps, we found that elevation explained 19 % -the largest detected proportion -of the spatial variation in carbon stocks Biogeosciences, 9, 2683-2696, 2012  throughout the study region (Fig. 5).Correlations between elevation and biomass are not new (e.g., Aplet et al., 1998), but it was surprising here because the study area is relatively flat, with over 85 % of the forested area within a narrow band of elevation (100-300 m).PV fractional cover accounted for 9 % of the regional variation in ACD, and bare soil fraction just 3 %.Slope, aspect and the TRI each explained less than 2 % of the ACD variation at the whole-region scale.
Within the forested classes, we observed higher ACD with increasing PV, peaking with a PV fraction of 84-96 %, above which we observed a significant decline in ACD (Fig. 5).Further review revealed that the high-PV/low-ACD condition corresponded to swampland, agriculture, or other secondary vegetation.Through iterative evaluation of ACD distributions within PV sub-strata, we established thresholds for further segmentation by PV fraction into three discrete classes: 80-83 %, 84-96 %, and 97-100 % (Fig. 2c).Within each of these PV fraction classes, the distribution of carbon stocks was also found to change with elevation (Fig. 6).ACD was most sensitive to elevation in the PV class representing the lowest overall canopy cover (80-83 %), with a pronounced leftward skew of the ACD distributions as elevation increased from less than 100 m to more than 500 m.
The mid-PV class (84-96 %) accounted for more than 70 % of forested area across all elevations and a majority of aboveground carbon storage within the LiDAR coverage.Within this class, we observed that variation in the TRI was locally associated with ACD (Fig. S4), a pattern not apparent at the full regional scale (Fig. 5).To minimize ACD variation within each mapping class for subsequent regional extrapolation, we iteratively identified a TRI threshold to bisect each class at mid-PV values (Fig. 2c).Pixels with low TRI (0-5) generally corresponded to flatter riparian corridors, often with lower-ACD forests than their higher TRI (5+) counterparts (Fig. S4).We found no other sub-regional or local-scale patterns relating ACD to other potential factors (e.g., aspect, soil fraction).
We emphasize that other explanatory factors -such as soil type or distance from infrastructure -were considered for use in these models, but the available data were found to lack the detail needed to assign additional spatial variance in ACD to them.This is a common problem when scaling up field or airborne data in many tropical regions, and here we limited our analysis to variables readily derived from globally-available satellite imagery.Nonetheless, we also recognize that some regions do have much more detailed data available, which has been used for analyses of environmental controls over carbon stocks (Asner et al., 2011).

Regional carbon stocks -comparing approaches
Results of the environmental analysis guided the upscaling of LiDAR ACD estimates to the regional level, first with the stratification-based approach.Here, we iteratively derived segmentation thresholds for each variable to minimize ACD variation within each stratum (Maniatis and Mollicone, 2010).Fractional PV cover, elevation and the TRI explained the largest proportion of the ACD variation, with PV and elevation operating at the regional level and TRI at the local scale as described earlier.While these factors emerged as important controls over ACD, we also incorporated catchment boundaries to account for additional variation observed in the LiDAR data that was expressed irrespective of terrain or fractional cover (Fig. 2c) (see Supplement).
The final stratification-based ACD map contains classes, each with a median ACD value computed from distributions developed with airborne LiDAR sampling, and thus representing the highest probability carbon value in each class (see Table S3).The stratification-based map indicated lower ACD levels in the northwest where elevation is higher (Fig. 7a).However, within this area, heavy deforestation and degradation have also driven ACD as low as 0-6 Mg C ha −1 in many areas.Near the center of the study region, sandstone plateaus are capped by short-statured vegetation with low ACD.Carbon stocks increase to the southeast, where elevations are lowest and annual rainfall reaches 3000 mm yr −1 (Tian et al., 2000).Median ACD in these wetter forests reaches 130 Mg C ha −1 , and over many stretches, forest cover is uninterrupted, suggesting low to no human use.Naturally low ACD features also exist in the far south, such as in topographic depressions that are often inundated by water (Duivenvoorden and Duque, 2010).Finally, the map highlights the mediating role of catchment and localized terrain controls on carbon stocks, including the suppression of biomass in active floodplains and riparian corridors to the southeast.
The regression method required relatively little processing, and it yielded an regional ACD pattern similar to that derived with the stratification approach (Fig. 7b).We determined that only PV and elevation influenced the fit of the regression model at the scale of the entire study area (Sup-plement).It is important to note that the TRI did not add to the strength of the model.In this region, the TRI operates at a localized scale, and these controls are difficult to incorporate into regional regression models.In addition to being unable to capture localized terrain effects, the regional regression approach, by definition, does not incorporate subregional features such as the active floodplain areas containing early successional vegetation with low carbon stocks.
Total regional aboveground carbon storage was estimated at 1.468 Pg and 1.454 Pg using the stratification and regression approaches, respectively.Similarly, mean ACD among forested classes was 103.8 and 102.8 Mg C ha −1 .Notwithstanding an overall systematic bias, we would expect the two maps shown in Fig. 7 to produce similar regional-scale carbon stock values, especially here because each map incorporates the same non-forest results derived through CLASlite    analysis.In general, spatially explicit carbon mapping methods will yield similar regionally-integrated results as long as there are no large systematic biases, such as those observed when generic ACD values are applied to low resolution forest cover maps (Gibbs et al., 2007;Goetz et al., 2009;Avitabile et al., 2011).We have found this to be a major drawback of using Tier-I mapping guidelines provided by the IPCC (IPCC, 2006), which can have large systematic biases in comparison to high-resolution results (Asner et al., 2010(Asner et al., , 2011)).
Despite the regionally-integrated similarities between upscaling approaches, there were systematic differences -some exceeding 15 Mg C ha −1 -at the sub-regional level (Fig. 8).Stratification-based mapping resulted in lower ACD in two catchments and in the active floodplains (red in Fig. 8), again because stratification and LiDAR sampling suggested that sub-regional scale variation was important.Conversely, the stratification approach yielded higher ACD values throughout the majority of upland forests, owing to the availability of the TRI to resolve local topographic variability.We   8. Difference in aboveground carbon density calculated using regional stratification and regression analysis (from Fig. 7).believe this may be associated with either higher soil fertility or increased drainage (see also Asner et al., 2010).Another difference arises in the central portion of the study region, where small residual artifacts within and along the edges of the satellite imagery are expressed in the regression map.An advantage of the stratification approach is its ability to diminish such effects because mapped carbon values are associated with each stratum, not with each satellite image.Overall, an advantage of quantitatively comparing both upscaling methods is that the areas of maximum disagreement can be targeted for additional airborne LiDAR sampling.Finally, the regional comparisons provided here serve as a starting point for additional comparisons using other methods, including the global approaches (Saatchi et al., 2011;Baccini et al., 2012).

Sources of uncertainty in carbon mapping
We identified and propagated four sources of error throughout the study (Table 1).Measurement errors included LiDAR height errors and allometric errors in field-estimated ACD.LiDAR height errors are influenced by uncertainty in sensor position and orientation as well as laser ranging.These errors are low at the top of the canopy (∼ 0.15 m by Asner et al., 2007), but may reach 1 m at the ground in Amazon forests with relatively poor GPS networks.Thus we estimated Li-DAR height errors to be 5 %, assuming a mean forest canopy height of 20 m.
Allometric errors might influence the application of the universal LiDAR calibration equation to a new region (below we separately consider errors in its predictions at the pixel level).The universal model was calibrated using field-based ACD estimates determined by the application of allometric regression equations; for most large trees and for most re-gions, a pan-tropical model by Chave et al. (2005) was used.Thus, the extent to which this pan-tropical model is not representative of the Colombian Amazon may introduce a bias in our LiDAR-based estimates.Chave et al. (2005) considered errors in model application of the pan-tropical model at the regional scale (e.g., its accuracy for a "new" region); for all sites, their average reported error is 9.4 %.Note that each measurement error considered is a source of potential bias; these errors are directional and do not differ across spatial resolutions.
We identified two sources of prediction error, including the application of the universal LiDAR equation to airborne data, and the extrapolation of LiDAR-scale ACD estimates to the regional level.The universal LiDAR equation predicted the ACD of field validation plots with an uncertainty of 18.1 Mg C ha −1 (Fig. 3), or about 20 % at the 0.28 ha resolution of the field plots.While we relied on a limited number of plots, 20 % error at 0.28 ha spatial resolution is well within the range of values previously reported for tropical LiDAR studies.Furthermore, as we detail in the Supplement, LiDAR calibration errors can be modeled predictably with spatial resolution.Briefly, Mascaro et al. (2011b) used a large inventory plot with mapped trees to demonstrate that CAO Li-DAR calibration errors scaled according to the inverse square root of the plot area between 30 m and 1 ha spatial resolution (i.e., from ∼ 35 to 10 %, respectively).Mascaro et al. further showed that at resolutions below 1 ha, nearly half of the error is caused by disagreement between LiDAR and field measurements concerning which trees (or portions of trees) are considered inside or outside of a calibration plot.This portion of the error is essentially "sub-tree" -existing at a spatial resolution smaller than tree crowns -and thus can be safely omitted for carbon mapping purposes.Following this empirical analysis, we estimate calibration errors of 21.5 % and 10.0 % at 30 m and 1 ha resolution, respectively (Table 1).
Finally, we evaluated pixel-level uncertainty associated with extending LiDAR results to the regional level.To do so, we produced regional ACD maps based on a 75 % of the Li-DAR coverage (following the same methodology), and evaluated their performance in the remaining 25 % extent.Pixellevel errors were 29.0 Mg C ha −1 (32.9 % of the mean carbon stock) using the stratification approach, and 32.6 Mg C ha −1 (37.0 %) using the regression approach.At 1-ha spatial resolution, these errors declined to 24.3 % and 28.7 % for stratification and regression approaches, respectively.
We assumed that all four sources of error were independent, and therefore propagated the errors as the square root of the sum of all squared errors (Table 1).Overall, we estimated a 28.3 % pixel-level error for ACD determinations across the study corridor at 1 ha resolution using the stratification method, and 32.2 % using the regression approach.A side-by-side comparison of LiDAR-derived ACD and regional ACD derived from land-cover stratification indicates good agreement between LiDAR-observed and regionally-predicted ACD using the stratification approach (Fig. S5).Finally, we contrasted mean LiDAR-observed and regionally-predicted ACD in the large sections of each flight polygon set aside for determination of pixel-level errors.We found that both approaches predicted ACD well, although a limited bias was observed using the regression approach (Fig. S6); this is consistent with its slightly weaker performance overall.Collectively, we contend that the environmental controls identified through iterative statistical analysis and incorporated into regional stratification effectively represent the ACD variation found throughout the region.Residuals for the stratification-based approach indicated very low spatial autocorrelation of errors: Moran's I ranged from 0.05 to 0.08 depending on spatial resolution (Table 2).By contrast, regression-based upscaling yielded higher spatial autocorrelation; this clustering of errors suggests that modest variation in ACD was not explained by the regression model.These results, along with higher pixel-level errors for the regression approach, suggest that stratification better captured landscape-scale controls on ACD due to a combination of (Eq. 1) ingesting additional possible controls over ACD (i.e., the inclusion of localized TRI variation in the stratification approach which was not retained in the regression analysis), and (Eq.2) the degree to which unknown localized controls on ACD are captured by stratification but remain absent from the regression inputs.

Conclusions
We demonstrate that high-resolution mapping of tropical forest carbon stocks assisted by airborne LiDAR can be accomplished with limited field calibration data and limited preexisting knowledge of the study region.In place of previous calibration models that relied on exhaustive inventories of tree diameters, heights, and wood densities, we used a new universal LiDAR model -adjusted by rapid field-based assessment of basal area and a regional wood density constant -Table 2. Moran's I values for residuals produced by stratificationand regression-based upscaling approaches.At 30-m resolution, the analysis was conducted on a subset of points representing 0.3 % of the total data set (n = 15 664), and 2 % at 1-ha resolution (n = 8611).All values are significant at P < 0.001.

Stratification Regression
30 m 1 ha 30 m 1 ha Moran's I 0.047 0.083 0.125 0.186 to produce a calibration for which we have equal confidence with the traditional approach.The universal approach provides cost-effective calibration and use of LiDAR in remote and difficult-to-access areas such as the Colombian Amazon.
We also found that, even with limited foreknowledge of carbon stock variation in a region, a systematic analysis of environmental controls on LiDAR-scale ACD can provide an effective means to upscale LiDAR measurements throughout a region.Through this process, we gained new insight into ecological drivers of ACD variation across the Colombian Amazon, revealing previously unknown variation mediated by elevation, terrain ruggedness, and canopy fractional cover.We found strengths and weaknesses in both stratification and regression-based approaches to upscaling.Stratificationbased mapping provides a means to dissect a region by potential environmental controls, sample the resulting strata with airborne LiDAR, and apply the results to minimize ACD variance per mapping class.On the other hand, regression approaches are less laborious and capture the general trends, but will miss the landscape features that are not expressed consistently throughout a region.Clearly, applying both approaches provides analytical leverage to identify areas of uncertainty for additional investigation.High-resolution approaches that report pixel-scale uncertainties will provide the most confidence in the effort to monitor changes in tropical www.biogeosciences.net/9/2683/2012/Biogeosciences, 9, 2683-2696, 2012 forest carbon stocks.This improved confidence will allow resource managers and decision-makers to more rapidly and effectively implement actions that better utilize and conserve forests in remote tropical regions.

Figure 1 .
Figure 1.A map showing the fractional cover of photosynthetic vegetation (PV; green), non-710

Fig. 2 .
Fig. 2. Decision trees used to stratify the study region with multi-temporal Landsat imagery analyzed with CLASlite, and a digital elevation model (DEM) derived from the NASA Shuttle Radar Topography Mission (SRTM).Preliminary stratification (i.e., prior to LiDAR flights) utilized (a) deforestation, degradation, and regrowth partitioning provided by CLASlite as well as (b) discrete SRTM DEM classes to further partition intact forest.The final stratification further partitioned intact forest by (c) fractional photosynthetic vegetation (PV) cover, terrain ruggedness index (TRI), and catchments.Components encircled in blue were later considered as inputs to the regression technique.

Fig. 3 .
Fig. 3. Regressions comparing the universal LiDAR equation of Asner et al. (2012b) with traditional validation field plots measured within the study region.

Figure 4 .
Figure 4. Example LiDAR-based maps of aboveground carbon density (ACD; left side) and 3

Fig. 4 .
Fig. 4. Example LiDAR-based maps of aboveground carbon density (ACD; left side) and ground elevation (right side) for (a) heavily deforested lands in the northwest region; (b) alluvial and riparian forests in the central region; and (c) terra firme forests in the southeast region of the study.Units are Mg C ha −1 for ACD and meters for elevation, and both are shown at 30 m spatial resolution.ACD values are shown in actual units; elevation data are stretched for maximum contrast, but generally varies by 200 m or less in each image.

Figure 6 .
Figure 6.Example frequency distributions of LiDAR-derived aboveground carbon density by

Fig. 7 .
Fig. 7. Aboveground carbon density (ACD) across the study region calculated using (a) regional stratification with elevation, fractional cover of photosynthetic vegetation and terrain ruggedness; and (b) regression analysis with elevation and fractional cover of photosynthetic vegetation.

Figure 8 .
Figure 8. Difference in aboveground carbon density calculated using regional stratification and 65 Fig.8.Difference in aboveground carbon density calculated using regional stratification and regression analysis (from Fig.7).

Table 1 .
Sources of error in high-resolution mapping of aboveground carbon density in the Colombian Amazon.Errors are provided at resolutions of 30 m and 1 ha using the stratification-and regression-based upscaling approaches.