LEAF AREA INDEX RETRIEVED FROM THERMAL HYPERSPECTRAL DATA

Leaf area index (LAI) is an important essential biodiversity variable due to its role in many terrestrial ecosystem processes such as evapotranspiration, energy balance, and gas exchanges as well as plant growth potential. A novel approach presented here is the retrieval of LAI using thermal infrared (8–14 μm, TIR) measurements. Here, we evaluate LAI retrieval using TIR hyperspectral data. Canopy emissivity spectral measurements were recorded under controlled laboratory conditions using a MIDAC (M4401-F) illuminator Fourier Transform Infrared spectrometer for two plant species during which LAI was destructively measured. The accuracy of retrieval for LAI was then assessed using partial least square regression (PLSR) and narrow band index calculated in the form of normalized difference index from all possible combinations of wavebands. The obtained accuracy from the PLSR for LAI retrieval was relatively higher than narrow-band vegetation index (0.54<R2<0.74). The results demonstrated that LAI may successfully be estimated from hyperspectral thermal data. The study highlights the potential of hyperspectral thermal data for retrieval of vegetation biophysical variables at the canopy level for the first time.


INTRODUCTION
Leaf area index (LAI) is a principal component of biogeochemical cycles in ecosystems (Bréda 2003).Among vegetation biophysical properties, LAI is of particular interest, as it exhibits significant control on the transpiration, respiration, and gas exchanges (e.g., uptake of CO2 and H2O by the canopy) between terrestrial ecosystems and atmosphere.Previous studies have shown the importance of LAI in ecological and remote sensing studies.LAI is a key input for climate and large-scale ecosystem models and also is a key structural characteristic of forest ecosystems (Chen et al. 1997;Myneni et al. 1997;Wang et al. 2004;Zheng and Moskal 2009).In the last decades, LAI has been successfully retrieved using hyperspectral data in the visible/near-infrared (0.35-1.0 µm, VNIR) and short-wave infrared (1.0-2.5 µm, SWIR) regions (Zheng and Moskal 2009).Despite the broadly recognized importance of LAI across ecological research, to our knowledge, LAI has not estimated from thermal infrared (8-14 µm, TIR) hyperspectral data.TIR hyperspectral data deserves the same exploration and development of methods, as hyperspectral data in the VNIR and SWIR regions.
Recently, Ullah (2013) showed that TIR hyperspectral data is supplementary to other remote sensing data and has the potential to explain the biochemical characteristics of vegetation (e.g.water content) at leaf level.The primary absorption features associated with water and cellulose, as important vegetation components, are only observable in the mid-wave infrared (3-5 µm, MIR) and TIR regions (Fabre et al. 2011;Gerber et al. 2011;Ribeiro da Luz 2006).In addition, previously, it has been presumed that plants are opaque and featureless in the TIR region.However, recent studies have uncovered that various plant species display distinct emissivity spectra and have recognizable spectral features in the TIR region (Ribeiro da Luz and Crowley 2010; Ullah et al. 2012).There are a number of reasons why to date limited attention was directed at using TIR hyperspectral remote sensing data for vegetation studies.In the TIR region, spectra result from the emissivity of surfaces, rather than from reflectance.Only a few instruments measure TIR emissivity spectra at high spectral resolution, and only very few are able to do this for complex surfaces such as canopies.Therefore, few studies have been conducted on TIR hyperspectral data for vegetation studies at canopy level (Ribeiro da Luz and Crowley 2010; Sepulcre-Cantó et al. 2006).More studies are required to assess and understand TIR remotely measured spectra from vegetation particularly at canopy level.
LAI has been predicted in many studies using hyperspectral data from the VNIR and SWIR regions using univariate (e.g.vegetation indices) and multivariate (e.g.partial least squares regression (PLSR)) models.This has led to the development of new vegetation indices which have further improved LAI prediction (Asner and Martin 2008;Baret and Guyot 1991;Darvishzadeh et al. 2009;Eriksson et al. 2006;Gao et al. 2000;Haboudane et al. 2004;Koetz et al. 2005).To our knowledge, neither the relation between LAI and univariate methods nor the relation between LAI and multivariate methods in the TIR region has yet been demonstrated.Therefore, the main objective of this study was to investigate the retrieval of LAI in the TIR region, using narrow band indices and PLSR.

LAI measurements
In the present study, two plant species were selected: Azalea japonica (n=10) and Ficus benjamina (n=6).The leaves from the plants were harvested in several steps and LAI (m 2 m −2 ) was calculated using the measured surface area of the leaves and the corresponding ground area of the canopy.The dataset included 16 plant specimens, which yielded 60 LAI estimates from the destructive sampling.

Canopy thermal infrared radiance measurement
The radiance spectral measurements were made using a portable MIDAC illuminator Fourier Transform Infrared (FTIR) spectrometer (Model M4401-F; MIDAC Corporation, CA, USA).Measurements were made at nadir position above the samples.To create optimal measurement conditions and reduce any possible source of errors due to the changes in atmospheric conditions or temperature, the measurements were made under controlled laboratory conditions in which the walls, ceiling, and floor were coated with a black material (Avis Aqua Blackboard Black) and black plastic of known emissivity.We reduced the laboratory temperature to 10˚C in order to generate a suitable thermal contrast with the plants, which were at a higher room temperature.In addition, plants were kept outside the laboratory at an ambient room temperature of 20°C and each one was briefly introduced to the laboratory in order to make the thermal measurements.Measurements were made with a fixed vertical distance between sensor and sample (60 cm).In this experiment, the background soil was covered with black plastic of known emissivity to minimize possible effects of soil.The radiance spectra of the plant canopies were measured between wavelengths of 2.5-20 μm at a resolution of two cm −1 .The emissivity spectra of plant canopies at each LAI value were obtained using a series of FTIR measurements performed in the following order: radiance measurements of the hot blackbody, radiance measurement of the cold blackbody, radiance measurement of the sample (i.e., the plant canopy with specific LAI value), and finally, radiance measurements of a highly diffuse reflecting gold plate (Infragold®).For instrument radiance calibration, two blackbodies (one hot, one cold) were used.The cold blackbody temperature was set just below the ambient temperature, at 5°C (Korb et al. 1996).The hot blackbody temperature was set above the sample temperature (Hori et al. 2006;Salvaggio and Miller 2001), at 30°C.Details of the radiometric calibration applied with the measurements of these blackbodies' radiances can be found in (Hook and Kahle 1996).A diffuse reflecting gold plate with an emissivity of ~ 0.04 was used to collect downwelling radiance (DWR) in order to determine the quantitative influence of the laboratory background selfemission (Eisele et al. 2015).The measurement series were taken within five minutes to minimize the possible drift of the instrument, physiological changes in the plants, and fluctuations in laboratory temperature.After each set of measurements, the canopy was rotated 90 degrees clockwise.The final corresponding canopy emissivity spectra of each sample (for a specific LAI value), was then calculated from the average of four sets of measurements (covering 360 degrees).In total, 240 (4*60) canopy radiance measurements were obtained for two plant species.Spectral emissivity of the plants was calculated from their absolute radiance using the following equation (Korb et al. 1996), Where ε sam (λ) = Directional emissivity of the sample at the wavelength λ L sam (λ) = Spectral radiance from the sample T sam = The true physical temperature of the sample B(λ, T sam ) = Planck function at the wavelength λ and sample temperature L DWR (λ) = Total spectral DWR from the hemisphere above the sample.
Tsam must accurately represent the emission temperature.Therefore, after measuring canopy temperature before and after each measurement, the blackbody fit method was used to calculate the sample emissivity value.Details about the blackbody fit method can be found in Kahle and Alley (1992) and Salvaggio and Miller (2001).The emissivity of 279 wavebands between 8 μm and 14 μm was selected for analysis; measurements outside this range had a very low signal strength.A Savitzky-Golay filter with a frame size of 15 data points and second-degree polynomial was used to reduce the noise of the canopy emissivity spectra (Savitzky and Golay 1964).

Narrow band Index
Generally, vegetation indices which are calculated from different combinations of bands are divided into two main categories: orthogonal and ratio indices (Baret and Guyot 1991).Of these vegetation indices, the normalized difference index used frequently to derive LAI from optical remote sensing data (Huete et al. 2002).Therefore, in this study, this widely used index was calculated from all the possible paired waveband combinations, using canopy emissivity spectra of pooled data and each species between 8 μm and 14 μm.To evaluate the strength of the relation between proposed index and LAI, the coefficients of determination (R 2 ) of all possible two-waveband combinations of vegetation index and LAI were used.From all these combinations, we selected the narrow band indices that performed best at estimating LAI, as having the maximum R 2 with LAI (Mutanga and Skidmore 2004).The naming conventions, acronyms, and calculations for the considered index is listed in Table 1.In this study, the concept of calculating normalized difference index is based upon the contrast in the emissivity between two different spectra wavebands (R λ1 and R λ2 ).
Table 1.The spectral index used in this paper.

Name Acronym TIR Equation
Vegetation indices can be related to LAI through linear or exponential regression models (Schlerf et al. 2005).Some studies have confirmed that vegetation indices have decreasing sensitivity to rising LAI values when optical multispectral and hyperspectral data are used (Baret and Guyot 1991;Broge and Mortensen 2002;Chen and Cihlar 1996;Turner et al. 1999).Darvishzadeh et al. (2009) showed that linear and exponential regression models are similarly accurate at estimating LAI from hyperspectral data.In our study, a linear regression model was used to model the relationship between narrow band index as the predictor variable and LAI.In addition, a cross-validation procedure also known as the "leave one out method" (Duda and Hart 1973), used to validate the regression models.The estimated LAI was assessed using the crossvalidated coefficient of determination (R 2 CV), and crossvalidated root mean squared error (RMSECV) for the best performing narrow band index.

Partial least squares regression
Hyperspectral data has a high dimensionality and high collinearity between adjacent wavebands.PLSR is a popular multivariate statistical technique for transforming the wavebands to new orthogonal factors, i.e. components, which eliminates this collinearity.Multivariate methods include models such as PLSR which as so-called "full spectrum methods" employ all available spectral data simultaneously (Atzberger et al. 2010) and reduce the effects of multicollinearity as a common problem inherent to the hyperspectral dataset (Mirzaie et al. 2014).
PLSR has been broadly applied in remote sensing vegetation analysis (Cho et al. 2007;Kooistra et al. 2004).Here, PLSR analysis was used to determine the relative contribution of LAI to the 279-waveband canopy emissivity spectra for all the sampled species.The canopy emissivity spectra as independent variables were mean-centered before performing the PLSR analysis.To avoid over-fitting and prevent collinearity, an optimum number of factors was determined, using a crossvalidation procedure.The criterion for adding an extra factor to the model was that it decreased RMSEcv by > 2% (Geladi and Kowalski 1986).The accuracy of the models was evaluated using RMSECV.All PLSR analyses were carried out using the TOMCAT toolbox 1.01 within MATLAB (Daszykowski et al. 2007).

Narrow band indices and leaf area index
The descriptive statistics of our experimental setup showed LAI ranges between 0.60 (m 2 m −2 ) and 8.36 (m 2 m −2 ).Together with the variation in LAI values, the measured canopy emissivity exhibited a large variability (Table 2).
Table 2. Summary statistics of the leaf area index (LAI) measurements for four plant species (n=60).Desired narrow band index (i.e.ND) was calculated from the measured canopy emissivity spectra, using all possible twowaveband combinations.The coefficients of determination (R 2 ) between ND and LAI was computed and tabulated in Table 3. Table 3 demonstrates the maximum R 2 between the ND and the LAI, for each species and pooled data, as well as the band position of the best performing narrow band index.The maximum R 2 value obtained between LAI and ND for Azalea japonica and Ficus benjamina were 0.65 and 0.48 respectively.

Species name LAI (m 2 m -2 ) Mean Max Min Sample size
The results showed that the maximum R 2 with LAI across species was yielded by wavebands in the 9.9-10.6µm range in combination with the bands in the 10.2-11.9µm range.In general, LAI was estimated with appropriate accuracy.An illustration of these results is depicted in the 2-D correlation plot (Figure 4).The meeting point of each pair of wavebands in the 2-D plots corresponded to the R 2 value of LAI and the ND calculated from the emissivity values in those twowavebands.3.

Partial least squares regression and leaf area index
The PLSR regression using TIR spectra yielded the maximum R 2 0.68 and 0.54, in Azalea japonica and Ficus benjamina, respectively.The optimal number of factors modelled using PLSR ranged from three to five.(Table 6).For both species, the PLSR models in the TIR region yielded better predictions of LAI in comparison with the narrow band index (i.e.resulted in a higher R 2 and lower relative RMSECV).For instance, the improvement was noticed for Azalea japonica where R 2 enhanced from 0.48 to 0.54 and RMSECV reduced from 0.43 to 0.30.Cross-validation estimation of LAI using the entire emissivity spectra in PLSR models is depicted in Figure 7.

CONCLUSIONS
It has long been thought that biophysical or biochemical properties of plants cannot be retrieved using emissivity spectra in the TIR region.This may partly be contributed to the complexity of measuring vegetation biophysical and biochemical variables, as well as the lack of suitable instruments to measure emissivity in the TIR region.Our results demonstrate that biophysical properties of vegetation at the canopy level, such as LAI, can be retrieved using canopy emissivity spectra in the TIR region.Univariate and multivariate techniques can be applied to retrieve LAI values using emissivity spectra.These results are comparable to the results of previous studies in the optical domain for LAI estimation (Darvishzadeh et al. 2008a;Schlerf et al. 2005).
The results of our experiments demonstrate that vegetation indices enable to accurately predict LAI for single species in the TIR domain rather than pooled data including two plant species.It can be attributed to differences in canopy structure (i.e., leaf size and orientation), and probably by the biochemical concentrations (i.e., cellulose and cutin) (Buitrago et al. 2016) of the leaves of these species.As can be seen from the Table 3, comparing the best performing narrow band index location among species showed that they were located in the slightly different spectral region.This can again be explained by differences in the biochemical concentration of leaves and canopy structural properties of the species studied which could effect on emissivity spectra and radiance scattering in the TIR region (Ribeiro da Luz and Crowley 2007;Salisbury 1986).When the PLSR was applied, the retrieved LAI had a relatively high accuracy compared with results obtained using narrow band index (Table 3).The results highlight that vegetation narrow band index utilize a limited number of spectral data bands (i.e.here only information from two-waveband combinations) from the massive spectral content of the total emissivity information available, whereas the PLSR technique employed all available spectral data.Our results are in a line with the Darvishzadeh et al. (2008b) finding which showed the PLSR technique improved the prediction of LAI compared with narrow band indices in the optimal domain.In all the models considered, the relationship between estimated and measured LAI was found to be linear.This study revealed that LAI could be estimated reasonably accurately in the TIR region.Therefore, we can conclude that LAI can be successfully be retrieved from TIR hyperspectral data, even at relatively high values of LAI.This study demonstrates TIR hyperspectral data can be expected to complement other remote sensing data and also it proof of concept for measuring LAI using the TIR, and brings a following interesting challenge: how to upscale the results from the laboratory (plant canopy) to field conditions.

Figure 5 .
Figure 5. Scatterplot of measured versus predicted leaf area index for the best narrow band index calculated from normalized difference index (ND) for Ficus benjamina (a), Azalea japonica (b), and pooled data (c).The optimum wavebands are those reported in Table3.

Figure 7 .
Scatterplot of measured versus estimated leaf area index using entire emissivity spectra in the partial least squares regression model from TIR region: Ficus benjamina (a), Azalea japonica (b), and pooled data (c).

Table 3 .
Band position and the coefficients of determination (R 2 ) values between the best performing narrow band index and LAI.R 2 CV is the cross-validated coefficient of determination between estimated and measured LAI; RMSECV is the relative root mean square error.

Table 6 .
The performance of partial least squares regression for estimating leaf area index.