A look-ahead approach to maximizing present value of genetic gains in genomic selection

Abstract Look-ahead selection is a sophisticated yet effective algorithm for genomic selection, which optimizes not only the selection of breeding parents but also mating strategy and resource allocation by anticipating the implications of crosses in a prespecified future target generation. Simulation results using maize datasets have suggested that look-ahead selection is able to significantly accelerate genetic gain in the target generation while maintaining genetic diversity. In this paper, we propose a new algorithm to address the limitations of look-ahead selection, including the difficulty in specifying a meaningful deadline in a continuous breeding process and slow growth of genetic gain in early generations. This new algorithm uses the present value of genetic gains as the breeding objective, converting genetic gains realized in different generations to the current generation using a discount rate, similar to using the interest rate to measure the time value of cash flows incurred at different time points. By using the look-ahead techniques to anticipate the future gametes and thus present value of future genetic gains, this algorithm yields a better trade-off between short-term and long-term benefits. Results from simulation experiments showed that the new algorithm can achieve higher genetic gains in early generations and a continuously growing trajectory as opposed to the look-ahead selection algorithm, which features a slow progress in early generations and a growth spike right before the deadline.


Introduction
The core of plant breeding is in the selection of breeding parents to improve traits of interest, such as yield, tolerance to environmental stress, resistance to pests, most of which are quantitatively inherited (Wricke and Weber 2010). Traditional selection strategies used to focus on observable phenotype or a handful of assisting markers related to the desired traits. However, these methods are not applicable to polygenic traits, i.e. traits consisting of many small-effect alleles for which effects are scattered and difficult to determine reliably . With the development of high-throughput genotyping and single nucleotide polymorphism (SNP) effect estimation (Heffner et al. 2009;Jannink et al. 2010) shedding lights on the field, genomic prediction emerged as a technique for linking genomic information of quantitative traits to phenotypic values. Genomic selection (GS) is an approach to exploiting genomic markers to cater to novel breeding programs and evaluation. In this technique, the genomic estimated breeding value (GEBV), i.e. the sum of the estimated marker effects for a specific individual becomes a popular criterion to evaluate the breeding potential for certain traits without relying on individual's phenotype (Goddard 2009). Improvements in genomic prediction models for complex patterns, such as genotype by environment interactions (G Â E), have received great attention. Much of the work on GS has been on the design and execution of field trials (Heslot et al. 2015;Crossa et al. 2017;Hickey et al. 2017). What appears to be missing in the body of literature is computational algorithms that use genomic prediction to intelligently select individuals or groups worthy of breeding. These algorithms should be able to not only select optimal breeding parents but also strategically determine the number of crosses to make and progenies to produce under resource and time constraints.
Conventional genomic selection (CGS) selects the individuals with the highest GEBVs as breeding parents (Meuwissen et al. 2001), which are assumed to be most likely to produce superior offspring. CGS has been widely adopted in both plant and animal breeding practices due to its simplicity and effectiveness. However, this truncation approach often leads to loss of genetic diversity after only a few breeding cycles. Several more sophisticated selection algorithms have been proposed to address such limitations. Weighted genomic selection modifies the SNP effects with weights to render the inheritance of favorable alleles with low population frequencies (Heffner et al. 2010). Optimal haploid value (OHV) evaluates a breeding parent not by its own genetic value but by the genetic value of the best gamete that it can produce in the immediate next generation (Daetwyler et al. 2015). OHV also aggregates adjacent SNP markers as recombination blocks distributed across chromosomes to accelerate the computation. Optimal population value (OPV) introduces the concept of group selection and selects a group of breeding parents that possess favorable alleles in complementary loci, thus can produce the best progeny in the long term (Goiffon et al. 2017). Look-ahead selection (LAS) extends the concept of OHV and OPV by anticipating the implications of crosses in the current generation on the progeny produced in a prespecific future generation. By selecting optimal crosses to maximize the long-term performance, LAS has been found to be able to not only accelerate genetic gains but also preserve genetic diversity (Moeinizade et al. 2019).
Building on the prior work, we propose an algorithm to overcome 2 limitations of LAS, which are the challenges to specify an appropriate deadline in the context of continuous genetic improvement and the weak performance before the deadline. The new algorithm borrows the concept of present value (PV) from the field of finance and uses it to define a new breeding objective, which converts genetic gains in different future generations over a planning window back to the current time using a discount rate, and the trade-off between short-term and long-term performances can be adjusted using the parameters of window length and discount rate. Details of this method are developed in the next section.

Nomenclature
Here, we define some of the parameters and variables used for modeling GS.

Formulation of GS
GS benefits from high-density markers used in whole-genome prediction models, which pave the way for effect estimation of quantitative trait locus for the traits of interest. Here, we use SNPs data, which are a common type of genetic variation among the population. With assumptions of linear additive SNP effects and appropriate high-dimensional point estimation methods, one can model the quantitative relationship between b and the GEBV of individual i as v i ¼ l þ P l b l P m G l;m;i ; 8i 2 f1; . . . ; Ng, where l is overall mean (Desta and Ortiz 2014).
With the aforesaid definitions, the goal of GS is to select S pairs of breeding parents to achieve specific breeding goals. The general program can be formulated as max x f ðx; GÞ (1) s:t: x i 2 f0; 1g i 2 f1; . . . ; Ng: Here, • Objective function f ðÁÞ represents the genetic breeding value or other appropriately defined breeding objectives. For example, CGS uses the total GEBVs of all selected breeding parents as its objective function: f ðx; GÞ ¼ P i x i v i . • Decision variables x ¼ ðx i Þ for all i 2 f1; . . . ; Ng are a binary variable that indicates whether individual i is selected ðx i ¼ 1Þ or not ðx i ¼ 0Þ, as shown in constraint (3).

Review of the look-ahead selection algorithm
LAS presents an efficient framework for selecting breeding parents that maximize genetic gain of future progeny in a userdefined deadline generation. It takes into account not only parent selection but also mating, time management, and resource allocation (number of progeny from each cross; Moeinizade et al. 2019). The formulation of LAS is shown as follows s:t: Pr½gðx; Y; G; b; r; T À tÞ 5 u 5 1 À c x i ; y i;j 2 f0; 1g i; j 2 f1; . . . ; Ng: Here, • decision variable x i indicates whether individual i is selected as a breeding parent (x i ¼ 1) or not (x i ¼ 0); • decision variable y i;j indicates whether individual i is mated with j (y i;j ¼ 1) or not (y i;j ¼ 0); • parameter t 2 f1; . . . ; Tg represents the current generation number; • u is the GEBV of the best progeny in the final generation that has a probability of occurrence of at least 1 À c; • parameter c is a risk tolerance parameter, a larger (smaller) value of which will incentivize the model to maximize the performance in more (less) optimistic scenarios; and • gðÁÞ is the GEBV of a random progeny in the final generation T, created using the breeding decisions x and Y through a lookahead simulation proposed in Moeinizade et al. (2020).
The introduction of the additional decision variable Y allows the model to further accelerate genetic responses by optimizing mating strategies of the selected individuals (Toro and Varona 2010;Wang et al. 2018). Constraint (5) defines u, which is the c quantile of g, the GEBV of a random progeny in the final generation T. This constraint also interprets the objective function (4), which is to maximize the genetic gain in the best possible scenario with a probability of at least 1 À c. Constraint (6) ensures that only selected individuals be mated with each other. Constraints (7) and (8) make sure that a total of S crosses are made and that the mating is symmetric. Constraint (9) requires all decision variables to be binary.

N
Number of the individuals in a population, a scalar L Number of SNPs of an individual, a scalar T Number of generations in a breeding program, a scalar S Number of breeding parents to be selected, a scalar G Genotype of a population, a binary matrix G 2 B LÂ2ÂN , with element G l;m;i indicating whether the allele in the first (m ¼ 1) or second (m ¼ 2) chromosome of diploid individual i at locus l is a major allele (G l;m;i ¼ 1) or a minor allele (G l;m;i ¼ 0) b SNP effect, a vector b 2 R L , with b l being the allele effect for locus l r Recombination frequencies, a vector r 2 R LÀ1 , with r l being the recombination frequency between loci l and l þ 1 v GEBVs, a vector v 2 R N , with v i being the GEBV of individual i

Motivation for improvement
Despite the effectiveness of the LAS algorithm in accelerating genetic gain while preserving genetic diversity, it has 2 major limitations that motivated the design of the proposed algorithm. First, it is challenging to determine an appropriate value for the breeding deadline, since the goal of breeding projects is usually to achieve continuous genetic improvements. To address this limitation, we replace the fixed deadline with a rolling horizon, which restarts its planning timeline over an adjustable time interval. The interval acts as a sliding window to smooth the variability in future offspring and consider the dependence along the breeding process.
Second, the LAS algorithm focuses exclusively on maximizing genetic gain in the terminal generation without considering earlier performance, potentially resulting in unacceptably low short-term genetic gains. In the proposed model, we redefine the breeding objective as the PV of genetic gains over the rolling horizon, which takes a user-defined discount rate, k, accounting for the time value of genetic gains and the market value of successful release of new commercial lines from the breeding program. A genetic gain achieved in generation s would be only 1 ð1þkÞ s times as valuable as the same genetic gain would have in generation 0. A larger k assigns a higher time value of genetic gain, putting higher weights on shorter-term performances.
As a result, the new method is expected to achieve higher genetic gains in earlier generations albeit at the cost of a weaker performance in the final generation, as illustrated in Fig. 1. This is achieved by maximizing the PV of genetic gains in all generations, with earlier gains having higher values than later ones, similar with the financial concept of "time value of money." Similar with the compound interest rate, the discount rate can be used to adjust the trade-off between short-term and longer-term genetic gains.

PV-based look-ahead selection algorithm
We propose to use the PV of GEBVs as the new objective for maximization. In finance (Weitzman 1998; Zi zlavskỳ 2014), for a series of cash flows f s over certain time period 8s 2 f1; . . . ; Tg, the PV of these cash flows is the summation of their discounted values at time 0: P T s¼1 fs ð1þkÞ s , where k is the compound interest rate, indicating the time value of money. In the context of GS, if we use w s to denote the GEBV of a progeny in the sth future generation, W the number of generations that we look ahead to, and k the discount rate that indicates the "time value of genetic gains," then the new objective becomes P W s¼1 w s ð1þkÞ s . As such, the proposed method, which we refer to as PV-LAS, can be formulated as follows: s:t: Constraints ð6Þ; ð7Þ; ð8Þ; ð9Þ (11) Pr½gðx; Y; G; b; r; sÞ 5 w s 5 1 À c; s 2 f1; . . . ; Wg: Here, • parameter W is the length of the sliding window, indicating the number of generations to look ahead; • parameter k is the discount rate, indicating the time value of genetic gains. Specifically, the genetic gain w s in s generations is as valuable as w s ð1þkÞ s in the current generation; and • w s is the c quantile of a random progeny' GEBV in the sth future generation given the current decision variables x and Y.
The model introduces 2 additional parameters: the length of the sliding window W and the discount rate k. The window length defines the longest term that the breeders look ahead in the planning horizon, whereas the discount rate determines the trade-off between short-term and long-term genetic gains. A larger (smaller) k places a higher emphasis on shorter-term (longerterm) performances.
PV-LAS can be seen as an extension of LAS with a modified approach to looking ahead. LAS maximizes the performance at a predefined target deadline, which gradually reduces the planning horizon as the breeding process progresses toward the deadline. In contrast, a moving horizon is used in PV-LAS, where the planning horizon is always the length of the sliding window, W. As such, PV-LAS is more applicable to breeding programs with goals for continuous genetic improvement. The change in planning horizon requires a different formula to calculate the recombination frequency in a future generation. In particular, equation (13)  For all l 2 f1; . . . ; L À 1g . . . ; Wg:

Optimization framework
We propose a heuristic search approach to find optimal selection and mating decisions for PV-LAS by iteratively searching the solution space. The workflow of this heuristic algorithm is described as follows: • Input: G, b, r, W, k, S, K, and c. . If ðx;ŶÞ is feasible and has a better objective value than the incumbent, then update the incumbent solution ðx Ã; Y Ã Þ ðx;ŶÞ and its corresponding objective value f Ã fðx;ŶÞ. Repeat this step until no further improvement of the incumbent can be made.
• Output: the locally optimal x and Y.
From our computational experiences, the LAS algorithm almost always found selfing as the optimal breeding strategy for the last generation. Our explanation for this observation is that when the breeding goal is to maximize the genetic gain in the immediate next generation, the value of a breeding parent can be largely determined by its best gamete that can be produced (within certain risk tolerance). As such, selfing the top #1 breeding parent would be more likely to produce a better progeny than crossing the top #1 with top #2. While this strategy has produced satisfactory results for the case study, it is unnecessarily the optimal strategy for all breeding programs.

Illustrative example
We use a toy example to illustrate the difference between LAS and PV-LAS. Suppose we aim to make S ¼ 3 crosses from a group of N ¼ 8 diploid individuals, each genotyped with L ¼ 10 SNPs. The breeding deadline in LAS and the length of sliding window in PV-LAS were both set as T ¼ W ¼ 3. The arbitrarily simulated SNP effect b, recombination frequency and finalized selection and mating results for LAS and PV-LAS are showed in Fig. 2. For each of these 3 future generations, K ¼ 500 future gametes have been simulated. We used c ¼ 0:8. Figure 2 illustrated the solutions given by the 2 methods. Both methods included a cross between the second and the eighth individuals from the left, whereas the remaining 2 crosses were different. Figure 3 shows vertical histograms of GEBVs of 500 random progeny in 3 generations using LAS and PV-LAS, which demonstrated the major differences of these algorithms: LAS resulted in higher genetic gains in the final generation, whereas PV-LAS improved the growth in the first 2 generations with a slightly compromised performance in the third.

Data sources
We conducted computational experiments using the same data set as Moeinizade et al. (2019): the genotypic data G contains genotypes of 369 maize inbred lines consisting of L ¼ 1,406,757 SNPs distributed across 10 maize chromosomes (Goiffon et al. 2017). These SNPs were aggregated into 10,000 blocks to increase the computational speed. SNP effects were estimated on the basis of 369 shoot apical meristem phenotypes using the BayesB model (Meuwissen et al. 2001;Leiboff et al. 2015). Recombination rates were based on the genetic map developed from maize nested association mapping (Yu et al. 2008). We assume that the marker effectsb and recombination ratesr have been estimated reasonably accurately and they stay fixed in our simulations. As a caveat for this simplifying assumption, however, we point out that LAS and PV-LAS are more sensitive to the accuracy of allele effects and recombination frequencies than CGS, because the former 2 extract more information from such data to make more sophisticated decisions.

Simulation setting
The aim of the simulation was to evaluate the performance of 3 algorithms: CGS, LAS, and PV-LAS with respect to genetic improvement throughout multiple breeding generations. For a group of offspring G with a size of N in at the end of the breeding program, the following criteria were used for evaluating the 3 algorithms: • Mean of GEBV: 1 N P i P l P m ðG l;m;i b l Þ. This criterion measures the average performance of genetic gains.
• Lower potential of GEBV: 2 P l min m;i ðG l;m;i b l Þ. This criterion gives the theoretical lower bound of GEBV based on the remaining genetic diversity.
• Upper potential of GEBV: 2 P l max m;i ðG l;m;i b l Þ. This criterion gives the theoretical upper bound of GEBV based on the remaining genetic diversity.
• PV of GEBV: s . This criterion measures the PV of GEBVs over a period of T generations for a given discount rate k.
Each of the simulations consisted of 3 steps: initialization, selection, and reproduction. In the initialization step, N ¼ 200 individuals were randomly selected among the given 369 maize inbred lines. In the selection step, GS algorithms were used to select the optimal crosses. In the reproduction step, crosses from the previous step are made, each producing N ¼ 200 progeny, and then it goes back to the selection step for the next generation until the final generation T ¼ 10. Each algorithm was tested for 500 independent simulations. Parameters k ¼ 0:1 and c ¼ 0:8 were used for all experiments.

Performance comparison
Results for genetic gain comparison are shown in Fig. 4. The GEBV mean in Fig. 4b matched our expected trajectory in Fig. 1.  Figure 4, a and c showed how PV-LAS compromised genetic diversity, with respect to LAS, to achieve higher growth in early generations, although PV-LAS still largely outperformed CGS in both genetic gain and genetic diversity. Figure 4, d-f showed mean, lower potential, and higher potential of GEBV together for the 3 selection methods. Although LAS and PV-LAS use a more forward-looking selection strategy than CGS, they have narrower ranges between upper and lower bounds in the first 2 generations. This counter-intuitive result is because they only select parents whose favorable alleles can be aggregated within time and resource constraints, which means that some otherwise high-performing parents may be discarded if their favorable alleles require more time or resources than available to be integrated with the selected ones. On the other hand, CGS will select all high-performing parents without anticipation of their future progeny, which may lead to higher genetic diversity in the first generations but inevitable loss of genetic gain and diversity in the longer term. Figure 5 compares the empirical cumulative distribution functions of the PVs of GEBVs resulted from different selection methods. We observe that PV-LAS is on the right-hand side of LAS in  almost all quantiles, which indicates the stochastic dominance of PV-LAS over LAS in terms of the PV of genetic gains. This was expected because PV-LAS was designed to optimize the PV of GEBVs.

Sensitivity analysis on PV-LAS parameters
Sensitivity analysis was performed to test the influence of window length W and discount rate k on the performance of PV-LAS. Figure 6 shows the average (over 500 independent experiments) differences in genetic gains between the benchmark case of W ¼ 1 and other window lengths. When W ¼ 1, the model focused on genetic gain in the immediate next generation; as a result, GEBVs jumped to a plateau after the first selection and lost diversity and potential for future growth. The figure shows that the effect of W is not monotonic and that a balance between shortterm and long-term growth requires a window length that is neither too short nor too long. For this particular case study, W ¼ 3 achieved the best performance between generations 6 and 10, but it is unnecessarily optimal for other studies or datasets. In general, the sensitivity analysis should be done for each new dataset to identify the best set of parameters. Figure 7 shows the differences in genetic gains between the benchmark case of k ¼ 0 and other discount rates, with window length fixed at W ¼ 3. When k ¼ 0, the model focused on the nominal genetic gain and ignored its time value. Larger (smaller) k values put higher (lower) emphasis on the time value of genetic gain and led to higher (lower) growth in the short-term and weaker (stronger) performance in later generations.

Conclusion
The introduction of the PV concept to GS is central to this work. PV-LAS uses the PV of GEBVs over a certain window period as the breeding objective by discounting genetic gain in future generations back to the current time. As such it balances the short-term and long-term benefits of genetic growth and provides a  continuous growth trajectory. At the same time, PV-LAS makes a moderate compromise on genetic diversity in order to achieve a higher growth in early generations.
Computational results demonstrated the effectiveness of the new algorithm. Optimal window length W and discount rate k for specific data sets can be determined using sensitive analyses.
Several research directions are worthy of future investigation. First, the simulations in this paper were based on the assumption that the estimates of additive effects of SNPs and recombination frequencies are reasonably accurate. Analysis should be conducted to assess the effects of inaccurate estimations and how to mitigate such effects. Second, we can extend PV-LAS to accommodate both additive and nonadditive effects, such as   7. Differences in genetic gains between the benchmark case of k ¼ 0 and other discount rates for PV-LAS, averaged over 500 independent experiments. Window length is fixed at W ¼ 3.
dominance and epistasis effects. Third, PV-LAS does not explicitly address G Â E, which can often be substantial and create problems in finding consistently superior genotypes, leading to reduced heritability and overall genetic gain. Future research should focus on the convergence of crop modeling and machine learning approaches to explore more advanced strategies to address G Â E in the breeding process.

Data availability
The datasets used the computational experiments were derived from sources in the public domain as described in Data Sources.