The prediction of phenotypic traits using high-density genomic data has many applications such as the selection of plants and animals of commercial interest; and it is expected to play an increasing role in medical diagnostics. Statistical models used for this task are usually tested using cross-validation, which implicitly assumes that new individuals (whose phenotypes we would like to predict) originate from the same population the genomic prediction model is trained on. In this paper we propose an approach based on clustering and resampling to investigate the effect of increasing genetic distance between training and target populations when predicting quantitative traits. This is important for plant and animal genetics, where genomic selection programs rely on the precision of predictions in future rounds of breeding. Therefore, estimating how quickly predictive accuracy decays is important in deciding which training population to use and how often the model has to be recalibrated. We find that the correlation between true and predicted values decays approximately linearly with respect to either F ST or mean kinship between the training and the target populations. We illustrate this relationship using simulations and a collection of data sets from mice, wheat and human genetics.

The availability of increasing amounts of genomic data is making the use of statistical models to predict traits of interest a mainstay of many applications in life sciences. Applications range from medical diagnostics for common and rare diseases to breeding characteristics such as disease resistance in plants and animals of commercial interest. We explored an implicit assumption of how such prediction models are often assessed: that the individuals whose traits we would like to predict originate from the same population as those that are used to train the models. This is commonly not the case, especially in the case of plants and animals that are parts of selection programs. To study this problem we proposed a model-agnostic approach to infer the accuracy of prediction models as a function of two common measures of genetic distance. Using data from plant, animal and human genetics, we find that accuracy decays approximately linearly in either of those measures. Quantifying this decay has fundamental applications in all branches of genetics, as it measures how studies generalise to different populations.

This approach is valuable in addressing several key questions in the implementation of genomic selection programs, such as: How often (e.g., in terms of future generations) will the genomic prediction model have to be re-estimated to maintain a minimum required accuracy in the predictions of the phenotypes? How should we structure our training population to maximise that accuracy? Which new, distantly related individuals would be beneficial to introduce in a selection program for the purpose of maintaining a sufficient level of genetic variability?

In this paper we will investigate the relationship between genetic distance and predictive accuracy in the prediction of quantitative traits. We will simulate training and target samples with varying genetic distances by splitting the training population into a sequence of pairs of subsets with increasing genetic differentiation. We will measure predictive accuracy with Pearson’s correlation, which we will estimate by performing genomic prediction from one subset to the other in each pair. Among various measures of relatedness available in the literature, we will consider mean kinship and F ST , although we will only focus on the latter. We will then study the mean Pearson’s correlation as a function of genetic distance, which we will refer to as the “decay curve” of the former over the latter.

Most of these issues have been explored in the literature, and have been tackled in various ways either from a methodological perspective or by producing larger data sets and more accurate phenotyping. However, the extent to which predictive models generalise from the populations used to train them to distantly related target populations appears not to have been widely investigated (two exceptions are [ 7 , 13 ]). The accuracy of prediction models is often evaluated in a general setting using cross-validation with random splits, which implicitly assumes that test individuals are drawn from the same population as the training sample; in that case accuracy to predict phenotypes is only bounded by heritability, although unaccounted “missing heritability” is common [ 14 , 15 ]. However, this assumption is violated in many practical applications, such as genomic selection, that require predictions of individuals that are genetically distinct from the training sample: for instance, causal variants may differ in both frequency and effect size between different ancestry groups (in humans, e.g. [ 16 ] for lactose persistence), subspecies (in plants and animals, e.g. [ 17 ] for rice) or even families [ 18 ]. In such cases cross-validation with random splits may overestimate predictive accuracy due to the mismatch between model validation and the prediction problem of interest [ 19 , 20 ] even when population structure is taken into account [ 21 ]. The more distantly the target population is related to the training population, the lower the average predictive accuracy of a genomic model; this has been demonstrated on both simulated and real dairy cattle data [ 20 , 22 , 23 ].

For the WHEAT data, we construct decay curves for grain yield, height, flowering time and grain protein content using the French wheat varieties as the training population. UK and German varieties are the target populations, for which we estimate . Furthermore, we also construct a second decay curve for yield using the varieties registered before 1990 as the training population, as in [ 4 ]. Varieties registered between 1990 and 1999, and those registered after 2000, are used as target populations.

Finally, we estimate the decay curves for some of the phenotypes available in the WHEAT and MICE data. For both data sets we also produce and average 40 values of using hold-out cross-validation. In hold-out cross-validation we repeatedly split the data at random into training and target subsamples whose sizes are fixed to be the same as those arising from clustering in step 1 of the decay curve estimation. Then we fit an elastic net model on the training subsamples and predict the phenotypes in the target subsamples to estimates . Ideally, the decay curve should cross the area in which the points cluster.

We explore cross-population predictions using the HUMAN data and simulated phenotypes. Similarly to the above, we pick 5, 20, 100, 2000, 10000 and 50000 causal variants at random among those with minor allele frequency > 5% and we assign them normally-distributed effects such that h 2 ≈ 0.55. The same effect sizes are used for all populations. We then use individuals from Asia as the training population to estimate the decay curves. Those from other continents are the target populations for which we are assessing predictive accuracy, and we compute their and the corresponding predictive correlations . We use the points as terms of comparison to assess the quality of the curve, which should be close to them or at least cross the respective 95% confidence intervals.

We then repeat this simulation after adding the varieties available at the end of the second round of selection to the training population while considering the scenario with 200 and 1000 causal variants. The size of the training population is thus increased to 800 varieties, allowing us to explore the effects of a larger sample size and of considering new varieties from the breeding program to update the genomic prediction models when their predictive accuracy is no longer acceptable. In the following, we refer to this second population as the “augmented population” as opposed to the “original population” including only the 200 varieties described in steps 1 and 2 above.

The size of the training (n TR ) and target (n TA ) subsamples is determined by k-means. For the data used in this paper, k-means splits the training populations in two subsamples of comparable size; but we may require a smaller n TA ≪ n TR to estimate and the while at the same time a larger n TR is needed to fit the genomic prediction model. In that case, we increase n TR by moving individuals from the target subsample while keeping the between the two as large as possible. The impact on the estimated is likely to be small, because its precision depends more on the number of markers than on n TR and n TA [ 40 ]. The estimated and might be inflated because we are altering the subsets, even when does not change appreciably. Its variance, which can be approximated as in [ 49 ], decreases linearly in n TA except that can be compensated by generating more pairs of subsamples for each value of m.

As an alternative approach, we also consider estimating the decay rate of by linear regression of the against the ; we will denote the resulting predictive accuracy estimates with . For any set value of , we compare the at that with the corresponding value from the decay curve estimated by averaging all the for which . Assuming that the decay curve is in fact a straight line reduces the number of subsamples that we need to generate, enforces smoothness and makes it possible to compute for values of F ST larger than . On the other hand, the estimated will be increasingly unreliable as , because the regression line will provide negative instead of converging asymptotically to zero. We also regress the against the to investigate whether they have a stronger linear relationship than the with the , as suggested in [ 22 ] using simulated genotypes and phenotypes mimicking a dairy cattle population.

The pair of subsets produced by k-means corresponds to m = 0, hence the notation , and we increase m by steps of 2 to 20 until the between the subsamples is at most 0.005. We choose the stepping for each data set to be sufficiently small to cover the interval as uniformly as possible. The larger m is, the smaller we can expect to be. We repeat step 3(a) and 3(b) 40 times for each m to achieve the precision needed for an acceptably smooth curve.

We consider 376 wheat varieties from the TriticeaeGenome project, described in [ 4 ]. Varieties collected from those registered in France (210 varieties), Germany (90 varieties) and the UK (75 varieties) between 1946 and 2007 were genotyped using a combination of 2712 predominantly DArT markers. Several traits were recorded; in this paper we will focus on grain yield, height, flowering time, and grain protein content. Genotype-environment interactions were accounted for by an incomplete block design over trial fields in different countries, to prevent genomic prediction being biased by the country of registration of each variety. As in [ 4 ], we also group varieties in three groups based on their year of registration: pre-1990 (103 varieties), 1990 to 1999 (120 varieties), and post-1999 (153 varieties).

At the population level, the divergence between two populations due to drift, environmental adaptation, or artificial selection is commonly measured with F ST . Several estimators are available in the literature, and reviewed in [ 39 ]. In this paper we will adopt the estimator from [ 40 ], which is obtained by maximising the Beta-Binomial likelihood of the allele frequencies as a function of F ST . then describes how far the target population has diverged from the training population, which translates to “how far” a genomic prediction model will be required to predict. In terms of kinship, we know from the literature that the mean kinship coefficient between two individuals in different populations is inversely related to [ 41 ]: kinship can be interpreted as the probability that two alleles are identical by descent, which is inversely related to F ST which is a mean inbreeding coefficient. Intuitively, the fact that individuals in the two populations are closely related implies that the latter have not diverged much from the former: if is large, the marker profiles (and therefore the corresponding allele frequencies) will on average be similar. As a result, any clustering method that uses the Euclidean distance to partition a population into subsets will maximise their F ST by minimising . The simulations and data analyses below confirm experimentally that and are highly correlated, which makes them equivalent in building the decay curves; thus we will report results only for (see Section C in S1 Text ).

Another possible choice is the prediction error variance (PEV). It is commonly used in conjunction with GBLUP because, for that model, it can be estimated (for small samples) or approximated (for large samples) in closed form from Henderson’s mixed model equations [ 34 ]. In the general case no closed form estimate is available, but PEV can still be derived from Pearson’s correlation [ 35 ] for any kind of model as both carry the same information: (3) For consistency with our previous work [ 31 ] and with [ 4 ], whose results we partially replicate below, we will only consider predictive correlation in the following.

Predictive accuracy is often measured by the Pearson correlation ( ) between the predicted and observed phenotypes. When we use the fitted values from the training population as the predicted phenotypes, and assuming that the model is correctly specified, coincides with the proportion of genetic variance of the trait explained by the model and therefore , the heritability of the trait. (An incorrect model may lead to overfitting, and in that case .) When using cross-validation with random splits, and typically the difference will be noticeable ( ). However, may still overestimate the actual predictive accuracy in practical applications where target individuals for prediction are more different from the training population than the test samples generated using cross-validation [ 14 ]. This problem may be addressed by the use of alternative model validation schemes that mirror more closely the prediction task of interest; for instance, by simulating progeny of the training population to assess predictive accuracy for a genomic selection program. This approach is known as forward prediction and is common in animal breeding [ 19 , 33 ].

A baseline model for genomic prediction of quantitative traits is the genomic BLUP (GBLUP; [ 24 , 25 ]), which is usually written as (1) where g is a vector of genetic random effects, Z is a design matrix that can be used to indicate the same genotype exposed to different environments, K is a kinship matrix and ε is the error term. Many of its properties are available in closed form thanks to its simple definition and normality assumptions, including closed form expressions of and upper bounds on predictive accuracy that take into account possible model misspecification [ 15 ]. Other common choices are additive linear regression models of the form (2) where y is the trait of interest; X are the markers (such as SNP allele counts coded as 0, 1 and 2 with 1 the heterozygote); β are the marker effects; and ε are independent, normally-distributed errors with variance . Depending on the choice of the prior distribution for β , we can obtain different models from the literature such as BayesA and BayesB [ 25 ], ridge regression [ 26 ], the LASSO [ 27 ] or the elastic net [ 28 ]. The model in Eq (1) is equivalent to that in Eq (2) if the kinship matrix K is computed from the markers X and has the form X X T and β ∼ N(0, VAR( β )) [ 29 , 30 ]. In the remainder of the paper we will focus on the elastic net, which we have found to outperform other predictive models on real-world data [ 31 ]. This has been recently confirmed in [ 32 ].

The cross-population prediction simulation based on the HUMAN data ( Fig 3 ) generated results consistent with those above. As before, the number of causal variants appears to influence the behaviour of the decay curve: while the decrease linearly for 20, 100 and 2000 casual variants, they converge to 0.65 for 5 causal variants. However, unlike in the genomic selection simulation, the quality of the estimated decay curve does not appear to degrade as the number of causal variants decreases. This difference may depend on the lack of a systematic selection pressure in the current simulation, which made the decay curve overestimate predictive correlation when considering 10 variants in the previous simulation. Finally, as in the analysis of the MICE data, the linear approximation to the decay curve provides a way to extend the reach of the decay curve to estimate predictive correlations for distantly related populations (AMERICA, AFRICA, OCEANIA). Again we observe some loss in precision (see Table 2 ), but the extension still crosses the 95% confidence intervals of those 14 times out of 18 (78%).

The decay curves fitted on the augmented training populations (800 varieties, now including those available at the end of the second round of selection, Fig 2 ) fit the first four generations well ( for the first two, for the third and the fourth). As before, the only exception is the first generation in the simulation with 1000 variants, with an absolute difference of 0.09. However, the decay curves are also able to capture the long-range decay rates through their linear approximations. When considering 200 causal variants, for generations 5 to 7 and ≈0.10 for generations 8 and 9; and for generations 4 to 9 when considering 1000 causal variants. This can be attributed to the increased sample size of the training population, which both improves the goodness of fit of the estimated decay curve; and makes the decay rate of the closer to linear, thus making it possible for the to approximate it well over a large range of F ST values. To investigate this phenomenon, we gradually increased the initial training population to 4000 varieties through random mating and we observed that for such a large sample size indeed decreases linearly as a function of F ST . We conjecture that this is due to a combination of the higher values observed for and their slower rate of decay, which prevents the latter from gradually decreasing as is still far from zero after 10 generations. In addition, we note that increasing the number of causal variants has a similar effect; with 200 and 1000 causal variants indeed decreases with an approximately linear trend, which is not the case with 10 and 50 causal variants.

The decay curves from the genomic selection simulation on the original training population (200 varieties), shown in blue in Fig 1 , span two rounds of selection and three generations. When considering 200 or 1000 causal variants, the curve overlaps the mean behaviour of the simulated data points (shown in green) almost perfectly: the difference between the generation means and the decay curve is ⩽ 0.06 for the first three generations, with the exception of the first generation in the simulation with 1000 variants ( ). As the number of causal variants decreases (50, 10), the decay curve increasingly overestimates , although the difference remains ⩽0.10 for the first two generations; and both show a slower decay than the . This appears to be due to a few alleles of large effect becoming fixed by the selection, leading to a rapid decrease of without a corresponding rapid increase in .

Furthermore, the decay curves for the phenotypes in the WHEAT data confirm two additional considerations originally made in [ 4 ]. Firstly, [ 4 ] noted that the distribution of the Ppd-D1a gene, which is a major driver of this flowering time, varies substantially with the country of registration and thus cross-country predictions are not reliable. Fig B.1 in S1 Text shows that the decay curve vastly overestimates the predictive correlation for both Germany and the UK. Splitting the WHEAT data in two halves that contain equal proportions of both alleles of Ppd-D1a and that are genetically closer overall ( ), we obtain a decay curve that fits the predictive correlations reported in the original paper ( , ). Secondly, we also split the data according to their year of registration and use the oldest varieties (pre-1990) as a training sample for predicting yield. Again the decay curve crosses the 95% confidence intervals for the predictive correlations reported in [ 4 ] and the correlations themselves are within 0.05 of the average from the decay curve both for 1990-1999 ( , , ) and post-2000 ( , , ) varieties.

Several interesting points arise from the analysis of the real phenotypes in the WHEAT and MICE data, shown in Table 2 and in Figs B.1, B.2 and B.3 in S1 Text . Firstly, cross-validation always produces pairs of subsamples with and high that are located at the left end of the decay curve. The average is 0.006 for the WHEAT data and 0.001 for the MICE data, and the difference between the average and the corresponding is ≪ 0.02 10 times out of 12 (83%, see Table B.4 in S1 Text ). The spread of the is also similar to that of the . Secondly, we note that in the WHEAT data all decay curves but that for flowering time cross the 95% confidence intervals for the cross-country predictive correlations for Germany and UK reported in [ 4 ]. Even in the MICE data, in which all families are near the end or beyond the reach of the decay curves, the latter (or their linear approximations) cross the 95% confidence intervals for the 18 times out of 24 (75%). However, we also note that those intervals are wide due to the limited sizes of those populations.

The range of the predictive correlations around the decay curves varies between 0.05 and 0.10, and it is constant over the range of observed for each curve. It does not appear to be related to either the size of the training subsample or the number of causal variants. This is apparent in particular from the genomic selection simulation, in which both are jointly set to different combinations of values. Similarly, there seems to be no relationship between the spread and the magnitude of the predictive correlations ( ). This amount of variability is comparable to that of other studies (e.g., the range of the is smaller than that in the cross-validated correlations in [ 32 ]) once we take into account that the are individual predictions and are not averaged over multiple repetitions. Furthermore, subsampling further reduces the size of the training subpopulations; and fitting the elastic net requires a search over a grid of values for its two tuning parameters, which may get stuck in local optima.

In all the simulations and the real-world data analyses the from the decay curve is close to the linear interpolation ; considering all the reference populations in Table 2 and the generation means in Tables A.1 and A.2 in S1 Text , 41 times out of 47 (87%). Both estimates of predictive correlation are close to the respective reference values and ; the difference (in absolute value) is ≪ 0.05 39 times (41%) and ≪ 0.10 69 times (73%) out of 94. The proportion of small differences increases when considering only target populations that fall within the span of the decay curve: 23 out of 44 (52%) are ≪ 0.05 and 38 are ≪ 0.10 (84%). This is expected because the decay curve is already an extrapolation from the training population, so extending it further with the linear interpolation reduces its precision. Regressing against the does not produce a stronger linear relationship than that represented by (p = 0.784, see Section D in S1 Text ).

Simulation of quantitative traits with 5 (top left), 20 (top right), 100 (middle left), 2000 (middle right), 10000 (bottom left) and 50000 (bottom right) causal variants from the Asian individuals in the HUMAN data. The blue circles are the used to build the curve, and the red point is . The blue line is the mean decay trend, with a shaded 95% confidence interval, and the dashed blue line is the linear interpolation provided by the . The red squares labelled EUROPE, MIDDLE EAST, AMERICA, AFRICA and OCEANIA correspond to the for the individuals from those continents, and the red brackets are the respective 95% confidence intervals.

The decay curves from the simulations are shown in Figs 1 , 2 and 3 , and the corresponding predictive correlations are reported in Tables 1 and 2 and S1 Text . The predictive correlations for the WHEAT and MICE data sets are reported in Table 2 , and the decay curves are shown in Figs 1 , 2 and 3 and S1 Text . A summary of the different predictive correlations defined in the Methods and discussed here is provided in Table 1 .

Discussion

Being able to assess the predictive accuracy is important in many applications, and will assist in the development of new models and in the choice of training populations. A number of papers have discussed various aspects of the relationship between training and target populations in genomic prediction, and of characterising predictive accuracy given some combination of genotypes and pedigree information. For instance, [51] discusses how to choose which individuals to include in the training population to maximise prediction accuracy for a given target population using the coefficient of determination. [52] separates the contributions of linkage disequilibrium, co-segregation and additive genetic relationships to predictive accuracy, which can help in setting expectations about the possible performance of prediction. [53] and [22] link predictive accuracy to kinship in a simulation study of dairy cattle breeding; and [54] investigates the impact of population size, population structure and replication in a simulated biparental maize populations. The approach we take in this paper is different in a few, important ways. Firstly, we choose to avoid the parametric assumptions underlying GBLUP and the corresponding approximations based on Henderson’s equations that provide closed-form results on predictive accuracy in the literature. It has been noted in our previous work [31] and in the literature (e.g. [32]) that in some settings GBLUP may not be competitive for genomic prediction; hence we prefer to use models with better predictive accuracy such as the elastic net for which the parametric assumptions do not hold. Our model-agnostic approach is beneficial also because decay curves can then be constructed for current and future competitive models, since the only requirement of our approach is that they must be able to produce an estimate of predictive correlation. Secondly, we demonstrate that the decay curves estimated with the proposed approach are accurate in different settings and on human, plant and animal real-world data sets. This complements previous work that often used synthetic genotypes and analysed predictive accuracy in a single domain, such as forward simulation studies on dairy cattle data. Finally, we recognise that the target population whose phenotypes we would like to predict may not be available or even known when training the model. In plant and animal selection programs, one or more future rounds of crossings may not yet have been performed; in human genetics, prediction may be required into different demographic groups for which no training data are available. Therefore, we are often limited to extrapolating a to estimate the we would observe if the target population were available. Prior information on values is available for many species such as humans [39, 43]; and can be used to extract the corresponding from a decay curve.

We observe that the decay rate of is approximately linear in for most of the curves, suggesting that regressing the against the is a viable estimation approach. This has the advantage of being computationally cheaper than producing a smooth curve with LOESS since it requires fewer points and thus fewer genomic prediction models to be fitted. In fact, if we assume that the decay rate is linear we could also estimate it as the slope of the line passing through and for a single, small value of m. It should be noted, however, that several factors can cause departures from linearity, including the number of causal variants underlying the trait, the use of small training populations and the confounding effect of exogenous factors. In the case of the MICE data, for instance, predictions may be influenced by cage effects; in the case of the WHEAT data, environmental and seasonal effects might not be perfectly captured and removed by the trials’ experimental design. We also note that the decay curves for traits with small heritabilities will almost never be linear, because converges asymptotically to zero. Unlike the results reported in [22], we do not find a statistically significant difference between the strength of the linear relationship between and and that between the respective squares. There may be several reasons for this discrepancy; the simulation study in [22] was markedly different from the analyses presented in this paper, since it used simulated genotypes to generate the population structure typical of dairy cattle and since it used GBLUP as a genomic prediction model.

We also observe that when , both and are, as expected, similar to the obtained by applying cross-validation to the training populations selected from the WHEAT and MICE data. This suggests that indeed is an accurate measure of predictive accuracy only when the target individuals for prediction are drawn from the same population as the training sample, as previously argued by [14] and [19], among others.

Some limitations of the proposed approach are also apparent from the results presented in the previous section. The most important of these limitations appears to be that in the context of a breeding program the performance of the decay curve depends on the polygenic nature of the trait being predicted, as we can see by comparing the panels in Fig 1. This can be explained by the fact that causal variants underlying less polygenic, highly and moderately heritable traits will necessarily have some individually large effects. As each of those variants approaches fixation due to selection pressure, allele frequencies in key areas of the genome will depart from those in the training population and the accuracy of any genomic prediction model will rapidly decrease [21]. However, these selection effects are genomically local and so have little impact on . A similar effect has been observed for flowering time in the WHEAT data. [4] notes that the Ppd-D1a gene is a major driver of early flowering, but it is nearly monomorphic in one allele in French wheat varieties and nearly monomorphic in the other allele in Germany and the UK. As a result, even though the for those countries are as small as 0.031 and 0.042, widely overestimates in both cases. A possible solution would be to compute only on the relevant regions of the genome or, if their precise location is unknown, on the relevant chromosomes; or to weight to promote genomic regions of interest.

On the other hand, in the case of more polygenic traits a larger portion of the genome will be in linkage disequilibrium with at least one causal variant, and their effects will be individually small. Therefore, will increase more quickly in response to selection pressure and changes in predictive accuracy will be smoother, thus allowing to track them more easily. Indeed, in the WHEAT data the genomic prediction model for flowering time has a much smaller number of non-zero coefficients (28) compared to yield (91), height (286) and grain protein content (121). Similarly, in the MICE data the model fitted on F010 to predict weight has only 168 non-zero coefficients while others range from 212 to 1169 non-zero coefficients. By contrast, all models fitted for predicting weight, which correspond to curves that well approximate other families’ , have between 1128 and 2288 non-zero coefficients.

The simulation on the HUMAN data suggests different considerations apply to outbred species. Having some large-effect causal variants does not necessarily result in low quality decay curves; on the contrary, if we assume that the trait is controlled by the same causal variants in the training and target populations it is possible to have a good level of agreement between the and the . Intuitively, we expect strong effects to carry well across populations and thus does not decrease beyond a certain F ST . However, this will mean that the curves will not be linear and will underestimate (see Fig 3, top left panel). We also note that effect sizes are the same in all the populations, which may make our estimates of predictive accuracy optimistic.

Another important consideration is that since the decay curve is extrapolated from the training population, its precision decreases as F ST increases, as can be seen from both simulations and by comparing the WHEAT and MICE data. Predictions will be poor in practice if the target and the training populations are too genetically distinct; an example are rice subspecies [17], which have been subject to intensive inbreeding. The trait to be predicted must have a common genetic basis across training and target populations. However, the availability of denser genomic data and of larger samples may improve both predictive accuracy and the precision of the decay curve for large F ST . Furthermore, the range of the decay curve in terms of F ST depends on the amount of genetic variability present in the training population; the more homogeneous it is, the more unlikely that k-means clustering will be able to split it in two subsets with high . One solution is to assume the decay is linear and use instead of to estimate ; but as we noted above this is only possible if . If , the decay curve estimated with LOESS from can converge asymptotically to zero as increases; but the linear regression used to estimate will continue to decrease until . Another possible solution is to try to increase by moving observations between the two subsets, but improvements are marginal at best and there is a risk of inflating .

Even with such limitations, estimating a decay curve for predictive correlation has many possible uses. In the context of plant and animal breeding, it is a useful tool to answer many key questions in planning genomic selection programs. Firstly, different training populations (in terms of allele frequencies, sample size, presence of different families, etc.) can be compared to choose that which results in the slowest decay rate. Secondly, the decay curve can be used to decide when genomic prediction can no longer be assumed to be accurate enough for selection purposes, and thus how often the model should be re-trained on a new set of phenotypes. Unlike genotyping costs, phenotyping costs for productivity traits have not decreased over the years. Furthermore, the rate of phenotypic improvements (i.e. selection cycle time) can be severely reduced by the need of performing progeny tests. Therefore, limiting phenotyping to once every few generations can reduce the cost and effort of running a breeding program. The presence of close ancestors in the training population suggests that decay curves are most likely reliable for this purpose, as we have shown both in the simulations and in predicting newer wheat varieties from older ones in the WHEAT data.

The other major application of decay curves is estimating the predictive accuracy of a model for target populations that, while not direct descendants of the training population, are assumed not to have strongly diverged and thus to have comparable genetic architectures. Some examples of such settings are the cross-country predictions for the WHEAT data, the cross-family predictions for the MICE data and across human populations. In human genetics, decay curves could be used to study the accuracy of predictions and help predict the success of interventions of poorly-studied populations. In plant and animal breeding, on the other hand, it is common to incorporate distantly related samples in selection programs to maintain a sufficient level of genetic variability. Decay curves can provide an indication of how accurately the phenotypes for such samples are estimated, since the model has not been trained to predict them well and they are not as closely related as the individuals in the program.