- Open Access
A sample design for globally consistent biomass estimation using lidar data from the Geoscience Laser Altimeter System (GLAS)
© Healey et al.; licensee BioMed Central Ltd. 2012
- Received: 27 April 2012
- Accepted: 3 September 2012
- Published: 31 October 2012
Lidar height data collected by the Geosciences Laser Altimeter System (GLAS) from 2002 to 2008 has the potential to form the basis of a globally consistent sample-based inventory of forest biomass. GLAS lidar return data were collected globally in spatially discrete full waveform “shots,” which have been shown to be strongly correlated with aboveground forest biomass. Relationships observed at spatially coincident field plots may be used to model biomass at all GLAS shots, and well-established methods of model-based inference may then be used to estimate biomass and variance for specific spatial domains. However, the spatial pattern of GLAS acquisition is neither random across the surface of the earth nor is it identifiable with any particular systematic design. Undefined sample properties therefore hinder the use of GLAS in global forest sampling.
We propose a method of identifying a subset of the GLAS data which can justifiably be treated as a simple random sample in model-based biomass estimation. The relatively uniform spatial distribution and locally arbitrary positioning of the resulting sample is similar to the design used by the US national forest inventory (NFI). We demonstrated model-based estimation using a sample of GLAS data in the US state of California, where our estimate of biomass (211 Mg/hectare) was within the 1.4% standard error of the design-based estimate supplied by the US NFI. The standard error of the GLAS-based estimate was significantly higher than the NFI estimate, although the cost of the GLAS estimate (excluding costs for the satellite itself) was almost nothing, compared to at least US$ 10.5 million for the NFI estimate.
Global application of model-based estimation using GLAS, while demanding significant consolidation of training data, would improve inter-comparability of international biomass estimates by imposing consistent methods and a globally coherent sample frame. The methods presented here constitute a globally extensible approach for generating a simple random sample from the global GLAS dataset, enabling its use in forest inventory activities.
- Forest monitoring
- Remote sensing
Methods are needed to monitor the magnitude and spatial distribution of global forest carbon storage, an important component of the global carbon cycle. Initiatives such as REDD (United Nations Collaborative Programmed on Reducing Emissions from Deforestation and Degradation in Developing Countries) depend upon accurate, precise, and consistent national-level reporting of forest carbon storage. Traditionally, estimates of carbon storage in the context of international monitoring have come from field-based inventories . In such inventories, well-developed principles of sample design support straightforward derivation of estimates and uncertainties. However, many countries do not have national forest inventories, and among those that do, important differences in methods and definitions can exist.
Satellite-based forest monitoring may offer observations which are more consistent across space and time, and potential contributions of remotely sensed estimation of carbon stored in biomass are widely recognized [2, 3]. However, barriers to broad acceptance of remotely sensed biomass estimates exist. Widely available satellite data, particularly from optical sensors such as Landsat and MODIS, may be relatively insensitive to different levels of biomass under closed forest canopies (e.g. [4, 5]). More importantly, while credible efforts have been made to empirically propagate errors through the process of summing pixel-level biomass predictions at the national level (e.g. ), acceptance of such approaches lags behind more formal estimation methods.
Several efforts to move beyond these limitations have centered around the use of lidar (light detection and ranging), not in wall-to-wall mapping (which can be relatively expensive) but as a vehicle for forest sampling . Lidar instruments measure characteristics of laser pulses as they return off of objects at different heights above the earth’s surface. The actively generated signals used by lidar typically penetrate deeper into the forest canopy than the passive signals used by optical sensors, and strong relationships are often found between lidar return data and forest structure parameters such as biomass and volume [8, 9].
While local- to regional-scale lidar monitoring missions are typically flown with instruments mounted on fixed-wing aircraft, globally consistent monitoring may be best achieved with spaceborne lidar. To date, the only widely available source of spaceborne lidar has been the GLAS (Geosciences Laser Altimeter System) instrument on NASA’s ICESat (Ice, Cloud, and land Elevation) satellite, which gathered data from 2002 to 2008. GLAS’ “full waveform” measurements are based upon time variation in the intensities of returned laser pulses, which resolve elliptical areas approximately 65 meters in diameter. GLAS measurements (“shots”) have been shown to be strongly correlated with biomass , and earlier problems with data quality on steeper slopes have been addressed to the point where such measurements can now be used in vegetation monitoring .
In this paper, we propose a method for identifying a subset of GLAS shots which can be treated as a simple random sample. We then demonstrate the use of such a sample over the U.S. state of California with a model-based estimator similar to that used by Stähl et al. . Model-based estimation, described below, allows us to predict, instead of measure, biomass at each sample point using relationships derived from a separate set of co-located ground and lidar measurements. Variance estimators used in this process take into account the uncertainty associated with the models used.
The sample design we describe is similar to that used by the U.S. national forest inventory (NFI), the Forest Inventory and Analysis program (FIA) administered by the U.S. Forest Service. Prior to a move toward a national sampling framework in the late 1990s, FIA plots were distributed and measured in slightly different ways in different regions of the country. The move to a nationally coherent sampling frame was accomplished by superimposing a hexagonal grid over the entire country, with the area of each grid cell equal to the nominal area represented by each FIA sample . In cells where one existing plot fell, that plot was kept. In those with more than one plot, only one was selected at random for retention. In those with no existing plots, a plot was established in a random location.
Establishment of this semi-systematic, equal-area sample frame, therefore, allowed FIA to accommodate existing measurement locations while drawing a sample which was spatially balanced across the country and yet was random with respect to forest conditions . The sample design we propose for GLAS follows a similar approach. One and only one GLAS shot is retained in each cell of an equal-area (but not equal-shape) tessellation of the area labeled as “forest” in a global land cover map. This tessellation is created following a fractal-based approach, using simple geometric rules to create equal-area clusters . Since retroactively “adding” GLAS measurements (the last of which were collected in 2008) is not possible, tessellation cell size (and, inversely, sample number) is limited by the constraint that each equal-area cell must contain at least one GLAS shot.
Given the elimination of all GLAS shots except one in every tessellation cell under this approach, it is of practical interest to know the precision of resulting biomass estimates. The precision (i.e. standard error) of model-based estimates of biomass in California using the GLAS sample will be compared to design-based estimates derived from FIA’s sample of more than 5500 field measurements in the state.
The R2 value the quadratic-only model was 0.87, although this figure should be viewed with the understanding that the R2 for a no-intercept model is calculated differently and represents a different aspect of correlation than R2 calculated for models which include an intercept term. To compare the two models, one can use a conditional R2, namely , where the full model is the intercept model and reduced model is the no-intercept model. The value of this R2 is 0.001 in our case. Thus, for the sake of comparison with relationships observed by others between GLAS and aboveground tree biomass, the fit of the model was quite similar to that of a quadratic-plus-intercept model, for which the r2 value was 0.64.
It should be noted that seven values (less than 4%) of the S1 sample exceeded the largest value in the model-building S2 dataset shown in Figure 4 (specifically, these were values of: 45, 46, 48, 50, 52, 54, and 60 meters). Ideally, the model-building dataset should span the entire range of the values to be modeled. However, given the small percentage of Lorey’s heights in S1 not represented in S2, we assume that the model is valid for the entire population. We likewise assume no spatial autocorrelation among S1 samples.
Model-based estimation using the sample design we describe provides a transparent method for estimating biomass for particular spatial domains. This sample design, in which one arbitrarily located sample point is drawn from equal-area sample units distributed across the landscape, is similar to that used by FIA, and our estimate of biomass density in the state of California closely matched FIA’s design-based estimate. The standard error of our estimate (approximately 9.8% of the estimate) was substantially larger than that of the FIA estimate (1.4%) and that derived through model-based estimation by Andersen et al.  using specifically acquired airborne lidar data (2011; 8%). However, the cost of the FIA estimate was approximately $10.5 million (using a commonly used valuation of $2000 per plot), and the lidar acquisition alone in Anderson et al.’s much smaller study area cost $60,000. While NASA’s investment in the GLAS mission was considerable, future use of GLAS data in the process described here represents an almost no-cost option for providing consistent, moderate-precision biomass estimates across the globe.
A primary advantage of the model-based inference used here is the capacity to apply models developed in areas of rich inventory data to GLAS shots informing estimates in ecologically similar areas where field data are sparse. For example, Nelson et al.  used relationships observed in a limited area of co-located biomass/GLAS observations to estimate biomass for the entire Canadian province of Quebec, following a modified model-based approach. However, the validity of inference in model-based approaches depends upon how well the stipulated models accord with the population of interest . The question of how well the model “applies” to the population of interest is a critical consideration in the application of our approach, whether the model was developed in situ or from a spatially remote but perhaps ecologically similar area. Since our model was created from an arbitrary subset of FIA’s presumably unbiased ground sample, there is a compelling argument that the model is appropriate for the forests of California.
The degree to which this model may apply beyond California remains an open question. Saatchi et al.  noted regional differences in the relationship between biomass and GLAS Lorey’s height in their pan-tropical study. Data collected to support biomass estimation using the global GLAS dataset would at least have to span major ecological systems. The need for broadly collected ground measurements in the composition of our S2 sample highlights the fact that there will always be demand for up-to-date ground data. Model-based estimation may spatially extend the value of available field data, but models must ultimately be grounded in actual observations that are relevant to the domain of interest.
The consolidation of ground data needed to support a global GLAS-based biomass inventory would require significant international cooperation and, as illustrated by our results, would likely not improve the precision of biomass estimates available in countries with established National Forest Inventories (NFIs). NFIs typically rely upon a denser sample than is available from GLAS, and do not have to account for model variance, which in our example made up approximately 44% of the total variance.
However, a GLAS-based biomass inventory would represent an internationally coherent basis for comparison among countries, especially those without established NFIs. Even moderate-precision biomass estimates would represent an improvement in many countries , and consistent sample design and estimation methods would remove an important source of uncertainty in international monitoring. The ICESat-2 mission, due in 2016, may provide an opportunity to update any GLAS-based biomass monitoring system. Although the scanning sensor on the ICESat-2 platform will provide continuous sample lines instead of discrete waveform returns, similar acquisition patterns from airborne lidar instruments have been discretized and used in model-based estimation approaches [16, 21].
An important variable not considered in this paper is how the area of forest is determined. As stated earlier, the domain of our estimation was the area in California mapped as “forest” by the MOD12Q1 global land cover product. However, significant disagreements can exist among land cover maps , due both to varying definitions and alternative mapping methods. Use of different maps may result in different S1 sample sizes, varying biomass density estimates, and different overall carbon estimates as density values are multiplied by mapped forest totals. Bearing in mind that the forest cover map used in this methodology functions as a proxy for the true distribution of forest, it is important to choose a map which best serves analytical needs. For international inventory purposes, it is reasonable to use a globally consistent product such as MOD12Q1.
In view of international efforts to increase or preserve forest carbon storage, the global GLAS height dataset presents an opportunity to establish how forest biomass was distributed internationally in 2005 (the mid-point of the GLAS mission). GLAS data were acquired in spatial patterns difficult to associate with either a systematic or random process. The sample design presented in this paper allows identification of a subset of GLAS data which may be used as a simple random sample to estimate biomass, perhaps globally, with consistent measures of uncertainty under a model-based estimation framework.
● The methods presented here constitute a globally extensible approach for generating a simple random sample from the global GLAS dataset. The properties of the sample collected by GLAS have hitherto not been strictly identifiable with any particular design.
● Model-based estimation, following Stähl et al. (2011), based upon GLAS data in the state of California produced an estimate of biomass density (biomass/hectare) almost identical to the estimate derived from the design-based NFI.
● Global application of model-based estimation using GLAS, while demanding significant consolidation of training data, would improve inter-comparability of international biomass estimates by imposing consistent methods and a globally coherent sample frame.
GLAS shots acquired in the following collections were intersected with the global MOD12Q1(v004) MODIS land cover product, subset for the state of California: L3B, L3C, L3D, L3E, L3F, L3G, L3H, and L3I. Shots were kept if they fell over one of five forest classes (“evergreen needle leaf”, “evergreen broadleaf”, “deciduous needle lead”, “deciduous broadleaf”, and “mixed”). This area became the domain over which average biomass density (tones/hectare) was to be estimated.
Shots were filtered only on the basis of quality flags due in many cases to clouds or other atmospheric anomalies. Topographic correction was applied following Lefsky et al. . Full-waveform signatures were processed to a crown-weighted height metric called “Lorey’s height” . Lorey’s height, used recently in a global tropical biomass mapping project , was the GLAS derivative upon which subsequent modeling and estimation were based.
A sample design for GLAS data
The choice of a particular statistical estimator does not necessarily imply any particular sample design . The model-based approach to inference that we describe in the next section has been employed with airborne lidar data, often using sample designs which consider strips of lidar measurements as systematic cluster samples (e.g. [16, 21]). However, as illustrated above, the irregular positioning of GLAS ground tracks poses difficulty in defining the terms under which the sample can be considered representative of the population. The primary contribution of this paper, which we describe in this section, is a means of identifying a subset of GLAS data which can be treated as a simple random sample in the estimation process.
Assign an ordinal number to each pixel in the forest map representing the domain of interest. The MODIS product referenced above was re-sampled from its native 1-kilometer resolution to 230 meters so that processing would occur at a scale closer to the field of view of the GLAS shots (approximately 70 meters). Re-sampling to 230 meters produced over 1.6 million pixels in the California study area, which, given subsequent operations, was near local computing limits. Next, a space-filling curve  was applied through the center point of each “forest” pixel to generate an ordered list of pixel locations. This fractal-based ordering process (described in detail in ) involved the generation of a self-similar line (Piano curve) that folded in upon itself as it occupied the set of pixel centers found on the landscape.
Align GLAS-based Lorey’s heights with spatially correspondent pixels on the ordinal number line. The GIS coverage of GLAS shots was spatially intersected with the ordered network of forested pixel centers in a combination of a GIS and Microsoft Access processes, and the Lorey’s heights were added to approximately 102,000 of the 1.6 million locations represented on the number line. In cases where multiple GLAS points fell within a single pixel, one was chosen at random as representative.
Divide the ordered number line into equal-length segments, such that there is at least one Lorey’s height measurement associated with each segment. A script was written using the open-source R statistical programming language , in which the ordered list of forested pixel centers (i.e., the number line) was iteratively broken into equal segments of varying length and tested for the condition of containing at least one pixel associated with a GLAS shot. This was accomplished by transforming the line into a matrix of n columns made up of equal contiguous line segments of lengths l, with the total length of the number line equal to n l plus a remainder, which was ignored. A matrix was considered a viable solution if each column contained at least one pixel center point with an associated GLAS measurement (see step 2).
Matrices representing different segment lengths were tested, starting with the shortest possible segment satisfying the requirement of ≥1 GLAS shot per segment (i.e. one half the length of the longest gap between GLAS shots on the number line) and working upward until a viable solution was found. Since the location of the first pixel represented in the number line was arbitrary, all possible segmentation starting points were tested for every tested segmentation length, “looping” the end of the number line to the beginning. R code for each of these operations has been uploaded to the Journal archive ( Additional File 1: “gap_finding_and_segment_sampling_R_code.pdf”).
For segments associated with more than one GLAS shot, choose one at random for the sample.
Similar to FIA’s sample design, this process assures a relatively uniform spatial distribution of plots but allows locally random positioning of measurements. Following FIA’s precedent , this sample is treated subsequently as a simple random sample.
Model-based inference depends upon fundamentally different assumptions than the design-based methods used by most field-based inventories, including FIA’s (for detailed description of the difference between model- and design-based inference, see ). Unlike design-based estimation, model-based methods treat observations as realizations of a random process (model).
The model-based approach we follow is similar to that of Stähl et al. . We make use of two samples; sample S1 is the “application sample” developed in the steps above, for which modeled Lorey’s heights are the only data available; and, sample S2, which is composed of co-located field and GLAS measurements which can be used to build and assess biomass models to be applied at all S1 plots. In this study, the S2 sample was not a subsample of S1; S2 was made up of the 35 single-condition, forested FIA ground plots in California which had plot centers falling within 120 meters of the center of a GLAS shot and which did not fall along condition boundaries, as determined by visual inspection of high-resolution National Land Cover Database maps . Care was taken to avoid condition boundaries to minimize mismatch in the forest measured by the satellite and forest measured in the field.
The predicted value of biomass Ŷ was constructed using maximum likelihood estimates based on the S 2 sample, denoted by the parameter estimates . The parameter estimates were constructed using linear model package, lm, in the R programming language . By standard theory for linear models , Ŷ is an unbiased estimator of the expectation of Y and there is an unbiased estimator of the variance covariance matrix of the parameter estimates, and an unbiased estimate of V(Ŷ), the variance of Ŷ.
where Ŷ i is the predicted biomass value for the ith element of S 1 and n is the number of elements in S 1. The estimator is an unbiased estimator of the population mean of the expected value of the Ŷ i with respect to the sampling distribution; that is: . The bias in the estimate is , where e i is the value of unknown error, with the expected value of this bias equal to zero.
where is the standard estimate with respect to simple random sampling of the variance of the mean. Given the form of the model used in this case , the double sum collapses to because there is only one term in the sum. R code for the model-building and estimation processes is given in Additional File 2 (1-modelBuilding_and_biomass_estimation_R_code.pdf).
This work was supported by the NASA Carbon Monitoring System and by FIA. The authors additionally would like to thank James Men love for expert help in interpreting FIA biomass data as well as Karen Waddell and Jock Blackard for their help with the FIA database. Gretchen Moisen also provided valuable discussion and feedback.
- Brown S: Estimating biomass and biomass change of tropical forests. A primer. In Book Estimating biomass and biomass change of tropical forests. A primer. 134th edition. Rome, Italy: Food and Agriculture Organization of the United Nations (FAO); 1997.Google Scholar
- Gibbs HK, Brown S, Niles JO, Foley JA: Monitoring and estimating tropical forest carbon stocks: making REDD a reality. Environ Res Lett 2007, 2: 13.Google Scholar
- Goetz SJ, Baccini A, Laporte NT, Johns T, Walker W, Kellndorfer J, Houghton RA, Sun M: Mapping and monitoring carbon stocks with satellite observations: a comparison of methods. Carbon Balance Manag 2009, 4: 2. 10.1186/1750-0680-4-2View ArticleGoogle Scholar
- Houghton RA, Butman D, Bunn AG, Krankina ON, Schlesinger P, Stone TA: Mapping Russian forest biomass with data from satellites and forest inventories. Environ Res Lett 2007, 2: 7.View ArticleGoogle Scholar
- Powell SL, Cohen WB, Healey SP, Kennedy RE, Moisen GG, Pierce KB, Ohmann JL: Quantification of live aboveground forest biomass dynamics with Landsat time-series and field inventory data: A comparison of empirical modeling approaches. Remote Sens Environ 2010, 114: 1053–1068. 10.1016/j.rse.2009.12.018View ArticleGoogle Scholar
- Saatchi SS, Harris NL, Brown S, Lefsky M, Mitchard ETA, Salas W, Zutta BR, Buermann W, Lewis SL, Hagen S, et al.: Benchmark map of forest carbon stocks in tropical regions across three continents. Proc Natl Acad Sci U S A 2011, 108: 9899–9904. 10.1073/pnas.1019576108View ArticleGoogle Scholar
- Wulder MA, Masek JG, Cohen WB, Loveland TR, Woodcock CE: Opening the archive: How free data has enabled the science and monitoring promise of Landsat. Remote Sens Environ 2012, 122: 2–10.View ArticleGoogle Scholar
- Dubayah RO, Drake JB: Lidar remote sensing for forestry. J For 2000, 98: 44–46.Google Scholar
- Miller ME, Lefsky M, Pang Y: Optimization of Geoscience Laser Altimeter System waveform metrics to support vegetation measurements. Remote Sens Environ 2011, 115: 298–305. 10.1016/j.rse.2010.09.002View ArticleGoogle Scholar
- Baccini A, Laporte NT, Goetz SJ, Sun M, Dong H: A first map of tropical Africa's above-ground biomass derived from satellite imagery. Environ Res Lett 2008, 3: 9.View ArticleGoogle Scholar
- Lefsky MA, Keller M, Pang Y, de Camargo PB, Hunter MO: Revised method for forest canopy height estimation from Geoscience Laser Altimeter System waveforms. J Appl Remote Sens 2007, 1: 18.Google Scholar
- Schutz BE, Zwally HJ, Shuman CA, Hancock D, DiMarzio JP: Overview of the ICESat mission. Geophys Res Lett 2005,32(L21S01):1–4.Google Scholar
- Wulder MA, White JC, Nelson RF, Naesset E, Ørka HO, Coops NC, Hilker T, Bater CW, Gobakken T: Lidar sampling for large-area forest characterization: A review. Remote Sens Environ 2012, 121: 196–209.View ArticleGoogle Scholar
- Nelson R, Boudreau J, Gregoire TG, Margolis H, Naesset E, Gobakken T, Ståhl G: Estimating Quebec provincial forest resources using ICESat/GLAS. Can J Forest Resour 2009, 39: 862–881. 10.1139/X09-002View ArticleGoogle Scholar
- Homer C, Dewitz J, Fry J, Coan M, Hossain N, Larson C, Herold N, McKerrow A, VanDriel JN, Wickham J: Completion of the 2001 National Land Cover Database for the Conterminous United States. Photogramm Eng Remote Sens 2007, 73: 337–341.Google Scholar
- Ståhl G, Holm S, Gregoire TG, Gobakken T, Naesset E, Nelson R: Model-based inference for biomass estimation in a LiDAR sample survey in Hedmark County, Norway. Can J Forest Resour 2011, 41: 96–107. 10.1139/X10-161View ArticleGoogle Scholar
- Reams GA, Smith WD, Hansen MH, Bechtold WA, Roesch FA, Moisen GG: The Forest Inventory and Analysis Sampling Frame. In Book The Forest Inventory and Analysis Sampling Frame. volth edition. USA: Forest Service Southern Research Station; 2005. General Technical Report SRS-80 General Technical Report SRS-80Google Scholar
- Bechtold WA, Patterson PL: The Enhanced Forest Inventory and Analysis Program - National Sampling Design and Estimation Procedures. Asheville, NC: USDA Forest Service, Southern Research Station; 2005.Google Scholar
- Lister A, Scott C: Use of space-filling curves to select sample locations in natural resource monitoring studies. Environ Monit Assess 2009, 149: 71–80. 10.1007/s10661-008-0184-yView ArticleGoogle Scholar
- Miles PD: Forest Inventory EVALIDator web-application version 4.01beta. Book Forest Inventory EVALIDator web-application version 4.01beta volth edition. 2011. Available only on internet: . City: USDA Forest Service Northern Research Station http://fiatools.fs.fed.us/Evalidator4/tmattribute.jsp Available only on internet: . City: USDA Forest Service Northern Research StationGoogle Scholar
- Andersen H-E, Strunk J, Temesgen H: Using airborne light detection and ranging as a sampling tool for estimating forest biomass resources in the Upper Tanana Valley of Interior Alaska. West J Appl For 2011, 26: 157–164.Google Scholar
- Gregoire TG: Design-based and model-based inference in survey sampling: appreciating the difference. Can J Forest Resour 1998, 28: 1429–1447. 10.1139/x98-166View ArticleGoogle Scholar
- Fritz S, See L: Identifying and quantifying uncertainty and spatial disagreement in the comparison of Global Land Cover for different applications. Glob Chang Biol 2008, 14: 1057–1075. 10.1111/j.1365-2486.2007.01519.xView ArticleGoogle Scholar
- Pang Y, Lefsky M, Andersen H-E, Miller ME, Sherrill K: Validation of the ICEsat vegetation product using crown-area-weighted mean height derived using crown delineation with discrete return lidar data. Can J Remote Sens 2008, 34: S471-S484. 10.5589/m08-074View ArticleGoogle Scholar
- R Development Core Team: R: A Language and Environment for Statistical Computing. Vienna: R Foundation for Statistical Computing; 2011.Google Scholar
- Neter J, Kutner MH, Nachtsheim CJ, Wasserman W: Applied Linear Statistical Models. 4th edition. Chicago: Irwin; 1996.Google Scholar
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.