Combined Use of Low-Cost Remote Sensing Techniques and δ13C to Assess Bread Wheat Grain Yield under Different Water and Nitrogen Conditions

Vegetation indices and canopy temperature are the most usual remote sensing approaches to assess cereal performance. Understanding the relationships of these parameters and yield may help design more efficient strategies to monitor crop performance. We present an evaluation of vegetation indices (derived from RGB images and multispectral data) and water status traits (through the canopy temperature, stomatal conductance and carbon isotopic composition) measured during the reproductive stage for genotype phenotyping in a study of four wheat genotypes growing under different water and nitrogen regimes in north Algeria. Differences among the cultivars were reported through the vegetation indices, but not with the water status traits. Both approximations correlated significantly with grain yield (GY), reporting stronger correlations under support irrigation and N-fertilization than the rainfed or the no N-fertilization conditions. For N-fertilized trials (irrigated or rainfed) water status parameters were the main factors predicting relative GY performance, while in the absence of N-fertilization, the green canopy area (assessed through GGA) was the main factor negatively correlated with GY. Regression models for GY estimation were generated using data from three consecutive growing seasons. The results highlighted the usefulness of vegetation indices derived from RGB images predicting GY.


Introduction
Bread wheat is one of the most cultivated herbaceous crops in the Mediterranean region [1], with water stress and low nitrogen fertility being the main constraints limiting productivity [2]. These limitations are likely to increase in the future because climatic change is expected to decrease precipitation and increase evapotranspiration in the Mediterranean region [3]. Increasing productivity in these semi-arid environments depends on the efficiency of crop management [2] and breeding [4], where efficient and affordable methodologies to monitor crop performance, or to assess phenotypic variability for breeding, are needed. Remote sensing techniques at the canopy level have become valuable tools for precision agriculture and high throughput phenotyping [5][6][7]. Thus, both spectral and thermal approaches have been proposed as potential solutions to identify crop N status and water stress across large areas [8,9]. In this way, these techniques can help farmers to practice more thermometry has been proposed as a low-cost approach in crop management to enable scheduling of support irrigation [41], to assess spatial soil heterogeneity [42], or to evaluate genotypic performance to drought [43]. However, potential interaction effects between N fertilization and water regime should be considered. In particular, how does N fertilization affect the water by temperature relationships and even how does water affect the interaction of N fertilization with vegetation indices. This is not trivial because haying-off, which is the negative effect of nitrogen fertilization on productivity caused by an imbalance between transpired biomass and the available water, is regarded as a potential problem for wheat cultivation in Mediterranean regions [44]. As a consequence, unexpected relationships between remote sensing readings and grain yield may occur.
Similarly, other physiological characteristics related to plant water status, such as stable carbon isotope composition (δ 13 C; frequently measured as discrimination from the surrounding air, ∆ 13 C) are also often used for evaluating genotypic performance under water stress [40,44] or even to monitor spatial variability and water status [45]. The natural 13 C abundance in plant matter provides time-integrated information of the effects of water stress on the photosynthetic carbon assimilation of C3 species, including wheat [46][47][48]. Conditions inducing stomatal closure (e.g., water deficit or salinity) restrict the CO 2 supply to carboxylation sites, which then decreases the ∆ 13 C (or increases the δ 13 C) of plant matter [47,49]. Under Mediterranean conditions the δ 13 C of mature kernels is better correlated with grain yield than the δ 13 C of other plant parts [50]. The costs of these analyses have decreased throughout the years, making their analysis increasingly feasible.
The objective of this study was to assess the grain yield performance of wheat under a range of water and fertilization conditions in the Mediterranean, using different low-cost remote sensing approaches to assess canopy green biomass (NDVI and vegetation indices derived from conventional RGB images), and characteristics associated with plant water status, (canopy temperature depression), together with additional traits informing on the water status (δ 13 C of mature grains and the stomatal conductance of the flag leaf). The novelty of the study centers on (i) testing how different low-cost, user-friendly remote sensing techniques may contribute to site-specific wheat management and eventually to the prediction of yield across seasons; and (ii) how interactions between growing conditions (water regime and N fertilization) may affect the predictive strength of these techniques. Moreover, to better explore the potential usefulness of our study for wheat phenotyping we developed a conceptual model of how the combination of these different traits explains genotypic variability in grain yield under different combinations of water regimes and nitrogen fertilization.

Plant Material and Growing Conditions
Field trials were conducted during the 2014-2015 crop season at Bir Ould Khelifa in the area of Khemis Miliana commune, approximately 230 km to the south west of Algiers (Algeria) and with geographical coordinates 36 • 11 50.23" N and 2 • 13 17.69" E. This commune receives an average rainfall of between 400 and 450 mm, and it is characterized by clay-silty fertile soils, with high organic matter content and high levels of total and mineral nitrogen (Table S1). Monthly total accumulated rainfall and temperatures for the study region for the 2014-2015 crop season are presented in Table S2. Four bread wheat (Triticum aestivum L.) genotypes were planted on 14 December 2014. The wheat genotypes were "Ain-Abid" and "Arz" (modern varieties) and "Wifak" and "Maaouna" (local varieties). The experimental design was a split-split-plot, with the main-plot factor being water regime, the subplot factor was the N amount and the sub-subplot factor was genotype ( Figure S1). A total of 108 plots (four genotypes, three replicates per genotype, three water regimes, and three nitrogen fertilization treatments) each with a size of 10 m × 1.2 m and six rows, 20 cm apart, were studied. The three water regimes consisted of one rainfed and two support irrigation treatments of 30 mm (SI-1, a single amount of supplemental irrigation) and 60 mm (SI-2, a double amount of supplemental irrigation) aimed at providing water in the typical range for the agronomic practices in this area. Supplementary irrigation Agronomy 2019, 9, 285 4 of 23 was applied with sprinklers around the beginning of stem elongation. For SI-1, one irrigation (30 mm) was applied at the beginning of stem elongation (31 Zadoks stage), whereas for SI-2, a second irrigation was also delivered at heading (58 Zadoks stage). Nitrogen fertilization was applied using urea fertilizer, and treatments consisted of no fertilizer (N0) and 60 kg ha −1 and 120 kg ha −1 of nitrogen fertilizer (N60, N120 respectively). Application of N fertilization was achieved at two growing stages, tillering and jointing (26 and 31 Zadoks stages). Plots were harvested with a sickle after physiological maturity and grain yield was estimated. Thousand kernel weight and the number of kernels per square meter (Kernel m −2 ) were evaluated.

Vegetation Indices
Remote sensing measurements were performed on a sunny day around anthesis (26 March 2015) between 10 h and 15 h solar time. The NDVI was determined with a portable spectroradiometer (GreenSeeker handheld crop sensor, Trimble, Sunnyvale, CA, USA). The NDVI is formulated using the following equation: (NIR − R)/(NIR + R), where R is the reflectance in the red band and NIR is the reflectance in the near-infrared band. The distance between the sensor and the plots was kept constant at around 50-60 cm above and perpendicular to the canopy. Additionally, one conventional digital picture was taken per plot, holding the camera about 80 cm above the plant canopy in a zenithal plane and focusing near the center of each plot. The images were acquired with an Olympus E-M10 camera (Olympus Corporation, Tokyo, Japan), using a 14 mm lens, triggered at a speed of 1/125 seconds with the aperture programmed in automatic mode. The size of the images was 4608 × 3456 pixels stored in JPG format using RGB color standard [51]. Pictures were analyzed with the free-access BreedPix 0.2 software, now integrated within the CerealScanner plugin (https://integrativecropecophysiology.com/software-development/cerealscanner/), from the Mediterranean Crop Ecophysiology Group, University of Barcelona [22], which was developed for digital image processing. This software quickly provides digital values from different color properties. The vegetation indices measured were the green area (GA) and the greener area (GGA). GA is the portion (as a %) of pixels with 60 < Hue < 120 from the total amount of pixels, whereas greener area is formulated as the % of pixels with 80 < Hue < 120 [31,52]. GGA is designed to capture the active photosynthetic area excluding senescent leaves. In addition, the leaf chlorophyll content of five flag leaf blades per plot was measured using a Minolta SPAD-502 portable meter (Spectrum Technologies Inc., Plainfield, IL, USA).

Canopy Temperature Measurements
Canopy temperature (CT) was measured at noon (12 h-14 h), on the same day around anthesis as the vegetation indices, using an infrared thermometer (PhotoTempTM MXSTMTD, Raytek®, Santa Cruz, CA, USA). Measurements were taken above the plants, pointing towards the canopy at a distance of about one meter and having the sun towards the rear. The air temperature was measured simultaneously for each plot with a temperature humidity meter (Testo 177-H1 Logger, Testo, Lenzkirch, Germany) and employed for the calculation of the canopy temperature depression (CTD) as the difference between the ambient and the canopy temperature.

Stomatal Conductance
Stomatal conductance was measured on the flag leaves on the same days as the remote sensing traits. Two measurements per plot were taken around noon (12 h-14 h), using a Decagon SC-1 Leaf Porometer (Decagon Devices Inc., Pullman, WA, USA).

Stable Carbon Isotope Composition
Carbon isotope composition was analyzed in mature grains using an Elemental Analyzer (Flash 1112 EA; ThermoFinnigan, Bremen, Germany) coupled with an Isotope Ratio Mass Spectrometer (Delta C IRMS, ThermoFinnigan, Bremen, Germany) operating in continuous flow mode to determine the Agronomy 2019, 9, 285 5 of 23 stable carbon ( 13 C/ 12 C) isotope ratios. Samples of about 0.7 mg of dry matter and reference materials were weighed into tin capsules, sealed, and then loaded into an automatic sampler (ThermoFinnigan, Bremen, Germany) prior to EAIRMS analysis. Measurements were carried out at the Scientific Facilities of the University of Barcelona. The 13 C/ 12 C ratios were expressed in δ notation determined by: δ 13 C = ( 13 C/ 12 C) sample /( 13 C/ 12 C) standard − 1 [53], where sample refers to plant material and standard to Pee Dee Belemnite (PDB) calcium carbonate. International isotope secondary standards of known 13 C/ 12 C ratios (IAEA CH 7 polyethylene foil, IAEA CH 6 sucrose, and USGS 40 L-glutamic acid) were used for calibration to a precision of 0.1% .

Statistical Analysis
Data were subjected to factorial ANOVA to test the effects of the growing conditions (water regime and nitrogen fertilization), genotype, and their interaction. Mean comparisons were performed using Tukey's honestly significant difference (HSD) test. Pearson correlation coefficients between grain yield and all different traits were calculated. Multiple linear regression analysis (stepwise) was used to analyze grain yield under different growing conditions. Data were analyzed using IBM SPSS Statistics 24 (SPSS Inc., Chicago, IL, USA). Figures were created using Sigma-Plot 11.0 for Windows (Systat Software Inc., Point Richmond, CA, USA).
The performance of the different remote sensing traits in predicting yield performance across the seasons was evaluated using data from the present study, together with data already published by our team from the two previous seasons (2012-2013 and 2013-2014) related to grain yield and different remote sensing traits measured in a set of genotypes grown under different water regimes [54]. The equations of the linear relationships between these traits and grain yield (GY) determined in the first crop season (2012-2013) were further tested during the two following crop seasons 2013-2014 and 2014-2015 (the latter of these conducted during the present study). Predicted and measured grain yields were expressed as relative values; grain yields were normalized with regard to the highest yield combination (genotype and growing condition), then the means of the three replicates per genotype were calculated. Finally, we performed path analyses to quantify the relative contributions of direct and indirect effects of water status (δ 13 C, g s and CTD) and vegetation indices (GA, GGA) on grain yield. This methodology offers the possibility of building associations between variables on the basis of prior knowledge. Mechanisms that play potential roles in grain yield variation and involving traits that exhibited genotypic differences have been proposed, as detailed in the conceptual model displayed in Figure S2. This model was aimed at understanding grain yield responses to genotypic differences under different levels of nitrogen fertilization (N0, N60, N120) and under different water regimes (supplementary irrigation and rainfed). A model with a comparative fit index (CFI) [42] with values > 0.9 was taken as indicative of a good fit. Data were analyzed using IBM, SPSS, Amos 21 (SPSS Inc., Chicago, IL, USA).

Irrigation and Fertilization Effects on Grain Yield
Both irrigation and nitrogen fertilization significantly affected the grain yield (GY) and the agronomic yield components ( Table 1). The doubled amount of supplemental irrigation (SI-2 with all nitrogen fertilization combined) and the highest nitrogen fertilization N120 (120 kg N ha −1 of fertilizer with all irrigations combined) were the growing conditions that exhibited the highest grain yield. Significant interactions only existed for grain yield between genotype and water regime (p = 0.049).

Vegetation Indices
The water regime significantly affected the NDVI (p < 0.001), as well as the green area (GA) (p < 0.001) and the greener area (GGA) (p < 0.001) indices ( Table 2). The values of these three vegetation indices were highest under SI-2 compared to the single amount of supplemental irrigation (SI-1) and rainfed conditions. However, no difference was observed in vegetation indices across fertilization treatments except for GGA (p = 0.005). Moreover, leaf chlorophyll content (LC) slightly increased ( Table 2) under SI-1 and rainfed conditions compared to SI-2 and under N60 and N120 compared to N0. Interactions were not significant, whatever the combination (genotypes, water regimes and fertilization levels) or variables considered.

Canopy Temperature Depression, Stable Carbon Isotope Composition and Stomatal Conductance
Water regime affected significantly the canopy temperature depression (CTD) (p < 0.001), the stomatal conductance (g s ) (p < 0.001) and the stable carbon isotope composition (δ 13 C) (p < 0.001) of mature grains. Rainfed conditions decreased g s and CTD, whereas δ 13 C increased compared to support irrigation conditions (Table 3). Nevertheless, fertilization treatments did not affect any of these three parameters (Table 3). No significant interactions were observed except for the δ 13 C between the water regime and nitrogen fertilization (p = 0.022). Table 3. Mean values of genotypes, water regimes, and nitrogen fertilization levels for canopy temperature depression (CTD), stomatal conductance (g s ) of the flag leaves and the stable carbon isotope composition (δ 13 C) of the mature grains. Means followed by different letters are significantly different (p < 0.05) according to Tukey's honestly significant difference (HSD) test.

Genotypic Effect on Grain Yield, Vegetation Indices and Water Status Traits
The genotypic effect was significant (p < 0.05) for grain yield (GY) and (p < 0.001) for vegetation indices (Tables 1 and 2) under the growing conditions analyzed together, whereas only the δ 13 C (in the case of water status traits) was significantly (p < 0.001) different between genotypes ( Table 3). The genotypic difference was also examined within each of the nine growing conditions, resulting from the combination of the three water regimes and the three nitrogen fertilization levels (Table 4). GY only showed a genotypic effect under SI-2 combined with nitrogen fertilization (either N120 or N60). However, under the SI-1 and rainfed conditions, GY did not show genotypic differences, regardless of the N fertilization regime. The vegetation indices, LC content, δ 13 C and g s also showed genotypic effects in some specific growing conditions while CTD did not (Table 4). Table 4. Genotype effect on the Normalized Difference Vegetation Index (NDVI) and the Green Area and the Greener Area (GA and GGA) vegetation indices, the leaf chlorophyll content (LC), the stomatal conductance (g s ) of the flag leaf, the canopy temperature depression (CTD), the stable carbon isotope composition (δ 13 C) of the mature grains and the grain yield (GY). Significance levels-ns, not significant; * p < 0.05; ** p < 0.01 and *** p < 0.001.

Relationships of the Grain Yield with the Vegetation Indices
Combining both irrigation (SI-2 and SI-1) conditions, the two RGB vegetation indices (GA and GGA) were positively correlated with GY at N120 and only GA at N60, while the NDVI was not correlated with GY ( Figure 1). Moreover, under rainfed conditions GA, GGA and the NDVI were correlated with GY at N60 (Figure 1). In the absence of nitrogen fertilization, no correlation was found between GY and the different vegetation indices (Table S3).

Relationships of the Grain Yield with the Canopy Temperature Depression and the δ 13 C
Grain yield was negatively correlated with the δ 13 C of mature grains under irrigation at both N120 and N60, and positively correlated with the δ 13 C of rainfed conditions at N120 (Figure 1). Moreover, GY was positively correlated with CTD under either irrigation (SI-1 and SI-2 combined) Figure 1. Relationships of grain yield (GY) with NDVI, GA, and GGA, δ 13 C of mature grains and the canopy temperature depression (CTD) under (left column) irrigation with N120, (central column) irrigation with N60, and (central column) rainfed conditions with N60. The irrigation data combines the two support-irrigation regimes (SI-1, SI-2). Significance levels-ns, not significant; * p < 0.05 and ** p < 0.01. Abbreviations for variables and growing conditions are as defined in Tables 1 and 2. 3.6. Relationships of the Grain Yield with the Canopy Temperature Depression and the δ 13 C Grain yield was negatively correlated with the δ 13 C of mature grains under irrigation at both N120 and N60, and positively correlated with the δ 13 C of rainfed conditions at N120 (Figure 1). Moreover, GY was positively correlated with CTD under either irrigation (SI-1 and SI-2 combined) or rainfed conditions at N120, while under N60 no association was found ( Figure 1 and Table S3). In the absence of N fertilizer, no correlations were found (Table S3).

Relationships of Vegetation and Water Status Indices
The NDVI correlated highly and positively with both the GA and GGA vegetation indices when all genotypes, growing conditions, and replicates analyzed were combined (Table 5). Moreover, the CTD was also significantly (p < 0.01) correlated with the other water status parameters; negatively with δ 13 C and positively with g s . Likewise, both vegetation and water status indices were significantly (p < 0.01) correlated; the NDVI, GA and GGA were positively associated with CTD and g s and negatively correlated with δ 13 C ( Table 5). Table 5. Correlation coefficients of the relationships between different indices used to measure crop biomass and water status. Parameters were measured at anthesis. Treatments and genotypes were analyzed together. Significance levels-** p < 0.010; *** p < 0.000. Abbreviations of variables as in Tables 1 and 3.

Vegetation Indices
Correlation

Grain Yield Estimation Using Vegetation Indices and Canopy Temperature
A general model using linear regression of the NDVI, GA and CTD with GY in the 2012-2013 crop season was developed in this study to estimate GY in the two successive crop seasons (2013-2014 and 2014-2015). Only the NDVI, GA and CTD were involved in this model because they were the three variables included in the stepwise models to explain the difference in GY in the present study (see Section 3.9 below) and were significantly (p < 0.01) correlated with GY in the first crop season [54]. With all growing conditions and genotypes analyzed together, the predicted and measured grain yields were positively and highly significantly (p < 0.01) correlated in 2013-2014 and 2014-2015 ( Figure 2) using any of the three parameters studied alone (NDVI, GA, CTD). However, the range of normalized values predicted was smaller, in general, using the NDVI (at 2014-2015 crop season) than either of the other two indices. Predicted and measured GY for the 2013-2014 crop season was also positively and highly correlated using the NDVI, GA and CTD ( Figure 3) under irrigation (SI-2 and SI-1 combined) with N fertilization (N120 and N60 combined) and without N fertilization (N0). For the 2014-2015 crop season, predicted and measured grain yields were also highly significantly (p < 0.01) correlated using either the NDVI, GA or CTD (Figure 3), but only under irrigation (SI-2 and SI-1 combined) with N fertilizer (N120 and N60 combined) and with no correlation under N0. In the absence of irrigation (rainfed conditions) and regardless of the N fertilization conditions, we did not find any association between the predicted and the measured grain yields (Figure 3). Grain yield estimation was also examined within each of the nine growing conditions resulting from the combination of the three water regimes and the three nitrogen fertilization levels ( Figure S3). Under any of the two irrigation conditions with and without N fertilizer and the rainfed conditions with N fertilizer, the predicted GY was positively correlated with the measured GY in both 2013-2014 and 2014-2015 using at least one of the NDVI, GA or CTD parameters ( Figure S3). Under rainfed conditions and N0 the predicted and measured GYs were not correlated in either crop season ( Figure S3).
conditions resulting from the combination of the three water regimes and the three nitrogen fertilization levels ( Figure S3). Under any of the two irrigation conditions with and without N fertilizer and the rainfed conditions with N fertilizer, the predicted GY was positively correlated with the measured GY in both 2013-2014 and 2014-2015 using at least one of the NDVI, GA or CTD parameters ( Figure S3). Under rainfed conditions and N0 the predicted and measured GYs were not correlated in either crop season ( Figure S3).   The irrigation with N values combine both irrigation regimes and the N120 and N60 treatments, irrigation without N combines the two irrigation regimes without N fertilization, and the rainfed conditions with N combine the two levels of N fertilization. Significance levels-ns, not significant; **p < 0.01 and ***p < 0.001. Abbreviations of variables are as in Tables 1, 2 and 3.

The Combined Effect of Remote Sensing Indices and Physiological Traits in Explaining Grain Yield
Stepwise regressions were performed for the irrigation conditions combined (SI-1 and SI-2) and the rainfed conditions, and considered both groups under each of the three N fertilization levels using grain yield (GY) as the dependent variable and the vegetation indices (NDVI, GA, GGA and LC) and water status traits (CTD, gs, and δ 13 C) as independent variables (Table 6). Except for the rainfed conditions with N60, the first trait selected by the model to explain GY was related to the plant water The irrigation with N values combine both irrigation regimes and the N120 and N60 treatments, irrigation without N combines the two irrigation regimes without N fertilization, and the rainfed conditions with N combine the two levels of N fertilization. Significance levels-ns, not significant; ** p < 0.01 and *** p < 0.001. Abbreviations of variables are as in Tables 1-3.

The Combined Effect of Remote Sensing Indices and Physiological Traits in Explaining Grain Yield
Stepwise regressions were performed for the irrigation conditions combined (SI-1 and SI-2) and the rainfed conditions, and considered both groups under each of the three N fertilization levels using grain yield (GY) as the dependent variable and the vegetation indices (NDVI, GA, GGA and LC) and water status traits (CTD, g s , and δ 13 C) as independent variables (Table 6). Except for the rainfed conditions with N60, the first trait selected by the model to explain GY was related to the plant water status, while under rainfed conditions the NDVI was the first and only trait chosen by the model (Table 6). Furthermore, the genotypic differences in grain yield within each of the nine growing conditions (resulting from the combination of the three water regimes and the three nitrogen fertilization levels) were assessed through a stepwise model having GY as the dependent variable and any of the vegetation indices and water status traits as independent variables. The model identified at least one trait positively correlated with GY in only four of the nine growing conditions (Table S4).
Additionally, a conceptual model based on a path analysis was proposed ( Figure S2) that separated direct acclimation responses in grain yield related to water status traits and vegetation indices (through GA and GGA). The three water status traits were included in the model because they represent different scales: Temporal (δ 13 C), individual organ (g s ) and canopy (CTD). Concerning the vegetation indices, the NDVI was discarded because GA (whole photosynthetic biomass) and GGA (non-senescent biomass) already tracked the same parameter. The final objective of the model was to dissect how these physiological traits may have directly or indirectly assessed GY performance within different growing conditions. The four path models proposed ( Figure 4) provided an acceptable fit to the data (CFI > 0.9 in all cases). Under irrigation and with and without N fertilizer, g s had a strong and negative association with δ 13 C and strong and positive associations with CTD and GA. Significant paths corresponding to a direct (negative) association of δ 13 C with GY were also observed in irrigation conditions without N fertilizer and without, while CTD had a positive association with GY only under irrigation with N fertilizer. A direct positive and strong association of GA with GGA was observed under irrigation. GGA was in turn strongly and negatively associated with GY only under irrigation without N fertilizer. Under rainfed conditions, GY was not associated with δ 13 C, but was positively associated with CTD and GGA with N fertilizer, and negatively associated with GGA without fertilizer. The irrigation with N values combine both the SI-1 and SI-2 irrigation regimes and the N120 and N60 nitrogen treatments, irrigation without N combines the two irrigation regimes (SI-1 and SI-2) without N fertilization, and the rainfed conditions with N combine the two levels of N fertilization (N60 and N120). Physiological parameters included in the model are: The stomatal conductance (gs), the stable carbon isotope composition (δ 13 C) of mature grains, the Relative Green Area (GA) and the Relative Greener Area (GGA), and the canopy temperature depression (CTD). The width of the arrows is proportional to the path coefficient values. Dashed lines indicate negative relationships. CFIs with values > 0.9 are taken as indicative of a good fit. Significance levels-*p < 0.05; **p < 0.01 and ***p < 0.001.
Path analysis was also examined within each of the nine growing conditions (resulting from the combination of the three water regimes and the three nitrogen fertilization levels ( Figure S4). Water irrigation conditions with and without N fertilization indicated the dependence of grain yield on the water status traits δ 13 C and CTD (and the latter also under rainfed conditions with N fertilizer). In the absence of N fertilization under both irrigation and rainfed conditions, GY seemed strongly and negatively dependent on the vegetation index GGA.

Discussion
The grain yields achieved, even under the best growing conditions (SI-2 with N120) are below 3 tonnes per hectare, and this clearly corresponds with moderate-low yielding conditions in the Mediterranean [44,54]. The vegetation indices tested in this study, generated either spectroradiometrically (NDVI) or derived from RGB images, performed well when assessing differences in water conditions. The efficacy of these indices in capturing differences in growth and senescence in response to water regime has been reported for wheat already [32,33,54]. The choice of The irrigation with N values combine both the SI-1 and SI-2 irrigation regimes and the N120 and N60 nitrogen treatments, irrigation without N combines the two irrigation regimes (SI-1 and SI-2) without N fertilization, and the rainfed conditions with N combine the two levels of N fertilization (N60 and N120). Physiological parameters included in the model are: The stomatal conductance (g s ), the stable carbon isotope composition (δ 13 C) of mature grains, the Relative Green Area (GA) and the Relative Greener Area (GGA), and the canopy temperature depression (CTD). The width of the arrows is proportional to the path coefficient values. Dashed lines indicate negative relationships. CFIs with values > 0.9 are taken as indicative of a good fit. Significance levels-* p < 0.05; ** p < 0.01 and *** p < 0.001.
Path analysis was also examined within each of the nine growing conditions (resulting from the combination of the three water regimes and the three nitrogen fertilization levels ( Figure S4). Water irrigation conditions with and without N fertilization indicated the dependence of grain yield on the water status traits δ 13 C and CTD (and the latter also under rainfed conditions with N fertilizer). In the absence of N fertilization under both irrigation and rainfed conditions, GY seemed strongly and negatively dependent on the vegetation index GGA.

Discussion
The grain yields achieved, even under the best growing conditions (SI-2 with N120) are below 3 tonnes per hectare, and this clearly corresponds with moderate-low yielding conditions in the Mediterranean [44,54]. The vegetation indices tested in this study, generated either spectroradiometrically (NDVI) or derived from RGB images, performed well when assessing differences in water conditions. The efficacy of these indices in capturing differences in growth and senescence in response to water regime has been reported for wheat already [32,33,54]. The choice of anthesis as the phenological stage for the remote sensing measurements was decided based on the results of Yousfi et al. (2016) [54] under similar agro-ecological conditions as in the present study, together with an additional study under Mediterranean conditions where NDVI and RGB indices were measured periodically during the crop cycle [22], with robust correlations between those indices measured at anthesis and grain yield being reported. While the argument may be valid that any action on N management at anthesis is probably too late to significantly affect yield, particularly for fully irrigated wheat, it may still positively affect grain quality. More importantly, scheduling irrigation at anthesis may be fully relevant for cereals under Mediterranean conditions where drought increases progressively during the reproductive stage of the crop. Thus, following previous reports in wheat [54][55][56], a significant association between the NDVI, GA and GGA vegetation indices and the water status parameters (δ 13 C, g s and CTD) was also found in this study. The NDVI, GA and GGA were positively associated with CTD and g s and negatively associated with δ 13 C. In this context, Lopes and Reynolds (2012) [57] reported that the relationship observed between chlorophyll retention or 'stay-green' (assessed via the NDVI) and canopy temperature would confirm the functionality of stay-green in terms of gas exchange and would explain a better capacity to use water by the stay-green genotypes under stressful environments related to low fertilization and the lack of water. Our data confirm the close association between vegetation indices and water status parameters and identified the canopy greenness as a good indicator of crop water status and irrigation management.

Vegetation Indices and Nitrogen Fertilization
Digital images have been used to evaluate the nitrogen status of crops [58,59]. Our results showed that GGA was the only vegetation index exhibiting a significant difference between N treatments, with lower values found under N0 compared with N60 and N120. The absence of N fertilizer limits plant growth, whereas it may accelerate plant senescence during the reproductive stage of the crop, therefore decreasing GGA compared to plants fertilized with nitrogen. In this context, [31] described the greener area index (GGA) derived from RGB images is a good parameter for capturing active photosynthetic area and plant senescence because it is formulated with green pixels alone. Furthermore, the NDVI failed to assess differences under different N treatments. Digital pictures provide information that is not currently acquired through spectral reflectance measurements, such as the portion of yellow leaves in wheat growing under field conditions [31,32,52]. In the case of GA, this index, which takes into account yellow/green pixels, is less stringent in terms of excluding non-senescent parts of the plant. This may explain why GA measured during anthesis was not affected by N fertilization.

Canopy Temperature and Water Status in Wheat
Many studies have recognized canopy temperature depression (CTD) as an indicator of overall plant water status [36,60] and a potential tool for irrigation management [9,61]. In our study, CTD measured with an infrared thermometer was lower (and even negative) under rainfed conditions compared to support irrigation. A priori, a higher CTD indicates a greater capacity for transpiration, for taking up water from the soil, and therefore for maintaining a better plant water status [36]. In the case of the rainfed trials, the fact that the leaf temperature was higher than the air temperature indicates that rainfed plants were subjected to severe water stress that closed the stomata. In fact, the stomatal conductance measured in the rainfed plants was very low and one third of that measured under the best support-irrigation regime.
In addition, Gutierrez et al. (2010) [62] reported that the association between canopy temperature and the normalized difference water index confirmed that canopy temperature is a good indicator of hydration status. According to this, our results showed highly significant associations of CTD with δ 13 C (negative) and g s (positive). Furthermore, CTD seems to be a better indicator of the water status at the crop level than other traits related to water status, such as leaf g s [63]. In our study, CTD was strongly associated with δ 13 C (r = 0.84 ***) under SI-1 at N60, with the NDVI (r = 0.67 **) under rainfed conditions without N, and with grain yield (r = 0.66 **) under rainfed conditions at N120, while g s was not correlated with any of these parameters. These results confirm the close association between canopy temperature and other water status parameters and identify the canopy temperature as a good indicator of crop water status.

Relationship of Vegetation Indices and Water Status Traits with GY
As found in previous studies in wheat, the RGB canopy indices measured at flowering were strongly correlated with GY [31,32,64]. For wheat under Mediterranean conditions, the reproductive stage is usually the best period for crop monitoring, since the crop is exposed to increasing stress (drought) conditions during the last part of the crop cycle. Following on from this, the present study revealed a positive relationship between GA and GGA with grain yield under irrigation. Additionally, stepwise analysis reinforced the evidence for the usefulness of the RGB vegetation indices to assess GY. The GA vegetation index was chosen by the model as the first independent variable, explaining 66% of GY variability under SI-2 without N fertilizer. Moreover, various studies have reported that RGB-based indices may perform far better than the NDVI for GY prediction in wheat [32,56,65]. In our study, the NDVI failed to assess GY under irrigation. In contrast, the NDVI was correlated (positively) with GY under rainfed conditions and was also the only variable chosen by the stepwise model, explaining 63% of GY variation under rainfed conditions. In this context, Casadesus et al. 2007 [18] reported that the NDVI measured at anthesis in durum wheat correlated positively with GY under severe water stress conditions, but failed to correlate under well-watered conditions. Verhulst and Govaerts (2010) [66] have also reported that the NDVI has been correlated with long-term water stress. The reason for the low correlation of NDVI against GY under well-watered conditions is because plant canopies during anthesis are very dense and the measured NDVI values become saturated. NDVI is an index based on the strong contrast between the near infrared and the red band reflectance of a vegetation canopy, and this difference becomes wider as the canopy cover increases. Thus, NDVI works better with stress conditions where canopies are sparse and/or early senescence is present [27,31]. In any case, anthesis proved to be the correct phenological stage for remote sensing evaluations when crops under different levels of stress were compared, which may be the case for crops exposed to a range of different combinations of water and nitrogen fertilization conditions. Furthermore, the water status of plants can also be associated with grain yield. In our study, CTD (positively) and δ 13 C (negatively) correlated with GY and both parameters were chosen by the stepwise model as the first variables to explain GY variation under different irrigated and rainfed growing conditions. In this context, previous studies have shown that a higher CTD is associated with increased wheat yield under irrigated, hot environments [38,67], but also under dryland environments [68].
However, in our study, and regardless of the water regime (rainfed or irrigation), neither the vegetation indices (NDVI, GA and GGA) nor the water status indices (CTD and δ 13 C) were associated with GY in the absence of N fertilizer. A lack of variability in green biomass and grain yield associated with the lack of nitrogen fertilization might explain this outcome.

Phenotyping Parameters under Different Water and N Supplies
The results of our study have shown the usefulness of vegetation indices with low implementation costs as a means to identify genetic variability under different growing conditions in the field. From nine of the growing conditions studied (resulting from the combination of the three water regimes and the three nitrogen fertilization levels), GY was significantly different between genotypes under only two of the growing conditions (SI-2 at N120 and N60). In contrast, two of the vegetation indices (NDVI, GGA) measured at anthesis were able to distinguish between the genotypes growing under six of the nine growing conditions (including rainfed), regardless of the N fertilization conditions. In this context, multispectral ground-based portable spectroradiometric devices have been used in wheat phenotyping [24,27]. The conventional RGB images have also been proposed as a selection tool for cereal breeding [18][19][20].
The genotypic differences observed using vegetation indices possibly reflect differences in canopy stay green during the reproductive stage. Lopes and Reynolds (2012) [57] reported that stay-green is regarded as a key indicator of stress adaptation. Thus, our study revealed the usefulness of the vegetation indices to select the most tolerant genotypes in terms of retaining a greener biomass during the last part of the crop cycle. The three vegetation indices assayed were able to identify genotypic differences, even under the most severe growing conditions, such as rainfed with and without N fertilizer and where GY failed to detect differences among genotypes. Phenotyping wheat genotypes for water and N fertilization deficit at anthesis using these vegetation indices should permit the formulation of the best crosses between genotypes. However, by comparison, canopy temperature performed much worse as a phenotyping parameter in our study. It has been reported that CTD is a poor indicator of plant performance when the yield is highly dependent on limited amounts of soil-stored water [69,70]. Moreover, the canopy in these trials, particularly during the reproductive stages, frequently leaves areas of bare soil exposed that may affect the canopy temperature readings. Leaf chlorophyll content measured by a portable device was perfect for distinguishing among genotypes, regardless of the water status (irrigation or rainfed), but only when trials were provided with nitrogen fertilizer. Therefore, for the agronomic conditions of our study, the vegetation indices assessed at the canopy level performed better as phenotyping tools than canopy temperature and chlorophyll content measures.

Grain Yield Prediction across Crop Seasons Using Low-Cost Remote Sensing Techniques
The results of data combining the growing conditions, genotypes and replicates support the use of different affordable remote sensing techniques to estimate grain yield across crop seasons. However, in agreement with Clevers (1997) [71], estimates of crop growth and yield using crop growth models often lost accuracy as the growing conditions became more stressed. The loss of accuracy may be the consequence of a very narrow range of variability in grain yield associated with stressed growing conditions. Moreover, vegetation indices derived from RGB images performed comparatively better than the NDVI, probably because GA was less saturated than the NDVI. The application timing could have played a critical role here-i.e., saturated NDVI at anthesis is indeed not expected to perform well, while the saturation pattern of RGB indices is less evident. In fact, the acquisition of high-resolution RGB images is fast and its dependence on atmospheric conditions (e.g., sunny versus cloudy days) is minimal [22,32]. Therefore, the availability, cost and practicality of digital cameras make them an ideal tool for the management of crop water and fertilization status [19,22]. However, in agreement with previous reports, the relationships between yield and vegetation were site and season specific [20,21]. In our study models were not able to predict absolute yields, but rather relative differences in yield, which makes the approach unfeasible for yield forecasting, but it is still useful in terms of crop management and even phenotyping. The strong relationships of these vegetation and water status indices with grain yield expressed in relative units support the effectiveness of these low-cost indices in crop management. Nevertheless, whereas the evaluations in the three successive seasons were performed in the same region, crop management conditions (water and fertilization regimes) affected the performance of the models. Hence, to make yield predictions more holistic and effective across different environments, it is necessary to use more robust calibration; for example, incorporating site-year covariates or a multi sensor approach [21].

An Integrated Model to Predict Grain Yield That Combines Remote-Sensing Traits, Canopy Reflectance Measurements and Grain δ 13 C
In this study, we performed a path analysis to dissect how the vegetation indices (GA and GGA) and the water status (CTD, δ 13 C and g s ) indices directly or indirectly assessed GY performance within each of the growing conditions assayed. The following parameters may provide suitable coverage of the factors affecting GY performance under a given water regime and nitrogen fertilization supply: Vegetation indices as indicators of photosynthetic capacity (GA) and the effect of early senescence (GGA) on the canopy; CTD which informs about the current water status of the canopy; the δ 13 C in mature grains as a time-integrated indicator of photosynthetic and transpirative gas exchange of the crop; and the water status at the single organ level (assessed as the g s of the flag leaf). Under irrigation and nitrogen fertilization, both δ 13 C and CTD (indicators of photosynthetic and water status) had a direct association with GY. Better grain yield performance is associated with higher CTD and lower δ 13 C under supplementary irrigation. The association of δ 13 C (negative) and CTD (positive) with GY is probably due the higher stomatal conductance and transpiration, (which increases CTD), therefore increasing the photosynthetic capacity even at the expense of a lower water use efficiency (and thus δ 13 C) and the consequent increases in GY [47,72]. Under rainfed conditions with N fertilizer, the transpiration/water status (assessed through CTD) and photosynthetic potential (evaluated through vegetation indices) affected GY. Under rainfed conditions with N120, CTD had a positive effect on GY. We suggest that N fertilization only had a positive effect on grain yield providing that there was water available to maintain transpiration (higher CTD), stomatal conductance and thus photosynthesis in the available canopy (higher vegetation indices).
In the absence of nitrogen fertilization, and despite the water conditions (irrigated or rainfed), the total canopy area (evaluated through GA) has a positive effect on GY. However, an excess of young (not senescing) leaf area (evaluated through GAA) had a strong negative association with GY. Stay-green character may have a negative effect on GY in the absence of nitrogen fertilization because it limits the retranslocation of N to the inflorescences and ultimately affects grain filling. In the case of the rainfed crop fertilized with a limited amount of nitrogen (N60), the active canopy area (evaluated through GGA) had a positive effect on grain yield, which may indicate that the limitation is imposed by the amount of photosynthetic area rather than by the availability of N to reproductive tissues. Under conditions of high nitrogen fertilization (N120) and irrespective of the water regime (irrigated or rainfed), GY is not affected by the size of the canopy or even by its greenness (assessed through GA or GGA), but by the water status of the crop (evaluated through CTD and δ 13 C).

Conclusions
This study demonstrated the potential of low-cost RGB vegetation indices and the canopy temperature for the management of growing conditions (essentially the water and nitrogen regimes) under Mediterranean conditions. Although the models did not predict absolute yields, they are still useful in terms of crop management and even phenotyping. Nevertheless, grain yield estimation performs better under irrigation than under the low-yielding conditions of rainfed cultivation in the absence of nitrogen fertilization and this illustrates one of the potential limitations associated with the remote sensing-based yield predictions; they are affected by specific environmental conditions. Even though a multispectral vegetation index, such as the NDVI is a widely accepted approach to monitor changes in growth under different conditions, in this study we have shown that vegetation indices derived from conventional images like the GA and GAA indices provided a similar if not better prediction of grain yield and at a comparatively lower cost than the NDVI. The use of vegetation indices derived from RGB images to assess GY could be popularized in the near future via apps installed on mobile phones.
The vegetation indices have also proven their suitability for differentiating among genotypes. Furthermore, these traits may contribute through path analysis to develop physiological models for assessing wheat ideotypes best suited to different water and nitrogen regimes. The models showed that for nitrogen fertilized trials, and regardless of the water regime imposed, the water status parameters were the main factors determining GY performance. Moreover, a larger green area at anthesis may also contribute to a larger yield. In the absence of nitrogen fertilization, a large greener canopy area (assessed through GGA) at anthesis is a factor that negatively affects grain yield.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-4395/9/6/285/s1, Figure S1: Scheme detailing the different plots of the experimental design, Figure S2: Conceptual model of the path analyses quantifying the relative strengths of the direct and indirect relationships of the different physiological traits and grain yield. Physiological parameters included in the model are: The stomatal conductance (g s ), the stable carbon isotope composition (δ 13 C) of mature grains, the Relative Green Area (GA) and the Relative Greener Area (GGA) indices calculated from digital pictures, and the canopy temperature depression (CTD) measured with an infrared thermometer, Figure S3: Relationships of the measured versus the predicted grain yields of wheat achieved during two successive crop seasons (2013-2014 and 2014-2015). Predicted grain yield values were calculated using the linear relationships of the grain yield with the two vegetation indices, the Normalized Difference Vegetation Index (NDVI) measured with a portable spectroradiometer, and the relative Green Area index (GA) calculated from digital images, and the canopy temperature depression (CTD) measured with an infrared thermometer. All variables were evaluated during the 2012-2013 crop season. Measurements were performed in the same region as the present study. Grain yield and vegetation index data of the two first crop seasons have been reported in Yousfi et al. (2016). For each crop season, the nine different growing conditions (resulting from the combination of the three water regimes and the three nitrogen fertilization levels) were analysed. The different combinations of nitrogen fertilization and water regimes for the 2013-2014 crop season were as follows: Supplementary irrigation (SI-2) with nitrogen fertilization N120 and N60 and without N fertilization; Supplementary irrigation (SI-1) with nitrogen fertilization N120 and N60 and without N fertilization; Rainfed with nitrogen fertilization N120 and N60 and without N fertilization. In addition, for the 2014-2015 crop season the combinations were as follows: Supplementary irrigation (SI-2) with nitrogen fertilization N120 and N60 and without N fertilization; Supplementary irrigation (SI-1) with nitrogen fertilization N120 and N60 and without N fertilization; Rainfed with nitrogen fertilization N120 and N60 and without N fertilization. Significance levels-ns, not significant; * p < 0.05; ** p < 0.01 and *** p < 0.001, Figure S4: Path analyses of four wheat genotypes grown under different combinations of nitrogen fertilization and water regimes. The different combinations of nitrogen fertilization and water regimes are as follows: (A) SI-2 with high fertilization (N120); (B) SI-2 with medium nitrogen fertilization (N60); (C) SI-2 without nitrogen fertilization; (D) SI-I with high nitrogen fertilization (N120); (E) SI-I without nitrogen fertilization (N60); (F) Rainfed with high fertilization (N120); (G) Rainfed with medium nitrogen fertilization (N60); (H) Rainfed without nitrogen fertilization. Physiological parameters included in the model are: The stomatal conductance (g s ), the stable carbon isotope composition (δ 13 C) of mature grains, the Relative Green Area (GA) and the Relative Greener Area (GGA) indices calculated from digital pictures, and the canopy temperature depression (CTD) measured with an infrared thermometer. The width of the arrows is proportional to the path coefficient values. Dashed lines indicate negative relationships. Overall fit statistics for each path model (chi-squared, the probability and comparative fit index, CFI) are shown at the bottom left of each panel. CFIs with values > 0.9 were taken as indicative of a good fit. Significance levels-* p < 0.05; ** p < 0.01 and *** p < 0.001; Table S1: Soil chemical characteristics at different depths, Table S2: Monthly total accumulated rainfall (PP), minimum air temperature (T min), maximum air temperature (T max) and average air temperature (T aver) for the 2014-2015 crop season. Values were collected at the meteorological station of Khemis Miliana (Algeria), Table S3: Correlation coefficients of the linear relationships of grain yield (GY) with NDVI, GA, GGA, CTD and δ 13 C under irrigation without N fertilization and rainfed conditions with N120, N60 and without N, Table S4. Multiple linear regressions (stepwise) across genotypes and replicates explaining grain yield (GY) variation as a dependent variable and the NDVI, GA, GGA, LC, CTD, g s and δ 13 C as independent variables. Data were analyzed within each level of nitrogen fertilization and water regime. Significance levels-* p < 0.05; ** p < 0.01 and *** p < 0.001. Abbreviations for variables and growing conditions as defined in Tables 1 and 2