Aim of study
: To test the use of LiDAR data from a single acquisition in order to estimate volume over_bark variations in a 5-yr period of
Pinus radiata
D. Don.
Area of study
: Province of Bizkaia in the Autonomous Community of the Basque Country (Spain).
Material and methods
: Two field plot measurements were made in 2011 and 2015 and two wood volume models (one for each year) were fitted using the metric variables of the 2012 LiDAR points cloud. The models were applied to a 26.59 m raster covering the study area and the increase in volume at each pixel was calculated by subtraction.
Main results
: The increase in estimated wood volume, when added to the volume of timber extracted in the area during the 5-yr period under consideration, yielded an average increase of 13.74 m
3
ha
-1
yr
-1
, which corresponds to the average growth of the
P. radiata
in that area. The harvest area estimated using this procedure largely coincides with the actual harvest area in the same period. The value of
R2
(85%) of the wood volume model for 2011 is similar to that obtained in other studies. However, as expected, the one obtained for the wood volume model for 2015 (80%) is significantly lower.
Research highlights
: The increase in wood volume can be estimated using a single LiDAR flight and field data from the 5-yr period provided that data from plots subjected to this kind of harvest is included in the models.
Key words:LiDAR forest inventorywood volume increasewood volume LiDAR modelAbbreviations used:CCF (Canopy Cover Fraction)CV (Coefficient of Variation)Dm (Mean Diameter)
G (basal area, m
2
ha
-1
)
Hm (mean height, m)Ho (top height, m)IFN4 (4th National Forest Inventory)LiDAR (Light Detection and Ranging)MAE (Mean Absolute Error)ME (Mean Error)RMSE (Root Mean Square Error)
Vob (stem volume over bark, m
3
ha
-1
)
Authors´ contributions:
MAV conceived and designed the experiments and contributed to writing, data production, and analysis. EM wrote the introduction and reviewed the work. FR revised the statistical calculations and fitted models.
Gobierno Vasco (Projects SAI10/147- SPE10UN90); Vicerrectorado de Investigación de la Universidad del País Vasco-Euskal Herriko Unibertsitatea (UPV/EHU) (Project NUPV14/11).
Competing interests:
The authors have declared that no competing interests exist.
Introduction
Over recent years, there has been increasing interest in the evaluation of the forestry variables, such as biomass and timber volume, especially within the context of the Kioto protocol (
Lindner & Karjalainen, 2007
;
Zhang
et al.
, 2007
). The quantification of these forestry parameters reveals the structure, operation and dynamics of forestry ecosystems and their assessment as an energy resource. For this reason, it is essential to determine these parameters in forestry management (
Torres
et al.
, 2010
;
Pan
et al.
, 2011
). In most European countries, the increase in these variables is greater than the amount harvested on an annual basis, converting forests into major carbon and biodiversity reserves (
Mateos
et al.
, 2016
). One of the key indicators for assessing forestry conditions is the estimation of forestry parameters, such as timber volume, but the usual methods for estimating this parameter require prolonged fieldwork (
Song
et al.
, 2016
). A large part of the forestry bibliography relating to this subject focuses on adjustment of timber volume equations by means of regression for geographical regions and specific tree species (
Nilsson, 1996
;
Hyyppä
et al.
, 2001
;
Popescu
et al.
, 2003
;
Holmgren, 2004
;
Zianis
et al.
, 2005
;
Dalponte
et al.
, 2009
;
Latifi
et al.
, 2010
;
Magnussen
et al.
, 2010
). The high cost of sampling, especially when this is done in inaccessible areas, limit in many cases the establishment of sufficient plots to determine real variability (
Popescu, 2007
).
The link between the parameters obtained in forestry inventories and LiDAR (
Light Detection And Ranging
) data through regression models offers an attractive approach to estimating forestry variables on a regional, continental or global scale (
Li & Weiskittel, 2010
;
Næsset, 2011
;
Babcock
et al.
, 2015
). Over recent decades, LiDAR technology has become an extremely powerful tool for estimating the biophysical properties of the vegetation. This is an active, non-destructive remote sensing method that offers a better alternative for the study of the structural variables of the canopy due to the penetrating capacity of the laser pulse through gaps in the vegetation (
Leiterer
et al.
, 2015
;
Giannico
et al.
, 2016
;
Koenig & Höfle, 2016
). The recording of multiple returns allows us to obtain sufficient three-dimensional, morphological information at a broad range of space and time scales of forestry variables such as the height of trees (Dean
et al.
, 2006;
Popescu & Zhao, 2008
), the basal area and the diameter of the crown (
Giannico
et al.
, 2016
), volume and aboveground biomass (
Popescu, 2007
;
Jochem
et al.
, 2010
;
Zolkos
et al.
, 2013
), and the density and classification of groups of forest species (
Brandtberg, 2007
;
Zhang
et al.
, 2016
). Moreover, the use of LiDAR data offers the possibility of reducing the number of required sampling plots (
Maltamo
et al.
, 2007
;
Montaghi
et al.
, 2013
) and the accuracy obtained with LiDAR is greater than that reached with traditional inventory techniques (
Condes
et al.
, 2013
;
Watt
et al.
, 2013
;
Maltamo
et al.
, 2014
).
Now that LiDAR technology has demonstrated its efficiency and accuracy in estimating most forestry variables, today the challenge is to try to take advantage of this tool in order to predict the annual growth of the plant structure at a large scale. The ease with which data can be updated allows us to readjust the estimations and quickly assess changes or disturbances in forest stands by estimating their development over time, or establish comparisons between different areas with a simple update of the forestry cartography (
Riaño
et al.
, 2004
). The bibliography contains numerous studies that used bi-temporal data from two LiDAR flights to estimate plant growth (
Yu
et al.
, 2004
,
2006
;
Næsset & Gobakken, 2005
;
Hopkinson
et al.
, 2008
;
Vepakomma
et al.
, 2011
;
Fekety
et al.
, 2014
;
Cao
et al.
, 2016
;
Song
et al.
, 2016
). Some researchers have estimated the growth of trees in forests on the basis of the differences in height derived from LiDAR data in flights separated by a 2-yr interval (
Næsset & Gobakken, 2005
;
Yu
et al.
, 2006
;
Vepakomma
et al.
, 2011
;
Cao
et al.
, 2016
). These studies allowed the growth in different forestry parameters to be calculated, such as the basal area or volume (
Næsset & Gobakken, 2005
) as well as estimations of the vertical and lateral growth of trees (
Vepakomma
et al.
, 2011
). Other researchers have used multi-temporal LiDAR datasets to estimate the growth of forests (
Yu
et al
., 2006
;
Hopkinson
et al
., 2008
;
Fekety
et al.
, 2014
;
Cao
et al
., 2016
;
Song
et al
., 2016
). However, we have not found any previous studies that use the data from a single LiDAR flight to estimate the increase in timber volume over a period of several years by fitting models that relate the timber volume at different moments in time with the metrics of the heights of a single LiDAR points cloud obtained at the beginning of the period. In view of the high cost of the flight, the optimisation of the long-term use of LiDAR data gathered during a campaign may represent a considerable economic advantage.
Material and methods
This study assesses the application of data from the last LiDAR flight made over the Autonomous Community of the Basque Country in 2012 to evaluate the growth of volume over bark (Vob) during 5 years on stands of
Pinus radiata
D. Don (PR) through the comparison of estimated volumes by applying fitted models to relate separately the field data obtained with 5 years of difference with the metrics of the LiDAR points cloud at the beginning of the period.
The area based approach (ABA) for forest inventory applications using LiDAR has been utilised in this study, following the methodology applied by other researchers (
Hopkinson
et al.
, 2008
;
Jochem
et al.
, 2010
;
Yu
et al.
, 2010
). In order to get a better understanding of how LiDAR and the set of field data can be used to map Vob, we used two sets of field data gathered in 2011 and 2015 respectively. This involves analysing and validating the proposed Vob models, fitting the corresponding models for both sets of data in accordance with the metric variables of the LiDAR points cloud obtained in the 2012 LiDAR flight. The results will allow us to assess whether it is possible to estimate –with an acceptable margin of error– the growth in the timber volume based on field data obtained over a 5-yr period and the metrics of a single LiDAR flight made during this time could optimise the available resources and would represent a large decrease in the cost compared to traditional estimation methods.
Study area
The study area lies in the Province of Bizkaia in the Autonomous Community of the Basque Country (ACBC) in the north of Spain. Its center is the municipality of Muxika and data from the surrounding municipalities between the coordinates 2º 30’ – 2º 54’ west and 43º 06’ – 43º 24’ north are also used. The average height is 230 m. The mean annual rainfall is 1277 mm and the mean annual temperature is 13.8 ºC. Frosts are infrequent and no physiological drought is apparent. With regard to the silviculture applied in the area, the rotation cycle is, on average, 35 yrs, varying between 30 and 40 according to the site quality and the price of wood at the time when harvest takes place. Clear cutting are applied followed by repopulation with initial densities of 1300-1600 trees ha
-1
. Between the 7th and 9th year, non-commercial thinning takes place, eliminating 500 trees ha
-1
and pruning up to 2.5 m on the trees that have not been eliminated. Between the 12th and 15th year, the first commercial thinning is carried out, leaving between 500 and 700 trees ha
-1
. Between the 17th and 20th year, the second commercial thinning is applied, leaving between 400 and 500 trees ha
-1
. The third thinning is carried out between the 22nd and 25th year, leaving between 250 and 300 trees for the clear cutting. The total production is on average 550 m
3
ha
-1
in one rotation, and therefore the average growth is 14.85 m
3
ha
-1
yr
-1
. Volumes and growths of 3153 ha occupied by
P. radiata
in the municipality of Muxika were estimated. The area in which data was gathered and results were obtained can be seen in Fig. 1.
The area in which the study was carried out. The area in which field data were gathered is shown in green. The
area in which the study was conducted and the results were obtained is shown in red. The points correspond to the location
of the IFN4 plots selected for the study.
2011 Field data
Two sets of data were used according to the year in which these were gathered. In order to fit the wood volume model from 2011, data was used from the Fourth National Forestry Inventory (IFN4) of the Province of Bizkaia (
http://www.euskadi.eus/inventario-forestal-2011/web01-a3estbin/es/
). Initially, all the plots included in the IFN4 and located within the boundaries of the municipality of Muxika and neighboring municipalities were selected. From these plots, 111 located within the area occupied by
P. radiata
were selected, according to the 2010 Forestry Map of the Basque Country (
http://www.euskadi.eus/mapa-forestal-del-pais-vasco-ano-2010/web01-a2hestat/es/
). From these 111 plots, 13 were eliminated as they had no data or did not have
P. radiata
as their main species. From the remaining 98, those that showed a clear incompatibility between the field data and the vegetation visible in the aerial photography were eliminated. This occurred because the coordinates of the centers of the plots were taken with GPS navigators that did not allow differential correction and had a potential error of more than 12 m. For the same reason, those plots in which the difference between the dominant height of the field trees and the 95th percentile of the height of the points of its LiDAR points cut was greater than 12 m were eliminated. Furthermore, those plots with a basal area (G) greater than 70 m
2
ha
-1
were also eliminated as these were not considered possible for this species in this area. Following the application of these filters, only 8 additional plots were discarded. Therefore, 90 plots were used finally in this study. Their location can be seen in Fig. 1.
IFN4 plots are of variable radius, so that a tree will or will not be measured in accordance with its diameter and the distance to the center of the plot. For each selected plot, the number of trees per hectare (N, trees ha
-1
), the dominant height (Ho, m) defined as the mean height of the 100 thickest trees per hectare, the mean height (Hm, m), the mean diameter (Dm, cm), the basal area (G, m
2
ha
-1
) and the wood volume over bark (Vob, m
3
ha
-1
) were calculated. For the last of these calculations, the Excel add-in CUBIFOR (
Rodríguez
et al.
, 2008
) was used. With the aim of reducing data-gathering costs for subsequent years, a Vob tariff was adjusted in accordance with G and Ho. In this way, in order to calculate volumes in measured plots it is sufficient to measure both variables quickly and economically. The Vob estimated with this tariff will be considered as the field Vob for the entire study.
2015 Field data
In order to be able to assess the possibility of using LiDAR data from a single flight when estimating the forestry resources in subsequent years, 55 field plots were measured during the month of December, 2015. So as to prevent any lack of accuracy in locating the center of the plot from affecting the estimation of the final results, those areas with a sufficient level of uniformity were selected so that a displacement of up to 10 m in the center of the plot would have only a slight effect on the final result.
With the aim of locating uniform areas, two influential variables were selected, the first being Canopy Cover Fraction (CCF), as a percentage of the surface covered by the canopy of trees, and the second being the height of the vegetation. The data gathering area was divided with a 10 m-sided grid and in each 100 m
2
square and the CCF was estimated as a percentage of first returns, higher than 3 m over the total of first returns of the LiDAR points contained in it. For each grid, the typical deviation and the average of its CCF value and those of the 8 surrounding squares were calculated and those squares in which the coefficient of variation (CV) was not greater than 30% were selected if the CCF was lower than 70% and those in which the CV was greater than 10% if the CCF was greater than 70%. To estimate the height of the vegetation, the data gathering area was divided with a 10 m-sided grid and in each grid cell, the height of the 95th percentile of the heights of the LiDAR points contained in it was calculated. For each grid cell, the mean and typical deviation of its value and that of the 8 surrounding cells was calculated. Those grid cells in which the CV was not greater than 20% were selected if the 95th percentile was less than 8 m; those grid cells in which the CV was not greater than 10% were selected if the 95th percentile was between 8 and 12 m, and those in which the CV was not greater than 5% were selected if the 95th percentile was greater than 12 m. Only those grid cells that met both conditions were candidates for containing a plot center for data gathering. In accordance with these conditions, 55 points were found and located in the field with a GPS navigator. From each point, the dominant height was measured as the average of the thickest three trees within a radius of 9.77 m from the point that corresponds to a plot with a surface area of 300 m
2
. The basal area was also measured by relascope inventory from the center of the plot. The measurement of these two variables allows the field work to be done cheaply and quickly. Later, the Vob was calculated for each plot applying the Vob tariff in accordance with Ho and G calculated with the plots from the IFN4. The values and statistics of the variables used in the measurements of the plots from the IFN4 (2011) and the new plots (2015) are shown in Table 1.
Summary statistics of the stand variables of Pinus radiata used to develop a model to predict stem
volume in plantations in Muxika (northern Spain)
LiDAR Data
The LiDAR data used came from the flight over the Autonomous Community of the Basque Country in July, 2012, with a mean density of 1 pulse m
-2
and up to 4 returns per pulse.
For the entire data-gathering area, those
las
format files that contained a selected plot from the IFN4 or some of the plots measured in the field in 2015 and all the sheets contained partially or completely in the municipality of Muxika were selected. All the selected
las
format files were normalized by subtracting from the height Z of each point the corresponding height of the 1 m × 1 m Digital Terrain Model (DTM) of the grid cell published by the Basque Government (
ftp://ftp.geo.euskadi.net/lidar/
). Once the points clouds had been normalized, the corresponding cylinders were clipped to the projection of each one of the plots, the center of which was that of the IFN4 plots and radius 15 m. In this way, 90 sets of LiDAR points coinciding with the 90 IFN4 plots were obtained. For example, the cut of plot 559 can be seen in Fig. 2.
Normalized points cloud cut from plot 559 of
the IFN4. The color of the points depends on their height
above the ground.
In order to make this clip, the PolyClipData tool of the FUSION program (
http://forsys.sefs.uw.edu/fusion/fusionlatest.html
) was used. Subsequently, the statistical parameters that describe each of the points clouds corresponding to plots were calculated. The results of the LiDAR points cloud statistics were compared to the results of the field measurements and volumes for each plot. In this way a table was obtained in which each row contains the field data and Vob of each IFN4 plot and the statistics corresponding to its LiDAR points cut. The selected plots and their data were separated in two groups of 60 and 30 plots respectively. The first of these was used to fit the model that relates the Vob with the statistics of the LiDAR points cloud and the second was used to validate this model. In order to separate the two groups, a random number between 0 and 100 was generated for each plot. Those plots with assigned numbers of under 65 were used to fit the model; those plots with a number over 65 were kept in order to validate the model. The same classification was used for adjusting the Vob tariff in accordance with G and Ho. The location of the plots of the model and of the validation group can be seen in Fig. 3.
Location of the IFN4 plots used to fit the model
(green points) and those used in the validation (red points).
The same procedure was followed using the LiDAR data for 2011 and the plots measured in the field in 2015. Fifty-five cuts of LiDAR points were obtained with their statistics and 55 Vob data calculated for each plot measured in 2015. Forty plots were selected to fit the model and the validation group was formed by the 15 plots. It should be pointed out that none of the plots measured in 2015 was located in areas in which clear cutting had been performed during the 5-yr period 2011-2015. This meant that the fitted model did not consider the possibility of volume 0. Some plots were indeed located in areas that had been subjected to thinning during the 5-yr period.
Clear cutting and thinning data in Muxika during the period 2011-2015
Clear cutting and thinning data in the municipality of Muxika between 2011 and 2015 was gathered in order to be able to use this to estimate the growth in wood volume during the period 2011-2015. Data on harvest provided by the Provincial Council of Bizkaia (the administrative office that manages the forest stands of the Province of Bizkaia), can be seen in Table 2.
Clear cutting and thinnig data in Muxika between 2011 and 2015.
In order to adjust the Vob tariff per hectare in accordance with the stand variables G and Ho, 90 plots of the IFN4 were used. Of these, 60 were used to fit the model and 30 to validate this. Plots were assigned to the model group or the validation group according to the random number assigned to them. In order to build the tariff, the STATGRAPHICS (
http://www.statgraphics.com
) was used. Using the regression model selection option, Vob was used as a dependent variable and Ho, Ho
2
, G, G
2
and Ho*G were used as independent variables (Table 3). The models were examined to determine whether all the terms should be retained in the final regression equations. Once the stepwise regression had been applied, the Mallows’ Cp selection method was used as the best selection technique (
Mallows, 1973
). Only models whose parameters were significant at a given level (in our case, 5%) were taken into account.
Stem volume models evaluated.
The best model was chosen with the following indicators: adjusted determination coefficient (
R2
), root mean square error (RMSE) and mean absolute error (MAE). The heterocedasticity was examined visually, by plotting residuals as a function of predicted values and using the Breusch-Pagan test (
Breusch & Pagan, 1979
). Moreover, the non-autocorrelation condition was checked using the Durbin-Watson test (
Durbin & Watson, 1971
). The estimations of the different fitted models were analyzed by means of numerical and graphical analysis. Once the model had been fitted, the Vob in each plot from the IFN4 was calculated in accordance with Ho and G. This value was the one considered for each plot as a Vob value measured in the field and was the dependent variable in the stem volume models in accordance with LiDAR statistics for 2011 and 2015.
Fitting the stem volume model based on 2011 LiDAR data (Vob_2011)
Linear models were used to establish empirical relationships between field measurements and LiDAR variables. In order to fit the Vob for 2011, the Vob calculated in each IFN4 plot based on the model fitted with the independent variables G and Ho was used as the dependent variable. In the STATGRAPHICS regression model selection procedure, all statistics relating to the height of the points obtained with the FUSION CloudMetrics order were included as independent variables. In order to fit the model, 60 plots from the IFN4 were used. The other 30 plots were used to validate the model. The highest model efficiency in accordance with the LiDAR metrics was obtained by ordering the proposed combinations according to the best Mallows Cp. The same goodness-of-fit test utilized for the models in the previous section was used. Likewise, for the selected models, a residuals analysis was made for the selected models, with the aim of identifying evidence of normality deviations, lack of adjustment and/or heterocedasticity.
Fitting the 2015 stem volume model (Vob_2015)
In order to fit the Vob for 2015, the Vob was used as dependent variable, calculated in each one of the 55 plots measured in the field in 2015 with the previously adjusted tariff of G and Ho. In order to maintain the uniformity of the results obtained, the model was fitted with the same independent variable as the Vob_2011 model, in other words, the 70 percentile of the height of the points cloud was used. In order to fit the model, 40 plots were used. The other 15 plots were used to validate the model.
Application of LiDAR data to Vob models
In order to apply the wood volume models obtained, the area occupied by
P. radiata
in Muxika (3,153 ha) was separated and a grid made up of 26.59-m grid cells was superimposed. Each grid cell had a surface area equal to a circular plot with a 15-m radius. Therefore, each grid cell was of 706.86 m
2
. The cells that intersected the polygons corresponding to the surface areas of
P. radiata
according to the forestry map of the Basque Country, were singled out from the vector layer. In this way, a 26.59-m grid mask was obtained with the location of
P. radiata
forests in the municipality of Muxika. Parallel to this, the GridMetrics procedure of FUSION was applied to all the LiDAR files that intersected the municipality of Muxika. In order to apply GridMetrics, a grid cell size equal to that of the 26.59-m grid was used. The work process diagram can be seen in Fig. 4.
Work process diagram
ResultsVob tariff per hectare
The lowest Mallows Cp value was 0.51. This corresponded to the model (3) that had a single independent variable: Ho*G (Table 3). This is a model with a high
R2
(0.99), although it presents heteroscedasticity (
p
<0.05) for the Breusch-Pagan test. This occurs because the higher the volume value, the greater the absolute error.
Vob_2011
The lowest Cp value was 3.28 and this corresponded to the model that had the value of the height of the percentile 70 (h
70
) of the points cloud (Eq. [4], Table 4) as single independent variable. This model had a good
R2
(0.85) but was heterocedastic with a
p
< 0.05 in the Breusch Pagan test. Once again, this occurs because the higher the volume value, the greater the absolute error.
Stem volume models for 2011 based on LiDAR data
Vob_2015
The stem volume model for 2015 based on LiDAR data is shown in Table 5.
Stem volume models for 2015 based on LiDAR data
The models obtained allowed us to estimate the total volume over bark with a different degree of accuracy, although in all cases it is quite high. The best regression model corresponds to the Vob model based on stand variables (H
0
, G). This model presents the highest value of the adjusted coefficient of determination
(R2adj
= 99), followed by the volume estimation model in accordance with LiDAR metrics (Vob_2011,
R2adj
= 85) and lastly, the model based on the field data obtained in 2015 (Vob_2015,
R2adj
= 80) (Table 6). Likewise, the lowest values of RMSE, MAE and ME corresponded to the best model (Vob) obtaining in this the RMSE values (13.61–12.70), MAE (10.04–9.86) and ME = (0–3.49) obtained in the adjustment plots and the validation, respectively (Table 6). Models Vob_2011 and Vob_2015 exhibited very similar values in RMSE (~ 60 m
3
ha
-1
) and MAE (46 m
3
ha
-1
). Residual plots for Vob prediction indicated that the proposed model was appropriate, as can be seem from the predicted versus the observed values for Vob in Fig. 5.
Chart showing models for stem volume and waste Vob model (a), Vob_2011 model (b) and Vob_2015 model (c)
Fits of the models
Each of the adjusted equations showed good fit statistics. All the parameters were significant at a confidence level of 95% (Table 7). In the graphical evaluation of the assumptions of normality and uniformity of the waste of the three fitted models, no evidence was found of any non-compliance with the assumptions.
The mean values of volume with bark of the three models are shown in Table 8.
Parameter estimates, standard errors of the fitted models
Mean stem volume in the models analyzed
Results of applying the Vob models to the entire surface area of Muxika
By applying the volume per hectare of wood with bark models to the entire surface area of Muxika for 2011 and 2015, the results shown in Table 9 were obtained.
The results of the distribution of the Vob for 2011 and 2015 and the difference between both can be seen in Fig. 6.
Distribution of the Vob in 2011 (a), distribution of the Vob in 2015 (b), and
difference 2015 Vob minus 2011 Vob (c).
Application of the Vob models to the entire surface area of Muxika
Annual increase in volume over bark
The difference between the Vob estimated for 2015 and that estimated for 2011 is 81,603 m
3
for the entire surface area of P. radiata in the municipality of Muxika. Therefore, the annual mean growth of the volume over bark (Vob) obtained during the 5-yr period was 5.14 m
3
ha
-1
yr
-1
(Table 9).
The forestry harvests data shows that the wood volume extracted during the period 2011-2015 was 135,624 m
3
(Table 2) and therefore, on average, 8.60 m
3
ha
-1
were extracted. The estimated mean growth for the period 2011 - 2015 is 13.74 m
3
ha
-1
yr
-1
, this being the value obtained by adding the estimated mean growth in the aforementioned period (Vob = 5.14 m
3
ha
-1
yr
-1
) and the mean volume extracted 8.6 m
3
ha
-1
.
As can be seen in Fig. 7, the fitted models predict the greatest increases for Vob of between 0 and 150 m
3
ha
-1
. The increases diminish gradually for Vob values close to 300 m
3
ha
-1
and become negative above that value.
Vob increase versus Vob_2011
Discussion
Both the
R2
and the RMSE of the model fitted for the Vob of 2011 are similar to those published by other authors under similar conditions (
Packalen
et al.
, 2011
;
Stone
et al.
, 2011
;
Tonolli
et al.
, 2011
;
Watt
et al.
, 2013
). The values of
R2
and RMSE of the Vob model for 2015 are more unfavorable than those obtained under similar working conditions. This is due to the fact that in the preparation of the Vob_2011 model, the field data and the LiDAR data were almost contemporary,
i.e.
, they were taken just a few months apart. However, for the fit of the Vob_2015 model, the field data is more than four years older than the LiDAR data. This time difference causes a reduction in the ability to predict the fitted model, even more so if we take into consideration that
P. radiata
is a fast-growing species and is normally subjected to intensive forestry regimes that bring about considerable changes within a short time span.
With regard to adapting the results obtained, the estimated mean growth being 13.74 m
3
ha
-1
yr
-1
, the volume will have increased 68.7 m
3
ha
-1
on average during the period 2011-2015. Therefore, in those places where harvest was more than 68.7 m
3
ha
-1
, the growth has failed to offset the volume extracted and as a result will have less volume in 2015 than in 2011. These are the grid cells with a negative growth value. On the other hand, the mean volume extracted in a thinning operation according to Table 2 is 39.2 m
3
ha
-1
. It is for this reason that the areas that were subjected to thinning between 2011 and 2015 must correspond to cells with a volume difference of less than 29.5 m
3
ha
-1
as the 68.7 m
3
of growth in the 5 yrs, minus the 39.2 m
3
extracted in thinning operations will have allowed this amount to accumulate on average in these areas. Those grid cells with a volume difference of under that amount add up to 1,552 ha. This is the surface area in which, according to the fitted models, some kind of harvest took place during 2011 and 2015 and coincides with those contributed as real registered harvest: 1068.1 ha. As shown in Table 1, the minimum Vob of the plots measured in the field data taken in 2015 was 90.84 m
3
ha
-1
and, therefore, no plot in which clear cutting had taken place was measured. Consequently, the 2015 model cannot detect the clear cutting that took place during this period and assimilates these to intensive thinning, and therefore the surface area of harvest estimated by the models includes both clear cutting and thinning. This circumstance must be borne in mind for subsequent applications of the method and those plots subjected to clear cutting during the period of use of LiDAR data. It can be seen, moreover, that the areas in which the models predict greater harvest intensity correspond to those that had a higher volume in 2011. These areas are the ones in which third thinnings and clear cuttings (being the harvest of greater weight) are normally carried out. Both the assimilation of clear cutting to intensive thinning and the fact that the areas that were actually thinned or clearcut after July 2012 are captured in the raster models, which are based on a Lidar dataset from July 2012, can be explained when analyzing the diagram resulting from a comparison of Vob_2011 and the growth between 2011 and 2015 (Fig. 7).
One of the aspects that must be clarified in subsequent studies is the development in growth in annual periods between both ends of the 5-yr interval. In this study it has not been possible to deal with this question due to a lack of field data in this intermediate period. Moreover, due to a lack of funding, locating the IFN plots was not possible for the second gathering of field data and for this reason the samples for 2011 and 2015 were different. Even though this circumstance might reduce the validity of the results obtained, it was decided to begin the study to verify whether this methodology could be used in subsequent studies with higher funding levels. Once sufficient funds are available, field data must be taken on an annual basis and in the same plot sample.
With regard to the multi-annual use and estimation of growth using a single LiDAR flight, no published references have been found concerning similar methods. This is the main contribution of this paper: the use of LiDAR data in successive years through the updating of the mass conditions with fast, reliable and cost-effective annual field data gathering. In this way, the investment involved in LiDAR data gathering is justifiable as this can be amortized in several years.
Fitted models that relate the descriptive variables of the points cloud of a LiDAR flight with the wood volume field data obtained in field plots with a 5-yr interval, provides a good approximation to the growth values produced in masses of
P. radiata
in this region (13.74 m
3
ha
-1
yr
-1
), as well as an estimation of the surface area subjected to harvest in this period. This allows the multi-annual use of LiDAR data gathered in a campaign with the consequent decrease in costs. To ensure the highest levels of accuracy in the estimation of the growth and volume of wood during the years following the gathering of LiDAR data, it is necessary to include in the sample of plots to be measured those in which both thinning and clear cutting operations have been carried out.
Acknowledgements
Our sincere thanks to Hazi Fundazioa for providing essential data for this study.
ReferencesBabcockC
,
FinleyAO
,
BradfordJB
,
KolkaR
,
BirdseyR
,
RyanMG
,
2015
.
169
:
113
-
127
.
https://doi.org/10.1016/j.rse.2015.07.028BrandtbergT
,
2007
.
Classifying individual tree species under leaf-off and leaf-onconditions using airborne LiDAR. ISPRS J Photogramm.61
:
325
-
340
.
https://doi.org/10.1016/j.isprsjprs.2006.10.006BreuschTS
,
PaganAR
,
1979
.
A simple test for heteroscedaticity and random coefficient variation.47
:
1287
-
1294
.
https://doi.org/10.2307/1911963CaoL
,
CoopsNC
,
InnesJL
,
SheppardSR
,
FuL
,
RuanH
,
SheG
,
2016
.
Estimation of forest biomass dynamics in subtropical forests using multi-temporal airborne LiDAR data.178
:
158
-
171
.
https://doi.org/10.1016/j.rse.2016.03.012Condes RuizS
,
Fernández-LandaA
,
RodríguezF
,
2013
.
Influencia del inventario de campo en el error de muestreo obtenido en un inventario con tecnología LidarDalponteM
,
CoopsN
,
BruzzoneL
,
GianelleD
,
2009
.
Analysis on the use of multiple returns LiDAR data for the estimation of tree stems volume.2
:
310
-
318
.
https://doi.org/10.1109/JSTARS.2009.2037523DeanTJ
,
CaoQV
,
RobertsSD
,
EvansDL
,
2009
.
Measuring heights to crownbase and crown median with LiDAR in a mature even-aged loblolly pine stand.257
:
126
-
133
.
https://doi.org/10.1016/j.foreco.2008.08.024DurbinJ
,
WatsonGS
,
1971
.
Testing for serial correlation in least squares regression.58
:
1
-
19
.
https://doi.org/10.2307/2334313FeketyPA
,
FalkowskiMJ
,
HudakAT
,
2014
.
Temporal transferability of LiDAR-based imputation of forest inventory attributes.45
(
4
):
422
-
435
.
https://doi.org/10.1139/cjfr-2014-0405GiannicoV
,
LafortezzaR
,
JohnR
,
SanesiG
,
PesolaL
,
ChenJ
,
2016
.
Estimating stand volume and above-ground biomass of urban forests using LiDAR84339https://doi.org/10.3390/rs8040339HolmgrenJ
,
2004
.
Prediction of tree height basal area and stem volume in forest stands using airborne laser scanning.19
:
543
-
553
.
https://doi.org/10.1080/02827580410019472HopkinsonC
,
ChasmerL
,
HallRJ
,
2008
.
The uncertainty in conifer plantation growth prediction from multi-temporal lidar datasets.112
(
3
):
1168
-
1180
.
https://doi.org/10.1016/j.rse.2007.07.020HyyppäJ
,
KelleO
,
LehikoinenM
,
InkinenM
,
2001
.
A segmentation-based method to retrieve stem volume estimates from 3-D tree height models produced by laser scanners.39
(
5
):
969
-
975
.
https://doi.org/10.1109/36.921414JochemA
,
HollausM
,
RutzingerM
,
HöfleB
,
2010
.
Estimation of aboveground biomass in alpine forests: A semi-empirical approach considering canopy transparency derived from airborne LiDAR data.11
(1):
278
-
295
.
https://doi.org/10.3390/s110100278KoenigK
,
HöfleB
,
2016
.
Full-waveform airborne laser scanning in vegetation studies-A review of point cloud and waveform features for tree species classification79198https://doi.org/10.3390/f7090198LatifiH
,
NothdurftA
,
KochB
,
2010
.
Non-parametric prediction and mapping of standing wood volume and biomass in a temperate forest: application of multiple optical/LiDAR-derived predictors.83
(
4
):
395
-
407
.
https://doi.org/10.1093/forestry/cpq022LeitererR
,
FurrerR
,
SchaepmanME
,
MorsdorfF
,
2015
.
Forest canopy-structure characterization: A data-driven approach.358
:
48
-
61
.
https://doi.org/10.1016/j.foreco.2015.09.003LiR
,
WeiskittelAR
,
2010
.
Comparison of model forms for estimating stem taper and volume in the primary conifer species of the North American Acadian Region.673302https://doi.org/10.1051/forest/2009109LindnerM
,
KarjalainenT
,
2007
.
Carbon inventory methods and carbon mitigation potentials of forests in Europe: a short review of recent progress.126
:
149
-
156
.
https://doi.org/10.1007/s10342-006-0161-3MagnussenS
,
NæssetE
,
GobakkenT
,
2010
.
Reliability of LiDAR derived predictors of forest inventory attributes: A case study with Norway spruce.1144700
-
712
.
https://doi.org/10.1016/j.rse.2009.11.007MallowsCL
,
1973
.
Some comments on CP.15
(
4
):
661
-
675
.
MaltamoM
,
PackalénP
,
PeuhkurinenJ
,
SuvantoA
,
PesonenA
,
HyyppäJ
,
2007
.
Experiences and possibilities of ALS based forest inventory in Finland270279MaltamoM
,
NæssetE
,
VauhkonenJ
,
2014
.
Forestry applications of airborne laser scanning. Concepts and case studies27460MateosE
,
GarridoF
,
OrmaetxeaL
,
2016
.
Assessment of biomass energy potential and forest carbon stocks in Biscay (Spain).7
(
4
):
75
-
89
.
https://doi.org/10.3390/f7040075MontaghiA
,
CoronaP
,
DalponteM
,
GianelleD
,
ChiriciG
,
OlssonH
,
2013
.
Airborne laser scanning of forest resources: an overview of research in Italy as a commentary case study.23
:
288
-
300
.
https://doi.org/10.1016/j.jag.2012.10.002NæssetE
,
GobakkenT
,
2005
.
Estimating forest growth using canopy metricsderived from airborne laser scanner data.96
:
453
-
465
.
https://doi.org/10.1016/j.rse.2005.04.001NæssetE
,
2011
.
Estimating above-ground biomass in young forests with airborne laser scanning.32
:
473
-
501
.
https://doi.org/10.1080/01431160903474970NilssonM
,
1996
.
Estimation of tree heights and stand volume using an airborne lidar system.56
:
1
-
7
.
https://doi.org/10.1016/0034-4257(95)00224-3PackalénP
,
MehtätaloL
,
MaltamoM
,
2011
.
ALS-based estimation of plot volume and site index in a eucalyptus plantation with a nonlinear mixed-effect model that accounts for the clone effect.68
(
6
):
1085
-
1092
.
https://doi.org/10.1007/s13595-011-0124-9PanYBirdseyRFangJHoughtonRKauppiPKruzWPhillipsOShvidenkoALewisSA large and persistent carbon sink in the world's forests333988993https://doi.org/10.1126/science.1201609PopescuSC
,
2007
.
Estimating biomass of individual pine trees using airborne Lidar.31
:
646
-
655
.
https://doi.org/10.1016/j.biombioe.2007.06.022PopescuSC
,
Zhao KA
,
2008
.
Voxel-based LiDAR method for estimating crownbase height for deciduous and pine trees.112
:
767
-
781
.
https://doi.org/10.1016/j.rse.2007.06.011PopescuSC
,
WynneRH
,
NelsonRF
,
2003
.
Measuring individual tree crown diameter with lidar and assessing its influence on estimating forest volume and biomass.29
(
5
):
564
-
577
.
https://doi.org/10.5589/m03-027RiañoD
,
ChuviecoE
,
CondésS
,
González-MatesanzJ
,
UstinSL
,
2004
.
Generation of crown bulk density for Pinus sylvestris L from lidar.92
(
3
):
345
-
352
.
https://doi.org/10.1016/j.rse.2003.12.014RodríguezF
,
BrotoM
,
LizarraldeI
,
2008
.
CubiFOR: Herramienta para cubicar clasificar productos y calcular biomasa y CO2 en masas forestales de Castilla y León.95
:
33
-
39
.
SongY
,
ImanishiJ
,
SasakiT
,
IokiK
,
MorimotoY
,
2016
.
Estimation of broad-leaved canopy growth in the urban forested area using multi-temporal airborne LiDAR datasets.16
:
142
-
149
.
https://doi.org/10.1016/j.ufug.2016.02.007StoneC
,
PenmanTD
,
TurnerR
,
2011
.
Determining an optimal model for processing lidar data at the plot level: results for a Pinus radiata plantation in New South Wales Australia.41
:
191
-
205
.
TonolliS
,
DalponteM
,
VescovoL
,
RodeghieroM
,
BruzzoneL
,
GianelleD
,
2011
.
Mapping and modeling forest tree volume using forest inventory and airborne laser scanning.130
(
4
):
569
-
577
.
https://doi.org/10.1007/s10342-010-0445-5TorresAB
,
MarchantR
,
LovettJC
,
SmartJC
,
TipperR
,
2010
.
Analysis of the carbon sequestration costs of afforestation and reforestation agroforestry practices and the use of cost curves to evaluate their potential for implementation of climate change mitigation.69
(
3
):
469
-
477
.
https://doi.org/10.1016/j.ecolecon.2009.09.007VepakommaU
,
St-OngeB
,
KneeshawD
,
2011
.
Response of a boreal forest to canopy opening: assessing vertical and lateral tree growth with multi-temporal LiDAR data.21
:
99
-
121
.
https://doi.org/10.1890/09-0896.1WattM
,
MeredithA
,
WattP
,
GunnA
,
2013
.
Use of LiDAR to estimate stand characteristics for thinning operations in young Douglas-fir plantations.43
(
18
):
1
-
10
.
https://doi.org/10.1186/1179-5395-43-1YuX
,
HyyppäJ
,
KaartinenH
,
MaltamoM
,
2004
.
Automatic detection of harvested trees and determination of forest growth using airborne laser scanning.90
:
451
-
462
.
https://doi.org/10.1016/j.rse.2004.02.001YuX
,
HyyppäJ
,
KukkoA
,
MaltamoM
,
KaartinenH
,
2006
.
Change detection techniques for canopy height growth measurements using airborne laser scanner data Photogramm.72
:
1339
-
1348
.
https://doi.org/10.14358/PERS.72.12.1339YuX
,
HyyppäJ
,
HolopainenM
,
VastarantaM
,
2010
.
Comparison of area-based and individual tree-based methods for predicting plot-level forest attributes.6
:
1481
-
1495
.
https://doi.org/10.3390/rs2061481ZhangJ
,
GeY
,
ChangJ
,
JiangB
,
JiangH
,
PengC
,
ZhuJ
,
YuanW
,
QiL
,
YuS
,
2007
.
Carbon storage by ecological service forests in Zhejiang Province subtropical China.245
:
64
-
75
.
ZhangZ
,
KazakovaA
,
MoskalLM
,
StyersDM
,
2016
.
Object-based tree species classification in urban ecosystems using LiDAR and hyperspectral data76122https://doi.org/10.3390/f7060122ZianisD
,
MuukkonenP
,
MäkipääR
,
MencucciniM
,
2005
.
Biomass and stem volume equations for tree species in Europe463ZolkosSG
,
GoetzSJ
,
DubayahRA
,
2013
.
Meta-analysis of terrestrial aboveground biomass estimation using LiDAR remote sensing.128
:
289
-
298
.
https://doi.org/10.1016/j.rse.2012.10.017