Articles | Volume 15, issue 12
Research article
22 Jun 2018
Research article |  | 22 Jun 2018

Estimating aboveground carbon density and its uncertainty in Borneo's structurally complex tropical forests using airborne laser scanning

Tommaso Jucker, Gregory P. Asner, Michele Dalponte, Philip G. Brodrick, Christopher D. Philipson, Nicholas R. Vaughn, Yit Arn Teh, Craig Brelsford, David F. R. P. Burslem, Nicolas J. Deere, Robert M. Ewers, Jakub Kvasnica, Simon L. Lewis, Yadvinder Malhi, Sol Milne, Reuben Nilus, Marion Pfeifer, Oliver L. Phillips, Lan Qie, Nathan Renneboog, Glen Reynolds, Terhi Riutta, Matthew J. Struebig, Martin Svátek, Edgar C. Turner, and David A. Coomes

Borneo contains some of the world's most biodiverse and carbon-dense tropical forest, but this 750 000 km2 island has lost 62 % of its old-growth forests within the last 40 years. Efforts to protect and restore the remaining forests of Borneo hinge on recognizing the ecosystem services they provide, including their ability to store and sequester carbon. Airborne laser scanning (ALS) is a remote sensing technology that allows forest structural properties to be captured in great detail across vast geographic areas. In recent years ALS has been integrated into statewide assessments of forest carbon in Neotropical and African regions, but not yet in Asia. For this to happen new regional models need to be developed for estimating carbon stocks from ALS in tropical Asia, as the forests of this region are structurally and compositionally distinct from those found elsewhere in the tropics. By combining ALS imagery with data from 173 permanent forest plots spanning the lowland rainforests of Sabah on the island of Borneo, we develop a simple yet general model for estimating forest carbon stocks using ALS-derived canopy height and canopy cover as input metrics. An advanced feature of this new model is the propagation of uncertainty in both ALS- and ground-based data, allowing uncertainty in hectare-scale estimates of carbon stocks to be quantified robustly. We show that the model effectively captures variation in aboveground carbon stocks across extreme disturbance gradients spanning tall dipterocarp forests and heavily logged regions and clearly outperforms existing ALS-based models calibrated for the tropics, as well as currently available satellite-derived products. Our model provides a simple, generalized and effective approach for mapping forest carbon stocks in Borneo and underpins ongoing efforts to safeguard and facilitate the restoration of its unique tropical forests.

1 Introduction

Forests are an important part of the global carbon cycle (Pan et al., 2011), storing and sequestering more carbon than any other ecosystem (Gibbs et al., 2007). Estimates of tropical deforestation rates vary, but roughly 61 300 km2 of forest were lost each year between 2000 and 2012, and an additional 30 % were degraded by logging or fire (Asner et al., 2009; Hansen et al., 2013). Forest degradation and deforestation result in about 1–2 billion tonnes of carbon being released to the atmosphere each year, which equates to about 10 % of global emissions (Baccini et al., 2012). Even if nations decarbonize their energy supply chains within agreed schedules, a rise of 2 C in mean annual temperature is unavoidable unless 300 million hectares of degraded tropical forests are protected and land unsuitable for agriculture is reforested (Houghton et al., 2015). Signatories to the Paris Agreement, brokered at COP21 in 2015, are now committed to reducing emissions from tropical deforestation and forest degradation (i.e. REDD+; Agrawal et al., 2011), whilst recognizing that these forests also harbour rich biodiversity and support livelihoods for around a billion people (Vira et al., 2015).

The accurate monitoring of forest carbon stocks underpins these initiatives to generate carbon credits through REDD+ and similar forest conservation and climate change mitigation programmes (Agrawal et al., 2011). Airborne laser scanning (ALS) has shown particular promise in this regard because it generates high-resolution maps of forest structure from which aboveground carbon density (ACD) can be estimated (Asner et al., 2010; Lefsky et al., 1999; Nelson et al., 1988; Popescu et al., 2011; Wulder et al., 2012). The principle of ALS is that laser pulses are emitted downwards from an aircraft, and a sensor records the time it takes for individual beams to strike a surface (e.g. leaves, branches or the ground) and bounce back to the emitting source, thereby precisely measuring the distance between the object and the airborne platform. Divergence of the beam means it is wider than leaves and allows for penetration into the canopy, resulting in a 3-D point cloud that captures the vertical structure of the forest. By far the most common approach to using ALS data for estimating forest carbon stocks involves developing statistical models relating ACD estimates obtained from permanent field plots to summary statistics derived from the ALS point cloud, such as the mean height of returns or their skew (Zolkos et al., 2013). These “area-based” approaches were first used for mapping the structural attributes of complex multi-layered forests in the early 2000s (Drake et al., 2002; Lefsky et al., 2002) and have since been applied to carbon mapping in several tropical regions (Asner et al., 2010, 2014; Jubanski et al., 2013; Laurin et al., 2014; Réjou-Méchain et al., 2015).

Figure 1Panel (a) shows the location of the Sepilok and Kuamut Forest Reserve, the Danum Valley and Maliau Basin Conservation Area, and the SAFE landscape within Sabah (Malaysia). Green shading in the background represents forest cover at 30 m resolution in the year 2000 (Hansen et al., 2013). In panel (b), the relationship between field-measured aboveground carbon density and ALS-derived top-of-canopy height found across the study sites (coloured circles, n=173) is compared to measurements taken mostly in the Neotropics (Asner and Mascaro, 2014; grey circles, n=754).


This paper develops a statistical model for mapping forest carbon and its uncertainty in South-east Asian forests. We work with ALS and plot data collected in the Malaysian state of Sabah on the north-eastern end of the island of Borneo (Fig. 1), which is an important test bed for international efforts to protect and restore tropical forests. Borneo lost around 62 % of its old-growth forest in just 40 years as a result of heavy logging and the subsequent establishment of oil palm and forestry plantations (Gaveau et al., 2014, 2016). Sabah lost its forests at an even faster rate in this period (Osman et al., 2012), and because these forests are amongst the most carbon dense in the tropics, carbon loss has been considerable (Carlson et al., 2012a, b; Slik et al., 2010). In response to past and ongoing forest losses, the Sabah state government has recently taken a number of concrete steps towards becoming a regional leader in forest conservation and sustainable management. Among these was commissioning a new high-resolution wall-to-wall carbon map for the entire state, which will inform future forest conservation and restoration efforts across the region. Here we develop the ALS-based model that underpins this new carbon map (Asner et al., 2018).

The approach we take builds on the work of Asner and Mascaro (2014), who proposed a general model for estimating ACD (in Mg C ha−1) in tropical forests using a single ALS metric – the mean top-of-canopy height (TCH, in m) – and minimal field data inputs. The method relates ACD to TCH, stand basal area (BA; in m2 ha−1) and the community-weighted mean wood density (WD; in g cm−3) over a prescribed area of forest, such as 1 ha, as follows:

(1) ACD General = 3.836 × TCH 0.281 × BA 0.972 × WD 1.376 .

Asner and Mascaro (2014) demonstrated that tropical forests from 14 regions differ greatly in structure. Remarkably, they found that a generalized power-law relationship could be fitted that transcended these contrasting forests types, once regional differences in structure were incorporated as sub-models relating BA and WD to TCH. However, this general model may generate systematic errors in ACD estimates if applied to regions outside the calibration range, and Asner and Mascaro (2014) make clear that regional models should be obtained where possible. Since South-east Asian rainforests were not among the 14 regions used to calibrate the general model and are phylogenetically and structurally distinct from Neotropical and Afrotropical forests (Banin et al., 2012), new regional models are needed before Borneo's forest carbon stocks can be surveyed using ALS. Central to the robust estimation of ACD using ALS data is identifying a metric which captures variation in basal area among stands. Asner and Mascaro's (2014) power-law model rests on an assumption that basal area is closely related to top-of-canopy height, an assumption supported in some studies, but not in others (Coomes et al., 2017; Duncanson et al., 2015; Spriggs, 2015). The dominance of Asian lowland rainforests by dipterocarp species makes them structurally unique (Banin et al., 2012; Feldpausch et al., 2011; Ghazoul, 2016) and gives rise to greater aboveground carbon densities than anywhere else in the tropics (Avitabile et al., 2016; Sullivan et al., 2017), highlighting the need for new ALS-based carbon estimation models for this region.

Here we develop a regional model for estimating ACD from ALS data that underpins ongoing efforts to map Sabah's forest carbon stocks at high resolution to inform conservation and management decisions for one of the world's most threatened biodiversity hotspots (Asner et al., 2018; Nunes et al., 2017). Building on the work of Asner and Mascaro (2014), we combine ALS data with estimates of ACD from a total of 173 permanent forests plots spanning the major lowland dipterocarp forest types and disturbance gradients found in Borneo to derive a simple yet general equation for predicting carbon stocks from ALS metrics at hectare resolution. As part of this approach we also develop a novel framework for propagating uncertainty in both ALS- and ground-based data, allowing uncertainty in hectare-scale estimates of carbon stocks to be quantified robustly. To assess the accuracy of this new model, we then benchmark it against existing ALS-derived equations of ACD developed for the tropics (Asner and Mascaro, 2014) and satellite-based carbon maps of the region (Avitabile et al., 2016; Pfeifer et al., 2016).

Table 1Summary of permanent forest plot data collected at each study site and description of which ALS sensor was used at each location. Plot size is in hectares, while minimum stem diameter thresholds (Dmin) are given in centimetres.

a Mean plot size after applying slope correction (see Sect. 2.2.2 for further details).
b Plots established as part of the SAFE experiment, and those located along riparian buffer zones in the SAFE landscape were aggregated into spatial blocks prior to statistical analyses (n=27 with a mean plot size of 0.5 ha; see Sect. 2.2.4 for further details).

Download Print Version | Download XLSX

2 Materials and methods

2.1 Study region

The study was conducted in Sabah, a Malaysian state in northern Borneo (Fig. 1a). Mean daily temperature is 26.7 C and annual rainfall is 2600–3000 mm (Walsh and Newbery, 1999). Severe droughts linked to El Niño events occur about once every 10 years (Malhi and Wright, 2004; Walsh and Newbery, 1999). Sabah supports a wide range of forest types, including dipterocarp forests in the lowlands that are among the tallest in the tropics (Fig. 1b; Banin et al., 2012).

2.2 Permanent forest plot data

We compiled permanent forest plot data from five forested landscapes across Sabah (Fig. 1a): Sepilok Forest Reserve, Kuamut Forest Reserve, Danum Valley Conservation Area, the Stability of Altered Forest Ecosystems (SAFE) experimental forest fragmentation landscape (Ewers et al., 2011) and Maliau Basin Conservation Area. Here we provide a brief description of the permanent plot data collected at each site, which are summarized in Table 1. Additional details are provided in Supplement S1.

2.2.1 Sepilok Forest Reserve

The reserve is a protected area encompassing a remnant of coastal lowland old-growth tropical rainforest (Fox, 1973) and is characterized by three strongly contrasting soil types that give rise to forests that are structurally and functionally very different (Dent et al., 2006; DeWalt et al., 2006; Nilus et al., 2011): alluvial dipterocarp forest in the valleys (hereafter alluvial forests), sandstone hill dipterocarp forest on dissected hillsides and crests (hereafter sandstone forests), and heath forest on podzols associated with the dip slopes of cuestas (hereafter heath forests). We used data from nine permanent 4 ha forest plots situated in the reserve, three in each forest type. These were first established in 2000–2001 and were most recently re-censused in 2013–2015. All stems with a diameter (D, in cm)  10 cm were recorded and identified to species (or closest taxonomic unit). Tree height (H, in m) was measured for a subset of trees (n=718) using a laser range finder. For the purposes of this analysis, each 4 ha plot was subdivided into 1 ha subplots, giving a total of 36 plots 1 ha in size. The corners of the plots were geolocated using a Geneq SXBlue II global positioning system (GPS) unit, which uses satellite-based augmentation to perform differential correction and is capable of a positional accuracy of less than 2 m (95 % confidence intervals).

2.2.2 Kuamut Forest Reserve

The reserve is a former logging area that is now being developed as a restoration project. Selective logging during the past 30 years has left large tracts of forest in a generally degraded condition, although the extent of this disturbance varies across the landscape. Floristically and topographically the Kuamut reserve is broadly similar to Danum Valley – with which it shares a western border – and predominantly consists of lowland dipterocarp forests. Within the forest reserve, 39 circular plots with a radius of 30 m were established in 2015–2016 spanning a range of forest successional stages, including young secondary forests characterized by the presence of species with low wood densities (e.g. Macaranga spp.). Coordinates for the plot centres were taken using a Garmin GPSMAP 64s device with an accuracy of ±10 m (95 % confidence intervals). Within each plot, all stems with D≥10cm were recorded and identified to species (or the closest taxonomic unit), and H was measured using a laser range finder. Because the radius of the plots was measured along the slope of the terrain (as opposed to a horizontally projected distance), we slope-corrected the area of each plot by multiplying by cos (θ), where θ is the average slope of the plot in degrees as calculated from the digital elevation model obtained from the ALS data. The average plot size after applying this correction factor was 0.265 ha (6 % less than if no slope correction had been applied).

2.2.3 Danum Valley Conservation Area

The site encompasses the largest remaining tract of primary lowland dipterocarp forest in Sabah. Within the protected area, we obtained data from a 50 ha permanent forest plot which was established in 2010 as part of the Centre for Tropical Forest Science (CTFS) ForestGEO network (Anderson-Teixeira et al., 2015). Here we focus on 45 ha of this plot for which all stems with D≥1cm have been mapped and taxonomically identified (mapping of the remaining 5 ha of forest was ongoing as of January 2017). For the purposes of this study, we subdivided the mapped area into 45 1 ha plots, the coordinates of which were recorded using the Geneq SXBlue II GPS. In addition to the 50 ha CTSF plot, we also secured data from 20 circular plots with a 30 m radius that were established across the protected area by the Carnegie Airborne Observatory (CAO) in 2017. These plots were surveyed following the same protocols as those described previously for the plots at Kuamut in Sect. 2.2.2.

2.2.4 SAFE landscape and Maliau Basin Conservation Area

Plot data from three sources were acquired from the SAFE landscape and the Maliau Basin Conservation Area: research plots established through the SAFE project, plots used to monitor riparian buffer zones and plots from the Global Ecosystem Monitoring (GEM) network (, last access: 20 June 2018). As part of the SAFE project, 166 plots 25×25m in size were established in forested areas (Ewers et al., 2011; Pfeifer et al., 2016). Plots are organized in blocks which span a land-use intensity gradient ranging from twice-logged forests that are currently in the early stages of secondary succession within the SAFE landscape to relatively undisturbed old-growth forests at Maliau Basin (Ewers et al., 2011; Struebig et al., 2013). Plots were surveyed in 2010, at which time all stems with D≥10cm were recorded and plot coordinates were taken using a Garmin GPSMAP60 device (accurate to within ±10 m; 95 % confidence intervals). Of these 166 plots, 38 were re-surveyed in 2014, at which time all stems with D≥1cm were recorded and tree heights were measured using a laser range finder. Using these same protocols, a further 48 plots were established in 2014 along riparian buffer zones in the SAFE landscape. As with the SAFE project plots, riparian plots are also spatially clustered into blocks. The small size of the SAFE and riparian plots (0.0625 ha) makes them prone to high uncertainty when modelling carbon stocks from ALS (Réjou-Méchain et al., 2014), especially given the relatively low positional accuracy of the GPS coordinates. To minimize this source of error, we chose to aggregate individual plots into blocks for all subsequent analyses (n=27, with a mean size of 0.5 ha). Lastly, we obtained data from six GEM plots – four within the SAFE landscape and two at Maliau Basin. The GEM plots are 1 ha in size and were established in 2014. All stems with D≥10cm were mapped, measured for height using a laser range finder and taxonomically identified. The corners of the plots were georeferenced using the Geneq SXBlue II GPS.

2.3 Estimating aboveground carbon density and its uncertainty

Across the five study sites we compiled a total of 173 plots that together cover a cumulative area of 116.1 ha of forest. For each of these plots we calculated aboveground carbon density (ACD, in Mg C ha−1) following the approach outlined in the BIOMASS package in R (R Core Development Team, 2016; Réjou-Méchain et al., 2017). This provides a workflow to not only quantify ACD, but also propagate uncertainty in ACD estimates arising from both field measurement errors and uncertainty in allometric models. The first step is to estimate the aboveground biomass (AGB, in kg) of individual trees using Chave et al.'s (2014) pantropical biomass equation: AGB=0.067×D2×H×WD0.976. For trees with no height measurement in the field, H was estimated using a locally calibrated HD allometric equation, while wood density (WD, in g cm−3) values were obtained from the global wood density database (Chave et al., 2009; Zanne et al., 2009; see Supplement S1 for additional details on both H and WD estimation).

In addition to quantifying AGB, Réjou-Méchain et al.'s (2017) workflow uses Monte Carlo simulations to propagate uncertainty in biomass estimates due to (i) measurement errors in D (following Chave et al.'s (2004) approach, in which 95 % of stems are assumed to contain small measurement errors that are in proportion to D, while the remaining 5 % is assigned a gross measurement error of 4.6 cm), (ii) uncertainty in HD allometries, (iii) uncertainty in WD estimates arising from incomplete taxonomic identification and/or coverage of the global wood density database, and (iv) uncertainty in the AGB equation itself. Using this approach, we generated 100 estimates of AGB for each recorded tree. ACD was then quantified by summing the AGB of all trees within a plot, dividing the total by the area of the plot and applying a carbon content conversion factor of 0.47 (Martin and Thomas, 2011). By repeating this across all simulated values of AGB, we obtained 100 estimates of ACD for each of the 173 plots that reflect the uncertainty in stand-level carbon stocks (note that a preliminary analysis showed that 100 iterations were sufficient to robustly capture mean and standard deviation values of plot-level ACD, while also allowing for efficient computing times). As a last step, we used data from 45 plots in Danum Valley – where all stems with D≥1cm were measured – to develop a correction factor that compensates for the carbon stocks of stems with D<10cm that were not recorded (Phillips et al., 1998; see Eq. S2 in Supplement S1).

Stand basal area and wood density estimation

In addition to estimating ACD for each plot, we also calculated basal area (BA, in m2 ha−1) and the community-weighted mean WD, as well as their uncertainties. BA was quantified by summing π×D/22 across all stems within a plot and then applying a correction factor that accounts for stems with D<10cm that were not measured (see Eq. S3 in Supplement S1). In the case of BA, uncertainty arises from measurement errors in D, which were propagated through following the approach of Chave et al. (2004) described in Sect. 2.3. The community-weighted mean WD of each plot was quantified as BAij×WDi, where BAij is the relative basal area of species i in plot j, and WDi is the mean wood density of species i. Uncertainty in plot-level WD reflects incomplete taxonomic information and/or lack of coverage in the global wood density database.

2.4 Airborne laser scanning data

ALS data covering the permanent forest plots described in Sect. 2.2 were acquired through two independent surveys, the first undertaken by NERC's Airborne Research Facility (ARF) in November of 2014 and the second by the Carnegie Airborne Observatory (CAO) in April of 2016. Table 1 specifies which plots were overflown with which system. NERC ARF operated a Leica ALS50-II lidar sensor flown on a Dornier 228-201 at an elevation of 1400–2400 m a.s.l. (depending on the study site) and a flight speed of 120–140 knots. The sensor emits pulses at a frequency of 120 kHz, has a field of view of 12 and a footprint of about 40 cm. The average point density was 7.3 points m−2. The Leica ALS50-II lidar sensor records both discrete point and full-waveform ALS, but for the purposes of this study only the discrete-return data, with up to four returns recorded per pulse, were used. Accurate georeferencing of the ALS point cloud was ensured by incorporating data from a Leica base station running in the study area concurrently to the flight. The ALS data were preprocessed by NERC's Data Analysis Node and delivered in LAS format. All further processing was undertaken using LAStools software (, last access: 20 June 2018). The CAO campaign was conducted using the CAO-3 system, a detailed description of which can be found in Asner et al. (2012). Briefly, CAO-3 is a custom-designed, dual-laser, full-waveform system that was operated in discrete-return collection mode for this project. The aircraft was flown at 3600 m a.s.l. at a flight speed of 120–140 knots. The ALS system was set to a field of view of 34 (after 2 cut-off from each edge) and a combined-channel pulse frequency of 200 kHz. The ALS pulse footprint at 3600 m a.s.l. was approximately 1.8 m. With adjacent flight-line overlap, these settings yielded approximately 2.0 points m−2. Despite differences in the acquisition parameters of the two surveys which can influence canopy metrics derived from ALS data (Gobakken and Næsset, 2008; Roussel et al., 2017), a comparison of regions of overlap between the flight campaigns showed strong agreement between data obtained from the two sensors (Supplement S2).

2.4.1 Airborne laser scanning metrics

ALS point clouds derived from both surveys were classified into ground and non-ground points, and a digital elevation model (DEM) was fitted to the ground returns to produce a raster at 1 m resolution. The DEM was then subtracted from the elevations of all non-ground returns to produce a normalized point cloud, from which a canopy height model (CHM) was constructed by averaging the first returns. Finally, any gaps in the raster of the CHM were filled by averaging neighbouring cells. From the CHMs we calculated two metrics for each of the permanent field plots: top-of-canopy height (TCH, in m) and canopy cover at 20 m aboveground (Cover20). TCH is the mean height of the pixels which make up the surface of the CHM. Canopy cover is defined as the proportion of area occupied by crowns at a given height aboveground (i.e. 1 – gap fraction). Cover20 was calculated by creating a plane horizontal to the ground in the CHM at a height of 20 m aboveground, counting the number of pixels for which the CHM lies above the plane and then dividing this number by the total number of pixels in the plot. A height of 20 m aboveground was chosen as previous work showed this to be the optimal height for estimating plot-level BA in an old-growth lowland dipterocarp forest in Sabah (Coomes et al., 2017).

2.4.2 Accounting for geopositional uncertainty

Plot coordinates obtained using a GPS are inevitably associated with a certain degree of error, particularly when working under dense forest canopies. However, this source of uncertainty is generally overlooked when attempting to relate field estimates of ACD to ALS metrics. To account for geopositional uncertainty, we introduced normally distributed random errors in the plot coordinates. These errors were assumed to be proportional to the operational accuracy of the GPS unit used to geolocate a given plot: ±2 m for plots recorded with the Geneq SXBlue II GPS and ±10 m for those geolocated using either the Garmin GPSMAP60 or Garmin GPSMAP 64s devices. This process was iterated 100 times, and at each step we calculated TCH and Cover20 across all plots. Note that for plots from the SAFE project and those situated along riparian buffer zones, ALS metrics were calculated for each individual 0.0625 ha plot before being aggregated into blocks (as was done for the field data).

2.5 Modelling aboveground carbon density and associated uncertainty

We started by using data from the 173 field plots to fit a regional form of Asner and Mascaro's (2014) model, in which ACD is expressed as the following function of ALS-derived TCH and field-based estimates of BA and WD:

(2) ACD = ρ 0 × TCH ρ 1 × BA ρ 2 × WD ρ 3 ,

where ρ0−3 represent constants to be estimated from empirical data. In order to apply Eq. (2) to areas where field data are not available, the next step is to develop sub-models to estimate BA and WD from ALS metrics. Of particular importance in this regard is the accurate and unbiased estimation of BA, which correlates very strongly with ACD (Pearson's correlation coefficient (ρ)=0.93 across the 173 plots). Asner and Mascaro (2014) found that a single ALS metric – TCH – could be used to reliably estimate both BA and WD across a range of tropical forest regions. However, recent work suggests this may not always be the case (Duncanson et al., 2015; Spriggs, 2015). In particular, Coomes et al. (2017) showed that ALS metrics that capture information about canopy cover at a given height aboveground – such as Cover20 – were better suited to estimating BA. Here we compared these two approaches to test whether Cover20 can prove a useful metric to distinguish between forests with similar TCH but substantially different BA.

2.5.1 Basal area sub-models

Asner and Mascaro (2014) modelled BA as the following function of TCH:

(3) BA = ρ 0 × TCH .

We compared the goodness-of-fit of Eq. (3) to a model that additionally incorporates Cover20 as a predictor of BA. Doing so, however, requires accounting for the fact that TCH and Cover20 are correlated. To avoid issues of collinearity (Dormann et al., 2013), we therefore first modelled the relationship between Cover20 and TCH using logistic regression and used the residuals of this model to identify plots that have higher or lower than expected Cover20 for a given TCH.

(4) ln Cover 20 1 - Cover 20 = ρ 0 + ρ 1 × ln TCH

Predicted values of canopy cover (Cover^20) can be obtained from Eq. (4) as follows:

(5) Cover ^ 20 = 1 1 + e - ρ 0 × TCH - ρ 1 .

From this, we calculated the residual cover (Coverresid) for each of the 173 field plots as Cover20-Cover^20 and then modelled BA as the following non-linear function of TCH and Coverresid:

(6) BA = ρ 0 × TCH ρ 1 × 1 + ρ 2 × Cover resid .

Equation (6) was chosen after careful comparison with alternative functional forms. This included modelling BA directly as a function of Cover20, without including TCH in the regression. We discarded this last option as BA estimates were found to be highly sensitive to small variations in canopy cover when Cover20 approaches 1.

2.5.2 Wood density sub-models

Following Asner and Mascaro (2014), we modelled WD as a power-law function of TCH:

(7) WD = ρ 0 × TCH ρ 1 .

The expectation is that, because the proportion of densely wooded species tends to increase during forest succession (Slik et al., 2008), taller forests should on average have higher stand-level WD values. While this explicitly ignores the well-known fact that WD is also influenced by environmental factors that have nothing to do with disturbance (e.g. soils or climate; Quesada et al., 2012), we chose to fit a single function for all sites, as from an operational standpoint applying forest-type-specific equations would require information on the spatial distribution of these forest types across the landscape (something which may not necessarily be available, particularly for the tropics). For comparison, we also tested whether replacing TCH with Cover20 would improve the fit of the WD model.

2.5.3 Error propagation and model validation

Just as deriving accurate estimates of ACD is critical to producing robust and useful maps of forest carbon stocks, so too is the ability to place a degree of confidence on the mean predicted values obtained from a given model (Réjou-Méchain et al., 2017). In order to fully propagate uncertainty in ALS-derived estimates of ACD, as well as robustly assessing model performance, we developed the following approach based on leave-one-out cross validation: (i) of the 173 field plots, 1 was set aside for validation, while the rest were used to calibrate models; (ii) the calibration dataset was used to fit both the regional ACD model (Eq. 2) and the BA and WD sub-models (Eqs. 3, 6–7); and (iii) the fitted models were used to generate predictions of BA, WD and ACD for the validation plot previously set aside. In each case, Monte Carlo simulations were used to incorporate model uncertainty in the predicted values. For Eqs. (4) and (6), parameter estimates were obtained using the L-BFGS-B non-linear optimization routine implemented in Python (Morales and Nocedal, 2011). For power-law models fit to log–log-transformed data (i.e. Eqs. 2 and 7), we applied the Baskerville (1972) correction factor by multiplying predicted values by exp (σ2∕2), where σ is the estimated standard deviation of the residuals (also known as the residual standard error); (iv) model fitting and prediction steps (ii)–(iii) were repeated 100 times across all estimates of ACD, BA, WD, TCH and Cover20 that had previously been generated for each field plot. This allowed us to fully propagate uncertainty in ACD arising from field measurement errors, allometric models and geopositional errors; (v) lastly, steps (i)–(iv) were repeated for all 173 field plots.

Once predictions of ACD had been generated for all 173 plots, we assessed model performance by comparing predicted and observed ACD values (ACDpred and ACDobs, respectively) on the basis of root mean square error (RMSE) calculated as 1Ni=1NACDobs-ACDpred2 and relative systematic error (or bias), which we calculated as 1Ni=1NACDpred-ACDobsACDobs×100 (Chave et al., 2014). Additionally, we tested how plot-level errors (calculated for each individual plot as ACDobs-ACDpredACDobs×100) varied as a function of forest carbon stocks and in relation to plot size (Réjou-Méchain et al., 2014).

The modelling and error propagation framework described above was chosen after a thorough comparison with a number of alternative approaches. The objective of this comparison was to identify the approach that would yield the lowest degree of systematic bias in the predicted values of ACD, as we consider this to be a critical requirement of any carbon estimation model, particularly if – as is the case here – that model is to underpin the generation of a carbon map designed to inform management and conservation policies (Asner et al., 2018). Of the two alternative approaches we tested, the first relied on fitting a combination of ordinary and non-linear least squares regression models to parameterize the equations presented above. As with the modelling routine described above, this approach did not account for potential spatial autocorrelation in the residuals of the models, which could result in a slight underestimation of the true uncertainty in the fitted parameter values. We contrasted this approach with one that used generalized and non-linear least squares regression that explicitly account for spatial dependencies in the data. Both these approaches underperformed compared to the routine described above, as they substantially overestimated ACD values in low-carbon-density forests and underestimated ACD in carbon-rich ones (see Supplement S3 for details). This tendency to introduce a systematic bias in the ACD predictions was particularly evident in the case of the spatially explicit models (see Fig. S4b in the Supplement). In light of this we opted for the approach presented here, even though we acknowledge that it may slightly underestimate uncertainty in modelled ACD values due to spatial non-independence in the data.

2.6 Comparison with satellite-derived estimates of aboveground carbon density

We compared the accuracy of ACD estimates obtained from ALS with those of two existing carbon maps that cover the study area. The first of these is a carbon map of the SAFE landscape and Maliau Basin derived from RapidEye satellite imagery (Pfeifer et al., 2016). The map has a resolution of 25×25m and makes use of textural and intensity information from four wavebands to model forest biomass (which we converted to carbon by applying a conversion factor of 0.47; Martin and Thomas, 2011). The second is a recently published consensus map of pantropical forest carbon stocks at 1 km resolution (Avitabile et al., 2016). It makes use of field data and high-resolution locally calibrated carbon maps to refine estimates from existing pantropical datasets obtained through satellite observations (Baccini et al., 2012; Saatchi et al., 2011).

Figure 2Relationship between field-estimated and modelled aboveground carbon density (ACD). Panel (a) shows the fit of the regionally calibrated ACD model (Eq. 8 in Sect. 3), which incorporates field-estimated basal area (BA) and wood density (WD), while panel (b) corresponds to Asner and Mascaro's (2014) general ACD model (Eq. 1 in Sect. 1). Panels (c)(d) illustrate the predictive accuracy of the regionally calibrated ACD model when field-measured BA and WD values are replaced with estimates derived from airborne laser scanning. In panel (c) BA and WD were estimated from top-of-canopy height (TCH) using Eqs. (9) and (11), respectively. In contrast, ACD estimates in panel (d) were obtained by modelling BA as a function of both TCH and canopy cover at 20 m aboveground following Eq. (10). In all panels, predicted ACD values are based on leave-one-out cross validation. Dashed lines correspond to a 1 : 1 relationship. Error bars are standard deviations and the RMSE of each comparison is printed in the bottom right-hand corner of the panels.


Figure 3Relationship between ALS-derived canopy cover at 20 m aboveground and top-of-canopy height. Panel (a) shows the distribution of the field plots with a line of best fit passing through the data, with error bars corresponding to standard deviations. Panel (b) illustrates how estimates of aboveground carbon density (ACD; obtained using Eq. 8 with Eqs. 10 and 11 as inputs) vary as a function of the two ALS metrics for the range of values observed across the forests of Sabah.


Figure 4Relationship between field-measured basal area and (a) top-of-canopy height and (b) canopy cover at 20 m aboveground as measured through airborne laser scanning. Error bars correspond to standard deviations.


Figure 5Relationship between community-weighted mean wood density (from field measurements) and top-of-canopy height (from airborne laser scanning). Error bars correspond to standard deviations.


To assess the accuracy of the two satellite products, we extracted ACD values from both carbon maps for all overlapping field plots and then compared field and satellite-derived estimates of ACD on the basis of RMSE and bias. For consistency with previous analyses, ACD values for SAFE project plots and those in riparian buffer zones were extracted at the individual plot level (i.e. 0.0625 ha scale) before being aggregated into the same blocks used for ALS model generation. In the case of Avitabile et al. (2016), we acknowledge that because of the large difference in resolution between the map and the field plots, comparisons between the two need to be made with care. This is particularly true when only a limited number of field plots are located within a given 1 km2 grid cell. To at least partially account for these difference in resolution when assessing agreement between Avitabile et al.'s (2016) map and the field data, we first averaged ACD values from field plots that fell within the same 1 km2 grid cell. We then compared satellite- and plot-based estimates of ACD for (i) all grid cells within which field plots were sampled, regardless of their number and size (n=135), and for (ii) a subset of grid cells for which at least five plots covering a cumulative area  1 ha were sampled in the field (n=8). The expectation is that grid cells for which a greater number of large plots have been surveyed should show closer alignment between satellite- and plot-based estimates of ACD.

3 Results

The regional model of ACD – parameterized using field estimates of wood density and basal area and ALS estimates of canopy height – was

(8) ACD Regional = 0.567 × TCH 0.554 × BA 1.081 × WD 0.186 .

The model had an RMSE of 19.0 Mg C ha−1 and a bias of 0.6 % (Fig. 2a; see Supplement S4 for confidence intervals on parameter estimates for all models reported here). The regional ACD model fit the data better than Asner and Mascaro's (2014) general model (i.e. Eq. 1 in Sect. 1), which had an RMSE of 32.0 Mg C ha−1 and tended to systematically underestimate ACD values (bias =7.1 %; Fig. 2b).

3.1 Basal area sub-models

When modelling BA in relation to TCH, we found the best-fit model to be

(9) BA = 1.112 × TCH .

In comparison, when BA was expressed as a function of both TCH and Coverresid we obtained the following model:

(10) BA = 1.287 × TCH 0.987 × 1 + 1.983 × Cover resid ,

where Coverresid=Cover20-11+e12.431×TCH-4.061 (Fig. 3). Of the two sub-models used to predict BA, Eq. (10) proved the better fit to the data (RMSE = 9.3 and 6.6 m2 ha−1, respectively; see Supplement S5), reflecting the fact that in our case BA was more closely related to canopy cover than TCH (Fig. 4).

3.2 Wood density sub-model

When modelling WD as a function of TCH, we found the best-fit model to be

(11) WD = 0.385 × TCH 0.097 .

Across the plot network WD showed a general tendency to increase with TCH (Fig. 5; RMSE of 0.056 g cm−3). However, the relationship was weak and Eq. (11) did not capture variation in WD equally well across the different forest types (see Supplement S5). In particular, heath forests at Sepilok, which have very high WD despite being much shorter than surrounding lowland dipterocarp forests (0.64 against 0.55 g cm−3), were poorly captured by the WD sub-model. We found no evidence to suggest that replacing TCH with canopy cover at 20 m aboveground would improve the accuracy of these estimates (see Supplement S5).

Figure 6Model errors (calculated for each individual plot as ACDobs-ACDpred/ACDobs×100) in relation to (a) field-estimated aboveground carbon density (ACD) and (b) plot size. Curves (±95 % shaded confidence intervals) were obtained by fitting linear models to log–log-transformed data. Black lines correspond to the regionally calibrated ACD model (Eq. 8 in Sect. 3). Orange lines show model errors when basal area (BA) was estimated from top-of-canopy height (TCH) using Eq. (9). In contrast, blue lines show model errors when BA was expressed as a function of both TCH and canopy cover at 20 m aboveground following Eq. (10). Vertical dashed lines along the horizontal axis show the distribution of the data (in panel b plot size values were jittered to avoid overlapping lines).


Figure 7Comparison between field-estimated aboveground carbon density (ACD) and satellite-derived estimates of ACD reported in (a) Pfeifer et al. (2016) and (b) Avitabile et al. (2016). In panel (b) large black points correspond to grid cells in Avitabile et al.'s (2016) pantropical biomass map for which at least five plots covering a cumulative area  1 ha were sampled in the field. By contrast, grid cells for which comparisons are based on less than five plots are depicted by small grey circles. Error bars correspond to standard deviations, while the RMSE of the satellite estimates is printed in the bottom right-hand corner of the panels (note that for panel b the RMSE in grey is that calculated across all plots, whereas that in black is based only on the subset of grid cells for which at least five plots covering a cumulative area  1 ha were sampled in the field). For comparison with ACD estimates obtained from airborne laser scanning, a kernel density plot fit to the points in Fig. 2d is displayed in the background.


3.3 Estimating aboveground carbon density from airborne laser scanning

When field-based estimates of BA and WD were replaced with ones derived from TCH using Eqs. (9) and (11), the regional ACD model generated unbiased estimates of ACD (bias =1.8 %). However, the accuracy of the model decreased substantially (RMSE = 48.1 Mg C ha−1; Fig. 2c). In particular, the average plot-level error was 21 % and remained relatively constant across the range of ACD values observed in the field data (yellow line in Fig. 6a). In contrast, when the combination of TCH and Cover20 was used to estimate BA through Eq. (10), we obtained more accurate estimates of ACD (RMSE = 39.3 Mg C ha−1, bias = 5.3 %; Fig. 2d). Moreover, in this instance plot-level errors showed a clear tendency to decrease in large and high-carbon-density plots (blue line in Fig. 6a), declining from an average 25.0 % at 0.1 ha scale to 19.5 % at 0.25 ha, 16.2 % at 0.5 ha and 13.4 % at 1 ha (blue line in Fig. 6b).

3.4 Comparison with satellite-derived estimates of aboveground carbon density

When compared to ALS-derived estimates of ACD, both satellite-based carbon maps of the study area showed much poorer agreement with field data (Fig. 7). Pfeifer et al.'s (2016) map covering the SAFE landscape and Maliau Basin systematically underestimated ACD (bias =36.9 %) and had an RMSE of 77.8 Mg C ha−1 (Fig. 7a). By contrast, Avitabile et al.'s (2016) pantropical map tended to overestimate carbon stocks. When we compared field and satellite estimates of ACD across all grid cells for which data were available we found that carbon stocks were overestimated by 111.2 % on average, with an RMSE of 100.1 Mg C ha−1 (grey circles in Fig. 7b). As expected, limiting this comparison to grid cells for which at least five plots covering a cumulative area  1 ha were sampled led to greater agreement between field and satellite estimates of ACD (large black circles in Fig. 7b). Yet the accuracy of the satellite-derived estimates of ACD remained much lower than that derived from ALS data (RMSE = 82.8 Mg C ha−1; bias = 59.3 %).

4 Discussion

We developed an area-based model for estimating aboveground carbon stock from ALS data that can be applied to mapping the lowland tropical forests of Borneo. We found that adding a canopy cover term to estimate BA to Asner and Mascaro's (2014) general model substantially improved its goodness-of-fit (Fig. 2c–d), as it allowed us to capture variation in stand basal area much more effectively compared to models parameterized solely using plot-averaged TCH. In this process, we also implemented an error propagation approach that allows various sources of uncertainty in ACD estimates to be incorporated into carbon mapping efforts. In the following sections we place our approach in the context of ongoing efforts to use remotely sensed data to monitor forest carbon stocks, starting with ALS-based approaches and then comparing these to satellite-based modelling. Finally, we end by discussing the implication of this work for the conservation of Borneo's forests.

4.1 Including canopy cover in the Asner and Mascaro (2014) carbon model

We found that incorporating a measure of canopy cover at 20 m aboveground in the Asner and Mascaro (2014) model improves its goodness-of-fit substantially without compromising its generality. Asner and Mascaro's (2014) model is grounded in forest and tree geometry, drawing its basis from allometric equations for estimating tree aboveground biomass such as that of Chave et al. (2014), in which a tree's biomass is expressed as a multiplicative function of its diameter, height and wood density: AGB=ρ0×WD×D2×Hρ1. By analogy, the carbon stock within a plot is related to the product of mean wood density, total basal area and top-of-canopy height (each raised to a power). Deriving this power-law function from knowledge of the tree size distribution and tree–biomass relationship is far from straightforward mathematically (Spriggs, 2015; Vincent et al., 2014), but this analogy seems to hold up well in a practical sense. When fitted to data from 14 forest types spanning aridity gradients in the Neotropics and Madagascar, Asner and Mascaro (2014) found that a single relationship applied to all forests types once regional differences in structure were incorporated as sub-models relating BA and WD to TCH. However, the model's fit depends critically on there being a close relationship between BA and TCH, as BA and ACD tend to be tightly coupled (ρ=0.93 in our case). Whilst that held true for the 14 forest types previously studied, in Bornean forest we found that the BA sub-models could be improved considerably by including canopy cover as an explanatory variable, particularly when it came to estimating BA in densely packed stands. This makes intuitive sense if one considers an open forest comprised of just a few trees; the crown area of each tree scales with its basal area, so the gap fraction at the ground level of a plot is negatively related to the basal area of its trees (Singh et al., 2016). A similar principle applies in denser forests, but in forests with multiple tiers formed by overlapping canopies such as those that occur in Borneo, the best-fitting relationship between gap fraction and basal area is no longer at ground level, but is instead further up the canopy (Coomes et al., 2017). Meyer et al. (2018) recently came to a very similar conclusion, showing that the cumulative crown area of emergent trees estimated at a height of 27 m aboveground using ALS data was strongly and linearly related to ACD across a diverse range of Neotropical forest types.

The functional form used to model BA in relation to TCH and residual forest cover (i.e. Eq. 10 presented above) was selected for two reasons: first, for a plot with average canopy cover for a given TCH, the model reduces to the classic model of Asner and Mascaro (2014), making comparisons straightforward. Second, simpler functional forms (e.g. ones relating BA directly to Cover20) were found to have very similar goodness-of-fit, but predicted unrealistically high ACD estimates for a small fraction of pixels when applied to mapping carbon across the landscape. This study is the first to formally introduce canopy cover into the modelling framework of Asner and Mascaro (2014), but several other studies have concluded that gap fraction is an important variable to include in multiple regression models of forest biomass (Colgan et al., 2012; Meyer et al., 2018; Ni-Meister et al., 2010; Pflugmacher et al., 2012; Singh et al., 2016). Regional calibration of the Asner and Mascaro (2014) model was necessary for the lowland forests of South-east Asia because dominance by dipterocarp species makes them structurally unique (Ghazoul, 2016; Jucker et al., 2018): trees in the region grow tall but have narrow stems for their height (Banin et al., 2012; Feldpausch et al., 2011), creating forests that have among the greatest carbon densities of any in the tropics (Avitabile et al., 2016; Sullivan et al., 2017).

The structural complexity and heterogeneity of Sabah's forests is one reason why even though accounting for canopy cover substantially improved the accuracy of our model (particularly in the case of tall, densely packed stands), a certain degree of error remains in the ACD estimates (Fig. 6). This error also reflects an inevitable trade-off between striving for generality and attempting to maximize accuracy when modelling ACD using ALS. In this regard, our modelling framework differs from the multiple-regression-with-model-selection approach that is typically adopted for estimating the ACD of tropical forests using ALS data (Chen et al., 2015; Clark et al., 2011; D'Oliveira et al., 2012; Drake et al., 2002; Hansen et al., 2015; Ioki et al., 2014; Jubanski et al., 2013; Réjou-Méchain et al., 2015; Singh et al., 2016). These studies, which build on 2 decades of research in temperate and boreal forests (Lefsky et al., 1999; Nelson et al., 1988; Popescu et al., 2011; Wulder et al., 2012), typically calculate between 5 and 25 summary statistics from the height distribution of ALS returns and explore the performance of models constructed using various combinations of those summary statistics as explanatory variables. The “best-supported” model is then selected from the list of competing models on offer by comparing relative performance using evaluation statistics such as R2, RMSE or AIC.

There is no doubt that selecting regression models in this way provides a solid basis for making model-assisted inferences about regional carbon stocks and their uncertainty (Ene et al., 2012; Gregoire et al., 2016). However, a well-recognized problem is that models tend to be idiosyncratic by virtue of local fine-tuning, so they cannot be applied more widely than the region for which they were calibrated and cannot be compared very easily with other studies. For example, it comes as no surprise that almost all publications identify mean height or some metric of upper-canopy height (e.g. 90th or 99th percentile of the height distribution) as being the strongest determinant of biomass. But different choices of height metric make these models difficult to compare. Other studies have included variance terms or measures of laser penetration to the lower canopy in an effort to improve goodness-of-fit. For instance, a combination of 75th quantile and variance of return heights proved effective in modelling ACD of selectively logged forests in Brazil (D'Oliveira et al., 2012). Similarly, a model developed for lower montane forests in Sabah included the proportion of last returns within 12 m of the ground (Ioki et al., 2014), while the proportions of returns in various height tiers were selected for the ALS carbon mapping of sub-montane forest in Tanzania (Hansen et al., 2015). Working with Asner and Mascaro's (2014) power-law model may well sacrifice goodness-of-fit compared with these locally tuned multiple regression models. However, it provides a systematic framework for the ALS modelling of forest carbon stocks, which we expect will prove hugely valuable for calibrating and validating the next generation of satellite sensors being designed specifically to monitor the world's forests.

4.2 Quantifying and propagating uncertainty

One of the most important applications of ACD estimation models is to infer carbon stocks within regions of interest. Carbon stock estimation has traditionally been achieved by networks of inventory plots designed to provide unbiased estimates of timber volumes within an acceptable level of uncertainty using well-established design-based approaches (Särndal et al., 1992). Forest inventories are increasingly supported by the collection of cost-effective auxiliary variables, such as ALS-estimated forest height and cover, that increase the precision of carbon stock estimation when used to construct regression models, which are in turn used to estimate carbon across areas where the auxiliary variables have been measured (e.g. McRoberts et al., 2013). But just as producing maps of our best estimates of carbon stocks across landscapes is critical to informing conservation and management strategies, so is the ability to provide robust estimates of the uncertainty associated with these products (Réjou-Méchain et al., 2017).

Assessing the degree of confidence which we place on a given estimate of ACD requires uncertainty to be quantified and propagated through all processes involved in the calculation of plot-level carbon stocks and statistical model fitting (Chen et al., 2015). Our Monte Carlo framework allows field measurement errors, geopositional errors and model uncertainty to be propagated in a straightforward and robust manner (Yanai et al., 2010). Our approach uses Réjou-Méchain et al.'s (2017) framework as a starting point for propagating errors associated with field measurements (e.g. stem diameter recording, wood density estimation) and allometric models (e.g. height–diameter relationships, tree biomass estimation) into plot-level estimates of ACD. We then combine these sources of uncertainty with those associated with co-location errors between field and ALS data and propagate these through the regression models we develop to estimate ACD from ALS metrics. This approach, which is fundamentally different to estimating uncertainty by comparing model predictions to validation field plots, is not widely used within the remote sensing community (e.g. Gonzalez et al., 2010) despite being the more appropriate technique for error propagation when there is uncertainty in field measurements (Chen et al., 2015).

Nevertheless, several sources of potential bias remain. Community-weighted wood density is only weakly related to ALS metrics and is estimated with large errors (Fig. 5). The fact that wood density cannot be measured remotely is well recognized, and the assumptions used to map wood density from limited field data have major implications for carbon maps produced by satellites (Mitchard et al., 2014). For Borneo, it may prove necessary to develop separate wood density sub-models for estimating carbon in heath forests versus other lowland forest types (see Fig. 5). Height estimation is another source of potential bias (Rutishauser et al., 2013): four published height–diameter curves for Sabah show similar fits for small trees (< 50 cm diameter) but diverge for large trees (Coomes et al., 2017), which contain most of the biomass (Bastin et al., 2015). Terrestrial laser scanning is likely to address this issue in the coming years, providing not only new and improved allometries for estimating tree height, but also much more robust field reference estimates of ACD from which to calibrate ALS-based models of forest carbon stocks (Calders et al., 2015; Gonzalez de Tanago Menaca et al., 2017). As this transition happens, careful consideration will also need to be given to differences in acquisition parameters among ALS campaigns and how these in turn influence ACD estimates derived from ALS metrics. While we found strong agreement between canopy metrics derived from the two airborne campaigns (Supplement S2), previous work has highlighted how decreasing ALS point density and changing footprint size can impact the retrieval of canopy parameters (Gobakken and Næsset, 2008; Roussel et al., 2017). In this regard new approaches designed to explicitly correct for differences among ALS flight specifications (e.g. Roussel et al., 2017) offer great promise for minimizing this source of bias.

Lastly, another key issue influencing uncertainty in ACD estimates derived from ALS data is the size of the field plots used to calibrate and validate prediction models. As a rule of thumb, the smaller the field plots the poorer the fit between field estimates of ACD and ALS-derived canopy metrics (Asner and Mascaro, 2014; Ruiz et al., 2014; Watt et al., 2013). Aside from the fact that small plots inevitably capture a greater degree of heterogeneity in ACD compared to larger ones (leading to more noise around regression lines), they are also much more likely to suffer from errors associated with poor alignment between airborne and field data, as well as exhibiting strong edge effects (e.g. large trees whose crowns straddle the plot boundary). As expected, for our best-fitting model of ACD we found that plot-level errors tended to decrease with plot size (blue curve in Fig. 6b), going from 25.0 % at 0.1 ha scale to 13.4 % for 1 ha plots. This result is remarkably consistent with previous theoretical and empirical work conducted across the tropics, which reported mean errors of around 25–30 % at 0.1 ha scale and approximately 10–15 % at 1 ha resolution (Asner and Mascaro, 2014; Zolkos et al., 2013). These results have led to the general consensus that 1 ha plots should become the standard for calibrating against ALS data. That being said, because there is a trade-off between the number of plots one can establish and their size, working with 1 ha plots inevitably comes at the cost of replication and representativeness. As such, in some cases it may be preferable to sacrifice some precision (e.g. by working with 0.25 ha plots, which in our case had a mean error of 19.5 %) in order to gain a better representation of the wider landscape – so long as uncertainty in ACD is fully propagate throughout.

4.3 Comparison with satellite-derived maps

Our results show that when compared to independent field data, existing satellite products systematically underestimate or overestimate ACD (depending on the product; Fig. 7). While directly comparing satellite-derived estimates with independent field data in not entirely straightforward, particularly when the resolution of the map is much coarser than that of the field plots (Réjou-Méchain et al., 2014), as is the case with Avitabile et al. (2016), it does appear that ALS is able to provide much more robust and accurate estimates of ACD and its heterogeneity within the landscape than what is possible with current satellite sensors. This is unsurprising given that in contrast to optical imagery, which only captures the surface of the canopy, ALS data provide high-resolution information on the 3-D structure of canopies, which directly relates to ACD. However, ALS data are limited in their temporal and spatial coverage due to high operational costs. Consequently, there is a growing need to focus on fusing ALS-derived maps of ACD with satellite data to advance our ability to map forest carbon stocks across large spatial scales and through time (e.g. Asner et al., 2018). In this regard, NASA plans to start making high-resolution laser-ranging observations from the international space station in 2018 as part of the GEDI mission, while the ESA BIOMASS mission will use P-band synthetic aperture radar to monitor forests from space from 2021. Pantropical monitoring of forest carbon using data from a combination of space-borne sensors is fast approaching, and regional carbon equations derived from ALS data such as the one we develop here will be critical to calibrate and validate these efforts.

5 Conclusions

Since the 1970s Borneo has lost more than 60 % of its old-growth forests, the majority of which have been replaced by large-scale industrial palm oil plantations (Gaveau et al., 2014, 2016). Nowhere else has this drastic transformation of the landscape been more evident than in the Malaysian state of Sabah, where forest clearing rates have been among the highest across the entire region (Osman et al., 2012). Certification bodies such as the Roundtable on Sustainable Palm Oil (RSPO) have responded to criticisms by adopting policies that prohibit planting on land designated as high conservation value (HCV) and have recently proposed to supplement the HCV approach with high carbon stock (HCS) assessments that would restrict the expansion of palm oil plantations onto carbon-dense forests. Yet enforcing these policies requires an accurate and spatially detailed understating of how carbon stocks are distributed cross the entire state, something which is currently lacking. With the view of halting the further deforestation of carbon-dense old-growth forests and generating the necessary knowledge to better manage its forests into the future, in 2016 the Sabah state government commissioned CAO to deliver a high-resolution ALS-based carbon map of the entire state (Asner et al., 2018). The regional carbon model we develop here underpins this initiative (Asner et al., 2018; Nunes et al., 2017) and more generally will contribute to ongoing efforts to use remote sensing tools to provide solutions for identifying and managing the more than 500 million ha of tropical lands that are currently degraded (Lamb et al., 2005).

Data availability

The data supporting the results of this paper have been archived on the NERC Open Research Archive website (, last access: 20 June 2018).


The supplement related to this article is available online at:

Author contributions

DAC and YAT coordinated the NERC airborne campaign, while GPA led the CAO airborne surveys of Sabah. TJ and DAC designed the study, with input from GPA and PGB; TJ, MD, PGB and NRV processed the airborne imagery, while other authors contributed field data; TJ analysed the data, with input from DAC, GPA, PGB and CDP; TJ wrote the first draft of the paper, with all other authors contributing to revisions.

Competing interests

The authors declare that they have no conflict of interest.


This study was funded by the UK Natural Environment Research Council's (NERC) Human Modified Tropical Forests research programme (grant numbers NE/K016377/1 and NE/K016407/1 awarded to the BALI and LOMBOK consortiums, respectively). We are grateful to NERC's Airborne Research Facility and Data Analysis Node for conducting the survey and preprocessing the airborne data and to Abdullah Ghani for manning the GPS base station. David A. Coomes was supported in part by an International Academic Fellowship from the Leverhulme Trust. The Carnegie Airborne Observatory portion of the study was supported by the UN Development Programme, the Avatar Alliance Foundation, the Roundtable on Sustainable Palm Oil, the World Wildlife Fund and the Rainforest Trust. The Carnegie Airborne Observatory is made possible by grants and donations to Gregory P. Asner from the Avatar Alliance Foundation, the Margaret A. Cargill Foundation, the David and Lucile Packard Foundation, the Gordon and Betty Moore Foundation, the Grantham Foundation for the Protection of the Environment, the W. M. Keck Foundation, the John D. and Catherine T. MacArthur Foundation, the Andrew Mellon Foundation, Mary Anne Nyburg Baker and G. Leonard Baker Jr., and William R. Hearst III. The SAFE project was supported by the Sime Darby Foundation. We acknowledge the SAFE management team, Maliau Basin Management Committee, Danum Valley Management Committee, South East Asia Rainforest Research Partnership, Sabah Foundation, Benta Wawasan, the State Secretary, the Sabah Chief Minister's Departments, the Sabah Forestry Department, the Sabah Biodiversity Centre and the Economic Planning Unit for their support, access to the field sites and permission to carry out fieldwork in Sabah. Field data collection at Sepilok was supported by an ERC Advanced Grant (291585, T-FORCES) awarded to Oliver L. Phillips, who is also a Royal Society Wolfson Research Merit Award holder. Martin Svátek was funded through a grant from the Ministry of Education, Youth and Sports of the Czech Republic (grant number INGO II LG15051), and Jakub Kvasnica was funded through an IGA grant (grant number LDF_VP_2015038). We are grateful to the many field assistants who contributed to data collection.

Edited by: Nobuhito Ohte
Reviewed by: two anonymous referees


Agrawal, A., Nepstad, D., and Chhatre, A.: Reducing emissions from deforestation and forest degradation, Annu. Rev. Environ. Resour., 36, 373–396,, 2011. 

Anderson-Teixeira, K. J., Davies, S. J., Bennett, A. C., Gonzalez-Akre, E. B., Muller-Landau, H. C., Joseph Wright, S., Abu Salim, K., Almeyda Zambrano, A. M., Alonso, A., Baltzer, J. L., Basset, Y., Bourg, N. A., Broadbent, E. N., Brockelman, W. Y., Bunyavejchewin, S., Burslem, D. F. R. P., Butt, N., Cao, M., Cardenas, D., Chuyong, G. B., Clay, K., Cordell, S., Dattaraja, H. S., Deng, X., Detto, M., Du, X., Duque, A., Erikson, D. L., Ewango, C. E. N., Fischer, G. A., Fletcher, C., Foster, R. B., Giardina, C. P., Gilbert, G. S., Gunatilleke, N., Gunatilleke, S., Hao, Z., Hargrove, W. W., Hart, T. B., Hau, B. C. H., He, F., Hoffman, F. M., Howe, R. W., Hubbell, S. P., Inman-Narahari, F. M., Jansen, P. A., Jiang, M., Johnson, D. J., Kanzaki, M., Kassim, A. R., Kenfack, D., Kibet, S., Kinnaird, M. F., Korte, L., Kral, K., Kumar, J., Larson, A. J., Li, Y., Li, X., Liu, S., Lum, S. K. Y., Lutz, J. A., Ma, K., Maddalena, D. M., Makana, J.-R., Malhi, Y., Marthews, T., Mat Serudin, R., McMahon, S. M., McShea, W. J., Memiaghe, H. R., Mi, X., Mizuno, T., Morecroft, M., Myers, J. A., Novotny, V., de Oliveira, A. A., Ong, P. S., Orwig, D. A., Ostertag, R., den Ouden, J., Parker, G. G., Phillips, R. P., Sack, L., Sainge, M. N., Sang, W., Sri-ngernyuang, K., Sukumar, R., Sun, I.-F., Sungpalee, W., Suresh, H. S., Tan, S., Thomas, S. C., Thomas, D. W., Thompson, J., Turner, B. L., Uriarte, M., Valencia, R., Vallejo, M. I., Vicentini, A., Vrška, T., Wang, X., Wang, X., Weiblen, G., Wolf, A., Xu, H., Yap, S., and Zimmerman, J.: CTFS-ForestGEO: a worldwide network monitoring forests in an era of global change, Glob. Change Biol., 21, 528–549,, 2015. 

Asner, G. P. and Mascaro, J.: Mapping tropical forest carbon: calibrating plot estimates to a simple LiDAR metric, Remote Sens. Environ., 140, 614–624,, 2014. 

Asner, G. P., Rudel, T. K., Aide, T. M., DeFries, R., and Emerson, R.: A contemporary assessment of change in humid tropical forests, Conserv. Biol., 23, 1386–1395,, 2009. 

Asner, G. P., Powell, G. V. N., Mascaro, J., Knapp, D. E., Clark, J. K., Jacobson, J., Kennedy-Bowdoin, T., Balaji, A., Paez-Acosta, G., Victoria, E., Secada, L., Valqui, M., and Hughes, R. F.: High-resolution forest carbon stocks and emissions in the Amazon, P. Natl. Acad. Sci. USA, 107, 16738–16742,, 2010. 

Asner, G. P., Knapp, D. E., Boardman, J., Green, R. O., Kennedy-Bowdoin, T., Eastwood, M., Martin, R. E., Anderson, C., and Field, C. B.: Carnegie Airborne Observatory-2: increasing science data dimensionality via high-fidelity multi-sensor fusion, Remote Sens. Environ., 124, 454–465,, 2012. 

Asner, G. P., Knapp, D. E., Martin, R. E., Tupayachi, R., Anderson, C. B., Mascaro, J., Sinca, F., Chadwick, K. D., Higgins, M., Farfan, W., Llactayo, W., and Silman, M. R.: Targeted carbon conservation at national scales with high-resolution monitoring, P. Natl. Acad. Sci. USA, 111, E5016–E5022,, 2014. 

Asner, G. P., Brodrick, P. G., Philipson, C., Vaughn, N. R., Martin, R. E., Knapp, D. E., Heckler, J., Evans, L. J., Jucker, T., Goossens, B., Stark, D. J., Reynolds, G., Ong, R., Renneboog, N., Kugan, F., and Coomes, D. A.: Mapped aboveground carbon stocks to advance forest conservation and recovery in Malaysian Borneo, Biol. Conserv., 217, 289–310,, 2018. 

Avitabile, V., Herold, M., Heuvelink, G. B. M., Lewis, S. L., Phillips, O. L., Asner, G. P., Armston, J., Asthon, P., Banin, L. F., Bayol, N., Berry, N., Boeckx, P., de Jong, B., DeVries, B., Girardin, C., Kearsley, E., Lindsell, J. A., Lopez-Gonzalez, G., Lucas, R., Malhi, Y., Morel, A., Mitchard, E., Nagy, L., Qie, L., Quinones, M., Ryan, C. M., Slik, F., Sunderland, T., Vaglio Laurin, G., Valentini, R., Verbeeck, H., Wijaya, A., and Willcock, S.: An integrated pan-tropical biomass map using multiple reference datasets, Glob. Change Biol., 22, 1406–1420,, 2016. 

Baccini, A., Goetz, S. J., Walker, W. S., Laporte, N. T., Sun, M., Sulla-Menashe, D., Hackler, J., Beck, P. S. A., Dubayah, R., Friedl, M. A., Samanta, S., and Houghton, R. A.: Estimated carbon dioxide emissions from tropical deforestation improved by carbon-density maps, Nat. Clim. Change, 2, 182–185,, 2012. 

Banin, L., Feldpausch, T. R., Phillips, O. L., Baker, T. R., Lloyd, J., Affum-Baffoe, K., Arets, E. J. M. M., Berry, N. J., Bradford, M., Brienen, R. J. W., Davies, S., Drescher, M., Higuchi, N., Hilbert, D. W., Hladik, A., Iida, Y., Salim, K. A., Kassim, A. R., King, D. A., Lopez-Gonzalez, G., Metcalfe, D., Nilus, R., Peh, K. S. H., Reitsma, J. M., Sonké, B., Taedoumg, H., Tan, S., White, L., Wöll, H., and Lewis, S. L.: What controls tropical forest architecture? Testing environmental, structural and floristic drivers, Glob. Ecol. Biogeogr., 21, 1179–1190,, 2012. 

Baskerville, G. L.: Use of logarithmic regression in the estimation of plant biomass, Can. J. For. Res., 2, 49–53, 1972. 

Bastin, J.-F., Barbier, N., Réjou-Méchain, M., Fayolle, A., Gourlet-Fleury, S., Maniatis, D., de Haulleville, T., Baya, F., Beeckman, H., Beina, D., Couteron, P., Chuyong, G., Dauby, G., Doucet, J.-L., Droissart, V., Dufrêne, M., Ewango, C., Gillet, J. F., Gonmadje, C. H., Hart, T., Kavali, T., Kenfack, D., Libalah, M., Malhi, Y., Makana, J.-R., Pélissier, R., Ploton, P., Serckx, A., Sonké, B., Stevart, T., Thomas, D. W., De Cannière, C., and Bogaert, J.: Seeing Central African forests through their largest trees, Sci. Rep., 5, 13156,, 2015. 

Calders, K., Newnham, G., Burt, A., Murphy, S., Raumonen, P., Herold, M., Culvenor, D., Avitabile, V., Disney, M., Armston, J., and Kaasalainen, M.: Nondestructive estimates of above-ground biomass using terrestrial laser scanning, Methods Ecol. Evol., 6, 198–208,, 2015. 

Carlson, K. M., Curran, L. M., Asner, G. P., Pittman, A. M., Trigg, S. N., and Marion Adeney, J.: Carbon emissions from forest conversion by Kalimantan oil palm plantations, Nat. Clim. Change, 3, 283–287,, 2012a. 

Carlson, K. M., Curran, L. M., Ratnasari, D., Pittman, A. M., Soares-Filho, B. S., Asner, G. P., Trigg, S. N., Gaveau, D. a, Lawrence, D., and Rodrigues, H. O.: Committed carbon emissions, deforestation, and community land conversion from oil palm plantation expansion in West Kalimantan, Indonesia., P. Natl. Acad. Sci. USA, 109, 1–6,, 2012b. 

Chave, J., Condit, R., Aguilar, S., Hernandez, A., Lao, S., and Perez, R.: Error propagation and scaling for tropical forest biomass estimates, Philos. T. R. Soc. B, 359, 409–20,, 2004. 

Chave, J., Coomes, D. A., Jansen, S., Lewis, S. L., Swenson, N. G., and Zanne, A. E.: Towards a worldwide wood economics spectrum, Ecol. Lett., 12, 351–366,, 2009. 

Chave, J., Réjou-Méchain, M., Búrquez, A., Chidumayo, E., Colgan, M. S., Delitti, W. B. C., Duque, A., Eid, T., Fearnside, P. M., Goodman, R. C., Henry, M., Martínez-Yrízar, A., Mugasha, W. A., Muller-Landau, H. C., Mencuccini, M., Nelson, B. W., Ngomanda, A., Nogueira, E. M., Ortiz-Malavassi, E., Pélissier, R., Ploton, P., Ryan, C. M., Saldarriaga, J. G., and Vieilledent, G.: Improved allometric models to estimate the aboveground biomass of tropical trees, Glob. Change Biol., 20, 3177–3190,, 2014. 

Chen, Q., Vaglio Laurin, G., and Valentini, R.: Uncertainty of remotely sensed aboveground biomass over an African tropical forest: Propagating errors from trees to plots to pixels, Remote Sens. Environ., 160, 134–143,, 2015. 

Clark, M. L., Roberts, D. A., Ewel, J. J., and Clark, D. B.: Estimation of tropical rain forest aboveground biomass with small-footprint lidar and hyperspectral sensors, Remote Sens. Environ., 115, 2931–2942,, 2011. 

Colgan, M. S., Asner, G. P., Levick, S. R., Martin, R. E., and Chadwick, O. A.: Topo-edaphic controls over woody plant biomass in South African savannas, Biogeosciences, 9, 1809–1821,, 2012. 

Coomes, D. A., Dalponte, M., Jucker, T., Asner, G. P., Banin, L. F., Burslem, D. F. R. P., Lewis, S. L., Nilus, R., Phillips, O. L., Phuag, M.-H., Qiee, L., Phua, M.-H., and Qie, L.: Area-based vs tree-centric approaches to mapping forest carbon in Southeast Asian forests with airborne laser scanning data, Remote Sens. Environ., 194, 77–88,, 2017. 

Dent, D. H., Bagchi, R., Robinson, D., Majalap-Lee, N., and Burslem, D. F. R. P.: Nutrient fluxes via litterfall and leaf litter decomposition vary across a gradient of soil nutrient supply in a lowland tropical rain forest, Plant Soil, 288, 197–215,, 2006. 

DeWalt, S. J. S. S. J., Ickes, K., Nilus, R., Harms, K. E. K., and Burslem, D. D. F. R. P.: Liana habitat associations and community structure in a Bornean lowland tropical forest, Plant Ecol., 186, 203–216,, 2006. 

D'Oliveira, M. V. N., Reutebuch, S. E., McGaughey, R. J., and Andersen, H.-E.: Estimating forest biomass and identifying low-intensity logging areas using airborne scanning lidar in Antimary State Forest, Acre State, Western Brazilian Amazon, Remote Sens. Environ., 124, 479–491,, 2012. 

Dormann, C. F., Elith, J., Bacher, S., Buchmann, C., Carl, G., Carré, G., Marquéz, J. R. G., Gruber, B., Lafourcade, B., Leitão, P. J., Münkemüller, T., McClean, C., Osborne, P. E., Reineking, B., Schröder, B., Skidmore, A. K., Zurell, D., and Lautenbach, S.: Collinearity: A review of methods to deal with it and a simulation study evaluating their performance, Ecography, 36, 27–46,, 2013. 

Drake, J. B., Dubayah, R. O., Clark, D. B., Knox, R. G., Blair, J. B. B., Hofton, M. A., Chazdon, R. L., Weishampel, J. F., and Prince, S. D.: Estimation of tropical forest structural characteristics using large-footprint lidar, Remote Sens. Environ., 79, 305–319,, 2002. 

Duncanson, L. I., Dubayah, R. O., Cook, B. D., Rosette, J., and Parker, G.: The importance of spatial detail: Assessing the utility of individual crown information and scaling approaches for lidar-based biomass density estimation, Remote Sens. Environ., 168, 102–112,, 2015. 

Ene, L. T., Næsset, E., Gobakken, T., Gregoire, T. G., Ståhl, G., and Nelson, R.: Assessing the accuracy of regional LiDAR-based biomass estimation using a simulation approach, Remote Sens. Environ., 123, 579–592,, 2012. 

Ewers, R. M., Didham, R. K., Fahrig, L., Ferraz, G., Hector, A., Holt, R. D., and Turner, E. C.: A large-scale forest fragmentation experiment: the Stability of Altered Forest Ecosystems Project, Philos. T. R. Soc. B, 366, 3292–3302,, 2011. 

Feldpausch, T. R., Banin, L., Phillips, O. L., Baker, T. R., Lewis, S. L., Quesada, C. A., Affum-Baffoe, K., Arets, E. J. M. M., Berry, N. J., Bird, M., Brondizio, E. S., de Camargo, P., Chave, J., Djagbletey, G., Domingues, T. F., Drescher, M., Fearnside, P. M., França, M. B., Fyllas, N. M., Lopez-Gonzalez, G., Hladik, A., Higuchi, N., Hunter, M. O., Iida, Y., Salim, K. A., Kassim, A. R., Keller, M., Kemp, J., King, D. A., Lovett, J. C., Marimon, B. S., Marimon-Junior, B. H., Lenza, E., Marshall, A. R., Metcalfe, D. J., Mitchard, E. T. A., Moran, E. F., Nelson, B. W., Nilus, R., Nogueira, E. M., Palace, M., Patiño, S., Peh, K. S.-H., Raventos, M. T., Reitsma, J. M., Saiz, G., Schrodt, F., Sonké, B., Taedoumg, H. E., Tan, S., White, L., Wöll, H., and Lloyd, J.: Height-diameter allometry of tropical forest trees, Biogeosciences, 8, 1081–1106,, 2011. 

Fox, J. E. D.: A handbook to Kabili-Sepilok Forest Reserve, Sabah Forest Record No. 9, Borneo Literature Bureau, Kuching, Sarawak, Malaysia, 1973. 

Gaveau, D. L. A., Sloan, S., Molidena, E., Yaen, H., Sheil, D., Abram, N. K., Ancrenaz, M., Nasi, R., Quinones, M., Wielaard, N., and Meijaard, E.: Four decades of forest persistence, clearance and logging on Borneo, PLoS One, 9, e101654,, 2014. 

Gaveau, D. L. A., Sheil, D., Husnayaen, Salim, M. A., Arjasakusuma, S., Ancrenaz, M., Pacheco, P., and Meijaard, E.: Rapid conversions and avoided deforestation: examining four decades of industrial plantation expansion in Borneo, Sci. Rep., 6, 32017,, 2016. 

Ghazoul, J.: Dipterocarp Biology, Ecology, and Conservation, Oxford University Press, Oxford, UK, 2016. 

Gibbs, H. K., Brown, S., Niles, J. O., and Foley, J. A.: Monitoring and estimating tropical forest carbon stocks: making REDD a reality, Environ. Res. Lett., 2, 045023,, 2007. 

Gobakken, T. and Næsset, E.: Assessing effects of laser point density, ground sampling intensity, and field sample plot size on biophysical stand properties derived from airborne laser scanner data, Can. J. Forest Res., 38, 1095–1109,, 2008. 

Gonzalez, P., Asner, G. P., Battles, J. J., Lefsky, M. A., Waring, K. M., and Palace, M.: Forest carbon densities and uncertainties from Lidar, QuickBird, and field measurements in California, Remote Sens. Environ., 114, 1561–1575,, 2010. 

Gonzalez de Tanago Menaca, J., Lau, A., Bartholomeusm, H., Herold, M., Avitabile, V., Raumonen, P., Martius, C., Goodman, R., Disney, M., Manuri, S., Burt, A., and Calders, K.: Estimation of above-ground biomass of large tropical trees with Terrestrial LiDAR, Methods Ecol. Evol., 8, 223–234,, 2017. 

Gregoire, T. G., Næsset, E., McRoberts, R. E., Ståhl, G., Andersen, H.-E., Gobakken, T., Ene, L., and Nelson, R.: Statistical rigor in LiDAR-assisted estimation of aboveground forest biomass, Remote Sens. Environ., 173, 98–108,, 2016. 

Hansen, E., Gobakken, T., Bollandsås, O., Zahabu, E., and Næsset, E.: Modeling aboveground biomass in dense tropical submontane rainforest using airborne laser scanner data, Remote Sens., 7, 788–807,, 2015. 

Hansen, M. C., Potapov, P. V, Moore, R., Hancher, M., Turubanova, S. A., Tyukavina, A., Thau, D., Stehman, S. V, Goetz, S. J., Loveland, T. R., Kommareddy, A., Egorov, A., Chini, L., Justice, C. O., and Townshend, J. R. G.: High-resolution global maps of 21st-century forest cover change, Science, 342, 850–853,, 2013. 

Houghton, R. A., Byers, B., and Nassikas, A. A.: A role for tropical forests in stabilizing atmospheric CO2, Nat. Clim. Change, 5, 1022–1023,, 2015. 

Ioki, K., Tsuyuki, S., Hirata, Y., Phua, M.-H., Wong, W. V. C., Ling, Z.-Y., Saito, H., and Takao, G.: Estimating above-ground biomass of tropical rainforest of different degradation levels in Northern Borneo using airborne LiDAR, Forest Ecol. Manag., 328, 335–341,, 2014. 

Jubanski, J., Ballhorn, U., Kronseder, K., Franke, J., and Siegert, F.: Detection of large above-ground biomass variability in lowland forest ecosystems by airborne LiDAR, Biogeosciences, 10, 3917–3930,, 2013. 

Jucker, T., Bongalov, B., Burslem, D. F. R. P., Nilus, R., Dalponte, M., Lewis, S. L., Phillips, O. L., Qie, L., and Coomes, D. A.: Topography shapes the structure, composition and function of tropical forest landscapes, Ecol. Lett., 21, 989–1000,, 2018. 

Lamb, D., Erskine, P. D., and Parrotta, J. A.: Restoration of degraded tropical forest landscapes, Science, 310, 1628–1632, 2005. 

Laurin, G. V., Cheung-Wai Chan, J., Chen, Q., Lindsell, J. A., Coomes, D. A., Guerriero, L., Del Frate, F., Miglietta, F., and Valentini, R.: Biodiversity mapping in a tropical West African forest with airborne hyperspectral data, PLoS One, 9, e97910,, 2014. 

Lefsky, M., Cohen, W., Acker, S., and Parker, G.: Lidar remote sensing of the canopy structure and biophysical properties of Douglas-fir western hemlock forests, Remote Sens. Environ., 70, 339–361, 1999. 

Lefsky, M. A., Cohen, W. B., Parker, G. G., and Harding, D. J.: Lidar remote sensing for ecosystem studies, BioScience, 52, 19–30,[0019:LRSFES]2.0.CO;2, 2002. 

Malhi, Y. and Wright, J.: Spatial patterns and recent trends in the climate of tropical rainforest regions, Philos. T. R. Soc. B, 359, 311–329,, 2004. 

Martin, A. R. and Thomas, S. C.: A reassessment of carbon content in tropical trees, PLoS One, 6, e23533,, 2011. 

McRoberts, R. E., Næsset, E., and Gobakken, T.: Inference for lidar-assisted estimation of forest growing stock volume, Remote Sens. Environ., 128, 268–275,, 2013. 

Meyer, V., Saatchi, S., Clark, D. B., Keller, M., Vincent, G., Ferraz, A., Espírito-Santo, F., d'Oliveira, M. V. N., Kaki, D., and Chave, J.: Canopy area of large trees explains aboveground biomass variations across neotropical forest landscapes, Biogeosciences, 15, 3377–3390,, 2018. 

Mitchard, E. T. A., Feldpausch, T. R., Brienen, R. J. W., Lopez-Gonzalez, G., Monteagudo, A., Baker, T. R., Lewis, S. L., Lloyd, J., Quesada, C. A., Gloor, M., ter Steege, H., Meir, P., Alvarez, E., Araujo-Murakami, A., Aragão, L. E. O. C., Arroyo, L., Aymard, G., Banki, O., Bonal, D., Brown, S., Brown, F. I., Cerón, C. E., Chama Moscoso, V., Chave, J., Comiskey, J. A., Cornejo, F., Corrales Medina, M., Da Costa, L., Costa, F. R. C., Di Fiore, A., Domingues, T. F., Erwin, T. L., Frederickson, T., Higuchi, N., Honorio Coronado, E. N., Killeen, T. J., Laurance, W. F., Levis, C., Magnusson, W. E., Marimon, B. S., Marimon Junior, B. H., Mendoza Polo, I., Mishra, P., Nascimento, M. T., Neill, D., Núñez Vargas, M. P., Palacios, W. A., Parada, A., Pardo Molina, G., Peña-Claros, M., Pitman, N., Peres, C. A., Poorter, L., Prieto, A., Ramirez-Angulo, H., Restrepo Correa, Z., Roopsind, A., Roucoux, K. H., Rudas, A., Salomão, R. P., Schietti, J., Silveira, M., de Souza, P. F., Steininger, M. K., Stropp, J., Terborgh, J., Thomas, R., Toledo, M., Torres-Lezama, A., van Andel, T. R., van der Heijden, G. M. F., Vieira, I. C. G., Vieira, S., Vilanova-Torre, E., Vos, V. A., Wang, O., Zartman, C. E., Malhi, Y., and Phillips, O. L.: Markedly divergent estimates of Amazon forest carbon density from ground plots and satellites, Glob. Ecol. Biogeogr., 23, 935–946,, 2014. 

Morales, J. L. and Nocedal, J.: Remark on Algorithm 778: L-BFGS-B, FORTRAN routines for large scale bound constrained optimization, ACM Trans. Math. Softw., 38, available at: (last access: 20 June 2018), 2011. 

Nelson, R., Krabill, W., and Tonelli, J.: Estimating forest biomass and volume using airborne laser data, Remote Sens. Environ., 24, 247–267, 1988. 

Nilus, R., Maycock, C., Majalap-Lee, N., and Burslem, D.: Nutrient limitation of tree seedling growth in three soil types found in Sabah, J. Trop. For. Sci., 23, 133–142, 2011. 

Ni-Meister, W., Lee, S., Strahler, A. H., Woodcock, C. E., Schaaf, C., Yao, T., Ranson, K. J., Sun, G., and Blair, J. B.: Assessing general relationships between aboveground biomass and vegetation structure parameters for improved carbon estimate from lidar remote sensing, J. Geophys. Res.-Biogeosci., 115, G00E11,, 2010. 

Nunes, M. H., Ewers, R. M., Turner, E. C., and Coomes, D. A.: Mapping aboveground carbon in oil palm plantations using LiDAR: A comparison of tree-centric versus area-based approaches, Remote Sens., 9, 1–13,, 2017. 

Osman, R., Phua, M.-H., Ling, Z. Y., and Kamlun, K. U.: Monitoring of deforestation rate and trend in Sabah between 1990 and 2008 using multitemporal landsat data, J. For. Environ. Sci., 28, 144–151,, 2012. 

Pan, Y., Birdsey, R. A., Fang, J., Houghton, R. A., Kauppi, P. E., Kurz, W. A., Phillips, O. L., Shvidenko, A., Lewis, S. L., Canadell, J. G., Ciais, P., Jackson, R. B., Pacala, S. W., McGuire, A. D., Piao, S., Rautiainen, A., Sitch, S., Hayes, D., Canadell, J. G., Khatiwala, S., Primeau, F., Hall, T., Quéré, C. Le, Dixon, R. K., Kauppi, P. E., Kurz, W. A., Stinson, G., Rampley, G. J., Dymond, C. C., Neilson, E. T., Stinson, G., Birdsey, R. A., Pregitzer, K., Lucier, A., Kauppi, P. E., Pan, Y., Pan, Y., Birdsey, R. A., Hom, J., McCullough, K., van Mantgem, P. J., Breshears, D. D., Ciais, P., Fang, J., Chen, A., Peng, C., Zhao, S., Ci, L., Lewis, S. L., Phillips, O. L., Gloor, M., Lewis, S. L., Lloyd, J., Sitch, S., Mitchard, E. T. A., Laurance, W. F., Houghton, R. A., Friedlingstein, P., Tarnocai, C., Hooijer, A., Page, S. E., Rieley, J. O., Banks, C. J., McGuire, A. D., Goodale, C. L., Sarmiento, J. L., Schulze, E. D., Pacala, S. W., Phillips, O. L., Metsaranta, J. M., Kurz, W. A., Neilson, E. T., Stinson, G., Zhao, M., Running, S. W., and Houghton, R. A.: A large and persistent carbon sink in the world's forests, Science, 333, 988–993,, 2011. 

Pfeifer, M., Kor, L., Nilus, R., Turner, E., Cusack, J., Lysenko, I., Khoo, M., Chey, V. K., Chung, A. C., and Ewers, R. M.: Mapping the structure of Borneo's tropical forests across a degradation gradient, Remote Sens. Environ., 176, 84–97,, 2016. 

Pflugmacher, D., Cohen, W. B., and Kennedy, R. E.: Using Landsat-derived disturbance history (1972–2010) to predict current forest structure, Remote Sens. Environ., 122, 146–165,, 2012. 

Phillips, O. L., Malhi, Y., Higuchi, N., Laurance, W. F., Núñez, P. V, Vásquez, R. M., Laurance, S. G., Ferreira, L. V., Stern, M., Brown, S., and Grace, J.: Changes in the carbon balance of tropical forests: Evidence from long-term plots, Science, 282, 439–442, 1998. 

Popescu, S. C., Zhao, K., Neuenschwander, A., and Lin, C.: Satellite lidar vs. small footprint airborne lidar: Comparing the accuracy of aboveground biomass estimates and forest structure metrics at footprint level, Remote Sens. Environ., 115, 1–12,, 2011. 

Quesada, C. A., Phillips, O. L., Schwarz, M., Czimczik, C. I., Baker, T. R., Patiño, S., Fyllas, N. M., Hodnett, M. G., Herrera, R., Almeida, S., Alvarez Dávila, E., Arneth, A., Arroyo, L., Chao, K. J., Dezzeo, N., Erwin, T., di Fiore, A., Higuchi, N., Honorio Coronado, E., Jimenez, E. M., Killeen, T., Lezama, A. T., Lloyd, G., López-González, G., Luizão, F. J., Malhi, Y., Monteagudo, A., Neill, D. A., Núñez Vargas, P., Paiva, R., Peacock, J., Peñuela, M. C., Peña Cruz, A., Pitman, N., Priante Filho, N., Prieto, A., Ramírez, H., Rudas, A., Salomão, R., Santos, A. J. B., Schmerler, J., Silva, N., Silveira, M., Vásquez, R., Vieira, I., Terborgh, J., and Lloyd, J.: Basin-wide variations in Amazon forest structure and function are mediated by both soils and climate, Biogeosciences, 9, 2203–2246,, 2012. 

R Core Development Team: R: A language and environment for statistical computing, R Foundation for Statistical Computing, Vienna, Austria, 2016. 

Réjou-Méchain, M., Muller-Landau, H. C., Detto, M., Thomas, S. C., Le Toan, T., Saatchi, S. S., Barreto-Silva, J. S., Bourg, N. A., Bunyavejchewin, S., Butt, N., Brockelman, W. Y., Cao, M., Cárdenas, D., Chiang, J.-M., Chuyong, G. B., Clay, K., Condit, R., Dattaraja, H. S., Davies, S. J., Duque, A., Esufali, S., Ewango, C., Fernando, R. H. S., Fletcher, C. D., Gunatilleke, I. A. U. N., Hao, Z., Harms, K. E., Hart, T. B., Hérault, B., Howe, R. W., Hubbell, S. P., Johnson, D. J., Kenfack, D., Larson, A. J., Lin, L., Lin, Y., Lutz, J. A., Makana, J.-R., Malhi, Y., Marthews, T. R., McEwan, R. W., McMahon, S. M., McShea, W. J., Muscarella, R., Nathalang, A., Noor, N. S. M., Nytch, C. J., Oliveira, A. A., Phillips, R. P., Pongpattananurak, N., Punchi-Manage, R., Salim, R., Schurman, J., Sukumar, R., Suresh, H. S., Suwanvecho, U., Thomas, D. W., Thompson, J., Uríarte, M., Valencia, R., Vicentini, A., Wolf, A. T., Yap, S., Yuan, Z., Zartman, C. E., Zimmerman, J. K., and Chave, J.: Local spatial structure of forest biomass and its consequences for remote sensing of carbon stocks, Biogeosciences, 11, 6827–6840,, 2014. 

Réjou-Méchain, M., Tymen, B., Blanc, L., Fauset, S., Feldpausch, T. R., Monteagudo, A., Phillips, O. L., Richard, H., and Chave, J.: Using repeated small-footprint LiDAR acquisitions to infer spatial and temporal variations of a high-biomass Neotropical forest, Remote Sens. Environ., 169, 93–101,, 2015. 

Réjou-Méchain, M., Tanguy, A., Piponiot, C., Chave, J., and Hérault, B.: BIOMASS: An R Package for estimating aboveground biomass and its uncertainty in tropical forests, Methods Ecol. Evol., 8, 1163–1167,, 2017. 

Roussel, J.-R., Caspersen, J., Béland, M., Thomas, S., and Achim, A.: Removing bias from LiDAR-based estimates of canopy height: Accounting for the effects of pulse density and footprint size, Remote Sens. Environ., 198, 1–16,, 2017. 

Ruiz, L., Hermosilla, T., Mauro, F., and Godino, M.: Analysis of the influence of plot size and LiDAR density on forest structure attribute estimates, Forests, 5, 936–951,, 2014. 

Rutishauser, E., Noor'an, F., Laumonier, Y., Halperin, J., Rufi'ie, Hergoualc'h, K., and Verchot, L.: Generic allometric models including height best estimate forest biomass and carbon stocks in Indonesia, Forest Ecol. Manag., 307, 219–225,, 2013. 

Saatchi, S. S., Harris, N. L., Brown, S., Lefsky, M., Mitchard, E. T. A., Salas, W., Zutta, B. R., Buermann, W., Lewis, S. L., Hagen, S., Petrova, S., White, L., Silman, M., and Morel, A.: Benchmark map of forest carbon stocks in tropical regions across three continents, P. Natl. Acad. Sci. USA, 108, 9899–9904,, 2011. 

Särndal, C.-E., Swensson, B., and Wretman, J. H.: Model assisted survey sampling, Springer-Verlag, New York, USA, 1992. 

Singh, M., Evans, D., Coomes, D. A., Friess, D. A., Suy Tan, B., and Samean Nin, C.: Incorporating canopy cover for airborne-derived assessments of forest biomass in the tropical forests of Cambodia, PLoS One, 11, e0154307,, 2016. 

Slik, J. W. F., Aiba, S. I., Brearley, F. Q., Cannon, C. H., Forshed, O., Kitayama, K., Nagamasu, H., Nilus, R., Payne, J., Paoli, G., Poulsen, A. D., Raes, N., Sheil, D., Sidiyasa, K., Suzuki, E., and van Valkenburg, J. L. C. H.: Environmental correlates of tree biomass, basal area, wood specific gravity and stem density gradients in Borneo's tropical forests, Glob. Ecol. Biogeogr., 19, 50–60,, 2010. 

Slik, J. W. F. F., Bernard, C. S., Breman, F. C., Van Beek, M., Salim, A., and Sheil, D.: Wood density as a conservation tool: quantification of disturbance and identification of conservation-priority areas in tropical forests, Conserv. Biol., 22, 1299–1308,, 2008. 

Spriggs, R.: Robust methods for estimating forest stand characteristics across landscapes using airborne LiDAR, University of Cambridge, Cambridge, UK, 2015. 

Struebig, M. J., Turner, A., Giles, E., Lasmana, F., Tollington, S., Bernard, H., and Bell, D.: Quantifying the biodiversity value of repeatedly logged rainforests: gradient and comparative approaches from Borneo, Adv. Ecol. Res., 48, 183–224,, 2013. 

Sullivan, M. J. P., Talbot, J., Lewis, S. L., Phillips, O. L., Qie, L., Begne, S. K., Chave, J., Cuni-Sanchez, A., Hubau, W., Lopez-Gonzalez, G., Miles, L., Monteagudo-Mendoza, A., Sonké, B., Sunderland, T., ter Steege, H., White, L. J. T., Affum-Baffoe, K., Aiba, S., de Almeida, E. C., de Oliveira, E. A., Alvarez-Loayza, P., Dávila, E. Á., Andrade, A., Aragão, L. E. O. C., Ashton, P., Aymard C., G. A., Baker, T. R., Balinga, M., Banin, L. F., Baraloto, C., Bastin, J.-F., Berry, N., Bogaert, J., Bonal, D., Bongers, F., Brienen, R., Camargo, J. L. C., Cerón, C., Moscoso, V. C., Chezeaux, E., Clark, C. J., Pacheco, Á. C., Comiskey, J. A., Valverde, F. C., Coronado, E. N. H., Dargie, G., Davies, S. J., De Canniere, C., Djuikouo K., M. N., Doucet, J.-L., Erwin, T. L., Espejo, J. S., Ewango, C. E. N., Fauset, S., Feldpausch, T. R., Herrera, R., Gilpin, M., Gloor, E., Hall, J. S., Harris, D. J., Hart, T. B., Kartawinata, K., Kho, L. K., Kitayama, K., Laurance, S. G. W., Laurance, W. F., Leal, M. E., Lovejoy, T., Lovett, J. C., Lukasu, F. M., Makana, J.-R., Malhi, Y., Maracahipes, L., Marimon, B. S., Junior, B. H. M., Marshall, A. R., Morandi, P. S., Mukendi, J. T., Mukinzi, J., Nilus, R., Vargas, P. N., Camacho, N. C. P., Pardo, G., Peña-Claros, M., Pétronelli, P., Pickavance, G. C., Poulsen, A. D., Poulsen, J. R., Primack, R. B., Priyadi, H., Quesada, C. A., Reitsma, J., Réjou-Méchain, M., Restrepo, Z., Rutishauser, E., Salim, K. A., Salomão, R. P., Samsoedin, I., Sheil, D., Sierra, R., Silveira, M., Slik, J. W. F., Steel, L., Taedoumg, H., Tan, S., Terborgh, J. W., Thomas, S. C., Toledo, M., Umunay, P. M., Gamarra, L. V., Vieira, I. C. G., Vos, V. A., Wang, O., Willcock, S., and Zemagho, L.: Diversity and carbon storage across the tropical forest biome, Sci. Rep., 7, 39102,, 2017. 

Vincent, G., Sabatier, D., and Rutishauser, E.: Revisiting a universal airborne light detection and ranging approach for tropical forest carbon mapping: scaling-up from tree to stand to landscape, Oecologia, 175, 439–443,, 2014. 

Vira, B., Christoph, W., and Mansourian, S.: Forests, Trees and Landscapes for Food Security and Nutrition – A Global Assessment Report, IUFRO World Series Volume 33, Vienna, 2015. 

Walsh, R. P. D. and Newbery, D. M.: The ecoclimatology of Danum, Sabah, in the context of the world's rainforest regions, with particular reference to dry periods and their impact, Philos. T. R. Soc. B, 354, 1869–1883,, 1999.  

Watt, M. S., Adams, T., Gonzalez Aracil, S., Marshall, H., and Watt, P.: The influence of LiDAR pulse density and plot size on the accuracy of New Zealand plantation stand volume equations, New Zeal. J. For. Sci., 43, 15,, 2013. 

Wulder, M. A., White, J. C., Nelson, R. F., Næsset, E., Ørka, H. O., Coops, N. C., Hilker, T., Bater, C. W., and Gobakken, T.: Lidar sampling for large-area forest characterization: A review, Remote Sens. Environ., 121, 196–209,, 2012. 

Yanai, R. D., Battles, J. J., Richardson, A. D., Blodgett, C. A., Wood, D. M., and Rastetter, E. B.: Estimating uncertainty in ecosystem budget calculations, Ecosystems, 13, 239–248,, 2010. 

Zanne, A. E., Lopez-Gonzalez, G., Coomes, D. A., Ilic, J., Jansen, S., Lewis, S. L., Miller, R. B., Swenson, N. G., Wiemann, M. C., and Chave, J.: Global wood density database, Dryad Digit. Repos.,, 2009. 

Zolkos, S. G., Goetz, S. J., and Dubayah, R.: A meta-analysis of terrestrial aboveground biomass estimation using lidar remote sensing, Remote Sens. Environ., 128, 289–298,, 2013. 

Short summary
Efforts to protect tropical forests hinge on recognizing the ecosystem services they provide, including their ability to store carbon. Airborne laser scanning (ALS) captures information on the 3-D structure of forests, allowing carbon stocks to be mapped. By combining ALS with data from 173 field plots on the island of Borneo, we develop a simple yet general model for estimating forest carbon stocks from the air. Our model underpins ongoing efforts to restore Borneo's unique tropical forests.
Final-revised paper