Scientia Forestalis, volume 42, n. 104
p.591-604, dezembro de 2014

Mapping aboveground carbon stocks using LiDAR data in Eucalyptus spp. plantations in the state of São Paulo, Brazil

Mapeamento de estoques de carbono acima do solo utilizando dados LiDAR em  plantações de Eucalyptus spp no estado de São Paulo, Brasil

Carlos Alberto Silva1
Carine Klauberg2
Samuel de Pádua Chaves e Carvalho3
Andrew T. Hudak4
Luiz Carlos Estraviz Rodriguez5

1MSc Forest Resources. Departament of Forest Sciences, ESALQ/USP - Av. Pádua Dias,11 – 13418-900 - Piracicaba, SP, Brasil. E-mail: carlos_engflorestal@outlook.com.
2Post-doctoral Research Associate. USDA - US Forest Service - Rocky Mountain Research Station, RMRS, 1221 South Main Street, Moscow, Idaho, USA – 83843. E-mail: carine_klauberg@hotmail.com.
3Associated Professor. Departament of Forest Sciences, UFMT - Universidade Federal de Mato Grosso - 78060-900 - Cuiabá-MT. E-mail: sam.padua@gmail.com
4Research Forester. ISDA - US Forest Service - Rocky Mountain Research Station, RMRS, 1221 South Main Street, Moscow, Idaho, USA – 83843. E-mail: ahudak@fs.fed.us
5Associated Professor. Departament of Forest Sciences, ESALQ/USP - Av. Pádua Dias,11 – 13418-900- Piracicaba, SP. E-mail: lcer@usp.br

Recebido em 22/01/2014 - Aceito para publicação em 31/07/2014

Resumo

As plantações florestais de rápido crescimento fornecem um significativo baixo custo de sequestro de carbono para redução de gases de efeito estufa. O objetivo deste estudo foi avaliar a utilização de dados LiDAR (Light Detection And Ranging) aerotransportado para estimativa do estoque de carbono acima do solo (AGC) em plantações de Eucalyptus spp. Parâmetros biométricos tais como altura (Ht) e diâmetro à altura do peito (DAP) das árvores foram coletadas em parcelas de inventários convencionais. Os modelos de regressão linear múltipla, com o intuito de estimar o estoque de carbono total acima do solo (AGCt), em toras comerciais (AGCc) e em resíduos de colheita (AGCr), foram desenvolvidos utilizando métricas derivadas da nuvem de pontos LiDAR. Os melhores modelos de um conjunto de seis modelos foram selecionados com base no Critério de informação de Akaike corrigido (AICc) e avaliados segundo a Raiz Quadrada do Erro Médio (RMSE) e Coeficiente de determinação ajustado (R²-adj). Os três melhores modelos para as estimativas do estoque de AGC foram AGCt: R²-adj= 0,81; RMSE = 7,70 Mg.ha-1; AGCc: R²-adj= 0,83; RMSE = 5,26 Mg.ha-1; AGCr: R²-adj = 0,71; RMSE = 2,67 Mg.ha-1. Este estudo mostrou que métricas derivadas a nuvem de pontos LiDAR podem ser usadas para estimar o estoque de AGC em plantações Eucalyptus spp. no Brasil com precisão. Concluímos que existe um bom potencial para monitorar o crescimento e a fixação de carbono em plantações de Eucalyptus spp. com o uso da tecnologia LiDAR.
Palavras-chave: Perfilhamento a Laser aerotrasportado-ALS; métricas LiDAR, estoque de carbono; plantações de rápido crescimento.

Abstract

Fast growing plantation forests provide a low-cost means to sequester carbon for greenhouse gas abatement. The aim of this study was to evaluate airborne LiDAR (Light Detection And Ranging) to predict aboveground carbon (AGC) stocks in Eucalyptus spp. plantations. Biometric parameters (tree height (Ht) and diameter at breast height (DBH)) were collected from conventional forest inventory sample plots. Regression models predicting total aboveground carbon (AGCt), aboveground carbon in commercial logs (AGCc), and aboveground carbon in harvest residuals (AGCr) from LiDAR-derived canopy structure metrics were developed and evaluated for predictive power and parsimony. The best models from a family of six models were selected based on corrected Akaike Information Criterion (AICc) and assessed by the root mean square error (RMSE) and coefficient of determination (R²-adj). The best three models to estimate AGC stocks were AGCt: R²-adj = 0.81, RMSE = 7.70 Mg.ha-1; AGCc: R²-adj = 0.83, RMSE = 5.26 Mg.ha-1; AGCr: R²-adj = 0.71, RMSE = 2.67 Mg.ha-1. This study showed that LiDAR canopy structure metrics can be used to predict AGC stocks in Eucalyptus spp. plantations in Brazil with high accuracy. We conclude that there is good potential to monitor growth and carbon sequestration in Eucalyptus spp. plantations using LiDAR.
Keywords: Airborne Laser Scanning ALS; LiDAR metrics, C stock; fast growing plantation.


INTRODUCTION

The international scientific community has issued warnings over steadily increasing concentrations of greenhouse gases in the atmosphere as a driver of climate change. According to the Intergovernmental Panel on Climate Change (IPCC), the atmospheric concentrations of the carbon dioxide (CO2) have increased to levels unprecedented in at least the last 800,000 years; concentrations of CO2 have increased by 40% since pre-industrial times, primarily from fossil fuel emissions and secondarily from net land use change emissions (IPCC, 2013). The forest sector contributes 17.4 percent of all greenhouse gases anthropogenic sources; most of this due to deforestation and forest degradation (IPCC, 2007). Therefore, quantifying the substantial roles of forests as C stores, as sources of C emissions and as C sinks has become key to understanding the global C cycle (FAO, 2010).

Eucalyptus plantations in Brazil cover an area of 5.10 million hectares, accounting for 76.6% of the country's total reforested areas (ABRAF, 2013). These plantations, when well established, produce between 100 and 400 Mg.ha-1 at the peak of maximum mean annual increment (ECOAR, 2003). Thus, timber management as well as quantifying total C in these forests has received special attention for the possibility of production regulation and the issuance of C credits for the role they can play as C sequestration reservoirs (YU, 2004).

Measurements of aboveground C (AGC) stocks in Brazilian forest plantations are currently limited by budgets and time, making it impractical to efficiently implement an extensive inventory (ZONETE, 2009). Remote sensing techniques, using both active and passive optical sensors, have been presented as a viable alternative for the estimation of forest stand parameters in planted and natural forests (ASNER et al., 2011; HUDAK et al., 2006, 2012; NÆSSET, 2002, 2004; NÆSSET; GOBAKKEN, 2008; SAATCHI et al., 2011).

Among the techniques of remote sensing currently available, Light Detection and Ranging (LiDAR) scanning assessments have emerged as the most prominent in the forestry sector (EVANS et al., 2009; NÆSSET, 2002, 2004; NÆSSET; GOBAKKEN, 2008). LiDAR is an active remote sensing technology that uses emitted laser pulses to measure distances between objects (JENSEN, 2007). It has proven useful for predictions of forest structure attributes at larger spatial scales than can be measured practically in the field (HUDAK et al., 2009). LiDAR also provides a good alternative for measuring C stocks in forests, since it may enable precise data collection at high spatial resolution and in a relatively short time compared to conventional methods (EVANS et al., 2009, HUDAK et al., 2012).   

In Brazil, the use of airborne LiDAR technology is relatively new; on the other hand, a few studies have already been developed and published such as Rodriguez et al. (2010) and Zonete et al. (2010), that used LiDAR data to predict diameter at breast height (DBH), height (Ht), and basal area (BA) of Eucalyptus spp. plantation in Bahia. Furthermore, Macedo (2009) and Zandoná (2006) used LiDAR data to model individual tree and stand volume in Pinus spp. and Eucalyptus spp. plantations, respectively. Therefore, these works have demonstrated the great potential of airborne technology to predict relevant biometric parameters in forest plantations. However, more studies are needed to improve the accuracy of predictions for operational implementation by the Brazilian forestry sector.

The aim of this study was to evaluate the use of LiDAR to predict aboveground C (AGC) stocks in Eucalyptus spp plantations located in São Paulo state, Brazil. The specific objectives were to: (i) select the best LiDAR metrics to build AGC models; (ii) select the best models to predict AGC; (iii) estimate AGC stock in the total biomass (trunk and canopy) (AGCt), commercial logs (AGCc) and in the harvest residuals (AGCr) (bark, leaves and branches); and (iv) apply the predictive models selected across the landscape to map the spatial distribution of C stocks of Eucalyptus spp. at the stand level for the benefit of forest managers.


METHOD


Study area

The study was based on data collected in a set of temporary and permanent sample plots installed for the purpose of annual forest inventory in industrial Eucalyptus plantations managed by Fibria Celulose S/A, a pulp company located in  São Paulo state, Brazil. The sample plots were randomly distributed among eight farms described in Table 1. These farms are located in the Paraíba Paulista Valley, in the Eastern part of the state of São Paulo, near the municipalities of Jacarei, São Luiz do Paraitinga and Paraibuna (Latitude between 23°30' and 23°00' S, Longitude between 46°00' and 45°00' W) (Figure 1). The climate of the region is characterized as humid temperate, with dry winters and hot summers (Cwa) (KÖPPEN; GEIGER, 1928). The annual mean precipitation is around 1200 – 1232 mm and the annual mean temperature ranges from 17.1º C in the coldest month (July) to 23.9° C in the hottest month (February). The topography in the selected plantations is complex, ranging from mildly to very hilly and from 577.52 m to 1310.23 m in elevation. The soils of the region are predominantly Oxisols, and depending on the topographic variations may also include Inceptisols and Ultisols; however, all are clayey.


Figure 1. Location of study areas in Sao Paulo state, Brazil.
Figura 1. Localização da área de estudo no estado de São Paulo, Brasil.

The sampled trees in the sample plots are hybrid clones of two Eucalyptus species, Eucalyptus grandis and Eucalyptus urophylla, with ages varying from two to eight years and planted predominantly in a 3m x 2m grid, resulting in a density of 1,667 trees.ha-1.

Table 1. Main characteristics of the farms comprising the study - in São Paulo state.
Tabela 1. Principais características das fazendas em estudo - Estado de São Paulo-SP.
Eucalyptus plantation farms City Area (ha) Age
F987 Jacarei 39.53 2.3
F986 Jacarei 94.16 3.3
F849 São Luiz of Paraitinga 138.96 4.7
F950 Paraibuna 86.72 5.5
F184 São Luiz of Paraitinga 58.34 5.9
F166 São Luiz of Paraitinga 84.35 6.1
F948 Paraibuna 79.33 6.8
F634 Paraibuna 84.80 8.0


Field data collection

For the inventory of C stock in the field, 136 circular plots (400 m²) were installed on eight farms.  All plots were geo-referenced with a geodetic GPS with differential correction capability (Trimble Pro-XR). The projected coordinate system used was UTM SIRGAS 2000, zone 23 S. All trees were measured for diameters at breast height (DBH) and a subsample (15%) of trees for tree heights. For trees in the plot that were not directly measured for Ht, the inventory team of Fibria S/A company provided predictions calculated from hypsometric equations. The mean values of DBH and Ht measured in plots in the farms included in the study are shown in Table 2.

Table 2. Descriptive values of the biometric parameters of the network of plots inventoried in the study area (Fall 2012).
Tabela 2. Descrição dos parâmetros biométricos da rede de parcelas inventariadas na área de estudo (Outono 2012).
Eucalyptus plantation farms DBH (cm) Ht (m) Number of plots
s s
F987 12.73 0.76 18.60 1.45 14
F986 14.59 0.95 24.17 1.38 21
F849 15.26 1.20 25.24 2.35 26
F950 8.82 1.57 11.65 2.55 17
F184 14.14 0.60 22.16 1.54 14
F166 14.57 0.70 23.74 1.19 17
F948 13.91 0.63 22.69 1.45 12
F634 13.70 1.20 23.29 1.68 15
DBH= diameter breast height; Ht= tree hieght; = average; s = standard deviation.

The AGCt, AGCc and AGCr (kg. tree-1) were obtained through allometric models, employing as independent variables the logarithm of DBH and Ht, and as dependent variables the logarithms of the AGCt, AGCc and AGCr (SILVA, 2013). The AGC models had adjusted coefficients of determination (R²-adj) of 0.97 and relative root mean squared errors (RMSE) ranged from 1.44 to 4.57 Mg.ha-1 (Table 3).

Table 3. Ordinary least square regression models predicting AGC stocks.
Tabela 3. Modelos de regressão de mínimos quadrados prevendo estoques de AGC.
Predictive models of AGC in tree R²-adj RMSE (kg.tree-1) RMSE (%)
ln(AGCt)     = -2.87 + 1.95*ln(DBH) + 0.44*ln(Ht) 0.97 4.57 12.38
ln(AGCc)    = -3.89 + 1.72*ln(DBH) + 0.83*ln(Ht) 0.97 2.73 11.04
ln(AGCr)     = -2.61 + 2.49*ln(DBH) + 0.47*ln(Ht) 0.97 1.44 11.24
R²-adj = adjusted coefficient of determination; RMSE= root-mean-square error.

The summed C content of all trees in the sample plot was multiplied by the plot area (0.04 ha) to calculate the C stored in the sample plot in Mg.ha-1. The summary statistics of AGC stocks (Mg.ha-1) measured in the farms evaluated are presented in Table 4.

Table 4. Predicted AGC stocks (Mg.ha-1) in the farms studied.
Tabela 4. Estoque de AGC predito (Mg.ha-1) nas fazendas estudadas.
Eucalyptus plantation farms AGCt AGCc AGCr Age Number of Plots
s s s
F987 14.31 5.29 8.32 3.31 6.44 1.94 2.3 14
F986 45.45 6.19 28.86 4.36 16.67 1.88 3.3 21
F849 62.71 6.69 41.72 4.93 20.69 1.94 4.7 26
F950 63.67 9.96 42.11 6.66 21.27 3.32 5.5 17
F184 64.59 6.81 43.81 4.81 20.37 2.03 5.9 14
F166 69.76 9.77 47.67 6.89 21.65 2.91 6.1 17
F948 60.42 13.12 41.19 9.01 18.86 4.24 6.8 12
F634 72.54 9.06 50.00 6.75 22.11 2.66 8.0 15
AGCt= Aboveground carbon total; AGCc= Aboveground carbon in commercial logs; AGCr= Aboveground carbon in harvesting residuals ; = average; s = standard deviation.


LiDAR data processing

LiDAR data were acquired using a Riegl ALS (model LMS-Q680I) mounted on a Piper Seneca II aircraft. LiDAR flight parameters are shown in Table 5.

Table 5. LiDAR flight acquisition parameters (Fall 2012).
Tabela 5. Parâmetros de aquisição dos dadosLiDAR (Outono de 2012).
Parameter Value
Density of the laser pulse 10 m²
Dynamic range 12bits
Speed 57 m/s (205.20 km.h-1)
Flying altitude average 422.94 m
Sweep angle 45º
Horizontal Precision 0.1-0.15m (1.0 sigma)
Frequency sweep 400kHz
IMU / GPS Applainix 510
IMU=Inertial measurement unit; GPS=Global Positioning System.

FUSION software (v. 3.3) (McGaughey, 2013) was used to process the LiDAR data and generate four main products: the digital terrain model (DTM), the digital surface model (DSM), the canopy height model (CHM), and the LiDAR metrics used in this study.

Initially the catalog function was used to evaluate the quality of the LiDAR data set. Ground returns were separated from the vegetation returns using groundfilter followed by gridsurfacecreate to generate the 1-meter-resolution digital terrain model (DTM), which represents the bare earth elevation. The digital surface model (DSM) which represents the earth's surface and includes the trees and other objects on it, was created using canopymodel also at 1 m resolution. The clipdata function was applied to normalize heights and to make sure that the z coordinate for each point corresponded to the height above ground and not the elevation of the point. After heights were normalized, the canopy height model (CHM), which represents the height of the forest, was created at 1 m resolution using canopymodel. The polyclipdata function was used to select only the LiDAR points falling within the field sample plots of interest. The cloudmetrics function was used to generate the LiDAR metrics within the sample plots, while the gridmetrics function was used to generate the same LiDAR metrics in across the LiDAR surveys at 5 m resolution.

LiDAR first returns provide a more stable distribution of object heights than the other returns (Bater et al., 2011). For instance, Kim et al. (2009) used the first returns to distinguish between live and dead standing tree biomass, while Ørka et al. (2009) and Korpela et al. (2010) used the first returns to classify individual tree species by LiDAR intensity metrics. Therefore, we used just the LiDAR first returns to generate the LiDAR metrics in this study. The twenty LiDAR metrics used in this paper are presented in table 6.

Table 6. LiDAR-derived vegetation height and canopy cover metrics.
Tabela 6. Métricas de altura e cobertura de copa derivadas do LiDAR.
Category LiDAR metric Acronym
Height Maximum Height hmax
Mean height hmean
Standard deviation of mean height hsd
Coefficient of variation of height hcv
Mode of height hmod
5th percentile of height H5
10th percentile of height h10
20th percentile of height h20
25th percentile of height h25
30th percentile of height h30
40th percentile of height h40
50th percentile of height h50
60th percentile of height h60
70th percentile of height h70
75th percentile of height h75
80th percentile of height h80
90th percentile of height h90
95th percentile of height h95
99th percentile of height h99
Cover Percentage of first returns above 2 m Cdens


Regression modeling

Initially, Pearson's correlation (r) was used to identify highly correlated predictor variables (r>0.9); redundant predictors were subsequently excluded. Secondly, we determined the best subsets using the regsubsets function in the "leaps" package in R (R Development Core Team, 2010), which selects the highest correlated subset of predictive variables through exhaustive search. To determine the best subset the “regsubsets” function uses the Mallows Cp statistic, which compares the error sum of squares for a reduced model to the mean square error of the full model (MALLOWS, 1973).

Third, after selecting the best eight variables, we used the “lm” linear model function in R to define the prospective models. Finally, the corrected Akaike information criterion (AIC) (AKAIKE, 1973; 1974) was calculated to measure the relative quality of each model and to rank them accordingly by minimum AIC. However, since we had small sized samples (n/p < 40, where n is number of samples and p is number of parameters of the model), we used the corrected information criterion (AICc) to rank the models and select the best of them (HURVICH; TSAI, 1989).

Residuals from all prospective models were analyzed graphically and tested for normality using the Shapiro-Wilk test (SHAPIRO; WILK, 1965) and for heteroscedasticity using the Breusch-Pagan test (BREUSCH; PAGAN, 1979). Statistics used to evaluate the alternative regression models were: (i) R²-adj and (ii) RMSE, both absolute and relative (%). The performance of the models was determined by means of a leave-one-out cross-validation strategy. For this paper, we defined an acceptable accuracy as a relative RMSE below 15%.


RESULTS AND DISCUSSION

Pearson's correlation (r) applied to the 21 candidate LiDAR metrics determined that 15 were highly correlated (r>0.9). The six remaining metrics were not so highly correlated, and therefore were used to build prospective models by best subset analysis. Even though the number of metrics was reduced to six, they still represented the top height (e.g. h99) and canopy density (Cdens) metrics of the LiDAR data. The six LiDAR metrics selected were hcv, hstd, h5, h10, h99 and Cdens. Some of them were positively correlated, such as hcv and hstd, while others were negatively correlated, such as Cdens and hcv (Figure 2).  


Figure 2. Pearson's correlation among selected LiDAR metrics.
Figura 2. Correlação de Person para as métricas LiDAR.

The best subsets of LiDAR metrics resulted in six models for predicting C content with the number of independent variables increasing from 1 to 6. Figure 3 shows box-and-whisker plots of the six selected LiDAR metrics.


Figure 3. Box and whisker plots of the six selected LiDAR metrics.
Figura 3. Diagrama de caixa das métricas LiDAR selecionadas.

Table 7 shows the competing AGC models and diagnostic statistics for predicting the three AGC response variables. For AGCt, the six best subset model R²-adj ranged from 0.79 to 0.81. The most common metric was h99 (present in all six models), followed by h10 (present in five models) and hcv (present in three models). These variables explained 79-81% of the variability in AGCt. In addition, the AICc ranged from 966.24 to 953.64, absolute RMSE ranged from 7.70 to 8.25 Mg.ha-1, and relative RMSE ranged from 13.43 to 14.40%.  

The AGCc model R²-adj ranged from 0.82 with only one independent variable (h99) to 0.83 using six independent variables (h5, h10, h99, hcv, hstd and Cdens). The variables that appeared most in the composition of these models were the same as were selected in the AGCt models. The AICc ranged from 848.06 to (three predictors) 856.51 (one predictor). Absolute RMSE ranged from 5.26 to 5.51 (Mg.ha-1) and relative RMSE from 13.75 to 14.40%.  

For the AGCr, R²-adj ranged from 0.67 for the model with only one independent variable (h99) to 0.73 for models with six explanatory variables h5, h10, h99, hcv, hstd and Cdens). The independent variables that appeared in most models were h99, h10, and Hcv. The AICc ranged from 664.03 (one predictor) to 678.46 (six predictors). Absolute RMSE ranged from 2.67 to 2.86 (Mg.ha-1); relative RMSE ranged from 14.22 to 15.23%.

Table 7. AGC models created using LiDAR-derived vegetation height and canopy cover metrics.
Tabela 7. Modelos de AGC criado a partir das métricas de altura e cobertura de copas derivadas do LiDAR.
AGC # Models R²-adj AICc RMSE (Mg.ha-1) RMSE (%)
AGCt 1 = -46.76 + 3.98h99 0.79 966.24 8.25 14.40
2 = -47.20 + 0.66h10 + 3.51h99 0.80 958.98 7.97 13.91
3 =-55.73 + 1.51h10 +2.81h99 + 66.59h10 0.81 955.04 7.80 13.60
4 =-42.92 - 0.21h5 + 1.56h10 +1.93h99 + 4.66hstd 0.81 953.64 7.70 13.43
5 = -72.69-0.16h5 + 1.61h10 + 2.83h99 + 65.89hcv + 0.18Cdens 0.81 959.43 7.78 13.57
6 = -43.90 + 0.20h5 + 1.47h10 +1.49h99 -72.06hcv + 8.18hstd+ 0.15Cdens 0.81 955.59 7.65 13.75
AGCc 1 = -34.49 + 2.90h99 0.82 856.51 5.51 14.40
2 = -37.76 + 0.41h99 + 2.60h99 0.82 850.45 5.35 13.97
3 =-42.72 + 0.91h10 + 2.20h99+ 38.68hcv 0.83 848.06 5.26 13.75
4 =-42.44 – 36.66h5 - 0.92h10+ 2.21h99 + 36.66h10 0.83 850.21 5.26 13.74
5 = -45.51 + 0.04h5 + 0.93h10 + 1.64h99 + 2.93hstd + 0.11Cdens 0.83 850.91 5.20 13.59
6 = -38.31 + 0.14h5 + 0.87h10 +1.48h99 -31.97hcv + 4.36hstd + 0.09Cdens 0.83 848.61 5.20 13.58
AGCr 1 = -7.85 + 1.02h99 0.67 678.46 2.86 15.23
2 = -8.00 + 0.24h10 + 0.84h99 0.69 670.40 2.76 14.67
3 =-11.49 + 0.59h10 + 0.56h99 + 27.18hcv 0.71 664.03 2.67 14.22
4 =-10.54 – 0.11h5 + 0.62h10 + 0.60h99 + 20.25hstd 0.71 665.54 2.67 14.19
5 = -13.49 + 0.04h5 + 0.65h10 + 0.18h99+ 2.01hcv + 0.08Cdens 0.72 670.85 2.60 14.54
6 = -3.87 + 0.05h5 + 0.57h10 - 0.04h99 -42.67hcv + 3.91hstd + 0.05Cdens 0.73 664.58 2.57 13.69
R²-adj = adjusted coefficient of determination; AICc= corrected Akaike information criterion; RMSE= root mean square error.

All the aforementioned AGC model residuals (Table 7) exhibited both normality and homogeneity of variances when evaluated by Shapiro-Wilk and Breusch-Pagan tests. The three best models with the lowest AICc statistic identified for predicting AGCt, AGCc and AGCr had by definition the lowest AICc statistics and were composed of either three or four independent variables. (Table 7).

The LiDAR metrics such as height percentiles and canopy density cover included in these models are frequently observed in studies that use LiDAR data for forest inventory (HUDAK et al., 2006, 2012; NÆSSET, 2002, 2004; NÆSSET; GOBAKKEN, 2008). The best three models for predicting AGCt, AGCc and AGCr explained 81, 83 and 71% of the variation, respectively. In addition, all three AGC models selected had a low RMSE between the predicted and observed values, indicating good model accuracy. The model accuracies were also favorably compared with other studies predicting forest biomass C (GARCÍA et al., 2010; HUDAK et al., 2012; PATENAUDE et al., 2004).

The LiDAR canopy height and density metrics selected in this work are also commonly used to predict AGC stocks in forests. For instance, Hudak et al. (2012) used LiDAR mean height and other metrics to quantify aboveground forest C pools and ?uxes from repeated LiDAR surveys in Moscow Mountain in Northern Idaho, USA, and Stephens et al. (2012) used LiDAR height percentiles and canopy density for estimation of C stocks in New Zealand planted forests. LiDAR height and canopy density metrics, and in some cases intensity metrics, also have been used to predict stand biomass (e.g. GARCIA et al., 2010; NAESSET, 2002). LIDAR metrics can also be useful to predict stand height, volume and basal area (RODRIGUEZ et al., 2010; ZONETE et al., 2010) in Brazilian Eucalyptus plantations.   

In this paper, equivalence plots of observed versus predicted AGC via leave-one-out cross-validation indicated that predicted and observed C were similar (Figure 4A1, B1 and C1). The models accuracy is confirmed by the low relative value of RMSE between observed and predicted AGC by cross-validation. The results from the equivalence tests support the stability of the chosen models. In addition, the patterns of residuals also confirm that the models do not violate the assumption of homogeneity of variances (Figure 4. A2, B2 and C2).


Figure 4. Equivalence plots of leave-one-out cross-validation of best models for: A) predicted versus observed AGCt (A1) and residuals (A2); B) predicted versus observed AGCc (B1) and residuals (B2); and C) predicted versus observed AGCr (C1) and residuals (C2). The black line indicates the line of best fit of a simple linear model of observations regressed on predictions. The gray shaded bar defines the region of similarity in the intercept. If the error bar about the line of best fit falls within the shaded gray region, then the intercept of the linear model does not significantly differ from its range of expected values. If the black error bar falls between the dotted lines, then the slope of the linear model does not significantly differ from its range of expected values.
Figura 4. Gráfico de equivalência para a validação cruzada “deixa um de fora” dos melhores modelos: A) Valores de AGCt preditos versus observados (A1) e resíduos (A2); B) Valores de AGCc preditos versus observados(B1) e resíduos (B2); C) Valores de AGCr predito versus observados (C1) e resíduos (C2). A barra cinza sombreada define a região de semelhança para o intercepto. Se a barra de erro preta estiver dentro da região sombreada cinza, o intercepto do modelo linear não difere significativamente do seu intervalo de valores esperados. Se a barra de erro preta estiver entre as linhas pontilhadas, em seguida, o termo o coeficiente angular do modelo linear não difere significativamente do seu intervalo de valores esperados.

AGC was slightly over-predicted for the F987, F986, F166, F634 farms (Figure 5). These even-aged plantations do not show height or other structural variability, so a relatively small number of field plots equitably distributed across the range of age and height variability was sufficient to develop robust models with the LiDAR metrics.   


Figure 5. Observed versus predicted AGC by the selected models for each farm in the study.
Figura 5. AGC predito versus observado de acordo com os modelos selecionados para cada fazenda em estudo.

AGC maps are useful to evaluate the forest homogeneity, because they show the distribution of the C content across of the whole area. The three best AGC models were applied to predict AGC across the landscape, and we created 24 maps showing the C stored at the stand level (Figures 6 and 7).


Figure 6. Predicted AGC stocks (AGCt, AGCc and AGCr) for Eucalyptus plantations located in F987, F986, F849 and F950 farms. Map resolution is 5 m.
Figura 6. Estoque de AGC predito (AGCt, AGCceAGCr) para as plantações de Eucalyptus spp. localizadas nas fazendas F987, F986, F849 e F950. O mapa tem resolução de 5 m.


Figure 7. Predicted AGC stocks (AGCt, AGCc and AGCr) for Eucalyptus plantations located in F184, F166, F948 and F634 farms. Map resolution is 5 m.
Figura 7. Estoque de AGC predito (AGCt, AGCc e AGCr) para as plantações de Eucalyptus spp. localizadas nas fazendas F184, F166, F948 e F634. O mapa tem resolução de 5 m.


CONCLUSIONS

Best subsets regression models predicting AGCt, AGCc and AGCr from LiDAR metrics produced R2-adj of 0.81, 0.83 and 0.71, and relative RMSE of 7.70, 5.26 and 2.67 Mg.ha-1, respectively. Therefore, we conclude that LiDAR metrics can be used to predict aboveground C stocks in fast growing Eucalyptus plantations in Brazil with acceptable accuracy. High spatial resolution maps of AGC stocks and the component pools can be created from LiDAR data once robust models are developed. Our results point to the need to continue advancing and promoting the application of LiDAR data for inventory of forest plantations of Eucalyptus spp.


ACKNOWLEDGEMENTS

We thank the US Forest Service (USFS) Rocky Mountain Research Station at the Forestry Sciences Laboratory in Moscow, Idaho and USFS International Programs for their support. This study was also supported by FAPESP (Process nº 2012/03176-0 MSc assistantship). LiDAR data collections were funded by FIBRIA, a Brazilian pulp and paper company.


REFERENCES

ABRAF - Brazilian Association of Forest Plantation Producers. Statistical Yearbook 2013 years: based ABRAF 2012. 7ed. Brasília: ABRAF, 2013. 142 p.

AKAIKE, H. Information theory and an extension of the maximum likelihood principle. In: PETROV, B.N.; CSAKE, F. Proceedings of the 2nd International Symposium on Information Theory. Budapest: Akademiai Kiado, 1973. p. 267–281.

AKAIKE, H. A new look at the statistical model identification. IEEE Transactions on Automatic Control, Piscataway, v. 19, n. 6, p. 716–723, 1974.

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.; HUGHES, R. F. High-resolution forest carbon stocks and emissions in the Amazon. . Frontiers in Ecology and the Environment, New York, v. 9, n. 8, p. 434-439, 2011. 

BATER, W. C.; WULDER, M. A.; COOPS, N. C.; NELSON, R. F.; HILKER, T.; NASSET, E. Stability of Sample-Based Scanning-LiDAR-Derived Vegetation Metrics for Forest Monitoring. Transactions on Geosciences and Remote Sensing, Piscataway, v. 49, n. 6, p. 2385–2392, 2011.

BREUSCH, T. S.; PAGAN, A. R. A simple test for heteroscedasticity and random coefficient variation. Econometrica, Philadelphia, v. 47, n. 5, p. 1287–1294, 1979.

ECOAR – INSTITUTO ECOAR PARA CIDADANIA. Efeito estufa. São Paulo: Ecoar, 2003. 5 p.

EVANS J. S.; HUDAK A. T.; FAUX R.; SMITH A. M. Discrete Return LiDAR in Natural Resources: Recommendations for Project Planning, Data Processing, and Deliverables. Remote Sensing, Basel, v.1, n.4, p.776-794, 2009.

GARCÍA, M.; RIAÑO, D.; CHUVIECO, E.; DANSON, F. M. Estimating biomass carbon stocks for a Mediterranean forest in central Spain using LiDAR height and intensity data. Remote Sensing of Environment, New York, v. 114, n. 4, p. 816-830, 2010.

IPCC - INTERGOVERNMENTAL PANEL ON CLIMATE CHANGE. Climate change 2007: Geneva, World Meteorological Organization, 2007. 520 p.

IPCC - INTERGOVERNMENTAL PANEL ON CLIMATE CHANGE. Working Group I Contribution to the IPCC Fifth Assessment Report Climate Change 2013: The Physical Science Basis Summary for Policymakers. Geneva, World Meteorological Organization. 2013. 380 p.

HUDAK, A. T.; STRAND, E. K.; VIERLING, L. A.; BYRNE, J. C.; EITEL, J. U. H.; MARTINUZZI, S.; FALKOWSKI, M. J. Quantifying aboveground forest carbon pools and fluxes from repeat LiDAR surveys. Remote Sensing of Environment, New York, v.82, p.397–416, 2012.

HUDAK, A. T.; EVANS, J. S.; SMITH, A. M. S. Review: LiDAR utility for natural resource managers. Remote Sensing, Basel, v. 1, n. 4, p. 934–951, 2009.

HUDAK, A. T.; CROOKSTON, N. L.; EVANS, J. S.; FALKOWSKI, M. J.; SMITH, A. M. S.; GESSLER, P. E. Regression modeling and mapping of coniferous forest basal area and tree density from discrete-return LiDAR and multispectral satellite data. Canadian Journal of Remote Sensing, Kanata, v. 32, p.126–138, 2006.

HURVICH, C. M.; TSAI, C. L. Regression and time series model selection in small samples. Biometrika, New York, v. 76, n. 2, p. 297–307, 1989.

JENSEN, J. R. Remote Sensing of the Environment: An Earth Resource Perspective, 2 ed. Upper Saddle River: Prentice-Hall, 592 p. 2007.

KIM, Y.; YANG, Z.; COHEN, W.B.; LAUVER, C.L.; VANKAT, J.L. Distinguishing between live and dead standing tree biomass on the North Rim of Grand Canyon National Park, USA using small-footprint lidar data. Remote Sensing of Environment, New York, v.113, n.11, p.2499-2510, 2009.

KÖPPEN, W.; GEIGER, R. Klimate der Erde. Gotha: Verlag Justus Perthes, 1928. Wall-map 150cmx200cm.

KORPELA, I.; ØRKA, H.O.; MALTAMO, M.; TOKOLA, T.; HYYPPA, J. Tree species classification using airborne LiDAR- effects of stand and tree parameters, downsizing of training set, intensity normalization, and sensor type. Silva Fennica, Helsinki, v.44, n.2, p.319-339, 2010.

MACEDO, R. D. C. Estimativa volumétrica de povoamento clonal de Eucalyptus spp. através de laser scanner aerotransportado. 2009. 143 p. Dissertação (Mestrado em Sensoriamento Remoto) - Instituto Nacional de Pesquisas Espaciais, São José dos Campos, 2009.

MALLOWS, C. L. Some comments on Cp. Technometrics, Milwaukee, v. 15, n. 4, p. 661–675, 1973.

MCGAUGHEY, R. J. FUSION/LDV: software for LiDAR Data Analysis and Visualization [Computer program]. Washington: USDA, Forest Service Pacific Northwest Research Station, 2013. 150 p. Available at: < http://http://forsys.cfr.washington.edu/fusion/ FUSION_manual.pdf >. Access: May. 20, 2013.

NÆSSET, E. Estimation of above- and below-ground biomass in boreal forest ecosystems.International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences, Freiburg, v. 36, p. 145-148, 2004

NÆSSET, E. Predicting forest stand characteristics with airborne scanning laser using a practical two-stage procedure and field data. Remote Sensing of Environment, New York, v.80, p. 88-99, 2002.

NÆSSET, E.; GOBAKKEN, T. Estimation of above and below-ground biomass across regions of the boreal forest zone using airborne laser. Remote Sensing of Environment, Basel, v.112, n. 6, p. 3079-3090, 2008.

ØRKA, H. O.; NÆSSET, E.; BOLLANDSA’S, O. M. Classifying species of individual trees by intensity and structure features derived from airborne laser scanner data. Remote Sensing of Environment, Basel, v. 113, n. 6, p.1163-1174, 2009.

PATENAUDE, G. L.; HILL, R. A.; MILNE, R.; GAVEAU, D. L. A.; BRIGGS, B. B. J.; DAWSON, T. P. Quantifying forest above ground carbon content using LiDAR remote sensing. Remote Sensing of Environment, Basel, v. 93, n. 3, p. 368-380, 2004.

R DEVELOPMENT CORE TEAM. R: A language and environment for statistical computing. Vienna: R Foundation for Statistical Computing, 2013. Available at: < http://www.R-project.org/ >. Access: Sep. 20, 2013. 

RODRIGUEZ, L. C. E.; POLIZEL, J. L.; FERRAZ, S. F. B.; ZONETE, M. F.; FERREIRA, M. Z. Inventário florestal com tecnologia laser aerotransportada de plantios de Eucalyptus spp. no Brasil. Ambiência, Guarapuava, v. 6, p. 67-80, 2010.

SAATCHI, S.; HARRIS, N. L.; BROWN, S.; LEFSKY, M.; MITCHARD, E. T.; SALAS, W.; ZUTTA, B. R.; BUERMANN, W.; LEWIS, S. L.; HAGEN, S.; PETROVA, S.; SHAPIRO, S. S.; WILK, M. B. An analysis of variance test for normality (complete samples). Biometrika, London, n. 52, v. 3/4, p. 591-611, 1965.

SILVA, C. A. Carbono na parte aérea de plantios de Eucalyptus spp. – em nível de árvore por amostragem destrutiva e para talhões inteiros após o ajuste de métricas LiDAR. 2013. 154 p. Dissertação (Mestrado em Recursos Florestais) - Escola Superior de Agricultura “Luiz de Queiroz”, Universidade de São Paulo, Piracicaba, 2013.

STEPHENS, P. R.; KIMBERLEY, M. O.; BEETS, P. N.; PAUL, T. S. H.; SEARLES, N.; BELL, A.; BRACK, C.; BROADLEY, J. Airborne scanning LiDAR in a double sampling forest carbon inventory. Remote Sensing of Environment, New York, v. 117, p. 348-357, 2012.

YU, C. M. Sequestro florestal do carbono no Brasil: dimensões políticas socioeconômicas e ecológicas. São Paulo: Annablume; IEB, 2004. 280 p.

ZANDONÁ, D. F. Potencial do uso de dados laser scanner aerotransportado para estimativa de variáveis dendrométricas. 2006. 92 p. Dissertação (Mestrado em Engenharia Florestal) - Universidade Federal do Paraná, Curitiba, 2006.

ZONETE, M. F.; RODRIGUEZ, L. C. E.; PACKALÉN, P. Estimação de parâmetros biométricos de plantios clonais de eucalipto no sul da Bahia: uma aplicação da tecnologia laser aerotransportada. Scientia Florestalis, Piracicaba, v. 38, n. 86, p. 225-235, 2010.

ZONETE, M. F. Análise do uso da tecnologia laser aerotransportado para inventários florestais em plantios clonais de Eucalyptus spp. no sul da Bahia. 2009. 95 p. Dissertação (Mestrado em Recursos Florestais) - Escola Superior de Agricultura “Luiz de Queiroz”, Universidade de São Paulo, Piracicaba, 2009.