A Multi-Year Study on Rice Morphological Parameter Estimation with X-Band Polsar Data

Rice fields have been monitored with spaceborne Synthetic Aperture Radar (SAR) systems for decades. SAR is an essential source of data and allows for the estimation of plant properties such as canopy height, leaf area index, phenological phase, and yield. However, the information on detailed plant morphology in meter-scale resolution is necessary for the development of better management practices. This letter presents the results of the procedure that estimates the stalk height, leaf length and leaf width of rice fields from a copolar X-band TerraSAR-X time series data based on a priori phenological phase. The methodology includes a computationally efficient stochastic inversion algorithm of a metamodel that mimics a radiative transfer theory-driven electromagnetic scattering (EM) model. The EM model and its metamodel are employed to simulate the backscattering intensities from flooded rice fields based on their simplified physical structures. The results of the inversion procedure are found to be accurate for cultivation seasons from 2013 to 2015 with root mean square errors less than 13.5 cm for stalk height, 7 cm for leaf length, and 4 mm for leaf width parameters. The results of this research provided new perspectives on the use of EM models and computationally efficient metamodels for agriculture management practices.


Introduction
Rice is the main source of food and income for several highly populated countries.The increasing population and limited arable lands bring out the need for higher yields that depends on the development of better management practices.The traditional method monitoring by visual inspection is not possible for kilometer-square areas.For such large scales, remote sensing based methods are good alternatives.Among different data sources, information provided by Synthetic Aperture Radar (SAR) is advantageous with its sensitivity to geometric and dielectric properties of the objects and its availability in all light and weather conditions.Therefore, SAR is a valuable tool for monitoring the rice fields [1].
Understanding the phenological evolution of rice fields is important to develop effective management strategies.In the literature, several SAR data based algorithms have been developed to investigate the phenological evolution of rice fields including their canopy height, growth stage, leaf area index and yields.In SAR based rice monitoring, one way of monitoring the phenological cycle is the use of different SAR data analysis techniques including temporal trend analysis [2,3], interferometric analysis [4,5] and polarimetric interferometry analysis [6,7].The other approach employs EM models to simulate the polarimetric parameters from plant properties [5,8,9].This letter presents the inversion results of the metamodel-driven EM model for the estimations of stalk height, leaf length, and leaf width parameters from copolar X-band TerraSAR-X time series data, by following the methodology given in Figure 1, detailed and presented in [10] for a single year and growth stage parameter, specifically BBCH.Unlike [10], this study focuses on the feasibility of the proposed approach for morphology parameter estimation under different conditions.The chosen EM model [11] uses the biophysical properties of rice plants to simulate the backscattering intensities ( σo ) with Monte-Carlo simulations.The computation costs related to the multi-dimensional algorithm and the Monte-Carlo simulations were reduced by introducing metamodels, which mimics the EM model after trained for once.The stability of the inversion was tested over a dataset of three different cultivation periods (2013)(2014)(2015), which includes data from rice fields under different agricultural and environmental conditions.This paper has four sections.It starts with the proposed methodology in Section 2 by providing the EM model, its metamodel, and the inversion procedure.Sections 3 and 4 present the SAR and ground data followed by the inversion analysis results.The overall summary is provided in Section 5.

Theoretical EM Model and Its Metamodel
In this study, the EM model [11], M(ξ), is employed to simulate the backscattering intensities, σo , from a set of rice plant morphology parameters, ξ, through Monte Carlo simulations for varying scatterer locations.The ξ set includes stalk dimensions (height and diameter), leaf dimensions (length and width), panicle dimensions (length and width), and their structural densities.The simulations are done for a unit area A, having randomly placed non-overlapping cylindrical stalks with a specific height and diameter.Each stalk is modeled to have elliptical leaves with a fixed length and width.In this study, we assumed flooded ground and plant components with fixed complex dielectric constants for the complete growth cycle by relying on the sensitivity analysis of the EM model [12].The M(ξ), given in (1), provides the relation between an incident Ēi and a scattered wave Ēs through the coherent sum of four different scattering mechanisms (S n ), as shown in Figure 2.
Figure 2. The four different scattering mechanisms considered within the EM model [11].S 1 : Direct scattering from the scatterers, S 2 : Scattering from the canopy followed by reflection from the ground, S 3 : Reflection from the ground followed by scattering from the canopy, S 4 : Reflection from the ground followed by scattering from the canopy and followed by reflection from the ground.
In the EM model given in (1), q and p subscripts correspond to transmitted and received horizontal (H) and vertical (V) linear polarization channels.The parameters k and r represent the free-space wavenumber and the distance between the sensor and the target, respectively.The σ0 qq for the different polarimetric channels, qq, are approximated from the ratio between E s q and E i q .The M(ξ) has high computation cost due to its multi-dimensional algorithm and Monte Carlo simulations.The computation costs of the algorithm were reduced using sparse Polynomial Chaos Expansion (PCE) metamodels.The PCE metamodels are spectral variance decompositions of the original model with low training cost and wide coverage in the parameter domain [13].For the chosen M(ξ), its PCE metamodel, PCE EM (ξ), is developed from the (2), given below.
In ( 2), a j ∈ R is a set of scalar coefficients and the Ψ j (ξ) ∈ R form a polynomial orthonormal basis [14].For practical reasons and to avoid over fitting conditions, the metamodels were limited to N (=20) expansions.In this study, the PCE EM metamodel, Y(ξ), was implemented the UQLab toolbox [15] with Legendre polynomial family with the uniform [−1,+1] input distributions.Details of the metamodel implementation for the EM model can be found in [10].

Probabilistic Particle Swarm Optimization
The inversions of multi-dimensional EM models are ill-posed problems with the higher number of unknowns compared to the number of equations.For such problems, optimization algorithms are used to reach the optimum solution in the parameter space using different constraints.However, the optimization of multi-dimensional EM model inversions may result in multiple optimum solutions since different inputs can lead to similar outputs.The presence of multiple optimum solutions prevents the use of deterministic optimization algorithms that focus on a single solution.For the existence of multiple solutions, stochastic optimization algorithms can be considered as an alternative.In stochastic optimization, the procedure is initiated several times to obtain all local solutions in a given parameter space based on the defined set of rules.
In this study, the Particle Swarm Optimization (PSO) algorithm [16] is utilized, which is based on updating the position of the particles, i.e., possible solutions, until they converge to an optimum solution within a parameter space.In each iteration, the locations of the particles are updated according to the position of the particle with the best position.The iterations continue until the particles converge to a solution that agrees with the defined constraints.For the estimation of rice morphology parameters from copolar X-band SAR backscattering intensities, the PCE EM metamodel is inverted with the stochastic PSO algorithm.
The fitness function for the PSO is given in (3) for HH and VV polarimetric channels and an arbitrary i th iteration.The fitness function is defined to minimize the difference between measured σ o HH,VV and estimated σo HH,VV values.The consistency of the solutions for different polarimetric channels is provided by considering the same input vector for both HH and VV polarimetric channels.
Optimization problems need constraints to simplify the problem by reducing the dimensional complexity.In PCE EM metamodel inversion, three plant morphology dependent constraints were established that are based on the natural limitations excerpted from the available rice morphology data.
Positivity constraint ensures positive and real morphological estimations for all iterations.Min-Max constraint limits the morphological estimations based on the phenological phase boundaries.Phenological phases are defined by the International Rice Research Institute (IRRI) [17].
The scale divides the growth cycle to five major phases.Figure 3 presents the boundaries obtained from ground data, for the chosen morphological parameters.
Natural limitations provide non-linear relationships among rice morphology parameters during their development.These relations eliminates the solutions which might be impossible for a healthy plant, such as a 10-cm-tall stalk having two 50-cm-length leaves.This condition restricts the parameter space with a convex hull.Convex hull defines non-linear boundaries in a parameter space according to the experimental ground measurements that involve the agronomically possible solutions.The stochastic PSO optimization provides distributions for each morphology parameter.For the accuracy analysis, the mean value of the resulting stochastic distributions is assigned as the estimated dimensions of the rice morphological parameters.For a single stochastic PSO optimization, a total number of 200 iterations was found to be sufficient for the optimization convergence which changes less than 0.1% regarding the mean of the estimated values.

The Ipsala Test Site and Ground Data
The selected test area, Ipsala, is located in the North-West part of Turkey with its center at 37°7 ′ 53 ′′ N and 6°19 ′ 32 ′′ N coordinates.Ipsala is one of the biggest rice cultivation sites in Turkey with approximate acreage of 190 square kilometers.Based on the knowledge gathered from the Trakya Agriculture Research Institute (TARI), rice cultivation is done in the area by local farmers between May and September.
As shown in Figure 4, field campaigns were conducted ±5 days of SAR acquisitions to have representative rice morphology parameters.In each year, the test fields were selected by the expertise of the TARI researchers.To monitor the evolution of the morphology parameters, at each field following parameters were measured: above water stalk height, stalk diameter, leaf length, leaf width, the number of plants per m 2 , the number of tillers per plant, and the number of leaves per tiller.In Figure 3 the evolution of some morphology parameters is depicted for each IRRI phase using a box-whisker plot.For each parameter, a quasi-linear increase is observed until IRRI-3.From IRRI-3 on, a decrease is noticed in the stalk diameter due to reduced water content.

SAR Data
During the study, the Ipsala test site is monitored with the data acquired from the TerraSAR-X satellite with an average incidence angle of 31°.The TerraSAR-X satellite has the central frequency of 9.65 GHz and temporal resolution of 11 days.The acquisition dates are given in Figure 4.The data were delivered in single look complex format and were spatially and temporally co-registered by bilinear interpolation.

Results and Discussion
In this section, we present the stochastic inversion results of the same PCE EM metamodel for the estimation of stalk height, leaf length and leaf width parameters from X-band co-polar SAR data on a dataset spanned over three cultivation periods.For the analysis, the noise in the TerraSAR-X data was reduced using 13 × 13 boxcar smoothing windows, resulting in 33 m × 25 m spatial resolution.The estimation accuracies of the chosen parameters are evaluated based on their correlation against the ground measured values.The results are reported with their coefficient of determination (R 2 ) and the root-mean-square error (RMSE) values.
Table 1 provides the input parameters for the PCE EM metamodel evolutions that were assumed constant.The details of the EM model, its metamodel, their simulation accuracies and growth phase based global sensitivity analysis results of the PCE EM metamodel are provided in [12].The chosen EM model and the proposed stochastic inversion approach consider the simplified plant morphology.As observed in Figure 3, the evolution of chosen morphological parameters are related to each other.In the proposed method these morphology based relations provide the base of the optimization constraints.As the stochastic inversion algorithm provides multi-dimensional parameter distributions, the mean values of the estimations are used for the calculation of R 2 and RMSE values.
The stochastic inversion results for stalk height, leaf length, and leaf width parameters are presented in Figure 5 and Table 2.The R 2 and RMSE were calculated against the estimated biophysical parameters are given for the cultivation periods from 2013 to 2015.The accuracy analyses were conducted using a total number of 93 different σ o values (32 from 2013, 25 from 2014, and 36 from 2015) measured from different fields having different morphologies, agricultural practices and environmental conditions.The global sensitivity analysis of the PCE EM metamodel was previously discussed in [12].The analysis results emphasize the importance of stalk height, leaf length and structural density parameters (number of plant components in a unit area) throughout the growth cycle.The following deductions can be made from the stochastic inversion results.Stalk Height estimations had the highest accuracy (R 2 ≥ 0.86) considering the entire dataset of three cultivation periods.As presented in Figure 5, the stalk height estimation results are slightly over-estimated for rice canopies shorter than 0.6 m and under-estimated for the taller canopies.The performance of the algorithm for the 2015 dataset was calculated to be lower compared to the other years with an RMSE value of 13.5 cm.This situation can be related to the variance of the σ o values and the PCE EM metamodel simulations for the corresponding morphology parameter ranges in the parameter space.On the other hand, the estimation bias between the measured and estimated values shows an increasing spread with increasing stalk height.The variation in the bias can be interpreted by the presence of plants with varying physical structures and similar backscattering behaviors at later phases of the growth cycle.The RMSE values were calculated to be less than 13.5 cm, which are acceptable for the stalk height estimations calculated from copolar SAR data.
Leaf Length estimation accuracy is calculated to be lower than the stalk height estimation accuracy for the dataset of three cultivation periods.As shown in Figure 5, the leaf lengths are mainly over-estimated when they are shorter than 50 cm and under-estimated when they become longer.The stochastic inversion results of the PCE EM metamodel shows highly accurate results (R 2 ≥ 0.78 and RMSE ≤ 7 cm) for the evaluations from the complete dataset.From three different years, the lowest accuracy was obtained from the data acquired during 2015, as it was observed for the stalk height.Regarding the estimation bias, it is noted that the errors tend to increase with increasing leaf length.Similar to the interpretation provided for stalk height estimations, the increasing error with increasing leaf length can be explained by the presence of higher variance in the parameter space at later growth phases.
Leaf Width is morphologically related to the leaf length due to natural growth limitations.The accuracy analysis on leaf width estimations presented acceptable values with R 2 values higher than 0.61 and RMSE values lower than 3 mm.Similar to stalk height and leaf length parameters, the lowest accuracy was again obtained from the copolar SAR measurements of 2015.The estimation bias of the analysis is observed to vary between ±5 mm for the complete growth cycles of three years.
The stochastic inversion of the PCE EM metamodel provided successful estimations for stalk height and leaf length parameters for the cultivation periods of three years.However, from the chosen rice morphology parameters, leaf width estimations had lower accuracies with higher relative RMSE values.This situation is supported by the global sensitivity analysis results, which states the importance of stalk height and leaf length parameters on the PCE EM metamodel simulations, while mentions the lower effect of the leaf width [12].
The accuracy analysis exhibited in this study combines 93 different rice plant morphology and SAR measurements collected from the cultivation periods of 2013 to 2015.Concerning the presence of various environmental factors and agricultural practices, the results are considered to be representative for the wet-cultivated and broadcast seeded rice fields that are located in Ipsala (Turkey).Therefore, the results are encouraging for the development of new management practices, which can use the estimated morphology parameters to interpret the yield and the detailed growth stage.

Overview and Summary
In this research, a stochastic inversion method has been presented to invert a multi-dimensional PCE EM metamodel, trained from an EM model [11].Apart from the previous studies that focuses on the estimation of the stalk height, this study extends the estimations to the stalk height, leaf length and leaf width parameters of rice plants by considering the natural growth limitations.This study was conducted over three cultivation cycles for validation purposes.
We have tested the proposed inversion algorithm on a three-year cultivation period of ground, and copolar high spatial resolution and high frequency SAR datasets.The data were acquired by TerraSAR-X over broadcast seeded and flooded rice fields.For method validation we trained the PCE EM metamodel by employing the EM model [11], ground measurements and field average SAR backscattering intensities.
As a result, we obtained significant correlations between the estimated and measured values of stalk height, leaf length and leaf width parameters using the proposed PCE EM metamodel inversion scheme.From the analysis, we calculated RMSE values less than 13.5 cm for stalk height, 6.9 cm for leaf length and 34 mm for leaf width parameters.The results pointed out that the use of spaceborne X-band PolSAR data is powerful for the development of new agriculture monitoring practices.We should mention that the overall performance of the proposed approach mainly relies on the accuracy of the EM model and the PCE EM metamodel, which shows the importance of EM model selection.
The presented algorithm has presented several advantages such as coupling the agronomical growth rules with the EM models for efficiency, inversion with stochastic optimization to handle environmental variability and most importantly being computationally efficient with the use of metamodels that are trained for once.In the presented study, the requirement of a priori knowledge of the growth phase and plant morphology information can be seen as a disadvantage.However, considering the importance of rice as the main source or food and income for several highly populated countries, related information can be found in the crop databases.Besides, the literature currently covers several different methodologies to determine the growth phase of rice plants [3,18,19].
In the future, our studies are going forward to improve the approach by substituting the a priori growth phase information with canopy height estimations that can be directly calculated from Pol-InSAR data.In addition, it is planned to extend the applicability of the inversion procedure to monitor other major crops.

Figure 1 .
Figure 1.Block-diagram of the proposed stochastic inversion of metamodel-driven electromagnetic scattering (EM) model.

Figure 3 .
Figure 3. Temporal variations of the rice crop biophysical collected between 2013-2015 from the regularly conducted ground campaigns in Ipsala (Turkey).Due to technical limitations, stalk diameter data was not measured in 2015.The measurements are grouped according to the IRRI growth phases as a Box-and-Whisker plot.Box presents the information for the quartiles while the whiskers present minimum and maximum values for IRRI growth stages IRRI-1 [Early vegetative], IRRI-2 [Late vegetative], IRRI-3 [Early reproductive], IRRI-4 [Late reproductive], IRRI-5 [Maturative].

Figure 4 .
Figure 4.The ground and Synthetic Aperture Radar (SAR) data collection dates for the cultivation period between 2013 and 2015.On the right side, the VV channel PolSAR intensity image is also provided from the 30 June 2013.In the image, rice fields can be detected with their higher backscattering intensities.

Figure 5 .
Figure 5.The measured values of stalk height, leaf length and leaf width parameters versus their estimations given in a correlation scatter plot for years 2013 to 2015 with a reference line.Symbols are colored with respect to the year: 2013 (◆), 2014 (◆) and 2015 (◆).

Table 1 .
Input parameters that are assumed constant during the EM model and PCE EM metamodel simulations.

Table 2 .
Accuracy analysis of the stochastic inversion of PCE EM metamodel given with the calculated R 2 and root-mean-square error (RMSE) values for the years 2013, 2014, 2015 and the complete dataset.