Improving the efficiency of multi-location field trials with complete and incomplete relationship information

The increasingly cost-efficient availability of ‘omics’ data has led to the development of a rich framework for predicting the performance of non-phenotyped selection candidates in recent years. The improvement of phenotypic analyses by using pedigree and/or genomic relationship data has however received much less attention, albeit it has shown large potential for increasing the efficiency of early generation yield trials in some breeding programs. The aim of this study was accordingly to assess the possibility to enhance phenotypic analyses of multi-location field trials with complete relationship information as well as when merely incomplete pedigree and/or genomic relationship information is available for a set of selection candidates. For his purpose, four winter bread wheat trial series conducted in Eastern and Western Europe were used to determine the experimental efficiency and accuracy of different resource allocations with a varying degree of relationship information. The results showed that modelling relationship between the selection candidates in the analyses of multi-location trial series was up to 20% more efficient than employing routine analyses, where genotypes are assumed to be unrelated. The observed decrease in efficiency and accuracy when reducing the testing capacities was furthermore less pronounced when modelling relationship information, even in cases when merely partial pedigree and/or genomic information was available for the phenotypic analyses. Exploiting complete and incomplete relationship information in both preliminary yield trials and multi-location trial series has thus large potential to optimize resource allocations and increase the selection gain in programs that make use of various predictive breeding methods.


Introduction
The increasing availability of cost-efficient 'omics' data has led to the development of a rich framework for predicting genotype performance in recent years (Robertsen et al. 2019;Montesinos-López et al. 2021;Sneller et al. 2021;Bayer et al. 2021). The usage of genome-wide distributed markers for a socalled genomic selection has thereby gained an especially large popularity in many breeding programs 1 3 10 Page 2 of 11 Vol:. (1234567890) (Belamkar et al. 2018;Juliana et al. 2019;Haikka et al. 2020;Raffo et al. 2022). A major goal in the genomic selection framework is given by obtaining as accurate as possible predictions of non-phenotyped selection candidates in early generations for traits like grain yield (Tsai et al. 2020;Borrenpohl et al. 2020), disease resistance (Beukert et al. 2020;Moreno-Amores et al. 2020) as well as costly and laborious to phenotype quality traits (Schmidt et al. 2016;Lado et al. 2018). The improvement of phenotypic analyses by using pedigree and/or genomic relationship data has on the other hand received much less attention (Endelman et al. 2014;Terraillon et al. 2022), albeit it has shown large potential for increasing the efficiency of early generation observation and preliminary yield trials in some breeding programs (Michel et al. 2019;Tsai et al. 2020;Borrenpohl et al. 2020). These studies generally assumed that a particular source of relationship information is fully covering the entire set of selection candidates, which is however not always the case in practice, for example when some of the tested lines were developed by another breeding program. The aim of this study was accordingly to take a step towards generalizing these previous results obtained for preliminary yield trials, and assess the possibility of enhancing phenotypic analyses of multi-location field trials when merely incomplete pedigree and/or genomic relationship information is available for a set of selection candidates.

Plant material and genotypic data
Four panels of 147-177 recombinant inbred and double haploid breeding lines developed in the winter bread wheat breeding program of Saatzucht Donau GesmbH & CoKG in Austria were analysed in this study. Each panel was phenotyped for grain yield in a different trial series each with four locations in Western Europe in 2015 (151 lines) and 2016 (150 lines) as well as in Eastern Europe in 2015 (177 lines) and 2016 (147 lines). Environmental means were available for each of the locations within the respective trial series (year-by-region combinations), so each line occurred with four replicates (total number of observations) within each of the trial series. The breeding lines were part of 277 different families with a size of 1-13 lines per family and a genealogy of 338 ancestors tracing back up to 8 generations. All lines were genotyped with the DArT genotyping-by-sequencing approach (Diversity Arrays Technology Pty Ltd 2020), and markers with more than 10% missing data and a minor allele frequency smaller than 5% were filtered out. Only one marker of identical marker pairs was furthermore retained for all subsequent analyses. A chromosome-wise imputation of missing data points with the missForest algorithm (Stekhoven and Bühlmann 2012) and after quality filtering the final marker dataset contained 1908 markers, which were used for investigating the population structure (Suppl. Fig. S1).

Phenotypic analysis
An across-trial analysis was conducted for each of the four trial series individually by using a linear mixed model of the form: where y jk are the observations for grain yield, μ is the grand mean, and l k is the effect of the kth location that was modelled as random. The effect of the jth line g j was modeled as random with ∼ N , σ 2 g to obtain an estimate of the genetic variance and best linear unbiased predictions (BLUP) of the lines' performances. The effect e jk that incorporated both the line-by-trial interaction variance and the residual effect was assumed to be random following a normal distribution with ∼ N , σ 2 e . The entry-mean heritability was subsequently estimated following the suggestion by Cullis et al. (2006): where σ 2 g is the genetic variance and VD the mean variance of a difference of two genotypic BLUPs. All phenotypic analyses were conducted with the package sommer (Covarrubias-Pazaran 2016) for the R statistical environment (R Core Team 2022).
Empirical assessment of the experimental efficiency Sets of 70 lines, each coming from a different family, were 50 times randomly sampled from each of (1) y jk = μ + g j + l k + e jk (2) h 2 = 1 − VD 2σ 2 g Page 3 of 11 10 Vol.: (0123456789) the four investigated trial series individually. This resulted in 50 unique and different sets per trial series (year-by-region combination) and a total of 200 sets across all trial series. These sets were subsequently analysed separately for assessing the efficiency of several experimental layouts without including relationship information as well as with including complete and incomplete relationship information into the analyses of the phenotypic data. The experimental designs comprised a fully orthogonal testing of the lines across all four locations of a given trials series as well as a reduction of the testing capacities to three or two locations. Furthermore, the merit of allocating lines to the locations according to an incomplete block design was tested by reducing the number of total observations per line from four to two. The percentage of genotyped lines within these sets was subsequently varied between 0 and 100%, whereas the number of lines with pedigree information was varied between 0, 30, 50, 80, and 100%. The corresponding proportions of lines were thereby sampled randomly, while the sampling of pedigreed and genotyped lines was additionally independent from each other. The lines were in this way allocated to groups for which both genomic and pedigree relationship information was available, only genomic or pedigree relationship information was available, and one group without relationship information (Fig. 1). Model [1] was again used to obtain the mean variance of a difference VD of two genotypic BLUPs for each of these scenarios, where the effect of the jth line g j followed in this case ∼ N , σ 2 g , where was computed as: with the genomic relationship matrix adj and the pedigree relationship matrix . The matrix 11 contained the pedigree relationship between non-genotyped lines, 22 the pedigree relationship between genotyped lines, while 12 and 21 modelled the pedigree relationship between genotyped and nongenotyped lines. All lines were assumed to have undergone seven cycles of selfing for the computation of , while non-pedigreed lines i.e., lines that were part of above-mentioned group without pedigree information were likewise included into the pedigree relationship matrix but possessed a covariance of zero with all other lines in . The genomic relationship matrix was computed following Endelman and Jannink (2012): where is a centered marker matrix of the j lines with W jl = Z jm + 1 − 2p m and m being the allele frequency at the mth marker locus. The genomic relationship matrix was moreover adjusted by solving: and setting adj = a + b as suggested by Christensen et al. (2012) before computing , in order to account for the impact of genetic trends across multiple generations and the reduction in genetic variance from the base population to the population of genotyped lines. It should be noticed that in the case all lines possess genotypic information reduces to , while in the case no genotypic information is available reduces to , and given no relationship is available reduces to the identity matrix .
The efficiency modelling complete and incomplete relationship information in the evaluated experimental designs for each sampled set of lines in a given trial series was determined by using the mean variance of a difference analogous to Piepho et al. (2006): where VD REF is the mean variance of a difference (squared standard error of a difference) of all pairwise comparisons among the genotypic BLUPs obtained from the analysis of a set of lines that was completely orthogonal tested in all four locations in a given trial series without including any relationship information i.e., ∼ N , σ 2 g . This reference value was compared with VD HBLUP , which is the mean variance of a difference of all pairwise comparisons among the genotypic BLUPs obtained from the multi-location analysis of the different experimental designs that

VD HBLUP
Euphytica (2023) 219:10 10 Page 4 of 11 Vol:. (1234567890) Fig. 1 The percentage of lines for which no relationship information was available (1st row), only pedigree (2nd row) or genomic relationship information (3rd row) was available, and the percentage of lines for which both pedigree and genomic relationship information was available (4th row) with different proportions of pedigreed and genotyped lines in the randomly sampled subsets Page 5 of 11 10 Vol.: (0123456789) included complete or incomplete pedigree and genomic relationship information as described above.
Simulation layout for assessing the prediction accuracy Analogous to the empirical study different sets of 70 lines, each coming from a different family, were 50 times randomly sampled each of the four investigated trial series individually in order to assess the prediction accuracy i.e., the correlation of the true genotypic value with the predicted genotypic value in a simulation study. True genotypic values of each line were derived by randomly sampling N QTL = 150 marker loci as causal variants of a quantitative inherited trait, for which a vector of effects was randomly sampled from a normal distribution with ∼ N(0, 1): where is the vector of effects of the causal loci, is the marker matrix of the investigated set of lines, and is the vector of their true genotypic values. The observed genotypic values were accordingly computed as: where the vector of error effects was randomly sampled from a normal distribution with zero mean and a variance equal to where σ 2 TGV is the variance of the true genotypic values, and h 2 an aspired repeatability of h 2 = 0.10 , h 2 = 0.30 and h 2 = 0.50 for a given location. Like in the empirical study, 2-4 locations with completely orthogonal testing or an allocation according to an incomplete block design were simulated and used assess the prediction accuracy. The prediction accuracy was in this case measured as r , , where is the predicted line performance when analysing the data with linear mixed models following Eq. (1) including the relationship matrix described in Eq. (3). The same proportions of complete and incomplete relationship information described above for the empirical study were tested in the simulations. The pedigree relationship matrix was like in the empirical study based on the original pedigree records, while the genomic relationship matrix was constructed with random samples of N SNP = 1500 markers that served as linked loci to the causal variants.
The genomic relationship matrix was computed with sommer (Covarrubias-Pazaran 2016), the pedigree relationship matrix was obtain with the package pedigreeTools (Vazquez et al. 2018), the combined relationship matrix was derived with the package AGHmatrix (Amadeu et al. 2016), and the incomplete bock designs were randomized with the package crossdes (Sailer 2022) for the R statistical environment (R Core Team 2022). All models for assessing the experimental efficiency and prediction accuracy were fitted with the R package sommer (Covarrubias-Pazaran 2016). An example dataset and accompanied R Code are available as supplemental material to illustrate the utilized models.

Results
The different trial series conducted in Eastern and Western Europe in 2015-2016 showed a substantial genotype-by-environment interaction exemplified by an average correlation of r = 0.18-0.27 between the series-specific locations. Nevertheless, a broad genetic variation was observed in each trial series for which the estimated entry-mean heritability varied between h 2 = 0.45 and h 2 = 0.58, which suggested that the dataset at hand was suitable for investigating the experimental efficiency with complete and incomplete relationship information as well as varying resource allocations (Table 1).
The empirical assessment of the experimental efficiency with different sets of randomly sampled lines revealed that a reduction in the number of test locations resulted on average in a 10-20% loss in efficiency in comparison to a completely orthogonal testing in four locations (Fig. 2). This reduction in efficiency was however much less pronounced if pedigree and/or genomic relationship was integrated into the phenotypic analysis of the data. Modelling the relationship between the lines by = , i.e., ∼ N , σ 2 g and ∼ N 0, σ 2 g , was in fact up to 20% more efficient than modelling lines independent in a routine analysis with , i.e., 10 Page 6 of 11 Vol:. (1234567890) ∼ N , σ 2 g , in the investigated resource allocations. A similar observation was made for cases in which merely incomplete genomic and/or pedigree relationship was available, where an increase in the availability of relationship information steadily increased the experimental efficiency. Interestingly, decreasing the testing capacities by up to one quarter appeared to be feasible in all investigated trial series without losing any efficiency in the case that most of the lines possessed genotypic data (Suppl. Figs. S2-S5).
The simulations showed furthermore that modelling relationship can result in a higher prediction accuracy (Fig. 3), which was in this study defined as the correlation between the true genotypic value with the predicted genotypic value obtained in the phenotypic analyses of the data. The results of the simulations generally followed the same pattern that has been observed for the experimental efficiency in the empirical investigations, with some advantage of genomic over pedigree relationship information. A gradual increase in prediction accuracy was accordingly observed when modelling pedigree and/or genomic relationship between the lines in comparison to a baseline analysis that assumed independence between the lines. The advantage of modelling relationship information in comparison to this baseline model diminished however with an increase in the repeatability of the individual trials at each location, as did the superiority of genomic over pedigree relationship information (Suppl. Figs. S6-S7).

Discussion
Modelling genetic relationship as an additional source of information in phenotypic analyses has shown promising results for determining the performance of genotypes both in simulation (Bauer et al. 2006;Möhring et al. 2014;Selle et al. 2019;Terraillon et al. 2022) and empirical studies (Moreau et al. 1999;Oakey et al. 2007;Endelman et al. 2014). Although the estimated experimental efficiency and accuracy was highest with complete genomic relationship information in the study at hand, a marked advantage was likewise observed when utilizing pedigree records for modelling relationship. The latter resulted in lower accuracies in comparison to a genomic prediction of non-phenotyped individuals (Auinger et al. 2016;Cericola et al. 2017) as the Mendelian sampling term i.e., segregation within families cannot be addressed in such a case, but pedigree best linear unbiased predictions can readily distinguish between family members if phenotypic observations are already available for them .
This issue renders pedigree records a relatively cost-efficient alternative to genomic data, depending on the strategy of a program that employs predictive breeding methods for an array of various target traits. Nevertheless, this assumes generally an ideal case where a particular source of relationship information is fully available for the entire set of selection candidates. However, pedigree and genomic data might not be available for every breeding line, especially when it was developed by another breeding program and is tested together with 'in-house' developed material in the framework of the breeders' exemption and Table 1 Mean, range, variance components and heritability of grain yield (dt ha −1 ) as well as the number of lines tested in each trial series of four locations in Eastern and Western Europe in 2015-2016 Genetic variance ( σ 2 g ), location variance ( σ 2 l ), residual variance ( σ 2 e ), and entry-mean heritability (h 2 ) † Environmental means were available for each line and the four locations within each trial series, i.e., each line occurred with a total number of four observations within each of the trial series bilateral germplasm exchange. Hence, some lines have to be assumed independent in the phenotypic analysis since no relationship information is available for them, even though it might be desirable to integrate such information into the analysis. Employing the single-step framework developed in animal breeding for combining different relationship matrices into a common matrix (Legarra et al. 2009;Christensen and Lund 2010) enabled to commonly rank lines with and without relationship information as well as Fig. 2 Average efficiency in the empirical study across all four investigated trial series, expressed relatively to a completely orthogonal testing in all four locations in these trial series without including any relationship information. The investigated experimental designs included a fully orthogonal testing of the lines across two to four locations as well as allo-cating lines to the locations according to an incomplete block design by reducing the number of total observations per line in a given trial series from four to two. The percentage of genotyped lines was additionally varied between 0 and 100%, while the number of lines with pedigree information was varied between 0, 30, 50, 80, and 100% 10 Page 8 of 11 Vol:. (1234567890) simultaneously improving the ranking between lines with different types i.e., pedigree and/or genomic relationship information in study at hand. This step towards generalization of modelling relationships in phenotypic analyses led to a considerable increase both in the experimental efficiency and accuracy. The simulations suggested the largest benefit for traits with a low to medium heritability like grain yield Fig. 3 Average prediction accuracy in the simulation study with a repeatability of h 2 = 0.10 at each trial location. The investigated experimental designs included a fully orthogonal testing of the lines across two to four locations as well as allocating lines to the locations according to an incomplete block design by reducing the number of total observations per line in a given trial series from four to two. The percentage of geno-typed lines was additionally varied between 0 and 100%, while the number of lines with pedigree information was varied between 0, 30, 50, 80, and 100%. The horizontal black lines correspond to the prediction accuracy of a completely orthogonal testing in all four locations of these trial series without including any relationship information Page 9 of 11 10 Vol.: (0123456789) or protein yield in winter bread wheat, whereas this advantage appeared to be rather marginal for traits with a high heritability like plant height or flowering date. This observation was in line with previous reports, which stated that the relative advantage of modelling relationship information in the phenotypic analysis depends on the heritability as well as the testing intensity of field trials (Bauer et al. 2006;Endelman et al. 2014;Terraillon et al. 2022).
The results of the empirical and simulation study showed moreover that a reduction in testing capacities appears to be feasible in multi-location trials when the pedigree and/or genomic data are integrated into the phenotypic analysis. Although resource allocations that follow an incomplete block (Montesinos-Lopez et al. 2022) or augmented design (Lell et al. 2021) have shown promise in the genomic prediction framework, the logistically easiest option for reducing the testing capacities is given by reducing the number of test locations. The feasibility of this option might however be strongly dependant on the target population of environments of a breeding program. Reducing the number of test locations might for example be readily feasible for regional breeding program with a strong focus on local adaptation, whereas it is probably less suitable for a breeding program that targets many different environments. Resource allocations with sparse testing strategies might thus be more suitable for the latter (Jarquin et al. 2020;Atanda et al. 2022), while including relationship information both into the randomization of specific trial designs as well as the subsequent analyses has shown large merit to increase the experimental efficiency in general (Cullis et al. 2020).

Conclusions
The usage of complete or incomplete relationship information has the potential to render phenotypic analyses more efficient. Such an application has been primarily suggested for preliminary yield trials (Endelman et al. 2014) and applied in practical breeding programs (Michel et al. 2019;Tsai et al. 2020;Borrenpohl et al. 2020), but can also have some merit for multi-location trials in advanced generations. A strong increase of efficiency and accuracy was accordingly observed within the four investigated trials series. The usage of relationship information might furthermore benefit the analysis across multiple trial series especially if their connectivity by commonly tested genotypes is low. Hence, exploiting complete and incomplete relationship information in both preliminary yield trials and multi-location trial series has a large potential to further optimize resource allocations and increase the yearly selection gain in programs that make use of various predictive breeding methods.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.