High-resolution spectral information enables phenotyping of leaf epicuticular wax in wheat

Abstract

Background

Epicuticular wax (EW) is the first line of defense in plants for protection against biotic and abiotic factors in the environment. In wheat, EW is associated with resilience to heat and drought stress, however, the current limitations on phenotyping EW restrict the integration of this secondary trait into wheat breeding pipelines. In this study we evaluated the use of light reflectance as a proxy for EW load and developed an efficient indirect method for the selection of genotypes with high EW density.

Results

Cuticular waxes affect the light that is reflected, absorbed and transmitted by plants. The narrow spectral regions statistically associated with EW overlap with bands linked to photosynthetic radiation (500 nm), carotenoid absorbance (400 nm) and water content (~ 900 nm) in plants. The narrow spectral indices developed predicted 65% (EWI-13) and 44% (EWI-1) of the variation in this trait utilizing single-leaf reflectance. However, the normalized difference indices EWI-4 and EWI-9 improved the phenotyping efficiency with canopy reflectance across all field experimental trials. Indirect selection for EW with EWI-4 and EWI-9 led to a selection efficiency of 70% compared to phenotyping with the chemical method. The regression model EWM-7 integrated eight narrow wavelengths and accurately predicted 71% of the variation in the EW load (mg·dm⁻²) with leaf reflectance, but under field conditions, a single-wavelength model consistently estimated EW with an average RMSE of 1.24 mg·dm⁻² utilizing ground and aerial canopy reflectance.

Conclusions

Overall, the indices EWI-1, EWI-13 and the model EWM-7 are reliable tools for indirect selection for EW based on leaf reflectance, and the indices EWI-4, EWI-9 and the model EWM-1 are reliable for selection based on canopy reflectance. However, further research is needed to define how the background effects and geometry of the canopy impact the accuracy of these phenotyping methods.

Background

Wheat is a major staple food and an important source of calories in developing countries [1]. More than 220 million ha of wheat is cultivated worldwide [2], and 600 million tons of wheat grain is produced each year [3]. The expected global population of 9 billion by 2050 will require an increase in wheat production of 60% to 100% [2, 4], but the current genetic gains of < 1% per year [5] will be insufficient to fulfill this expected demand. Annually, more than 600 million tons of wheat are harvested [6], but maintaining this production is already a challenge in the face of climate change. It is estimated that climate change can reduce global wheat production by 6% for every degree centigrade increase in the temperature [7]. Therefore, the development of wheat cultivars that are adapted to high temperatures and limited irrigation is crucial for responding to a changing climate and ensuring wheat production.

Developing wheat cultivars that are adapted to a wide range of environments, have resilience to abiotic stresses and high yield potential are priorities of the main public breeding programs [2]. Physiological trait-based improvements for tolerance to heat and drought stress rely on the favorable expression of morphological and physiological plant traits (PT) [8,9,10]. Independent conceptual models for grain yield (GY) under heat and drought have been proposed based on the following main drivers: light interception (LI), radiation use efficiency (RUE), partitioning of total assimilates [8], water use efficiency (WUE) and harvest index [11]. Each of these main drivers contains genetically determined PTs that can potentially lead to an additive genetic effect for resilience to heat and drought when combined through strategic crossing [12, 13]. Physiological traits such as canopy temperature (CT) are already utilized as selection criteria in breeding pipelines [5, 10], but key PTs such as epicuticular wax (EW) remain unexplored because of the expensive, subjective and laborious method for phenotyping [14].

EW is the outermost layer of leaves that is located on the top of the cutin matrix and intracuticular wax [15] and consist of hydrocarbon compounds [16, 17] derived from long chains of C₂₀ and C₃₀ fatty acids [18, 19]. The hydrophobic layer that creates the EW covers the aerial epidermis of plants maintaining the integrity of the plant against high UV radiation [20] and external environmental stresses such as insect infestation [21, 22] and pathogen infection [23]. This cuticle also minimizes the water loss via transpiration in wheat [18, 24] and reduces leaf temperature [25, 26]. Early studies estimated a decrease in the internal temperature of the plant by 0.7° C under simulated drought conditions in a glasshouse, saving 30 g of water per plant during the growth season and extending grain filling by 3 days [27].

Waxes and trichomes affect the interaction of the plant with the environment, particularly the reflection and absorbance of light. Surface waxes are very effective in reflecting excessive radiation in specific ranges of the spectrum, namely at 330 and 680 nm [28], but the main increases in radiation reflection occur at the photosynthetic active radiation (PAR) range to dissipate excess energy and avoid damage to the PSII reaction center [29, 30]. In wheat, increases in light reflectance of 12% to 35% were detected at the PAR (400 to 700 nm) range in wax covered genotypes [31]. Several studies have reported that EW and its constituents are an important protective barrier against harmful UV-B radiation [20, 28, 32,33,34,35], but these fluctuations in reflectance are species-specific and can range from < 10% in most species to 70% in others [21, 36].

Limitations on field phenotyping restrict our capacity to unravel complex morphological and physiological traits. Spectral technologies have the potential to increase the efficiency, precision and accuracy of phenotyping platforms. In breeding programs, high-precision phenotyping can enable the screening of segregating material, advanced lines and germplasm [5, 37]. Increasing the accuracy of phenotyping can provide more reliable estimates of heritability and variance components [38], facilitate gene discovery and enable prediction of complex traits with approaches such as genomic selection [39, 40]. The strong association of spectral secondary traits with GY [41, 42] highlights the potential of canopy reflectance to increase productivity in wheat. A detailed list of sensors and their applications for plant phenotyping is provided by [43].

Recent advances in technology have maximized the throughput of phenotyping measurements under field and controlled conditions [44,45,46]. RGB and hyperspectral sensors have enabled the rapid and noninvasive acquisition of spectral information. Spectral vegetation indices (SVI) are a quick, easy and inexpensive method of transforming light reflectance into simple indicators of photosynthetic and canopy variations. The simple ratio index (SR) [47] and the normalized difference vegetation index (NDVI) [48] are two of the first SVIs developed for detecting green vegetation. Both indices combine the percentage of reflectance at the wavelengths where plants absorb (~ 750 to 800 nm) and reflect (800 to 2500 nm) light. Several other SVIs have been built for sensing the water content of plants [49], photosynthetic radiation [50], carotenoid pigments [51], plant height [52], leaf area [53], and diseases [54].

In this study, the aim was to develop spectral methods to indirectly phenotype EW accumulated on the surface of leaves. This wax index will serve as a proxy to detect genotypes with a thick wax cover, in order to integrate the trait into breeding pipelines to enhance resilience to heat and drought stress. The goal is to facilitate frequent screening for EW at multiple field trial locations by replacing conventional sample-based methods. Additionally, these methods will support ongoing research on the underlaying physiological and genetic mechanisms of cuticular waxes. We conducted a set of theoretical studies with the following specific objectives: i) detect the wavelengths at which reflectance is affected by cuticular waxes, ii) develop spectral indices and models to detect wheat cultivars with high and low EW content, and iii) validate the spectral methods for phenotyping under field conditions.

Results

Light interactions associated with leaf EW

The differences in the light interactions detected after the removal of EW confirmed the role of the cuticle in avoiding and dissipating excess radiation (Fig. 1). Variations in the percentage of light absorbed (Fig. 1a), transmitted (Fig. 1b) and reflected (Fig. 1c) by leaves were detected when EW was partially eliminated. The removal of the cuticle increased the light absorbance in the visible range from 0.02 to 0.04%, with a subsequent decrease to 0% reflectance at 710 nm and 0.03% in the NIR. An increase in light transmission through the leaf from 0.01 to 0.06% in the visible region was also observed, with a significant increment of 0.13% in the red-edge (680 to 740 nm). In the NIR (740 to 980 nm), the transmittance also increased by approximately 0.06%. Light reflectance was most affected when the wax cuticle was removed. Its removal revealed that EW contributed from 0.05% to 0.15% of the increase in reflectance at various wavelengths. Further analysis enabled the estimation of both positive and negative variations in the percentage of light reflected by the unit (mg·dm²) of wax accumulated on top of the leaf surface. Figure 1d presents the slopes of the linear regression models individually fitted with the EW content as the independent variable and the percentage of light reflectance detected with the spectroradiometer as the dependent variable. From 424 to 450 nm, there was an increase of ~ 0.82% in reflectance, and from 544 to 575 nm the increase was 0.79%. The light reflectance in the 700 and 730 nm was not affected by the cuticular wax; however, there was a reduction of 0.77% from 713 to 720 nm, with the highest peak in 717 nm (− 0.8%), and a subsequent increase of ~ 1.5% from 756 to 825 nm.

Fig. 1

Variation in a absorbance, b transmittance and c reflectance derived by the removal of the EW coat with HPLC chloroform (CHCl3). These variations are presented as the difference of the spectral signature of the leaf after the wax coat was removed minus the the spectral response of the leaf with the wax coat in place. The red solid lines define the red edge and the dash line mark the end of the visible and start of the NIR region. The slope of the linear regression models (d) were fitted as Y = a + bX, where Y corresponds to the independent variable EW (mg·dm⁻²), X is the percentage of light reflectance at one nanometer resolution, a and b are the intercept and the slope of the fitted model, respectively. The statistical significance of the models was $P\le 0.05$ or less

The partial least square regression (PLSR) analysis reduced the total number of spectral bands and integrated only uncorrelated bands in the predictive model (Fig. 2). The correlation coefficients of the regression between the wavelengths and EW content are presented in Fig. 2a. Three main components enabled the maximum correlation between the wavelengths and the EW content and explained 97.34% of the variability of the trait. These three components combine the follow spectral regions: 424 to 448, 625, 660, 712 to 727, 775 to 835, 994 and 997 nm. Most of these wavelengths coincide with the regions detected in Fig. 1d. The most influential variables were 712 to 727 nm, where the transition from low reflectance in the visible region to high reflectance in the NIR wavelengths occurs. Overall, the selected latent variables or wavelengths enabled the prediction of EW content in the data subset for validation and lead to a positive linear association between the predicted and the actual values of EW (Fig. 2b).

Fig. 2

Correlation coefficients a of the three main partial least square components with the EW content (mg·dm²) and b association of the EW load predicted with the PLS’s model vs EW measured by the chemical method

Spectral indices for indirect phenotyping of EW

Spectral indices developed in previous studies for indirect phenotyping of morphological and physiological characteristics of the plant (Tables 1 and 2) were calculated and statistically associated with EW. The spectral resolution of the light reflectance extracted from the spectroradiometer was adjusted from 1 to 3 nm by averaging the percentage of light reflectance of every 3 bands. The broad spectral bands were calculated by estimating the average of reflectance w ithin the range (nm) corresponding to the blue, green, yellow, orange, red, red-edge and near infrared regions. The ranges in nm of every spectral region are included in Table 2. The narrow spectral indices PRI-1 (r = −0.57), CARI-2 (r = − 0.67), PSSR-b (r = − 0.5 7), PSSR-a (r = − 0.55 ), ARI-1 (r = − 0.52) , ARI-2 (r = − 0.58) and SIPI-2 (− 0.61) were significantly correlated (p < 0.001) with EW. These indices are effective to detect variations in carotenoids, chlorophyll and anthocyanins in plants [51, 55,56,57,58]. However, these indices were not able to predict more than 38% of the total variability of EW. Among the broad vegetation indices calculated, only RGRI (− 0.57) and ARI-1 (− 0.67) were strongly associated with EW.

Table 1 Narrow vegetation indices for phenotyping specific traits in plants

Table 2 Broad band vegetation indices for phenotyping specific traits in plants

The broad and narrow indices developed in this study are presented in Table 3. The selection of these indices was based on their high R² in the cross validation (LOOCV) and low root mean square error (RMSE) estimates in the bootstrapping analysis. EWI-13 and EWI-14 estimated 65% and 62% of the EW variation combining the wavelengths 625, 736 and 832 nm. The indices EWI-6, EWI-9 and EWI-12 integrated only two wavelengths, 658 and 712 nm; 670 and 718 nm; and 622 and 718 nm, respectively. The proportion of the variance in the EW explained by the EWI (R²) was as follow: EWI-6 = 0.52, EWI-9 = 0.51 and EWI-12 = 0.51. When the broad spectral bands blue, red and NIR were combined in the spectral indices, the prediction accuracy ranged from 31 to 44%. Specifically, EWI-1 estimated 44% of the variability with a RMSE of 1.19 mg·dm⁻². The slope of the linear models (B) in most cases was positive, except those for EWI-3, EWI-4, EWI-8, EWI-10 and EWI-11. The increase in EW content of 1 mg·dm⁻² caused wide variations in the values of the broad and narrow indices from 0.002 to 5.73.

Table 3 Coefficients of determination (R²) and root mean square error (RMSE in mg·dm⁻²) of the indices developed for phenotyping EW in leaves

Prediction accuracy of spectral indices for phenotyping of EW load

The EW content determined with the chemical method from samples collected in the field experimental trials of the mapping population ranged from 1.54 to 2.4 mg·dm⁻² (Table 4). The heritability estimate (h²) of EW ranged from 0.51 to 0.58 across all three trials. Overall, the CV of all trials was low, 6.9 in DS-1, 5.8 in DS-2 and 7.6 in DS-3. All fourteen indices included in Table 3 were calculated with ground spectral information. However, only the indices EWI-1, EWI-2, EWI-3, EWI-4, EWI-9 and EWI-13 were strongly associated with EW deposition when estimated with canopy reflectance and are the only indices included and discussed in Table 4.

Table 4 Mean, genetic variance (${{\varvec{\sigma}}}_{{\varvec{g}}}^{2}$), heritability estimate (${{\varvec{h}}}^{2}$), error variance (${{\varvec{\sigma}}}_{{\varvec{e}}}^{2}$) and coefficient of variation (CV in %) of EW content (mg·dm⁻²), EWI-1, EWI-2, EWI-3, EWI-4, EWI-9 and EWI-13

The phenotypic (r_p) and genotypic (r_g) correlations of the top performing indices estimated with the ground hyperspectral information and EW content are presented in Fig. 3. All correlations were statistically significant at P ≤ 0.01. As expected, genotypic correlations were higher than phenotypic correlations in all cases. According to these parameters, the index EWI-4 and EWI-9 better estimated the variation of EW. The average response of r_p and r_g were -0.51 and -0.55 for the EWI-4, and 0.33 and 0.48 for EWI-9.

Fig. 3

Phenotypic $({\sigma }_{p}$) and genotypic (${\sigma }_{g})$ correlation of the epicuticular wax indices (EWI) 1, 2, 3, 4, 9 and 13 with EW content measured with the chemical method (mg·dm⁻²). All six indices were statistically significant (P < 0.01) across all three sets (DS-1, 2 and 3)

Although the lack of variance of EW in the mapping population used in this study might limit the response for direct selection, the moderate h² of the trait would lead to genetic advances when selection is applied (Table 5). The genetic gain (GG) for EW with selection pressure of 10% is also included in Table 5. In DS-1, the GG was 0.65 mg·dm⁻²; in DS-2, it was 0.89 mg·dm⁻², and in DS-3 it was 0.59 mg·dm⁻². However, when the mean of the actual population was considered, the genetic advance with direct selection averaged 2.5%. The correlated response of the EW indices with the EW content derived an increase in EW. Improvement in the EWI-1, 2, 9 and 13 resulted in increases of 0.063, 0.053, 0.047 and 0.043 mg·dm⁻² of EW content, respectively. On the other hand, decreases of 0.047 and 0.063 mg·dm⁻² were calculated with a positive selection of the indices EWI-3 and 4. The efficiency of selection (RE) based on the secondary characters or indices (EWI) ranged from 46 to 78% in average. However, EWI-1, 2 and 4 in DS-2 were almost as efficient in selection as the direct selection of the trait with the chemical method with 112, 99 and 90% efficiency, respectively.

Table 5 Genetic gain (GG in mg·dm⁻²), genetic advance with respect to the mean (GAM in %) and response to direct selection (R) of EW

Multivariate regression models integrating narrow spectral bands to predict the EW load

The statistics of the multivariate models developed with the selected bands in the PLSR analysis are presented in Table 6. The final selection of the variables led to seven models in single and/or multiple combinations of eight wavelengths. The spectral response at 424 nm predicted almost 33% of the total variability of the trait in the validation set. However, when as many as seven spectral bands were incorporated in a model (424, 547, 574, 658, 712, 721, 775 and 817 nm), the accuracy increased by 38% (EWM-7 with R² = 0.71). The RMSE of the prediction was consistent across the models, ranging from 0.49 to 0.52 mg·dm².

Table 6 Statistics of regression models (EWM)

The RMSE of the seven EW models was estimated in the four experimental trials in which ground and aerial reflectance were collected. The square root of the residuals is presented in Fig. 4a. A considerable increase in the error of the prediction models calculated with ground and aerial hyperspectral information was observed across all trials in comparison to the RMSE estimated with reflectance from single leaves. The highest prediction accuracy was obtained with EWM-1, with an average RMSE of 1.4 mg·dm⁻² from the ground measurements and of 0.63 mg·dm⁻² from the aerial information. EWM-2 seems to accurately estimate EW load utilizing light reflectance in the same way as EWI-1. However, in cases as DS-1, the prediction accuracy with the EWM-2 led to a RMSE as high as 5.7 mg·dm⁻².

Fig. 4

Root Mean Square Error (RMSE in mg·dm⁻²) of the multivariate models (EWM) for predictions of epicuticular wax load utilizing the ground-based and aerial hyperspectral reflectance. The RMSE of prediction for the ground-based information is presented as the average response across the four sets of wheat inbreed lines evaluated. DS stands for data set 1 to 3 and Aerial corresponds to the hyperspectral information collected with the aircraft. (O-P)² is the square of the difference between the observed minus the predicted values

Discussion

In this study we evaluated the spectral response of leaves and derived and validated a set of indirect methods for phenotyping the trait. Furthermore, differences in light interaction derived by cuticular waxes and detected in this study coincide with results from studies conducted in Vitis vinifera [59], Leucadendron lanigerum [60] and Cotyledon orbiculata [61]. The increase of approximately 10% of light reflectance in the NIR region was considerably less than changes in reflectance previously detected on wheat ~ 15% [62] and in rosette succulent plants ~ 50% [31]. Additionally, the violet (r = 0.64) and blue (r = 0.63) spectral regions strongly correlated with EW, and it is in line with preliminary studies where waxes were reported as photoprotective mechanisms against short wavelength radiation [17, 63,64,65]. A significant number of research studies also reported that waxes enhance UV (~ 100 to 400 nm) reflectance [17, 28, 63, 66], but the analysis of these wavelengths is outside of the scope of this study due to limitations of the equipment utilized to collect the spectral information.

The wax cuticle enhances light reflectance almost by 1% per every unit of wax (mg·dm⁻²) accumulated on top of the leaf surface, but specifically in the PAR wavelengths where the absorption of photosynthetic pigments occurs [17]. In sorghum (Sorghum bicolor L.), a similar increase of 3% in reflectance by the cuticular leaf coat was reported, but the result was based on wavelengths of the spectrum from 400 to 1000 nm, without a detailed examination of specific narrow spectral bands [67]. Among the spectral regions associated with EW, the wavelength at 424 and 448 nm are linked to the absorption of light by carotenoids in plants [56], while the surrounding wavelengths at approximately 500 nm are associated with the dissipation of excess radiation and the efficiency of photosynthetic radiation [50, 68, 69]. Additional peaks of absorption of chlorophyll a and b in the 600 nm coincides with two main wavelengths linked to EW in 625 and 660 nm.

Several narrow and broad indices for phenotyping additional traits in plants were correlated with EW load, but the moderate to low correlation (Table 7) of these indices with EW limits any form of application for phenotyping. Although a reasonable r² value of 0.65 was estimated when three narrow spectral bands were integrated in the novel spectral indices (EWI-13 and EWI-14 in Table 2), the high cost of sensors required to acquire hyperspectral reflectance can limit the utilization of these indices. On the contrary, EWI-1 requires two main spectral ranges (blue and red) that can easily be extracted from RGB images. Broad-sense heritability was estimated for the EW indices and the EW measured by the chemical method (Table 4). In all three trials h² was considerably higher than in prelaminar published results [70]. The six EW indices enabled a more reliable and precise quantification of the proportion of the genetic variance of the trait by considerably decreasing the error variance. However, the coefficients of variation of the indices EWI-3, EWI-4, EWI-9 and EWI-13 in DS-1 were considerably large due to the dispersion of the data around the mean of the population.

Table 7 Correlation of the narrow and broad vegetation indices and the EW content (mg·dm⁻²)

The moderated h² of EW let to genetic gains of up to 3%, a reasonable advancement for a quantitative trait and superior to genetic gains of ~ 1% in grain yield [71, 72]. All four broad and two narrow indices presented in Table 5 positively improve EW content, except EWI-3 and EWI-4, for which negative selection is needed to increase the EW load on leaf surfaces. The efficiency of indirect selection with the spectral indices was highly dependent on the experimental trial and its coefficient of variation, as it is the case in the DS-2 where selection with the EWI-1 was 12% more efficient that the direct selection. Examining the residuals of the model against the EW measurements, we observed a shift towards an increase in the residuals of EWI-1 as the EW content increases, suggesting a potential restriction on utilizing this index for phenotyping genotypes where the EW is above 6 mg·dm⁻². However, this was not observed with the residuals of EWM-2. We suspect that implementing the EWI developed in this study with aerial spectral reflectance might lead to a low-quality phenotyping of EW and could potentially lead to confounding results.

Conclusions

EW is the outermost cuticle of leaves and directly affects light interactions, especially reflectance. This cuticle increases light reflectance at the visible and NIR regions by 0.5% and 1.6%, respectively. Integrating specific narrow wavelengths that are highly sensitive to variations in the EW load, we generated several spectral linear models and vegetation indices for predicting the EW content and detecting cultivars with low and high EW. The prediction accuracy of these phenotyping methods was dependent on the characteristics of the sensor utilized to capture the spectral information, as well as on the canopy architecture and the distance of the sensor from the ground. With light reflectance captured from either the adaxial or abaxial side of the leaf, the broad index EWI-1 and the narrow index EWI-13 can accurately estimate EW. However, for canopy reflectance, the indices EWI-4 and EWI-9 more accurately estimate the density of the cuticle and led to a similar genetic advance than that from direct selection for the trait. The multivariate regression model EWM-7 integrated eight wavelengths distributed across the visible and NIR spectra and estimated 71% of the variation of the trait from the reflectance of a single leaf. In contrast, with ground and aerial reflectance, EWM-1 and EWI-2 accurately estimated the EW content (mg·dm⁻²), but insensitivity to variation at EW values larger than 6 mg·dm⁻² was detected for EWI-1.

Methods

Plant material and culture

The first set of genotypes evaluated were twenty-four recombinant inbred lines (RILs) derived from a cross of the spring cultivars Halberd (tolerant to heat stress) and Len (susceptible to heat stress). The lines were grown in a completely randomized design (CRD) with four replications in a growth chamber programmed with intervals of twelve hours of light and dark. Plants were sown in nursery pots 0.185 m in height with a diameter of 0.162 m that were filled with peat moss. The plants were fertilized twice during the growing season with the standard fertilizer 20–20-20 (N-P₂O₅-K₂O).

Leaf radiometric measurements

The spectral response from 350 to 1050 nm was captured with a CI-710 miniature leaf spectrometer from CID Bio-Science in 3022 spectral channels. The equipment was calibrated every five minutes with an integrated BaSO4 white reference disk for 100% reflectance and a black panel for 0% reflectance. Ten readings of the light reflected by the flag leaf were obtained prior to collecting the leaf sample for wax quantification. The spectral range of the signatures were adjusted to 400–900 nm and the spectral resolution to 3 nm by averaging the spectral channels included every 3 nm. The last step was to estimate the average of the ten spectral readings, only the averaged signature was considered for further analysis.

Quantification of epicuticular wax

Leaf samples were collected after the light reflectance at the adaxial and abaxial sides of the leaf was recorded, approximately 10 days after pollination (DAP). Each sample consisted of six leaf punches of 0.01 m diameter and were collected in 2.0 ml (ml) glass vials. The EW was extracted by immersing the leaf punches in 1.5 ml of HPLC chloroform (CHCl₃) for 20 s, and the EW was quantified via the colorimetric method described by [14]. The optical density of every sample at 590 nm was measured with PHERAstar® spectrophotometer. A standard curve was developed to transform the readings of absorbance to milligrams (mg) per square decimeter (dm²) of EW.

Partial least square (PLSR)

A supervised multivariate model was built to predict EW (Y) in a training set of data by applying the partial least square regression (PLSR) approach. PLSR is a statistical method that combines the theoretical principles of multiple linear regression and principal component analysis (PCA) to address the situations where several highly correlated predictor variables and relatively fewer samples exist. This approach decomposes the response variables (X) into orthogonal scores (T), loadings (P) and the error (E) while simultaneously incorporating the information from the variables:

$$X = TP' + E$$

(1)

Two hundred spectral bands were integrated in the PLS analysis to identify a set of components that best estimated EW content. The RMSE (root mean square error) of the prediction was estimated with a leave-one-out cross-validation analysis (LOOCV) in a subset of the data with 66.7% of the observations. The EW content was predicted in the remaining 33% of the observations (validation set), integrating the optimum number of components detected in the PLSR model. The analysis was conducted with the plsr function included in the pls package in the statistical software R [73].

Narrow and broad empirical spectral indices for the indirect estimation of the EW content of leaves

Vegetation indices developed to phenotype the morphological and physiological characteristics of the plant (Tables 1 and 2) were calculated with the light reflectance. The correlation coefficients and the statistical significance of each of the VIs and the EW content were estimated with the cor function in the statistical software R.

Spectral indices were calculated by combining the light reflectance at different wavelengths. Eleven mathematical combinations of the spectral bands were calculated including the reflectance at every 3 nm. The adjustment of the spectral resolution was done by averaging the percentage of light reflectance every 3 wavelengths. Additional combinations were also calculated with the average light reflectance of the spectral range of the blue, green, yellow, red and NIR regions (Table 2). In each of the mathematical combinations one, two and up to three bands were integrated. The significance of the linear models and the coefficient of determination (R²) of every pairwise combination of the spectral bands was calculated with the lm function in a leave one out cross-validation (LOOCV) analysis. The format of the linear models was the follow:

$${y_i} = {\beta _0} + {\beta _1}\left( {{x_i}} \right)$$

(2)

where ${\mathrm{y}}_{\mathrm{i}}$ corresponds to the EW content (mg·dm⁻²), ${\mathrm{x}}_{\mathrm{i}}$ is the value of the spectral index, and ${\upbeta }_{0}$ and ${\upbeta }_{1}$ are the intercept and the slope of the model, respectively. A set of indices were selected based on the highest R² criteria in a LOOCV analysis. An inverse prediction approach was applied for assessing the predictive capability of the indices selected. These selected indices were fitted in a calibration data set (66% of the observations) using the same format as in Eq 2 but with ${y}_{i}$ as the spectral index and ${x}_{i}$ the EW content [74, 75]. The step after the models were fitted was to solve for ${x}_{i}$ with the estimates of the parameters calculated in the training data set with the follow equation: ${\mathrm{X}}_{\mathrm{i}}=\frac{{\mathrm{y}}_{\mathrm{i}}-{\upbeta }_{0}}{{\upbeta }_{\mathrm{i}}}$. Estimates of the RMSE were bootstrapped 1000 times for each of indices and in the remaining 34% of the observations with the estimated values of EW with the equation and the chemically measured EW content.

Stepwise regression (SWR)

The spectral bands statistically associated with the EW content (Figs. 1 and  2) were incorporated to build a multivariate model for prediction. These variables were included and/or removed based on the significance of the partial F-values. The final models were selected when the inclusion of more spectral bands was statistically not justifiable. This analysis was conducted with the PROC REG statement in the statistical analysis software SAS [84] in a random training set of the total data (60% of the observations). The prediction models were selected based on the Mallows’ Cp estimator. The estimate of the RMSE of each prediction model was calculated in the remaining 40% of the observations, and final models were choosen based on the lowest value of the RMSE.

Plant material and field experimental trials for validation

Two panels of spring wheat cultivars were evaluated during the agronomic cycle in 2013 at the Norman E. Borlaug Experimental Station (CENEB), Ciudad Obregon, Sonora in northwestern Mexico (27.20°N, 109.54°W, 38 masl). The panels were two sets of 114 and 216 landraces and product of interspecific hybridization with wild relatives. These experimental trials were established as an alpha-lattice design with two replications in a raised bed system with two rows per bed and were planted 80 days later than the normal planting date of wheat in the Yaqui Valley. Late panting allowed the genotypes to be exposed to average daily temperatures of 28 °C and maximum environmental temperatures of 39 °C during the heading and anthesis stages of the crop. There was an intern row spacing within each bed of 10 cm (cm), and a space between beds of 80 cm. In 2016, an additional panel of synthetic derived wheat lines (SDLs) was also evaluated in Bushland and College Station, Texas under non-irrigated conditions. The panel of SDLs were established in an alpha-lattice design with two replications and a plot size of 3.0 × 1.5 m.

Ground base radiometric measurements for the direct validation

The canopy reflectance was collected with a FieldSpec 4 Hi-Res spectroradiometer that captured the light reflected in 2151 continuous bands with a spectral resolution of 3 nm (nm) from 0.35 to 0.7 µm and 8 nm from 1.4 to 2.1 µm. The measurements were collected from 11 AM to 1 PM by placing the optic fiber of the spectroradiometer 40 cm (cm) above the canopy. The sensor was radiometrically calibrated with a white BaSO4 reference card for 100% reflectance and by blocking the light intercepted by the optic fiber for 0% reflectance. Ten readings per plot were captured and the average response of these signatures at a single wavelength was utilized in further analysis.

Airborne hyperspectral information

A set of aerial hyperspectral images were captured from the panel of wheat SDLs in College Station, Texas. The images were obtained with an Aisa KESTREL-10 hyperspectral camera, developed by SPECIM®, and mounted on a Cessna 355 II aircraft. An altitude of 5000 feet (ft) and a speed of 192 km per hour (km/h) were maintained through the flight of the aircraft. The camera captured 120 spectral bands with spectral and spatial resolutions of 5 nm and 0.25 m, respectively. For calibration, four 8 m by 8 m ground tarps with nominal reflectance values of 8%, 16%, 32% and 48% were laid out in the field and captured in the hyperspectral images. The exact percentage of reflectance of the tarps was captured with a Hand-held 2 spectroradiometer. The range of this spectroradiometer is from 0.325 to 1.075 µm, and the spectral resolution is 3 nm. The hyperspectral images were georeferenced and ensembled using the image analysis software ERDAS®. Digital counts (DCs) were extracted individually for each tarp and for individual plots with the software ENVI®. A linear regression model for a single spectral band was developed using light reflectance captured with the spectroradiometer from the tarps as the response variable (dependent) and the DCs as the independent variable. The linear equations were utilized for the estimation of the total canopy reflectance of each of the two hundred spectral bands in each plot.

Efficiency of indirect selection of EW with spectral information

The fourteen spectral indices (EWI) and the eleven regression models developed in this study (Table 4) were calculated with the ground based and aerial spectral information collected in the four experimental trials. Each of the indirect selection methods (spectral indices and models) was considered as an independent variable and subjected to an analysis of variance (ANOVA) for an alpha-lattice experimental design with the lmer function included in the package lme4 in the statistical software R. The variance components were extracted with the function varComp and estimates of the heritability in a broad sense (${\mathrm{h}}^{2})$ calculated according to the formula described by [85]: ${\mathrm{h}}^{2}=\frac{{\sigma }_{g}^{2}}{{\sigma }_{g}^{2}+({\sigma }_{e}^{2}/r)}$ where ${\sigma }_{g}^{2}$ corresponds to the genetic variance, ${\sigma }_{e}^{2}$ to the error variance and r is the number of replications in the experimental trial.

The statistical relationship of the target trait (phenotypic correlation) and the spectral method was calculated with the cor function of the stats package, while the genetic relationship of the traits (genotypic correlation) was estimated with the following equation: ${\upsigma }_{\mathrm{g}}=\frac{{\mathrm{COV}}_{\mathrm{XY}}}{\surd {\mathrm{Var}}_{\mathrm{x}}{\mathrm{Var}}_{\mathrm{y}}}$, where COV_XY corresponds to the covariance estimate of the EWI and EW content calculated with the chemical method, Var_x is the variance of the EWI and Var_y is the variance of EW [86]. The COV_XY was calculated with the cov function, and Var_x and Var_y with the varComp function in R. The genetic correlation of the traits helps us to understand the pleiotropic action of the genes controlling the trait and its indirect selection.

The genetic gain (GG), the genetic advance (GA), the genetic advance with respect to the mean (GAM), the expected response to selection (R), the correlated response to selection (CR), the relative efficiency of indirect selection (RE) were all calculated according to Falconer [87]. GG, GAM, R, CR and RE were estimated as follows:

$\mathrm{GG}={\mathrm{h}}^{2}*\mathrm{SDiff}$, where ${\mathrm{h}}^{2}$ is the estimate of the broad sense heritability of the trait, and SDiff is the selection differential of the trait (EW) with a selection pressure of 10% ($\mathrm{SDiff}={\overline{\mathrm{x}} }_{\mathrm{p}}-{\overline{\mathrm{x}} }_{\mathrm{S}}$).

$\mathrm{GA}=\mathrm{K}\left({\upsigma }_{\mathrm{p}}\right){\mathrm{h}}^{2}$ where K is the selection differential, ${\upsigma }_{\mathrm{p}}$ is the phenotypic standard deviation of every spectral index or prediction model, and ${h}^{2}$ corresponds to the broad sense heritability. The k was estimated for 10% selection intensity as $\mathrm{k}=\overline{{\mathrm{x}}_{\mathrm{p}}}-\overline{{\mathrm{x}}_{\mathrm{s}}}$, where $\overline{{\mathrm{x}}_{\mathrm{p}}}$ and $\overline{{\mathrm{x}}_{\mathrm{s}}}$, are the population mean and the mean of the selected individuals, respectively.

$\mathrm{GAM }\left(\mathrm{\%}\right)=\frac{\mathrm{GA}}{\overline{\mathrm{x}}}\mathrm{ x }100$, where $\overline{\mathrm{x}}$ is the grand mean of the specific character.

$\mathrm{R}={\mathrm{h}}_{\mathrm{x}}{\upsigma }_{\mathrm{x}}$, where ${\mathrm{h}}_{\mathrm{x}}$ is the square root of the heritability and ${\sigma }_{x}$ is the genotypic standard deviation.

$\mathrm{CR}={\mathrm{h}}_{\mathrm{x}}{\mathrm{r}}_{\mathrm{gx}}{\upsigma }_{\mathrm{gy}}$, where ${\mathrm{h}}_{\mathrm{x}}$ is the square root of the heritability for trait X (spectral index), ${\mathrm{r}}_{\mathrm{gx}}$ is the genetic correlation of the spectral index and EW, and ${\upsigma }_{\mathrm{gy}}$ is the genotypic standard deviation of trait Y (EW).

$\mathrm{RE}=\frac{\mathrm{CR}}{\mathrm{R}}$, where CR is the correlated response to selection and R is the expected response to selection for the trait.