Estimating Plant Traits of Grasslands from UAV-Acquired Hyperspectral Images: A Comparison of Statistical Approaches

Grassland ecosystems cover around 40% of the entire Earth’s surface. Therefore, it is necessary to guarantee good grassland management at field scale in order to improve its conservation and to achieve optimal growth. This study identified the most appropriate statistical strategy, between partial least squares regression (PLSR) and narrow vegetation indices, for estimating the structural and biochemical grassland traits from UAV-acquired hyperspectral images. Moreover, the influence of fertilizers on plant traits for grasslands was analyzed. Hyperspectral data were collected from an experimental field at the farm Haus Riswick, near Kleve in Germany, for two different flight campaigns in May and October. The collected image blocks were geometrically and radiometrically corrected for surface reflectance. Spectral signatures extracted for the plots were adopted to derive grassland traits by computing PLSR and the following narrow vegetation indices: the MERIS Terrestrial Chlorophyll Index (MTCI), the ratio of the Modified Chlorophyll Absorption in Reflectance and Optimized Soil-Adjusted Vegetation Index (MCARI/OSAVI) modified by Wu, the Rededge Chlorophyll Index (CIred-edge), and the Normalized Difference Red Edge (NDRE). PLSR showed promising results for estimating grassland structural traits and gave less satisfying OPEN ACCESS ISPRS Int. J. Geo-Inf. 2015, 4 2793 outcomes for the selected chemical traits (crude ash, crude fiber, crude protein, Na, K, metabolic energy). Established relations are not influenced by the type and the amount of fertilization, while they are affected by the grassland health status. PLSR is found to be the best strategy, among the approaches analyzed in this paper, for exploring structural and biochemical features of grasslands. Using UAV-based hyperspectral sensing allows for the highly detailed assessment of grassland experimental plots.


Introduction
Grassland covers roughly 40% of the total world land area and this extension corresponds to approximately 52.5 million km 2 [1]. Moreover, grassland ecosystems have extremely variable features, since they are strongly influenced by environmental conditions, such as topographic location, and by anthropogenic effects, such as (in)organic fertilizer application [2]. Therefore, efficient management plans of grasslands require insight into the health status and spatial variation in order to improve conservation and to achieve optimal growth. The planning of grassland management is even more important if investigated at the field scale, since the farmers need a guide for identifying the optimal time for fertilizer application and for predicting the ideal time of harvest [3]. Up until now, the schedule of grassland management was mainly organized according to qualitative information deducted from farmers' experiences. The integration of quantitative and spatial data into the already existing management plans could significantly improve the health status of grasslands and the feeding quality of the final harvest.
Currently, traditional field-based methods are usually preferred, even if their application is strongly restricted because they are time-consuming, destructive, and cannot be repeated often enough for high spatial resolution investigations over large areas [3,4]. Over the years, alternative techniques have been developed, and, nowadays, remote sensing is recognized as a suitable technique thanks to its ability to characterize the land surface in a fast and relatively cheap way. For these reasons, it has been widely used for estimating various biophysical and biochemical vegetation variables in grasslands, such as the Leaf Area Index (LAI) or chlorophyll [5][6][7][8][9][10][11]. The subsequent introduction of hyperspectral sensors has allowed researchers to improve the retrieval of grassland traits considerably [3][4][5]. Based on its narrow contiguous wavebands, hyperspectral sensors are more sensitive to vegetation variables since they provide a continuous reflectance spectrum of the vegetation target [4,[12][13][14][15][16]. Nevertheless, the hyperspectral images include hundreds of spectral bands, many of which are strongly correlated. Consequently, it has been necessary to select a subset of data in order to decrease the dimensionality of the dataset and reduce redundant information. Therefore, the selection of the subsets can be carried out by taking into account the sensitivity of the vegetation variables to the spectral bands [17].
Further improvements could be achieved by adapting the spatial resolution to the size of the object investigated. In this way, the image captures the spatial detail needed to resolve important patterns in the field. Given the fine spatial scale of grassland biophysical traits, Unmanned Aerial Vehicles (UAV) can provide the necessary fine spatial scales because they can be flown at a low altitude. In order to analyze biophysical traits of grassland plants, it is necessary to have the smallest possible resolution. Therefore, the application of Unmanned Aerial Vehicles (UAV) could provide a useful support, as they permit us to drastically reduce altitude in order to increase the spatial resolution and determine it based on the object size under examination and the aim of the research [18][19][20]. Furthermore, UAVs can provide improved access to more remote targets and can be deployed flexibly, eliminating problems due to clouds. For these reasons, the utility of hyperspectral cameras mounted on UAVs has been explored in a few papers [21]. Indeed, although the application of hyperspectral imaging and of UAVs is significantly widespread for quantitatively the Earth system, their combination is still rather limited.
Two different statistical approaches are commonly used to analyze the relationship between spectral measurements and vegetation traits: univariate and multivariate regression models. The first approach is mainly based on a regression model between a so-called spectral vegetation index, a limited set of combined spectral bands, and the vegetation trait of interest, while the multivariate approach is based on the application of regression modeling between all observed spectra variables and a vegetation trait [22].
Depending on the spectral resolution of the sensor applied, the calculation of vegetation indices can be based on both narrow or broad spectral bands, obtaining different values for the same index. Nevertheless, when comparing the results of the two different approaches for computing the same index, some studies have shown that the use of narrow-band indices results in improved statistical relationships between biochemical and biophysical traits and vegetation indices [15,23]. Moreover, earlier studies have shown that narrow-band information of hyperspectral images is important for providing information for the estimation of biochemical and biophysical vegetation variables by applying the vegetation indices [14,15,[23][24][25][26][27][28]. However, in this way, it is not possible to fully exploit the hyperspectral data potential of the large number of spectral bands available. For this reason, numerous researchers have focused their attention on Stepwise Multiple Linear Regression (SMLR), which takes advantage of almost all of the hyperspectral information available [7,[29][30][31][32][33][34]. However, SMLR suffers from multicollinearity problems and the extensive spectral overlap of individual biochemical properties [35]. Partial least squares regression (PLSR) has been recognized as an alternative technique to SMLR, and it has been used earlier in the spectral data analysis for wheat [8] and heterogeneous grasslands [36,37].
Several studies have dealt with validity and accuracy in the analysis of grasslands using the two statistical methods described. However, no studies have been published on the evaluation of these techniques for estimating grassland variables relevant for forage quality assessment, such as crude ash or metabolic energy. In addition, grassland traits such as sodium are also relevant to assess the taste for feeding cattle. Thus, the utility of the two statistical methods for estimating biochemical traits in grasslands from spectral datasets and the potential influence of the type and amount of fertilization on these grasslands needs investigation.
The overall objective of the present study is to investigate the utility of hyperspectral images acquired using an unmanned airborne platform for predicting vegetation traits in grasslands. Specific sub-objectives can be summarized as follows: 1. Compare two regression methods based on vegetation indices and the PLSR approach in order to choose the best strategy for predicting bio-physical and bio-chemical plant traits of grasslands; 2. Investigate the influence of the amount and the type of fertilization on grassland traits and the resulting spectral curve; 3. Evaluate the influence of the phenology of grasslands and the timing of spectral data collection on the established regression relations.

State of the Art of UAVs Applications
In the last few years, UAV applications in the civil field have become increasingly popular. Their use is becoming widespread in many scientific disciplines, such as earth observation and precision agriculture [38][39], in which high temporal frequency and detailed spatial resolution are necessary in order to improve the health status of soil and crops [40][41][42].
This has resulted in the development of suitable instrumentation for such purposes. With the recent technical electronic and optical improvements for aerial platforms and the devices mounted on them, the field has developed fast. Indeed, thanks to the evolution of the integrated circuits and radio-controlled systems, it is possible to adopt the use of UAVs [43]. They can fly at lower flight altitudes than traditional airborne platforms, allowing an increase in spatial resolution, reaching difficult-to-access areas and acquiring images in flexible data, eliminating cloud cover problems. Moreover, their large market penetration and continuous development have led to a drastic reduction in their cost [19].
Furthermore, the size and weight of electronic devices, used for capturing images, have been modified, and strongly decreased overall, so that they can be mounted on the UAVs [44]. Also, the introduction of small digital multispectral sensors has further allowed us to expand the areas investigated and to improve the results obtained [45], as shown in several previous studies for grassland monitoring.
After comparing and assessing the potential of several multispectral sensors and small RGB digital cameras, the authors of [46] has shown the limitations of RGB cameras and the superiority of multispectral sensors in grassland traits estimation. In the same way, the authors of [47] has shown that there is a strong relationship between vegetation indices, calculated from multispectral images, and grassland bio-physical parameters, useful for temporal changes and, consequently, future crop growth monitoring. Subsequently, the introduction and the miniaturization of hyperspectral sensors has allowed us to obtain further improvements in grassland analysis [21].

Study Area and Field Data Collection
The study area consisted of a grassland fertilization experiment located on the experimental farm Haus Riswick near Kleve in Germany (51°47′12.5″N, 6°10′08.7″E). The experiment consists of 60 grassland plots with different levels of organic and inorganic fertilizer ( Figure 1). Each plot had an area of 12 m 2 with a length of 8 m and a width of 1.5 m. The plots were split in four groups and for each of them different types and levels of fertilization were applied ( Figure 2). In particular, a comparison was made between an inorganic fertilization (calcium ammonium nitrate) with six different levels (0, 85, 115, 170, 230, and 340 kgN/ha) and an organic fertilization with three different levels (170, 230, and 340 kgN/ha). The choice of the fertilizer levels was planned based on the customs of German farmers.
In Germany, farmers normally are allowed to apply 170 kg· N· ha −1 in the form of organic manure and only the maximum of 230 kg· N· ha −1 under special conditions (no grazing). Therefore, these two levels were tested in the grassland experiment, with an additional inflated level of 340 kg· N· ha −1 . This high level of 340 kg· N· ha −1 is the maximal amount of inorganic nitrogen that can normally be taken up by grassland during a year in this countryside. To control the efficiency of the organic fertilizer, nitrogen-equivalent levels of inorganic nitrogen increasing from 0 up to 340 kg· N· ha −1 (by steps: 0, 85, 170, 230, 340 kg· N· ha −1 ) were tested. Moreover, in order to study the effect of slurry on a grassland experiment, different amounts of slurry were applied for one (2012)  After the flights, the grassland height (H) was determined in the field using a grass height meter (Eijkelkamp, The Netherlands). Every plot was completely harvested using the Haldrup C-65 plot combine (Haldrup, Denmark) and vegetation traits were analyzed based on the harvested grassland material. The selection of the grassland traits was based on their relevance for the analysis of grassland forage quality (including nutrition and taste) and to which extent this can already be determined in situ in the field. First of all, the vegetation fresh biomass (FB) was determined, after which the samples were dried for 24 h at 105 °C and the dry weight and dry matter (DM) were determined. Subsequently, the bio-chemical composition of the grassland was determined. Analysis of the collected material was carried out by Landwirtschaftliche Untersuchungs-und Forschungsanstalt Nordrhein-Westfalen (LUFA NRW) using the following methods of the Verband der Landwirtschaftlichen Untersuchungs-und Forschungsanstalten (VDLUFA) [49]. The following grassland traits were determined: (2) Figure 2. Map of different types and levels of fertilization applied on the experimental grassland plots. The 60 plots were split into four groups of 15 plots and for every group the treatments were applied as indicated in the table. The black line is a measurement rail which is also shown in Figure 1.
In order to take into account the influence of fertilizer type and amount, the datasets acquired during the two field work campaigns were integrated and three datasets were prepared: (1) all plots prepared with different levels of inorganic fertilizer for both May and October campaigns; (2) all the data of all plots arranged with both inorganic and organic fertilizer for both May and October campaigns; (3) this dataset takes into account all the May campaign data related to all plots arranged with different levels of inorganic and organic fertilizer.
The subsequent regression analysis was performed on all three datasets and then compared. Comparison of the datasets allows us to investigate the influence of fertilizer type (comparison of 1 and 2) and differences due to the growing season (comparison 2 and 3).

UAV and Its Equipment
The study area was surveyed using an octocopter UAV (Aerialtronics Altura AT8 v1A) equipped with a ground station for planning and control of the flight paths and all necessary hardware and software tools for its flight programming ( Figure 3). In order to acquire hyperspectral images, the Wageningen UR Hyperspectral Mapping System (HYMSY) was mounted under the octocopter [51]. HYMSY is a lightweight hyperspectral pushbroom system (2.0 kg) specifically designed for small UAVs and developed at Wageningen University and Research Centre. The sensor consists of a photogrammetric camera (Panasonic GX1), a global position system, and a hyperspectral camera, combining the PhotonFocus SM2-D1312 industrial camera and the Specim ImSpector V10 2/3 spectrograph. HYMSY operates over the wavelength range of 400-950 nm with a spectral resolution of 9 nm. After processing of the imagery, HYMSY can deliver a RGB ortho-mosaic, Digital Surface Model (DSM), and a hyperspectral dataset.

Hyperspectral Data Collection and Pre-Processing of Spectra
Two flights with the Unmanned Aerial Vehicle (UAV) were made over the experimental grassland plots on 15 May 2014 and 14 October 2014. For both days, flights were carried out under clear sky conditions and with the same flight pattern, altitude, and speed. The flights were prepared by defining a set of way-points which covered the grassland experiment with three flight lines over the three rows of experimental plots ( Figure 1).
For the mapping flights, the UAV altitude and speed were programmed to obtain a Ground Sampling Distance (GSD) of 200 mm for the hyperspectral data (cross-track) and 20 mm for the aerial photos. Flights were made at a height of 70 m with a flight speed of 5 m/s. The acquired raw dataset was subsequently radiometrically calibrated, converted to surface reflectance, and geometrically corrected. Hyperspectral processing was divided into the following steps: conversion of the digital numbers into radiance spectra using laboratory calibration; conversion of the radiances into reflectance factor spectra using Spectralon panel data; and georeferencing the scanlines using ReSe PARGE software, produced by ReSe Applications Schlä pfer and Remote Sensing Laboratories (RSL) of the University of Zurich [52]. Radiance was converted to reflectance in several steps. First, the scan lines were calibrated to the illumination conditions by taking ground-based HYMSY measurements of a 25% Spectralon reference panel acquired before and after the flight. These reference panel radiance spectra constituted two different datasets, which are separately analyzed. Both reference sets were used to transform all radiance values of the scan lines to reflectance values using a linear interpolation. This procedure allows us to linearly interpolate the radiance values through time in order to correct the irradiance changes, depending on the atmospheric conditions, during the UAV image acquisition. At this point, the two datasets of reflectance values are merged into one file in 16-bit ENVI BSQ format. The georectification and ortho-mosaicing of the aerial photos was carried out using PhotoScan Pro software (v1.0.0, Agisoft, St. Petersburg, Russia,) which involves the image orientation by extracting the tie point followed by the alignment of photos using block bundle adjustment [53][54][55]. In order to improve the georeferencing and tie point extraction, during the alignment step, it was essential to insert the coordinates of the ground control points (GCPs), measured with an RTK-GPS. More details on the calibration and processing procedures can be found in [52].
From the processed hyperspectral images of the two datasets, complete reflectance factor spectra were extracted for every plot of the grassland experiment by selecting the pixels within a region of interest (ROI). For every plot the ROI was defined according to the plot dimensions but taking into account a buffer of two pixels from the edge in order to exclude boundary effects. The selected spectra for each grassland plot at the two measurement dates were averaged, integrated, and analyzed in the following steps when computing the vegetation indices and the partial least squares regression technique. The average operation allowed us to minimize the noise in the reflectance spectra.

The Narrow-Band Vegetation Indices
The use of vegetation indices is the most common method to predict vegetation variables, drastically reducing the dimensionality of the hyperspectral data set. For this study, it was preferred to compute narrow-band vegetation indices and not broad vegetation indices, in order to take into account the sensitivity of vegetation indices to the narrow spectral bands [15,23].
Narrow-band vegetation indices were calculated from the average reflectance spectra of each plot. Although additional vegetation indices present in previous studies [5,[12][13][14][15][24][25][26][27][28]37] were analyzed, in this paper, only the indices that gave the best results are shown. In particular, the MERIS Terrestrial Chlorophyll Index (MTCI) [56], the ratio of the Modified Chlorophyll Absorption in Reflectance and Optimized Soil-Adjusted Vegetation Index (MCARI/OSAVI) modified by Wu [57], the Red-edge Chlorophyll Index (CIred-edge) [58,59] and the Normalized Difference Red Edge (NDRE) [60] were selected for the estimation of structural and biochemical traits. The equations and wavelengths of the selected four vegetation indices are shown in Table 1.
Linear regression was performed in order to analyze the relationship between grassland traits and the vegetation indices. The choice of the pattern is strongly influenced by the distribution of data and, consequently, by the type of relationship between vegetation index and trait. Consequently, the assessment of the quality of the relationship was achieved by calculating the coefficient of determination (R 2 ) and the root mean square error (RMSE). According to classification suggested by [27], we distinguish three classes: (1) strong correlation: R 2 > 0.7; (2) moderate correlation: 0.5 < R 2 < 0.7; and (3) weak correlation: R 2 < 0.5

Partial Least Squares Regression (PLSR)
Partial least squares regression (PLSR) is a useful statistical multivariate regression technique for explaining the relationship between hyperspectral data variables, the independent variables (X), and grassland traits, the dependent variables (Y). Indeed, its goal is to predict the dependent variables (Y), selecting a limited number of new orthogonal factors (T) and a set of specific loadings (P) able to simultaneously both X and Y, in order to reduce the high dimensionality of the input dataset [8,22,33,36,[61][62][63][64]. Thus, X is calculated by applying the following equation: T is called the "score matrix" and it is composed by the latent variables (LVs), while P is called the "loading matrix". It is important to underline that P is not orthogonal while T has to maximize the covariance between X and Y with the following property (Equation (4)): where I is the identify matrix. Finally, Y is estimated by applying the following equation: where B is a diagonal matrix, which has regression weights as diagonal elements, and C is the "weight matrix" of the dependent variables.
The goodness of fit of the PLSR model depends on the number of selected LVs. Typically, it increases until it reaches a certain number of LVs before declining again. For that reason, it is fundamental to select the optimal number of LVs for each model, which allow us to obtain the best quality [64][65][66]. In this study, the optimal number of LVs for each PLSR model per grassland trait was estimated based on the root mean square error (RMSE) for the leave-one-out cross-validation (LOOCV). The criterion to add an additional LV to the model was that it had to reduce the RMSE of LOO cross-validation with >2%. The RMSE of the LOO (RMSELOO) was computed using the following equation: where ̅ and yi represent the leave-one-out predictions and observed values of the traits under investigation, respectively. Indeed, the model was fine-tuned without taking into account the observations i. Subsequently, the values of these observation were predicted by applying this model.
Although recent studies have shown that the optimal use of PLSR can be achieved by splitting the original dataset into a training dataset and an independent test set [67], in this work, the use of LOOCV, applied to a single calibration dataset, was necessary because we did not have a large enough experimental dataset to divide the samples into separate training and validation subsets. The assessment of the quality of the model was estimated using R 2 between measured and predicted values from the LOO crossvalidation and the RMSELOO, while the usefulness and goodness-of-fit of calibration models were performed by residual predictive deviation (RPD), calculated as the ratio between the standard deviation of the reference dataset and the RMSE. On the contrary, there is not a specific statistical rule for setting the RPD thresholds and, for these reasons, in this paper, the commonly used classification in literature was adopted [51][52]. It involves three levels of classification: RPD < 3 is not useful because it means low predictive power; RPD > 3 is suitable for screening; RPD > 5 is suitable for detecting the accurate analysis (laboratory analysis).
All regression analysis was done using R (R_Core_Team 2013), by applying the pls package [68] to the three sub-datasets defined before.

Grassland Traits
As a first step, the correlation between the grassland traits under study was investigated for the analyzed traits of the May campaign (Table 2). Only the values for May are reported because the Pearson coefficients of the two analyses are nearly similar. The results show that the structural variables (height, fresh and dry biomass) are significantly correlated to each other because the Pearson correlation values are higher than 0.89. In contrast, the biochemical characteristics show a substantially lower correlation which is trait-dependent. For example, Na shows a low correlation with metabolic energy and crude fiber, a moderate correlation with K and crude ash, and a high correlation with crude protein. The correlation between structural traits and biochemical traits is moderate, except for Na and crude protein which show a low correlation. In addition, all traits are positively correlated except for metabolic energy which shows an inverse correlation with the other traits. Harvest statistics of the structural and biochemical traits related to the May harvest are shown in Table 3. It shows that grassland biochemical characteristics have low variability, with the exception of Na. Metabolic energy shows the lowest coefficient of variation (CV) while K, crude fiber, and crude ash are also characterized by a small variability.  Table 4 illustrates the harvest statistics of all traits from the October harvest. It shows that especially the variability of the biochemical grassland characteristics is relatively small, with the exception of Na. Analyzing more in depth, it is possible to see that the metabolic energy shows the lowest variability next to K, crude fiber, and crude ash. For the structural traits, grassland height has the lowest variation.
Comparing Tables 3 and 4, it is apparent that the range of the grassland traits between the two harvest moments shows different patterns for the structural and biochemical traits. For all three structural traits, the range of values for October is significantly lower compared to May. In contrast, the biochemical traits are comparable, with the exception Na. This is also confirmed by Figure 4 presenting the effect of the applied fertilizer treatments on the level and variability of grassland traits. The influence of the type and amount of fertilization is clearly expressed and differs depending on the structural and biochemical traits under research and on the growth status of the grassland. The May data show the increase of structural trait values with the increase of the inorganic fertilization level until 230 kgN/ha; after that threshold, they have a small decrease again ( Figure 4). The only biochemical feature with the same pattern of variation as the structural traits is crude fiber. On the contrary, the biochemical variables all increase their values with increasing amounts of inorganic fertilization. The metabolic energy pattern of variation deserves a separate mention. Indeed, its values decreased with the increasing inorganic fertilization until 230 kgN/ha, after which its value increased. Regarding the influence of organic fertilization treatments, all the structural and biochemical features have the same pattern, with the exception of metabolic energy, which shows a dip at 230 kgN/ha.
The pattern of variation of the structural and biochemical traits of the grassland from the October observation under the different treatments is substantially different. Also, all of the structural traits show the same pattern of variation: under the influence of inorganic fertilization, structural traits increased, with the exception of a dip at 230 KgN/ha. Under the organic manure fertilization, they show an increasing pattern. On the contrary, it is not possible to identify a common pattern of variation for all the biochemical characteristics. The inorganic fertilizer shows the same influence on crude ash and crude protein, a pattern of variation which shows a decrease up to 115 kgN/ha followed by an increase. However, when considering the influence of the organic manure on crude ash and crude protein, it is possible to see that the pattern of variation is clearly different, where crude protein peaks at 230 kgN/ha while crude ash decreases. Metabolic energy, like potassium, is more or less constant over all treatments, while sodium shows variation independent of the amount of fertilizer treatment.  The significance levels from applying a Student's t-test between organic and inorganic fertilized plots in May, in October, and between May and October data are shown in Table 5. These results show that both grassland structural and biochemical traits are strongly influenced by the plant-growing season but not by the type of fertilizer. Moreover, the plant-growing season also influences the result of the Student's t-test based on organic and inorganic fertilizer in May and October.

Hyperspectral data
For both the May and October data, the spectral variation within the acquired images is relatively small, especially in the visible part of the spectra, while it is more noticeable in the NIR part. The maximum and the minimum value of both average reflectance spectra are more or less the same and only in the NIR part is a small difference noticeable. Indeed, even if the red edge slopes do not show differences, the maximum value of the average reflectance spectra of May is slightly higher than the maximum value of the October profile, and this trend is confirmed in all of the NIR part expect for 950 nm, where the trend is reversed. In addition, in the NIR part, another difference can be observed: NIR reflectance increases with an increase in the wavelength in October, but it is spectrally flat in May. Also, the minimums of the average reflectance spectra are comparable; however, in the NIR region, the profile of October is slightly higher than the profile of May. Observing the spectral profile of the average of the three different organic levels for May and October surveys, the same trends can be recognized ( Figure 5). Moreover, increasing the concentration of organic manure influences the overall level of reflectance with a maximum increase between 170 and 230 kgN/ha for the May dataset, while in October, this is between 230 and 340 kgN/ha.  Figure 6 shows the correlogram between the average reflectance spectra of hyperspectral images and previously selected grassland traits for May. The trends of the two structural traits, height and fresh biomass, are highly comparable, and the NIR area is the most relevant region for these traits. In fact, the NIR area (700-950 nm) is sensitive to LAI, grassland structure, and canopy thickness of the vegetation and, therefore, the correlation with height and fresh biomass is relatively high. For the same reason, NIR is also the most important region for estimating metabolic energy, since it is inversely proportional to the amount of green vegetation and, thus, the correlation is high but negative. Instead, crude protein needs a separate mention because the most interesting region of the spectra is the VIS area (400-700 nm), which is more sensitive to the absorption of chlorophyll at the leaf level. It turns out that the correlation between grassland traits and wavelengths is significantly influenced by the growth status of the grassland. This is also confirmed by the correlogram for the October harvest, which has lower correlations (<0.6) for all grassland traits.

Narrow-Band Vegetation Indices
Based on the extracted average reflectance spectra for all the plots in the experiment (Figure 1), the selected narrow-band vegetation indices were calculated (Table 1). In order to investigate the influence of the fertilizer regime, the variation of the vegetation indices as a function of type and amount of fertilization was explored for both moments of harvest (Figure 7). The variation of the vegetation indices related to the different levels of inorganic fertilization is more or less the same within the same survey while it is different between the May and October campaigns. On the contrary, the variation of the vegetation indices under the influence of the organic manure is the same for all indices and for both campaigns. Indeed, for the May harvest, the value of the indices increases with an increase in the amount of fertilizer up to 170 kgN/ha, then it decreases. In contrast, all vegetation indices, computed from the October data as a function of inorganic fertilization, show a decreasing trend up to 115 kgN/ha and an increasing value after that point. Instead, the tendency of the vegetation indices as a function of the organic fertilization is totally different: it decreases at the value of 230 kgN/ha and then it increases again. The only exception is represented by MCARI/OSAVI, which shows a peak at 230 kgN/ha for the May campaign but not for the October survey. Table 6 summarizes the results of the linear relation between narrow-band vegetation indices and the grassland traits for two sub-datasets: (1) the integrated dataset of May and October and (2) the sub-dataset of May only. For both datasets, all organic and inorganic fertilization plots are included. Analyzing the R 2 and RMSE values resulting from the linear regression models, a lot of variation between the grassland traits can be observed. In general, the variation in structural traits (height, fresh and dry biomass) could be better explained by the vegetation indices than the biochemical traits. For example, the R 2 between height and MCARI/OSAVI was 0.599, while it was only 0.49 between crude fiber and NDRE. However, in general, R 2 values for the biochemical grassland traits are considerably lower (Table 6). For the structural traits, the highest correlation was observed for NDRE and MCARI/OSAVI for the May dataset. When the May and October datasets were combined, MCARI/OSAVI and CIred-edge produced the highest correlation. Overall, it can be stated that the relations are not influenced by the type of fertilization regime applied. On the contrary, the growth status of the grassland strongly on the R 2 value.   Nevertheless, considering the classification suggested by [28], it is evident that there is a moderate relation between the grassland traits and the vegetation indices in both subsets, even if the results obtained for the structural traits of grassland are sufficiently promising.
In reality, it is not possible to determine the best predictor, but the prediction ability of the vegetation indices depends on the feature investigated and the sub-dataset considered. Figure 8 illustrates the best results for some grassland traits, chosen as examples.
The different shades of color are indicative of the variations: the green plots correspond to areas with the lowest values, while the red points are related to the regions characterized by higher values.  Two examples of grassland trait maps are presented in Figures 9 and 10, showing the vegetation variation within plots based on the linear regression relation for the best narrow vegetation index on the dataset of May for height and metabolic energy, respectively. Height is quit homogeneous over all plots (>34 cm), however, with clearly lower values (<30 cm) for the non-fertilized plots (see arrows down in Figure 9), while other lower fertilized plots also show decreased height values. However, highly fertilized plots can also have low height values (see arrows up in Figure 9), which can be explained by the falling over of grassland stems when the biomass gets too heavy for the stems to support it. For metabolic energy, a larger spatial variation can be observed between the plots (Figure 10).

Partial Least Squares Regression (PLSR)
Partial least squares regression (PLSR) was used for modeling the relationships between grassland traits and the average reflectance spectra of the two integrated datasets of plots (Table 7) in order to detect its prediction capacity and the influence of the type and the amount of fertilizer.   On the basis of the classifications suggested by [28], it is possible to point out in which cases and with which grassland features a strong, moderate, or weak correlation is verifiable. The best results (R 2 > 0.7) were obtained on the integrated datasets of May and October, both considering the levels of inorganic fertilization and all treatments of organic and inorganic fertilizations. Height, fresh biomass, dry matter, crude protein, crude fiber, and metabolic energy show a strong correlation, while the others (crude ash, K, Na) show a weak correlation (R 2 < 0.5). The situation is different when the subject under investigation is the sub-dataset related to all levels of both organic and inorganic fertilizations of May. In this case, only height and fresh matter yield show a high correlation; crude fiber, Na, and metabolic energy have a low correlation, while dry matter, crude ash, crude protein, and K are characterized by a moderate correlation.  These results are also confirmed by the RPD. In fact, the best results were obtained in the same sub-dataset which showed a high correlation. In particular, in the integrated dataset of May and October, including all treatments of organic and inorganic fertilizations, height, dry matter, and fresh matter yield show a RPD suitable for screening; instead, all the other features are characterized by a low prediction value. The sub-dataset related to all levels of both organic and inorganic fertilization of the May dataset presents a RPD value suitable for laboratory analysis regarding K prediction. However, it is important to emphasize that the low RPD values are not really small because they are around 2.0 in all sub-datasets. Figure 11 illustrates the scatterplot of measurement and predicted values for the same four variables selected in Figure 8.
The PLSR model allowed us to identify and select the most important wavelengths for the variables under investigation. In particular, Figure 12 shows the wavelengths which mostly influenced predictions of height, fresh matter, crude protein, and metabolic energy. For each trait analyzed, the most interesting wavelengths are located in the first part (450-545 nm) of the visible part of the spectrum, except for metabolic energy, where 545 nm is the central wavelength of the region of interest. For crude protein, it is possible to identify another area of interest around 645 nm (545 nm-765 nm).

Discussion
In this paper, the best strategy for estimating structural and biochemical traits of grasslands from UAV-acquired hyperspectral images and the study of the type and the amount of the fertilization influence on their prediction capacity was explored. The application of nitrogen (N) to grassland is one of the major management tools which farmers have to influence yield and quality. Currently, there are knowledge gaps in how organic fertilizer application affects N balances in grasslands. Especially weather conditions can influence nitrogen mineralization of the organic N of slurry and thus influence the longterm efficiency of slurry application [69]. Therefore, to study this long-term effect of slurry, a grassland experiment was planned by applying different amounts of slurry for the years of 2012, 2013, and 2014. The level of fertilizer applied was chosen based on the customs of German farmers. Moreover, this data combination seemed especially eligible for studying the best statistical strategy for the estimation of structural and biochemical traits from UAV-acquired hyperspectral images as the different levels of N fertilization provided a wide, well-defined range of both biomass yield and N concentration at every flying date.
The hyperspectral images of an experimental grassland field at the farm Haus Riswick, near Kleve in Germany (Figure 1), were taken by applying an octocopter UAV, equipped with the Wageningen UR Hyperspectral Mapping System (HYMSY) (Figure 3) [51]. Two different flight campaigns were carried out, one on 15 May 2014 and the other on 14 October 2014. Following from the flight campaigns, the structural and biochemical grassland traits analysis was also carried out.
The results of the grassland traits analysis in both campaigns showed different values for each of the parameters under investigation, but similar statistical relationships (Tables 2-4). The Pearson correlation among the structural traits (height, fresh and dry matter yield) is significantly high (>0.89), while among the biochemical traits, the correlation is substantially different depending on the feature considered ( Table 2): for example, metabolic energy and crude fiber are negatively correlated (−0.98), while metabolic energy and crude protein show a very low correlation (0.04). Indeed, the metabolic energy has not been directly measured but it has been calculated by applying Equations 1 and 2, and surely this influences the correlation. On the contrary, the correlation between crude protein and crude fiber is quite low (0.10). Normally, when the analysis is carried out on a growing plant at subsequent dates, they are negatively correlated [70]. In this case, however, the grassland was only analyzed at one harvest date and, therefore, it is not possible to find a good correlation among them. This is also confirmed by the analysis of the N fertilization influence. Indeed, N fertilization increases crude protein concentration while normally it has little effect on the crude fiber concentration since crude fiber depends mainly on temperature, day length, and light intensity (Figure 4).
The use of a UAV, applied to adjust the spatial resolution to the size of the grassland field investigated, obtaining a GSD of 200 mm for the hyperspectral image, due to the reduced altitude (70 m) and speed (5 m/s) compared with the traditional quotas and speed of the airplane [18][19][20][21]. Moreover, it also permits us to acquire images in flexible dates, exploiting the same sunny illumination conditions and eliminating problems due to clouds, in such a way that the imagery of two different periods can be compared. In addition, the possibility to mount a hyperspectral camera on it has allowed us also to fully take advantage of all the potentials of remote sensing, processing hundreds of narrow contiguous spectral bands to which the vegetation is sensitive [4,[12][13][14][15][16]. Indeed, in general, quite homogeneous grassland vegetation cover has relatively small variations in reflectance properties of the canopy, which could be visible and analyzed only from a continuous reflectance spectrum. For these reasons, the field under investigation was analyzed by applying hyperspectral data. From this detailed spectral data, the variation in plant traits of grassland, caused by the type and the amount of fertilization ( Figure 5), correspond to a relative small variability in reflectance values ( Figure 6). Thus, in order to identify this variation, the narrow-band vegetation indices were calculated instead of broad vegetation indices, as shown in previous studies [14,15,[24][25][26][27][28]. Narrow-band vegetation indices improve the results in processing the hyperspectral images since they allow us to exploit all of the information contained in the spectrum [14,15,[24][25][26][27][28].
The selected narrow-band vegetation indices (Table 1) were computed separately on the three sub-datasets in which the original dataset was split in order to understand the influence of the type and the amount of fertilizer and the weight of the grassland growth status. Results showed that, although narrow vegetation indices detect the effects of the type and the amount of fertilization, their prediction ability is moderate. Indeed, the results of the two integrated sub-datasets, including both May and October data, but one characterized only by inorganic fertilization and the other by both organic and inorganic fertilization, are comparable. On the contrary, grassland health status influences the estimation ability of this technique (Table 6). This is perfectly in line with the results of Table 5. Indeed, the best results (higher correlation values) were obtained using a linear regression model on the data related to the more productive part of the season in May. This is justified by the distribution of the two datasets. Indeed, the distribution of data acquired during the two campaigns does not follow a specific linear relation, but it forms two separate clusters. For that reason, it is not possible to fix the best narrow vegetation index, since it depends on the grassland health status and the trait under analysis. In particular, the best estimator for grassland structural traits is MCARI/OSAVI [57], which shows a value > 0.5, while for the most biochemical features, it is CIred-edge [57,58] (Figure 8). Moreover, this approach shows promising results for estimating grassland structural traits, but not biochemical features, with the exception of metabolic energy. The UAVbased set-up allows the delivery of spatial maps showing the spatial distribution of grassland traits following the application of this statistical approach (Figures 9 and 10).
Like the narrow vegetation indices, the PLSR method was implemented on the same three sub-datasets, in order to analyze the consequence of the type and the amount of the fertilization and the influence of the grassland growing season. Also, in this case, results showed that the estimation ability of this approach is mainly influenced by the grassland growing season and not by the different type and the amount of fertilization. Indeed, the results related to the two different populations including both May and October, one with only the inorganic fertilization and the other with both organic and inorganic fertilization, are comparable. Instead, comparing one of them with the sub-dataset related to the combined data of the inorganic and organic fertilization data of the May survey, it is possible to see that the method performs better with the integrated sub-dataset of May and October, with the exception of crude ash, Na, and K. Indeed, in the integrated sub-dataset of May and October, including both inorganic and organic fertilization, the R 2 values are good (>0.7) for all characteristics, except for crude ash (0.4), Na (0.2), and K (0.3) (Table 7). Also, the value of the RPD is good (Table 7) while, in the sub-dataset of May including both inorganic and organic fertilization, the R 2 values are quite high for all traits, and regarding the RPD, only K has a value higher than 3 (Table 7). Probably, the higher coefficient of determination of the May-October group is related to a higher number of LVs, which could be a sign for overfitting. In addition, the leave-one-out cross-validation approach as adopted in this study cannot be treated as a completely independent assessment of prediction quality. In order to generate a more robust validation and representative error estimates, a larger dataset is necessary so that a training and independent test dataset can be used to test the PLSR approach, as suggested in recent studies [67]. Therefore, for future application of the described approach, an additional evaluation using an independent test set would be required to assess the use of PLSR in practice. Moreover, other future applications could contemplate the PLSR use for integrating different vegetation indices over the whole spectrum in order to evaluate their performance and to compare them with the performance of the vegetation indices computed on pure spectral reflectance values.
From the above results, the best strategy for detecting structural and biochemical features is hyperspectral remote sensing in combination with PLSR. The outcome from this study is in agreement with the results of previous research studies [3,13,37,38]. For example, [3] compared a vegetation index approach with PLSR for estimating fresh grassland biomass. In their study, PLSR gave a better retrieval of fresh biomass with a RMSECV of 2.10, slightly better compared to the findings from this study (Table  7). However, for biochemical grassland traits, this finding gives new opportunities to map grassland quality before harvest, which is one of the requirements for the future adoption of precision farming practices in grassland management [3].
In order to confirm the results of this paper and achieve a schedule of precision grassland management, the approach presented in this study should be tested on other fields of different size and geographic position. In addition, more additional points need to be investigated. For example, a possible future development could consider more observations over the growing season, since this study has taken into account only the beginning and the end of the growth.

Conclusions
The results of this study have shown that hyperspectral data acquired from an UAV platform can be adopted to characterize both structural and biochemical traits of grasslands. Therefore, using the resulting spatial information, effective management measures can be taken to maintain the health status of grassland ecosystems. Moreover, the continuous reflectance spectrum of hyperspectral images has shown a high potentiality to describe the small variations in reflectance properties in the canopy of the quite homogeneous grassland vegetation cover. For these reasons, hyperspectral images have been very useful for detecting the small variability in reflectance values caused by the different types and the amounts of fertilization.
This study has also demonstrated that both statistical methods analyzed present promising results for estimating the structural traits, but not for assessing biochemical features. Indeed, only the PLSR model has displayed a good performance in both cases. In addition, both statistical approaches were not influenced by the amount and the type of fertilization.
In conclusion, therefore, the PLSR model based on airborne hyperspectral imagery is the best strategy for estimating grassland traits, even if more studies are necessary to confirm the results obtained in this paper.