# Integration of genomic selection into winter-type bread wheat breeding schemes: a simulation pipeline including economic constraints

Sarah Ben-Sadoun, Aline Fugeray-Scarbel, Jérôme Auzanneau, Gilles Charmet, Stephane Lemarié, Sophie Bouchet

► To cite this version:

Sarah Ben-Sadoun, Aline Fugeray-Scarbel, Jérôme Auzanneau, Gilles Charmet, Stephane Lemarié, et al. Integration of genomic selection into winter-type bread wheat breeding schemes: a simulation pipeline including economic constraints. Crop Breeding, Genetics and Genomics, 2021, 3, 10.20900/cbgg20210008. hal-03359314

HAL Id: hal-03359314

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

Submitted on 30 Sep 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

# Integration of Genomic Selection into Winter-Type Bread Wheat Breeding Schemes: A Simulation Pipeline including Economic Constraints

Sarah Ben-Sadoun¹, Aline Fugeray-Scarbel², Jérôme Auzanneau³,
Gilles Charmet¹, Stéphane Lemarié², Sophie Bouchet¹,*

¹ Genetics Diversity and Ecophysiology of Cereals, UMR1095, INRAE-UCA, 5 chemin de Beaulieu 63000 Clermont-Ferrand, France

² University Grenoble Alpes, INRAE, CNRS, Grenoble INP, GAEL, 38000 Grenoble, France

³ Agri-Obtentions, Ferme de Gauvilliers, 78660 Orsonville, France

* Correspondence: Sophie Bouchet, Email: sophie.bouchet@inrae.fr; Tel.: +33-044-376-1509.

## ABSTRACT

*Background:* Relatively little effort has been made yet to optimize the allocation of resources when using genomic predictions to maximize the return to investment in terms of genetic gain per unit of time and cost.

*Methods:* We built a simulation pipeline in the R environment designed to become a decision tool to help breeders adjusting breeding schemes, according to their either short or long-term objectives. We used it to explore different scenarios in order to investigate under which conditions (at what step of the breeding program) genomic predictions could improve genetic gain. For a given budget per cycle, we compared 36 scenarios, varying strategies (phenotypic selection PS or genomic selection + phenotypic selection: GPS), trait heritability, relative selection rate at two key steps and genotyping cost. With GPS strategy, we also optimized mating using genomic predictions. The reference population is a 20 years historical data set from the INRAE-Agri-Obtentions bread wheat breeding program. We simulated 3 cycles of 5 years selection.

*Results:* We showed that GPS selection using mating optimization significantly improved genetic gain for all scenarios while GPS without mating optimization and PS had similar efficiency in terms of genetic gain. Our results also highlighted that the loss of genetic diversity over successive cycles was faster using GPS strategies. Those were more efficient to increase favourable allele frequency, rare alleles in particular.

**KEYWORDS:** breeding scheme; genomic prediction; mating optimization; economic constraint; bread wheat

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

## Open Access

Received: 01 December 2020

Accepted: 23 September 2021

Published: 29 September 2021

Copyright © 2021 by the author(s). Licensee Hapres, London, United Kingdom. This is an open access article distributed under the terms and conditions of [Creative Commons Attribution 4.0 International License](https://creativecommons.org/licenses/by/4.0/).
ABBREVIATIONS

DH, double haploid; GEBV, genomic estimated breeding values; GPS, genomic and phenotypic selection; GS, genomic selection; He, expected heterozygosity; INRAE, Institut National de la Recherche en Agriculture, Alimentation et Environnement; OCS, optimal contribution selection; OHV, optimal haploid value; PS, phenotypic selection; QTL, quantitative trait locus; SNK, Student-Newman-Keuls; SNP, single nucleotide polymorphism; SSD, single seed descent; TBV, true breeding value; UC, usefulness criterion; UCPC, usefulness criterion parental contribution

INTRODUCTION

The objective of bread wheat breeding programs is to develop new varieties that outperform current varieties in terms of yield, adaptation, resistance to biotic and abiotic stresses, and/or end use qualities. A great challenge in plant breeding is to improve the genetic gain per unit of time for a given investment. To meet this goal, the optimization of resource allocation appears to be a key point. In addition, breeders must find a compromise between short-term genetic gain and the conservation of genetic diversity within their germplasm in order to guarantee long-term genetic gain [1].

Furthermore, the exponential decrease of genotyping costs, improvement of computing tools, data storage capacities and algorithms' efficiency, and the development of new statistical methods have led to the development of a new powerful approach to optimize breeding schemes: genomic selection (GS). GS is a marker-based selection method that uses thousands to millions of molecular markers evenly spread along the genome to predict the genetic value of candidates to selection [2,3]. According to the breeder's equation [4], GS could improve genetic gain by (i) accelerating genetic gain by shortening the breeding cycle, replacing field evaluation with genotyping at juvenile stage [5] for long-cycle plants like trees in particular, implement recurrent selection with rapid cycles in cereals with up to 3 steps of GS per year in greenhouse for maize for instance [6,7], (ii) decreasing the budget allocated to field evaluation by optimizing the number of genotypes and replicates per environment in the experimental design [8] and by this way increasing the size of the breeding program (number of crosses or progenies per cross), (iii) increasing genetic variance by optimizing mating to cumulate favourable alleles [9–12], (iv) increasing accuracy of selection. For wheat in particular, theoretical estimates of genetic gain showed that GS can accelerate recurrent selection using off-season field for spring wheat or rapid cycles in greenhouses for winter wheat. The gain was higher when applying GS in F2 compared to F3 or F4 [13]. We can postulate that optimized mating should increase the gain further.

Several factors influence the accuracy of GS. These factors include trait
architecture and trait heritability [14], training set size and composition
[15,16], i.e., congruency between the allelic composition represented in the training population to estimate marker effects and the allelic composition of the candidates whose performance is to be predicted [17–23], marker density [23–26], and statistical model for estimation of marker effects [27]. For wheat in particular, it was shown that the size of the training population should be 50 if individuals are full-sibs, 100 if they are half-sibs and 1000 if they are unrelated [13].

Most of previous researches on GS in plant breeding focused on the prediction accuracy of unphenotyped lines. In bread wheat, studies evaluated the prediction accuracy of grain yield [28–32], traits linked to bread making quality [33–38], or disease [39–42]. At the breeding program scale, some simulation works showed an interest of GS compared to classical phenotypic evaluation in terms of genetic gain. For example, a study [43] showed that GS accuracies were high enough (GS accuracy twice as high as PS accuracy) to achieve greater gain from selection per unit time compared with phenotypic selection for bread wheat adult stem rust resistance. However, the two strategies they compared represented different budgets. In contrast, another study [44] simulated and compared several hybrid and line wheat breeding programs using GS for a fixed budget. This study showed that GS could be advantageous in terms of genetic gain for line but even more for hybrid breeding in wheat. Furthermore, the efficiency of PS and GS for grain yield assuming a single selection cycle and a given budget were compared using a biparental population of maize double haploid (DH) lines, as discussed by [45]. They showed that with large Genotype × Environment interactions and under limited resources, it was beneficial to use an index combining PS and GS to maximize genetic gain. They also noticed that DH price was a limiting factor for large genetic gain. But none of those studies evaluated the interest of mating optimization using genomic predictions. In maize, several studies showed its interest for long term genetic gain, in pre-breeding programs in particular [46–50].

Simulation studies actually enable the comparison of a wide range of scenarios that would not be possible to test experimentally. They also allow to evaluate an unlimited number of selection cycles (long-term selection) in a short amount of time with a cost limited to data storage and processing.

In this study we compared the genetic gain and the evolution of genetic diversity in two types of simulated breeding schemes: one called Phenotypic Selection (PS) with two steps of selection based on field trials, and one called Genomic and Phenotypic Selection (GPS) that combines a first step of genomic selection and a second step of selection based on field trials. We also evaluated the interest of mating optimization (GPSopt). We explored different scenarios in order to investigate under which conditions GPS and GPSopt were more cost-effective than PS. We compared scenarios for a given budget covering all breeding costs.
## MATERIALS AND METHODS

We developed a R pipeline to simulate and compare winter wheat breeding programs in terms of genetic gain and genetic diversity evolution for a given budget (Figure 1). The scripts of the pipeline are available using the following link: https://forgemia.inra.fr/umr-gdec/gps.

<figure><img src="image_3.png" /><figcaption>**Figure 1.** PS and GPS breeding schemes. PS: Phenotypic Selection. GPS: Genomic and Phenotypic Selection. $N_p$ and $N_c$: number of parents and crosses respectively. $N_2$, $N_3$ and $N_4$: number of progenies at the beginning of steps 2, 3 and 4 respectively. $\alpha_2$, $\alpha_3$ and $\alpha_4$: selection rate on steps 2, 3 and 4 respectively. $\alpha_{34}$: global selection rate on steps 3 and 4.</figcaption></figure>

## Data Set

The pipeline was tested using a real breeding population of 757 winter-type bread wheat lines developed by the Institut National de la Recherche en Agriculture, Alimentation et Environnement (INRAE, formerly INRA) and its subsidiary company Agri-Obtentions. These lines were selected between 2000 and 2013. Each line was genotyped using a 280K SNP array [51]. The datasets with the genotyping data are available in the INRAE Dataverse repository (https://data.inrae.fr/dataset.xhtml?persistentId=doi:10.15454/M8SAYH).
In order to limit computing time in simulation, we used a subset composed of 12,119 SNP evenly spread along the genetic map. Such marker number was previously shown as giving similar predictive ability as the full marker set [52].

## Trait Simulation

We considered that 100 QTLs control the traits of interest. We simulated 20 traits using random positions of the 100 QTLs, i.e., 20 random samples of 100 SNPs were assigned as QTL, marker effects drawn from a gaussian distribution:

$$ \beta \sim N(0,1) $$

The favourable allele was attributed at random to one of the two SNP alleles, so that coupling and repulsion associations also occur at random. The entry-mean heritability (h²) was set to either 0.2, 0.4 or 0.7. Phenotypic values of lines were obtained by adding a normally distributed noise to the genotypic values. In our study, the residual variance was 80% (h² = 0.2), 60% (h² = 0.4) or 30% (h² = 0.7) of phenotypic variance. We simulated phenotypes as:

$$ Y = TBV + N(0, \sqrt{\frac{(1-h^2) \cdot var(TBV)}{h^2}}) = Q\beta + N(0, \sqrt{\frac{(1-h^2) \cdot var(X\beta)}{h^2}}), \quad (1) $$

with Q the genotyping matrix of 100 QTLs and β the vector of marker effects.

## Simulation of the Breeding Programs

We compared two types of breeding schemes (Figure 1): one called Phenotypic Selection (PS) with two steps of selection based on field trials, and one called Genomic and Phenotypic Selection (GPS) that combines a first step of genomic selection and a second step of selection based on field trials. For the sake of simplicity, both breeding programs were designed using doubled haploids (DHs), instead of successive selfing, to reduce the time required for inbred development. We counted three years for cross, DH productions from F1 and seed multiplication, one year of PS (or GPS) selection and a last year of phenotypic selection, to fit breeder's requirement of having real data to apply candidates in official registration trials. For both PS and GPS approaches, a breeding cycle lasts five years. We simulated three successive cycles. For the first cycle, the training population is composed of 757 Nref lines that have been genotyped and virtually “phenotyped” for one trait (20 times). Then, we extend the training population at each cycle with the N4 individuals that are genotyped and “virtually” phenotyped to update the genomic prediction equation.
## Genomic Prediction Models

We start the simulation with $N_{ref}$ lines (757 in our example) that have been genotyped and virtually “phenotyped” for one quantitative trait, twenty times to avoid any bias in QTL positions. For GPS, the first vector of marker effects $\hat{\beta}_0$ is calculated using a ridge regression on those phenotypes. The database of phenotypes is incremented at each phenotyping step and the vector of marker effects $\hat{\beta}_k$ is updated at each cycle k.

$$ \hat{\beta} = (X'X + \lambda I)^{-1} X'y, \quad (2) $$

with y the vector of phenotypes, X the matrix of the 12,019 markers after removing the 100 markers sampled to simulate QTLs, $\lambda$ chosen to make X'X non-singular, using the R package rrBLUP.

We evaluated the accuracy of genomic predictions using the correlation between TBV and GEBV of lines at the third step.

At the first step of each cycle, $N_c$ crosses are produced between the $N_p$ best individuals according to phenotyping results at step 4 of previous cycle (or $N_{ref}$ individuals for the first cycle). In PS scheme, crosses are obtained by randomly mating these $N_p$ best individuals. Note that we could have chosen to cross the two best individuals according to phenotypes, or the first with the second, the second with the third. But none of those strategy is realistic. In a real breeding program, breeders would select parents that complement for different traits. More realistic mating scenarios in PS strategy will be tested when a multi-trait strategy will be implemented in the pipeline. For this paper, we decided to compare PS and GPS strategy using the same random mating strategy among the $N_p$ best individuals, and a GPS strategy with random or optimized mating. This second strategy is called GPSopt and optimizes the complementarity between parental alleles [53]. We calculate the value of a cross as the value of the individual that would have inherited the best chromosomes of its parents.

$$ UC_{ij} = \sum_{c=1}^{nchr} \max(GEBV_{ci}, GEBV_{cj}), \quad (3) $$

with $UC_{ij}$ the usefulness criterion of the cross between the i-th and the j-th parents, $i \neq j$, c the chromosome number, $GEBV = X\hat{\beta}$. This UC assumes chromosome being inherited without recombination, which is true for half of the chromatids in a single meiosis, as it is the case for doubled haploids from F1. The best cross is the cross that could produce the best possible gamete if the progeny size was unlimited. This is unrealistic but it was shown to be an acceptable approximation of cross value when using DH in wheat programs by [53].

To simulate progeny, each chromosome is either parental or recombined. If recombined, the number of cross-overs on each chromosome is sampled from a Poisson distribution, and cross-overs positions are distributed randomly along the genetic map [51].
At step 2, N₂ progeny are multiplied to implement one trial at the following step. We considered that the selection at this step was made on agronomic visual traits, excluding plants with major defaults, disease sensibility in particular. Since it cannot be focused on targeted traits (yield), we considered it had no impact on the distribution of genetic value of the N₃ selected progeny. So, the selection was considered random for the targeted trait.

At step 3, N₃ progeny are evaluated in one trial.

$$N_3 = N_2\alpha_2, \qquad (4)$$

with α₂ the selection rate at step 2 (α₂ = 0.2 in this study).

In PS scheme, these individuals are evaluated in the field and we consider a multiplication cost corresponding to five trials for next step. Selection is based on phenotypic data with a rate of selection α₃. In GPS scheme, there is no trial, the selection is based on GEBVs calculated as the cross product between the vector of marker effects and the matrix of genotypes of progenies, excluding markers to which QTLs were assigned. Only N₄ (N₃α₃) progeny are multiplied in nursery.

At step 4, N₄ progeny are evaluated in five trials (different environments, e.g., different locations).

$$N_4 = N_3\alpha_3 \qquad (5)$$

The last step is a phenotypic selection of Nₚ parents for the next cycle in both PS and GPS schemes with a selection rate α₄.

$$N_p = N_4\alpha_4 \qquad (6)$$

Instead of defining α₃ and α₄ independently, we fixed a relative selection rate (called λ) between steps 4 and 3 as follow:

$$\alpha_3 = \alpha_{34}^{\lambda} \qquad (7)$$

and

$$\alpha_4 = \alpha_{34}^{1-\lambda}, \qquad (8)$$

with:

$$\alpha_{34} = \alpha_3\alpha_4 = \frac{N_p}{N_2\alpha_2}, \qquad (9)$$

If λ = 0.5, α₃ = α₄. If λ > 0.5, α₃ < α₄. If λ < 0.5, α₃ > α₄.

## Costs Modelling

We defined Cₓ the cost of each operation X of the breeding program (Supplementary Table S1). The cost of an operation sometimes varies depending on the step where the operation was done. The step is then specified as an exponent. For example, we assumed that field evaluations were realized in one trial at step 2, five trials at step 3 and ten trials at step 4 (for registration), which explains the different costs of seed multiplication and evaluation. Note that those costs that include DNA extraction, genotyping and biometrician salaries are specific to INRAE.
Genotyping cost may decrease in the future. Some other genotyping platforms are more flexible and allow to genotype with high density the parents and low density the progeny that can be imputed [13].

Total cost (CT) of the PS scheme was defined as follow:

$$CT_{PS}(N_c, N_2, N_3, N_4) = K[N_c C_c + N_2(C_{DH} + C_M^2) + N_3(C_M^3 + C_P^3) + N_4(C_M^4 + C_P^4)], \quad (10)$$

with K the number of cycle (K = 3 in this study) and N₂ the progeny number at the beginning of step 2. In our study, we supposed that the number of progenies was the same for all crosses. Note however that the pipeline offers the possibility to make it proportional to cross value (Usefulness Criterion). Total cost depends on 4 parameters concerning progeny size at each step. The number of possible combinations for a same budget is very large. In order to simplify comparisons, as progeny sizes are dependant of each other for a fixed budget, we defined three parameters (N_c, N_p and λ) that are fixed for one scenario. The total cost CT and progeny sizes N₃ and N₄ can be expressed as functions of these parameters and N₂.

Thanks to previous equations, we have:

$$N_4 = \frac{N_p}{\alpha_4} = \frac{N_p}{\alpha^{1-\lambda}} = \frac{N_p}{(\frac{N_p}{N_2\alpha_2})^{1-\lambda}} = N_p^{\lambda}(N_2\alpha_2)^{1-\lambda}, \quad (11)$$

Therefore, the total cost of the PS scheme can be defined as follows:

$$CT_{PS}(N_p, N_c, \lambda, N_2) = K[N_c C_c + N_2(C_{DH} + C_M^2 + \alpha_2(C_M^3 + C_P^3)) + N_p^{\lambda}(N_2\alpha_2)^{1-\lambda}(C_M^4 + C_P^4)], \quad (12)$$

For a given total cost and set of parameters (N_c, N_p and λ) we search for the value of N₂ that solve this equation.

To define the total cost of the GPS scheme (CT<sub>GPS</sub>), we introduced the genotyping cost of the reference population composed of N<sub>ref</sub> lines. Note that cost of field evaluation at the third step is replaced by the cost of genotyping and that only N₄ lines that are selected at step 3 were multiplied for GPS (instead of N₃ for PS):

$$CT_{GPS}(N_p, N_c, \lambda, N_2) = N_{ref} C_G + K[N_c C_c + N_2(C_{DH} + C_M^2 + \alpha_2 C_G) + N_p^{\lambda}(N_2\alpha_2)^{1-\lambda}(C_M^3 + C_P^4)], \quad (13)$$

As genotyping cost is lower than phenotyping cost, the number of progenies is larger in GPS strategy. As we fixed the total cost, either the number of crosses or the number of progenies per cross will be different between PS and GPS schemes. For GPS and GPSoft schemes, we tested scenarios with the same number of crosses and a different number of progenies n2 per cross (called GPS.n2 and GPSoft.n2) and scenarios with the same number of progenies per cross and a different number of crosses (called GPS.Nc and GPSoft.Nc).

## Simulation of Different Scenarios

We evaluated the impact of several parameters on the final genetic gain. To do so, we simulated breeding programs with two total costs for 15 years
(three cycles of five years; CT ∈ {22.5 M€, 45 M€}, i.e., average CT per year of 1.5 or 3 M€), two genotyping cost (CG ∈ {37€, 10€}), three relative selection rate λ (λ ∈ {0.25,0.5,0.75}), and three levels of heritability of the trait (h² ∈ {0.2,0.4,0.7}). The description of each scenario is available in Supplementary Table S2.

It led to the evaluation of 36 scenarios for five different strategies (PS, GPS.n2: fixed number of crosses, GPS.Nc: fixed number of progenies per cross, GPSopt.n2: optimized mating design with a fixed number of crosses, GPSopt.Nc: optimized mating design with a fixed number of progenies per cross). For each strategy / scenario, we tested 20 simulated traits (100 QTLs randomly sampled for each simulated trait) to evaluate the variance due to different QTL positions. For each combination of strategy, scenario and trait, we ran the algorithm 10 times to estimate the variance due to mendelian sampling only. For each simulation, we performed three cycles of five years.

The algorithm used to run the simulations required several input data: including the simulated strategy (PS, GPS.Nc, GPS.n2, GPSopt.Nc, GPSopt.n2) trait and replication, a matrix with chromosome number and genetic position in columns 1 and 2 respectively, a matrix of genotypes with N<sub>ref</sub> rows, a vector of phenotypes of N<sub>ref</sub> length, one vector including true marker effects for the simulated trait, one scenario (combination of CT, C<sub>G</sub>, λ and h<sup>2</sup> values), the cost of each operation and the selection rate on step 2. The algorithm is illustrated in Supplementary Figure S1.

Evolution of the Genetic Gain

For each strategy and scenario, at the end of each cycle, we computed the genetic gain as the difference between the mean of the true breeding values (TBVs) of the N<sub>p</sub> (200) best lines and the reference (initial) data set. The average TBV of the 200 best parental lines was 4.85. This average TBV was higher than the average TBV of the parental population which was 0.33. We performed an analysis of variance (ANOVA) to compare the impact of the input parameters on the genetic gain.

We tested a model without interactions between parameters (14), models that take into account interactions between the strategy and the simulated trait (15), the strategy and the total cost of the breeding scheme (16), the strategy and the relative selection rate (λ) (17), the strategy and the trait heritability (18).

$$
g_{ijklm} = \mu_0 + S_i + \lambda_j + h^2_k + CT_l + T_m + \varepsilon_{ijklm}, \quad (14)
$$

$$
g_{ijklm} = \mu_0 + S_i + \lambda_j + h^2_k + CT_l + T_m + S_i \times T_m + \varepsilon_{ijklm}, \quad (15)
$$

$$
g_{ijklm} = \mu_0 + S_i + \lambda_j + h^2_k + CT_l + T_m + S_i \times CT_l + \varepsilon_{ijklm}, \quad (16)
$$

$$
g_{ijklm} = \mu_0 + S_i + \lambda_j + h^2_k + CT_l + T_m + S_i \times \lambda_j + \varepsilon_{ijklm}, \quad (17)
$$

$$
g_{ijklm} = \mu_0 + S_i + \lambda_j + h^2_k + CT_l + T_m + S_i \times h^2_k + \varepsilon_{ijklm}, \quad (18)
$$
with gijklmn is the genetic gain of the scenario with the i-th strategy (S), the j-th value of λ, the k-th value of h², the l-th value of CT, and for the m-th simulated trait (T).

For GPS strategies, we also studied the impact of genotyping cost on genetic gain:

$$g_{ijklmn} = \mu_0 + S_i + \lambda_j + h^2_k + CT_l + T_m + CG_n + \varepsilon_{ijklmn} \quad (19)$$

with CGn the n-th value of CG. F tests were considered significant at α < 0.05.

For the different strategies, genetic gains for pairwise scenarios were compared using the Student-Newman-Keuls (SNK) test from the R library Agricolae. Means were judged significantly different when P-values < 0.05.

## Evolution of Genetic Diversity

For each strategy and scenario, we analysed the evolution of genetic diversity over successive cycles. To do so, we measured the percentage of alleles present in both the reference population and the 200 best progenies of each cycle. We estimated the significance of the various factors on genetic diversity with the same models used for genetic gain above. We also measured the difference between the expected heterozygosity (He) in the initial population and at the third cycle for (i) all the markers but QTLs, (ii) markers assigned to QTLs and (iii) markers located at the vicinity of QTLs (positioned at less than 1, 5, 10 or 15 cM from a QTL), using the following formula:

$$He = \frac{\sum_{1}^{L}(1-\sum_{1}^{n} p_{i}^{2})}{L}, \quad (20)$$

with L the number of loci, n the number of alleles (2 in our case) for each locus, and pᵢ the frequency of the allele i.

For the different strategies, He differences at the end of the last cycle were compared using the SNK test.

We also compared the number of cumulated favourable alleles (difference between the third cycle and the initial population) between strategies and scenarios.

## Parental Contribution

For GPSopt strategies, a usefulness criterion was used to choose crosses that maximize the probability to get lines that cumulate a maximum of favourable alleles among the Np selected lines. Note that in our example, we did not constrain the contribution of any parents, but it is possible to fix a contribution threshold using this pipeline. For GPS and PS strategies, crosses were chosen at random among the Np selected lines.

We compared the distribution of the contributions of reference lines between strategies and scenarios based on the pedigree of the Np lines selected at the end of each cycle.

We defined the contribution of a reference line as the number of selected progenies presenting them in their pedigree divided by the total number of parental lines that contributed to the global progeny. Note that
all the reference lines or Np lines selected at the end of previous cycle did not contribute to crosses. Then we calculated the number of parents that contributed to 25% and 75% of the progeny. In addition, we used the Shannon index to evaluate the diversity of the parental lines that contributed using the following formula:

$$H' = - \sum_{i=1}^{N_{ref}} p_i \times \ln(p_i), \qquad (21)$$

with Nref the number of initial parental lines, and pi = ni/N with ni the contribution of parent i to the Np progenies and N the total number of lines that contributed to progeny. Note that the more unequal the parental contributions, the smaller the corresponding Shannon index.

In order to check if the distribution of favourable alleles was different among parents selected by UC or at random, we computed for each cross that contributed to selected lines the number of shared and specific favourable alleles among the two parents. We calculated it for each cycle.

# RESULTS

## Progeny Size under Different Strategies and Scenarios

When the budget CT was doubled, the number of progenies was 1.9 times larger at each step on average (for given values of Cg and λ). When the cost to genotype one line (Cg) decreased from 37 to 10 euros, the number of progenies at each step was only 1.1 times larger at each step on average (for given values of CT and λ). Note that for given values of CT and Cg, the higher the value of λ, the higher the population size N3 in step 3, the smaller the population size in step 4 N4 and the higher N3 + N4.

Supplementary Table S3 reports the number of progenies at each step, for each strategy (PS, GPS.Nc, GPS.n2, GPSopt.Nc, GPSopt.n2) and scenario (combination of CT, Cg and λ).

## Cost Distribution between Operations

We compiled the percentage of the total budget allocated to each operation (crosses, DH, multiplication, field experiment and genotyping) under different strategies (PS, GPS.Nc, GPS.n2, GPSopt.Nc, GPSopt.n2) and scenarios (combination of CT, Cg, λ). Results were similar for GPS.Nc and GPS.n2, and for GPSopt.Nc and GPSopt.n2, for that reason we only show results for PS, GPS.n2 and GPSopt.n2 in Table 1. For each strategy, the production of DHs was the major source of expenses. Indeed, for PS strategy between 30.8% and 46.5% of the global budget was allocated to DH production. This percentage is even higher for GPS schemes (between 37.3% and 72.6%). Multiplication steps required on average between 20% and 25% of the global budget. Since one field evaluation was more expensive at step 4 than at step 3 (more replicated plots), scenarios with a higher number of lines in the last step (λ < 0.5) used more budget for field evaluations.
For GPS strategies, genotyping required between 7.6% and 13.1% of the global budget when the cost of genotyping was 37€. It was obviously lower when the cost of genotyping was 10 € (ranging from 2.2% to 3.9%). But note that the number of tested lines only slightly increased (x1.1).

<table><caption>**Table 1.** Percentage of total budget allocated to each operation.</caption><thead><tr><th rowspan="2">Breeding scheme</th><th rowspan="2">Average annual CT</th><th rowspan="2">CG</th><th rowspan="2">λ</th><th colspan="4">% percentage of the budget allocated to</th><th rowspan="2">Genotyping</th></tr><tr><th>crosses</th><th>DH</th><th>multiplication</th><th>field experiment</th></tr></thead><tbody><tr><td rowspan="6">PS</td><td rowspan="3">1.5M€</td><td rowspan="3">---</td><td>0.25</td><td>0.5</td><td>30.8</td><td>22.5</td><td>46.2</td><td rowspan="6">0</td></tr><tr><td>0.5</td><td>0.5</td><td>40.8</td><td>24.5</td><td>34.2</td></tr><tr><td>0.75</td><td>0.5</td><td>45.4</td><td>25.3</td><td>28.7</td></tr><tr><td rowspan="3">3M€</td><td rowspan="3">---</td><td>0.25</td><td>0.3</td><td>32.9</td><td>23.0</td><td>43.8</td></tr><tr><td>0.5</td><td>0.3</td><td>42.8</td><td>24.9</td><td>32.0</td></tr><tr><td>0.75</td><td>0.3</td><td>46.5</td><td>25.6</td><td>27.6</td></tr><tr><td rowspan="10">GPS and GPSopt</td><td rowspan="4">1.5M€</td><td rowspan="4">10€</td><td>0.25</td><td>0.5</td><td>39.9</td><td>21.5</td><td>35.9</td><td>2.2</td></tr><tr><td>0.5</td><td>0.5</td><td>60.6</td><td>20.8</td><td>14.8</td><td>3.3</td></tr><tr><td>0.75</td><td>0.5</td><td>70.6</td><td>20.5</td><td>4.6</td><td>3.8</td></tr><tr><td>0.25</td><td>0.5</td><td>37.3</td><td>20.4</td><td>34.2</td><td>7.6</td></tr><tr><td rowspan="2">37€</td><td rowspan="2">0.5</td><td>0.5</td><td>55.1</td><td>19.1</td><td>14.1</td><td>11.2</td></tr><tr><td>0.75</td><td>0.5</td><td>63.6</td><td>18.6</td><td>4.5</td><td>12.8</td></tr><tr><td rowspan="4">3M€</td><td rowspan="4">10€</td><td>0.25</td><td>0.3</td><td>43.6</td><td>21.4</td><td>32.3</td><td>2.4</td></tr><tr><td>0.5</td><td>0.3</td><td>64.7</td><td>20.7</td><td>10.8</td><td>3.5</td></tr><tr><td>0.75</td><td>0.3</td><td>72.6</td><td>20.4</td><td>2.8</td><td>3.9</td></tr><tr><td>0.25</td><td>0.3</td><td>40.6</td><td>20.2</td><td>30.7</td><td>8.2</td></tr><tr><td rowspan="2">37€</td><td rowspan="2">0.5</td><td>0.5</td><td>58.7</td><td>18.9</td><td>10.3</td><td>11.8</td></tr><tr><td>0.75</td><td>0.3</td><td>65.4</td><td>18.5</td><td>2.7</td><td>13.1</td></tr></tbody></table>

DH: Doubled haploids. CT: Total cost. CG: Genotyping cost. λ: Relative selection rate. Results were similar for GPS.N<sub>c</sub> and GPS.n2, and for GPSopt.N<sub>c</sub> and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt).

## Contribution of Input Parameters to Final Genetic Gain

We evaluated the contribution of input parameters (i.e., strategy, total cost, trait heritability, QTL sampling, and λ) to genetic gain for the 200 lines selected at last step.

We found that each of the five factors had a significant effect on the final genetic gain (Supplementary Table S4). Trait heritability contributed the most. As expected, scenarios with higher heritabilities led to higher genetic gain (the average genetic gain was 13.32, 16.36 or 18.90 when the heritability was 0.2,0.4 or 0.7, respectively). The strategy was the second most significant parameter. The genetic gain for GPSopt.N<sub>c</sub> and GPSopt.n2 strategies was larger than the genetic gain for strategies using random mating (the average genetic gain was 18.0 for GPSopt, 15.0 for GPS and 15.3 for PS). The strategy and the trait heritability accounted for 16.7% and 40.1%
of the sum of squares, respectively. Note that doubling the total cost had a rather low impact on the final genetic gain (+ 1.07 on average). In addition, λ and QTL sampling had a low impact on the final genetic gain. We did the same analysis with the genetic gain obtained for various selection rates at last step (for the 10, 50 and 100 best lines) and we obtained the same conclusions (Supplementary Table S5). For GPS strategies, decreasing genotyping cost from 37€ to 10€ had no significant effect on the genetic gain (the average genetic gain was 16.4 for 37€, 16.5 for 10€).

<table><caption>**Table 2.** Contribution of input parameters to the genetic gain.</caption><thead><tr><th>Trait</th><th>Factor</th><th>F</th><th>P-value</th><th>% of SS</th><th>Factor level</th><th>Number of records</th><th>Means</th></tr></thead><tbody><tr><td rowspan="10">Genetic gain</td><td rowspan="3">Strategy</td><td rowspan="3">1909.9</td><td rowspan="3">&lt;2 × 10-16</td><td rowspan="3">16.7</td><td>PS</td><td>7200</td><td>15.3</td></tr><tr><td>GPS.n2</td><td>7200</td><td>15.0</td></tr><tr><td>GPSopt.n2</td><td>7200</td><td>18.0</td></tr><tr><td rowspan="3">λ</td><td rowspan="3">300.0</td><td rowspan="3">&lt;2 × 10-16</td><td rowspan="3">0.7</td><td>0.25</td><td>12,000</td><td>16.3</td></tr><tr><td>0.5</td><td>12,000</td><td>16.7</td></tr><tr><td>0.75</td><td>12,000</td><td>15.6</td></tr><tr><td rowspan="3">h<sup>2</sup></td><td rowspan="3">18,379.2</td><td rowspan="3">&lt;2 × 10-16</td><td rowspan="3">40.1</td><td>0.2</td><td>12,000</td><td>13.3</td></tr><tr><td>0.4</td><td>12,000</td><td>16.4</td></tr><tr><td>0.7</td><td>12,000</td><td>18.9</td></tr><tr><td rowspan="12">Average annual CT</td><td rowspan="12">1123.2</td><td rowspan="12">&lt;2 × 10-16</td><td rowspan="12">2.5</td><td>1.5M€</td><td>18,000</td><td>15.6</td></tr><tr><td>3M€</td><td>18,000</td><td>16.7</td></tr><tr><td>PS, 0.25</td><td>2400</td><td>14.6</td></tr><tr><td>PS, 0.5</td><td>2400</td><td>15.62</td></tr><tr><td>PS, 0.75</td><td>2400</td><td>15.8</td></tr><tr><td>GPS.n2, 0.25</td><td>2400</td><td>14.9</td></tr><tr><td>GPS.n2, 0.5</td><td>2400</td><td>15.4</td></tr><tr><td>GPS.n2, 0.75</td><td>2400</td><td>14.3</td></tr><tr><td>GPSopt.n2, 0.25</td><td>2400</td><td>18.6</td></tr><tr><td>GPSopt.n2, 0.5</td><td>2400</td><td>18.5</td></tr><tr><td>GPSopt.n2, 0.75</td><td>2400</td><td>16.9</td></tr><tr><td>Strategy * λ</td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr></tbody></table>

h²: trait heritability. λ: relative selection rate. CT: Total cost. CG: Genotyping cost (37€). % of SS: (Sum of squares) / (Total sum of squares). Number of records: number of records for each factor value. Means: average genetic gain for each factor value. Results were similar for GPS.Nc and GPS.n2, and for GPSopt.Nc and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt). QTL sampling accounted for less than 1%, for that reason we do not detail the mean of each of the 20 traits in this table.

Considering the interaction between strategy and one of the other input parameters, we found that strategy and λ had the most significant effect on final genetic gain (Table 2). The genetic gain for GPS strategies (GPS.Nc, GPS.n2, GPSopt.Nc and GPSopt.n2) was always larger than genetic gain for PS when λ = 0.25 (i.e., α₃ > α₄, N₃ is decreased and N₄ is increased compared to λ = 0.75, the selection is less stringent at step 3 where GPS is applied). The interaction between strategy and total cost and the
interaction between strategy and trait heritability had a significant effect on final genetic gain but they accounted for less than 1% of the sum of squares. Finally, the interaction between strategy and QTL sampling (i.e., the simulated trait) had no significant effect on final genetic gain.

## Comparison of the Evolution of Genetic Gain between Scenarios

As expected, we observed an increase of the 200 best lines average TBV along generations. Figure 2 illustrates the cumulative genetic gain at the end of each cycle with an average annual CT = 3M€ and CG = 37€. Results for other values of CT and CG are available in Supplementary Table S6.

Genetic gain was always significantly larger for GPSopt.Nc and GPSopt.n2 than for PS, GPS.Nc and GPS.n2 whatever the values of CT, λ, CG and h². The difference between GPSopt.Nc or GPSopt.n2 and PS was even larger when λ = 0.25 (α₃ > α₄). Indeed, the average differences between GPSopt and PS were +3.94 for λ = 0.25, +2.90 for λ = 0.5 and +1.46 for λ = 0.75 (Table 2). Genetic gain was also larger for GPS.Nc and GPS.n2, compared to PS when λ = 0.25. However, this difference was significant in only one third of the scenarios. When λ = 0.5, no significant difference was observed between PS, GPS.Nc and GPS.n2. In contrast, when λ = 0.75, genetic gain was significantly larger for PS compared to GPS. Nc and GPS.n2 (+1.55).

Note that GPS. Nc and GPSopt. Nc allowed to do 156 more crosses per cycle on average compared to GPS. n2 and GPSopt. n2. The number of progenies per cross was decreased from 412 to 333 for n3 and from 82 to 67 for n4 (with λ = 0.25, average annual CT = 3 M€ and CG = 37€). However, the genetic gain for GPS. Nc (GPSopt. Nc) and GPS.n2 (GPSopt. n2) was not significantly different whatever the scenario.

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Evolution of genetic gain. h²: trait heritability. λ: relative selection rate. Average annual total cost (CT) = 3M€ and genotyping cost (CG) = 37€. Results were similar for GPS.Nc and GPS.n2, and for GPSopt.Nc and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt).</figcaption></figure>
## Evolution of Correlations Between TBV and GEBV or Phenotypic Values

We calculated the Pearson correlation between TBV and GEBV (GS accuracy) and between TBV and phenotypic values (PS accuracy) for the progenies at step 3 (Table 3). As expected, accuracy increased with trait heritability. Whatever the values of trait heritability and $\lambda$, PS accuracy at step 3 slightly decreased over the cycles. In contrast, GS accuracy increased over cycles. For instance, for GPSopt strategies, accuracy at step 3 of third cycle was on average twice as high as the one at the first cycle. Although $\lambda$ had low impact on PS accuracy, it was significant on GPS and GPSopt. When $\lambda$ = 0.25 ($\alpha_3$ > $\alpha_4$), final accuracy was superior compared to scenarios with $\lambda$ = 0.75 (+ 0.12 and + 0.17 in GPS and GPSopt respectively). When $\lambda$ decreases, the number of progenies $N_4$ evaluated in the field at step 4 increases (Supplementary Table S3) and the number of lines added to the training population increases.

<table><caption>**Table 3.** PS and GS accuracies at step 3.</caption><thead><tr><th rowspan="3">h<sup>2</sup></th><th rowspan="3">strategy</th><th colspan="9">Λ</th></tr><tr><th colspan="3">0.25</th><th colspan="3">0.50</th><th colspan="3">0.75</th></tr><tr><th>Cycle 1</th><th>Cycle 2</th><th>Cycle 3</th><th>Cycle 1</th><th>Cycle 2</th><th>Cycle 3</th><th>Cycle 1</th><th>Cycle 2</th><th>Cycle 3</th></tr></thead><tbody><tr><td rowspan="3">0.2</td><td>PS</td><td>0.56</td><td>0.54</td><td>0.53</td><td>0.56</td><td>0.53</td><td>0.53</td><td>0.56</td><td>0.53</td><td>0.53</td></tr><tr><td>GPS</td><td>0.37</td><td>0.70</td><td>0.79</td><td>0.39</td><td>0.64</td><td>0.74</td><td>0.38</td><td>0.56</td><td>0.65</td></tr><tr><td>GPSopt</td><td>0.31</td><td>0.79</td><td>0.88</td><td>0.30</td><td>0.71</td><td>0.80</td><td>0.29</td><td>0.59</td><td>0.65</td></tr><tr><td rowspan="3">0.4</td><td>PS</td><td>0.73</td><td>0.70</td><td>0.69</td><td>0.73</td><td>0.69</td><td>0.69</td><td>0.73</td><td>0.69</td><td>0.69</td></tr><tr><td>GPS</td><td>0.42</td><td>0.76</td><td>0.86</td><td>0.42</td><td>0.70</td><td>0.80</td><td>0.42</td><td>0.63</td><td>0.73</td></tr><tr><td>GPSopt</td><td>0.37</td><td>0.85</td><td>0.92</td><td>0.37</td><td>0.80</td><td>0.86</td><td>0.38</td><td>0.70</td><td>0.76</td></tr><tr><td rowspan="3">0.7</td><td>PS</td><td>0.88</td><td>0.86</td><td>0.86</td><td>0.88</td><td>0.86</td><td>0.86</td><td>0.88</td><td>0.86</td><td>0.85</td></tr><tr><td>GPS</td><td>0.46</td><td>0.83</td><td>0.90</td><td>0.46</td><td>0.77</td><td>0.87</td><td>0.46</td><td>0.70</td><td>0.81</td></tr><tr><td>GPSopt</td><td>0.45</td><td>0.92</td><td>0.95</td><td>0.45</td><td>0.88</td><td>0.92</td><td>0.46</td><td>0.80</td><td>0.85</td></tr></tbody></table>

h²: trait heritability. $\lambda$: relative selection rate. Average annual total cost (CT) = 3 M€ and genotyping cost (CG) = 37€. Results were similar for GPS.Nc and GPS.n2, and for GPSopt.Nc and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt).

## Contribution of Input Parameters to Genetic Diversity Evolution

As expected, we observed a decrease of polymorphism rate (calculated as the percentage of alleles present in both the reference lines and the selected progenies) over cycles for all scenarios and strategies.

The five input factors had a significant effect on polymorphism rate (percentage of alleles) in the 200 best progenies of the last cycle compared
to the reference (initial) population (Supplementary Table S8). Figure 3 illustrates the results with an average annual CT = 3M€ and CG = 37€. Results for other values of CT and CG are available in Supplementary file. Each of the five factors had a significant effect. Strategy and λ accounted for 39.1% and 13.7% of the sum of squares, respectively. Trait heritability, annual total cost and QTL sampling (the simulated trait) accounted for less than 1% of the sum of squares. This decrease was more pronounced for GPS breeding schemes than for PS, when crosses were optimized in particular, and when λ = 0.75.

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

**Figure 3.** Evolution of polymorphism rate over cycles. h²: trait heritability. λ: relative selection rate. Annual total cost (CT) = 3 M€ and genotyping cost (CG) = 37€. The evolution of polymorphism rate was calculated as the percentage of alleles present in both the reference population and the selected progenies. Results were similar for GPS.Nc and GPS.n2, and for GPSopt. Nc and GPSopt. n2, for that reason we only show results for PS, GPS. n2 (GPS) and GPSopt. n2 (GPSopt).

We also estimated the evolution of the expected heterozygoty (He) over the cycles (Supplementary Table S7) and the difference He in the reference population and at the end of each cycle for the all set of markers except QTLs, for markers assigned to QTLs and for markers at the vicinity of QTLs (1, 5 and 10 cM) (Table 4). Initial He in the reference population (He ini) was 0.26. He decreased faster for GPSopt strategies (-0.14 < effect at the end of third cycle < -0.09) whatever the scenario compared to GPS strategies (-0.08 < effect < -0.03) and PS strategy (-0.03 < effect < -0.02). The difference was even larger when λ = 0.75 for GPS strategies (-0.14 < effect < -0.12 for GPSopt, -0.1 < effect < -0.08 for GPS). The loss of genetic diversity was stronger at QTL locations but not at the vicinity of QTLs.

The number of additional favourable alleles cumulated at the end of the third cycle was also significantly larger for GPSopt strategies. For GPSopt, it ranged from 15.57 for h² = 0.2 and λ = 0.75 to 22.39 for h² = 0.7 and λ = 0.25. For GPS, it ranged from 13.2 for h² = 0.2 and λ = 0.75 to 19.1 for h² = 0.7 and λ = 0.5. For PS, it ranged from 12.75 for h² = 0.2 and λ = 0.25 to 19.47 for h² = 0.7 and λ = 0.75.
<table><caption>**Table 4.** Evolution of genetic gain, number of favourable alleles and diversity.</caption><thead><tr><th>Strategy</th><th>h<sup>2</sup></th><th>λ</th><th>Gain<sup>a</sup></th><th>Nfav<sup>b</sup></th><th>He<sup>c</sup></th><th>He.QTL<sup>d</sup></th><th>He.QTL1<sup>e</sup></th><th>He.QTL5<sup>f</sup></th></tr></thead><tbody><tr><td>GPSopt</td><td>0.7</td><td>0.25</td><td>25.7 (a)</td><td>22.39 (a)</td><td>-0.09 (jk)</td><td>-0.14 (k)</td><td>-0.1 (j)</td><td>-0.1 (j)</td></tr><tr><td>GPSopt</td><td>0.7</td><td>0.5</td><td>25.6 (a)</td><td>22.5 (a)</td><td>-0.1 (n)</td><td>-0.16 (m)</td><td>-0.12 (k)</td><td>-0.11 (k)</td></tr><tr><td>GPSopt</td><td>0.7</td><td>0.75</td><td>24.3 (b)</td><td>21.27 (b)</td><td>-0.12 (p)</td><td>-0.17 (n)</td><td>-0.13 (n)</td><td>-0.13 (n)</td></tr><tr><td>GPSopt</td><td>0.4</td><td>0.25</td><td>23.3 (c)</td><td>19.94 (c)</td><td>-0.09 (kl)</td><td>-0.14 (j)</td><td>-0.1 (j)</td><td>-0.1 (j)</td></tr><tr><td>GPSopt</td><td>0.4</td><td>0.5</td><td>23.2 (c)</td><td>20.01 (c)</td><td>-0.11 (o)</td><td>-0.15 (lm)</td><td>-0.12 (l)</td><td>-0.12 (l)</td></tr><tr><td>PS</td><td>0.7</td><td>0.75</td><td>23 (cd)</td><td>19.47 (d)</td><td>-0.03 (de)</td><td>-0.09 (f)</td><td>-0.04 (e)</td><td>-0.04 (de)</td></tr><tr><td>PS</td><td>0.7</td><td>0.5</td><td>22.9 (cd)</td><td>19.35 (de)</td><td>-0.03 (de)</td><td>-0.09 (f)</td><td>-0.04 (e)</td><td>-0.04 (de)</td></tr><tr><td>GPS</td><td>0.7</td><td>0.5</td><td>22.6 (de)</td><td>19.1 (e)</td><td>-0.05 (f)</td><td>-0.1 (g)</td><td>-0.06 (f)</td><td>-0.06 (f)</td></tr><tr><td>GPS</td><td>0.7</td><td>0.25</td><td>22.3 (e)</td><td>18.75 (f)</td><td>-0.04 (e)</td><td>-0.09 (f)</td><td>-0.05 (e)</td><td>-0.04 (e)</td></tr><tr><td>PS</td><td>0.7</td><td>0.25</td><td>22.2 (e)</td><td>18.66 (f)</td><td>-0.03 (cd)</td><td>-0.08 (e)</td><td>-0.04 (cd)</td><td>-0.04 (cd)</td></tr><tr><td>GPSopt</td><td>0.4</td><td>0.75</td><td>21.6 (f)</td><td>18.47 (f)</td><td>-0.13 (q)</td><td>-0.17 (n)</td><td>-0.14 (o)</td><td>-0.14 (o)</td></tr><tr><td>GPS</td><td>0.7</td><td>0.75</td><td>21.4 (f)</td><td>18.08 (g)</td><td>-0.08 (i)</td><td>-0.12 (h)</td><td>-0.08 (h)</td><td>-0.08 (h)</td></tr><tr><td>PS</td><td>0.4</td><td>0.75</td><td>20.7 (g)</td><td>17.11 (hi)</td><td>-0.03 (bc)</td><td>-0.07 (d)</td><td>-0.04 (bc)</td><td>-0.03 (bc)</td></tr><tr><td>GPSopt</td><td>0.2</td><td>0.25</td><td>20.5 (g)</td><td>17.19 (h)</td><td>-0.1 (l)</td><td>-0.13 (i)</td><td>-0.1 (j)</td><td>-0.1 (j)</td></tr><tr><td>PS</td><td>0.4</td><td>0.5</td><td>20.5 (g)</td><td>16.98 (hi)</td><td>-0.03 (bc)</td><td>-0.07 (d)</td><td>-0.04 (bc)</td><td>-0.03 (bc)</td></tr><tr><td>GPSopt</td><td>0.2</td><td>0.5</td><td>20.4 (g)</td><td>17.28 (h)</td><td>-0.12 (p)</td><td>-0.15 (l)</td><td>-0.12 (m)</td><td>-0.12 (m)</td></tr><tr><td>GPS</td><td>0.4</td><td>0.5</td><td>20.2 (g)</td><td>16.82 (i)</td><td>-0.06 (g)</td><td>-0.1 (g)</td><td>-0.06 (g)</td><td>-0.06 (g)</td></tr><tr><td>GPS</td><td>0.4</td><td>0.25</td><td>19.7 (h)</td><td>16.25 (j)</td><td>-0.04 (e)</td><td>-0.07 (d)</td><td>-0.04 (de)</td><td>-0.04 (de)</td></tr><tr><td>PS</td><td>0.4</td><td>0.25</td><td>19.4 (h)</td><td>15.97 (jk)</td><td>-0.03 (b)</td><td>-0.06 (c)</td><td>-0.03 (b)</td><td>-0.03 (b)</td></tr><tr><td>GPS</td><td>0.4</td><td>0.75</td><td>18.9 (i)</td><td>15.66 (kl)</td><td>-0.09 (j)</td><td>-0.12 (h)</td><td>-0.1 (i)</td><td>-0.09 (i)</td></tr><tr><td>GPSopt</td><td>0.2</td><td>0.75</td><td>18.5 (i)</td><td>15.57 (l)</td><td>-0.14 (r)</td><td>-0.17 (n)</td><td>-0.15 (p)</td><td>-0.15 (p)</td></tr><tr><td>GPS</td><td>0.2</td><td>0.5</td><td>17.4 (j)</td><td>14.22 (m)</td><td>-0.06 (h)</td><td>-0.09 (f)</td><td>-0.07 (g)</td><td>-0.07 (g)</td></tr><tr><td>PS</td><td>0.2</td><td>0.75</td><td>17.3 (j)</td><td>13.99 (mn)</td><td>-0.02 (a)</td><td>-0.05 (a)</td><td>-0.02 (a)</td><td>-0.02 (a)</td></tr><tr><td>PS</td><td>0.2</td><td>0.5</td><td>17 (j)</td><td>13.86 (n)</td><td>-0.02 (a)</td><td>-0.05 (a)</td><td>-0.03 (a)</td><td>-0.02 (a)</td></tr><tr><td>GPS</td><td>0.2</td><td>0.25</td><td>16.4 (k)</td><td>13.32 (o)</td><td>-0.03 (cd)</td><td>-0.06 (b)</td><td>-0.04 (bc)</td><td>-0.03 (bc)</td></tr><tr><td>GPS</td><td>0.2</td><td>0.75</td><td>16.2 (kl)</td><td>13.2 (o)</td><td>-0.1 (m)</td><td>-0.12 (h)</td><td>-0.1 (j)</td><td>-0.1 (j)</td></tr><tr><td>PS</td><td>0.2</td><td>0.25</td><td>15.8 (l)</td><td>12.75 (p)</td><td>-0.02 (a)</td><td>-0.04 (a)</td><td>-0.02 (a)</td><td>-0.02 (a)</td></tr></tbody></table>

a: Genetic gain (TBV difference between the 200 selected lines at year 15 and the reference population); b: number of favourable alleles difference; c: diversity (He) difference along the genome (QTL excluded); d: diversity (He) difference at QTLs; e: diversity (He) difference at 1 cM interval around QTLs; f: diversity (He) difference at 5 cM interval around QTLs; Two different crossing strategies are compared, random crossing and optimized crossing using UC criterion (opt), three different levels of heritability (0.2,0.4,0.7) and three different levels of selection ratio λ between steps 3 and 4 (0.25,0.5,0.75). The letters in brackets (a-r) correspond to significant different groups according to SNK test (α = 10⁻⁵). Results were similar for GPS.Nc and GPS.n2, and for GPSopt.Nc and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt).

## Parental Contributions

The number of reference lines that contributed to progenies decreased over cycles (Supplementary Table S9). When average annual CT = 3 M€ and C<sub>G</sub> = 37€, on average 31 reference (founder) lines were used as parents to produce the 200 best lines of the last cycle in GPSopt strategies. In contrast, on average 53 founder lines were used as parents in other GPS strategies and 81 in PS strategy. We also noticed that the number of founder lines was slightly smaller for GPS strategies with 400 crosses
(GPS.N<sub>c</sub> and GPSopt.N<sub>c</sub>) compared to GPS strategies with more crosses (GPS.n2 and GPSopt.n2).

We also calculated the number of reference lines that contributed to 25% or 75% of progeny. These two indicators followed the same trends as the number of reference lines that contributed to the whole progeny (Supplementary Table S9). Indeed, only 1.4 reference lines contributed to 25% of progeny for GPSopt strategy on average compared to 3.2 and 4.8 different initial lines were needed in GPS and PS strategies, respectively. On average, 7 founders contributed to 75% of progeny, compared to 17 and 27.6 for GPS and PS strategies, respectively. These reveals that some parents contribute to a large number of selected progenies for GPSopt if no contribution threshold is applied. This is consistent with the Shannon index that was smaller for GPSopt approach at the end of the third cycle (2.4,3.3 or 3.8 in GPSopt, GPS or PS).

<table><caption>**Table 5.** Number of common and specific favourable alleles between parents.</caption><thead><tr><th>Strategy</th><th>h<sup>2</sup></th><th>com0<sup>a</sup></th><th>spe0<sup>b</sup></th><th>ef.com0<sup>c</sup></th><th>ef.spe0<sup>d</sup></th><th>fq.ef.com0<sup>e</sup></th><th>fq.ef.spe0<sup>f</sup></th><th>com5<sup>g</sup></th><th>com10<sup>h</sup></th><th>com15<sup>i</sup></th><th>spe5<sup>j</sup></th><th>spe10<sup>k</sup></th><th>spe15<sup>l</sup></th></tr></thead><tbody><tr><td>GPSopt</td><td>0.2</td><td>42 (a;41)</td><td>25 (b;25)</td><td>36 (a;34)</td><td>20 (b;20)</td><td>28 (a;27)</td><td>8 (b;9)</td><td>43 (b)</td><td>48 (b)</td><td>55 (c)</td><td>23 (d)</td><td>22 (e)</td><td>17 (e)</td></tr><tr><td>GPSopt</td><td>0.4</td><td>42 (a;41)</td><td>26 (ab;26)</td><td>37 (a;35)</td><td>20 (ab;21)</td><td>29 (a;28)</td><td>8 (b;9)</td><td>43 (ab)</td><td>49 (a)</td><td>57 (b)</td><td>24 (c)</td><td>23 (c)</td><td>18 (d)</td></tr><tr><td>GPSopt</td><td>0.7</td><td>43 (a;42)</td><td>26 (a;26)</td><td>38 (a;36)</td><td>21 (ab;21)</td><td>29 (a;28)</td><td>8 (b;9)</td><td>44 (a)</td><td>49 (a)</td><td>58 (a)</td><td>26 (a)</td><td>24 (a)</td><td>19 (c)</td></tr><tr><td>GPS</td><td>0.2</td><td>37 (b;37)</td><td>26 (ab;26)</td><td>30 (b;30)</td><td>21 (ab;21)</td><td>24 (b;25)</td><td>10 (a;10)</td><td>41 (d)</td><td>46 (d)</td><td>51 (f)</td><td>25 (b)</td><td>22 (d)</td><td>20 (b)</td></tr><tr><td>GPS</td><td>0.4</td><td>37 (b;37)</td><td>26 (ab;26)</td><td>30 (b;30)</td><td>21 (ab;21)</td><td>25 (b;25)</td><td>10 (a;10)</td><td>42 (c)</td><td>47 (c)</td><td>53 (e)</td><td>25 (b)</td><td>23 (c)</td><td>20 (b)</td></tr><tr><td>GPS</td><td>0.7</td><td>37 (b;37)</td><td>26 (ab;26)</td><td>30 (b;30)</td><td>21 (a;21)</td><td>24 (b;25)</td><td>10 (a;10)</td><td>43 (c)</td><td>48 (b)</td><td>54 (d)</td><td>26 (a)</td><td>24 (b)</td><td>21 (a)</td></tr></tbody></table>

We also observed that the number of favourable alleles present in both parents was larger in GPSopt compared to the other strategies (Table 5) and that this number increased over cycles. On the opposite, the number of favourable alleles specific to one parent was not different between optimized and random crossing strategies in our elite germplasm (reference population) but became significantly inferior over cycles as the number of favourable alleles and inbreeding increased. There was 5 more QTLs in common between parents at the first generation in GPSopt compared to GPS, 4 more at the third cycle (year 15). Finally, GPSopt was more efficient to increase the frequency of favourable alleles, especially

<footnote><sup>a</sup> number of common favourable alleles between parents of the 10 (400) best crosses according to UC at generation 0 (among genitors from the reference population); <sup>b</sup> number of favourable alleles that are specific to one of the parent of the reference population; <sup>c</sup> number of common favourable alleles multiplied by their effect; <sup>d</sup> number of specific favourable alleles multiplied by their effect; <sup>e</sup> number of common favourable alleles multiplied by their effect and frequency; <sup>f</sup> number of specific favourable alleles multiplied by their effect and frequency; <sup>g</sup> number of common favourable alleles between parents of the 200 lines selected at year 5; <sup>h</sup> number of specific favourable alleles between parents of the 200 lines selected at year 5; <sup>i</sup> number of common favourable alleles between parents of the 200 lines selected at year 10; <sup>j</sup> number of specific favourable alleles between parents of the 200 lines selected at year 10; <sup>k</sup> number of common favourable alleles between parents of the 200 lines selected at year 15; <sup>l</sup> number of specific favourable alleles between parents of the 200 lines selected at year 15. Two different crossing strategies are compared, random crossing and optimized crossing using UC criterion (opt), three different levels of heritability (0.2,0.4,0.7). The letters in brackets (a-f) correspond to significant different groups according to SNK test (α = 10<sup>-5</sup>). Results were similar for GPS.N<sub>c</sub> and GPS.n2, and for GPSopt.N<sub>c</sub> and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt).</footnote>
when those alleles were rare in the reference population, heritability was
high and $\lambda$ was low ($\alpha_3$ was high) (Figure 4).

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Increase of favourable alleles frequency. $h^2$: trait heritability. $\lambda$: relative selection rate. Annual total cost (CT) = 3 M€ and genotyping cost (CG) = 37€. Results were similar for GPS.N<sub>c</sub> and GPS.n2, and for GPSopt.N<sub>c</sub> and GPSopt.n2, for that reason we only show results for PS, GPS.n2 (GPS) and GPSopt.n2 (GPSopt).</figcaption></figure>

DISCUSSION

This study focused on two types of breeding schemes: Phenotypic
Selection (PS) with two steps of selection based on field trials, and Genomic
and Phenotypic Selection (GPS) that combines a first step of genomic
selection and a second step of selection based on field trials. We finalized
the cycle by phenotypic selection because plant breeders have to provide
field data for candidates when applying for official variety registration. In
order to compare PS and GPS schemes with a fixed budget, GPS could
either have the same number of crosses than PS (GPS.n2) or it could have
the same number of progenies per cross at the beginning of step 2 (GPS.Nc).
In addition, we explored two methods to choose couples to be crossed at
each cycle for GPS schemes: random mating of the best lines (GPS.Nc or
GPS.n2 strategies) or optimized mating (GPSopt.Nc or GPSopt.n2 strategies)
using the optimal haploid value (OHV) usefulness criterion [53]. We tested
several scenarios by varying cost of the breeding scheme, genotyping cost,
trait heritability and selection rates at each step of selection, in order to
investigate under which conditions GPS was more efficient than PS. We
evaluated each strategy by comparing the genetic gain of lines selected at
the end of the breeding program and the evolution of genetic diversity
over cycles. Note that computational time was three to five times longer
for the simulations of breeding schemes with optimized mating.

Comparison of Genetic Gain of Selected Lines at the End of the
Breeding Programs

The objective of bread wheat breeding programs is to develop new
varieties that outperform existing varieties, particularly those used as
control in official registration trials. To reach this goal, breeders have to
choose the strategy that leads to maximal True Breeding Value (TBV, with
the hope that the realized phenotype will also be superior. As expected,
the TBV of the 200 best lines increased over cycles for all strategies and
scenarios in our simulations.

We observed that trait heritability and strategy (PS, GPS, GPSoft) are
the two main factors affecting genetic gain which is consistent with the
breeder's equation [4]. To a lesser extent, relative selection rate (λ)
between genomic and phenotypic selection steps had a significant effect
on the genetic gain variance.

In this study, we simulated traits controlled by 100 QTLs randomly
sampled from 12,119 genomic markers. The next step will be to vary the
number of sampled QTLs, the distribution of their effects and epistasis.

We simulated breeding programs with two different budgets. Doubling
the global budget did not lead to a large increase of the genetic gain (+ 6%
on average). The genotyping cost had no significant effect on the genetic
gain variance. This may be due to DH expense that is so high that varying
the other costs has low impact on progenies number. This suggests that the
INRAE-AO breeding program and reference panel may not have reached
the necessary critical size to benefit fully from genomic predictions.

We however identified under which conditions genomic predictions
were the most interesting. We compared two strategies, genotyping
predictions were used (1) to select candidates for next steps (GPS.Nc and
GPS.n2 strategies), or (2) to select candidates and optimize crosses
(GPSopt.Nc and GPSopt.n2 strategies). In our germplasm, we found that
only strategies with optimized crosses were significantly superior to PS in
terms of final genetic gain. We showed that parents involved in crosses
that have led to the production of selected lines had a larger number of
favourable alleles when the crosses were optimized. In addition, we
highlighted the fact that some parents were involved in a large proportion
of crosses in GPSopt strategies, leading to efficient favourable allele
frequency increase, even rare alleles, but also global decrease of diversity
and a genetic gain slowdown. Strategies with optimized mating were more
advantageous when the selection rate was more intense at the step of
phenotypic selection (step 4) than at the step of genomic selection (step 3).
In this case, the number of progenies that were both genotyped and
phenotyped at each cycle in GPS and GPSopt strategies were larger. As
genomic prediction equations were updated at each cycle, the training
population was larger and the predictions were more accurate. This result
is consistent with previous studies that have shown that training
population size and composition are key factors affecting accuracy of
genomic predictions [15,16].

In our simulations GPS and PS breeding programs required the same
number of years. It would be interesting to simulate more realistic
breeding programs were genomic predictions results in accelerating the
breeding programs using rapid cycles for instance [44,54,55].
## Evolution of Genetic Diversity

Genetic diversity is a key parameter in plant breeding since it has an effect on long term genetic gain (according to breeder's equation [4]). However, it has been shown in both experimental study and by stochastic simulations [1,56] that GS accelerates the loss of diversity compared to phenotypic selection due to the rapid fixation of regions of the genome with an effect on the trait of interest. Our results were consistent with those studies. Indeed, we observed a faster decrease of polymorphism rate in GPS strategies compared to PS schemes whatever the simulated scenario. This loss of alleles was even faster for strategies with optimized crosses (GPSopt.N<sub>c</sub> and GPSopt.n2). In order to reduce this loss, we could place additional weight on low-frequency favourable alleles as recommended by [1]. Since the loss of diversity is particularly fast when crosses were optimized, a second option would be to define a maximum number of crosses for which each parent could contribute in order to avoid having too many lines with one (or two) common parent like the optimal contribution selection (OCS) strategies [49]. Note that favourable alleles were more efficiently fixed for GPS and GPSopt and the frequency of favourable alleles, the rare ones in particular increased faster. In this study, we considered that only lines from the previous generation could be mated. This assumes no introduction of neo-diversity in the breeding scheme, which is not realistic, but simulations are always a simplification of real life. Indeed, it is common use to introduce external bread wheat registered material at each cycle. In order to simulate more realistic breeding programs, it would be important to give the possibility to add parents from an external pool (for example lines selected by other breeders or coming from genetic resources collections).

## Resource Allocation

A great challenge in plant breeding is to improve and accelerate the genetic gain with a fixed budget. To meet this goal, breeders have to evaluate the best way to allocate resources. In this study, we found that production of doubled haploids (DHs), used to reduce the time required for inbred development [57], was the major expense in the INRAE-AO breeding program. Indeed, up to 46.5% and 72.6% of the global budget were allocated to DH production in PS and GPS schemes respectively. Therefore, it could be interesting to simulate breeding schemes based on other breeding methods than DH production, Single Seed Descent (SSD) in rapid cycles for instance. It would also double the number of efficient crossovers, increasing the probability to obtain beneficial recombinations. Another possibility to reduce the total cost of production of DHs without reducing the genetic gain would be to genotype F2 and produce DH only for F2 with highest GEBVs.

We also noticed that field evaluation accounted for a significant part of the budget (up to 46.2% for PS and up to 35.9% for GPS strategies). In
addition, the percentage of the global budget allocated to field evaluation was higher for the most interesting scenarios in terms of genetic gain, i.e., the scenarios with the higher number of lines at step 4. To reduce phenotypic cost, each line could be phenotyped in a relevant subset of trials. We could optimize the experimental design in order to decrease the number of lines observed in each environment without decreasing the number of alleles observed in each environment and the selection accuracy [8]. When some cheap traits correlated to the targeted traits exist, we could also improve resource allocation by optimizing phenotyping between target and secondary traits [37,38,58–61]. We could simulate traits that are controlled by QTLs with partly pleiotropic effects [62].

## Perspectives for the Improvement of the Pipeline

We highlighted the importance of optimizing crosses using genomic predictions. We used as a first test a usefulness criterion that is fast to calculate (OHV). It calculates the TBV of the progeny that would get the best chromosome from each parent. This is not realistic with limited progeny size, since the probability of obtaining this best OHV can be very low, and this criterion assumes no recombination. We will include in the second version of the pipeline some usefulness criterion that takes into account the actual recombination rates across the genome for a limited number of progenies [11,63]. We could also give the possibility to monitor long term diversity by estimating simultaneously genetic gain and parental contribution using OCS and usefulness criterion parental contribution (UCPC) methodologies [12,49]. This way, we will be able to adjust the number of progenies according to predicted cross values using genetic algorithms.

In this study, we replaced one step of phenotypic selection by genomic selection. However, genomic selection predictions could also be used at an earlier stage. Indeed, genomic prediction could be applied during the second step in order to reduce the budget allocated to multiplication. In this case, the number of genotyped lines would be higher. In addition, genomic predictions could be used to accelerate recurrent selection in a two or three-part breeding scheme composed of: (i) a population improvement component (pre-breeding) through recurrent selection, (ii) an optional bridging component if the performance gap between donors and elites is too large [50] and (iii) a traditional breeding program component [64].

In this study, we focused on single-trait selection. However, real breeding programs most often deal with simultaneous improvement of several traits that can be negatively or un-correlated with each other. In a multi-trait context, replacing one step of phenotyping by a step of genotyping would be even more cost effective. Simulations of breeding programs with several selection objectives would be more realistic. We will have to build a selection index [65] that depends on the relative economical weight of each selected traits. In addition, a multi-objective
optimization framework was proposed by [48] to define the best compromise.

## CONCLUSIONS

In this study, we showed that GPS selection using mating optimization in INRAE-AO breeding program significantly improved genetic gain over successive cycles for different levels of heritability, selection rates and genotyping cost, compared to PS. In contrast, GPS without mating optimization and PS had a similar efficiency in terms of genetic gain. In addition, the loss of genetic diversity was faster in GPS breeding programs. Our results also highlighted that GPS strategies were more efficient to increase favourable allele frequency, rare alleles in particular. The results may be different in programs with different genotyping and phenotyping costs.

The simulation pipeline we developed can help breeders to test the effectiveness of their breeding programs when changing parameters (genotyping and phenotyping costs, heritability of the trait, selection intensity, progeny size, number of crosses, cross value). For instance, we could test a scenario in which F2 are genotyped and only those with high genetic variance prediction are used for DH production. This would accelerate the cycles [13]. In the next version, the pipeline will also improve several trait simultaneously and be able to integrate external varieties at each cycle.

## SUPPLEMENTARY MATERIALS

The following supplementary materials are available online: https://doi.org/10.20900/cbgg20210008. Supplementary Table S1: Operation costs; Supplementary Table S2: Description of the 36 scenarios; Supplementary Table S3: Progeny size under different scenarios; Supplementary Table S4: Contribution of input parameters to genetic gain; Supplementary Table S5: Contribution of input parameters to genetic gain for the 10, 50 and 100 best lines; Supplementary Table S6: Evolution of genetic gain and polymorphism rate; Supplementary Table S7: Contribution of input parameters to polymorphism rate; Supplementary Table S8: Evolution of expected heterozygosity (He) over cycles; Supplementary Table S9: Evolution of parental contribution over cycles; Supplementary Figure S1. Description of the algorithm.

## DATA AVAILABILITY

The dataset analyzed in the study can be found in the INRAE Dataverse repository:

https://data.inrae.fr/dataset.xhtml?persistentId=doi:10.15454/M8SAYH.
## AUTHOR CONTRIBUTIONS

JA, SL, AFS, GC, and SB designed the simulated breeding programmes and defined costs of each part of the breeding programmes. SBS made the simulations. SBS and SB analyzed the simulations and wrote the manuscript. SL, AFS and GC helped improving the manuscript. All authors approved the final manuscript.

## CONFLICTS OF INTEREST

The authors declare that they have no conflicts of interest.

## FUNDING

Doctoral work of SBS was funded by a grant from the Auvergne Rhônes-Alpes region and from the European Regional Development Fund (FEDER). This work was supported by the Breed wheat project thanks to the funding from the French Government managed by the National Research Agency (ANR) in the framework of Investments for the Future (ANR-10-BTBR-03) France AgriMer and the French Fund to support Plant Breeding (FSOV).

## ACKNOWLEDGMENTS

The authors thank the genotyping platform GENTYANE at INRAE Clermont-Ferrand (gentyane.clermont.inra.fr) which has conducted the genotyping. The authors also thank their INRAE colleagues who managed the INRAE-Agri-Obtentions bread wheat breeding programme.

## REFERENCES

1. Jannink JL. Dynamics of long-term genomic selection. Genet Sel Evol. 2010;42(1):35. doi: 10.1186/1297-9686-42-35

2. Whittaker JC, Thompson R, Denham MC. Marker-assisted selection using ridge regression. Genet Res. 2000;75(2):249-52. doi: 10.1017/S0016672399004462

3. Meuwissen THE, Hayes BJ, Goddard ME. Prediction of Total Genetic Value Using Genome-Wide Dense Marker Maps. Genetics. 2001;157(4):1819-29.

4. Lush JL. Animal Breeding Plans. Ann Arbor (US): Collegiate Press; 1937.

5. Crossa J, Pérez-Rodríguez P, Cuevas J, Montesinos-López O, Jarquín D, de los Campos G, et al. Genomic Selection in Plant Breeding: Methods, Models, and Perspectives. Trends Plant Sci. 2017;22(11):961-75.

6. Lorenzana RE, Bernardo R. Accuracy of genotypic value predictions for marker-based selection in biparental plant populations. Theor Appl Genet. 2009;120(1):151-61. doi: 10.1007/s00122-009-1166-3

7. Combs E, Bernardo R. Accuracy of Genomewide Selection for Different Traits with Constant Population Size, Heritability, and Number of Markers. Plant Genome. 2013;6(1). doi: 10.3835/plantgenome2012.11.0030

8. Heslot N, Feoktistov V. Optimization of selective phenotyping and population design for genomic prediction. J Agric Biol Environ Stat. 2020;25(4):579-600.
9. Zhong S, Jannink JL. Using Quantitative Trait Loci Results to Discriminate Among Crosses on the Basis of Their Progeny Mean and Variance. Genetics. 2007;177(1):567-76. doi: 10.1534/genetics.107.075358

10. Akdemir D, Sánchez JI. Efficient Breeding by Genomic Mating. Front Genet. 2016;7:210. doi: 10.3389/fgene.2016.00210

11. Lehermeier C, Teyssèdre S, Schön CC. Genetic Gain Increases by Applying the Usefulness Criterion with Improved Variance Prediction in Selection of Crosses. Genetics. 2017;207(4):1651-61. doi: 10.1534/genetics.117.300403

12. Allier A, Moreau L, Charcosset A, Teyssèdre S, Lehermeier C. Usefulness Criterion and Post-selection Parental Contributions in Multi-parental Crosses: Application to Polygenic Trait Introgression. G3. 2019;9(5):1469-79. doi: 10.1534/g3.119.400129

13. Bassi FM, Bentley AR, Charmet G, Ortiz R, Crossa J. Breeding schemes for the implementation of genomic selection in wheat (Triticum spp.). Plant Sci. 2016;242:23-36. doi: 10.1016/j.plantsci.2015.08.021

14. Daetwyler HD, Villanueva B, Woolliams JA. Accuracy of Predicting the Genetic Risk of Disease Using a Genome-Wide Approach. PLoS One. 2008;3(10):e3395. doi: 10.1371/journal.pone.0003395

15. Jannink JL, Lorenz AJ, Iwata H. Genomic selection in plant breeding: from theory to practice. Brief Funct Genomics. 2010;9(2):166-77. doi: 10.1093/bfgp/elq001

16. Lorenz AJ, Chao S, Asoro FG, Heffner EL, Hayashi T, Iwata H, et al. Chapter Two—Genomic Selection in Plant Breeding: Knowledge and Prospects. Adv Agron. 2011;110:77-123. doi: 10.1016/B978-0-12-385531-2.00002-5

17. Habier D, Fernando RL, Dekkers JCM. The Impact of Genetic Relationship Information on Genome-Assisted Breeding Values. Genetics. 2007;177(4):2389-97. doi: 10.1534/genetics.107.081190

18. Charmet G, Storlie E, Oury FX, Laurent V, Beghin D, Chevarin L, et al. Genome-wide prediction of three important traits in bread wheat. Mol Breeding. 2014;34(4):1843-52. doi: 10.1007/s11032-014-0143-y

19. Crossa J, Pérez P, Hickey J, Burgueño J, Ornella L, Cerón-Rojas J, et al. Genomic prediction in CIMMYT maize and wheat breeding programs. Heredity. 2014;112(1):48-60. doi: 10.1038/hdy.2013.16

20. Isidro J, Jannink JL, Akdemir D, Poland J, Heslot N, Sorrells ME. Training set optimization under population structure in genomic selection. Theor Appl Genet. 2015;128(1):145-58. doi: 10.1007/s00122-014-2418-4

21. Rincent R, Laloë D, Nicolas S, Altmann T, Brunel D, Revilla P, et al. Maximizing the Reliability of Genomic Selection by Optimizing the Calibration Set of Reference Individuals: Comparison of Methods in Two Diverse Groups of Maize Inbreds (Zea mays L.). Genetics. 2012;192(2):715-28. doi: 10.1534/genetics.112.141473

22. Rincent R, Charcosset A, Moreau L. Predicting genomic selection efficiency to optimize calibration set and to assess prediction accuracy in highly structured populations. Theor Appl Genet. 2017;130(11):2231-47. doi: 10.1007/s00122-017-2956-7
23. Zaïm M, Kabbaj H, Kehel Z, Gorjanc G, Filali-Maltouf A, Belkadi B, et al. Combining QTL analysis and genomic predictions for four durum wheat populations under drought conditions. Front Genet. 2020;11:316. doi: 10.3389/fgene.2020.00316

24. Calus MP, Veerkamp RF. Accuracy of multi-trait genomic selection using different methods. Genet Sel Evol. 2011;43(1):26. doi: 10.1186/1297-9686-43-26

25. Solberg TR, Sonesson AK, Woollams JA, Meuwissen THE. Genomic selection using different marker types and densities. J Anim Sci. 2008;86(10):2447-54. doi: 10.2527/jas.2007-0010

26. Roos A, Schrooten C, Mullaart E, Beek SVD, Jong GD and Voskamp W. Genomic selection at CRV. Available from: https://www.researchgate.net/profile/Erik-Mullaart/publication/242042641_Genomic_selection_at_CRV/links/0c96052a8205001611000000/Genomic-selection-at-CRV.pdf. Accessed 2021 Sept 24.

27. Desta ZA, Ortiz R. Genomic selection: genome-wide prediction in plant improvement. Trends Plant Sci. 2014;19(9):592-601. doi: 10.1016/j.tplants.2014.05.006

28. Poland J, Endelman J, Dawson J, Rutkoski J, Wu S, Manes Y, et al. Genomic Selection in Wheat Breeding using Genotyping-by-Sequencing. Plant Genome. 2012;5(3):103-13. doi: 10.3835/plantgenome2012.06.0006

29. Lado B, Matus I, Rodríguez A, Inostroza L, Poland J, Belzile F, et al. Increased Genomic Prediction Accuracy in Wheat Breeding Through Spatial Adjustment of Field Trial Data. G3. 2013;3(12):2105-14. doi: 10.1534/g3.113.007807

30. Storlie E, Charmet G. Genomic Selection Accuracy using Historical Data Generated in a Wheat Breeding Program. Plant Genome. 2013;6(1). doi: 10.3835/plantgenome2013.01.0001

31. Zhao Y, Mette MF, Reif JC. Genomic selection in hybrid breeding. Plant Breed. 2015;134(1):1-10. doi: 10.1111/pbr.12231

32. Norman A, Taylor J, Tanaka E, Telfer P, Edwards J, Martinant JP, et al. Increased genomic prediction accuracy in wheat breeding using a large Australian panel. Theor Appl Genet. 2017;130(12):2543-55. doi: 10.1007/s00122-017-2975-4

33. Battenfield SD, Guzmán C, Gaynor RC, Singh RP, Peña RJ, Dreisigacker S, et al. Genomic Selection for Processing and End-Use Quality Traits in the CIMMYT Spring Bread Wheat Breeding Program. Plant Genome. 2016;9(2). doi: 10.3835/plantgenome2016.01.0005

34. Guzman C, Peña RJ, Singh R, Autrique E, Dreisigacker S, Crossa J, et al. Wheat quality improvement at CIMMYT and the use of genomic selection on it. Appl Transl Genomic. 2016;11:3-8. doi: 10.1016/j.atg.2016.10.004

35. Liu G, Zhao Y, Gowda M, Longin CFH, Reif JC, Mette MF. Predicting Hybrid Performances for Quality Traits through Genomic-Assisted Approaches in Central European Wheat. PLoS One. 2016;11(7):e0158635. doi: 10.1371/journal.pone.0158635

36. Hayes BJ, Panozzo J, Walker CK, Choy AL, Kant S, Wong D, et al. Accelerating wheat breeding for end-use quality with multi-trait genomic predictions incorporating near infrared and nuclear magnetic resonance-derived
phenotypes. *Theor Appl Genet.* 2017;130(12):2505-19. doi: 10.1007/s00122-017-2972-7

37. Lado B, Vázquez D, Quincke M, Silva P, Aguilar I, Gutiérrez L. Resource allocation optimization with multi-trait genomic prediction for bread wheat (*Triticum aestivum* L.) baking quality. *Theor Appl Genet.* 2018;131(12):2719-31. doi: 10.1007/s00122-018-3186-3

38. Michel S, Kummer C, Gallee M, Hellinger J, Ametz C, Akgöl B, et al. Improving the baking quality of bread wheat by genomic selection in early generations. *Theor Appl Genet.* 2018;131(2):477-93. doi: 10.1007/s00122-017-2998-x

39. Ornella L, Singh S, Perez P, Burgueño J, Singh R, Tapia E, et al. Genomic Prediction of Genetic Values for Resistance to Wheat Rusts. Plant Genome. 2012;5(3). doi: 10.3835/plantgenome2012.07.0017

40. Rutkoski J, Benson J, Jia Y, Brown-Guedira G, Jannink JL, Sorrells M. Evaluation of Genomic Prediction Methods for Fusarium Head Blight Resistance in Wheat. Plant Genome. 2012;5(2). doi: 10.3835/plantgenome2012.02.0001

41. Daetwyler HD, Bansal U, Bariana H, Hayden M, Hayes B. Genomic prediction for rust resistance in diverse wheat landraces. *Theor Appl Genet.* 2014;127:1795-803. doi: 10.1007/s00122-014-2341-8

42. Arruda MP, Brown PJ, Lipka AE, Krill AM, Thurber C, Kolb FL. Genomic Selection for Predicting Fusarium Head Blight Resistance in a Wheat Breeding Program. Plant Genome. 2015;8(3). doi: 10.3835/plantgenome2015.01.0003

43. Rutkoski JE, Poland JA, Singh RP, Huerta-Espino J, Bhavani S, Barbier H, et al. Genomic Selection for Quantitative Adult Plant Stem Rust Resistance in Wheat. Plant Genome. 2014;7(3). doi: 10.3835/plantgenome2014.02.0006

44. Longin CFH, Mi X, Würschum T. Genomic selection in wheat: optimum allocation of test resources and comparison of breeding strategies for line and hybrid breeding. *Theor Appl Genet.* 2015;128(7):1297-306. doi: 10.1007/s00122-015-2505-1

45. Riedelsheimer C, Melchinger AE. Optimizing the allocation of resources for genomic selection in one breeding cycle. *Theor Appl Genet.* 2013;126(11):2835-48. doi: 10.1007/s00122-013-2175-9

46. Gorjanc G, Jenko J, Hearne SJ, Hickey JM. Initiating maize pre-breeding programs using genomic selection to harness polygenic variation from landrace populations. BMC Genom. 2016;17(1):30. doi: 10.1186/s12864-015-2345-z

47. Gorjanc G, Hickey JM. AlphaMate: a program for optimizing selection, maintenance of diversity and mate allocation in breeding programs. *Bioinformatics.* 2018;34(19):3408-11. doi: 10.1093/bioinformatics/bty375

48. Akdemir D, Beavis W, Fritsche-Neto R, Singh AK, Isidro-Sánchez J. Multi-objective optimized genomic breeding strategies for sustainable food improvement. Heredity. 2019;122(5):672-83. doi: 10.1038/s41437-018-0147-1

49. Allier A, Lehermeier C, Charcosset A, Moreau L, Teyssèdre S. Improving Short- and Long-Term Genetic Gain by Accounting for Within-Family Variance in Optimal Cross-Selection. *Front Genet.* 2019;10:1006. doi: 10.3389/fgene.2019.01006
50. Allier A, Teyssèdre S, Lehermeier C, Charcosset A, Moreau L. Genomic prediction with a maize collaborative panel: identification of genetic resources to enrich elite breeding programs. Theor Appl Genet. 2020;133(1):201-15. doi: 10.1007/s00122-019-03451-9

51. Rimbert H, Darrier B, Navarro J, Kitt J, Choulet F, Leveugle M, et al. High throughput SNP discovery and genotyping in hexaploid wheat. PLoS One. 2018;13(1):e0186329. doi: 10.1371/journal.pone.0186329

52. Oury FX, Heumez E, Rolland B, Auzanneau J, Bérard P, Brancourt-Hulmel M, et al. Winter wheat (Triticum aestivum L) phenotypic data from the multiannual, multilocal field trials of the INRA Small Grain Cereals Network. Available from: https://data.inrae.fr/dataset.xhtml?persistentId=doi:10.15454/1.4489666216568333E12. Accessed 2021 Sept 24.

53. Daetwyler HD, Hayden MJ, Spangenberg GC, Hayes BJ. Selection on Optimal Haploid Value Increases Genetic Gain and Preserves More Genetic Diversity Relative to Genomic Selection. Genetics. 2015;200(4):1341-8. doi: 10.1534/genetics.115.178038

54. Heffner EL, Lorenz AJ, Jannink JL, Sorrells ME. Plant Breeding with Genomic Selection: Gain per Unit Time and Cost. Crop Sci. 2010;50(5):1681-90. doi: 10.2135/cropsci2009.11.0662

55. Bernardo R, Charcosset A. Usefulness of Gene Information in Marker-Assisted Recurrent Selection: A Simulation Appraisal. Crop Sci. 2006;46(2):614-21. doi: 10.2135/cropsci2005.05-0088

56. Lin Z, Cogan NOI, Pembleton LW, Spangenberg GC, Forster JW, Hayes BJ, et al. Genetic Gain and Inbreeding from Genomic Selection in a Simulated Commercial Breeding Program for Perennial Ryegrass. Plant Genome. 2016;9(1). doi: 10.3835/plantgenome2015.06.0046

57. Maluszynski M, Kasha K, Forster BP, Szarejko I. Doubled Haploid Production in Crop Plants: A Manual. New York(US): Springer Science & Business Media; 2013.

58. Ben-Sadoun S, Rincent R, Auzanneau J, Oury FX, Rolland B, Heumez E, et al. Economical optimization of a breeding scheme by selective phenotyping of the calibration set in a multi-trait context: application to bread making quality. Theor Appl Genet. 2020;133(7):2197-212. doi: 10.1007/s00122-020-03590-4

59. Rutkoski J, Poland J, Mondal S, Autrique E, Pérez LG, Crossa J, et al. Canopy Temperature and Vegetation Indices from High-Throughput Phenotyping Improve Accuracy of Pedigree and Genomic Selection for Grain Yield in Wheat. G3. 2016;6(9):2799-808. doi: 10.1534/g3.116.032888

60. Sun J, Rutkoski JE, Poland JA, Crossa J, Jannink JL, Sorrells ME. Multitrait, Random Regression, or Simple Repeatability Model in High-Throughput Phenotyping Data Improve Genomic Prediction for Wheat Grain Yield. Plant Genome. 2017;10(2). doi: 10.3835/plantgenome2016.11.0111

61. Crain J, Mondal S, Rutkoski J, Singh RP, Poland J. Combining High-Throughput Phenotyping and Genomic Information to Increase Prediction and Selection Accuracy in Wheat Breeding. Plant Genome. 2018;11(1). doi: 10.3835/plantgenome2017.05.0043
62. Jia Y, Jannink JL. Multiple-Trait Genomic Selection Methods Increase Genetic Value Prediction Accuracy. Genetics. 2012;192(4):1513-22. doi: 10.1534/genetics.112.144246

63. Müller D, Schopp P, Melchinger AE. Selection on Expected Maximum Haploid Breeding Values Can Increase Genetic Gain in Recurrent Genomic Selection. G3. 2018;8(4):1173-81. doi: 10.1534/g3.118.200091

64. Gaynor RC, Gorjanc G, Bentley AR, Ober ES, Howell P, Jackson R, et al. A Two-Part Strategy for Using Genomic Selection to Develop Inbred Lines. Crop Sci. 2017;57(5):2372-86. doi: 10.2135/cropsci2016.09.0742

65. Hazel LN. The Genetic Basis for Constructing Selection Indexes. Genetics. 1943;28(6):476-90.

How to cite this article:

Ben-Sadoun S, Fugeray-Scarbel A, Auzanneau J, Charmet G, Lemarié S, Bouchet S. Integration of genomic selection into winter-type bread wheat breeding schemes: a simulation pipeline including economic constraints. Crop Breed Genet Genom. 2021;3(4):e210008. https://doi.org/10.20900/cbgg20210008