**Remote Sensing-Based Biomass Estimation**

José Mauricio Galeana Pizaña, Juan Manuel Núñez Hernández and Nirani Corona Romero

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/61813

#### **Abstract**

Over the past two decades, one of the research topics in which many works have been done is spatial modeling of biomass through synergies between remote sensing, for‐ estry, and ecology. In order to identify satellite-derived indices that have correlation with forest structural parameters that are related with carbon storage inventories and forest monitoring, topics that are useful as environmental tools of public policies to focus areas with high environmental value. In this chapter, we present a review of dif‐ ferent models of spatial distribution of biomass and resources based on remote sens‐ ing that are widely used. We present a case study that explores the capability of canopy fraction cover and digital canopy height model (DCHM) for modeling the spatial distribution of the aboveground biomass of two forests, dominated by *Abies Religiosa* and *Pinus* spp., located in Central Mexico. It also presents a comparison of different spatial models and products, in order to know the methods that achieved the highest accuracy through root-mean-square error. Lastly, this chapter provides con‐ cluding remarks on the case study and its perspectives in remote sensing-based bio‐ mass estimation.

**Keywords:** Aboveground biomass, forest, remote sensing, modeling, resources

#### **1. Introduction**

Forest ecosystems are about 31% of the total land cover of the earth [1], being one of the most important ecosystems due to economic goods and environmental services they provide. One of these services is as an environmental regulator, reducing the concentration of carbon dioxide (greenhouse gas) from the atmosphere and transforming it into oxygen and biomass through photosynthesis, thereby playing an important role in the global carbon cycle [2–4].

© 2016 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Biomass is defined as the dry weight of both aboveground biomass (AGB) and belowground biomass (BGB) living mass of vegetation, such as wood, bark, branches, twigs, stumps, or roots as well as dead mass of litter associated with the soil [4–6]. According to this, it can be considered as a measure of forest structure and function. Thus, by knowing the spatial distribution of biomass, it is possible to calculate the net flow of terrestrial carbon, nutrient cycling, forest productivity, biomass energy, and carbon storage and sequestration by the forest, reducing the uncertainty of carbon emission and sequestration measures to support climate change modeling studies [5–8].

Since calculating field measures of BGB is difficult, most studies have focused on calculating AGB. The most accurate way to obtain AGB data is by using field measurements and allometric equations for individual trees; however, these techniques are difficult to implement because they are time consuming and labor intensive. Furthermore, forests are a complex and widely distributed ecosystem, which makes these techniques expensive to apply in large areas; therefore, they cannot provide the spatial distribution of biomass [4,6].

An alternative form to map and monitor spatial distribution of AGB is through the use of remote sensing-based techniques, because through them it is possible to obtain a continuous and repetitive collection of digital data from the same area with different spatial resolutions, covering large areas and reducing processing time and costs [5,6]. Due to the increasing availability of satellite imagery, several researches have been developed to prove the effec‐ tiveness of both imagery data provided by different sensors and diverse modeling approaches [4,6,9]. In order to estimate AGB, two main kinds of models have been used: direct and indirect. The direct models measure biomass throughout the relationship between spectral data response and biomass field measurements. For the indirect models, biomass is estimated from biophysical parameters or forest structural metrics [10].

This chapter reviews the main models used for estimating biomass and key resources used in remote sensing (Sections 1 and 2). The case study integrates some models and resources applied in forests dominated by *Abies Religiosa* and *Pinus* spp. located in Central Mexico. In the last section, concluding remarks are provided about the best biomass estimation models as well as their limitations.

#### **2. Biomass modeling**

Several factors can affect the remote sensing-based AGB estimation, such as insufficient sample data, atmospheric conditions, complex biophysical environments, scale of the study area, availability of software, spatial resolution of remotely sensed data, or mixed pixels, among others [6, 10]. In order to introduce different approaches that have been developed to reduce the uncertainties produced by these factors on the estimation of the spatial distribution of AGB, the most commonly used models will be described in this section.

#### **2.1. Field measurements and allometric equations**

Biomass is defined as the dry weight of both aboveground biomass (AGB) and belowground biomass (BGB) living mass of vegetation, such as wood, bark, branches, twigs, stumps, or roots as well as dead mass of litter associated with the soil [4–6]. According to this, it can be considered as a measure of forest structure and function. Thus, by knowing the spatial distribution of biomass, it is possible to calculate the net flow of terrestrial carbon, nutrient cycling, forest productivity, biomass energy, and carbon storage and sequestration by the forest, reducing the uncertainty of carbon emission and sequestration measures to support

Since calculating field measures of BGB is difficult, most studies have focused on calculating AGB. The most accurate way to obtain AGB data is by using field measurements and allometric equations for individual trees; however, these techniques are difficult to implement because they are time consuming and labor intensive. Furthermore, forests are a complex and widely distributed ecosystem, which makes these techniques expensive to apply in large areas;

An alternative form to map and monitor spatial distribution of AGB is through the use of remote sensing-based techniques, because through them it is possible to obtain a continuous and repetitive collection of digital data from the same area with different spatial resolutions, covering large areas and reducing processing time and costs [5,6]. Due to the increasing availability of satellite imagery, several researches have been developed to prove the effec‐ tiveness of both imagery data provided by different sensors and diverse modeling approaches [4,6,9]. In order to estimate AGB, two main kinds of models have been used: direct and indirect. The direct models measure biomass throughout the relationship between spectral data response and biomass field measurements. For the indirect models, biomass is estimated from

This chapter reviews the main models used for estimating biomass and key resources used in remote sensing (Sections 1 and 2). The case study integrates some models and resources applied in forests dominated by *Abies Religiosa* and *Pinus* spp. located in Central Mexico. In the last section, concluding remarks are provided about the best biomass estimation models

Several factors can affect the remote sensing-based AGB estimation, such as insufficient sample data, atmospheric conditions, complex biophysical environments, scale of the study area, availability of software, spatial resolution of remotely sensed data, or mixed pixels, among others [6, 10]. In order to introduce different approaches that have been developed to reduce the uncertainties produced by these factors on the estimation of the spatial distribution of AGB,

therefore, they cannot provide the spatial distribution of biomass [4,6].

biophysical parameters or forest structural metrics [10].

the most commonly used models will be described in this section.

as well as their limitations.

**2. Biomass modeling**

climate change modeling studies [5–8].

4 Environmental Applications of Remote Sensing

All spatial distribution AGB estimation models need high quality and representative field data to be implemented. Therefore, forest inventories are the most common approach to obtain detailed and periodic data. They are held for monitoring, modeling, and predicting several biophysical processes, such as stocking levels, harvests, diseases, and pests, among others. Therefore, they are generally implemented at several scales to obtain different structural parameters, either through the data aggregated from stand-level management inventories or by plots established through a statistical sampling design [3,4,6,11].

Since the grouping of data from stand-level management inventories tend to underestimate the forested areas and stock volume calculation, presently in most of the countries a statistical sampling design is used. In this procedure, the sampling plots are randomly selected from a population where each one of them has the same positive known probability to be chosen. They can be selected by different methods. One of them is a systematic sampling procedure based on the use of grids of randomly selected points in two dimensions, with a 0.5–20-km separation range. In other cases, plot locations are randomly selected by regular polygons created by a tessellation of large areas, which can be stratified when different sampling intensities are required, for example, when different land uses and covers exist [11].

In order to make a more efficient sampling, a plot generally conforms clusters, commonly of four plots but, in some cases, they can be as large as 18, with a huge variety of shapes and sizes. Circular plots are commonly used in Boreal and Temperate forests, whereas square and rectangular ones are used in rainforests. Plots can be combined with other sampling methods like transect for measuring deadwood or soil pits to assess soil carbon [11].

Commonly, the information gathered about a plot is location, number of trees, species, health, and site description, among others. In addition, individual tree dendrometric variables are considered, such as diameter at breast height (DBH), tree height, crown size, and canopy cover, DBH and tree height being the most commonly used parameters to derivate AGB through allometric equations [3,4,10].

The allometric biomass equation is the most common and accurate method to translate forest inventory reports of individual tree data to tree and stand biomass. It is a mathematical relation between total tree biomass (stem, branch, and foliage) and its DBH or both DBH and height, applying a least-squares regression of logarithmic equation. It can be both species- or site specific (e.g., *Pinus montezume* or *Abies religiosa*) or more generic (e.g., pine gender or tropical hardwoods); however, it has been observed that biomass equations at the site-specific level produce better results than generic equations [3,8,12].

The most accurate method to obtain allometric equations is a destructive process in which individual trees from a wide range of DBH, distributed in a local forest, must be felled and separated into boles, branches, and leaves. Then boles are cut into sections and weighted in the field as leaves and branches. After that, a thick disc sample must be cut from the base of each bole section and a subsample is extracted for the other side of the disk, both being dried at 105°C until a constant weight is reached to obtain the dry mass, to estimate the moisture content and the wood density. The total biomass therefore is the sum of the dry mass of the branches, the leaves, and the various sections of the stem [3,13].

Since it has been observed that growing plants maintain the weight proportion between different parts, it is possible to build a nonlinear mathematical model to relate biomass with DBH using these parameters [12]. One of the most common models used is the logarithmic equation (1):

$$\operatorname{Lnn}(B) = \operatorname{Lnn}(a) + b\operatorname{Lnn}(D) \tag{1}$$

where *B* is the biomass, *a* and *b* are scalar coefficients estimated by least-squares linear regression, and *D* is the DBH.

Using allometric equations, it is possible to compute the total AGB for a given area using biomass expansion factors or conversion tables; however, these approaches do not provide the spatial distribution of AGB [3,10]. Later in this chapter, different models created to obtain spatial distribution of AGB are explained.

#### **2.2. Regression models**

One of the most common methods to estimate biomass is the regression analysis, which is a statistical technique to investigate and model the relationship between variables. Traditionally, in the remote sensing approach, the regression analysis techniques applied to AGB estimation are based on the quantitative relationship between ground-based data and satellite informa‐ tion, such as spectral reflectance, radar, or light detection and ranging (LiDAR) data [6,14,15]. Models based on regression analysis are considered to be relatively easy to implement and can provide accurate results through their application at all spatial scales [16]. Generally, these methodologies consist of three major steps: biomass estimation based on fieldwork, establish‐ ment of regression model between field biomass and satellite information of corresponding pixels, and the use of regression models to generate a biomass image with the spatial predic‐ tion.

These remote sensing-based biomass estimation methods assume that the forest information, obtained by the sensors, is highly correlated with AGB. According to this, the keys for biomass estimation are the use of appropriate variables and the development of suitable estimation models for sufficient sample plots, using regression methods that aim for an efficient integra‐ tion of multisource data, necessary to get better biomass estimation. Spectral data, radar, and LiDAR have their own positive and negative characteristics and proper integration of them can improve biomass estimation accuracy [17].

After data integration, the correct use of regression methods for establishing biomass estima‐ tion models is also important. Many models have been developed based on multiple combi‐ nations of in situ tree parameters calculated through linear regression (LR) or nonlinear regression (NLR) models [18–20]. Multiple regression analysis may be the most frequently used approach for developing biomass estimation models [21–23]. In both cases, these parametric algorithms assume that the relationships between dependent (e.g., biomass) and independent (derived from remote sensing data) variables have explicit model structures that can be specified a priori by parameters [15]. Generally, the independent variables can be spectral bands, vegetation indices, textural images, LiDAR height, and synthetic aperture radar (SAR) backscatter; in some cases, the use of subpixel information offered better estimation results than per-pixel-based spectral signatures [24].

However, the linear regression approach has been known to mislead the prediction of the studied variable at values beyond a saturation point of the canopy reflectance [25]. Since biomass is usually nonlinearly related to remote sensing variables, nonlinear models such as power models [26], logistic regression models [27], and geographically weighted regression models [16] were often used to estimate biomass with more accuracy. Nonetheless, some estimation methods have been established as a nonparametric alternative to the use of regression approaches for biomass modeling: *k*-nearest neighbor (*k*-NN), artificial neural network (ANN), regression tree, random forest, support vector machine (SVM), and maximum entropy (MaxEnt) [15].

#### **2.3. Geostatistical models**

content and the wood density. The total biomass therefore is the sum of the dry mass of the

Since it has been observed that growing plants maintain the weight proportion between different parts, it is possible to build a nonlinear mathematical model to relate biomass with DBH using these parameters [12]. One of the most common models used is the logarithmic

where *B* is the biomass, *a* and *b* are scalar coefficients estimated by least-squares linear

Using allometric equations, it is possible to compute the total AGB for a given area using biomass expansion factors or conversion tables; however, these approaches do not provide the spatial distribution of AGB [3,10]. Later in this chapter, different models created to obtain

One of the most common methods to estimate biomass is the regression analysis, which is a statistical technique to investigate and model the relationship between variables. Traditionally, in the remote sensing approach, the regression analysis techniques applied to AGB estimation are based on the quantitative relationship between ground-based data and satellite informa‐ tion, such as spectral reflectance, radar, or light detection and ranging (LiDAR) data [6,14,15]. Models based on regression analysis are considered to be relatively easy to implement and can provide accurate results through their application at all spatial scales [16]. Generally, these methodologies consist of three major steps: biomass estimation based on fieldwork, establish‐ ment of regression model between field biomass and satellite information of corresponding pixels, and the use of regression models to generate a biomass image with the spatial predic‐

These remote sensing-based biomass estimation methods assume that the forest information, obtained by the sensors, is highly correlated with AGB. According to this, the keys for biomass estimation are the use of appropriate variables and the development of suitable estimation models for sufficient sample plots, using regression methods that aim for an efficient integra‐ tion of multisource data, necessary to get better biomass estimation. Spectral data, radar, and LiDAR have their own positive and negative characteristics and proper integration of them

After data integration, the correct use of regression methods for establishing biomass estima‐ tion models is also important. Many models have been developed based on multiple combi‐ nations of in situ tree parameters calculated through linear regression (LR) or nonlinear regression (NLR) models [18–20]. Multiple regression analysis may be the most frequently used approach for developing biomass estimation models [21–23]. In both cases, these

Ln Ln Ln (*B ab D* ) = + ( ) ( ) (1)

branches, the leaves, and the various sections of the stem [3,13].

equation (1):

regression, and *D* is the DBH.

6 Environmental Applications of Remote Sensing

**2.2. Regression models**

tion.

spatial distribution of AGB are explained.

can improve biomass estimation accuracy [17].

Some studies have used geostatistics as the main approach for biomass forest estimation in order to predict variables related to forest structure (e.g., tree height and volume) and aboveground biomass and carbon measurements in unsampled sites based on known values of adjacent spatial data as forest inventory sites [28,29]. Other recent works have explored the synergy between geostatistical models with remotely sensed data to improve estimations using remote sensing indices as spatial secondary variables [10, 16, 30–32].

Geostatistics was defined in the 1960s by Georges Matheron, who generalized a set of techni‐ ques developed by Krige (1951) in order to exploit the spatial correlation to make predictions in evaluating reserves of gold mines in South Africa. This generalization is detailed in his regionalized variable theory in 1970 [33]. The purpose of geostatistics is the estimation, prediction, and simulation of the values of a variable that is distributed through space [34]. This theory assumes that a variable measured in a spatial domain corresponds to a random variable *z*(*x*), assuming that the structure of the phenomenon having spatial correlation is considered a regionalized variable; therefore, a set of spatially distributed random variables will be a random function *Z*(*x*). This provides the theoretical basis for establishing the spatial structural characteristics of natural phenomena. Moreover, it can be used as a tool for calcu‐ lating the value of a variable in a certain position in space, knowing the values of that variable among adjacent positions in space, which is known as interpolation [35].

There are diverse methods of interpolations, which can be classified into two main groups: deterministic and geostatistical. The deterministic techniques are based directly on some properties of similarity of adjacent measured values (e.g., distance), which establish a set of mathematical formulas that determine the smoothness of the resultant surface interpolated. Examples of these are the inverse distance weighting (IDW), nearest neighbors, splines, and triangular irregular network (TIN) [35]. Geostatistical techniques studied spatial autocorrela‐ tion of the variables in order to fit a spatial dependence model to a set of random variables. This approach produces predictions and also generates an error surface concerning the uncertainty-associated analyzed model [33].

The spatial dependence is the spatial behavior of a phenomenon, derived from spatial patterns in terms of distance and similarity and/or contrast of a spatial unit or relative spatial location to other spatial units [36]. It has its basis in Tobler's first law of geography proposed in 1970 that says, "everything is related to everything else, but near things are more related than distant things" [37].

The kriging algorithms are one such example, which is mostly used in geosciences, ecology, and geomatics [35]. Kriging is a generic name for a family of generalized least-squares regression techniques, where the spatial structural characteristics are accomplished by the semivariogram function as a metric of the spatial autocorrelation [33].

All kriging estimators are variants of the following basic equation (2):

$$\hat{Z}\{\mathbf{x}\_0\} - \mu = \sum\_{i=1}^{n} \mathbb{X}\_i \Big[ Z(\mathbf{x}\_i) - \mu(\mathbf{x}\_0) \Big] \tag{2}$$

where *µ* is a known stationary mean, assumed to be constant over the whole domain and calculated as the average of the data. The parameter *λ* is the kriging weight; *n* is the number of sampled data points used to make the estimation, and *µ(x*0*)* is the mean of the samples within the search window [33].

Different studies have applied univariate and bivariate geostatistical interpolations in order to calculate the forest volume [28], aboveground biomass [10,16,29–31], and carbon in the aboveground biomass [32]. The most commonly used technique for univariate-based model‐ ing is kriging [28], whereas in bivariate-based modeling regression-kriging [10,30–32], cokriging [31], kriging with external drift [29], cokriging regression [31], and geographically weight regression (GWR) are used [16].

#### **2.4. Nonparametric models**

Similar to regression models, nonparametric algorithms are based on the use of different sensor data, for example, spectral, radar, and LiDAR [21,38], using many of these models in the forest attributes estimation [39–41]. They are a framework for creating complex nonlinear biomass models based on the use of remote sensing variables and as alternatives for the parametric approaches. Common nonparametric algorithms include *k*-nearest neighbor (*k*-NN), artificial neural network (ANN), random forest, support vector machine (SVM), and maximum entropy (Max Ent).

One of the most applied nonparametric methods is the nearest neighbor approach (NN). In the context of forest attribute estimation, the NN methods have been first introduced in the late 1980s [42]. In the NN methods, the value of the target variable at a certain location is predicted as a weighted average of the values of neighboring observations, with *k*-nearest neighbors spectrally using a weighting method [43]. Several methods have been offered to measure the distance from the target unit to the neighbors. In the NN approach, the choice of the *k* value, type of distance measure including Euclidean and Mahalanobis, along with weighed function are the critical factors influencing the estimation accuracy [44,45].

The ANN provides a more robust solution for complicated and nonlinear problems due to its properties [46]. The network commonly consists of one input layer, one or more hidden layers, and one output layer. Since it does not require the assumption that data have normal distri‐ bution and linear relationships between biomass and independent variables, the ANN can deal with different data through approximation, using various complex mathematical functions, with independent variables from different data sources such as remote sensing and ancillary data [15]. A detailed overview of ANN approach is provided in Ref. 47.

Regression tree and random forest are a family of tree-based models; in the first one, data are stratified into homogeneous subsets by decreasing the within-class entropy, whereas in the second one, a large number of regression trees are constructed by selecting random bootstrap samples from the discrete or continuous dataset. In fact, the random forest algorithm is now widely used for biomass estimation [48,49].

SVM is an important method to estimate forest biophysical parameters using remote sensing data [50, 51]. It is a statistical learning algorithm with the ability to use small training sample data to produce relatively higher estimation accuracy than other approaches like ANN [15]. Ref. 51 provides a detailed overview of the SVM approach used in remote sensing. The Max Ent approach is a general-purpose machine-learning method for predicting or inferring target probability distribution from incomplete information [52]. These kinds of nonparametric algorithms have become popular in biomass modeling when large representative field datasets exist for calibration [53,54].

### **3. Remote sensing products for biomass modeling**

#### **3.1. From optical sources**

This approach produces predictions and also generates an error surface concerning the

The spatial dependence is the spatial behavior of a phenomenon, derived from spatial patterns in terms of distance and similarity and/or contrast of a spatial unit or relative spatial location to other spatial units [36]. It has its basis in Tobler's first law of geography proposed in 1970 that says, "everything is related to everything else, but near things are more related than distant

The kriging algorithms are one such example, which is mostly used in geosciences, ecology, and geomatics [35]. Kriging is a generic name for a family of generalized least-squares regression techniques, where the spatial structural characteristics are accomplished by the

> m( ) ( )

*Z x Zx x* (2)

semivariogram function as a metric of the spatial autocorrelation [33]. All kriging estimators are variants of the following basic equation (2):

<sup>ˆ</sup> *<sup>n</sup>*

ml

=

*i*


*i i*

where *µ* is a known stationary mean, assumed to be constant over the whole domain and calculated as the average of the data. The parameter *λ* is the kriging weight; *n* is the number of sampled data points used to make the estimation, and *µ(x*0*)* is the mean of the samples within

Different studies have applied univariate and bivariate geostatistical interpolations in order to calculate the forest volume [28], aboveground biomass [10,16,29–31], and carbon in the aboveground biomass [32]. The most commonly used technique for univariate-based model‐ ing is kriging [28], whereas in bivariate-based modeling regression-kriging [10,30–32], cokriging [31], kriging with external drift [29], cokriging regression [31], and geographically

Similar to regression models, nonparametric algorithms are based on the use of different sensor data, for example, spectral, radar, and LiDAR [21,38], using many of these models in the forest attributes estimation [39–41]. They are a framework for creating complex nonlinear biomass models based on the use of remote sensing variables and as alternatives for the parametric approaches. Common nonparametric algorithms include *k*-nearest neighbor (*k*-NN), artificial neural network (ANN), random forest, support vector machine (SVM), and maximum entropy

One of the most applied nonparametric methods is the nearest neighbor approach (NN). In the context of forest attribute estimation, the NN methods have been first introduced in the late 1980s [42]. In the NN methods, the value of the target variable at a certain location is predicted as a weighted average of the values of neighboring observations, with *k*-nearest neighbors spectrally using a weighting method [43]. Several methods have been offered to

( )

uncertainty-associated analyzed model [33].

8 Environmental Applications of Remote Sensing

things" [37].

the search window [33].

**2.4. Nonparametric models**

(Max Ent).

weight regression (GWR) are used [16].

Optical sensors are those that detect electromagnetic radiation emitted or reflected from the earth, the main source of light being the sun. Among the passive sensors are photographic and optical-electronic sensors that combine similar photographic optics and electronic detection system (detectors and push scanning) and image spectrometers [55,56]. Optical remote sensing refers to methods and technologies that acquire information from the visible, near-infrared, shortwave infrared, and thermal infrared regions of the electromagnetic spectrum. They are called optical, because energy is directed through optical components such as lenses and mirrors.

#### **3.2. Spectral indices**

Remotely sensed spectral vegetation indices represent an integrative measure of both vegeta‐ tion photosynthetic activity and canopy structural variation that are widely used and have benefited numerous disciplines interested in the assessment of biomass estimation [57,58]. Different kinds of vegetation indices have been in use for a long time in AGB estimation and the number of publications is immense [6].

The key to the development of vegetation indices is the ability of the canopy of green vegetation to interact differently with certain portions of the electromagnetic spectrum. Since this contrast is particularly strong between red and near-infrared regions (NIR), it has been the focus of a large variety of attempts to build up quantitative indices of vegetation condition using remotely sensed imagery [59]. Theoretically, the ideal vegetation indices should be particularly sensitive to different vegetation covers, insensitive to soil brightness and color, and little affected by atmospheric effects [60]. However, in reality, different factors affect the reflectance of vegetation and consequently the vegetation index (e.g., atmospheric correction is essential when biomass is extracted from the vegetation indices as a final product).

According to Ref. 58, vegetation indices can be classified into (1) slope-based, (2) distancebased, and (3) transformation indices. The slope-based indices are simple arithmetic combi‐ nations that focus on the contrast between the spectral responses patterns of vegetation in the red and near-infrared portions of the electromagnetic spectrum. The most known of them are the ratio vegetation index (RVI) proposed by Birth and Mc Vey [61], normalized difference vegetation index (NDVI) introduced by Rouse et al. [62], and soil-adjusted vegetation index (SAVI) developed by Huete [63]. The vegetation indices show better sensitivity than individual spectral bands for the detection of biomass [64].

Distance-based indices measure the degree of vegetation from the soil background (known as a soil line) to the pixel with the highest content of vegetation in a perpendicular incremental distance. In this group of indices, the slope and intercept of the soil line have to be defined for each particular image; perpendicular vegetation index (PVI) introduced by Richardson and Weigand [65] cancels the effect of soil brightness in cases where vegetation is sparse and the pixels contain a mixture of green vegetation and soil background. The effect of the background soil is a major limiting factor in certain statistical analyses geared toward the quantitative assessment of AGB [59].

Transformation indices are transformations of the available spectral bands to form a new set of uncorrelated bands within which a vegetation index band can be defined. Tasseled cap (TC) is one of the most widely used indexes of this type and may apply to various remote sensing images with multiple resolutions [66–72].

In the specialized literature, many vegetation indices have been proposed, and depending on the complexity of the forest stand structure, indices vary in their relationships with biomass [23,60,73]. For forest sites with complex stand structures, vegetation indices including nearinfrared wavelength have weaker relationships with biomass than those including shortwave infrared wavelength. In contrast, for areas with poor soil conditions and relatively simple forest stand structure, near-infrared vegetation indices had a strong relationship with biomass, and finally the results of transformation indices showed stronger relationships with biomass independent of different biophysical conditions [15].

#### **3.3. Textural indices**

Different kinds of vegetation indices have been in use for a long time in AGB estimation and

The key to the development of vegetation indices is the ability of the canopy of green vegetation to interact differently with certain portions of the electromagnetic spectrum. Since this contrast is particularly strong between red and near-infrared regions (NIR), it has been the focus of a large variety of attempts to build up quantitative indices of vegetation condition using remotely sensed imagery [59]. Theoretically, the ideal vegetation indices should be particularly sensitive to different vegetation covers, insensitive to soil brightness and color, and little affected by atmospheric effects [60]. However, in reality, different factors affect the reflectance of vegetation and consequently the vegetation index (e.g., atmospheric correction is essential

According to Ref. 58, vegetation indices can be classified into (1) slope-based, (2) distancebased, and (3) transformation indices. The slope-based indices are simple arithmetic combi‐ nations that focus on the contrast between the spectral responses patterns of vegetation in the red and near-infrared portions of the electromagnetic spectrum. The most known of them are the ratio vegetation index (RVI) proposed by Birth and Mc Vey [61], normalized difference vegetation index (NDVI) introduced by Rouse et al. [62], and soil-adjusted vegetation index (SAVI) developed by Huete [63]. The vegetation indices show better sensitivity than individual

Distance-based indices measure the degree of vegetation from the soil background (known as a soil line) to the pixel with the highest content of vegetation in a perpendicular incremental distance. In this group of indices, the slope and intercept of the soil line have to be defined for each particular image; perpendicular vegetation index (PVI) introduced by Richardson and Weigand [65] cancels the effect of soil brightness in cases where vegetation is sparse and the pixels contain a mixture of green vegetation and soil background. The effect of the background soil is a major limiting factor in certain statistical analyses geared toward the quantitative

Transformation indices are transformations of the available spectral bands to form a new set of uncorrelated bands within which a vegetation index band can be defined. Tasseled cap (TC) is one of the most widely used indexes of this type and may apply to various remote sensing

In the specialized literature, many vegetation indices have been proposed, and depending on the complexity of the forest stand structure, indices vary in their relationships with biomass [23,60,73]. For forest sites with complex stand structures, vegetation indices including nearinfrared wavelength have weaker relationships with biomass than those including shortwave infrared wavelength. In contrast, for areas with poor soil conditions and relatively simple forest stand structure, near-infrared vegetation indices had a strong relationship with biomass, and finally the results of transformation indices showed stronger relationships with biomass

when biomass is extracted from the vegetation indices as a final product).

the number of publications is immense [6].

10 Environmental Applications of Remote Sensing

spectral bands for the detection of biomass [64].

assessment of AGB [59].

images with multiple resolutions [66–72].

independent of different biophysical conditions [15].

Some studies have used regression analyses between remote sensing textural indices and biomass data from sampling sites [30,32,74–77]. This framework has been applied to different optical and synthetic aperture radar-derived indices in order to use the textural parameters as continuous spatial variables to improve biomass estimations.

In 1973, Haralick proposed a statistical analysis based on a set of parameters according to spatial dependence of gray tones in an image, defined as second-order statistics [78]. Textures are intrinsic properties of surfaces and its importance lies in image–objects segmentation, because they are related to structural arrangements of the land surface and the connections between neighboring spatial objects [79]. The most common mathematical method used to measure texture parameters is the co-occurrence matrix (gray-level co-occurrence matrix) [78]. It describes the frequency of a gray level displayed, in a specific spatial relationship, to another gray value within a neighborhood represented by a mobile window or kernel [79]. Its con‐ struction is based on four steps: (1) Window or kernel size definition, (2) band selection input, (3) texture parameters selection, and (4) spatial dependency criteria.

Textures applied parameters are described below:


**Table 1.** Textures

P*ij* = P*r* (I*s* = *i* ∩ I*t* = *j)* = pixel probability of *s* is *i* and *t* is *j*, for separated pixels through one pixel distance in a relative displacement vector between neighboring pixels.

These texture parameters and the indices derived (e.g., ratios) have been applied in different satellite inputs, for example, on high spatial resolution optical-infrared images as SPOT-V images [30,32] and ALOS AVNIR-2 [74], medium spatial resolution resources as Landsat TM [76] and Landsat OLI [76], and SAR images as Jers-1 [76] and ALOS Palsar [32,75].

#### **3.4. Biophysical variables**

A forest biophysical parameter is a measure that simplifies the aboveground organization of plant materials [4]. In this context, several studies in remote sensing field have focused on determining canopy structural parameters such as leaf area index (LAI), canopy height or canopy fraction [80–86], the second one being the most commonly used for biomass estimation. Because canopy height is most commonly derived from active sensors [87], this part of the chapter is focused on canopy fraction.

The reflectance of forest canopy cover recorded by the instantaneous field of view (IFOV) of the sensor is a spectral mixture obtained from the interaction between electromagnetic radiation and both canopy and forest elements such as nonphotosynthetic vegetation (NPV; such as branches, stem, and litter), photosynthetic vegetation (PV; such as leaves), and others, such as bare soil and shadow. Therefore, image pixels are generally composed of more than one element, making the image interpretation difficult, which can result in a poor relation between AGB and spectral bands [88–91].

One of the most widely used remote sensing approaches to derive and extract fraction covers from mixed pixels is the spectral mixture analysis (SMA) [92]. It decomposes the mixed pixel using a collection of constituent spectra (end members) to obtain their areal proportions or abundances in a pixel, and therefore unmixing a multispectral image into fraction images of end members [88,93]; it can be linear or nonlinear. The linear model assumes a single interac‐ tion between each incident photon and the surface object, and therefore the mixed pixel is a linear combination of pure spectral signatures (end members) of the surface materials weight‐ ed by their area covered. The nonlinear model is the opposite of linear model, since electro‐ magnetic radiation can intercept more than one element of surface, with mixed pixels resulting from a multiple-scattered signal [94].

As green leaves scatter radiation at NIR spectra, vertical structure of vegetation commonly produces a multiple-scattered signal; however, nonlinear spectral mixture approaches are barely used because they require more specific information, such as scattering properties of end members, the illumination of the sensor, and certain geometrical parameters of the scene [91,93]. For these reasons, linear approaches have been largely implemented [92].

The linear mixture model is expressed in matrix form in the following equation (3):

$$p = Mf + \varepsilon \tag{3}$$

where *p* = *[p*1*...pn]T* is the mixed pixel, *f* = *[f*1*...fm]T* is the fractional end-member abundance, *M* is an *n* ×*m* matrix with *n* end-member spectra as column vectors, and *ε* is the residual not explained by the model. If *p* and *M* are known, it is possible to estimate it from ordinary leastsquares procedure [93,95]), with two common constraints: the full additive condition, which determines that the fraction must sum to one, and the nonnegative condition, which makes all abundances nonnegative, thus making the model physically meaningful [92].

SMA has been successfully applied in vegetation studies to estimate the land cover fractions of PV, NPV, bare soil, and shade [96–98] or mapping the fractional cover of coniferous species [99–101]. In biomass studies, it has been proved that using ASTER fraction images (green vegetation, soil, and shade) in regression models improve the AGB estimation in Mediterra‐ nean forests, being better than NDVI or tasseled cap components [88]. Ref. [15] uses the same Landsat Thematic Mapper (TM) fractions and TM spectral signatures to relate with AGB, finding that fractional images perform better for successional forest biomass estimation than TM spectral signatures, but not for primary biomass estimation. SMA has also been applied to remove subpixel atmospheric and soil reflectance contamination in order to improve dry biomass estimations, showing that unmixed vegetation indices are better [102] than those which are not unmixed. Through geometric-optical reflectance models, Peddle et al. also estimated areal fractions of sunlit canopy, sunlit background, and shadow at subpixel scales showing higher accuracy than NDVI [103]. The other case in which sunlit crowns, background, and shadow fractions were compared with seven different vegetation indices (NDVI, SR, MSR, RDVI, WDVI, GEMI, and NLI) and three different soil-adjusted vegetation indices (SAVI, SAVI-1, and SAVI-2) to estimate biomass, LAI, net primary productivity (NPP), DBH, stem density, and basal area was the study conducted by Peddle et al. In this study, the authors concluded that the SMA shadow fraction improves the results by about 20% compared to vegetation indices [104].

#### **3.5. From active sources**

P*ij* = P*r* (I*s* = *i* ∩ I*t* = *j)* = pixel probability of *s* is *i* and *t* is *j*, for separated pixels through one pixel

These texture parameters and the indices derived (e.g., ratios) have been applied in different satellite inputs, for example, on high spatial resolution optical-infrared images as SPOT-V images [30,32] and ALOS AVNIR-2 [74], medium spatial resolution resources as Landsat TM

A forest biophysical parameter is a measure that simplifies the aboveground organization of plant materials [4]. In this context, several studies in remote sensing field have focused on determining canopy structural parameters such as leaf area index (LAI), canopy height or canopy fraction [80–86], the second one being the most commonly used for biomass estimation. Because canopy height is most commonly derived from active sensors [87], this part of the

The reflectance of forest canopy cover recorded by the instantaneous field of view (IFOV) of the sensor is a spectral mixture obtained from the interaction between electromagnetic radiation and both canopy and forest elements such as nonphotosynthetic vegetation (NPV; such as branches, stem, and litter), photosynthetic vegetation (PV; such as leaves), and others, such as bare soil and shadow. Therefore, image pixels are generally composed of more than one element, making the image interpretation difficult, which can result in a poor relation

One of the most widely used remote sensing approaches to derive and extract fraction covers from mixed pixels is the spectral mixture analysis (SMA) [92]. It decomposes the mixed pixel using a collection of constituent spectra (end members) to obtain their areal proportions or abundances in a pixel, and therefore unmixing a multispectral image into fraction images of end members [88,93]; it can be linear or nonlinear. The linear model assumes a single interac‐ tion between each incident photon and the surface object, and therefore the mixed pixel is a linear combination of pure spectral signatures (end members) of the surface materials weight‐ ed by their area covered. The nonlinear model is the opposite of linear model, since electro‐ magnetic radiation can intercept more than one element of surface, with mixed pixels resulting

As green leaves scatter radiation at NIR spectra, vertical structure of vegetation commonly produces a multiple-scattered signal; however, nonlinear spectral mixture approaches are barely used because they require more specific information, such as scattering properties of end members, the illumination of the sensor, and certain geometrical parameters of the scene

[91,93]. For these reasons, linear approaches have been largely implemented [92].

The linear mixture model is expressed in matrix form in the following equation (3):

*p Mf* = +

e

(3)

[76] and Landsat OLI [76], and SAR images as Jers-1 [76] and ALOS Palsar [32,75].

distance in a relative displacement vector between neighboring pixels.

**3.4. Biophysical variables**

12 Environmental Applications of Remote Sensing

chapter is focused on canopy fraction.

between AGB and spectral bands [88–91].

from a multiple-scattered signal [94].

Active sensors are those that provide its own energy source in order to control the double operation of signal emission and reception of known characteristics. These sensors have the advantage of an operational capacity to produce information both at night and in the day, in addition to working in a region of the electromagnetic spectrum that makes them less sensitive to atmospheric conditions. Of these, radar and LiDAR systems [55,56] are the most known. In this section, we briefly describe each of them, pointing out the relevant examples of their application in biomass modeling.

#### **3.6. Radio detection and ranging**

RAdio detection and ranging (RADAR) is the system name of active sensors that work in microwave region of the electromagnetic spectrum. Their mechanism is performed through signal transmission–reception of a portion of energy that interacts with the surface, which is referred as backscattered, being a measure of strength and time delay of the returned signals [105].

Such energy is considered consistent or coherent (illumination beam has same wavelength and phase), which makes it possible to use different polarization schemes (orientation of emitted and detected electromagnetic fields) to generate images [106]. The spatial resolution of radar images is strongly dependent on the antenna length (aperture) of the receptor and sensor inclination angle.

Synthetic aperture radar is widely used in forest monitoring through remote sensing, which is able to generate high-resolution imagery by taking advantage of the movement of the aircraft or satellite platform. SAR simulates a long virtual antenna that comprises long coherent successive radar signals, transmitted and received by a small antenna, which simultaneously moves along a given flight path [105,106].

Since microwave energy can penetrate forest canopies, the backscattered energy of SAR systems is modulated or influenced by the structural parameters of trees (e.g., branches, leaves, and stems), which in turn depends on different ecological variables [6,107,108]. Analyses of these data have been used to determine forest state [109], forest types [110], biomass density [32,111], forest canopy height [112], forest fire degradation [49,69], deforestation [113], and forest soil moisture [114].

The sensor sensitivity to forest parameters is a function of the wavelength, for instance, bands X and C with wavelengths 2.4 and 7.5 cm, respectively, are more sensitive to backscatter of leaves, whereas bands L and P with wavelengths 15 and 100 cm, respectively, are associated with backscatter of branches and stems [55,106,115].

Three of the main approaches from SAR systems that are widely used include (1) SAR backscattering coefficient [108,111,116], (2) interferometric SAR data [32,84], and (3) polari‐ metric SAR data [49,117].

#### **3.7. SAR backscattering coefficient**

Forest components scatter energy transmitted by the SAR systems in all directions. A portion of energy recorded by radar is translated to a proportional ratio between density of energy scattered and density of energy transmitted from the surface targets per unit area. Backscatter coefficient (*σ°*) or sigma nought is the amount of radar cross section [106,108]. Generally, this magnitude is expressed as a logarithm through decibel units as the following linear-form equation (4):

$$
\sigma^{\circ}\_{\,\,(db)} = 10^{\ast} \log\_{10} \sigma^{\circ} \tag{4}
$$

The backscatter coefficient value is related to two variables of sensor and target parameters. Sensor characteristics are a function of wavelength, polarization, and incidence angle, whereas target characteristics are associated with roughness, geometry, and dielectric properties [106, 108].

Biomass modeling through this approach has been usually applied from simple regression models under the assumption of correlation between backscatter coefficient and aboveground biomass/carbon [32,108,111,116,118]. The results are different because they rely on saturation of the signal, which is a function of wavelength, polarization, and the characteristics of the vegetation cover as well as of the difficulties caused by the specific properties of the ground as slope and aspect.

Recently, some of them have combined spatial models with remotely sensed data to improve geostatistical estimations using backscatter coefficient as spatial secondary variables [32].

#### **3.8. Interferometry SAR data**

Such energy is considered consistent or coherent (illumination beam has same wavelength and phase), which makes it possible to use different polarization schemes (orientation of emitted and detected electromagnetic fields) to generate images [106]. The spatial resolution of radar images is strongly dependent on the antenna length (aperture) of the receptor and sensor

Synthetic aperture radar is widely used in forest monitoring through remote sensing, which is able to generate high-resolution imagery by taking advantage of the movement of the aircraft or satellite platform. SAR simulates a long virtual antenna that comprises long coherent successive radar signals, transmitted and received by a small antenna, which simultaneously

Since microwave energy can penetrate forest canopies, the backscattered energy of SAR systems is modulated or influenced by the structural parameters of trees (e.g., branches, leaves, and stems), which in turn depends on different ecological variables [6,107,108]. Analyses of these data have been used to determine forest state [109], forest types [110], biomass density [32,111], forest canopy height [112], forest fire degradation [49,69], deforestation [113], and

The sensor sensitivity to forest parameters is a function of the wavelength, for instance, bands X and C with wavelengths 2.4 and 7.5 cm, respectively, are more sensitive to backscatter of leaves, whereas bands L and P with wavelengths 15 and 100 cm, respectively, are associated

Three of the main approaches from SAR systems that are widely used include (1) SAR backscattering coefficient [108,111,116], (2) interferometric SAR data [32,84], and (3) polari‐

Forest components scatter energy transmitted by the SAR systems in all directions. A portion of energy recorded by radar is translated to a proportional ratio between density of energy scattered and density of energy transmitted from the surface targets per unit area. Backscatter coefficient (*σ°*) or sigma nought is the amount of radar cross section [106,108]. Generally, this magnitude is expressed as a logarithm through decibel units as the following linear-form

> s

The backscatter coefficient value is related to two variables of sensor and target parameters. Sensor characteristics are a function of wavelength, polarization, and incidence angle, whereas target characteristics are associated with roughness, geometry, and dielectric properties [106,

Biomass modeling through this approach has been usually applied from simple regression models under the assumption of correlation between backscatter coefficient and aboveground

°= ° 10 10\*log *db* (4)

( ) s

inclination angle.

14 Environmental Applications of Remote Sensing

forest soil moisture [114].

metric SAR data [49,117].

equation (4):

108].

**3.7. SAR backscattering coefficient**

moves along a given flight path [105,106].

with backscatter of branches and stems [55,106,115].

Interferometric synthetic aperture radar (Figure 1) is a framework containing diverse methods or techniques that use phase information derived by recording phase difference or state of vibration of the wave at the instant that is received by the radar between two SAR images (known as master and slave) acquired from different sensor positions [106], called the interferometric phase (Δ*Φ*Int). The interferometric phase can be written as equation (5):

**Figure 1.** Interferometric synthetic aperture radar.

$$
\Delta\Phi\text{INT} = \Phi\text{S} - \Phi M = \Phi\text{Topo} + \Phi\text{Mov} + \Phi\text{Atm} + \Phi\text{Noise} \tag{5}
$$

where *ΦS* and *ΦM* are the phase observations of slave image and phase observations of master image, respectively, *Φ*Topo is the topographic component, *Φ*Mov is the shift component, *Φ*Atm is the atmospheric component, and *Φ*Noise is the phase noise.

One parameter obtained by this approach is the coherence or correlation image as an indicator of the phase stability [105]. Its mathematical expression is (6)

$$\gamma = \frac{\left(\mathbf{S}\_1 \mathbf{S}\_2^\*\right)}{\sqrt{\left(\mathbf{S}\_1 \mathbf{S}\_1^\*\right) \left(\mathbf{S}\_2 \mathbf{S}\_2^\*\right)}}\tag{6}$$

where *γ* is the coherence, the brackets *()* denote the ensemble average, and *\** denotes the complex conjugate; *S*<sup>1</sup> denotes the complex value of the master single-look complex (SLC) for pixels *x*,*y* and *S*<sup>2</sup> \* indicates the complex conjugate value for pixels *x*,*y* in the second image known as slave. Coherence values are included between 0 and 1 and a coherence image quantifies the decorrelation between two SLC images as the loss of coherence [106]. Decorre‐ lation is the combination of impacts in the radar phase: (1) the baseline decorrelation due to changes in the acquisition geometry of the images (which increases as the distance between the satellite orbits at the moment of acquisition increases) and (2) the temporal decorrelation due to variations in reflectivity in the zone, which can be caused by rain, phonological changes, and agricultural factors [105].

The interferometric coherence in biomass modeling is used under the assumption that for forested areas, coherence diminishes with the increment in vegetation density, as the volu‐ metric scattering increases with movement (wind) and forest growth. Biomass modeling through this approach has been typically used from simple regression models under the assumption of correlation between interferometric coherence and aboveground biomass/ carbon [32,84,116], another approach is by combining methods such as regression-kriging [32] and classification algorithms [119]. In this case, results are related to baseline and temporal decorrelations, forest type, polarization, and sensor wavelength.

#### **3.9. Polarimetric SAR data**

Antennas of radar systems can be configured to transmit and receive electromagnetic radiation polarized either horizontally or vertically. The two most common bases of polarizations are horizontal linear or H, and vertical linear or V. When the energy transmitted is polarized in the same direction as the received, it is called as like-polarized and when the transmitted energy depolarizes in a direction orthogonal to the received system, it is known as crosspolarization [120]. The polarization schemes are HH (for horizontal transmit and horizontal receive), VV (for vertical transmit and vertical receive), HV (for horizontal transmit and vertical receive), and VH (for vertical transmit and horizontal receive).

A radar system can be configured in different levels of polarization complexity:

**•** Single polarized – HH, VV, HV, or VH


where *ΦS* and *ΦM* are the phase observations of slave image and phase observations of master image, respectively, *Φ*Topo is the topographic component, *Φ*Mov is the shift component,

One parameter obtained by this approach is the coherence or correlation image as an indicator

( )

*S S*

\* 1 2 \* \* 11 22

*SS SS* (6)

( )( )

where *γ* is the coherence, the brackets *()* denote the ensemble average, and *\** denotes the complex conjugate; *S*<sup>1</sup> denotes the complex value of the master single-look complex (SLC) for

known as slave. Coherence values are included between 0 and 1 and a coherence image quantifies the decorrelation between two SLC images as the loss of coherence [106]. Decorre‐ lation is the combination of impacts in the radar phase: (1) the baseline decorrelation due to changes in the acquisition geometry of the images (which increases as the distance between the satellite orbits at the moment of acquisition increases) and (2) the temporal decorrelation due to variations in reflectivity in the zone, which can be caused by rain, phonological changes,

The interferometric coherence in biomass modeling is used under the assumption that for forested areas, coherence diminishes with the increment in vegetation density, as the volu‐ metric scattering increases with movement (wind) and forest growth. Biomass modeling through this approach has been typically used from simple regression models under the assumption of correlation between interferometric coherence and aboveground biomass/ carbon [32,84,116], another approach is by combining methods such as regression-kriging [32] and classification algorithms [119]. In this case, results are related to baseline and temporal

Antennas of radar systems can be configured to transmit and receive electromagnetic radiation polarized either horizontally or vertically. The two most common bases of polarizations are horizontal linear or H, and vertical linear or V. When the energy transmitted is polarized in the same direction as the received, it is called as like-polarized and when the transmitted energy depolarizes in a direction orthogonal to the received system, it is known as crosspolarization [120]. The polarization schemes are HH (for horizontal transmit and horizontal receive), VV (for vertical transmit and vertical receive), HV (for horizontal transmit and vertical

indicates the complex conjugate value for pixels *x*,*y* in the second image

*Φ*Atm is the atmospheric component, and *Φ*Noise is the phase noise.

g=

decorrelations, forest type, polarization, and sensor wavelength.

receive), and VH (for vertical transmit and horizontal receive).

**•** Single polarized – HH, VV, HV, or VH

A radar system can be configured in different levels of polarization complexity:

of the phase stability [105]. Its mathematical expression is (6)

pixels *x*,*y* and *S*<sup>2</sup>

\*

16 Environmental Applications of Remote Sensing

and agricultural factors [105].

**3.9. Polarimetric SAR data**

A quadrature polarization or polarimetric radar uses these four polarizations in order to measure the magnitudes and relative phase difference between the polarization schemes or channels through an ellipse shape [106,120]. These kinds of radar systems promoted a new framework called polarimetry of synthetic aperture radar, which describes the surface through different combinations of polarization under the assumption that the interaction of electro‐ magnetic energy and elements of the land surface can change the polarization of a portion of the wavelength transmitted by the sensor, and therefore receive information of amplitude and relative phase of the same target in four channels of information, which is considered as a basis for description of scattering polarimetric of surface targets [106,120]. It is mathematically simplified in the so-called scattering matrix (equation (7)):

$$\mathbf{S} = \begin{bmatrix} \mathbf{S}\_{hh} & \mathbf{S}\_{hv} \\ \mathbf{S}\_{vh} & \mathbf{S}\_{vv} \end{bmatrix} \tag{7}$$

which describes different forms of the polarized electric fields between incident wave and the scatter wave in order to be the basis for diverse ways to analyze the scattering properties of a target (e.g., the covariance and coherency matrices) and diverse transformations as polarimet‐ ric decompositions and through the synergy between polarimetry and interferometry (polari‐ metric interferometric coherence) [120].

The use of polarimetry in biomass modeling is under the assumption of correlation between forest structural properties and polarimetric behavior. It may be construed through scattering mechanism analysis. This has been addressed mainly by polarimetric decompositions, such as Freeman Durden [49,121], eigenvector–eigenvalue [49], and Cloude and Pottier [31,121]. Biomass modeling through this framework has been usually performed from simple and multiple regression models [31,121] and nonparametric model random forest regression [49]. In this framework, results are related to forest type, spatial resolution, and sensor wavelength. 3.10. Light detection and ranging

Light detection and ranging is an active laser sensor, which emits pulses of polarized light or pulse echo, which can be calibrated within a narrow range of wavelength. The most commonly used wavelength is 1,064 nm (near-infrared), although it can range from ultraviolet to infrared range of electromagnetic spectrum (500–1500 nm) [122–124].

These laser scanners consist of a range finder, global positioning system (GPS), inertial measuring unit (IMU), and a clock capable of recording travel times to within 0.2 of a nano‐ second. The integration of these systems produce accurate measurements of the position and orientation of objects registered. These technologies allow us to measure elapsed time of pulse echo between laser transmitter and objects on the surface. The energy that interacts with surfaces is backscattered over different times exhibiting multiple laser pulse returns associated with distinct surface layers toward a laser scanner that can be mounted on an aerial or terrestrial platform [123,125,126].

LiDAR information is essentially a three-dimensional point cloud composed of simple derivative returns and multiple laser pulses, this type of LiDAR data is known as discrete LiDAR returns. In addition to three-dimensional information, most LiDAR systems record the intensity as a fraction of pulse energy reflected at that location [127].

The use of LiDAR in forest areas is mainly to analyze forest vertical structural metrics under the assumption that laser can be sensed from the top of the canopies, elements of different canopies, or even to the ground, which will be reflected in the number of returns. The depth of laser penetration depends on the density of canopies and density of point clouds, which vary from less than one point per square meter to several dozens, with vertical accuracies around 12.5 cm [123,127].

One of the most widely used products in forest analysis is the result of the processing of threedimensional point cloud, the canopy height model (CHM) (Figure 2). It is derived from the difference or subtraction between digital elevation models and digital terrain models, both datasets generally are a result of different interpolation methods, such as nearest neighbor, splines, inverse distance weighting, and kriging [127]. The first one is associated with first returns and the second one is related to the last returns. Other forest structural measurements are the fractional crown cover, crown area, crown diameter, basal area [38,125], and canopy volume [126], which are of key interest to the managers and represent information that is expensive and time consuming to collect in the field. When small-footprint LiDAR data are acquired at very high enough densities, individual tree crowns can readily be observed in the point clouds, processing algorithm for automated measuring and modeling of vegetation at individual tree crown segmentation (e.g., watershed segmentation) [95,128,129].

Biomass modeling through this approach has been usually used in simple and multiple regression models [38,125,126,130]. Other works have explored the use of learning machines [131,132]. Other approach is through combining methods that integrate LiDAR information with other sensors [22,31,124,133].

#### **4. Case study**

Mexico City is a continuous truss of multiple ecosystems, which is administratively divided into two large areas: urban land (41%) governed by the Urban Development Programs and Environmental Conservation Zone (ECZ; 61%) steered by the General Ecological Planning Program (Figure 3). The ECZ provides Mexico City with environmental services such as carbon capture, aquifer recharge, biodiversity, and scenic beauty. The zone is under anthropogenic pressure, including human settlements, land use changes, and extraction of natural resources, and therefore immediate action for conservation and appropriate resource management is necessary. This has led to deforestation, degradation, development of pest infestations, fires, and erosion. Models of the spatial variability of forest density are required in order to obtain

**Figure 2.** Canopy height model.

with distinct surface layers toward a laser scanner that can be mounted on an aerial or

LiDAR information is essentially a three-dimensional point cloud composed of simple derivative returns and multiple laser pulses, this type of LiDAR data is known as discrete LiDAR returns. In addition to three-dimensional information, most LiDAR systems record the

The use of LiDAR in forest areas is mainly to analyze forest vertical structural metrics under the assumption that laser can be sensed from the top of the canopies, elements of different canopies, or even to the ground, which will be reflected in the number of returns. The depth of laser penetration depends on the density of canopies and density of point clouds, which vary from less than one point per square meter to several dozens, with vertical accuracies

One of the most widely used products in forest analysis is the result of the processing of threedimensional point cloud, the canopy height model (CHM) (Figure 2). It is derived from the difference or subtraction between digital elevation models and digital terrain models, both datasets generally are a result of different interpolation methods, such as nearest neighbor, splines, inverse distance weighting, and kriging [127]. The first one is associated with first returns and the second one is related to the last returns. Other forest structural measurements are the fractional crown cover, crown area, crown diameter, basal area [38,125], and canopy volume [126], which are of key interest to the managers and represent information that is expensive and time consuming to collect in the field. When small-footprint LiDAR data are acquired at very high enough densities, individual tree crowns can readily be observed in the point clouds, processing algorithm for automated measuring and modeling of vegetation at

individual tree crown segmentation (e.g., watershed segmentation) [95,128,129].

Biomass modeling through this approach has been usually used in simple and multiple regression models [38,125,126,130]. Other works have explored the use of learning machines [131,132]. Other approach is through combining methods that integrate LiDAR information

Mexico City is a continuous truss of multiple ecosystems, which is administratively divided into two large areas: urban land (41%) governed by the Urban Development Programs and Environmental Conservation Zone (ECZ; 61%) steered by the General Ecological Planning Program (Figure 3). The ECZ provides Mexico City with environmental services such as carbon capture, aquifer recharge, biodiversity, and scenic beauty. The zone is under anthropogenic pressure, including human settlements, land use changes, and extraction of natural resources, and therefore immediate action for conservation and appropriate resource management is necessary. This has led to deforestation, degradation, development of pest infestations, fires, and erosion. Models of the spatial variability of forest density are required in order to obtain

intensity as a fraction of pulse energy reflected at that location [127].

terrestrial platform [123,125,126].

18 Environmental Applications of Remote Sensing

around 12.5 cm [123,127].

with other sensors [22,31,124,133].

**4. Case study**

an inventory of carbon storage, useful for public policies in areas with high environmental value in order to facilitate decision-making by reducing the complexity involved in integrating and interpreting values at a pixel level.

The study area lies within the ECZ of Mexico City (882 km2 ) and is covered by sacred fir or oyamel (*Abies religiosa*) and Mexican mountain pine (*Pinus hartwegii*) forests. Fir forests generally grow at 2,700–3,500 m above sea level. They are evergreen forests with heights of 20–40 m and their understory is densely shaded. This type of forest contains important elements, including alders (*Alnusfir mifolia*), white cedar (*Cupressus lindleyi*), oak (*Quercus laurina*), Mexican Douglas-fir (*Pseudotsuga macrolepis*), willows (*Salix oxylepis*), and black cherry (*Prunus serotina*). Pine forests (*Pinus hartwegii*) grow above 3,000 m, this species being most tolerant to the extreme environmental conditions of the high mountains. This pine develops along with Festuca and Muhlenbergia grasses [134].

**Figure 3.** Environmental Conservation Zone.

The present case study is a comparative analysis between regression-cokriging and multiple regression approaches using satellite-derived indices for modeling the aboveground biomass of forests in the Environmental Conservation Zone of Mexico City. The objectives are:


In order to achieve these objectives, a correlation analysis was performed between digital canopy height model (LiDAR data) and vegetation fraction cover (SPOT-5 data) and, on the other hand, ground biomass estimates at forest inventory sites. Then, the spatial autocorrela‐ tion was calculated for residuals in order to define the variables to be used in multiple models and regression-cokriging methods. Once models were obtained, the root mean square error (RMSE) was computed for each approach.

#### **4.1. Methods**

#### *4.1.1. Models of aboveground biomass carbon used three sources of data: Forest inventory data from in situ measurements*

Since aboveground carbon is the amount of carbon stored in aboveground biomass, compris‐ ing all living plant material above the soil (e.g., trunks, branches, and leaves) [3], the calculation of carbon stock from biomass consists of multiplying the total biomass by a conversion factor that represents the average carbon content in biomass. A common assumption is that biomass is around 50% carbon expressed in tons of dry matter per unit area [135]. Typically, the terms of measurement are density of biomass expressed as mass per unit area (e.g., t/ha). Here, forest inventory data were obtained from Mexico City Environmental and Land Planning Authority (PAOT, because of its name in Spanish) and were derived from sampling 283 plots during 2008–2010. Their sampling is based on the design of the National Forest and Soil Inventory of the National Forest Commission (CONAFOR, because of its name in Spanish). In it, each sampling conglomerate is composed of four circular secondary sampling plots in an inverted "Y" shape, each of which covers an area of 400 m2 and peripheral plots are at 45.15 m from the center of the conglomerate (Figure 4). Of the 283 plots, 155 were among pines, 86 in fir forest, 30 in mixed forest, 10 in scrub, and 2 in planted forests. Per-tree carbon was estimated from allometric carbon equations developed by Acosta-Mireles et al., Jiménez, and Avendaño-Hernández et al. for the species of the region [136–138]. Conversion of biomass carbon from conglomerate to hectares [139] used the "ratio of means" as shown in equation (8):

$$
\hat{R} = \frac{\overline{Y}}{\overline{X}} = \frac{\sum\_{i=1}^{n} Y\_i}{\sum\_{i=1}^{n} X\_i} \tag{8}
$$

where *Yi* is the total aboveground carbon in all plots of 400 m2 and *Xi* is the total area sampled in *i* plots.

**Figure 4.** Sampling plots and conglomerate form.

The present case study is a comparative analysis between regression-cokriging and multiple regression approaches using satellite-derived indices for modeling the aboveground biomass of forests in the Environmental Conservation Zone of Mexico City. The objectives are:

**1.** to identify satellite-derived indices that are better associated with aboveground biomass,

**2.** to quantify spatial patterns of the residuals derived from simple regression between

**3.** to determine whether spatial statistical methods improve the estimates of aboveground

In order to achieve these objectives, a correlation analysis was performed between digital canopy height model (LiDAR data) and vegetation fraction cover (SPOT-5 data) and, on the other hand, ground biomass estimates at forest inventory sites. Then, the spatial autocorrela‐ tion was calculated for residuals in order to define the variables to be used in multiple models and regression-cokriging methods. Once models were obtained, the root mean square error

*4.1.1. Models of aboveground biomass carbon used three sources of data: Forest inventory data from in*

Since aboveground carbon is the amount of carbon stored in aboveground biomass, compris‐ ing all living plant material above the soil (e.g., trunks, branches, and leaves) [3], the calculation of carbon stock from biomass consists of multiplying the total biomass by a conversion factor

either from LiDAR or fraction cover (SPOT-5) imagery

(RMSE) was computed for each approach.

**Figure 3.** Environmental Conservation Zone.

20 Environmental Applications of Remote Sensing

**4.1. Methods**

*situ measurements*

satellite indices and carbon values, using spatial autocorrelation

biomass carbon pools over nonspatial conventional regression methods

#### *4.1.1.1. SPOT image*

Four multispectral SPOT-5 HRG images were used: two from February 25, 2010 (zenith 51.72° and 52.03°, azimuth 136.96° and 136.43°) and two from March 28, 2010 (zenith 62.67° and 62.89°, azimuth 125.66° and 124.75°). These were radiometric and atmospherically corrected with the second simulation of a satellite signal in the solar spectrum (S-6) code through CLASlite software and orthorectified with the polynomial coefficients and geoid information based on Geocover 2000 of Landsat as reference images.

#### *4.1.2. LiDAR data*

LiDAR data used were acquired from the ALS50-II sensor flown by the National Institute for Statistics and Geography (INEGI, because of its name in Spanish) between November and December of 2007 over the entire Mexico Valley. The data had an average horizontal distance of 2.0 m, minimum point density of 0.433 points/m2 , and vertical root mean square error of 7.3 cm. These points are used as the basis for the generation of digital terrain model (DTM) and digital surface model (DSM) with a resolution of 5 m.

#### *4.1.2.1. Photosynthetic vegetation fraction cover*

Photosynthetic vegetation fraction cover was estimated throughout the Automated Monte Carlo Unmixing (AutoMCU) model. This model integrates spectral mixture analysis and spectral end-member libraries resulting from fieldwork (ground spectrometer) and highresolution hyperspectral information of Hyperion Sensor, in order to separate photosynthetic vegetation, non-photosynthetic vegetation, and bare substrate. The photosynthetic vegetation fraction cover was calculated by CLASlite v3.2 software (Figure 5) [140].

#### *4.1.2.2. Canopy height model*

The calculation of canopy height model used altitude values of different digital terrain model (DTM) and digital surface model (DSM) in order to extract differences between both models (Figure 6).

#### *4.1.2.3. Correlation and autocorrelation coefficients*

Correlation analysis based on multiple regressions explored statistical relationships between aboveground biomass carbon and satellite indices. The sampling points were randomly divided into 50% for model calibration and 50% for model verification. Residuals from regressions were retained and their spatial autocorrelation was quantified [141]. Moran's I index was used to identify the type and intensity of spatial pattern, measuring the degree of autocorrelation or dependence of a distribution. Moran's I index can be written as in the equation (9):

$$I = \frac{n}{\sum\_{i} \sum\_{j} \sigma\_{ij}} \ast \frac{\sum\_{i} \sum\_{j} \sigma\_{ij} \left(\mathbf{x}\_{i} - \overline{\mathbf{x}}\right) \left(\mathbf{x}\_{j} - \overline{\mathbf{x}}\right)}{\sum\_{i} \left(\mathbf{x}\_{i} - \overline{\mathbf{x}}\right)^{2}} \tag{9}$$

**Figure 5.** Photosynthetic vegetation fraction cover.

where *n* is the number of spatial units indexed by *i* and *j*, *x* is the variable of interest, and *ϖij* is the spatial weight matrix. Moran's I values near to 1 indicate clustering negative values near to −1 represent spatial dispersion, and a value of zero indicates randomness. Statistical significance was expressed in terms of the Z descriptor and confidence level 1−α. Construction of the spatial weight matrix was distance based, since the spatial representation units are points. This distance or spatial lag includes at least 12 samples as recommended by Isaaks and Srivsatava [33].

#### *4.1.3. Regression models*

*4.1.1.1. SPOT image*

22 Environmental Applications of Remote Sensing

*4.1.2. LiDAR data*

Geocover 2000 of Landsat as reference images.

of 2.0 m, minimum point density of 0.433 points/m2

*4.1.2.1. Photosynthetic vegetation fraction cover*

*4.1.2.3. Correlation and autocorrelation coefficients*

*4.1.2.2. Canopy height model*

(Figure 6).

equation (9):

digital surface model (DSM) with a resolution of 5 m.

Four multispectral SPOT-5 HRG images were used: two from February 25, 2010 (zenith 51.72° and 52.03°, azimuth 136.96° and 136.43°) and two from March 28, 2010 (zenith 62.67° and 62.89°, azimuth 125.66° and 124.75°). These were radiometric and atmospherically corrected with the second simulation of a satellite signal in the solar spectrum (S-6) code through CLASlite software and orthorectified with the polynomial coefficients and geoid information based on

LiDAR data used were acquired from the ALS50-II sensor flown by the National Institute for Statistics and Geography (INEGI, because of its name in Spanish) between November and December of 2007 over the entire Mexico Valley. The data had an average horizontal distance

cm. These points are used as the basis for the generation of digital terrain model (DTM) and

Photosynthetic vegetation fraction cover was estimated throughout the Automated Monte Carlo Unmixing (AutoMCU) model. This model integrates spectral mixture analysis and spectral end-member libraries resulting from fieldwork (ground spectrometer) and highresolution hyperspectral information of Hyperion Sensor, in order to separate photosynthetic vegetation, non-photosynthetic vegetation, and bare substrate. The photosynthetic vegetation

The calculation of canopy height model used altitude values of different digital terrain model (DTM) and digital surface model (DSM) in order to extract differences between both models

Correlation analysis based on multiple regressions explored statistical relationships between aboveground biomass carbon and satellite indices. The sampling points were randomly divided into 50% for model calibration and 50% for model verification. Residuals from regressions were retained and their spatial autocorrelation was quantified [141]. Moran's I index was used to identify the type and intensity of spatial pattern, measuring the degree of autocorrelation or dependence of a distribution. Moran's I index can be written as in the

v

å å å å <sup>å</sup> <sup>2</sup> \* *ij i <sup>j</sup> i j ij i j <sup>i</sup> <sup>i</sup>*

*x xx x <sup>n</sup> <sup>I</sup>*

v

=

( )( )


( )

*x x*


fraction cover was calculated by CLASlite v3.2 software (Figure 5) [140].

, and vertical root mean square error of 7.3

(9)

In addition to multiple regression, the present study compared models derived from regres‐ sion-cokriging. The regression-cokriging was calculated through the estimation of a simple linear regression approach between aboveground carbon and canopy height model (equation (10)) and the addition of interpolated layer via ordinary cokriging integrated by regression residuals and the secondary variable (photosynthetic vegetation fraction cover) [31, 142]. The predictions were carried out separately for drift and residuals, and were added together later as in the following equation (11):

$$
\hat{z}RCoK\left(\mathcal{S}\_{0}\right) = \hat{m}\left(\mathcal{S}\_{0}\right) + \hat{e}\left(\mathcal{S}\_{0}\right) \tag{10}
$$

**Figure 6.** Digital canopy height model.

$$\hat{\mathbf{z}}\mathbf{R}\mathbf{C}\mathbf{o}\mathbf{K}\{\mathbf{S}\_{0}\} = \sum\_{k=0}^{p} \hat{\boldsymbol{\beta}}\_{k}^{\ast} \cdot \boldsymbol{q}\_{k}\{\mathbf{S}\_{0}\} + \left(\sum\_{l=1}^{n} \boldsymbol{\sigma}\_{l} \boldsymbol{e}\{\mathbf{S}\boldsymbol{o}\} + \sum\_{j=1}^{n} \boldsymbol{\sigma}\_{2j} \boldsymbol{Z}\_{2}\{\mathbf{S}\_{i}\}\right) \tag{11}$$

where *βk* are the coefficients of the drift model, *qk* is the number of auxiliary variables, *ϖi(S*0*)* and *ϖ*<sup>2</sup> *<sup>j</sup>* are the weights determined by covariogram for regression residuals, and secondary variable and *e* are the regression residuals [31].

#### *4.1.3.1. Model accuracy*

The regression models were validated with data from the field sampling [143]. The RMSE criterion was used to determine which regression models have more precision in the estimation of stored carbon in the area, the RMSE can be written as equation (12):

$$\text{RMSE} = \sqrt{\sum\_{i=1}^{n} \frac{\left(Z\_{(si)} - z\_{(si)}\right)^2}{n}} \tag{12}$$

where *Z(si)* is the reference value, *z(si)* is the estimated value, and *n* is the total number of samples.

#### **4.2. Results**

#### *4.2.1. Correlation*

The degree of association between carbon stored and each index derived from remote sensing and multiple associations (Table 2) was the synergy between canopy height model and photosynthetic vegetation fraction cover with an *r* coefficient of 0.88. All these correlations were positive, indicating that, as stored carbon increases, vegetation indices increase.


**Table 2.** Correlation coefficients above-ground carbon and remote sensing indices.

#### *4.2.2. Moran's I index*

( )

variable and *e* are the regression residuals [31].

*4.1.3.1. Model accuracy*

**Figure 6.** Digital canopy height model.

24 Environmental Applications of Remote Sensing

samples.

ˆ ˆ \*

b

of stored carbon in the area, the RMSE can be written as equation (12):

RMSE

where *βk* are the coefficients of the drift model, *qk* is the number of auxiliary variables, *ϖi(S*0*)* and *ϖ*<sup>2</sup> *<sup>j</sup>* are the weights determined by covariogram for regression residuals, and secondary

The regression models were validated with data from the field sampling [143]. The RMSE criterion was used to determine which regression models have more precision in the estimation

( ( ) ( ) )

*<sup>n</sup> si si*

where *Z(si)* is the reference value, *z(si)* is the estimated value, and *n* is the total number of

*Z z n*


2

(12)

=

*i*

1

<sup>=</sup> å

= ++ ç ÷ è ø <sup>0</sup> å åå <sup>0</sup> 2 2 0 11

*k k i j i*

= ==

*k ij*

*p n n*

vv( ) ( ) ( )

*zRCoK S q S e So Z S* (11)

æ ö

Once the satellite-derived index, that is the most associated with carbon storage, was identified, Moran's I was calculated from regression residuals according to sampling sites in order to identify spatial autocorrelation and, hence, the information could be included as auxiliary information in the regression-cokriging model. A value of 0.31 (*z* = 2.96, α = 0.01) was obtained for the regression residuals between canopy height model and aboveground carbon. In this case, a low-positive spatial autocorrelation was present, with statistical significance, indicating that the neighboring spatial units presented near or close values and a slight trend toward plots.

#### *4.2.3. Spatial distribution*

To identify the spatial distribution of stored carbon, models were developed with the appli‐ cation of equations resulting from multiple regression and regression-cokriging spatial methods. By obtaining the covariograms for the theoretical fit, no significant variation was found between anisotropic and isotropic covariograms, therefore, the isotropic model was settled.

As the empirical covariogram showed strong spatial dependence, it did not present constant semivariances as a function of distance. In this case, the adjusted theoretical models ranged from 10,000 to 15,000 m, distances at which the observations were independent.

We used the best fit (according to root mean squares in prediction errors) for the ordinary cokriging interpolation to evaluate which is more effective in the predictions throughout the study area. Table 3 shows the parameters obtained.


**Table 3.** Covariogram indices

#### *4.2.4. Model accuracy*

Comparison of the models with the set of verification sites produced the RMSE in tC/ha (Table 4). The models based on regression-cokriging presented the least error. Figure 5 shows the spatial distribution obtained by regression-cokriging and multiple regressions for stored carbon. The delineation of forest types (fir and pine) was based on the map of land use and vegetation of PAOT [144].


**Table 4.** RMSE

**Figure 7.** Spatial distribution of the aboveground carbon storage in Environmental Conservation Zone, estimated by two regression models.

The results of this research indicated that consideration of spatial autocorrelation can improve estimates of carbon content in aboveground biomass, specifically using the regression cokriging method. This could be due to its sensitivity to local variations [142], since it is particularly developed to consider the adjustment of the spatial variance model in order to improve predictions obtained from global models. The present work uses spatial modeling of indices related to carbon for the purpose of exploring the spatial autocorrelation of auxiliary variables, which is useful for representing the way in which a phenomenon radiates through spatial units. Although the primary method used for the estimation of carbon stocks in Mexico is the stratify-and-multiply approach, which assigns a single value or a range of values to each vegetation type and then multiplies these values by the areas covered by the vegetation to estimate total carbon stock values, this investigation has demonstrated a more accurate, spatial-explicit, repeatable alternative.

According to Ref. [145], autocorrelation is "perhaps, after the average and variance, the most important property of any geographic variable and, unlike these, it is explicitly linked with spatial patterns." The present comparative analysis demonstrates the importance of the use of spatial methods to model carbon stored in the aboveground biomass, since these methods consider the spatial pattern of the data. The hypothesis of the homogeneity of the relationships between stored carbon and remote sensing indices sometimes does not consider the spatial heterogeneity of many factors affecting this relationship, such as geographic differences in orientation and climatic and soil conditions [28].

This analysis provided a synoptic mapping of aboveground biomass as a potentially valuable tool for environmental protection policies in the ECZ of Mexico City, one of the most important ecological reserves for the inhabitants of the Mexico Valley in the economic, cultural, and social sense, as well as for the volume and quality of the environmental services it provides.

### **5. Conclusion**

**Remote Sensing indices model Sill Range Nugget**

Comparison of the models with the set of verification sites produced the RMSE in tC/ha (Table 4). The models based on regression-cokriging presented the least error. Figure 5 shows the spatial distribution obtained by regression-cokriging and multiple regressions for stored carbon. The delineation of forest types (fir and pine) was based on the map of land use and

**Figure 7.** Spatial distribution of the aboveground carbon storage in Environmental Conservation Zone, estimated by

**Models RMSE (tC/ha)**

Multiple regression 34.1 Regression-Cokriging 28.6

Exponential 0.89 12,007 0

Photosynthetic Vegetation Fraction Cover and

26 Environmental Applications of Remote Sensing

regression residuals

**Table 3.** Covariogram indices

vegetation of PAOT [144].

**Table 4.** RMSE

two regression models.

*4.2.4. Model accuracy*

Remote sensing-derived indices play a major role in forest monitoring, because traditional methodologies derive their estimates of carbon content in the biomass through forest inven‐ tories and, for its implementation, they require much time and money and are generally limited to 10-year intervals. The information resulting from them is designed to present average timber volumes linked to administrative regions, which not represent the spatial variability and therefore it generate a bias in carbon measures..

This has led to great interest to estimate, map, and monitor the carbon stored in forests more precisely, enhancing the recognition of their role in the global carbon cycle, particularly in the mitigation of greenhouse gases. Through the estimation of carbon content, a base line for calculating the dynamics of this gas as a mitigation strategy can be established.

An overview of modeling options and remote sensing resources that have been used for monitoring and researching the forest biomass was presented. Techniques ranging from collecting georeferenced data in the field to the information extraction methods from satellite images and synergies between remote sensing and geostatistics were described. A case study was selected to illustrate the application of some of these techniques in the modeling of spatial distribution of aboveground carbon in Mexican coniferous forests.

Based on the study, the following conclusions can be drawn:


### **Author details**

José Mauricio Galeana Pizaña\* , Juan Manuel Núñez Hernández and Nirani Corona Romero

\*Address all correspondence to: mgaleana@centrogeo.edu.mx

Centro de Investigación en Geografía y Geomática "Ing. Jorge L. Tamayo", México DF, Mexico

#### **References**


[5] Zhu X, Liu D. Improving forest aboveground biomass estimation using seasonal Landsat NDVI time-series. ISPRS Journal of Photogrammetry and Remote Sensing. 2015;102:222–231. DOI: 10.1016/j.isprsjprs.2014.08.014

images and synergies between remote sensing and geostatistics were described. A case study was selected to illustrate the application of some of these techniques in the modeling of spatial

**1.** According to these results, the synergy between remote sensing and geostatistics has the potential to estimate forest biomass to improve estimations using remote sensing indices

**2.** Geospatial methods have a better modeling adjustment (e.g., RMSE) than conventional statistical methods as multiple regressions, because geospatial methods considerer local

**3.** Best Pearson coefficient between the two variables tested in the study is the digital canopy model, which resulted from LiDAR data. This kind of information is very expensive, so integrating multispectral information can be a way to capitalize on multitemporal study

Centro de Investigación en Geografía y Geomática "Ing. Jorge L. Tamayo", México DF,

[1] FAO. Global Forest Resources Assessment 2010: Main Report. Roma: Food and Agri‐

[2] Pavlic G, Chen W, Fernandes R, Cihlar J, Price DT, Latifovic R, Fraser R, Leblanc SG. Canada-wide maps of dominant tree species from remotely sensed and ground data.

[3] Brown S. Measuring carbon in forests: current status and future challenges. Environ‐

[4] Wulder M. Optical remote-sensing techniques for the assessment of forest inventory and biophysical parameters. Progress in Physical Geography. 1998;22:449–476. DOI:

Geocarto International. 2007;22:185–204. DOI: 10.1080/10106040701201798

mental Pollution. 2002;116:363–372. DOI: 10.1016/S0269-7491(01)00212-3

, Juan Manuel Núñez Hernández and Nirani Corona Romero

distribution of aboveground carbon in Mexican coniferous forests.

Based on the study, the following conclusions can be drawn:

\*Address all correspondence to: mgaleana@centrogeo.edu.mx

culture Organization of the United Nations. 2010. 340 p.

as spatial secondary variables.

28 Environmental Applications of Remote Sensing

spatial variations.

of biomass.

**Author details**

Mexico

**References**

José Mauricio Galeana Pizaña\*

10.1177/030913339802200402


[28] Maselli F, Chiesi M. Evaluation of statistical methods to estimate forest volume in a Mediterranean Region. IEEE Transaction on Geoscience and Remote Sensing. 2006;44:2239–2250. DOI: 10.1109/TGRS.2006.872074

[17] Kellndorfer JM, Walker WS, LaPoint E, Kirsch K, Bishop J, Fiske G. Statistical fusion of Lidar, InSAR, and optical remote sensing data for forest stand height characteriza‐ tion: a regional‐scale method based on LVIS, SRTM, Landsat ETM+, and ancillary da‐ ta sets. Journal of Geophysical Research: Biogeosciences. 2010;115:1-10. DOI:

[18] McRoberts RE, Westfall JA. Effects of uncertainty in model predictions of individual tree volume on large area volume estimates. Forest Science. 2014;60:34–42. DOI:

[19] Seidel D, Fleck S, Leuschner C, Hammett T. Review of ground-based methods to measure the distribution of biomass in forest canopies. Annals of Forest Science.

[20] Segura M, Kanninen M. Allometric models for tree volume and total aboveground biomass in a tropical humid forest in Costa Rica. Biotropica. 2005;37:2–8. DOI:

[21] Mitchard ETA, Saatchi SS, Lewis SL, Feldpausch TR, Woodhouse IH, Sonké B, Row‐ land C, Meir P. Measuring biomass changes due to woody encroachment and defor‐ estation/degradation in a forest–savanna boundary region of central Africa using multi-temporal L-band radar backscatter. Remote Sensing of Environment.

[22] Sun G, Ranson KJ, Guo Z, Zhang Z, Montesano P, Kimes D. Forest biomass mapping from lidar and radar synergies. Remote Sensing of Environment. 2011;115:2906–2916.

[23] Lu D. Aboveground biomass estimation using Landsat TM data in the Brazilian Am‐ azon. International Journal of Remote Sensing. 2005;26:2509–2525. DOI:

[24] Lu D, Batistella M, Moran E. Satellite estimation of aboveground biomass and im‐ pacts of forest stand structure. Photogrammetric Engineering & Remote Sensing.

[25] Boyd DS, Danson FM. Satellite remote sensing of forest resources: three decades of research development. Progress in Physical Geography. 2005;29:1–26. DOI:

[26] Chen Q, Laurin GV, Battles JJ, Saah D. Integration of airborne lidar and vegetation types derived from aerial photography for mapping aboveground live biomass. Re‐ mote Sensing of Environment. 2012;121:108–117. DOI: 10.1016/j.rse.2012.01.021 [27] McRoberts RE, Næsset E, Gobakken T. Inference for lidar-assisted estimation of for‐ est growing stock volume. Remote Sensing of Environment. 2013;128:268–275. DOI:

10.1029/2009JG000997

30 Environmental Applications of Remote Sensing

http://dx.doi.org/10.5849/forsci.12-141

10.1111/j.1744-7429.2005.02027.x

DOI: 10.1016/j.rse.2011.03.021

10.1080/01431160500142145

10.1191/0309133305pp432ra

10.1016/j.rse.2012.10.007

2011;68:225–244. DOI: 10.1007/s13595-011-0040-z

2011;115:2861–2873. DOI: 10.1016/j.rse.2010.02.022

2005;71:967–974. DOI: 10.14358/PERS.71.8.967


[51] Mountrakis G, Im J, Ogole C. Support vector machines in remote sensing: a review. ISPRS Journal of Photogrammetry and Remote Sensing. 2011;66:247–259. DOI: 10.1016/j.isprsjprs.2010.11.001

[40] McRoberts RE, GobakkenT, Næsset E. Post-stratified estimation of forest area and growing stock volume using lidar-based stratifications. Remote Sensing of Environ‐

[41] Saatchi S, Malhi Y, Zutta B, Buermann W, Anderson LO, Araujo AM, Phillips OL, Peacock J, ter Steege H, Lopez Gonzalez G, Baker T, Arroyo L, Almeida S, Higuchi N, Killen T, Monteagudo A, Neill D, Pitman N, Prieto A, Salomão R, Silva N, Vásquez Martínez R, Laurence W, Ramírez HA. Mapping landscape scale variations of forest structure, biomass, and productivity in Amazonia. Biogeosciences Discussions.

[42] Kilkki P, Päivinen R. Reference sample plots to combine field measurements and sat‐ ellite data in forest inventory. Department of Forest Mensuration and Management,

[43] Tomppo EO, Gagliano C, De Natale F, Katila M, McRoberts RE. Predicting categori‐ cal forest variables using an improved k-Nearest Neighbour estimator and Landsat imagery. Remote Sensing of Environment. 2009;113:500–517. DOI: 10.1016/j.rse.

[44] Latifi H, Nothdurft A, Koch B. Non-parametric prediction and mapping of standing timber volume and biomass in a temperate forest: application of multiple optical/ LiDAR-derived predictors. Forestry. 2010;83:395–407. DOI: 10.1093/forestry/cpq022

[45] Tomppo E, Olsson H, Ståhl G, Nilsson M, Hagner O, Katila M. Combining national forest inventory field plots and remote sensing data for forest databases. Remote

[46] Foody GM, Cutler ME, Mcmorrow J, Pelz D, Tangki H, Boyd DS, Douglas I. Map‐ ping the biomass of Bornean tropical rain forest from remotely sensed data. Global

[47] Mas JF, Flores JJ. The application of artificial neural networks to the analysis of re‐ motely sensed data. International Journal of Remote Sensing. 2008;29:617–663. DOI:

[48] Pflugmacher D, Cohen WB, Kennedy RE, Yang Z. Using Landsat-derived disturb‐ ance and recovery history and lidar to map forest biomass dynamics. Remote Sens‐

[49] Tanase M, Panciera R, Lowell K, Tian S, Kacker J, Walker J. Airborne multi-temporal L-band polarimetric SAR data for biomass estimation in semi-arid forests. Remote

[50] Marabel M, Alvarez-Taboada F. Spectroscopic determination of aboveground bio‐ mass in grasslands using spectral transformations, support vector machine and parti‐ al least squares regression. Sensors. 2013;13:10027–10051. DOI: 10.3390/s130810027

ing of Environment. 2014;151:124–137. DOI: 10.1016/j.rse.2013.05.033

Sensing of Environment. 2014;145:93–104. DOI: 10.1016/j.rse.2014.01.024

Sensing of Environment. 2008;112:1982–1999. DOI: 10.1016/j.rse.2007.03.032

ment. 2012;125:157–166. DOI: 10.1016/j.rse.2012.07.002

2009;6:5461–5505. DOI: 10.5194/bgd-6-5461-2009

Ecology and Biogeography. 2001;10:379–387.

10.1080/01431160701352154

2008.05.021.

32 Environmental Applications of Remote Sensing

University of Helsinki, Research Notes. 1987;19:210–215.


tions between regions. ISPRS Journal of Photogrammetry and Remote Sensing. 2012;70:66–77. DOI: 10.1016/j.isprsjprs.2012.03.011

[77] Dube T, Mutanga O. Investigating the robustness of the new Landsat-8 Operational Land Imager derived texture metrics in estimating plantation forest aboveground bi‐ omass in resource constrained areas. ISPRS Journal of Photogrammetry and Remote Sensing. 2015;108:12–32. DOI: 10.1016/j.isprsjprs.2015.06.002

[64] Asrar G, Fuchs M, Kanemasu ET, Hatfield JL. Estimating absorbed photosynthetic radiation and leaf area index from spectral reflectance in wheat. Agronomy Journal.

[65] Richardson AJ, Weigand CL. Distinguishing vegetation from soil background infor‐ mation. Photogrammetric Engineering and Remote Sensing. 1977;43:1541–1552. [66] Baig MHA, Zhang L, Shuai T, Tong Q. Derivation of a tasselled cap transformation based on Landsat 8 at-satellite reflectance. Remote Sensing Letters. 2014;5:423–431.

[67] Yarbrough LD, Easson G, Kuszmaul JS. Using at-sensor radiance and reflectance tas‐ seled cap transforms applied to change detection for the ASTER sensor. ASTER.

[68] Horne JH. A tasseled cap transformation for IKONOS images. In: Proceedings of the InASPRS 2003 Annual Conference, May 2003; Anchorage Alaska; 2003. pp. 60–70. [69] Huang C, Wylie B, Yang L, Homer C, Zylstra G. Derivation of a tasselled cap trans‐ formation based on Landsat 7 at-satellite reflectance. International Journal of Remote

[70] Zhang X, Schaaf CB, Friedl M, Strahler AH, Gao F, Hodges JC. MODIS tasseled cap transformation and its utility. In: Proceedings of the Geoscience and Remote Sensing Symposium, IGARSS'02. 2002 IEEE International (Vol. 2); 24–28 June 2002. Canada:

[71] Crist EP. A TM tasseled cap equivalent transformation for reflectance factor data. Re‐ mote Sensing of Environment. 1985;17:301–306. DOI: 10.1016/0034-4257(85)90102-6 [72] Kauth RJ, Thomas GS. The tasseled cap — A graphic description of the spectral-tem‐ poral development of agricultural crops as seen by Landsat. In: Proceedings of the Symposium on Machine Processing of Remotely Sensed Data; June 29–July 1 1976;

[73] McDonald AJ, Gemmell FM, Lewis PE. Investigation of the utility of spectral vegeta‐ tion indices for determining information on coniferous forests. Remote Sensing of

[74] Sarker L, Nichol J. Improved forest biomass estimates using ALOS AVNIR-2 texture indices. Remote Sensing of Environment. 2011;115:968–977. DOI: 10.1016/j.rse.

[75] Sarker L, Nichol J, Ahmad B, Busu I, Rahman AA. Potential of texture measurements of two-date dual polarization PALSAR data for the improvement of forest biomass estimation. ISPRS Journal of Photogrammetry and Remote Sensing. 2012;69:146–166.

[76] Cutler MEJ, Boyd DS, Foody GM, Vetrivel A. Estimating tropical forest biomass with a combination of SAR image texture and Landsat TM data: an assessment of predic‐

Environment. 1998;66:250–272. DOI: 10.1016/S0034-4257(98)00057-1

1984;76:300–306. DOI: 10.2134/agronj1984.00021962007600020029x

Sensing. 2002;23:1741–1748. DOI: 10.1080/01431160110106113

DOI: 10.1080/2150704X.2014.915434

2005;2:141–145.

34 Environmental Applications of Remote Sensing

IEEE; 2002. pp. 1063–1065.

2010.11.010

West Lafayette Indiana. 1976. pp. 41–51.

DOI: 10.1016/j.isprsjprs.2012.03.002


International Journal of Applied Earth Observation and Geoinformation. 2014;31:45– 56. DOI: 10.1016/j.jag.2014.03.005


[100] Vikhamar D, Solberg R. Subpixel mapping of snow cover in forests by optical remote sensing. Remote Sensing of Environment. 2002;84:69–82. DOI: 10.1016/ S0034-4257(02)00098-6

International Journal of Applied Earth Observation and Geoinformation. 2014;31:45–

[89] Cabacinha CD, de Castro SS. Relationship between floristic diversity and vegetation indices, forest structure and landscape metrics of fragments in Brazilian Cerrado. Forest Ecology and Management. 2009;257:2157–2165. DOI: 10.1016/j.foreco.

[90] Gilabert MAJ, González-Piqueras J, García-Haro J. Acerca de los índices de vegeta‐

[91] Roberts DA, Smith MO, Adams JB. Green vegetation, nonphotosynthetic vegetation and soils in AVIRIS data. Remote Sensing of Environment.

[92] Quintano C, Fernández-Manso A, Shimabukuro YE, Pereira G. Spectral unmixing. International Journal of Remote Sensing. 2012;33:5307–5340. DOI:

[93] Keshava N, Mustard JF. Spectral unmixing. IEEE. Signal Processing Magazine. 2002;

[94] Somers B, Asner GP, Tits L, Coppin P. Endmember variability in spectral mixture analysis: a review. Remote Sensing of Environment. 2011;115:1603–1616. DOI:

[95] Silván-Cárdenas JL, Wang L. Retrieval of subpixel Tamarix canopy cover from Land‐ sat data along the Forgotten River using linear and nonlinear mixture models. Re‐ mote Sensing of Environment. 2010; 114:1777-1790. DOI: 10.1016/j.rse.2010.04.003 [96] Okin GS, Clarke KD, Lewis MM. Comparison of methods for estimation of absolute vegetation and soil fractional cover using MODIS normalized BRDF-adjusted reflec‐ tance data. Remote Sensing of Environment. 2013;130:266–279. DOI: 10.1016/j.rse.

[97] Yang J, Weisberg PJ, Bristow NA. Landsat remote sensing approaches for monitoring long-term tree cover dynamics in semi-arid woodlands: comparison of vegetation in‐ dices and spectral mixture analysis. Remote Sensing of Environment. 2012;119:62–71.

[98] Souza Jr. C, Firestone L, Moreira Silva L, Roberts D. Mapping forest degradation in the Eastern Amazon from SPOT 4 through spectral mixture models. Remote Sensing

[99] Olthof I, Fraser RH. Mapping northern land cover fractions using Landsat ETM+. Re‐ mote Sensing of Environment. 2007;107:496–509. DOI: 10.1016/j.rse.2006.10.009

of Environment. 2003;87:494–506. DOI: 10.1016/j.rse.2002.08.002

56. DOI: 10.1016/j.jag.2014.03.005

10.1080/01431161.2012.661095

10.1016/j.rse.2011.03.003

DOI: 10.1016/j.rse.2011.12.004

ción. Revista de Teledetección. 1997;8:1–10.

1993;44:255−44:255−44:255−44:255−44:255−44:255−269

2009.02.030

36 Environmental Applications of Remote Sensing

19:44–57.

2012.11.021


[123] Andersen H, Reutebuch S, McGaughey R. A rigorous assessment of tree height meas‐ urements obtained using airborne LiDAR and conventional field methods. Canadian Journal of Remote Sensing. 2006;32:355–366.

[111] Cartus O, Santoro M, Kellndorfer J. Mapping forest aboveground biomass in the Northeastern United States with ALOS PALSAR dual-polarization L-band. Remote

[112] Karjalainen M, Kankare V, Vastaranta M, Holopainen M, Hyyppä J. Prediction of plot-level forest variables using TerraSAR-X stereo SAR data. Remote Sensing of En‐

[113] Schlund M, von Poncet F, Hoekman DH, Kuntz S, Schmullius C. Importance of bi‐ static SAR features from TanDEM-X for forest mapping and monitoring. Remote

[114] Bourgeau-Chavez L, Leblon B, Charbonneau F, Buckley J. Evaluation of polarimetric Radarsat-2 SAR data for development of soil moisture retrieval algorithms over a chronosequence of black spruce boreal forest. Remote Sensing of Environment.

[115] Balzter H. Forest mapping and monitoring with interferometric Synthetic Aperture Radar (InSAR). Progress in Physical Geography. 2001;25:159–177. DOI:

[116] Santoro M, Fransson J, Eriksson L, Magnusson M, Ulander L, Olsson H. Signatures of ALOS PALSAR L-band backscatter in Swedish forest. IEEE Transaction on Geosci‐ ence and Remote Sensing. 2009;47:4001–4019. DOI: 10.1109/TGRS.2009.2023906

[117] Le Toan T, Quegan S, Davidson M, Balzter H, Paillou P, Papathanassiou Plummer S, Rocca F, Saatchi S, Shugart H, Ulander L. The BIOMASS mission: mapping global forest biomass to better understand the terrestrial carbon cycle. Remote Sensing of

[118] Carreiras J, Vasconcelos M, Lucas R. Understanding the relationship between above‐ ground biomass and ALOS PALSAR data in the forest of Guinea-Bissau (West Afri‐ ca). Remote Sensing of Environment. 2012;121:426–442. DOI: 10.1016/j.rse.2012.02.012

[119] Tansey KJ, Luckman AJ, Skinner L, Balzter H, Strozzi T, Wagner W. Classification of forest volume resources using ERS tandem coherence and JERS backscatter data. In‐ ternational Journal of Remote Sensing. 2004;25:751–768. DOI:

[120] Lee JS, Pottier E. Polarimetric radar imaging: from basics to applications. Boca Raton:

[121] Gonçalves F, Santos J, Treuhaft R. Stem volume of tropical forest from polarimetric radar. International Journal of Remote Sensing. 2011;32:503–522. DOI:

[122] Lefsky MA, Cohen W, Parker G, Harding D. LiDAR remote sensing for ecosystem

Environment. 2011;115:2850–2860. DOI: 10.1016/j.rse.2011.03.020

CRC press; 2009. 391 p. DOI: 10.1201/9781420054989.fmatt

Sensing of Environment. 2014;151:16–26. DOI: 10.1016/j.rse.2013.08.024

Sensing of Environment. 2012;124:466–478. DOI: 10.1016/j.rse.2012.05.029

vironment. 2012;117: 338–347. DOI: 10.1016/j.rse.2011.10.008

2013;132:71–85. DOI: 10.1016/j.rse.2013.01.006

10.1177/030913330102500201

38 Environmental Applications of Remote Sensing

10.1080/0143116031000149970

10.1080/01431160903475217

studies. BioScience. 2002;52:19–30.


#### **Chapter 2**

## **Detection of Tree Crowns in Very High Spatial Resolution Images**

Marilia Ferreira Gomes and Philippe Maillard

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/62122

#### **Abstract**

[135] Houghton R. Aboveground forest biomass and the global carbon balance. Global

[136] Acosta-Mireles M, Vargas-Hernández J, Velásquez-Martínez A, Etchevers-Barra JD. Estimación de la biomasa aérea mediante el uso de relaciones alométricas en seis es‐

[137] Jiménez M. Ecuaciones alométricas para determinar biomasa y carbono para *Pinus hartwegii* en el Parque Nacional Izta- Popo [thesis]. Texcoco: Universidad Autónoma

[138] Avendaño-Hernández DM, Acosta-Mireles M, Carrillo-Anzures F, Etchevers-Barra JD. Estimación de biomasa y carbono en un bosque de Abies religiosa. Fitotecnia

[139] Šmelko Š, Merganič J. Some methodological aspects of the National Forest Inventory

[140] Asner GP, Heidebrecht KB. Spectral unmixing of vegetation, soil and dry carbon cov‐ er in arid regions: comparing multispectral and hyperspectral observations. Interna‐ tional Journal of Remote Sensing. 2002;23:3939–3958. DOI:

[141] Anselin L, Rey S, editors. Perspectives on Spatial Data Analysis. New York: Springer;

[142] Hengl T, Heuvelink GBM, Stein A, 2003. Comparison of kriging with external drift and regression-kriging. Technical note, ITC, WEB. 3 august 2015. Available online at

[143] Goovaerts P. Geostatistical approaches for incorporating elevation into the spatial in‐ terpolation of rainfall. Journal of Hydrology. 2000;228:113–129. DOI: 10.1016/

[144] GDF. Atlas geográfico del suelo de conservación del Distrito Federal. Ciudad de México: Gobierno del Distrito Federal, Secretaría del Medio Ambiente, Procuraduría

[145] Celemín JP. Autocorrelación espacial e indicadores locales de asociación espacial: Im‐ portancia, estructura y aplicación. Revista Universitaria de Geografía. 2009;18:11–31.

Ambiental y del Ordenamiento Territorial del Distrito Federal; 2012. 96 p.

https://www.itc.nl/library/Papers\_2003/misca/hengl\_comparison.pdf

and Monitoring in Slovakia. Journal of Forest Science. 2008;54:476–483.

Change Biology. 2005;11:945–958. DOI: 10.1111/j.1365-2486.2005.00955.x

pecies arbóreas en Oaxaca, México. Agrociencia. 2002;36:725–736.

Chapingo; 2009.

40 Environmental Applications of Remote Sensing

Mexicana. 2009;32:233–238.

10.1080/01431160110115960

S0022-1694(00)00144-X

2010. 290 p. DOI: 10.1007/978-3-642-01976-0

The requirements for advanced knowledge on forest resources have led researchers to de‐ velop efficient methods to provide detailed information about trees. Since 1999, orbital re‐ mote sensing has been providing very high resolution (VHR) image data. The new generation of satellite allows individual tree crowns to be visually identifiable. The in‐ crease in spatial resolution has also had a profound effect in image processing techniques and has motivated the development of new object-based procedures to extract informa‐ tion. Tree crown detection has become a major area of research in image analysis consid‐ ering the complex nature of trees in an uncontrolled environment. This chapter is subdivided into two parts. Part I offers an overview of the state of the art in computer detection of individual tree crowns in VHR images. Part II presents a new hybrid ap‐ proach developed by the authors that integrates geometrical-optical modeling (GOM), marked point processes (MPP), and template matching (TM) to individually detect tree crowns in VHR images. The method is presented for two different applications: isolated tree detection in an urban environment and automatic tree counting in orchards with an average performance rate of 82% for tree detection and above 90% for tree counting in orchards.

**Keywords:** Tree crown detection, VHR image, Template matching, Marked point process, Valley-following, Watershed segmentation, Local maxima, Region growing

#### **1. Introduction**

Inventories on forest communities are performed with the objective of providing support to the management and conservation activities in rural or urban forests or even in tree plantations. The traditional method of obtaining information on forest communities is to use systematic or random sampling or by sampling stands, so that the final parameters for the population are obtained on the basis of statistical extrapolation [1, 2]. Usually, the following parameters are

© 2016 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

determined for each tree included in the sampling: location, diameter at breast height (DBH), basal area (BA), height, identification of the species, crown size, and crown closure. Based on these measurements, other parameters such as volume of wood and biomass can be derived for the community stand. This renders the field survey techniques for forest inventories expensive, time consuming, and unsuited for large areas.

Remote sensing with high spatial resolution is a cost-effective and reliable way to obtain information about trees. It may be the only practical manner to assure sustainable management of forests with the necessary information, such as biochemical and biophysical data on the vegetation in a synoptic and repetitive manner for large areas and over long periods of time [3]. The tree crown is the basis of the data required for the inventory, for it allows to determine not only its size but also its position, crown closure, and, in some cases, the species. It also allows the derivation of parameters such as the density of the population, the health condition of the trees, the volume, the biomass, and the carbon sequestration rates [3–6]. This information is crucial to a series of applications such as the inventory and management of forested areas as well as in parks and urban forests. It can also be used for counting and monitoring trees in orchards or under power lines to prevent damage and accidents.

The study of individual trees with remote sensing started with the use of aerial photography with very high spatial resolution (scale greater than 1:10.000), driven mainly by the use of stereoscopy techniques. The task was performed by photointerpreters trained to recognize individual tree species, extract a series of measurements, or evaluate different types of damage [2]. The use of orbital optical remote sensing data for forest studies began in the 1970s, with the development of techniques to separate forested from non-forested areas [7]. The spatial resolution of these satellite images was the main limiting factor for more detailed studies about the forests, and as a result, the studies remained focused on the disturbances affecting forests (such as land clearing, burning, diseases, and pest) or to estimate some biophysical parameters of the vegetation [3,8]. It was only toward the end of the 1990s that orbital remote sensing began to provide very high resolution (VHR) data with a spatial resolution under 1 m, allowing the study of individual trees. Launched in 1999, Ikonos was the first of what is now a series of VHR satellites (Table 1), consolidating the use of orbital data for the study of individual trees. However, the increase in spatial resolution was not always accompanied by an increase in spectral resolution for VHR data which is often restricted to a single panchromatic band.

The increase in spatial resolution changed the focus of many remote sensing studies, which started to analyze not only classes of objects but also each object individually [9]. Branches and irregularities within the crowns became visible, and as a result, the spectral response of a tree is influenced by variations in the shape of the crown (differential illumination) and background effects. This causes an increase in the intra-class variance and often results in a reduced accuracy when using conventional pixel-based classification [10]. This had a significant effect on the image processing techniques for forest studies and generated the development of new forms of information extraction.

Within the study of individual objects, the automatic detection and delineation of tree crowns using remote sensing VHR imagery have attracted much attention from researchers in forestry and computer vision [4,7]. Researchers have developed several automatic and semiautomatic methods for extracting individual trees and their characteristics using digital aerial photos of various types and VHR satellite images [11]. The applications range from the identification of tree crowns to their delineation and classification and are often based on image segmentation algorithms and other advanced image processing and analysis techniques [9,12]. Most of these algorithms were specifically developed for the detection and delineation of trees in temperate forests based on the assumption that the trees are cone shaped and round (mostly coniferous) in the images, with the apex of the tree having the highest reflectance of the crown area [4].


\* Panchromatic (Px) and Multispectral (Mx) resolution at nadir.

determined for each tree included in the sampling: location, diameter at breast height (DBH), basal area (BA), height, identification of the species, crown size, and crown closure. Based on these measurements, other parameters such as volume of wood and biomass can be derived for the community stand. This renders the field survey techniques for forest inventories

Remote sensing with high spatial resolution is a cost-effective and reliable way to obtain information about trees. It may be the only practical manner to assure sustainable management of forests with the necessary information, such as biochemical and biophysical data on the vegetation in a synoptic and repetitive manner for large areas and over long periods of time [3]. The tree crown is the basis of the data required for the inventory, for it allows to determine not only its size but also its position, crown closure, and, in some cases, the species. It also allows the derivation of parameters such as the density of the population, the health condition of the trees, the volume, the biomass, and the carbon sequestration rates [3–6]. This information is crucial to a series of applications such as the inventory and management of forested areas as well as in parks and urban forests. It can also be used for counting and monitoring trees in

The study of individual trees with remote sensing started with the use of aerial photography with very high spatial resolution (scale greater than 1:10.000), driven mainly by the use of stereoscopy techniques. The task was performed by photointerpreters trained to recognize individual tree species, extract a series of measurements, or evaluate different types of damage [2]. The use of orbital optical remote sensing data for forest studies began in the 1970s, with the development of techniques to separate forested from non-forested areas [7]. The spatial resolution of these satellite images was the main limiting factor for more detailed studies about the forests, and as a result, the studies remained focused on the disturbances affecting forests (such as land clearing, burning, diseases, and pest) or to estimate some biophysical parameters of the vegetation [3,8]. It was only toward the end of the 1990s that orbital remote sensing began to provide very high resolution (VHR) data with a spatial resolution under 1 m, allowing the study of individual trees. Launched in 1999, Ikonos was the first of what is now a series of VHR satellites (Table 1), consolidating the use of orbital data for the study of individual trees. However, the increase in spatial resolution was not always accompanied by an increase in spectral resolution for VHR data which is often restricted to a single panchromatic band.

The increase in spatial resolution changed the focus of many remote sensing studies, which started to analyze not only classes of objects but also each object individually [9]. Branches and irregularities within the crowns became visible, and as a result, the spectral response of a tree is influenced by variations in the shape of the crown (differential illumination) and background effects. This causes an increase in the intra-class variance and often results in a reduced accuracy when using conventional pixel-based classification [10]. This had a significant effect on the image processing techniques for forest studies and generated the development of new

Within the study of individual objects, the automatic detection and delineation of tree crowns using remote sensing VHR imagery have attracted much attention from researchers in forestry and computer vision [4,7]. Researchers have developed several automatic and semiautomatic methods for extracting individual trees and their characteristics using digital aerial photos of

expensive, time consuming, and unsuited for large areas.

42 Environmental Applications of Remote Sensing

orchards or under power lines to prevent damage and accidents.

forms of information extraction.

**Table 1.** Very High Resolution satellites (1999-2006), with their spatial resolutions and spectral bands. Note that there is not any available satellite with VHR in multispectral bands.

The analysis of individual trees based on remote sensing images is a complex problem. Images of trees with varied crown size increase the difficulty of the analysis. What is detected as a single object may in fact represent a separate branch or even a group of trees [10]. Other sources of error are caused by the proximity between neighboring trees, trees located under other trees, trees in the shade, or trees that have a low spectral contrast with the background [13]. Conse‐ quently, high-level complex algorithms are necessary to exploit this contextual information [1].

This chapter provides an overview of the state of the art in individual tree crown detection based on optical VHR remote sensing data. An original method developed by the authors is also presented as an alternative approach to the problem of tree crown detection. In Part I, we present the main algorithms developed for the detection of individual trees, be it for tree identification or delineation. The principle of each approach is presented as well as its potential and limitations. Part II is dedicated to outlining the original MPP–TM approach, a hybrid method that combined two methods used in pattern recognition: marked point process and template matching. The results are shown for tree detection and delineation in an urban environment and for tree counting in orchards.

#### **2. Part I – Review of tree crown detection methods**

We present six of the main algorithms used in individual tree detection in high spatial resolution images. The algorithms are summarily described individually, but it should be noted that many approaches use hybrid methods for the detection and delineation of tree crowns. For instance, some authors might use one algorithm for detecting the trees and another to delineate them; some may even use one approach as a first approximation and another to fine-tune the results.

#### **2.1. Local maxima filtering**

Local maxima (LM) filtering is a technique used for identifying tree crowns in high spatial resolution imagery which is based on the recognition of the points with the greatest brightness within a search window that scans the entire image [4,14]. The search window, with a fixed size, defines which pixel has the greatest reflectance compared to all the other pixels inside the window. The pixels with the highest digital number are identified as possible tree locations. This method is adequate for trees which have the greatest reflectance at their top, surrounded by lower intensity pixels, and due to its concept, it is widely used for detecting conifers.

When the kernel window passes over the image, it does not take into account the presence of trees with different crown sizes, and the success of the LM tree recognition depends on the careful selection of the size of the search window. If it is too small, errors of commission occur by selecting nonexistent trees or multiple radiance peaks for an individual tree crown; if it is too big, the algorithm is likely to miss some trees (omission errors) [13].

The identification of trees by LM is affected by false bright pixels, which are not part of the brightest part of the crown. An effective method for dealing with the problem is to apply a Gaussian filter to the image. This allows the low-pass filter function to grant more weight to the crown center pixels (surrounded by much lower values) compared to those located toward the crown edge which might belong to other bright objects or noise. Applying a Gaussian filter directly affects the number of local maxima identified and causes the smoothing off of the brightness values on the tree crown edges [15].

In order to minimize the problem of the window size with LM, reference [13] used windows of varying sizes based on the assessment of the spatial structure of the image obtained by analyzing the local semi-variogram with different pixel lags and different window sizes. This results in a personalized window for each pixel, leading to greater accuracy when compared to using a single fixed window size. Reference [16] used LM to identify the centroid of eucalyptus trees in Australia. The search for the trees is carried out based on the maxima in four linear kernels pertaining to the four main directions (0°, 45°, 90°, and 135°) of the image and by summing the individual maxima found in each pass (Figure 1).

**Figure 1.** Examples of the surface produced by applying a LM kernel operator: (a) original image, (b) the local maxima appearing in the third dimension are associated with the presence of trees, (c) application of local maxima filter in four linear cumulative kernels (0°, 45°, 90°, and 135°).

#### **2.2. Template matching**

identification or delineation. The principle of each approach is presented as well as its potential and limitations. Part II is dedicated to outlining the original MPP–TM approach, a hybrid method that combined two methods used in pattern recognition: marked point process and template matching. The results are shown for tree detection and delineation in an urban

We present six of the main algorithms used in individual tree detection in high spatial resolution images. The algorithms are summarily described individually, but it should be noted that many approaches use hybrid methods for the detection and delineation of tree crowns. For instance, some authors might use one algorithm for detecting the trees and another to delineate them; some may even use one approach as a first approximation and another to

Local maxima (LM) filtering is a technique used for identifying tree crowns in high spatial resolution imagery which is based on the recognition of the points with the greatest brightness within a search window that scans the entire image [4,14]. The search window, with a fixed size, defines which pixel has the greatest reflectance compared to all the other pixels inside the window. The pixels with the highest digital number are identified as possible tree locations. This method is adequate for trees which have the greatest reflectance at their top, surrounded by lower intensity pixels, and due to its concept, it is widely used for detecting conifers.

When the kernel window passes over the image, it does not take into account the presence of trees with different crown sizes, and the success of the LM tree recognition depends on the careful selection of the size of the search window. If it is too small, errors of commission occur by selecting nonexistent trees or multiple radiance peaks for an individual tree crown; if it is

The identification of trees by LM is affected by false bright pixels, which are not part of the brightest part of the crown. An effective method for dealing with the problem is to apply a Gaussian filter to the image. This allows the low-pass filter function to grant more weight to the crown center pixels (surrounded by much lower values) compared to those located toward the crown edge which might belong to other bright objects or noise. Applying a Gaussian filter directly affects the number of local maxima identified and causes the smoothing off of the

In order to minimize the problem of the window size with LM, reference [13] used windows of varying sizes based on the assessment of the spatial structure of the image obtained by analyzing the local semi-variogram with different pixel lags and different window sizes. This results in a personalized window for each pixel, leading to greater accuracy when compared to using a single fixed window size. Reference [16] used LM to identify the centroid of

too big, the algorithm is likely to miss some trees (omission errors) [13].

brightness values on the tree crown edges [15].

environment and for tree counting in orchards.

44 Environmental Applications of Remote Sensing

fine-tune the results.

**2.1. Local maxima filtering**

**2. Part I – Review of tree crown detection methods**

Template matching (TM) is a technique used for object recognition widely cited in the specialized literature which uses quantitative descriptors, such as length, area, and texture to describe recurring patterns in an image [17,18]. Based on a synthetic model or a sample extracted from the image, the correlation coefficient between the model and the image is calculated in order to determine the strength of the match between the two matrices. The object is assumed to be located where the measurement of the match reaches a maximum [4].

For tree crown detection, the study of reference [19] was the first to propose an elliptical 3D model for tree crowns based on an ellipsoid of generalized revolution (Equation 1).

$$\frac{\mathbf{z}^u}{a^u} - \frac{\left(\mathbf{x}^2 + \mathbf{y}^2\right)^{u/2}}{b^u} - \mathbf{1} \tag{1}$$

where *z* is the vertical axis of the center of the tree crown in its origin, *a* is half the height of the ellipsoid, *b* is half the radius, and *n* is the parameter of the shape of the tree crown. Subse‐ quently, the model is illuminated using the acquisition parameters of the image (sun elevation and azimuth) and the characteristics of crown absorption and reflection of light in the chosen spectral band.

Because it is based on a physical model (rather than a complex mathematical concept), TM is considered a user-friendly method. Its limitation is mainly due to the need to use a library of models if many types of trees are present in the image, which may involve a complex phase for generating the models. Figure 2 shows examples of synthetic tree models and an application in a orchard.

References [20 and 21] used this technique to identify tree crowns in aerial images. Other researchers used this technique to recognize individual tree crowns, using templates made from small sub-images of the actual scene to identify the trees [22, 23]. Reference [24] proposed an improved version by generating separate models for trees and their shade in VHR images of unmanned aerial vehicles (UAV). The authors explored the relation between the tree and shade models separately and then joined them to generate a more powerful object detector.

**Figure 2.** Left: examples of synthetic tree models to different tree crown shape. Right: identification of trees in an apple orchard showing the model used at the upper right corner.

#### **2.3. Valley-following**

Valley-following (VF) is a crown delineation method which identifies the shaded areas between the trees. This methodology was initially described in reference [25] and makes an analogy with topographic data, where the shades of gray of the pixels represent local lows in the third dimension. In this analogy, the bright tree crowns would be the hills and the darker zones around the trees the valleys (Figure 3). This darker zone is the one which typically helps human interpreters to separate one tree crown from the other. In this approach the shaded areas are eliminated, making it possible to separate the trees in the image. This was not sufficient to separate all of the trees, so the authors developed an approach based on a series of rules (e.g., no discontinuity, checking directions, context, gap filling, etc.) to accurately describe the boundaries of each tree, one at a time [26].

This approach performed well in images with a combination of low solar elevation angle and conical trees. Conversely, the approach failed to produce good results when the canopy was composed of trees of very different sizes, or when the tree crowns were very large and have internal shadows. The latter case resulted in subdividing the individual tree into two or more parts. Smaller trees, in contrast, tended to be grouped together. Reference [27] found that this

**Figure 3.** Results of the valley-following applied on a forest image in Canada (source: Gougeon and Leckie, 2003; re‐ produced with permission from Natural Resources Canada).

approach causes many false positives (FPs) in open areas (clearings). As a solution, they suggested the exclusion of these areas by retaining only the high-value pixels in the normalized difference vegetation index (NDVI).

#### **2.4. Watershed**

Because it is based on a physical model (rather than a complex mathematical concept), TM is considered a user-friendly method. Its limitation is mainly due to the need to use a library of models if many types of trees are present in the image, which may involve a complex phase for generating the models. Figure 2 shows examples of synthetic tree models and an application

References [20 and 21] used this technique to identify tree crowns in aerial images. Other researchers used this technique to recognize individual tree crowns, using templates made from small sub-images of the actual scene to identify the trees [22, 23]. Reference [24] proposed an improved version by generating separate models for trees and their shade in VHR images of unmanned aerial vehicles (UAV). The authors explored the relation between the tree and shade models separately and then joined them to generate a more powerful object detector.

(a) (b)

orchard showing the model used at the upper right corner.

describe the boundaries of each tree, one at a time [26].

**2.3. Valley-following**

**Figure 2.** Left: examples of synthetic tree models to different tree crown shape. Right: identification of trees in an apple

Valley-following (VF) is a crown delineation method which identifies the shaded areas between the trees. This methodology was initially described in reference [25] and makes an analogy with topographic data, where the shades of gray of the pixels represent local lows in the third dimension. In this analogy, the bright tree crowns would be the hills and the darker zones around the trees the valleys (Figure 3). This darker zone is the one which typically helps human interpreters to separate one tree crown from the other. In this approach the shaded areas are eliminated, making it possible to separate the trees in the image. This was not sufficient to separate all of the trees, so the authors developed an approach based on a series of rules (e.g., no discontinuity, checking directions, context, gap filling, etc.) to accurately

This approach performed well in images with a combination of low solar elevation angle and conical trees. Conversely, the approach failed to produce good results when the canopy was composed of trees of very different sizes, or when the tree crowns were very large and have internal shadows. The latter case resulted in subdividing the individual tree into two or more parts. Smaller trees, in contrast, tended to be grouped together. Reference [27] found that this

in a orchard.

46 Environmental Applications of Remote Sensing

Like VF, the watershed segmentation (WS) is a technique related to thresholding that uses the gray levels in the images as if it were a topographic surface [28]. It is used not only for the delineation of individual tree crowns but also for generic segmentation of images. The watershed concept is based on a 3D image representation, with the third dimension being provided by the intensity of gray. The main objective of the watershed algorithm is to find the "drainage" divide lines. The "relief" in the image is inverted (high gray values become valleys) and progressively filled with a virtual liquid, and when the liquid is almost overflowing from one basin to another, a virtual dam is built, to create the watershed. These lines are considered the limits of each segment. The simplest approach to the construction of the dam is the use of morphological dilation of the minima, without merging the regions [17].

The images are usually preprocessed before the WS is applied. In fact, this segmentation is frequently applied to the gradient of an image, and not to the image itself. This is due to the relative homogeneity of the gray values of objects that do not provide sufficient contrast for an effective segmentation. In this formulation, the regional minimum value of the catchment basins usually correlates well with the lower gradient values that match the contours of the objects of interest [17,28]. The direct application of the WS algorithm generally leads to oversegmentation due to noise or other local irregularities of the gradient (Figure 4a). One of the approaches used to limit the number of regions is to use markers. The selection of markers can be based on simple procedures, intensity and connectivity between pixels, or even complex descriptors, such as size, format, location, relative distances, texture, and others. The use of markers provides prior knowledge to support the segmentation process [17].

The approaches that use WS for the delineation of the tree crowns normally use markers representing the center of the tree crown, to assist the segmentation process. For instance, reference [29] used WS to detect and delineate tree crowns in a VHR forest image in Canada but divided the approach into two phases, namely using LM to detect the crown and applying WS for the delineation. The LM image with the detected crowns was produced by using a Laplacian of Gaussian edge detection operator. The tree crowns were modeled based on their geometry and radiometry, resulting in an image of markers. This image then served to guide the WS in delineating the crowns. Reference [30] developed a bitemporal procedure for the automatic segmentation and reconciliation of groups of pixels (called blobs) within the forest using WS. By using two dates, they were able to increase the probability of properly defining the tree contours. Many problems were encountered in the segmentation process of the individual trees. For instance, trees with spread branches were sometimes split into two or more segments or contrarily by including several crowns in the same segment when trees were not sufficiently separated.

#### **2.5. Region growing**

Region growing (RG) is another segmentation technique that groups pixels or groups of pixels based on predefined growth criteria in an attempt to separate and recognize objects in the image [4]. Like WS, RG is used as a generic segmentation method and can be adapted for the delineation of individual tree crowns (Figure 4b). Starting with some seed pixels (which can be random if no other information is provided), the neighboring pixels are examined one by one and added to the growth region if their predefined properties are similar to those of the seeds (such as specific intervals of intensity or color) [17]. When no more pixels can be added or some predefined limit is reached (e.g., number of pixels), these pixels are labeled as belonging to the specific region of the seed pixel. Additional criteria can increase the power of an RG algorithm by introducing a higher concept like size and similarity between candidate pixels and the pixels selected or even the format of the region [17,31].

Reference [16] used RG integrated with LM to identify and delineate tree crowns in Australia. The LM method served to find the center of potential trees, which were then used as seeds for the RG. Reference [6] tested two different types of segmentation by RG, one by Brownian motion and the other by random walk, to detect conifers in a boreal forest. The methods were capable of detecting about 80% of the illuminated portion of the crowns, with a better per‐ formance found in larger crowns (Figure 5).

(a) (b)

(c)

**Figure 4.** Comparison between two segmentation algorithms: (a) watershed and (b) region growing on a WorldView-2 image (panchromatic band with a 50 cm resolution). The WS was applied to the gradient image without using markers and resulted in the over-segmentation of the tree crowns. In (b), the RG segmentation was performed within an objectoriented classification (GEOBIA) approach, where the correct delineation of the tree crowns is noteworthy (source: Gomes and Maillard, 2013).

#### **2.6. Marked point processes**

an effective segmentation. In this formulation, the regional minimum value of the catchment basins usually correlates well with the lower gradient values that match the contours of the objects of interest [17,28]. The direct application of the WS algorithm generally leads to oversegmentation due to noise or other local irregularities of the gradient (Figure 4a). One of the approaches used to limit the number of regions is to use markers. The selection of markers can be based on simple procedures, intensity and connectivity between pixels, or even complex descriptors, such as size, format, location, relative distances, texture, and others. The use of

The approaches that use WS for the delineation of the tree crowns normally use markers representing the center of the tree crown, to assist the segmentation process. For instance, reference [29] used WS to detect and delineate tree crowns in a VHR forest image in Canada but divided the approach into two phases, namely using LM to detect the crown and applying WS for the delineation. The LM image with the detected crowns was produced by using a Laplacian of Gaussian edge detection operator. The tree crowns were modeled based on their geometry and radiometry, resulting in an image of markers. This image then served to guide the WS in delineating the crowns. Reference [30] developed a bitemporal procedure for the automatic segmentation and reconciliation of groups of pixels (called blobs) within the forest using WS. By using two dates, they were able to increase the probability of properly defining the tree contours. Many problems were encountered in the segmentation process of the individual trees. For instance, trees with spread branches were sometimes split into two or more segments or contrarily by including several crowns in the same segment when trees were

Region growing (RG) is another segmentation technique that groups pixels or groups of pixels based on predefined growth criteria in an attempt to separate and recognize objects in the image [4]. Like WS, RG is used as a generic segmentation method and can be adapted for the delineation of individual tree crowns (Figure 4b). Starting with some seed pixels (which can be random if no other information is provided), the neighboring pixels are examined one by one and added to the growth region if their predefined properties are similar to those of the seeds (such as specific intervals of intensity or color) [17]. When no more pixels can be added or some predefined limit is reached (e.g., number of pixels), these pixels are labeled as belonging to the specific region of the seed pixel. Additional criteria can increase the power of an RG algorithm by introducing a higher concept like size and similarity between candidate

Reference [16] used RG integrated with LM to identify and delineate tree crowns in Australia. The LM method served to find the center of potential trees, which were then used as seeds for the RG. Reference [6] tested two different types of segmentation by RG, one by Brownian motion and the other by random walk, to detect conifers in a boreal forest. The methods were capable of detecting about 80% of the illuminated portion of the crowns, with a better per‐

pixels and the pixels selected or even the format of the region [17,31].

formance found in larger crowns (Figure 5).

markers provides prior knowledge to support the segmentation process [17].

not sufficiently separated.

48 Environmental Applications of Remote Sensing

**2.5. Region growing**

The marked point process refers to a probabilistic method which has been used in recent years for the recognition of objects in high spatial resolution imagery [5,11,32–35]. In an MPP, sets of random points in a given space (*x*, *y*) are provided with a mark which is complete and separable, allowing the definition of a topology (defined by the mark) and the attribution of a label. An image is considered a random model where the gray tones are the realization of a random point process [34]. This random configuration of gray levels in the images is then modeled based on geometric figures (ellipses, circles, rectangles, and lines), respecting certain geometric (nature of the objects) and radiometric constraints (type of image).

The laws of density and probability distinguish various types of point processes, which can be Poisson, Strauss, Markov, or Gibbs, among others. The Markov or Gibbs point processes have been used for the recognition of tree crowns by a number of authors [5,32,33]. These

**Figure 5.** Example of RG segmentation to delineate tree crowns in a boreal forest. (a) The original image (with a spatial resolution of 3 cm), (b) with the results using Brownian motion, (c) and random walk (Source: Erikson, 2004).

processes are defined by a density function using a form of energy expressed as a sum of the a priori energy and the local energy. The process seeks to minimize the global energy of the model, by iterating it with some optimization scheme (Markov random fields, algorithm of multiple births and deaths, and Monte Carlo-Monte Carlo simulations). (a) (b)

Reference [5] proposed two different models to serve as marks in an MPP, one in 2D for detection of trees in densely forested zones (Figure 6) and the other in 3D for scattered or isolated zones, based on aerial photos of high spatial resolution in the infrared band. The MPP was integrated with a reversible jump Markov Chain Monte Carlo in a simulated annealing

**Figure 6.** Example of the application of MPP for tree crown recognition on a poplar plantation. The original image is shown at the top and the results at the bottom (Source: Perrin, 2006).

method. Reference [32] used an MPP to automatically detect the tree crowns in high spatial resolution images, based on the modeling of the crowns as 2D circles on high-resolution spatial images. The method was successfully tested on mangrove forests and eucalyptus plantations.

#### **2.7. Discussion**

processes are defined by a density function using a form of energy expressed as a sum of the a priori energy and the local energy. The process seeks to minimize the global energy of the model, by iterating it with some optimization scheme (Markov random fields, algorithm of

(b)

(a)

(c)

**Figure 5.** Example of RG segmentation to delineate tree crowns in a boreal forest. (a) The original image (with a spatial resolution of 3 cm), (b) with the results using Brownian motion, (c) and random walk (Source: Erikson, 2004).

(a) (b)

50 Environmental Applications of Remote Sensing

Reference [5] proposed two different models to serve as marks in an MPP, one in 2D for detection of trees in densely forested zones (Figure 6) and the other in 3D for scattered or isolated zones, based on aerial photos of high spatial resolution in the infrared band. The MPP was integrated with a reversible jump Markov Chain Monte Carlo in a simulated annealing

multiple births and deaths, and Monte Carlo-Monte Carlo simulations).

In the previous section, we have presented some of the most common algorithms used in the detection of individual trees, be it for their identification, delineation, or both. Table 2 presents a summary of these principles through their main characteristics and limitations.

Trees may differ in shape, size, spectral properties, height, foliage type, and density, and their spatial context varies with illumination, ground type, and inclination. They can also be surrounded by many other objects, especially in an urban setting. As such, the task is not trivial and can become highly complex depending on the number of parameters involved. Converse‐ ly, in planted forest and orchards where trees have the same age and species, tree crown extraction can take advantage of their relative uniformity.



**Algorithm Usage Principle Researches Characteristics / Limitations**

Wulder *et al.* (2000)

Culvenor (2002) Pouliot (2002), Wang *et al.* (2004)

Pollock (1996) Larsen (1997) Larsen and Rudemo (1998) Quackenbush *et al.* (2000) Erikson (2004) Hung *et al.* (2012)

Gougeon (1995,

1999) Leckie and Gougeon (1998) Gougeon and Leckie (2003) Erickson (2004) Gougeon and Leckie (2006)

Appropriate for conifers, with a conical shape and high reflectance point at the

Results are affected by the spatial distribution of trees, variation of tree crowns size, search window size (increased omission errors in larger windows and commission errors in

Enables analysis of the tree crown from its spectral, textural and structural

Allows neighborhood analysis of the tree crown by considering its shadow.

Needs a template library, making it unpractical in complex forests. Recognition errors increase with irregularity of the tree crowns.

Easier to detect larger trees than smaller

Appropriate for trees with conical shape that create shadow areas between

Most successful to delineate populations of the same age without intertwined tree

Best performance for images in mid-low

Performance reduced when trees are asymmetrical, of from different species, with different tree crown sizes or when shadows of trees protrude over each

Tendency to group smaller trees together and split larger trees into

Performance reduced in very dense

top of the tree. Simple method to use.

smaller windows).

characteristics.

ones.

environments.

individuals.

crowns.

other.

solar elevation angle.

multiple segments.

User-friendly method.

Identification of brightest points locally as corresponding to the apex of a treetop within a search window.

Quantitative descriptors used to describe patterns. Calculate the correlation between the image and the model. Model may be a sample extracted from the image or not.

Derives from an analogy with a topographical surface, programmed to identify the shaded portion between the tree crowns (valleys).

**Local Maxima (LM)**

**Template Matching (TM)**

**(VF)**

Identification of tree crown

52 Environmental Applications of Remote Sensing

Pattern recognition

Delineation of tree crown

**Table 2.** Summary of Local Maxima, Template Matching, Valley Following, Watershed, Region Growing and Marked Point Processes algorithms used to individual tree crown detection. The principles, main researches and main characteristics and limitations are presented.

Reference [11] compared six different algorithms (valley-following, region growing, template matching, scale-space theory, marked point processes, and Markov random fields) in six different aerial images, ranging from a homogeneous plantation and an area with isolated tree crowns to an extremely dense deciduous forest type. The authors found that none of the algorithms can by itself reach a high rate of success in all of the tested images and concluded that there is no single optimum algorithm for all types of images and forests. They also emphasized that for complex types of forests, monoscopic images are insufficient for a consistent detection of tree crowns, even for human interpreters.

### **3. Part II – A hybrid approach integrating marked point process and template matching**

As shown in our brief review, many methods have been developed for trees in temperate forest environments. In an exploratory research [23], three algorithms in urban tropical environments were tested: region growing, watershed, and template matching. Better results were generally obtained by combining region-growing segmentation and geographic object-based image analysis (GEOBIA) for classification. Although highly effective, the approach requires much parameter setting and experience and is not especially dedicated to the problem of tree crown detection.

Studies that use marked point processes have triggered our attention and made us consider that they could benefit from using marks modeled from 3D objects in a different approach than from that developed by reference [5]. We propose to use a geometrical-optical tree model in a manner resembling that of template matching that uses some form of correlation between image and model to identify candidate pixels. An MPP taking advantage of a geometrical optical 3D model and measurements of similarity to seek tree crowns could represent a significant improvement to using simpler marks. Considering such a hypothesis, we devel‐ oped an algorithm for tree crown detection that combines elements from MPP, TM, and tree crown geometrical-optical modeling for the automatic detection and (simplified) delineation of trees in VHR satellite imagery. We have named our algorithm MPP–TM.

In our approach, the TM did not scan the whole image like it was initially conceived but rather uses an MPP approach to select random locations within the image. Additionally, the 3D marks receive a random diameter between a predetermined range depending on the type of envi‐ ronment. The geometrical-optical model includes both the sunlit and shaded areas of the crown and a portion of the projected shadow to allow a better match between model and image. Some statistical and spectral parameters were also included in the model-matching phase.

MPP-based algorithm for pattern recognition usually alternates between phases of birth and death during which the objects are created (placed) and destroyed when they do not comply with the matching rules. This is also a characteristic of MPP–TM, but we have somewhat deviated from the original concept where the destruction phase also incorporated a random process.

The following subsections are devoted to describe the construction process of the 3D geomet‐ rical-optical model and the functioning of the algorithm.

#### **3.1. Description of MPP–TM approach**

#### *3.1.1. A geometrical-optical 3D tree crown model*

The parameters that determine the radiance pattern of a tree crown are direct and indirect radiation, shape of tree, branch pattern, leaf reflectance, multiple reflectances within the canopy, etc. [36]. In creating a valid 3D geometrical-optical model, we have chosen a simplified version in which the crown is represented by a dome of varying skewness, a Lambertian reflectance model with ambient light, and a projected shadow on the ground (or on another tree). Equations 2 and 3 give the formulation of our model in which each pixel is treated as a singular surface.

$$\cos(\theta\_i) = \left(\cos(\theta\_s)\cos(\theta\_n) + \sin(\theta\_s)\sin(\theta\_n)\cos(\varphi\_s - \varphi\_n)\right) \tag{2}$$

where *θ<sup>i</sup>* is the local solar incidence angle, *θ*s is the solar zenith angle, *θ*n is the slope of the object surface, φs is the solar azimuth and *φ*n is the aspect of the object surface.

$$\mathbf{L} = L\_M \left( \frac{\cos(\theta\_s)}{\cos(\theta\_i)} \right) + \mathbf{amb} \tag{3}$$

where *LM* is the maximum reflectance of the model, "amb" represents the diffuse ambient lighting. The geometrical-optical model is adjusted according to the specific illumination parameters of the image, and the size of the trees present on the scene. Figure 7 shows two examples of tree models with similar reflectance but different solar elevations.

**Figure 7.** Illustration of the geometrical-optical model of tree crown as seen in the same sun azimuth (32°) but in two different solar elevation angles: (left) 20° and (right) 45°

A parameter of projected shadow clipping has also been added to account for the fact that it was not beneficial to use the whole shadow in situations where it was projected onto another tree and not on the ground. The height of the tree also affects the size of the shadow so that it did not appear wise to set the height to a fixed value. To illustrate this, Figure 8 shows a comparison between the tree model and an actual tree from the image both with whole and clipped shadows.

#### *3.1.2. Algorithm description*

**3. Part II – A hybrid approach integrating marked point process and**

As shown in our brief review, many methods have been developed for trees in temperate forest environments. In an exploratory research [23], three algorithms in urban tropical environments were tested: region growing, watershed, and template matching. Better results were generally obtained by combining region-growing segmentation and geographic object-based image analysis (GEOBIA) for classification. Although highly effective, the approach requires much parameter setting and experience and is not especially dedicated to the problem of tree crown

Studies that use marked point processes have triggered our attention and made us consider that they could benefit from using marks modeled from 3D objects in a different approach than from that developed by reference [5]. We propose to use a geometrical-optical tree model in a manner resembling that of template matching that uses some form of correlation between image and model to identify candidate pixels. An MPP taking advantage of a geometrical optical 3D model and measurements of similarity to seek tree crowns could represent a significant improvement to using simpler marks. Considering such a hypothesis, we devel‐ oped an algorithm for tree crown detection that combines elements from MPP, TM, and tree crown geometrical-optical modeling for the automatic detection and (simplified) delineation

In our approach, the TM did not scan the whole image like it was initially conceived but rather uses an MPP approach to select random locations within the image. Additionally, the 3D marks receive a random diameter between a predetermined range depending on the type of envi‐ ronment. The geometrical-optical model includes both the sunlit and shaded areas of the crown and a portion of the projected shadow to allow a better match between model and image. Some

MPP-based algorithm for pattern recognition usually alternates between phases of birth and death during which the objects are created (placed) and destroyed when they do not comply with the matching rules. This is also a characteristic of MPP–TM, but we have somewhat deviated from the original concept where the destruction phase also incorporated a random

The following subsections are devoted to describe the construction process of the 3D geomet‐

The parameters that determine the radiance pattern of a tree crown are direct and indirect radiation, shape of tree, branch pattern, leaf reflectance, multiple reflectances within the canopy, etc. [36]. In creating a valid 3D geometrical-optical model, we have chosen a simplified version in which the crown is represented by a dome of varying skewness, a Lambertian

statistical and spectral parameters were also included in the model-matching phase.

of trees in VHR satellite imagery. We have named our algorithm MPP–TM.

rical-optical model and the functioning of the algorithm.

**3.1. Description of MPP–TM approach**

*3.1.1. A geometrical-optical 3D tree crown model*

**template matching**

54 Environmental Applications of Remote Sensing

detection.

process.

According to reference [32], using MPP to extract objects consists in searching for the "best" possible object configuration in a scene, the one that will respect a certain number of properties

**Figure 8.** Comparison between an isolated tree from the (a) WorldView-2 image and (b) the geometrical-optical 3D model. A clipping factor of about 80% was applied to the same images in (c) and (d) to enable the use of only a portion of the shadow in cases where that shadow is not projected on the ground but on another object.

both of the objects being sought and the radiometric properties of the image. In our algorithm, the "best" configuration be it geometric or radiometric is given by the model.

The process consists in alternating phases of birth and death. The MPP starts with a birth phase during which tree crowns represented by circles of varying size (a randomized interval) are inserted on a matrix of equal size to the image being processed. Tree crowns are only inserted where no other crowns are present. Once all the circles have been inserted (determined by a density parameter *N*c), a similarity (*Sm*) value between the image and a version of the model fitted on each circle is computed and stored in a list along with the parameters of the model. A routine then sorts the list by decreasing values of *Sm*. During the death phase, the circles that do not comply with the acceptance restrictions (*Sm* and minimum and maximum standard deviation threshold) are successively deleted from the matrix and the list. At the end of each death phase, the overall parameters of the pixel distribution of the remaining crowns are updated. The crowns that have been found are definitively kept but are re-thrown in the bundle of crowns of the next iterations. If after an iteration 10 crowns are kept out of 100, then the next iteration will randomly place 90 more crowns, and the new set of 100 crowns are evaluated and sorted for the next iteration. All tree crowns are considered "found" when one of the three possibilities is encountered: 1) the number of trees found is equal to the number given by the density parameter, 2) one of the interruption criteria has been attained, or 3) the maximum number of iteration has been reached.

The *Sm* value is computed as the subtraction of two parameters: cross-correlation and the normalized absolute difference as defined by the following relation (Equation 4):

$$\text{Sm} = \gamma \cdot \alpha \,\text{ND} \tag{4}$$

where *γ* is the cross-correlation between image and model, ND is the normalized sum of absolute differences between them, and *α* is a constant weight factor (normally approx. 0.5). The cross-correlation and absolute difference are calculated as follows (Equations 5 and 6).

$$\chi(\mathbf{x}, y) = \frac{\sum\_{i} \sum\_{j} \left[ w\_{(i,j)} - \overline{w} \right] \sum\_{i} \sum\_{j} \left[ f\_{(i,j)} - \overline{f} \right]}{\left\{ \sum\_{i} \sum\_{j} \left[ w\_{(i,j)} - \overline{w} \right]^2 \sum\_{i} \sum\_{j} \left[ f\_{(i,j)} - \overline{f} \right]^2 \right\}^{1/2}} \tag{5}$$

In Equation 5, the cross-correlation is calculated between the model matrix (*w*(*i*, *<sup>j</sup>*)) and the portion of the image that corresponds to the circle of the same radius ( *f* (*i*, *<sup>j</sup>*)). In other words, two matrices of same dimensions are always compared. *w*¯ and *<sup>f</sup>* ¯ are their respective means. The values of gamma range between −1 and 1. The same logic is used in Equation 6 which computes a normalized difference value between the same two matrices.

both of the objects being sought and the radiometric properties of the image. In our algorithm,

**Figure 8.** Comparison between an isolated tree from the (a) WorldView-2 image and (b) the geometrical-optical 3D model. A clipping factor of about 80% was applied to the same images in (c) and (d) to enable the use of only a portion

(c) (d)

(a) (b)

The process consists in alternating phases of birth and death. The MPP starts with a birth phase during which tree crowns represented by circles of varying size (a randomized interval) are inserted on a matrix of equal size to the image being processed. Tree crowns are only inserted where no other crowns are present. Once all the circles have been inserted (determined by a density parameter *N*c), a similarity (*Sm*) value between the image and a version of the model fitted on each circle is computed and stored in a list along with the parameters of the model. A routine then sorts the list by decreasing values of *Sm*. During the death phase, the circles that do not comply with the acceptance restrictions (*Sm* and minimum and maximum standard deviation threshold) are successively deleted from the matrix and the list. At the end of each death phase, the overall parameters of the pixel distribution of the remaining crowns are updated. The crowns that have been found are definitively kept but are re-thrown in the bundle of crowns of the next iterations. If after an iteration 10 crowns are kept out of 100, then the next iteration will randomly place 90 more crowns, and the new set of 100 crowns are evaluated and sorted for the next iteration. All tree crowns are considered "found" when one of the three possibilities is encountered: 1) the number of trees found is equal to the number given by the

the "best" configuration be it geometric or radiometric is given by the model.

of the shadow in cases where that shadow is not projected on the ground but on another object.

56 Environmental Applications of Remote Sensing

$$\text{NDD} = \frac{\sum w - \sum f}{\sum w + \sum f} \tag{6}$$

In the death phase, tree crowns are kept if their similarity is larger or equal to a pre-set threshold. Because we found that such a threshold represented a weak element in our algorithm, we implemented a strategy by which it needs not be predetermined with a fixed value but rather adjusts itself as the number of iterations grows. The threshold is set very high at the beginning but then starts to decay when a certain number of iterations do not find any "new" tree crown (typically 100 iterations). Additionally, if more than a certain amount of iterations (say 1000) still does not add any new tree crown, then the process is stopped. Ultimately, it will be stopped if the maximum number of iterations is reached. A flowchart of our algorithm is presented in Figure 9 and schematically described in Table 3.

Figure 10a shows the state of the crown matrix after a single birth phase with 163 circles of random radius (between 3 and 15 m) and randomly located within the image matrix. After the death phase, using a similarity threshold of 0.98, only one tree crown was kept (Figure 10b).

**Figure 9.** Flowchart of the MPP–TM algorithm.

#### *3.1.3. A modified approach for orchards*

Because trees in orchards are often individually distinguishable and have similar shape and size, they are perfect candidates for TM with a 3D geometrical-optical model. By using a GOM, the effects of varying illumination (sun elevation and azimuth) become an advantage rather than an obstacle especially when the background is homogeneous. In terms of data, VHR image data such as a large proportion of Google Earth images have sufficient resolution for identi‐ **1. Task:** Tree crown detection in Very High Resolution images.

#### **2. Set parameters:**

**a. 3D model:** maximum reflectance, ambient light, sun elevation, sun azimuth, tree shape, clip factor.

**b. Descriptors of the objects:** minimum and maximum radius, minimum and maximum standard

deviation (δ) threshold,maximum and minimum similarity(sm) threshold, trees density.

**c. Change the process:** maximum iterations for decrease similarity

**d. Interruption of the process**: total iterations, maximum iterations without find new trees.

#### **3. Approach to tree crown detection:**

a. While the number of searched trees is not achieved or some of the interruption process (total iterations or minimum threshold for similarity).

#### **4. Starts the birth phase:**


d. If not:

i. Fill area with circle of radius *r*


#### **5. Starts the death phase:**

*3.1.3. A modified approach for orchards*

**Figure 9.** Flowchart of the MPP–TM algorithm.

58 Environmental Applications of Remote Sensing

Because trees in orchards are often individually distinguishable and have similar shape and size, they are perfect candidates for TM with a 3D geometrical-optical model. By using a GOM, the effects of varying illumination (sun elevation and azimuth) become an advantage rather than an obstacle especially when the background is homogeneous. In terms of data, VHR image data such as a large proportion of Google Earth images have sufficient resolution for identi‐

a. Input parameters: birth image matrix; crown statistics (*Sm* sorted); *Sm* threshold; tree models catalogue with radius between maximum and minimum radius

b. While smcrown < smthreshold:


**7. Update number of crowns eliminated for next birth phase**

**8. When the process finish:** reports, graphs and image with individual tree crowns

**Table 3.** Description of MPP-TM algorithm.

fying orchard trees. In this case, however, illumination parameters are not readily available and must be determined.

The objective of this modified approach is to introduce an adaptation of the algorithm described earlier to detect and count trees in orchards of different types. Because it was aimed at a more regional or even global application, Google Earth images were used in an attempt to simulate a generic operational framework. The modified approach uses a similarity measurement between the GOM and the image to calculate the probability of being the center of tree and then places trees in nonoverlapping positions (unless some overlapping is allowed). The algorithm also incorporates a module to determine the illumination parameters from a sample.

(c) (d) **Figure 10.** Illustration of the (a) birth and (b) death phases of the MPP–TM algorithm. In this example, of the 163 ran‐ domly positioned crowns, only one had a similarity value larger than the threshold of 0.98.

The algorithm is based on three principles. First, it assumes that the trees have a dome-like shape approximated with a GOM and the right illumination parameters. Second, there is little or no overlapping between trees, and third, the pixel with the highest similarity represents the most likely central position of the tree. (a) (b)

The GOM is a simple dome model for which the height is estimated at 1.5 times the diameter of the crown, and to simplify the problem we have assumed a unique diameter for all trees in the orchard (this can easily be modified to incorporate a range of diameters). The algorithm responsible for the detection of trees are best explained through a list of steps.

**Step 1.** Get user parameters: percent overlapping allowed minimum similarity value, tree diameter, and coordinate of sample tree. These parameters cannot be estimated automatically and are entered by the user. The illumination parameters can optionally be entered by the user, else they will be estimated by the program using the tree sample.

**Step 2.** If sun elevation and azimuth are not provided by the user, the parameters are auto‐ matically estimated by the program using the coordinates of a single-tree sample. The program then computes the similarity between the sample and all possibilities of illuminations param‐ eters in steps of 10 degrees.

**Step 3.** Calculate the similarity value for each pixel.

**Step 4.** Sort pixels by decreasing similarity and store coordinates. If the value is lower than the minimum allowed, the pixel is not stored.

**Step 5.** Place a temporary tree "stamp" (flat template) at the next pixel location with highest similarity value.

**Step 6.** Verify if space is already occupied by a tree. If some overlapping is allowed, make sure that the number of nonzero pixel is smaller than the percentage of overlapping allowed. An output image is created to receive a permanent "stamp" of the GOM shape with the value of similarity associated.

**Step 7.** Validate the results. Validation is performed by estimating the overall number of trees using the density of a representative sample and comparing with the number of trees found.

#### **3.2. Testing the MPP–TM approach**

#### *3.2.1. Urban trees*

The algorithm is based on three principles. First, it assumes that the trees have a dome-like shape approximated with a GOM and the right illumination parameters. Second, there is little or no overlapping between trees, and third, the pixel with the highest similarity represents the

(a) (b)

(a) (b) (c) (d)

**Figure 10.** Illustration of the (a) birth and (b) death phases of the MPP–TM algorithm. In this example, of the 163 ran‐

The GOM is a simple dome model for which the height is estimated at 1.5 times the diameter of the crown, and to simplify the problem we have assumed a unique diameter for all trees in the orchard (this can easily be modified to incorporate a range of diameters). The algorithm

**Step 1.** Get user parameters: percent overlapping allowed minimum similarity value, tree diameter, and coordinate of sample tree. These parameters cannot be estimated automatically and are entered by the user. The illumination parameters can optionally be entered by the user,

**Step 2.** If sun elevation and azimuth are not provided by the user, the parameters are auto‐ matically estimated by the program using the coordinates of a single-tree sample. The program then computes the similarity between the sample and all possibilities of illuminations param‐

**Step 4.** Sort pixels by decreasing similarity and store coordinates. If the value is lower than the

**Step 5.** Place a temporary tree "stamp" (flat template) at the next pixel location with highest

**Step 6.** Verify if space is already occupied by a tree. If some overlapping is allowed, make sure that the number of nonzero pixel is smaller than the percentage of overlapping allowed. An

responsible for the detection of trees are best explained through a list of steps.

domly positioned crowns, only one had a similarity value larger than the threshold of 0.98.

else they will be estimated by the program using the tree sample.

**Step 3.** Calculate the similarity value for each pixel.

minimum allowed, the pixel is not stored.

most likely central position of the tree.

60 Environmental Applications of Remote Sensing

eters in steps of 10 degrees.

similarity value.

Urban trees play an important role in the welfare and quality of life in cities. They contribute to improving air and water quality, mitigate the carbon dioxide and other pollutants, moderate the microclimate and air temperature, help control soil erosion, reduce the flow of rainwater, and provide biodiversity [37–39]. A good knowledge of the species planted in cities and their health contributes to the inventory and management of these trees. To fulfill their role in the urban environment, trees need to be looked after through maintenance practices such as pruning and monitoring them for pests and diseases.

A WorldView-2 (WV-2) image of the campus of the Universidade Federal de Minas Gerais (UFMG) (and surroundings) in Belo Horizonte, Brazil, was used as our test data (Figure 11). The scene was already orthorectified and radiometrically corrected. Although WV-2 offers nine different spectral bands, only the panchromatic band (*λ* ≈ 450–800 nm) with a ground resolution of 50 cm was used since all other bands have a ground resolution of 2 m.

Three WV-2 sub-images were selected to test the performance of MPP–TM algorithm (Figure 10). These images were chosen from different contexts with both isolated and grouped trees and with other objects present in the scene. A wide variety of crown radii is also present in these images. The first two images (Figure 12a and b) are from the university campus of UFMG, and the last is from an urban park (Figure 12c).

To assess the quality of the results produced by MPP–TM, validation was done by comparing our results with a visual interpretation of the trees in the image. For these, only tree counting was used as validation. For the crown counting validation, we considered the following situations: 1) true positives (TP) for found trees, 2) false positives (FP) when a detected object is not a tree, and 3) false negatives (FN) for trees not encountered. The success score was computed as follows (Equation 7):

$$\mathbf{A} = \left(\frac{\mathbf{T}\mathbf{D} - \mathbf{F}\mathbf{P}}{\mathbf{N} + \mathbf{F}\mathbf{N}}\right) \times 100\tag{7}$$

where TD represents the total detected trees and *N* is the total number of trees.

These results are shown below with their respective overall similarity and standard deviation graphs (Figure 13). The validation results are presented in Table 4.

**Figure 11.** Location of study area. The image on the right is a WorldView-2 false color composite


**Table 4.** Validation of the MPP-TM results with the three WV-2 images.

In the two images of the campus, the program was able to find 73% and 93% of the trees, respectively, with very few errors in isolated trees (Figure 13a and b). The presence of other objects (buildings, streets, and sidewalks) did not hinder the identification of trees and few false positives (3 and 8, respectively) were found. In both images, MPP–TM was able to find most grouped trees, but the crown diameter was often slightly off. It should be noted that some cases are even difficult to correctly identify and delineate visually. Mostly, the errors came from dividing a single crown into two, or including two different crowns as a single object.

The WV-2 image 3 is from a protected urban park area with predominantly isolated trees and relative homogeneous crown size of about 6 m (Figure 13 c). A total of 161 objects were detected with only 5 false positives and 20 false negatives for an overall success of 80%. Although most deciduous trees were selected, the crown size was often incorrect but given the highly irregular

**Figure 12.** Sub-images selected from the WV-2 image.

shape of many of these trees, this was somewhat expected, and similar problems have been reported by reference [4].

The behavior of the overall similarity during the iterations tend to increase as the image is progressively occupied by found trees and this is why the overall similarity increases. The standard deviation, however, is very different for each image and is mostly related to the amount of contrast in the original image. Images with highly contrasting objects (e.g., building tops) will tend to show a progressively decreasing standard deviation. Images of low contrast will tend to see it increasing as the trees are progressively added because of the double illumination nature of the trees.

#### *3.2.2. Orchards*

**Image Number of Trees Trees Detected False Positive False Negative Overall Accuracy**

**WV image 1** 47 43 3 8 72.73 **WV image 2** 50 59 8 5 92.73 **WV image 3** 175 161 5 20 80.00

**Figure 11.** Location of study area. The image on the right is a WorldView-2 false color composite

In the two images of the campus, the program was able to find 73% and 93% of the trees, respectively, with very few errors in isolated trees (Figure 13a and b). The presence of other objects (buildings, streets, and sidewalks) did not hinder the identification of trees and few false positives (3 and 8, respectively) were found. In both images, MPP–TM was able to find most grouped trees, but the crown diameter was often slightly off. It should be noted that some cases are even difficult to correctly identify and delineate visually. Mostly, the errors came from dividing a single crown into two, or including two different crowns as a single object.

The WV-2 image 3 is from a protected urban park area with predominantly isolated trees and relative homogeneous crown size of about 6 m (Figure 13 c). A total of 161 objects were detected with only 5 false positives and 20 false negatives for an overall success of 80%. Although most deciduous trees were selected, the crown size was often incorrect but given the highly irregular

**Table 4.** Validation of the MPP-TM results with the three WV-2 images.

62 Environmental Applications of Remote Sensing

**(%)**

Orchards are collections of individual trees often arranged regularly for which the MPP–TM algorithm could easily be adapted. Tree counting in orchards can be very useful for inventory and management purposes. For instance, the European Union (EU) Common Agricultural Policy (CAP) regulations (EC 73/2009) provide support for permanent crops such as hazelnuts, almonds, walnuts, and fruits in general [40–42]. Eligible orchards need to have a certain size and tree density depending on the type of crop. It has been estimated that orchard fruit production represents approximately 3– 4% of the total arable land [43], so the task of esti‐ mating fruit production needs tools for counting trees in a timely fashion. Furthermore, the task can take advantage of the near-global high-resolution image cover provided by Google (Google Earth and Google Map) and other Internet-based image services.

**Figure 13.** MPP–TM results obtained with the three WV-2 image 1–3 (left) and their graphs of global similarity (center) and standard deviation (right). The yellow circles correspond to correctly identified trees (true positive or TP), the ob‐ jects marked with a yellow "A" are false negatives (FN) and the objects marked with a yellow "B" are false positives (FP).

Orchards are plantation of trees of the same species and often of the same age. Consequently, trees of orchards usually have similar size and shape and are regularly spaced. Image proc‐ essing can easily be adapted to such a task providing VHR images are available. To illustrate the adapted MPP–TM algorithm (which no longer is a real MPP), we have tested over three different types of orchards: a mango plantation in Brazil near Juazeiro, a walnut plantation in France near Grenoble, and an olive plantation in Italy near Bracciano. The three images were directly extracted from Google Earth and had a relatively bad quality as they appeared to have been enhanced for sharpness. To validate the results, we have asked three geography students to manually interpret and mark the trees belonging to orchards for the three test images, and we have evaluated the results in the following way:


Orchards are plantation of trees of the same species and often of the same age. Consequently, trees of orchards usually have similar size and shape and are regularly spaced. Image proc‐ essing can easily be adapted to such a task providing VHR images are available. To illustrate the adapted MPP–TM algorithm (which no longer is a real MPP), we have tested over three different types of orchards: a mango plantation in Brazil near Juazeiro, a walnut plantation in France near Grenoble, and an olive plantation in Italy near Bracciano. The three images were directly extracted from Google Earth and had a relatively bad quality as they appeared to have been enhanced for sharpness. To validate the results, we have asked three geography students

(FP).

64 Environmental Applications of Remote Sensing

**Figure 13.** MPP–TM results obtained with the three WV-2 image 1–3 (left) and their graphs of global similarity (center) and standard deviation (right). The yellow circles correspond to correctly identified trees (true positive or TP), the ob‐ jects marked with a yellow "A" are false negatives (FN) and the objects marked with a yellow "B" are false positives

To be fair, the interpreters were told not to mark the trees that seem too small or too big for the orchards. In addition, valid trees that were found by the algorithm but did not pertain to an orchard were not computed as false positive. As a further improvement, restricting the search within the boundaries of the orchards would increase the accuracy and enable the similarity parameter to be relaxed. The addition of other spectral bands should also improve the results.


**Table 5.** Results of the tree counting algorithm for the three regions (France, Italy and Brazil).

Table 5 shows an overview of the results for the three test images, and Figure 14 shows the graphical results. The top row shows the original images, the center row shows the results of the tree identification (as well as false positives and negatives), and the bottom row displays a detailed section of the image on which the results were overlaid. The Grenoble test image (Figure 14 left column) was characterized by densely arranged walnut trees, which have a large round crown so that the model was well correlated with trees on the image, but the fact that the trees are close to one another produced a relatively large number of "miss" (103). This forced to relax the similarity threshold and caused a few false positives (69). In the case of the Bracciano image (Figure 14 center column), the olive trees are more ill- shaped than the walnut trees, and the relaxation of the similarity threshold caused a large number of false positives, especially in the nearby forested areas. Conversely, very few trees were missed. Finally, the last test image from Juazeiro (Figure 14 right column) is populated by mango trees that, like the walnut trees, have large round crowns. Still, the algorithm produced a fair amount of both false positives and false negatives mainly because of the variation of tree crown size and the particular situation of the dirt road at the top of the image that created a pattern of light and shade similar to the trees (approximately one-third of the false positives came from that road). The three very different images still produced similar accuracy results between 90 and 93%.

**Figure 14.** Illustration of the results of the tree counting for the three test images: Grenoble (left column), Briacciano (center column), and Juazeiro (right column). The empty circles represent the trees that were found, "*x*" represents the false positives and the black circles represent the false negatives.

#### **4. Final considerations**

The detection of individual tree crown in images of very high resolution is a growing and challenging field of research within the remote sensing community. In addition to the struc‐ tural complexity of the forest, many other factors such as the characteristics of the scene (topography, illumination, and other environmental variables) and forest type (season and biodiversity) make the task difficult. To reference [16], the ability to achieve individual tree crown delineation of all trees in a forest was recognized as an unrealistic expectation.

In an effort to provide the reader with an overview of the current state of the research in tree crown detection, Part I presented a brief review of some of the most common computerized techniques for detecting and delineating trees in optical VHR images. Part II describes the concepts and implementation of a novel approach based on two mathematical/pattern recognition concepts integrated to improve performance. MPP–TM was developed based on concepts from marked point processes and template matching for the former to take advantage of a mark built from a geometrical-optical model.

MPP–TM was highly effective in finding trees in urban environment with images from the WorldView-2 satellite (ground resolution of 50 cm). A total of 263 trees out of 272 were found (96%), and taking false positives into account, a success rate over 90% was still achieved. The algorithm was also adapted for a tree counting application such as is often needed in large orchards. To count trees in orchards, the approach works very well when the trees are easily distinguishable. Results from three datasets of different crops show an average success better than 90%. Out of 5806 trees, 5537 were found excluding all false positives.

The growing availability of VHR images from commercial satellites or even from web mapping services opens a wide field of applications especially that VHR multispectral images are becoming increasingly common. Multi-temporal studies will further strengthen these appli‐ cations for monitoring purposes.

Finally, we should mention that Lidar (light detection and ranging) data are also becoming widely available, and its integration with VHR images promises to further improve the results of tree detection algorithm. By adding a third dimension to the images, Lidar reduces the probability of errors by strengthening the evidence around the digital representation of trees.

### **Acknowledgements**

We are grateful to François Gougeon, Donald Leckie, Guillaume Perrin and Mats Erikson for having kindly provided the rights of reproduction of their figures.

### **Author details**

**Figure 14.** Illustration of the results of the tree counting for the three test images: Grenoble (left column), Briacciano (center column), and Juazeiro (right column). The empty circles represent the trees that were found, "*x*" represents the

The detection of individual tree crown in images of very high resolution is a growing and challenging field of research within the remote sensing community. In addition to the struc‐ tural complexity of the forest, many other factors such as the characteristics of the scene (topography, illumination, and other environmental variables) and forest type (season and

false positives and the black circles represent the false negatives.

**4. Final considerations**

66 Environmental Applications of Remote Sensing

Marilia Ferreira Gomes1,2\* and Philippe Maillard1

\*Address all correspondence to: mariliafgomes@yahoo.com

1 Geography Department, Geosciences Institute, Federal University of Minas Gerais, Belo Horizonte, Brazil

2 Cartography Department, National Institute of Land Reform, Belo Horizonte, Brazil

#### **References**


age Analysis - Spatial Concepts for Knowledge-Driven Remote Sensing Applications. 1st ed. Berlin: Springer; 2008. p. 1.4. DOI: 10.1007/978-3-540-77058-9

[13] Wulder M, Niemann KO, Goodenough DG. Local maximum filtering for the extrac‐ tion of tree locations and basal area from high spatial resolution imagery. Remote Sensing of Environment. 2000;73(1):103–114. DOI: 10.1016/S0034-4257(00)00101-2

**References**

68 Environmental Applications of Remote Sensing

[1] Pollock RJ. Model-based approach to automatically locating tree crowns in high spa‐ tial resolution images. In: Jacky Desachy, editor. Proceedings SPIE 2315: Image and Signal Processing for Remote Sensing; 30 December 1994; Rome. International Soci‐

[2] Avery TE, Berlin GL. Fundamentals of remote sensing and airphoto interpretation.

[3] Shao G, Reynolds KM, editors. Computer applications in sustainable forest manage‐ ment: Including perspectives on collaboration and integration. 1st ed. Dordrecht:

[4] Ke Y, Quackenbush LJ. A review of methods for automatic individual tree-crown de‐ tection and delineation from passive remote sensing. International Journal of Remote

[5] Perrin G. Etude du couvert forestier par processus ponctuels marqués [thesis]. Paris: Ecole Centrale Paris; 2006. 170 p. Available from: https://tel.archives-ouvertes.fr/

[6] Erikson M. Segmentation and classification of individual tree crowns in high spatial resolution aerial images [thesis]. Uppsala: Swedish University of Agricultural; 2004.

[7] Gougeon FA, Leckie DG. Forest information extraction from high spatial resolution images using an individual tree crown approach. 1st ed. Victoria: Canadian Forest

[8] Franklin SE. Remote sensing for sustainable forest management. 1st ed. New York:

[9] Blaschke T. Object based image analysis for remote sensing. ISPRS Journal of Photo‐ grammetry and Remote Sensing. 2010;65(1): 2–16. DOI: 10.1016/j.isprsjprs.2009.06.004

[10] Pu R, Landry S. A comparative analysis of high spatial resolution IKONOS and WorldView-2 imagery for mapping urban tree species. Remote Sensing of Environ‐

[11] Larsen M, Eriksson M, Descombes X, Perrin G, Brandtberg T, Gougeon F. Compari‐ son of six individual tree crown detection algorithms evaluated under varying forest conditions. International Journal of Remote Sensing. 2011;32(20):5827–5852. DOI:

[12] Hay GJ, Castilla G. Geographic Object-Based Image Analysis (GEOBIA): A new name for a new discipline?. In: Blaschke T, Lang S, Hay G, editors. Object-Based Im‐

ment. 2012;124:516–533. DOI: 10.1016/j.rse.2012.06.011

45 p. Available from: http://pub.epsilon.slu.se/676/[Accessed: 2015-04-12]

ety for Optics and Photonics; 1994. p. 526–537. DOI: 10.1117/12.196753

Springer Netherlands; 2006. 277 p. DOI: 10.1007/978-1-4020-4387-1

Sensing. 2011;32(17):4725–4747. DOI: 10.1080/01431161.2010.494184

5th ed. New Jersey: Prentice Hall; 1992. 472 p.

tel-00109074/ [Accessed: 2015-06-08]

Service; 2003. 27 p.

CRC Press; 2001. 424 p.

10.1080/01431161.2010.507790


[35] Baddeley AJ, Van Lieshout MNM. Stochastic geometry models in high-level vision. Journal of Applied Statistics. 1993;20(5-6):231–256. DOI: 10.1080/02664769300000065

[24] Hung C, Bryson M, Sukkarieh S. Multi-class predictive template for tree crown de‐ tection. ISPRS Journal of Photogrammetry and Remote Sensing. 2012;68:170–183.

[25] Gougeon F. A crown-following approach to the automatic delineation of individual tree crowns in high spatial resolution aerial images. Canadian Journal of Remote

[26] Leckie DG, Gougeon FA. An assessment of both visual and automated tree counting and species identification with high spatial resolution multispectral imagery. In: Hill DA, Leckie DG, editors. International Forum on Automated Interpretation of High Resolution Digital Imagery for Foresty; 10–12 February 1998; Victoria. Victoria: Cana‐

[27] Gougeon FA, Leckie DG. The individual tree crown approach applied to Ikonos im‐ ages of a coniferous plantation area. Photogrammetric Engineering & Remote Sens‐

[28] Szeliski R. Computer vision: algorithms and applications. 1st ed. London: Springer-

[29] Wang L, Gong P, Biging GS. Individual tree-crown delineation and treetop detection in high-spatial-resolution aerial imagery. Photogrammetric Engineering & Remote

[30] Lamar WR, McGraw JB, Warner TA. Multitemporal censusing of a population of eastern hemlock (Tsuga canadensis L.) from remotely sensed imagery using an auto‐ mated segmentation and reconciliation procedure. Remote Sensing of Environment.

[31] Bunting P, Lucas RM. The delineation of tree crowns in Australian mixed species for‐ ests using hyperspectral Compact Airborne Spectrographic Imager (CASI) data. Re‐ mote Sensing of Environment. 2006;101(2):230–248. DOI: 10.1016/j.rse.2005.12.015

[32] Zhou J. Application de l'identification d'objets sur images à l'étude de canopées de peuplements forestiers tropicaux: cas des plantations d'Eucalyptus et des mangroves [thesis]. Montpellier: Université Montpellier II - Sciences et Techniques du Langue‐

doc; 2012. 191 p. Available from: https://tel.archives-ouvertes.fr/tel-00763706/

[33] Descombes X. Méthodes stochastiques en analyse d'image: des champs de Markov aux processus ponctuels marqués [thesis]. Nice: Université de Nice Sophia-Antipolis; 2004. 225 p. Available from: https://tel.archives-ouvertes.fr/tel-00506084/[Accessed:

[34] Ortner M. Processus ponctuels marqués pour l'extraction automatique de caricatures de bâtiments à partir de modèles numériques d'élévation [thesis]. Nice: Université de Nice Sophia-Antipolis; 2004. 250 p. Available from: https://tel.archives-ouvertes.fr/

ing. 2006;72(11):1287–1297. DOI: 10.14358/PERS.72.11.1287

Sensing. 2004;70(3):351–357. DOI: 10.14358/PERS.70.3.351

Verlag; 2010. 812 p. DOI: 10.1007/978-1-84882-935-0

2005;94(1):133–143. DOI: 10.1016/j.rse.2004.09.003

2015-09-20]

tel-00189803 [Accessed: 2015-09-01]

DOI: 10.1016/j.isprsjprs.2012.01.009

dian Forest Service; 1998. p. 141–154.

Sensing. 1995;21(3):274–284.

70 Environmental Applications of Remote Sensing


**Remote Sensing Based Glacier Studies**
