Empirical and Physical Estimation of Canopy Water Content from CHRIS/PROBA Data

Abstract

Efficient monitoring of Canopy Water Content (CWC) is a central feature in vegetation studies. The potential of hyperspectral high spatial resolution CHRIS/PROBA satellite data for the retrieval of CWC was here investigated using empirical and physical based approaches. Special attention was paid to the spectral band selection, inversion technique and training process. Performances were evaluated with ground measurements from the SEN3EXP field campaign over a range of crops. Results showed that the optimal band selection includes four spectral bands: one centered about 970 nm absorption feature which is sensible to C_w, and three bands in green, red and near infrared to estimate LAI and compensate from leaf- and canopy-level effects. A simple neural network with a single hidden layer of five tangent sigmoid transfer functions trained over PROSAIL radiative transfer simulations showed benefits in the retrieval performances compared with a look up table inversion approach (root mean square error of 0.16 kg/m² vs. 0.22 kg/m²). The neural network inversion approach showed a good agreement and performances similar to an empirical up-scaling approach based on a multivariate iteratively re-weighted least squares algorithm, demonstrating the applicability of radiative transfer model inversion methods to CHRIS/PROBA for high spatial resolution monitoring of CWC.

1. Introduction

Canopy Water Content (CWC), defined as the mass of water per unit ground area, is a key biophysical parameter in agricultural and forestry applications. It is required for monitoring drought conditions [1], assessing forest fire susceptibility [2,3], or improving soil moisture retrievals [4], among many others applications. High spatial resolution sensor may allow detecting the within-field spatial variability of CWC useful in agricultural management [5].

Several studies have demonstrated the existing link between canopy reflectance and CWC. The estimation of CWC from optical remote sensing generally relies on the water absorption features centered at 970 nm, 1,200 nm, 1,450 nm and 1,950 nm [6]. Clevers et al. [7,8] demonstrated the potential of using water absorption band at 970 nm for the retrieval of CWC. Peñuelas et al. [9] focused on the 950–970 nm region and defined the so-called water band index (WI) as the ratio between the reflectance at 970 nm and the one at 900 nm (as a reference wavelength). However, the absorption band of atmospheric water vapor at about 940 nm should be considered when using satellite data around 970 nm. Interference with atmospheric absorption can be avoided by using contiguous spectral coverage in the near infrared region with hyperspectral data. There are many studies in the literature showing the applicability of hyperspectral data for the estimation of CWC concurrent with water vapor from absorption feature at 970 nm [7,10–12]. However, very few studies have been devoted to investigate the potential of using CHRIS/PROBA satellite data for water content retrieval.

CHRIS/PROBA is a hyperspectral satellite sensor with high spatial resolution (34 m), which combines a multi-viewing capability in five nominal angles (−55°, −36°, 0°, +36°, +55°) and sixty two bands from 405 to 1,005 nm. Estimating CWC from CHRIS/PROBA appears challenging due to the limited spectral range of the sensor (405–1,005 nm) and the weaker water absorption feature at 970 nm as compared to the main water absorption bands (at >1,200 nm), which are out of CHRIS spectral range. Stagakis et al. [13] and Sykioti et al. [14] analyzed a two-year period of CHRIS/PROBA observations and assessed the potential of 900–1,000 nm spectral region for leaf water potential monitoring over a Mediterranean ecosystem dominated by the semi-deciduous shrub Phlomis fruticosa. Sykioti et al. [14] found that in the 900–1,000 nm spectral area, both CHRIS and leaf spectra showed a “weak” (band depth < 0.1) but constant in time absorption feature centered around 970 nm with significant correlation (R² = 0.81) with leaf water potential, only for a viewing zenith angle of −36°. Stagakis et al. [13] analyzed different CHRIS band combinations and found that the spectral indices with significant correlation with leaf water potential were actually linked to chlorophyll variations while specific indices designed for water retrieval such as WI [9] showed poor correlation. These findings may be interpreted based on the particularities of the studied species physiology as authors suggested [13] but also the intense covariance between leaf water content and other biochemical and structural leaf- and canopy-level parameters [15]. Several combinations of leaf and canopy variables may provide very similar reflectance responses and radiometric information may be not sufficient to identify a unique solution of variables: ill-posed nature of the inverse problem in remote sensing [16]. Regularization of the inversion problem to get more stable and accurate solutions requires introducing prior information on the distribution of the variables [16]. Since leaf water content (C_w) and the leaf area index (LAI) are strongly correlated [17], different combinations of C_w and LAI values may have very similar values of C_w × LAI product. Combal et al. [16] showed that C_w × LAI is more accurately estimated than C_w. Recognizing the challenge of estimating leaf water content, we propose the retrieval at the canopy level of CWC corresponding to C_w × LAI product. This choice is consistent with previous studies [18] including the proposed algorithm for CWC retrieval for the Sentinel 2 mission [19].

A wide variety of methods have been proposed to estimate CWC from optical remote sensing. Much of the effort has been to develop multispectral indices and establish empirical relationships based on regression analysis with water content [4,5,7–9,12]. These relationships are, however, side and species specific and their applicability is limited to calibration conditions. Recent research efforts focus on the application of radiative transfer models and model inversion techniques for leaf and canopy water content estimation from satellite reflectance imagery [18]. To invert canopy reflectance models, a range of techniques have been proposed from numerical inversion procedures [18,20,21] to look up tables [22,23] or artificial neural networks [24,25]. Physical approaches are, however, limited by several aspects linked not only to the realism of the radiative transfer model but also to the inversion methods themselves [26].

Recently, Trombetti et al. [25] used neural networks trained over radiative transfer simulations for monitoring CWC for the continental USA on a monthly basis from MODIS data. Neural networks have been increasingly used for large datasets processing due to their (i) computational efficiency and ability to (ii) learn complex pattern, taking into account any nonlinear complex relationship between the variables, (iii) generalize in noisy environments, which makes networks robust estimators in the presence of incomplete or imprecise data and (iv) incorporate a priori knowledge and realistic physical constraints into the analysis [27]. Verger et al. [28] analyzed the optimal modalities for radiative transfer-neural network estimation of LAI, fraction of green vegetation cover (FCOVER) and fraction of absorbed photosynthetically active radiation (FAPAR) with due attention to the impact of the architecture of the network, the training dataset and the inversion method to extract the solution. This study extends our previous research by evaluating if the principles of neural network inversion established in Verger et al. [28] for LAI, FAPAR and FCOVER variables are still valid for CWC.

This investigation studies the potential of CHRIS/PROBA satellite data for the retrieval of CWC at high spatial resolution. Special attention was paid to the spectral band selection and inversion technique. The applicability of empirical transfer functions and physical radiative transfer inversion based on neural networks and look up tables was investigated. The PROSAIL [29] radiative transfer model is used. The impact of the training dataset and the influence of the inversion modalities for optimal extraction of the solution are analyzed. The performance of the estimates is assessed through the comparison with ground data acquired in a cropland area during the Sentinel-3 Experiment (SEN3EXP) field campaign [30].

2. Study Area and Data

2.1. Study Area

The Barrax agricultural site (39°30′N, 2°6′W) is located in La Mancha, a plateau 700 m above sea level in south-east Spain. The area is characterized by a large variety of uniform land use units of different crops and dry bare soils. The dominant cultivation pattern is approximately 65% dry land with winter/spring cereals and bare soils/fallow/land. One third of the land (35%) is irrigated, comprising alfalfa, maize, potatoes, sunflower, onion, garlic, sugar beet and vineyard, among other crops.

2.2. Field Campaign

The SEN3EXP experimental campaign took place in Barrax between 20 June 2009 and 24 June 2009 as a part of the preparatory activities for the European Space Agency (ESA) Sentinel-3 mission [30]. Spaceborne and airborne data combined with co-located in situ information were acquired for algorithm prototype and ground segment processor development. In situ vegetation measurements included characterization of Leaf Area Index (LAI), Fresh Weight (FW), Dry Weight (DW) and leaf Area (A). A total number of 38 statistically representative Elementary Sampling Units (ESUs) of 20 × 20 m² were sampled, corresponding to 12 different fields and 7 vegetation classes: vineyards, fruits, sunflowers, alfalfa, corn, garlic and onion. Additionally, the location of 11 bare soil control areas was identified. The position of ESUs is shown in Figure 1.

LAI was estimated from the directional variation of the gap fraction measured by using upward looking digital hemispheric photography [31]. Between 12 and 15 photos per ESU were acquired. Images were processed using the CAN-EYE software [32]. Photography processing provides Plant Area Index (PAI) rather than LAI since all green or non green vegetation elements were accounted for. However, most of the crops in the study area were fully green and measured PAI is a good approximation of the LAI. For the sake of simplicity, the measured PAI will be termed LAI.

Destructive samples were conducted to estimate water content. Samples included approximately 95% of leaves and 5% of stems. The area of samples (A) varies from 0.01 to 0.25 m². Samples were weighted before and after being dried at 105 °C during 24 h to retrieve FW and DW, respectively. Further details on the field protocols and data acquisition are provided in Camacho et al. [33]. Leaf water content (C_w), being the amount of water per unit leaf area, was estimated from FW, DW and A:

$${\text{C}}_{\text{w}}=(\text{FW}-\text{DW})/\text{A}$$

(1)

CWC is defined as the product of leaf water content (C_w) and the leaf area per unit of ground area (i.e., LAI):

$$\text{CWC}={\text{C}}_{\text{w}}\times \text{LAI}$$

(2)

Note that the measured CWC refers to the total mass of liquid water content in foliage and stems per ground area (kg/m²) since both C_w and LAI measurements account for all the vegetation elements. In situ CWC values at Barrax area (Figure 2) range from 0.05 kg/m² for vineyards and fruit trees to 1.2 kg/m² for garlic crops.

2.3. Satellite Data

CHRIS/PROBA high spatial resolution hyperspectral imagery were acquired over Barrax on 19 June 2009 at 10:07 AM, under partly cloudy conditions [30]. A quasi-nadir image was here considered with 13.31° view zenith angle, 138.08° observation azimuth angle and 30° solar zenith angle. Data corresponds to mode 1 of acquisition characterized by 62 spectral bands ranging from 405 nm to 1,005 nm with a nadir ground sampling distance of 34 m. The image of 34 m pixel size was re-sampled to 17 m and projected to UTM 30 North WGS-84 using 50 ground control points ensuring about 13 m geometrical accuracy [34]. The pixel size (17 m) was selected in order to preserve spectral signatures and avoid resampling artifacts (e.g., aliasing) during the orthorectification of the image (the orientation of the satellite image is different to the UTM grid) by nearest neighbor interpolation [35]. Note that the pixel size of the CHRIS/PROBA re-sampled data approximately coincides with the size of the ground ESUs (20 × 20 m²). The BEAM [36] toolbox specific for CHRIS/PROBA was used for atmospheric correction.

3. Canopy Water Content Modeling and Retrieval Approaches

3.1. PROSAIL Model

The widely used PROSAIL radiative transfer model [29], which is a combination of the leaf PROSPECT [37] and the canopy SAIL [38,39] reflectance models, was selected in this study. The bidirectional top of the canopy reflectance for the viewing and illumination geometry of CHRIS/PROBA acquisition is simulated by SAIL model as a function of three structural parameters: the LAI, the average leaf angle inclination (ALA), and the hot spot size parameter (H) [40]. The leaf reflectance and transmittance are simulated by PROSPECT as a function of four structural and biochemical parameters: the mesophyll structure parameter (N), leaf chlorophyll concentration (C_ab), leaf water content (C_w) and the leaf dry matter content (C_m). The soil spectral reflectance was here simulated using the soil reflectance spectra measured in Barrax multiplied by a brightness coefficient B_s allowing to represent the variability induced by soil roughness and moisture [41].

Prior information available on Barrax site was taken into account for defining the canopy and leaf parameters of PROSAIL model and reducing the size of the parameter space. The ranges of model input variables (

Table 1. Statistics of the input variables of PROSAIL model used for the simulations.

	Model Variables	N class	Min	Max	Mean	Std. Dev.
Leaf	N	4	1.0	2.5	1.5	1.0
	Cab (μg/cm2)	4	20	50	33	5
	Cw (g/cm2)	4	0	0.08	0.03	0.02
	Cm (g/cm2)	4	0.003	0.02	0.01	0.003
Canopy	LAI (m2/m2)	6	0	8	2	2
	ALA (°)	4	30	80	50	10
	Hot	1	0.001	1	0.1	0.3
Soil	Bs	4	0.7	2.3	1.4	0.3

) were established according to the ground measurements from the SEN3EXP campaign including C_ab, C_w, C_m, LAI, ALA and B_s. The range of the other variables (N and H) was set according to the values reported in the literature in studies over the Barrax study area [28,42].

Model parameters were randomly sampled within their ranges (Table 1) applying a stratified sampling scheme. The whole range of variation of each variable was split into a small number of classes (Table 1) to ensure that values from each class were combined with values from each other variable class [28]. This results in a number of around 100,000 parameter combinations. To reduce the number of parameter combination, the final training dataset was made of 10,000 cases randomly selected. This size was found to be a good compromise between computer resource requirements and the accuracy of the estimates [28,42].

3.2. Spectral Sensitivity Analysis and Band Selection

Selection of a reduced set of optimal spectral bands may be preferable for operational perspectives [42]. Selection of spectral bands is based on the information content of the bands but also on technical aspects in CHRIS/PROBA data processing in order to avoid noisy bands. A spectral sensitivity analysis using PROSAIL model was first conducted for selecting, among the 62 bands of CHRIS sensor, the optimal spectral regions for CWC retrieval. Second, an uncertainty assessment of CHRIS reflectances based on the spectral stability to signal noise was done to select the best band inside each spectral region.

• The PROSAIL model was run in forward mode and the canopy reflectance was computed by fixing the different parameters to their mean values (Table 1) and varying the two components of CWC, i.e., leaf water content, C_w, and LAI, between their ranges of variation in the study area. The sensitivity of canopy reflectance to C_w and LAI variations in the optical spectra (400–2,400 nm) is illustrated in Figure 3. The maximum sensitivity to C_w (Figure 3a) is observed in the near infrared (NIR) and middle infrared with water absorption features centered around 970 nm, 1,200 nm, 1,450 nm and 1,950 nm [24]. Only the spectral information around 970 nm can be exploited for the estimation of C_w from CHRIS data due to the limited spectral range of the sensor (405–1,005 nm). The last four spectral bands (from band b59 to b62) of CHRIS ranging from 960 to 1,005 nm appear to be the optimal spectral domain for C_w retrieval. In the case of LAI (Figure 3b), the sensitivity analysis confirms that red (bands b21–b25 ranging from 627 to 677 nm) and NIR (bands b41–b52 ranging from 773 to 891 nm) are the domains of major interest for the estimation of LAI. Several studies have demonstrated that the combination of two bands in NIR (high sensitivity to LAI) and red (used as a reference band to minimize the influence of soil background) spectral domains is optimal for LAI retrieval [28]. An additional band in the green region (bands b11–b14 ranging from 526 to 566 nm) provides complementary information and a wider range of reflectance sensitivity to LAI (Figure 3b).
• A coefficient of variability (CV), defined here as the ratio between the standard deviation of reflectance measurements and their mean value expressed in percentage, was computed over a homogeneous surface as an indicator of the stability of reflectance measurements to signal noise. Results (Figure 4) show that first CHRIS bands located in the blue spectral domain are significantly affected by residual atmospheric effects (non-physical negative values of reflectances and high CVs). Residual noise also affects the CHRIS bands in the red domain.

According to the two complementary criteria (i.e., spectral sensitivity to the C_w and LAI variables, and stability to signal noise) bands b13 (centered at 551.9 nm), b22 (641.8 nm), b50 (869.1 nm) and b61 (988.3 nm) were finally selected.

3.3. Neural Network Inversion Approach

Neural networks have been broadly employed in remote sensing to map land variables from reflectance model inversion [28,43,44]. The reflectance model is used for generating a synthetic dataset, which is used to train the networks and establish non-linear relationships between the input reflectances and the corresponding canopy variables to be retrieved. To represent the relationship between input and output variables, networks are composed on n layers of neurons that transform the spectral signal into biophysical variables by synaptic weight and bias added to the neurons. Neural networks are mainly defined by their architecture, the training dataset and the learning algorithm:

• The architecture of the networks was made of one input layer with as much neurons as the number of inputs. The number of hidden layers and the number of neurons per layer was empirically defined by selecting the optimal values. Based on a comparison of different activation functions in literature [28,45,46] three combinations were tested by considering hyperbolic tangent in hidden layers, and hyperbolic tangent, linear or saturated linear in the output layer.
• The training process consists in adjusting the networks coefficients by minimizing a cost function using a back propagation algorithm [47]. The selected cost function was here defined as the root mean square error between the targeted variable in the simulation dataset and the network output. The Levenberg-Marquardt minimization algorithm was used because of its efficient convergence performance. For generating the training dataset, truncated Gaussian distributions that mimic the actual distribution of the radiative transfer model input variables based on prior knowledge of the study area (Table 1) were considered. Verger et al. [28] demonstrated that for neural networks training such Gaussian distributions of model parameters performs better than uniform distributions for which no prior information is exploited. Including moderate uncertainties in the reflectance simulations used in the training process improves the flexibility of the neural networks in cases where simulations slightly depart from observations [28]. Three different Gaussian-white noise levels (2, 4 and 6%) were added to the simulated reflectances in order to include instrumental noise and radiometric and atmospheric uncertainties.
• The training datasets were randomly split into three subsets [44]: one half of the simulated cases were used to train the network, one fourth to avoid over-specializations during the training process and one fourth to test the performance of the network and select the solution. The solution was finally extracted by training 10 parallel networks to select the one providing the best performance over the test dataset.

The several modalities used within the neural networks inversion were investigated by evaluating the performance of CWC estimates through the comparison with ground measurements (Section 4.1). Attention was paid to the uncertainties in the simulated reflectances of the training dataset, architecture of the neural networks and type of the activation functions. Rather than analyzing all the combinations, two of the three analyzed variables were fixed to the best set of modalities from which each modality of the third one was sequentially investigated. The overall performance of CWC estimates was quantified by the root mean square error (RMSE) decomposed into accuracy (B) and precision (S) components. The accuracy (B) was measured as the mean value of the differences between products and ground measurements while the precision (S) was computed as the standard deviation of estimates around the best linear fit also reflected by the correlation coefficient (R²).

3.4. Look up Table Inversion Approach

The look up table (LUT) is conceptually the simplest technique to solve the inversion of radiative transfer model. The LUT is built in advance by using the reflectance model in forward mode. The inversion consists in searching the parameter combination in the LUT that minimizes a distance between the simulated and the measured reflectances [42]. To build the LUT, the PROSAIL model was used with the parameter combinations specified in Table 1. To find the solution to the inverse problem, the distance criterion (cost function) was defined as the RMSE between the satellite reflectances and modeled spectra in the LUT.

Several LUTs were tested (Section 4.2) focusing on the influence of the distribution law (Gaussian vs. uniform distributions) as well as the number of simulations, which determines the size of the tables.

3.5. Empirical Transfer Function Approach

Empirical transfer functions relating ground measured CWC and satellite reflectances were established based on a multivariate ordinary least squares (OLS) algorithm [48,49]. The OLS assumes that the dependent variable Y (in situ CWC measurements) is related to the independent variable X_j, j=1, 2..., q, (satellite reflectance values) through the following functional relationship:

$${\text{Y}}_{\text{i}}={\beta }_0+\sum _{\text{j}=1}^{\text{q}}{\beta }_j\cdot {\text{X}}_{\text{ij}}+\epsilon \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}(\text{i}=1,\dots \text{m})$$

(3)

where m is the number of observations and β_j are the parameters of the multiple linear regression to be estimated. The considered multivariate OLS method minimizes the sum of squared deviations of the observed data values away from the fitted response value using an iteratively re-weighted least squares (IRLS) method [49]. IRLS includes a weight factor to adjust the influence of each response value on the model estimates. Samples with weight values lower than 0.7 usually correspond to ESUs located near the field border or having misregistration or experimental errors [49].

IRLS allows testing different spectral bands and selecting the optimal combination based on the following error metrics: (i) minimum weighted RMSE (RW), (ii) minimum cross-validation RMSE (RC) and (iii) minimum number of points with weight lower than 0.7. RC has been calculated from the leave-one-out method by using sequentially a single observation for the validation, and the remaining observations for the training [49]. The RW gives the mean prediction error assumed by the model for all the observations. The RC provides a more reliable evaluation of the model performance since it refers to the prediction error for data in the validation dataset which are not included in the training process to derive the model [50].

To establish optimal empirical transfer functions between ground CWC data and satellite measurements, different band combinations and the widely used normalized difference vegetation index (NDVI) [51] were tested (Section 4.3). NDVI was selected because previous studies have shown its interest for water content retrieval [52–54].

4. Results

4.1. Optimal Modalities of Neural Network Inversion

Table 2. RMSE between networks’ estimates and ground measurements of CWC (kg/m²) for the different modalities considered in the inversion scheme: level of uncertainty between simulated and observed reflectances, architecture of the networks (number of hidden layers and number of neurons per layer) and type of activation functions. The minimum and mean RMSE for three replications are indicated. The optimal modalities are highlighted in bold.

			Min	Mean
Noise level		No noise	0.17	0.19
		2%	0.21	0.23
		4%	0.16	0.19
		6%	0.17	0.19
Architecture	1 hidden layer	2	0.16	0.20
		5	0.16	0.19
		8	0.23	0.25
		11	0.22	0.24
	2 hidden layers	5-2	0.18	0.19
		8-5	0.20	0.21
Activation functions		Tang-Tang	0.16	0.19
		Tang-Lin	0.21	0.23
		Tang-Slin	0.20	0.22

summarizes the RMSE between networks’ estimates and ground measurements for the different modalities of neural networks inversion. Results showed that adding a noise level of 4% to the simulated reflectances slightly improved the performance of neural networks (RMSE of 0.16 kg/m²) although marginal differences are observed as compared to the initial case with non noise being added. Regarding the influence of the architecture, results seems to indicate that a relatively simple network with a single hidden layer of five neurons with 26 synaptic and 6 bias coefficients to be adjusted was performing the best (Table 2). Hyperbolic tangent as activation function in the network’s hidden and output layers provided the minimum RMSE over the ground dataset.

4.2. Optimal Modalities of Look up Table Inversion

Results (

Table 3. RMSE between look up table estimates and ground data for different number of simulations and distributions functions. The optimal modality is highlighted in bold.

		Number of Simulations
		600	3,600	7,200
Distribution functions	Uniform	0.32	0.25	0.25
	Gaussian	0.25	0.22	0.24

) showed that look up tables built with Gaussian distributions focusing on the most frequent cases performed better than those using uniform distributions for which no prior information was exploited. Performances for uniform distributions improve with the number of simulation since it results in a better representation of the space of variables. After 3,600 simulations, no variation in the performances is noticed. Similarly, for Gaussian distribution performances initially improve with the number of simulation but a slight degradation is observed after enough cases (more than 3,600) are considered in the training process. Note however that the differences are not significant and may be partially dependent on the data being used for the validation.

4.3. Optimal Modalities of Empirical Approach

Based on the error metrics (Section 3.5), the best band combination for establishing an empirical transfer function between CHRIS data and ground measurements corresponds to the one using simultaneously the 4 bands selected in Section 3.2 according to the sensibility analysis and data noise. That is, bands centered at 551.9, 641.8, 869.1 and 988.3 nm in green, red and NIR spectral domains respectively. Selection results from a compromise between error values (RC = 0.18 kg/m², RW = 0.15 kg/m²) and number of points (only 3) with associated weight lower than 0.7 (Figure 5a). High correlation (R² = 0.89) between empirical transfer function estimates and in situ data is found (Figure 5b).

4.4. Comparison of Empirical and Physical Approaches

This section is dedicated to the evaluation of retrieval performances and comparison of physical and the empirical approaches. For the validation of the inversion results, the CWC estimates resulting from inverting the PROSAIL canopy reflectance model with the optimal neural network (i.e., considering a moderate uncertainty of 4% in the input reflectances and a simple network of one hidden layer with five neurons and hyperbolic tangent as activation functions) were compared with the available ground measurements for different crop types (Figure 6a). For comparison purposes, the CWC estimates from the optimal look up table approach (i.e., 3,600 entries generated with Gaussian distributions of PROSAIL model parameters) were also validated (Figure 6b). Results (Figure 6,

Table 4. Performances of empirical transfer function (ETF), neural networks (NNT) and look up table inversion approaches as evaluated over the ground measurements. In addition, performances of NNT and LUT at map level as compared with the ETF map. The statistics are: correlation coefficient (R²), root mean square error (RMSE), bias (B), standard deviation (S) and computation time (Time).

	Ground Measurements				ETF Map
	R2	RMSE	Bias	S	R2	RMSE	Bias	S
ETF	0.89	0.16	0.02	0.15	−	−	−	−
NNT	0.82	0.16	−0.04	0.15	0.92	0.14	−0.05	0.13
LUT	0.64	0.22	−0.06	0.21	0.89	0.16	−0.01	0.15

) show that the neural network estimates are in better agreement with ground measurements than the look up table estimates: improvement in the overall performances (0.16 kg/m² as compared to 0.22 kg/m² in terms of RMSE) due to the better precision (higher R²: 0.82 vs. 0.64, and lower S: 0.15 vs. 0.21 kg/m² values) and accuracy (lower bias: −0.04 vs. −0.06 kg/m²) mostly observed for garlic, sunflower and vineyard crops. The performances of neural networks in terms of RMSE (0.16 kg/m²) are similar as those of the empirical transfer function approach although the latter shows a slightly better correlation with ground measurements (R²: 0.89 vs. 0.82) (cf. Figures 5b and 6a, Table 4).

Note that differences between satellite estimates and ground measurements are partially associated to the limitations of remote sensing inversion techniques but also to the uncertainties of in situ measurements. Further, the validation is limited to the reduced number of samples with in situ data (38 ESUs). To complement the direct comparison with ground measurements at ESU level, different approaches for CWC retrieval were compared at map level over the entire CHRIS/PROBA image (Figure 7). CWC map derived from the neural networks approach shows maximum values of 1.5 kg/m² over potatoes and pea crops and zero values over bare soil as expected (Figure 7). White areas correspond to pixels where CHRIS image shows negative reflectance values due to atmospheric correction errors. Comparison of neural networks and empirical transfer function maps shows an overall good agreement with differences in the range of −0.5 to 0.5 kg/m² (not shown for brevity). Statistics of the comparison (Table 4) indicate high correlation (R² = 0.92) and an overall agreement of 0.14 kg/m² in terms of RMSE with practically no bias (B = −0.05 kg/m²) and most of RMSE contribution due to the random differences (S = 0.13 kg/m²). The comparison of look up table maps with empirical transfer functions shows a slight degradation in the performance as compared to neural networks approach (R²: 0.89 vs. 0.92; RMSE: 0.16 vs. 0.14 kg/m²; Table 4).

5. Discussion

Results presented above seem to demonstrate the capacity of CHRIS/PROBA satellite sensor for canopy water content (CWC) monitoring. This is a significant achievement since up to know the potential of using CHRIS/PROBA for CWC retrieval was not studied. To the best of our knowledge, the very few studies in the literature using CHRIS/PROBA had focused on leaf level water content (C_w) but not at canopy level for the estimation of CWC, defined as C_w × LAI product. CHRIS with a maximum wavelength at 1005 nm may only exploit water absorption features centered around 970 nm which is a weaker absorption as compared to main absorption bands in the middle infrared part of the spectrum (at 1,200, 1,450 and 1,950 nm, Figure 3a) where maximum sensitivity to leaf water content (C_w) exists although out of CHRIS spectral range. Spectral sensitivity analysis based on model simulations (Figure 3) and the uncertainty assessment of CHRIS/PROBA bands (Figure 4) suggested that, in absence of middle infrared information, inclusion of bands in green, red and near infrared spectral regions may be optimal for CWC retrieval. Multivariate analysis based on an iteratively re-weighted least square algorithm confirmed that the optimal empirical transfer function relating ground measurements and satellite images makes use of four spectral bands centered at 551.9 nm, 641.8 nm, 869.1 nm and 988.3 nm (Figure 5a). The 4-band relationship outperforms other relationships based on 3 or 2 band combinations including the widely used NDVI (Figure 5a). The band centered at 988.3 nm captures liquid water absorption and enhances C_w sensibility (Figure 3a) but cannot be used alone to retrieve CWC because two other leaf parameters (internal structure (N) and dry matter (C_m)) [15] and canopy level parameters (leaf inclination (ALA) and LAI (Figure 3b)) also influence the reflectance in this band. The three selected bands centered at 551.9 nm, 641.8 nm and 869.1 nm allow optimal estimation of LAI (Figure 3b) and may compensate from leaf- and canopy-level effects. However, band selection may be dependent on the quality of the data being used and the specific conditions of our study. Further confrontation with other existing methods for band selection would be required for extending the applicability of the proposed approach to other datasets and study areas.

Note that the viewing and illumination geometry was fixed in the inversion and the influence of angular effects on CWC retrieval was not here investigated. Further studies should explore the multi-angular capabilities CHRIS/PROBA to improve structure features characterization. Forthcoming studies should also focus on the analysis of soil moisture effects on CWC retrieval based on available soil water content measurements [55] and model simulations.

Our results confirm several previous studies showing the interest of PROSAIL for CWC retrieval over crops and natural vegetation [7,8,15,18,21]. PROSAIL constitutes a good compromise between the realism, complexity of the model, accuracy and computation time requirements [56]. However, this study was carried out over a range of crops for which the turbid medium canopy radiative transfer hypothesis of PROSAIL model was approximately met. It is expected that the performance of PROSAIL may degrade for canopy architectures showing regular or clumped distributions. More realistic model representations may be necessary in such cases.

Results indicate that inversion performances improve through the use of a priori information on model parameters. Gaussian distributions that mimic the actual distribution of the radiative transfer model input variables performed better than uniform distributions for which no prior information was exploited. This validates a posteriori the choice implicitly made in Richter et al. [42] for look up tables and agrees with findings of Verger et al. [28] for neural network inversion methods. Other principles proposed by Verger et al. [28] regarding the architecture and the training process of neural networks were also verified for CWC retrieval.

Comparison of the two considered radiative transfer inversion approaches with ground CWC measurements and the up-scaled map resulting from the application of the empirical transfer function indicates that neural networks outperforms look up table (Table 4). Neural network approach provides similar performances as empirical transfer functions in terms of RMSE as evaluated through ground measurements (RMSE around 0.16 kg/m²) although the empirical approach provides an slightly improvement in terms of correlation (R² of 0.82 for neural network estimates as compared to 0.89 for the empirical approach) at expenses of being more site specific and ground based calibration demanding. The performance of the neural network estimates in terms of the coefficient of correlation with ground CWC data (R² = 0.82) obtained in this study is similar to results in the literature using other approaches such as the best derivative at 950.5 nm (R² = 0.80) reported by Clevers et al. [8]. Note, however, that the obtained results are very particular to the conditions of the study and an extensive validation and comparison with other existing approaches over a wider range of vegetation conditions should be conducted.

6. Conclusions

This study focused on the estimation of canopy water content from CHRIS/PROBA reflectances. Empirical and physical approaches, respectively based on up-scaling ground measurements and the inversion of the PROSAIL radiative transfer model, were assessed. Emphasis was put on the spectral band selection and the optimal modalities of inversion. Performances were evaluated based on the comparison with ground measurements acquired over Barrax cropland area during SEN3EXP experiment. Results from our analysis led to the following conclusions:

CHRIS/PROBA satellite sensor allows CWC monitoring with an overall performance of around 0.16 kg/m² in terms of RMSE. This is a very noticeable result since the main water absorption bands (>1,200 nm) are out of the limited spectral range of CHRIS/PROBA satellite sensor with maximum wavelength at 1005 nm. The proposed approach makes use of four spectral bands for estimating the two components (C_w and LAI) of CWC: one band centered around 970 nm water absorption feature which is sensible to C_w; and green, red and near infrared bands for LAI retrieval and minimizing leaf- and canopy-level effects.

Principles proposed by Verger et al. [28] for LAI, FAPAR and FCOVER variables were here verified for CWC retrieval. Gaussian distributions that mimic the actual distribution of the radiative transfer model input variables performed better than uniform distributions for which no prior information was exploited. Regarding the network’s architecture, a relatively simple one with a single hidden layer of five tangent sigmoid transfer functions outperformed more complex networks. Inclusion of moderate (4%) additive Gaussian noise in the training process increased the adaptive capacity of the system to input and model uncertainty and resulted in an improvement of the overall accuracy as evaluated with ground measurements.

Neural networks outperformed look up tables for CWC estimation when evaluated over the available ground measurements (0.16 kg/m² as compared to 0.22 kg/m² in terms of RMSE) due to the better precision (higher R²: 0.82 vs. 0.64, and lower S: 0.15 vs. 0.21 kg/m² values) and accuracy (lower bias: −0.04 vs. −0.06 kg/m²). In fact, radiative transfer-neural network inversion approach showed similar performances as empirical up-scaling transfer functions (RMSE around 0.16 kg/m²), the latter being more site specific and calibration demanding. Comparison at map level over the entire CHRIS/PROBA image showed an overall good agreement (0.14 kg/m² in terms of RMSE with practically no bias (B = −0.05 kg/m²) and correlation of R² = 0.92) between neural network estimates and up-scaling empirical approach. This suggests the interest of the proposed physical approach that could potentially be used for operational monitoring of CWC. The forthcoming Sentinel-2 sensors that provide both a frequent revisit capacity and a decametric spatial resolution with a wider range of spectral sensitivity to the water content as compared to CHRIS/PROBA sensor open new perspectives for CWC monitoring.