# Dissecting Bread Wheat Heterosis through the Integration of Agronomic and Physiological Traits

Kevin Gimenez, Pierre Blanc, Odile Argillier, Jean-Baptiste Pierre, Jacques

Le Gouis, Etienne Paux

► To cite this version:

Kevin Gimenez, Pierre Blanc, Odile Argillier, Jean-Baptiste Pierre, Jacques Le Gouis, et al. Dissecting Bread Wheat Heterosis through the Integration of Agronomic and Physiological Traits. Biology, 2021, 10 (9), pp.907. 10.3390/biology10090907. [https://hal.inrae.fr/hal-03438973](https://hal.inrae.fr/hal-03438973)

HAL Id: hal-03438973

https://hal.inrae.fr/hal-03438973

Submitted on 23 Nov 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Article

# Dissecting Bread Wheat Heterosis through the Integration of Agronomic and Physiological Traits

Kevin Gimenez <fnref n="1" />,<fnref n="2" /> ID, Pierre Blanc <fnref n="2" />, Odile Argillier <fnref n="2" />, Jean-Baptiste Pierre <fnref n="2" />, Jacques Le Gouis <fnref n="1" /> ID
and Etienne Paux <fnref n="1" />,*,† ID

<figure><img src="image_2.png" /></figure>

check for
updates

**Citation:** Gimenez, K.; Blanc, P.; Argillier, O.; Pierre, J.-B.; Le Gouis, J.; Paux, E. Dissecting Bread Wheat Heterosis through the Integration of Agronomic and Physiological Traits. *Biology* **2021**, *10*, 907. https://doi.org/10.3390/biology10090907

Academic Editors: Pierre Devaux and Pierre Sourdille

Received: 15 August 2021

Accepted: 6 September 2021

Published: 13 September 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_3.png" /></figure>

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

**Simple Summary:** To meet the challenge of feeding almost 10 billion people by 2050, wheat yield has to double by 2050. However, over the past 20 years, yield increase has slowed down and even stagnated in the main producing countries. Similar to what has been observed in maize, hybrids have been suggested as a solution to overcome yield stagnation in wheat. However, wheat heterosis, i.e., the fact that a progeny surpasses the performances of its parents, is still limited and poorly understood. To better characterize this phenomenon, we developed and phenotyped for physiological and agronomic traits 91 hybrids and their nineteen female and sixteen male parents. We showed that hybrids had a longer grain filling phase that led to bigger grains and an increased thousand kernel weight. This resulted in a better yield for 86% of hybrids compared to the average yield of their parents. In addition, hybrids appeared to be less affected by the negative correlation between protein content and yield compared to pure lines. These results shed light on the physiological bases underlying yield heterosis in wheat, paving new ways to breed for better wheat hybrids that can help to meet agriculture's challenges.

**Abstract:** To meet the challenge of feeding almost 10 billion people by 2050, wheat yield has to double by 2050. However, over the past 20 years, yield increase has slowed down and even stagnated in the main producing countries. Following the example of maize, hybrids have been suggested as a solution to overcome yield stagnation in wheat. However, wheat heterosis is still limited and poorly understood. Gaining a better understanding of hybrid vigor holds the key to breed for better varieties. To this aim, we have developed and phenotyped for physiological and agronomic traits an incomplete factorial design consisting of 91 hybrids and their nineteen female and sixteen male parents. Monitoring the plant development with normalized difference vegetation index revealed that 89% of the hybrids including the five higher yielding hybrids had a longer grain filling phase with a delayed senescence that results in larger grain size. This average increase of 7.7% in thousand kernel weight translated to a positive mid-parent heterosis for grain yield for 86% of hybrids. In addition, hybrids displayed a positive grain protein deviation leading to a +4.7% heterosis in protein yield. These results shed light on the physiological bases underlying yield heterosis in wheat, paving new ways to breed for better wheat hybrids.

**Keywords:** hybrids; wheat; heterosis; yield factor; grain protein deviation; grain filling; senescence

## 1. Introduction

With 216 million hectares and an annual production of 765 million tons, bread wheat (*Triticum aestivum* L.) is one of the most important crops worldwide and the staple food for one third of the world population [1]. It is also a major renewable resource for feed and industrial raw materials. Today's agriculture has to face an unprecedented challenge: to keep pace with the human demand in an environmentally and socially sustainable

<footnote n="1">¹ Université Clermont Auvergne, INRAE, Genetics, Diversity and Ecophysiology of Cereals, 63000 Clermont-Ferrand, France; kevin.gimenez@inrae.fr (K.G.); jacques.le-gouis@inrae.fr (J.L.G.)</footnote>

<footnote n="2">² Syngenta France SAS, 28000 Chartres, France; pierre.blanc@syngenta.com (P.B.); odile.argillier@syngenta.com (O.A.); jean-baptiste.pierre@syngenta.com (J.-B.P.)</footnote>

<footnote>* Correspondence: etienne.paux@inrae.fr</footnote>

<footnote>† Current address: VetAgro Sup, Campus Agronomique de Clermont, 63370 Lempdes, France.</footnote>
manner [2]. To feed a world population that is expected to reach 9.7 billion by 2050 [3], wheat yield should increase by 2.4% per year over the next 30 years [4]. This goal would be achievable under the assumption of favorable growing conditions but is unlikely under climate change that affects not only the average yield but also its stability [5-7]. Wheat yields have been multiplied by 2.5 between 1960 and 1995 as a result of genetic progress as well as the development of cultivation techniques that have been applied to crops since the "Green Revolution" [5,8,9]. However, yield increase has slowed down to an annual rate of approximately 1% per year over the past 20 years, and even stagnating in the main producing countries such as Europe [10]. Indeed, with its high yielding wheat production, the EU27 is the first world wheat producer and its production contributes significantly (22%) to the world supply. France alone ranks fifth in the world and first in the EU for both production and export. However, like in many other countries, French annual yield has been stagnating since the end of the 1990s and the genetic progress was only able to compensate for more unfavorable growing conditions [5,11]. Thus, there is an urgent need to accelerate genetic progress for yield potential as well as to improve tolerance to biotic and abiotic stresses that are expected to increase in frequency and intensity as a consequence of climate change [12].

Yield increase slow-down has also been observed in other crops such as rice and soybean [4]. By contrast, maize is doing better, yet not sufficiently, with an annual rate of +1.6% since the beginning of the 21st century. This current success can be at least partially attributed to the development of hybrid varieties since the 1930s that led to a 5.5-fold yield increase over 75 years [13]. Even though this annual yield gain results from several complementary factors including agronomic practices and development of new varieties, the genetic gain for maize was the main factor for +100 kg/ha/year [14]. By contrast, average wheat genetic gain in Europe is about +40 kg/ha/year [15]. Consequently, the use of hybrid varieties appears as a possible way to boost yield increase in wheat. A study performed on wheat commercial breeding programs during 21 years in North America reported a genetic gain for hybrids of +1.55% per year compared to 0.84% for pure lines [16]. However, despite having been a matter of interest since the 1970s, the deployment of hybrid wheat has long been hampered by two interconnected factors. The first factor is that hybrid wheat production has long relied on the use of chemical hybridizing agents (CHA), that are difficult to use, and also a low rate of cross pollination, both difficulties resulting in high cost for seed production [17]. This leads to the second factor which is the yield heterosis needed to compensate additional cost for hybrid seed production. In the 1980s and 1990s, a 10% yield heterosis was often observed for some combinations of hybrid wheat in experimental trials [18,19] and up to 13.5% for the best hybrids after breeding [20]. However, very few hybrids overpassed the +13% threshold estimated for profitability [21] and the interest decreased and few acres were grown with hybrid wheat varieties [17]. During the 2010s, the advent of the wheat genome sequencing projects allowed for a better knowledge of the genes involved in the hybrid production (cytoplasmic male sterility and fertility restauration) leading to a new generation of hybrids with high heterosis [22-26]. Additional interest was also observed in hybrids that showed a better grain protein deviation and a higher yield stability [23,27].

Heterosis for grain yield has mainly been associated with thousand kernel weight [28,29], suggesting that hybrids and pure lines may differ for the dynamic of canopy senescence. Normalized Difference Vegetation Index (NDVI) has been proposed as a way to rapidly characterize senescence kinetics on a large number of genotypes [30,31]. Measured on pure lines, NDVI has highlighted the importance of specific stages in yield determination [30-33]. Integrative approaches using NDVI measurements revealed particularly that stay green traits resulting in a delayed senescence could be linked to increased yield [34-37].

The objective of our study was to (i) evaluate heterosis for agronomic and physiological traits and (ii) assess the potential of physiological traits estimated from NDVI to explain heterosis. Our study is based on an incomplete factorial design including 19 females, 16 males, and 92 F1 combinations, which were grown at three locations in northern France.
## 2. Materials and Methods

### 2.1. Plant Material and Genotyping

A set of 136 wheat varieties were selected within the genetic pool from Syngenta breeding company comprising 43 lines carrying fertility restorer (Rf) genes, hereafter referred to as 'males', and 93 lines carrying a cytoplasmic male sterility (CMS) derived from *Triticum timopheevii*, referred to as 'females' [38]. Genotyping was performed by Trait Genetics (Gatersleben, Germany) using the Illumina® iSelect 15K wheat SNP array [39]. Heterozygous data were considered as missing data. Monomorphic markers and markers with more than 20% missing data were removed from the analysis. No minor allele frequency filtering was applied. SNPs were pruned for linkage disequilibrium ($r^2$ = 1) using the plink software [40]. Eventually, a subset of 2966 SNPs was selected to describe the diversity within the genetic panel.

Pearson coefficient correlations were calculated with the R *cor* function for each male x female combination (use = "pairwise(complete.obs") [41]. Dissimilarity matrices were built using the "as.dist (1-correlation)" function. Hierarchical clustering was performed using the "hclust" function (method = "average") and genetic groups were identified using function "cutree". Dendrogram representing the genetic distances between lines of the panel were produced using function "dendro_data_k" of package ggdendro [42]. Principal Component Analyses (PCAs) were conducted with the centered and scaled genotype matrix using package FactoMineR [43]. Individuals were projected onto the two first axes using the ggplot2 package [44].

The 35 parents selected from the genetic panel were genotyped with an improved Axiom array based on the TaBW280K SNP chip [45] and composed of 410K SNP markers which was conducted on the Affymetrix GeneTitan system according to the procedure described by Affymetrix (Axiom® 2.0 Assay Manual Workflow User Guide Rev3, Santa Clara, CA, USA). Allele calling was carried out using a modified version of the Affymetrix proprietary software packages named Affymetrix Power Tools (APT) and SNPolisher™ (http://www.affymetrix.com/estore/partners_programs/programs/developer/tools/devnettools.affx (accessed on 22 February 2021)). The objective was to take into account the specificities of the hexaploid wheat genome. For all SNPs, HomRO and HomFLD were calculated. The HomFLD filter was set to 3.6 (http://media.affymetrix.com/support/developer/downloads/Tools/SNPolisher_User_Guide.pdf (accessed on 22 February 2021)). As a first step, all the probesets were processed with a mild inbred penalty equal to 4 on all the samples. As a second step, the SNPs failing the QC criteria ("Other" and "NoMinorHom") were reprocessed using an inbred penalty of 16. Probesets classified as Off-Target Variants (OTVs) by SNPolisher were analyzed with OTV_caller in the two steps.

Crosses between 19 females with cytoplasmic male sterility [38] and 16 males carrying fertility restorer genes were conducted in an incomplete factorial design and 92 different hybrid combinations were selected for further analyses (Table S2).

### 2.2. Field Experiments

The 92 hybrids, 16 males, and 19 maintainer lines, i.e., isogenic to CMS female lines except for the CMS gene [46], were tested in an augmented design with four replicated checks (Rubisko, RGT Cesario, Tenor and LG Absalon) and a non-replicated check (Chevignon) resulting in a design composed of 160 entries in eight blocks. Because of an insufficient amount of seeds, one hybrid (FEM16 x MA25) was not sown in Arvillers. In addition, due to technical problems, FEM47 female and MA21 male lines were not sown in any location. These two lines were involved in five crosses each and were replaced by the 'Hyking' and 'Hypodrom' commercial hybrid varieties.

Previous studies performed in five locations in northern France demonstrated that no difference was observed on hybrid yield between a normal and a 15%-reduced sowing density (data not shown). For this reason, hybrids were sown at an 85% density relative to their parents.
Field evaluations were conducted during season 2019/2020 in three locations in France (Table 1). Treatments and fertilization were managed according to the local agricultural practices. In Moinville-la-Jeulin, with a clayey loam type soil, 190 kg N/ha were applied in three applications (respectively 50, 90, and 50 kg N/ha) as well as 40 kg/ha of phosphorous. Two herbicides were applied in autumn. Three fungicide treatments were sprayed in April and May and two insecticides were performed in May. In Arvillers, with a silt limon type soil, 220 kg N/ha in three applications, two herbicides, and three fungicides were applied. In Pomacle, with a chalk type soil, 190 kg N/ha in three applications, one herbicide, and one insecticide before winter and five fungicide treatments were applied. The three nitrogen fertilization were applied during three specific stages: tillering (Z = 26), stem elongation (Z = 30), and booting (Z = 41).

<table><caption>**Table 1.** Key characteristics of the three field trial sites in France where bread wheat experiments were conducted.</caption><thead><tr><th>Location</th><th>Coordinates</th><th>Plot Size (m²)</th><th>Sowing Date</th><th>Harvest Date</th><th>Cycle Duration (Days)</th><th>Average Air Temperature (°C)</th><th>Sum of Temperatures (°C Days)</th><th>Number of Days <0 °C</th><th>Number of Days >25 °C</th><th>Cumulative Rainfall (mm)</th></tr></thead><tbody><tr><td>Moinville-la-Jeulin</td><td>Lat: 48.38, Lon: 1.7</td><td>12.4</td><td>25/10/2019</td><td>11/07/2020</td><td>260</td><td>10.4</td><td>2709</td><td>35</td><td>13</td><td>347</td></tr><tr><td>Arvillers</td><td>Lat: 49.73, Lon: 2.65</td><td>10.6</td><td>29/10/2019</td><td>29/07/2020</td><td>274</td><td>10.4</td><td>2854</td><td>39</td><td>16</td><td>382</td></tr><tr><td>Pomacle</td><td>Lat: 49.33, Lon: 4.15</td><td>12.0</td><td>12/10/2019</td><td>19/07/2020</td><td>281</td><td>10.5</td><td>2953</td><td>47</td><td>21</td><td>430</td></tr></tbody></table>

## 2.3. Phenotyping Measurements

For all 92 hybrids, 35 parents, and the checks in all locations, eleven agronomic traits were measured for each plot (Table S3). Grain protein content (GPC) was determined by near infrared spectroscopy using an NIR XDS analyzer (FOSS; Hillerød, Denmark). Two hundred grains were scanned with the MARViN system and analyzed with the MARViN software [47] to estimate grain dimensions (area, width, and length). GPC, thousand kernel weight (TKW), specific weight (SW), and seed dimension data were missing for MA25 and six hybrids in Arvillers (FEM18 x MA22, FEM18 x MA24, FEM09 x MA21, FEM70 x MA25, FEM09 x MA28 and FEM18 x MA08) and two hybrids in Pomacle (FEM09 x MA28 and FEM18 x MA08).

Normalized Difference Vegetation Index (NDVI) was monitored from sowing to harvest with a GreenSeeker® crop sensing system (Trimble, CA, USA). For each plot, between 60 and 80 measurements were performed 50 cm above the canopy. The number of timepoints was on average 13 for Pomacle and 14 for Moinville-la-Jeulin and Arvillers.

## 2.4. Modelling Plant Dynamics

For each genotype, NDVI values were used to establish a plant development model. The developmental cycle was divided in three phases: the “growing phase”, from seedling to maximum NDVI (Nmax); the “flowering phase”, from Nmax to beginning of senescence; and the “senescence phase” from beginning of senescence to end of senescence (Figure 1). NDVI measurement dates were converted into thermal time and expressed in degree days (°C days) with heading date as a reference point (0 °C days) [48]. Thermal time was calculated with a base temperature of 0 °C. For each genotype, three equations were established to fit the three development phases:

$$ \text{Growing phase: } f(\text{TT}) = a_1 \cdot \text{TT} + b_1 \quad (1) $$

$$ \text{Flowering phase: } g(\text{TT}) = a_2 \cdot \text{TT} + b_2 \quad (2) $$

$$ \text{Senescence phase: } h(\text{TT}) = \text{NDVI}_{\text{final}} + \frac{\text{NDVI}_{\text{amplitude}}}{1 + e^{a_3 \cdot \text{TT} + b_3}} \quad (3) $$

where $a_1$ and $a_2$ are the slopes, and $b_1$ and $b_2$ the intercepts. $NDVI_{final}$ is the lowest NDVI value before harvesting and $NDVI_{amplitude}$ corresponds to the difference between the NDVI at the beginning of senescence and $NDVI_{final}$. The exponential term of Equation
(3) is composed of a linear model with $a_3$ being the slope and $b_3$ the intercept. Thermal time (TT) was centered on the heading date (0 °C days). Functions for Equations (1) and (2) were calculated with the “LINEST” Microsoft Excel function. Parameters $a_3$ and $b_3$ were determined with the Excel “Solver” add-in using the nonlinear least-squares minimization method.

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Dynamic changes of wheat vegetation index as a function of thermal time centered on heading date. Example for one wheat line of NDVI modelling based on three phases. NDVI was monitored (blue dots) from sowing to harvest. f, g, and h curves represent the functions calculated in Equations (1)–(3). Descriptions of modeled traits are given in Table S3.</figcaption></figure>

Several physiological traits linked with stay green were calculated from these functions: slopes, thermal time when 10%, 50%, 90%, and 99% (respectively TFN90, TFN50, TFN10, and TFN1) of NDVI$_{amplitude}$ was reached [34]. The beginning of the senescence when 90% of the NDVI remained (TFN90) was calculated as $-\frac{\ln(9)-b_3}{a_3}$, the mid-senescence (TFN50) was calculated as $-\frac{b_3}{a_3}$, the inflexion point corresponding to the end of the rapid phase of senescence (TFN10) was calculated as $\frac{\ln(9)-b_3}{a_3}$, and the end of the senescence (TFN1) was calculated as $\frac{\ln(99)-b_3}{a_3}$. Area under the curve of the growing phase (GPA) was calculated from sowing to thermal time for Nmax. Area under the curve of the “flowering phase” (FPA) was calculated between thermal time for Nmax and thermal time when H(x = 0) where H is the primitive of Equation (3). Area under the curve of the “senescence phase” (SPA) was calculated between thermal time when H(x = 0) and TFN1. We called “area of the declining phase” (DPA) the area between thermal time for Nmax to TNF1 by summing FPA and SPA.

## 2.5. Statistical Analysis

All traits were analyzed for each environment separately with the following fixed ANOVA model for estimating adjusted means:

$$Y_{ij} = \mu + \mu_i + \beta_j + \varepsilon_{ij} \\ \varepsilon_{ij} \sim N(0, \sigma_r), \quad \sum_{1}^{8} \beta_j = 0, \quad \sum_{1}^{133} \alpha_i = 0 \quad (4)$$

where $Y_{ij}$ is the phenotypic trait for the ith genotype in the jth block and $\varepsilon_{ij}$ the residual error. $\mu$ corresponds to the intercept term of the model; $\alpha_i$ is the genotypic effect of the ith genotype; $\beta_j$ is the jth block effect. The best linear unbiased estimators (BLUE) were calculated in R with function “emmeans” from the package “emmeans” for each genotype [49].
All traits were analyzed for each environment separately with the following mixed ANOVA model for calculating heritability:

$$ 
\begin{gathered}
Y_{ij} = \mu + \mu_i + \beta_j + \varepsilon_{ij} \\
\varepsilon_{ij} \sim N(0, \sigma_r), \quad \sum_{1}^{8} \beta_j = 0, \quad \alpha_i \sim N(0, \sigma_G)
\end{gathered}
\qquad (5) $$

where $Y_{ij}$ is the phenotypic trait for the $i$th genotype in the $j$th block and $\varepsilon_{ij}$ its residual error. $\mu$ corresponds to the intercept term of the model; $\alpha_i$ is the random genotypic effect of the $i$th genotype; $\beta_j$ is the $j$th block fixed effect.

Broad sense heritability was calculated as the ratio of genotypic variance to phenotypic variance with the following formula [50]:

$$ h^2 = \frac{\sigma_G^2}{\sigma_G^2 + \frac{\sigma_r^2}{n}} \qquad (6) $$

where $\sigma_G^2$ corresponds to the genotypic variance, $\sigma_r^2$ refers to the residual variance, and $n$ is the average number of repetitions.

The traits were also analyzed for estimating adjusted means by including the location effect in a fixed model called “combined environments”:

$$ 
\begin{gathered}
Y_{ijk} = \mu + \mu_i + \beta_{jk} + L_k + \gamma_{ik} + \varepsilon_{ijk} \\
\varepsilon_{ijk} \sim N(0, \sigma_r), \quad \sum_{1}^{133} \alpha_i = 0, \quad \sum_{1}^{8} \beta_{jk} = 0, \quad \sum_{1}^{3} L_k = 0
\end{gathered}
\qquad (7) $$

where $Y_{ijk}$ is the phenotypic trait for the $i$th genotype in the $j$th block at $k$th location and $\varepsilon_{ijk}$ the residual error. $\mu$ corresponds to the intercept term of the model, $\alpha_i$ is the genotypic effect of the $i$th genotype, $\beta_{jk}$ is the $j$th bloc effect at $k$th location, $L_k$ is the location effect of the $k$th location, $\gamma_{ik}$ is the interaction effect for the $i$th genotype at the $k$th location. The best linear unbiased estimators (BLUE) were used in the linear models for estimating the genotypic effect for each genotype.

All traits were analyzed for combined environments with the following mixed ANOVA model for calculating heritability:

$$ 
\begin{gathered}
Y_{ijk} = \mu + \mu_i + \beta_{jk} + L_k + \gamma_{ik} + \varepsilon_{ijk} \\
\varepsilon_{ijk} \sim N(0, \sigma_r), \quad \sum_{1}^{8} \beta_{jk} = 0, \quad \alpha_i \sim N(0, \sigma_G), \quad \sum_{1}^{3} L_k = 0
\end{gathered}
\qquad (8) $$

where $Y_{ijk}$ is the phenotypic trait for the $i$th genotype in the $j$th block at $k$th location and $\varepsilon_{ijk}$ the residual error. $\mu$ corresponds to the intercept term of the model, $\alpha_i$ is the genotypic effect of the $i$th genotype, $\beta_{jk}$ is the $j$th bloc effect at $k$th location, $L_k$ is the location effect of the $k$th location, $\gamma_{ik}$ is the interaction effect for the $i$th genotype at the $k$th location.

Broad sense heritability was calculated as the ratio of genotypic variance to phenotypic variance with the following formulas [50]:

$$ h^2 = \frac{\sigma_G^2}{\sigma_G^2 + \frac{\sigma_r^2}{nl} + \frac{\sigma_r^2}{l}} \qquad (9) $$

where $\sigma_G^2$ corresponds to the genotypic variance, $\sigma_r^2$ refers to the variance of the interaction, $\sigma_r^2$ refers to the residual variance, $n$ is the average number of repetitions, and $l$ is the number of locations.

Pearson coefficient correlations were calculated with the R cor.test function. Correlations matrices were plotted in R using the ggpairs function from the GGally package [51]. Traits were ordered with the following order: measured, estimated, then modelled. Linear models were estimated using the function 'lm'. For the different performed statistical tests, threshold of significative p-value was fixed to $\alpha = 0.05$.
## 2.6. Heterosis

Heterosis was calculated as the percentage of either mid-parent or best parent heterosis as [52]:

$$ \text{Mid-parent heterosis: } Het_i = \frac{HYB_i - mean(M_i + F_i)}{mean(M_i + F_i)} \times 100 \quad (10) $$

and

$$ \text{Best parent heterosis: } Bet_i = \frac{HYB_i - \max(M_i, F_i)}{\max(M_i, F_i)} \times 100 \quad (11) $$

where $HYB_i$ is the value for the hybrid $i$, $mean(M_i + F_i)$ is the average between the female and male parents for hybrid $i$, $max(M_i, F_i)$ is the highest value between the female and male parents for the hybrid $i$. No heterosis was calculated for hybrids for which the value for at least one parent was missing.

# 3. Results

## 3.1. Parents and Hybrids Selection

In order to produce a panel of genetically diverse wheat hybrids, we selected a set of 136 wheat lines comprising 43 genotypes carrying fertility restorer genes (Rf), hereafter referred to as 'males', and 93 lines carrying a cytoplasmic male sterility (CMS) derived from *Triticum timopheevii* [38], hereafter referred to as 'females'. These 136 genotypes were genotyped with an Illumina iSelect 15K SNP array [39] and a set of 2966 SNPs was selected to describe the diversity within the genetic panel. The dendrogram revealed four different groups and this clustering was consistent with the projection of the genotypic data onto a principal component analysis (Figure 2A,C). Indeed, the first axis discriminated the males and females while the second axis discriminated females into two groups. Males were divided in two groups on the third axis (data not shown). Eventually, 16 males and 19 females originating from different geographic areas were selected to maximize the genetic diversity within this panel. They were subsequently genotyped with the TaBW410K SNP array and 196,305 PolyHighResolution and 18,480 Off-Target Variant SNPs were selected to analyze the genetic diversity [45]. The dendrogram and the projection of the two first PCA axes of these 35 parent lines showed that diversity was conserved and that four main genetic groups were represented (Figure 2B,D).

Out of the 304 possible female x male combinations, 92 hybrids (30%) were selected for further analyses (Tables S1 and S2). These 92 combinations were chosen in order to have at least three hybrids from each male and two from each female. For males, the average number of hybrid combinations was 5.8 (from 3 to 13) while it was 4.8 (from 2 to 9) for females.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Genetic diversity of the wheat lines. (A) Dendrogram showing phylogenetic relationships between 136 wheat lines revealed by 2966 SNP markers (43 males carrying fertility restorer genes in green and orange, 93 females with cytoplasmic male sterility in red and blue). (B) Dendrogram showing phylogenetic relationships between 16 males (in green and orange) and 19 females (in red and blue) wheat lines revealed by 214,785 SNP markers. (C) Projection of the first and second axes of the principal component analysis based on the dissimilarity matrix of 136 wheat lines calculated with 2966 SNPs (males in green and orange, females in red and blue). (D) Projection of the first and second axes of the principal component analysis based on the dissimilarity matrix of 16 males (in green and orange) and 19 females (in red and blue) wheat lines calculated with 214,785 SNPs.</figcaption></figure>

3.2. *Phenotyping and Modelling*

In order to evaluate the physiological and agronomic performances of hybrids, the 92
aforementioned combinations as well as their 35 parents and six checks were sown in three
different locations in northern France: Moinville-la-Jeulin, Arvillers, and Pomacle. These three
locations were chosen in major wheat production areas and close to the official field trial sites
for wheat varieties registration. Average temperatures were similar for the three locations and
differences between sowing and harvesting dates did not exceed three weeks. The number
of days with high temperature (above 25 °C) was higher in Pomacle (n = 21) than in Arvillers
(n = 16) and Moinville-la-Jeulin (n = 13) (Table 1). To characterize the dynamics of the wheat
canopy, NDVI was monitored during the whole plant cycle, from emergence to harvest. NDVI
data were used to model plant development with the aim to compare similar cycle phases
between genotypes. A first phase named "growing phase" was defined with a linear model
from sowing date to the NDVI maximum considered as the maximum of plant development.
A slightly decreasing plateau modelled with a linear regression was observed after maximum
development and was called "flowering phase". A "senescence phase" was then defined by
a sigmoid curve till end of senescence (Figure 1). The "declining phase" included both the
"senescence phase" and the "flowering phase". Due to a lack of measured NDVI points, models
were not estimated for five hybrids in Arvillers (FEM18 x MA24, FEM09 x MA21, FEM70 x
MA25, FEM09 x MA28 and FEM18 x MA08) and two in Pomacle (FEM09 x MA28 and FEM18 x
MA08). For the remaining 466 'genotype x location' combinations, Equations (1)–(3) fitted with
measured points with respectively an average R of 0.98, 0.91 and 1.00 and a standard deviation
of 0.01, 0.11, and 0.00.
### 3.3. Heritability and Genetic Variation

To estimate the proportion of the variability of each trait accounted for by genetic variance, we calculated heritability in each location as well in the three locations combined, hereafter referred to as 'combined environments'. Yield heritability was moderate in Moinville-la-Jeulin ($h^2$ = 0.57) and Arvillers ($h^2$ = 0.66), and high for Pomacle ($h^2$ = 0.80) and in the combined environments ($h^2$ = 0.69). GPC heritability was moderate to high in all three locations (0.66 to 0.81). Thousand Kernel Weight (TKW) and specific weight (SW) had a high heritability in all environments (0.78–0.97). A similar trend was observed for traits controlled by major genes including plant height and heading date (0.84–0.99). Heritability for the different areas under the curve corresponding to plant phases were environment-dependent and varied from low to high values (0.35–0.92). Regarding the senescence indicators (TFN), their heritability was variable for each individual location (0.35–0.78) except for TFN50 and TFN10 in Moinville-la-Jeulin (0.87 and 0.90) and moderate to high in the three combined environments (0.62–0.82) (Figure 4 and Figures S4–S6).

The fact that heritabilities were high enough for a majority of traits allowed us to calculate adjusted means for each trait in each of the three locations and in the combined environments in order to use these means for trait comparisons and heterosis. The average yield was 8.5 t/ha in the three locations, ranging from 6.6 to 9.8 t/ha (Table 2). Grain protein content ranged from 10.5% to 13.7%, with an average of 11.9%. The average specific weight was above the threshold for superior quality required by the market (81.7 kg/hL), ranging from 78.6 to 84.2 kg/hL. TKW displayed a wide range, from 30.8 to 47.1 g. Similarly, a large diversity was observed for heading date and plant height with 20 days and 44 cm differences between extreme genotypes, respectively. The mean total area modelled from NDVI measurements was in the combined environments of 1314 °C days and variable according to the location with the lowest area in Arvillers (1100 °C days) and the highest in Pomacle (1573 °C days); Moinville-la-Jeulin was intermediary (1276 °C days). The duration of the growing phase was similar between environments (~800 °C days) but displayed a wide range between extreme genotypes (from 727 to 899). The main differences between locations were observed during the declining phase that covers both the flowering and senescence phases with less degree days in Moinville-la-Jeulin (466) compared to Pomacle and Arvillers (respectively 723 and 701) for an average of 627 °C days in the combined environments.

<table><caption>**Table 2.** Descriptive statistics of phenotypic traits measured on 92 wheat hybrids and their 35 parents studied in three locations in northern France. Average, minimal, and maximal values as well as heritabilities ($h^2$), genetic variances ($V_α$), and residual errors ($V_ε$) were calculated for all phenotypic traits.</caption><thead><tr><th>Trait</th><th>Unit</th><th>Mean</th><th>Min</th><th>Max</th><th>H<sup>2</sup></th><th>V<sub>α</sub></th><th>V<sub>γ</sub></th><th>V<sub>ε</sub></th></tr></thead><tbody><tr><td>YLD</td><td>t/ha</td><td>8.5</td><td>6.6</td><td>9.8</td><td>0.69</td><td>4.4 × 10<sup>-1</sup></td><td>3.7 × 10<sup>-1</sup></td><td>3.9 × 10<sup>-1</sup></td></tr><tr><td>GPC</td><td>%</td><td>11.9</td><td>10.5</td><td>13.7</td><td>0.81</td><td>4.0 × 10<sup>-1</sup></td><td>2.0 × 10<sup>-1</sup></td><td>3.0 × 10<sup>-1</sup></td></tr><tr><td>SW</td><td>kg/hL</td><td>81.7</td><td>78.6</td><td>84.2</td><td>0.87</td><td>1.1 × 10<sup>0</sup></td><td>6.1 × 10<sup>-1</sup></td><td>4.7 × 10<sup>-1</sup></td></tr><tr><td>HD</td><td>days</td><td>132.6</td><td>125.3</td><td>145.7</td><td>0.97</td><td>3.9 × 10<sup>0</sup></td><td>9.6 × 10<sup>-1</sup></td><td>6.6 × 10<sup>-1</sup></td></tr><tr><td>HT</td><td>cm</td><td>93.4</td><td>72.0</td><td>116.0</td><td>0.95</td><td>7.1 × 10<sup>0</sup></td><td>0.0 × 10<sup>0</sup></td><td>3.1 × 10<sup>0</sup></td></tr><tr><td>SEEDA</td><td>mm<sup>2</sup></td><td>16.0</td><td>12.6</td><td>18.4</td><td>0.94</td><td>8.9 × 10<sup>-1</sup></td><td>2.6 × 10<sup>-1</sup></td><td>3.3 × 10<sup>-1</sup></td></tr><tr><td>SEEDL</td><td>mm</td><td>6.3</td><td>5.5</td><td>7.0</td><td>0.95</td><td>2.2 × 10<sup>-1</sup></td><td>5.6 × 10<sup>-2</sup></td><td>6.7 × 10<sup>-2</sup></td></tr><tr><td>SEEDW</td><td>mm</td><td>3.5</td><td>3.2</td><td>3.7</td><td>0.87</td><td>9.0 × 10<sup>-2</sup></td><td>4.0 × 10<sup>-2</sup></td><td>4.8 × 10<sup>-2</sup></td></tr><tr><td>TKW</td><td>g</td><td>41.1</td><td>30.8</td><td>47.1</td><td>0.84</td><td>2.9 × 10<sup>0</sup></td><td>2.0 × 10<sup>0</sup></td><td>8.2 × 10<sup>-1</sup></td></tr><tr><td>GPY</td><td>t/ha</td><td>1.0</td><td>0.8</td><td>1.1</td><td>0.48</td><td>3.5 × 10<sup>-2</sup></td><td>4.6 × 10<sup>-2</sup></td><td>4.7 × 10<sup>-2</sup></td></tr><tr><td>GPSM</td><td>grains/m<sup>2</sup></td><td>20,614.9</td><td>16,514.2</td><td>24,633.9</td><td>0.76</td><td>1.3 × 10<sup>3</sup></td><td>9.0 × 10<sup>2</sup></td><td>9.5 × 10<sup>2</sup></td></tr><tr><td>Nmax</td><td></td><td>0.9</td><td>0.8</td><td>0.9</td><td>0.68</td><td>1.9 × 10<sup>-2</sup></td><td>1.7 × 10<sup>-2</sup></td><td>1.6 × 10<sup>-2</sup></td></tr><tr><td>GPA</td><td></td><td>834.8</td><td>727.8</td><td>899.0</td><td>0.88</td><td>3.2 × 10<sup>1</sup></td><td>1.2 × 10<sup>1</sup></td><td>1.8 × 10<sup>1</sup></td></tr><tr><td>FPA</td><td></td><td>339.7</td><td>257.4</td><td>467.3</td><td>0.87</td><td>4.3 × 10<sup>1</sup></td><td>1.9 × 10<sup>1</sup></td><td>2.3 × 10<sup>1</sup></td></tr><tr><td>SPA</td><td></td><td>287.7</td><td>230.2</td><td>342.2</td><td>0.35</td><td>1.0 × 10<sup>1</sup></td><td>1.7 × 10<sup>1</sup></td><td>1.8 × 10<sup>1</sup></td></tr><tr><td>DPA</td><td></td><td>627.4</td><td>530.7</td><td>732.0</td><td>0.83</td><td>3.5 × 10<sup>1</sup></td><td>2.0 × 10<sup>1</sup></td><td>2.1 × 10<sup>1</sup></td></tr><tr><td>TA</td><td></td><td>1314.7</td><td>1138.4</td><td>1443.8</td><td>0.75</td><td>4.3 × 10<sup>1</sup></td><td>3.2 × 10<sup>1</sup></td><td>3.0 × 10<sup>1</sup></td></tr><tr><td>TFN90</td><td>°C days</td><td>567.2</td><td>477.4</td><td>628.7</td><td>0.71</td><td>2.6 × 10<sup>1</sup></td><td>2.7 × 10<sup>1</sup></td><td>1.4 × 10<sup>1</sup></td></tr><tr><td>TFN50</td><td>°C days</td><td>688.5</td><td>592.2</td><td>742.5</td><td>0.82</td><td>2.8 × 10<sup>1</sup></td><td>1.7 × 10<sup>1</sup></td><td>1.7 × 10<sup>1</sup></td></tr><tr><td>TFN10</td><td>°C days</td><td>809.8</td><td>704.0</td><td>889.2</td><td>0.76</td><td>3.0 × 10<sup>1</sup></td><td>1.7 × 10<sup>1</sup></td><td>2.6 × 10<sup>1</sup></td></tr><tr><td>TFN1</td><td>°C days</td><td>942.2</td><td>820.3</td><td>1062.6</td><td>0.62</td><td>3.2 × 10<sup>1</sup></td><td>1.6 × 10<sup>1</sup></td><td>4.4 × 10<sup>1</sup></td></tr></tbody></table>

YLD, grain yield; GPC, grain protein content; SW, specific weight; HD, heading date; HT, plant height; SEEDA, seed area; SEEDL, seed length; SEEDW, seed width; TKW, thousand kernel weight; GPY, grain protein yield; GPSM, grains per square meter; Nmax, maximum measured NDVI; GPA, area of the growing phase; FPA, area of the flowering phase; SPA, area of the senescence phase; DPA, area of the declining phase; TA, total NDVI area; TFN90, 90% of the NDVI amplitude remains; TFN50, 50% of the NDVI amplitude remains; TFN10, 10% of the NDVI amplitude remains; TFN1, 1% of the NDVI amplitude remains.
3.4. *Correlation between Agronomic and Physiological Traits*

We then calculated correlations between all variables that were either measured or
calculated for all locations separately as well as in combined environments (Figure 3 and
Figures S1–S3).

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Correlations between adjusted means of phenotypic traits in the combined environments for 92 bread wheat hybrids and their 35 parents. Lower panel represents the linear regression plot. The diagonal shows histogram representing the distribution of each trait. The upper triangle shows Pearson correlations. Positive correlations are colored in blue, negative correlations in red. Significance of the correlation is indicated as either non-significant (n.s.), significant at a *p*-value of 0.05 (*), 0.01 (**), or 0.001 (***) . YLD, grain yield; GPC, grain protein content; GPD, grain protein deviation; SW, specific weight; HT, plant height; HD, heading date; TKW, thousand kernel weight; GPY, grain protein yield; GPSM, grains per square meter; GPA, area of the growing phase; SPA, area of the senescence phase; DPA, area of the declining phase area; TA, total NDVI area; Nmax, maximum measured NDVI; TFN90, 90% of the NDVI amplitude remains; TFN10, 10% of the NDVI amplitude remains.</figcaption></figure>
There was a strong negative correlation between grain yield and grain protein content (R = -0.57, *p*-value = 2.0 × 10⁻¹⁰). As expected, we also found a negative correlation between TKW and the number of grains per square meter (R = -0.63, *p*-value = 2.7 × 10⁻¹¹). Finally, stay-green traits were strongly and negatively correlated with heading date (-0.76 to -0.80). A linear model showed that the onset of senescence (TFN90) was significantly correlated to heading date (*p*-value < 2.2 × 10⁻¹⁶).

### 3.5. *Hybrid Performance and Heterosis*

To evaluate the performances of hybrids, we first compared hybrid values to the average values of both parents (mid-parent heterosis).

An average mid-parent heterosis of +5.0% was observed for yield in all environments, with 86% of hybrids performing better than the average of their parents (Table 3, Figure 4A). The maximum was observed in Moinville-la-Jeulin for hybrid FEM85 x MA25 and in Pomacle for hybrid FEM05 x MA25 that both displayed a +23% heterosis. Following this yield increase, hybrids had an average negative heterosis for grain protein content (-0.7%). However, the grain protein deviation was significatively higher for hybrids compared to parent lines leading to +0.21 percentage points increase for hybrids in protein content (*p*-value = 0.04) (Figure 5A). Consequently, the overall protein yield was higher for hybrids than for their mid-parent (+4.7%).

<table><caption>**Table 3.** Mid-parent and best-parent heterosis of wheat F₁ hybrids studied in three locations. For each trait, the number of hybrids tested, the average heterosis, the *p*-value of the Student test between hybrids and parent lines, the minimal, and maximal values, as well as the heterosis coefficient of variation are indicated.</caption><thead><tr><th rowspan="2">Trait</th><th colspan="4">Mid-Parent Heterosis</th><th colspan="4">Best-Parent Heterosis</th></tr><tr><th>Number of Hybrids</th><th>Mean (%)</th><th><i>p</i>-Value</th><th>Min (%)</th><th>Max (%)</th><th>σ/μ</th><th>Mean (%)</th><th><i>p</i>-Value</th><th>Min (%)</th><th>Max (%)</th><th>σ/μ</th></tr></thead><tbody><tr><td>YLD</td><td>80</td><td>5.0</td><td>3.2 × 10<sup>-12</sup></td><td>-16.4</td><td>14.6</td><td>1.1</td><td>0.9</td><td>n.s</td><td>-17.7</td><td>11.6</td><td>6.6</td></tr><tr><td>GPC</td><td>75</td><td>-0.7</td><td>2.2 × 10<sup>-2</sup></td><td>-5.7</td><td>4.4</td><td>-3.7</td><td>-3.6</td><td>2.2 × 10<sup>-15</sup></td><td>-11.8</td><td>2.7</td><td>-0.9</td></tr><tr><td>SW</td><td>75</td><td>0.8</td><td>4.5 × 10<sup>-15</sup></td><td>-1.0</td><td>2.2</td><td>0.9</td><td>-0.2</td><td>3.8 × 10<sup>-2</sup></td><td>-2.6</td><td>1.5</td><td>-4.1</td></tr><tr><td>HD</td><td>80</td><td>-2.3</td><td>3.7 × 10<sup>-34</sup></td><td>-4.4</td><td>0.1</td><td>-0.4</td><td>-5.2</td><td>4.5 × 10<sup>-35</sup></td><td>-9.8</td><td>-0.7</td><td>-0.4</td></tr><tr><td>HT</td><td>80</td><td>6.7</td><td>1.3 × 10<sup>-29</sup></td><td>-4.7</td><td>14.1</td><td>0.5</td><td>-1.5</td><td>1.4 × 10<sup>-2</sup></td><td>-13.2</td><td>9.6</td><td>-3.6</td></tr><tr><td>SEEDA</td><td>74</td><td>5.6</td><td>5.0 × 10<sup>-23</sup></td><td>-6.9</td><td>13.8</td><td>0.6</td><td>1.5</td><td>3.6 × 10<sup>-3</sup></td><td>-13.4</td><td>9.8</td><td>2.8</td></tr><tr><td>SEEDL</td><td>74</td><td>2.0</td><td>2.6 × 10<sup>-18</sup></td><td>-2.1</td><td>5.1</td><td>0.7</td><td>-0.7</td><td>1.5 × 10<sup>-2</sup></td><td>-7.2</td><td>3.6</td><td>-3.4</td></tr><tr><td>SEEDW</td><td>74</td><td>3.3</td><td>5.8 × 10<sup>-23</sup></td><td>-5.1</td><td>7.8</td><td>0.6</td><td>1.5</td><td>1.3 × 10<sup>-6</sup></td><td>-6.6</td><td>7.3</td><td>1.6</td></tr><tr><td>TKW</td><td>75</td><td>7.7</td><td>1.2 × 10<sup>-26</sup></td><td>-6.6</td><td>16.2</td><td>0.5</td><td>3.6</td><td>5.2 × 10<sup>-6</sup></td><td>-14.7</td><td>12.6</td><td>1.4</td></tr><tr><td>GPY</td><td>75</td><td>4.7</td><td>2.4 × 10<sup>-12</sup></td><td>-9.7</td><td>14.3</td><td>1.0</td><td>1.2</td><td>4.9 × 10<sup>-2</sup></td><td>-12.1</td><td>13.5</td><td>4.3</td></tr><tr><td>GPSM</td><td>66</td><td>-3.1</td><td>7.3 × 10<sup>-7</sup></td><td>-14.4</td><td>8.3</td><td>-1.5</td><td>-7.0</td><td>5.1 × 10<sup>-19</sup></td><td>-18.5</td><td>2.8</td><td>-0.6</td></tr><tr><td>Nmax</td><td>75</td><td>1.2</td><td>1.2 × 10<sup>-4</sup></td><td>-10.7</td><td>6.6</td><td>2.1</td><td>-0.6</td><td>n.s</td><td>-12.6</td><td>3.9</td><td>-4.5</td></tr><tr><td>GPA</td><td>75</td><td>1.2</td><td>4.0 × 10<sup>-5</sup></td><td>-11.7</td><td>6.6</td><td>2.0</td><td>-1.7</td><td>3.3 × 10<sup>-6</sup></td><td>-11.9</td><td>5.1</td><td>-1.7</td></tr><tr><td>FPA</td><td>75</td><td>-4.6</td><td>4.3 × 10<sup>-6</sup></td><td>-21.8</td><td>17.4</td><td>-1.7</td><td>-13.7</td><td>1.2 × 10<sup>-17</sup></td><td>-33.6</td><td>10.6</td><td>-0.8</td></tr><tr><td>SPA</td><td>75</td><td>3.1</td><td>3.1 × 10<sup>-5</sup></td><td>-10.5</td><td>17.2</td><td>1.9</td><td>-2.8</td><td>2.5 × 10<sup>-4</sup></td><td>-16.2</td><td>12.2</td><td>-2.2</td></tr><tr><td>DPA</td><td>75</td><td>-1.3</td><td>4.3 × 10<sup>-3</sup></td><td>-10.7</td><td>10.9</td><td>-2.9</td><td>-5.5</td><td>2.1 × 10<sup>-14</sup></td><td>-18.2</td><td>6.0</td><td>-0.9</td></tr><tr><td>TA</td><td>75</td><td>0.1</td><td>n.s</td><td>-11.4</td><td>5.5</td><td>26.8</td><td>-2.8</td><td>1.6 × 10<sup>-10</sup></td><td>-13.7</td><td>5.0</td><td>-1.2</td></tr><tr><td>TFN90</td><td>75</td><td>5.8</td><td>7.9 × 10<sup>-15</sup></td><td>-5.6</td><td>20.7</td><td>0.9</td><td>0.1</td><td>n.s</td><td>-13.8</td><td>11.5</td><td>50.1</td></tr><tr><td>TFN50</td><td>75</td><td>4.3</td><td>1.5 × 10<sup>-18</sup></td><td>-2.7</td><td>12.1</td><td>0.7</td><td>-0.7</td><td>n.s</td><td>-11.1</td><td>6.5</td><td>-5.8</td></tr><tr><td>TFN10</td><td>75</td><td>3.4</td><td>6.3 × 10<sup>-16</sup></td><td>-1.5</td><td>11.8</td><td>0.8</td><td>-1.3</td><td>2.6 × 10<sup>-3</sup></td><td>-9.8</td><td>6.1</td><td>-2.8</td></tr><tr><td>TFN1</td><td>75</td><td>2.6</td><td>2.9 × 10<sup>-9</sup></td><td>-4.7</td><td>11.6</td><td>1.3</td><td>-1.8</td><td>1.6 × 10<sup>-4</sup></td><td>-11.5</td><td>8.1</td><td>-2.2</td></tr></tbody></table>

YLD, grain yield; GPC, grain protein content; SW, specific weight; HD, heading date; HT, plant height; SEEDA, seed area; SEEDL, seed length; SEEDW, seed width; TKW, thousand kernel weight; GPY, grain protein yield; GPSM, grains per square meter; Nmax, maximum NDVI measured; GPA, area of the growing phase; FPA, area of the flowering phase; SPA, area of the senescence phase; DPA, area of the declining phase; TA, total NDVI area; TFN90, 90% of the NDVI amplitude remains; TFN50, 50% of the NDVI amplitude remains; TFN10, 10% of the NDVI amplitude remains; TFN1, 1% of the NDVI amplitude remains. n.s: non-significant.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Density distribution of mid-parent (A) and best-parent (B) heterosis for major traits in the combined environment. Dashed lines correspond to either the mid-parent or the best-parent values. Average heterosis and standard deviation are represented in black dots and bars respectively. DPA, area of the declining phase; GPA, area of the growing phase; GPC, grain protein content; GPY, grain protein yield; HD, heading date; HT, plant height; SW, specific weight; TA, total NDVI area; TFN10, 10% of the NDVI amplitude remains; TKW, thousand kernel weight; YLD, grain yield.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Negative correlations (A) between grain yield (YLD) and grain protein content (GPC), (B) between grains number per square meter (GPSM) and thousand kernel weight (TKW), (C) between heading date (HD) and beginning of the senescence (TFN90), and (D) between heading date (HD) and grain yield (YLD) for 92 bread wheat hybrids and their 35 parents studied in three locations. Parent genotypes are represented as blue triangles and hybrids as red circles. The blue and red lines correspond to linear regression models for parents and hybrids, respectively. Shaded areas represent the 95% confidence intervals of the linear regression models.</figcaption></figure>

Despite a 15% reduction in sowing density, hybrids had only 3% less grains per square meter at harvest compared to their parents. This decrease was largely compensated by an average +7.7% heterosis for TKW, with 97% of hybrids performing better than their mid-parent in combined environments (Figure 4A, Figures S4A, S5A and S6A). A linear model between TKW and grain per square meter showed that for a given number of grains per square meter, TKW was higher for hybrids compared to pure lines (p-value = 9.5 × 10⁻⁸) (Figure 5B, Figures S7B, S8B and S9B). No significant difference was observed between pure lines and hybrids for the slopes of the negative correlation (p-value = 0.2).

The analysis of the plant development dynamics revealed that hybrids grew faster than their parents until the maximal NDVI value (+1.2% heterosis for the area of the growing phase) and had a delayed beginning of senescence (+5.8% for TFN90). In addition, the rate of senescence was lower for hybrids compared to pure lines with a positive heterosis observed for the area of the senescence phase (+3.1%).

Finally, hybrids displayed a negative heterosis for heading date (−2.3%) that translated into heading data that were approximatively 30 °C days earlier.

These tendencies were also observed when only the five higher yielding hybrids were analyzed. Indeed, grain yield heterosis and TKW heterosis was more pronounced for the best hybrids and they were earlier than their parents (Table S11). However, their earliness is less pronounced and having a later heading date than the whole set of hybrids allows them to increase their yield compared to the parental lines (Figure 5D, Figures S7D, S8D and S9D).

Another method to estimate heterosis is to compare hybrids to the best parent line. For this best-parent heterosis, sometimes referred to as heterobeltiosis [53], we did not observe any general tendency for grain yield, grain protein yield, plant height, and TFN90 but rather positive or negative heterosis depending on the location (Table 3 and Tables S8–S10). However, a large variance was observed and for each trait, a number of hybrids over-passing their best parent. For yield, 54% of hybrids were better than both parents with up to +12% heterosis (Figure 4B, Figures S4B, S5B and S6B). Other agronomic and physiological traits were on average lower for
hybrids than for the best-parent except TKW and seed area, with respectively +3.6% and +1.5% for combined environments. A significant negative best-parent heterosis was observed in each location for the number of grains per square meter (−7%). The total NDVI area under the curve was smaller for hybrids compared to the best parent with a shorter growing phase due to an earlier heading date and with a shorter declining phase.

## 4. Discussion

Grain yield, grain protein content, and yield stability are three of the most important traits for farmers and the most important targets for breeders. However, not only has wheat yield been stagnating and fluctuating in the main producing countries over the past 20 years, but also grain yield and grain protein content are negatively correlated [54]. Hybrids have been proposed as a way to overcome these limitations by taking advantage of the underexploited wheat heterosis. Indeed, heterosis, i.e., the superiority of the hybrid performance compared to parental lines, has been demonstrated in several species such as maize, sunflower, etc. [13,55,56]. By contrast, the difficulties linked to the use of chemical hybridization agent as well as their limited superiority compared to pure lines have long hampered the development of wheat hybrids. The advent of cytoplasmic male sterile wheat varieties was a game changer and opened the way to the development of better wheat hybrids. However, further insight into wheat heterosis is still needed and, to this aim, we developed an incomplete factorial design composed of 92 hybrids from 35 parent lines and phenotyped them for agronomic traits including yield and yield components as well as for physiological traits.

A positive mid-parent heterosis for grain yield was observed for the vast majority of our hybrids (86%), with an average increase of 5% and up to 23%. In addition, 54% of our combinations displayed a higher grain yield than their best parent, 5% being over the +13% threshold of hybrid profitability compared to pure lines [21]. This result is consistent with previous studies conducted on small or large panels that observed an overall grain yield heterosis for the majority of hybrids [23,29,57–60] and confirms that inbreeding depression can be efficiently decreased in a self-pollinating species to optimize yield.

As expected from the grain protein concentration-grain yield negative relationship [54], this positive yield heterosis was accompanied with a slightly negative heterosis for grain protein content. However, this very limited negative heterosis (−0.7%) suggested that hybrids were less impacted by the negative correlation. The linear regression between these two traits confirmed this hypothesis. Indeed, while still showing a negative correlation, hybrids showed a positive grain protein deviation (GPD), i.e., the residual from the regression between grain yield and grain protein content [61]. This GPD translated in a +0.21 percentage point increase in protein content, consistent with a previous study conducted on a large panel composed of 1604 hybrids [23].

Overall, this led to a positive heterosis for final protein yield and suggested hybrids broke down the negative balance between high yield and nitrogen metabolism involving less protein in grains [62] resulting in a better nitrogen use efficiency of hybrid genotypes compared to pure lines under the same level of nitrogen fertilization. It seems due to a better nitrogen absorption from the soil for hybrids while the remobilization is similar between hybrids and pure lines [58]. In pure lines, grain protein deviation is highly dependent on the nitrogen absorption after anthesis [63].

To decipher the bases of grain yield heterosis, we analyzed yield component traits. TKW appeared as one of the main explanatory factors and an average positive heterosis of 7.7% was observed, consistent with the results of Askander and collaborators on a small dataset composed of 12 hybrids [28]. However, TKW is known to be negatively correlated with the number of grains per square meter. In our study, hybrids were sown at reduced density compared to pure lines that might have caused this increase in TKW. However, while the reduction in sowing density was 15%, hybrids had only 3% less grains per square meter at harvest compared to their parents, demonstrating that hybrids can compensate for reduced sowing density with increased tillering and spike density, as
suggested by Milan et al. [64]. In addition, we also showed that, for a given number of
grains per square meter, hybrids had a mean increase of nearly +5% in TKW compared to
pure lines, suggesting a more efficient accumulation of resources in the grain.

Besides measuring agronomic traits, we also used NDVI measurements, a non-destructive and heritable method which provides information on leaf area and radiation intercepted by the canopy [65,66]. Monitoring NDVI all along the plant development allowed us to model the different phases and to identify "stay green" varieties, i.e., those with a delayed senescence, a trait that has been shown to impact positively yield in pure lines [34]. We observed a positive heterosis for TFN10, associated with a delayed senescence, that translates into a stay-green phenotype and, subsequently, a longer grain filling time. Nevertheless, despite a normalization of plant cycle in degree days from heading date, stay green traits remained highly impacted by heading date, as reported by Montazeaud et al. and Bogard et al. [35,37]. However, linear models with heading date as covariable showed that hybrids were more susceptible to delayed senescence than the parental lines (Figure 5C, Figures S7C, S8C and S9C). In other words, not only were the hybrids earlier at heading date, but they also finished their senescence later. This finding was supported by the results of several studies that reported stay green QTLs that did not colocalize with genes involved in flowering time or vernalization [67–71]. Combined with a faster grain filling [72,73], the longer filling duration of hybrids is likely responsible for the accumulation of more resources in the grain, resulting in increased TKW and protein yield, and ultimately in a higher yield. Besides its interest for yield increase, early heading date might also help hybrids to escape environmental stresses, such as drought or high temperature, that are expected to become more and more frequent in the coming years, as a result of climate change [74]. The five higher yielding hybrids were intermediate, i.e., later than the average heading date of hybrids and earlier than their parents. They combine the advantages to escape from stressful environments with a longer plant cycle. In this regard, it might be interesting to evaluate this panel over several years to test this hypothesis and assess the effect of heterosis on yield stability.

Finally, while the selection of parental lines relied on SNP-based genetic distances, it
would be interesting to characterize the 16 male and 19 female lines for genomic structural
variations (SVs). Indeed, it is now widely admitted that SVs play a role in hybrid vigor
through the combination of complementary gene and allele contents [75,76]. While in maize,
hybrids result from crosses between different heterotic groups, no apparent structure was
observed in the registered wheat varieties [77]. Despite this lack of heterotic groups, many
wheat hybrids out-perform their parents. Combining additive-dominance effects with copy
variants, as well as better knowledge of the best pairing of fertility restoration genes, has
the potential to improve the prediction of the best-performing hybrids.

**5. Conclusions**

Hybrid wheat produced with chemical hybridization agent has been widely studied over the past decades. Unfortunately, with this technique, the production of hybrids is fastidious and difficult. The new generation of hybrids produced by cytoplasmic male sterility has helped to decrease these difficulties and to increase the number of combinations tested. Our study focused on heterosis in 92 CMS hybrid combinations compared to their 35 parents. It appears that a majority of hybrids have a mid-parent heterosis for grain yield with a very limited decrease in grain protein content thanks to a positive GPD. Grain yield heterosis is made possible thanks to a positive heterosis on the thousand kernel weight due to a longer filling time for hybrids.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/10.3390/biology1009090/s1, Table S1: Best Linear Unbiased Estimators (BLUEs) for each environment. Table S2: correspondence between the line names of the wheat genotypes and the code names. Table S3: Incomplete factorial bread wheat F1 hybrid mating design: Females (FEM; in blue and red) in rows and males (MA; in orange and green) in columns. “X” indicates hybrids that have been tested in this study. Colors corresponds to genetic groups identified in Figure 2. Table S4: Description of
phenotypic traits. Table S5: Descriptive statistics of phenotypic traits in Moinville-la-Jeulin. Average, minimal, and maximal values as well as heritabilities (H²), genetic variances (Vα), and residual errors (Vε) were calculated for all phenotypic traits. Table S6: Descriptive statistics of phenotypic traits in Arvillers. Average, minimal, and maximal values as well as heritabilities (H²), genetic variances (Vα), and residual errors (Vε) were calculated for all phenotypic traits. Table S7: Descriptive statistics of phenotypic traits in Pomacle. Average, minimal, and maximal values as well as heritabilities (H²), genetic variances (Vα) and residual errors (Vε) were calculated for all phenotypic traits. Table S8: Mid-parent and best-parent heterosis of F1 hybrids in Moinville-la-Jeulin. For each trait, the number of hybrids tested, the average heterosis, the *p*-value of the Student test between hybrids and parent lines, the minimal and maximal value, as well as the heterosis coefficient of variation are indicated. Table S9: Mid-parent and Best-parent heterosis of F1 hybrids in Arvillers. For each trait, the number of hybrids tested, the average heterosis, the *p*-value of the Student test between hybrids and parent lines, the minimal and maximal value, as well as the heterosis coefficient of variation are indicated. Table S10: Mid-parent and Best-parent heterosis of F1 hybrids in Pomacle. For each trait, the number of hybrids tested, the average heterosis, the *p*-value of the Student test between hybrids and parent lines, the minimal and maximal value, as well as the heterosis coefficient of variation are indicated. Table S11: Mid-parent and Best-parent heterosis of the five best F1 hybrids for yield in the three locations.. Figure S1. Correlations between adjusted means of phenotypic traits in Moinville-la-Jeulin. Lower panel represents the linear regression plot. The diagonal shows histogram representing the distribution of each trait. The upper triangle shows Pearson correlations. Positive correlations are colored in blue, negative correlations in red. Significance of the correlation is indicated as either non-significant (n.s.), significant at a *p*-value of 0.05 (*), 0.01 (**), or 0.001 (***) .Figure S2. Correlations between adjusted means of phenotypic traits in Arvillers. Lower panel represents the linear regression plot. The diagonal shows histogram representing the distribution of each trait. The upper triangle shows Pearson correlations. Positive correlations are colored in blue, negative correlations in red. Significance of the correlation is indicated as either non-significant (n.s.), significant at a *p*-value of 0.05 (*), 0.01 (**), or 0.001 (***) .Figure S3. Correlations between adjusted means of phenotypic traits in Pomacle. Lower panel represents the linear regression plot. The diagonal shows histogram representing the distribution of each trait. The upper triangle shows Pearson correlations. Positive correlations are colored in blue, negative correlations in red. Significance of the correlation is indicated as either non-significant (n.s.), significant at a *p*-value of 0.05 (*), 0.01 (**), or 0.001 (***) .Figure S4. Density distribution of mid-parent (A) and best-parent (B) heterosis for major traits in Moinville-la-Jeulin. Dashed lines correspond to either the mid-parent or the best-parent values. Average heterosis and standard deviation are represented in black dots and bars respectively. Figure S5. Density distribution of mid-parent (A) and best-parent (B) heterosis for major traits in Arvillers. Dashed lines correspond to either the mid-parent or the best-parent values. Average heterosis and standard deviation are represented in black dots and bars respectively. Figure S6. Density distribution of mid-parent (A) and best-parent (B) heterosis for major traits in Pomacle. Dashed lines correspond to either the mid-parent or the best-parent values. Average heterosis and standard deviation are represented in black dots and bars respectively. Figure S7. Negative correlations (A) between grain yield (YLD) and grain protein content (GPC), (B) between grains number per square meter (GPSM) and thousand kernel weight (TKW), (C) between heading date (HD) and beginning of the senescence (TFN90) and (D) between heading date (HD) and grain yield (YLD) for 92 bread wheat hybrids and their 35 parents studied in Moinville-la-Jeulin. Parent genotypes are represented as blue triangles and hybrids as red circles. The blue and red lines correspond to linear regression models for parents and hybrids, respectively. Shaded areas represent the 95% confidence intervals of the linear regression models. Figure S8. Negative correlations (A) between grain yield (YLD) and grain protein content (GPC), (B) between grains number per square meter (GPSM) and thousand kernel weight (TKW), (C) between heading date (HD) and beginning of the senescence (TFN90) and (D) between heading date (HD) and grain yield (YLD) for 92 bread wheat hybrids and their 35 parents experimentally studied in Arvillers. Parent genotypes are represented as blue triangles and hybrids as red circles. The blue and red lines correspond to linear regression models for parents and hybrids, respectively. Shaded areas represent the 95% confidence intervals of the linear regression models. Figure S9. Negative correlations (A) between grain yield (YLD) and grain protein content (GPC), (B) between grain number per square meter (GPSM) and thousand kernel weight (TKW), (C) between heading date (HD) and beginning of the senescence (TFN90), and (D) between heading date (HD) and grain yield (YLD) for 92 bread wheat hybrids and their 35 parents
experimentally studied in Pomacle. Parent genotypes are represented as blue triangles and hybrids as red circles. The blue and red lines correspond to linear regression models for parents and hybrids, respectively. Shaded areas represent the 95% confidence intervals of the linear regression models.

**Author Contributions:** Conceptualization, K.G., P.B., O.A., J.-B.P., J.L.G., and E.P.; formal analysis, K.G.; funding acquisition, J.-B.P., J.L.G., and E.P.; investigation, K.G.; methodology, K.G., P.B., O.A., J.-B.P., J.L.G., and E.P.; project administration, J.-B.P., J.L.G., and E.P.; resources, P.B., O.A., and J.-B.P.; software, K.G.; supervision, P.B., O.A., J.L.G., and E.P.; writing—original draft, K.G.; writing—review and editing, O.A., J.-B.P., J.L.G., and E.P. All authors have read and agreed to the published version of the manuscript.

**Funding:** The research leading to these results was funded by Syngenta, INRAE, and ANRT (Association Nationale de la Recherche et de la Technologie) under research grant #2019/0358.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available in the Supplementary Materials.

**Acknowledgments:** The authors would like to thank David Guiard and his field study teams for the significant contribution to the study. The authors thanks Vincent Allard for his expertise on the wheat physiological properties. The authors are also grateful to the ANRT (Association Nationale de la Recherche et de la Technologie) which supported the PhD thesis (2019/0358).

**Conflicts of Interest:** Odile Argillier, Pierre Blanc, and Jean-Baptiste Pierre are from company "Syngenta France SAS". The authors declare no conflict of interest.

## References

1. FAOSTAT. 2019. Available online: http://www.fao.org/ (accessed on 11 March 2021).

2. Godfray, H.C.J.; Beddington, J.R.; Crute, I.R.; Haddad, L.; Lawrence, D.; Muir, J.F.; Pretty, J.; Robinson, S.; Thomas, S.M.; Toulmin, C. Food Security: The Challenge of Feeding 9 Billion People. *Science* **2010**, *327*, 812–818. [CrossRef]

3. United Nations Department of Economic and Social Affairs Population Division. *World Population Prospects 2019: Highlights; Statistical Papers-United Nations (Ser. A);* Population and Vital Statistics Report; United Nations: New York, NY, USA, 2019; ISBN 9789210042352.

4. Ray, D.K.; Mueller, N.D.; West, P.C.; Foley, J.A. Yield Trends Are Insufficient to Double Global Crop Production by 2050. *PLoS ONE* **2013**, *8*, e66428. [CrossRef]

5. Brisson, N.; Gate, P.; Gouache, D.; Charmet, G.; Oury, F.-X.; Huard, F. Why are wheat yields stagnating in Europe? A comprehensive data analysis for France. *Field Crop. Res.* **2010**, *119*, 201–212. [CrossRef]

6. Porter, J.R.; Semenov, M.A. Crop responses to climatic variation. *Philos. Trans. R. Soc. B Biol. Sci.* **2005**, *360*, 2021–2035. [CrossRef] [PubMed]

7. Tester, M.; Langridge, P. Breeding Technologies to Increase Crop Production in a Changing World. *Science* **2010**, *327*, 818–822. [CrossRef]

8. Charmet, G. Wheat domestication: Lessons for the future. *Comptes Rendus Biologies* **2011**, *334*, 212–220. [CrossRef]

9. This week in Nature-Editorials Wheat Lag. *Nature* **2014**, *507*, 399–400. [CrossRef]

10. Eurostat Agricultural Production—Crops. Available online: https://ec.europa.eu/eurostat/fr/data/database (accessed on 23 April 2021).

11. Oury, F.-X.; Godin, C.; Mailliard, A.; Chassin, A.; Gardet, O.; Giraud, A.; Heumez, E.; Morlais, J.-Y.; Rolland, B.; Rousset, M.; et al. A study of genetic progress due to selection reveals a negative effect of climate change on bread wheat yield in France. *Eur. J. Agron.* **2012**, *40*, 28–38. [CrossRef]

12. IPCC. *Climate Change 2014: Synthesis Report; Contribution of Working Groups I, II and III to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change; Core Writing Team, Pachauri, R.K., Meyer, L.A., Eds.;* IPCC: Geneva, Switzerland, 2014; 151p.

13. Troyer, A.F. Adaptedness and heterosis in corn and mule hybrids. *Crop Sci.* **2006**, *46*, 528–543. [CrossRef]

14. Duvick, D.N. The Contribution of Breeding to Yield Advances in maize (*Zea mays* L.). *Adv. Agron.* **2005**, *86*, 83–145.

15. Le Gouis, J.; Oury, F.X.; Charmet, G. How changes in climate and agricultural practices influenced wheat production in Western Europe. *J. Cereal Sci.* **2020**, *93*, 102960. [CrossRef]

16. Koemel, J.E.; Guenzi, A.C.; Carver, B.F.; Payton, M.E.; Morgan, G.H.; Smith, E.L. Hybrid and Pureline Hard Winter Wheat Yield and Stability. *Crop Sci.* **2004**, *44*, 107–113. [CrossRef]

17. Longin, C.F.H.; Mühleisen, J.; Maurer, H.P.; Zhang, H.; Gowda, M.; Reif, J.C. Hybrid breeding in autogamous cereals. *Theor. Appl. Genet.* **2012**, *125*, 1087–1096. [CrossRef] [PubMed]
18. Borghi, B.; Perenzin, M.; Nash, R.J. Agronomic and qualitative characteristics of ten bread wheat hybrids produced using a chemical hybridizing agent. *Euphytica* **1987**, *39*, 185–194. [CrossRef]

19. Morgan, C.L.; Chowdhury, S.; Austin, R.B.; Ford, M.A. An evaluation of f1 hybrid winter wheat genotypes produced using a chemical hybridizing agent. *J. Agric. Sci.* **1989**, *112*, 143–149. [CrossRef]

20. Bruns, R.; Peterson, C.J. Yield and stability factors associated with hybrid wheat. *Euphytica* **1998**, *100*, 23–27. [CrossRef]

21. Brears, T.; Hydon, A.; Bingham, J. An assessment of the feasibility of producing F1 and F2 hybrids for the UK. In *Proceedings of the 7th International Wheat Genetics Symposium*; Miller, T.E., Koebner, R.N.D., Eds.; Institute of Plant Science Research, Cambridge Laboratory: Trumpington, Cambridge, UK, 1988; pp. 1157–1162.

22. Jatoi, W.A.; Baloch, M.J.; Khan, N.U.; Munir, M.; Khakwani, A.A.; Vessar, N.F.; Panhwar, S.A.; Gul, S. Heterosis for yield and physiological traits in wheat under water stress conditions. *J. Anim. Plant Sci.* **2014**, *24*, 252–261.

23. Thorwarth, P.; Piepho, H.P.; Zhao, Y.; Ebmeyer, E.; Schacht, J.; Schachschneider, R.; Kazman, E.; Reif, J.C.; Würschum, T.; Longin, C.F.H. Higher grain yield and higher grain protein deviation underline the potential of hybrid wheat for a sustainable agriculture. *Plant Breed.* **2018**, *137*, 326–337. [CrossRef]

24. Paux, E.; Sourdille, P.; Salse, J.; Saintenac, C.; Choulet, F.; Leroy, P.; Korol, A.; Michalak, M.; Kianian, S.; Spielmeyer, W.; et al. A physical map of the 1-gigabase bread wheat chromosome 3B. *Science* **2008**, *322*, 101–104. [CrossRef]

25. Melonek, J.; Duarte, J.; Martin, J.; Beuf, L.; Murigneux, A.; Varenne, P.; Comadran, J.; Specel, S.; Levadoux, S.; Bernath-Levin, K.; et al. The genetic basis of cytoplasmic male sterility and fertility restoration in wheat. *Nat. Commun.* **2021**, *12*, 1–14. [CrossRef]

26. Gupta, P.K.; Balyan, H.S.; Gahlaut, V.; Saripalli, G.; Pal, B.; Basnet, B.R.; Joshi, A.K. Hybrid wheat: Past, present and future. *Theor. Appl. Genet.* **2019**, *132*, 2463–2483. [CrossRef] [PubMed]

27. Mühleisen, J.; Piepho, H.P.; Maurer, H.P.; Longin, C.F.H.; Reif, J.C. Yield stability of hybrids versus lines in wheat, barley, and triticale. *Theor. Appl. Genet.* **2014**, *127*, 309–316. [CrossRef]

28. Askander, H.S.; Salih, M.M.; Altaweel, M.S. heterosis and combining ability for yield and its related traits in bread wheat (*Triticum aestivum* L.). *Plant Cell Biotechnol. Mol. Biol.* **2021**, *22*, 46–53.

29. Oury, F.; Koenig, J.; Bérard, P.; Rousset, M. Une comparaison entre des blés hybrides produits par voie chimique et leurs parents: Niveaux d'hétérosité et élaboration du rendement. *Agronomie* **1990**, *10*, 291–304. [CrossRef]

30. Babar, M.A.; Reynolds, M.P.; Van Ginkel, M.; Klatt, A.R.; Raun, W.R.; Stone, M.L. Spectral reflectance indices as a potential indirect selection criteria for wheat yield under irrigation. *Crop Sci.* **2006**, *46*, 578–588. [CrossRef]

31. Marti, J.; Bort, J.; Slafer, G.A.; Araus, J.L. Can wheat yield be assessed by early measurements of Normalized Difference Vegetation Index? *Ann. Appl. Biol.* **2007**, *150*, 253–257. [CrossRef]

32. Hernandez, J.; Lobos, G.A.; Matus, I.; del Pozo, A.; Silva, P.; Galleguillos, M. Using ridge regression models to estimate grain yield from field spectral data in bread wheat (*Triticum aestivum* L.) grown under three water regimes. *Remote Sens.* **2015**, *7*, 2109–2126. [CrossRef]

33. Gutierrez, M.; Reynolds, M.P.; Raun, W.R.; Stone, M.L.; Klatt, A.R. Spectral water indices for assessing yield in elite bread wheat genotypes under well-irrigated, water-stressed, and high-temperature conditions. *Crop Sci.* **2010**, *50*, 197–214. [CrossRef]

34. Christopher, J.T.; Veyradier, M.; Borrell, A.K.; Harvey, G.; Fletcher, S.; Chenu, K. Phenotyping novel stay-green traits to capture genetic variation in senescence dynamics. *Funct. Plant Biol.* **2014**, *41*, 1035–1048. [CrossRef]

35. Montazeaud, G.; Karatoğma, H.; Özturk, I.; Roumet, P.; Ecarnot, M.; Crossa, J.; Özer, E.; Özdemir, F.; Lopes, M.S. Predicting wheat maturity and stay-green parameters by modeling spectral reflectance measurements and their contribution to grain yield under rainfed conditions. *Field Crop. Res.* **2016**, *196*, 191–198. [CrossRef]

36. Gaju, O.; Allard, V.; Martre, P.; Snape, J.W.; Heumez, E.; LeGouis, J.; Moreau, D.; Bogard, M.; Griffiths, S.; Orford, S.; et al. Identification of traits to improve the nitrogen-use efficiency of wheat genotypes. *Field Crop. Res.* **2011**, *123*, 139–152. [CrossRef]

37. Bogard, M.; Jourdan, M.; Allard, V.; Martre, P.; Perretant, M.R.; Ravel, C.; Heumez, E.; Orford, S.; Snape, J.; Griffiths, S.; et al. Anthesis date mainly explained correlations between post-anthesis leaf senescence, grain yield, and grain protein concentration in a winter wheat population segregating for flowering time QTLs. *J. Exp. Bot.* **2011**, *62*, 3621–3636. [CrossRef]

38. Zhou, W.C.; Zhao, Y.H.; Zou, M.L.; Wang, S.W. Selection and genetic improvement of pollen fertility restorer lines with *Triticum timopheevii* cytoplasm in common wheat. *Plant Breed.* **1999**, *118*, 221–224. [CrossRef]

39. Muqaddasi, Q.H. 15k SNP chip data for spring and winter wheat. *Plant Genom. Phenomics Res. Data Repos.* **2017**. [CrossRef]

40. Purcell, S.; Neale, B.; Todd-Brown, K.; Thomas, L.; Ferreira, M.A.R.; Bender, D.; Maller, J.; Sklar, P.; De Bakker, P.I.W.; Daly, M.J.; et al. PLINK: A tool set for whole-genome association and population-based linkage analyses. *Am. J. Hum. Genet.* **2007**, *81*, 559–575. [CrossRef] [PubMed]

41. R Development Core Team. R: A Language and Environment for Statistical Computing. *R Found. Stat. Comput.* **2021**, *105*. Available online: https://www.R-project.org (accessed on 19 April 2021).

42. De Vries, A.; Ripley, B.D. Ggdendro: Create Dendrograms and Tree Diagrams Using "Ggplot2". 2020, 12. Available online: https://CRAN.R-project.org/package=ggdendro (accessed on 19 April 2021).

43. Lè, S.; Josse, J.; Husson, F. FactoMineR: An R Package for Multivariate Analysis. *J. Stat. Softw.* **2008**, *25*, 1–18. [CrossRef]

44. Wickham, H. Ggplot2: *Elegant Graphics for Data Analysis*; Use R! Springer International Publishing: Cham, Switzerland, 2016; Volume 174, ISBN 978-3-319-24275-0.

45. Rimbert, H.; Darrier, B.; Navarro, J.; Kitt, J.; Choulet, F.; Leveugle, M.; Duarte, J.; Rivière, N.; Eversole, K.; Le Gouis, J.; et al. High throughput SNP discovery and genotyping in hexaploid wheat. *PLoS ONE* **2018**, *13*, e0154756. [CrossRef]
46. Rakshit, S.; Bellundagi, A. Conventional Breeding Techniques in Sorghum. In *Breeding Sorghum for Diverse End Uses*; Elsevier: Amsterdam, The Netherlands, 2019; pp. 77–91.

47. Marvitech MARViN Software. Available online: https://www.marvitech.de/en/home-2/ (accessed on 21 January 2021).

48. Zadoks, J.C.; Chang, T.T.; Konzak, C.F. A decimal code for the growth stages of cereals. *Weed Res.* **1974**, *14*, 415–421. [CrossRef]

49. Lenth, R.V. Emmeans: Estimated Marginal Means, aka Least-Squares Means. 2021. Available online: https://CRAN.R-project.org/package=emmeans (accessed on 19 April 2021).

50. Nyquist, W.E.; Baker, R.J. Estimation of heritability and prediction of selection response in plant populations. *Crit. Rev. Plant. Sci.* **1991**, *10*, 235–322. [CrossRef]

51. Schloerke, B.; Cook, D.; Larmarange, J.; Briatte, F.; Marbach, M.; Thoen, E.; Elberg, A.; Crowley, J. GGally: Extension to "ggplot2". 2021. Available online: https://CRAN.R-project.org/package=GGally (accessed on 19 April 2021).

52. Zhang, Q.; Gao, Y.J.; Yang, S.H.; Ragab, R.A.; Maroof, M.A.S.; Li, Z.B. A diallel analysis of heterosis in elite hybrid rice based on RFLPs and microsatellites. *Theor. Appl. Genet.* **1994**, *89*, 185–192. [CrossRef] [PubMed]

53. Habib Ahmad, I.; Mohammad, F.; Siraj-Ud-Din; Hassan, G.; Gul, R. Evaluation of the heterotic and heterobeltiotic potential of wheat genotypes for improved yield. *Pak. J. Bot.* **2006**, *38*, 1159–1167.

54. Feil, B. The inverse yield-protein relationship in cereals: Possibilities and limitations for genetically improving the grain protein yield. *Trends Agron.* **1997**, *1*, e119.

55. Ahmed, M.A.; Hassan, T.H.A.; Zahran, H.A. Heterosis for seed, oil yield and quality of some different hybrids sunflower. *OCL* **2021**, *28*, 31. [CrossRef]

56. Tollenaar, M.; Ahmadzadeh, A.; Lee, E.A. Physiological basis of heterosis for grain yield in maize. *Crop Sci.* **2004**, *44*, 2086–2094. [CrossRef]

57. Zhao, Y.; Li, Z.; Liu, G.; Jiang, Y.; Maurer, H.P.; Würschum, T.; Mock, H.P.; Matros, A.; Ebmeyer, E.; Schachschneider, R.; et al. Genome-based establishment of a high-yielding heterotic pattern for hybrid wheat breeding. *Proc. Natl. Acad. Sci. USA* **2015**, *112*, 15624–15629. [CrossRef]

58. Oury, F.; Triboï, E.; Bérard, P.; Ollier, J.L.; Rousset, M. Étude des flux de carbone et d'azote chez des blés hybrides et leurs parents, pendant la période de remplissage des grains. *Agronomie* **1995**, *15*, 193–204. [CrossRef]

59. Le Gouis, J.; Beghin, D.; Heumez, E.; Pluchard, P. Diallel analysis of winter wheat at two nitrogen levels. *Crop Sci.* **2002**, *42*, 1129–1134. [CrossRef]

60. Oury, F.; Brabant, P.; Pluchard, P.; Bérard, P.; Rousset, M.; Gourdone, J.; Corti, H.; Meunier, J.; Blake, A. Étude multilocale de blés hybrides: Niveaux d'hétérosérosis et élaboration du rendement. *Agron. EDP Sci.* **1990**, *10*, 735–748. [CrossRef]

61. Monaghan, J.M.; Snape, J.W.; Chojecki, A.J.S.; Kettlewell, P.S. The use of grain protein deviation for identifying wheat cultivars with high grain protein concentration and yield. *Euphytica* **2001**, *122*, 309–317. [CrossRef]

62. Sultana, N.; Islam, S.; Juhasz, A.; Ma, W. Wheat leaf senescence and its regulatory gene network. *Crop J.* **2021**. [CrossRef]

63. Bogard, M.; Allard, V.; Brancourt-Hulmel, M.; Heumez, E.; MacHet, J.M.; Jeuffroy, M.H.; Gate, P.; Martre, P.; Le Gouis, J. Deviation from the grain protein concentration-grain yield negative relationship is highly correlated to post-anthesis N uptake in winter wheat. *J. Exp. Bot.* **2010**, *61*, 4303–4312. [CrossRef]

64. Milan, M.; Fogliatto, S.; Blandino, M.; Vidotto, F. Are wheat hybrids more affected by weed competition than conventional cultivars? *Agronomy* **2020**, *10*, 526. [CrossRef]

65. Pradhan, S.; Sehgal, V.K.; Bandyopadhyay, K.K.; Panigrahi, P.; Parihar, C.M.; Jat, S.L. Radiation interception, extinction coefficient and use efficiency of wheat crop at various irrigation and nitrogen levels in a semi-arid location. *Indian J. Plant Physiol.* **2018**, *23*, 416–425. [CrossRef]

66. Tan, C.W.; Zhang, P.P.; Zhou, X.X.; Wang, Z.X.; Xu, Z.Q.; Mao, W.; Li, W.X.; Huo, Z.Y.; Guo, W.S.; Yun, F. Quantitative monitoring of leaf area index in wheat of different plant types by integrating NDVI and Beer-Lambert law. *Sci. Rep.* **2020**, *10*, 1–10. [CrossRef] [PubMed]

67. Pinto, R.S.; Reynolds, M.P.; Mathews, K.L.; McIntyre, C.L.; Olivares-Villegas, J.J.; Chapman, S.C. Heat and drought adaptive QTL in a wheat population designed to minimize confounding agronomic effects. *Theor. Appl. Genet.* **2010**, *121*, 1001–1021. [CrossRef]

68. Olivares-Villegas, J.; Reynolds, M.; William, H.; Mcdonald, G.; Ribaut, J. Drought Adaptation Attributes and Associated Molecular Markers via BSA in the Seri / Babax Hexaploid Wheat (*Triticum aestivum*, L.) Population. *Funct. Plant. Biol.* **2008**, *34*, 1–3.

69. Hitz, K.; Clark, A.J.; Van Sanford, D.A. Identifying nitrogen-use efficient soft red winter wheat lines in high and low nitrogen environments. *Field Crop. Res.* **2017**, *200*, 1–9. [CrossRef]

70. Bennett, D.; Reynolds, M.; Mullan, D.; Izanloo, A.; Kuchel, H.; Langridge, P.; Schnurbusch, T. Detection of two major grain yield QTL in bread wheat (*Triticum aestivum* L.) under heat, drought and high yield potential environments. *Theor. Appl. Genet.* **2012**, *125*, 1473–1485. [CrossRef]

71. Christopher, M.; Paccapelo, V.; Kelly, A.; Macdonald, B.; Hickey, L.; Richard, C.; Verbyla, A.; Chenu, K.; Borrell, A.; Amin, A.; et al. QTL identified for stay-green in a multi-reference nested association mapping population of wheat exhibit context dependent expression and parent-specific alleles. *Field Crop. Res.* **2021**, *270*, 108181. [CrossRef]

72. Yashavanthakumar, K.J.; Baviskar, S.; Navathe, S. Impact of heat and drought stress on phenological development and yield in bread wheat. *Plant Physiol. Rep.* **2021**, *26*, 357–367. [CrossRef]

73. Oury, F.; Brabant, P.; Pluchard, P.; Bérard, P.; Rousset, M. Une étude sur la supériorité des blés hybrides au niveau des capacités de remplissage du grain: Résultats d'une expérimentation multilocale. *Agronomie* **1993**, *13*, 381–393. [CrossRef]
74. Lanning, S.P.; Kephart, K.; Carlson, G.R.; Eckhoff, J.E.; Stougaard, R.N.; Wichman, D.M.; Martin, J.M.; Talbert, L.E. Climatic change and agronomic performance of hard red spring wheat from 1950 to 2007. *Crop Sci.* **2010**, *50*, 835–841. [CrossRef]

75. Lyra, D.H.; Galli, G.; Alves, F.C.; Granato, Í.S.C.; Vidotti, M.S.; e Sousa, M.B.; Morosini, J.S.; Crossa, J.; Fritsche-Neto, R. Modeling copy number variation in the genomic prediction of maize hybrids. *Theor. Appl. Genet.* **2019**, *132*, 273–288. [CrossRef] [PubMed]

76. Sun, S.; Zhou, Y.; Chen, J.; Shi, J.; Zhao, H.; Zhao, H.; Song, W.; Zhang, M.; Cui, Y.; Dong, X.; et al. Extensive intraspecific gene order and gene structural variations between Mo17 and other maize genomes. *Nat. Genet.* **2018**, *50*, 1289–1295. [CrossRef]

77. Tyrka, M.; Mokrzycka, M.; Bakera, B.; Tyrka, D.; Szeliga, M.; Stojałowski, S.; Matysik, P.; Rokicki, M.; Rakoczy-Trojanowska, M.; Krajewski, P. Evaluation of genetic structure in European wheat cultivars and advanced breeding lines using high-density genotyping-by-sequencing approach. *BMC Genom.* **2021**, *22*, 1–17. [CrossRef]