# ABSTRACT:

Genome-wide selection (GWS) is based on a large number of markers widely distributed throughout the genome. Genome-wide selection provides for the estimation of the effect of each molecular marker on the phenotype, thereby allowing for the capture of all genes affecting the quantitative traits of interest. The main statistical tools applied to GWS are based on random regression or dimensionality reduction methods. In this study a new non-parametric method, called Delta-p was proposed, which was then compared to the Genomic Best Linear Unbiased Predictor (G-BLUP) method. Furthermore, a new selection index combining the genetic values obtained by the G-BLUP and Delta-p, named Delta-p/G-BLUP methods, was proposed. The efficiency of the proposed methods was evaluated through both simulation and real studies. The simulated data consisted of eight scenarios comprising a combination of two levels of heritability, two genetic architectures and two dominance status (absence and complete dominance). Each scenario was simulated ten times. All methods were applied to a real dataset of Asian rice (Oryza sativa) aiming to increase the efficiency of a current breeding program. The methods were compared as regards accuracy of prediction (simulation data) or predictive ability (real dataset), bias and recovery of the true genomic heritability. The results indicated that the proposed Delta-p/G-BLUP index outperformed the other methods in both prediction accuracy and predictive ability.

Keywords:
genomic prediction; selection index; genetic gain; asian rice

# Introduction

The evolution of sequencing and genotyping techniques has prompted a great advance in molecular genetics becoming a possibility using DNA information directly for the selection of genetically superior individuals. This framework is called genome-wide selection (GWS), which enables direct early marker-based selection that generally increases the genetic gain per unit of time (Meuwissen et al., 2001Meuwissen, T.H.E.; Hayes, B.J.; Goddard, M.E. 2001. Prediction of total genetic value using genome-wide dense marker maps. Genetics 157: 1819-1829.). Genome-wide selection is based on the estimation of molecular marker effects on the phenotypes, thus providing accurate additive and dominant genomic value predictions for individuals subject to selection, even for those without observed phenotypic values.

The Genomic Best Linear Unbiased Predictor (G-BLUP) (Goddard, 2009Goddard, M.E. 2009. Genomic selection: prediction of accuracy and maximization of long term response. Genetics 136: 345-357.; Habier et al., 2007Habier, D.; Fernando, R.L.; Dekkers, J.C.M. 2007. The impact of genetic relationship on genome-assisted breeding values. Genetics 117: 2389-2397.; Van Raden, 2008Van Raden, P.M. 2008. Efficient methods to compute genomic predictions. Journal of Dairy Science 91: 4414-4423.) has been extensively applied to GWS and recommended for additive and dominant genomic value predictions (Azevedo et al., 2015Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105.; De los Campos et al., 2013De los Campos, G.; Hickey, J.M.; Pong-Wong, R.; Daetwyler, H.D.; Callus, M.P.L. 2013. Whole genome regression and prediction methods applied to plant and animal breeding. Genetics 193: 327-345.; Gianola, 2013Gianola, D. 2013. Priors in whole-genome regression: the Bayesian alphabet returns. Genetics 194: 573-596.; Gianola et al., 2009Gianola, D.; De los Campos, G.; Hill, W.G.; Manfredi, E.; Fernando, R. 2009. Additive genetic variability and the Bayesian alphabet. Genetics 183: 347-363.). However, non-parametric and purely conceptual methods can also be used for this purpose. The Delta-p method was theoretically proposed by Resende (2015)Resende, M.D.V. 2015. Quantitative and Population Genetics = Genética Quantitativa e de Populações. Editora Suprema, Visconde do Rio Branco, MG, Brazil (in Portuguese). for application to genomic selection, but has not been implemented and evaluated so far. The method is based on the genetic distance between two subpopulations, using the concepts of changes in allele frequency due to selection and the genetic gain theory. It is also possible to combine the genomic values from the G-BLUP with the genomic values provided by Delta-p under a selection index framework called Delta-p/G-BLUP.

The aim of this study was to propose the Delta-p and Delta-p/G-BLUP index, and to compare it with the traditional G-BLUP method. The efficiency of these methods was evaluated through simulation and real studies. The simulated dataset consisted of eight scenarios comprising the combination of two levels of heritability, two genetic architectures and two dominance status (absence and complete dominance). The empirical dataset consisted of 370 rice accessions phenotyped for seven yield traits, which were genotyped for 44,100 SNP markers.

# Materials and Methods

## Simulated datasets

The data set was simulated using the Real Breeding software program (Viana, 2011Viana, J.M.S. 2011. Program for molecular and quantitative data analysis Real Breeding = Programa para análises de dados moleculares e quantitativos Real Breeding. Universidade Federal de Viçosa, Viçosa, MG, Brazil (in Portuguese).), and the generation of data has been described by Azevedo et al. (2015)Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105.. A population of 5,000 individuals from 100 families, generated from the random mating of two linkage equilibrium populations, was subjected to five generations of random mating without selection, mutation nor migration. The final population is an advanced generation composite in Hardy-Weinberg equilibrium and linkage disequilibrium (LD).

The prediction of genomic breeding value based on thousands of single nucleotide polymorphisms (SNPs) depends on LD between markers and QTL. Thus, the LD value (Δ) in a composite population can be calculated from $Δab=(1−2θab4)(pa1−pa2)(pb1−pb2)$ (Kempthorne, 1973Kempthorne, O. 1973. An Introduction to Genetic Statistics. Iowa State University Press, Ames, IA, USA.; Viana, 2004Viana, J.M.S. 2004. Quantitative genetics theory for non-inbred populations in linkage disequilibrium. Genetics and Molecular Biology 27: 594-601.; Viana et al., 2016Viana, J.M.S.; Piepho, H.P.; Silva, F.F. 2016. Quantitative genetics theory for genomic selection and efficiency of breeding value prediction in open-pollinated populations. Scientia Agricola 73: 243-251.), where a and b are two SNPs, or two QTLs, or one SNP and one QTL, θ the frequency of recombinant gametes, and p1 and p2 the allele frequencies in the parental populations 1 and 2, respectively.

Azevedo et al. (2015)Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105. simulated a total of 2,000 equidistant SNP markers separated by 0.10 centimorgans (cM) across the ten equally sized chromosomes. One hundred of the simulated markers were QTLs randomly distributed in the genome and 95 % of the genetic variation was explained by the markers, according to the expression proposed by Goddard et al. (2011)Goddard, M.E.; Hayes, B.J.; Meuwissen, T.H.E. 2011. Using the genomic relationship matrix to predict the accuracy of genomic selection. Journal of Animal Breeding and Genetics 128: 409-421.. An amount of 1,000 individuals from 20 full-sib families were genotyped and phenotyped. This simulation study provides a small effective population size (Ne = 39.22), which, according to Resende (2015)Resende, M.D.V. 2015. Quantitative and Population Genetics = Genética Quantitativa e de Populações. Editora Suprema, Visconde do Rio Branco, MG, Brazil (in Portuguese)., is the typical value in elite breeding populations of plant species.

The simulated phenotypes considered two types of genetic architecture, one with polygenic inheritance, wherein the 100 QTLs were assumed with limited additive effects on the phenotype, and the other with mixed inheritance, the 5 QTLs with effects accounting for 50 % of the genetic variability and the small additive effects were assigned to the remaining 95 QTLs. The genotypic values for homozygotes were calculated considering Gmax = 100(m + a) as the upper limit and Gmin = 100(ma), as the lower limit where m was the mean of genotypic values and a the homozygote genotypic value. The phenotypic values were computed by adding genotypic values to the true population mean, and error effects sampled from a normal distribution $N(0,σe2)$, where the error variance $σe2$ was computed according to the heritability levels. The additive and dominance effects were assumed to be independents and were normally distributed with zero mean and genetic variance, depending on the assumed genetic architecture and the desired heritability level. In the simulation study, it was also observed that SNPs had minor allele frequency (MAF) greater than 5 %.

## Scenarios

Eight different scenarios were used in the analyses: two genetic architectures (polygenic and mixed inheritance); dominance status (absence and complete), two broad-sense heritability levels (around 0.20 and 0.30) for absence of dominance and two broad-sense heritability levels (around 0.30 and 0.50) for complete dominance. These values were chosen because, according to Azevedo et al. (2015)Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105., it was expected in these cases, that the genomic selection would be superior to phenotypic selection. The description of the scenarios is presented in Table 1.

Table 1
Scenarios with the respective true additive heritabilities ($ha2$), due to dominance ($hd2$) and total ($hg2$), genetic architectures (traits controlled by genes of small effect - polygenic inheritance, and traits controlled by genes of small and greater effect - mixed inheritance) and dominance status (absence and complete).

## Real datasets

Asian rice (Oryza sativa) is one of the most widely consumed foods in the world (Subudhi et al., 2006Subudhi, P.K.; Sasaki, T.; Khush, G.S. 2006. Rice. v.1, p. 1-9. In: Kole, C., ed. Genome mapping and molecular breeding in plants. Springer, Berlin, Germany.). Thus, the high population growth rate has justified the interest of researchers in developing high-yielding rice varieties. The real/empirical dataset consisted of 370 rice accessions phenotyped for seven yield traits, which were genotyped for 44,100 SNPs markers. The dataset is publicly available and is part of two projects, the OryzaSNP Project and the OMAP Project (Ammiraju et al., 2006Ammiraju, J.S.S.; Luo, M.; Goicoechea, J.L.; Wang, W.; Kudrna, D.; Mueller, C.; Talag, J.; Kim, H.; Sisneros N.B.; Blackmon, B.; Fang, E.; Tomkins, J.B.; Brar, D.; Mackill, D.; MacCouch, S.; Kurata, N.; Lambert, G.; Galbraith, D.W.; Arumuganathan, K.; Rao, K.; Walling, J.G.; Gill, N.; Yu, Y.; Sanmiguel, P.; Soderlund, C.; Jackson, S.; Wing, R.A. 2006. The Oryza bacterial artificial chromosome library resource: construction and analysis of 12 deep-coverage large-insert BAC libraries that represent the 10 genome types of the genus Oryza. Genome Research 16: 140-147.; Zhao et al., 2011Zhao, K.; Tung, C.W.; Eizenga, G.C.; Wright, M.H.; Ali, M.L.; Price, A.H., Norton, J. G.; Islam, A.R.; Reynolds, A.; Mezey, J.; McClung, A.M.; Bustamante, C.D.; McClung, A.M 2011. Genome-wide association mapping reveals a rich genetic architecture of complex traits in Oryza sativa. Nature Communications 2: 467-477.), and is available at http://ricediversity.org/.

Plantations were supervised throughout the access phase, from May to Oct of the years 2006 and 2007. A complete block design with two replications was used, in which the planting lines had a length equal to 5 m. The plants were spaced 25 cm apart and 0.50 m between rows. Rice accessions were evaluated in the field at Stuttgart, Arkansas coordinates, Lat: 034°49.5805’ N; Long −091°54.8357’ W, 100 m. Further details can be found in Zhao et al. (2011)Zhao, K.; Tung, C.W.; Eizenga, G.C.; Wright, M.H.; Ali, M.L.; Price, A.H., Norton, J. G.; Islam, A.R.; Reynolds, A.; Mezey, J.; McClung, A.M.; Bustamante, C.D.; McClung, A.M 2011. Genome-wide association mapping reveals a rich genetic architecture of complex traits in Oryza sativa. Nature Communications 2: 467-477.. The quality control procedures were made considering a call rate of 70 % and minor allele frequency (MAF) less than 1 %. After the quality control of the genomic database the total was 36,901 SNPs markers.

The seven traits evaluated were: (i) Panicle number per plant, (ii) Plant height, (iii) Panicle length, (iv) Primary panicle branch number, (v) Seed number per panicle, (vi) Florets per panicle, (vii) Panicle fertility.

## G-BLUP method

The Genomic Best Linear Unbiased Predictor (G-BLUP) method is based on the following linear mixed model:

$y = 1 μ + Z a + Z d + e ,$

where, y is a vector of phenotypes (N × 1, being N the number of genotyped and phenotyped individuals); µ the general mean and 1 the vector with dimension (N × 1); a the vector of additive genomic values (N × 1) with incidence matrix Z (N × N), whose assumed distribution is $a~N(0,Gaσa2)$, where $σa2$ is the genetic additive variance and Ga(N × N) the additive genomic relationship matrix; d the vector of dominance genomic values (N × 1) with incidence matrix Z (N × N), whose assumed distribution is $d~N(0, Gdσd2)$, where $σd2$ is the dominance variance and Gd (N × N) the dominance genomic relationship matrix; and e the random residual vector, assumed as $e~N(0,Iσe2)$, where $σe2$ is the residual variance and I an identity matrix.

According to Vitezica et al. (2013)Vitezica, Z.G.; Varona, L.; Legarra, A. 2013. On the additive and dominance variance and covariance of individuals within the genomic selection scope. Genetics 195: 1223-1230., the genomic relationship matrices for the additive and dominance effects, Ga and Gd, are given, respectively, by:

$G a = W W ' ∑ i = 1 n ( 2 p i q i ) a n d G d = S S ' ∑ i = 1 n ( 2 p i q i ) 2 ;$

where pi and qi are the allele frequencies of the locus i. The elements of W and S are given by:

$W = { i f M M , t h e n 2 → 2 q i f M m , t h e n 1 → q − p a n d S = i f m m , t h e n 0 → 2 p { i f M M , t h e n 0 → 2 q 2 i f M m , t h e n 1 → 2 p q i f m m , t h e n 0 → − 2 p 2$

## Delta-p method

The Delta-p method uses the concepts of changes in allele frequency due to selection and genetic gain theory (contrast between the averages of two subpopulations). The training population is initially divided into two subpopulations: one with individuals above the general mean (subpopulation one, with a higher average genetic value u1) and the other with individuals below the general mean (subpopulation two, with lower average genetic value u2). The difference between the average genetic values of the two subpopulations (u1u2) is due to the frequency of favorable alleles (p) among the two populations. Thus, (u1u2) is explained by Δp = p1p2, where Δp is the difference of allelic frequencies between the two subpopulations. The Δp values are calculated for each marker i, so that the positive signals are allocated as favorable. Thus, the allelic substitution effects (αi) of these markers are taken as positive, whereas for those with negative signs of Δp their αi's are assigned as negative.

In this context, the incidence matrix W for the vector of additive effects of markers is reparametrized, making the markers genotypes as 0 (bb), 1 (Bb) and 2 (BB). Here, the B allele is the favorable allele, and the BB or bb allocation is dictated by the αi signal. Thus, the correct allocation in BB or bb is probabilistic. In general, it is expected that the highest number of allocation errors will occur in loci with very minor effects.

The Delta-p method requires the estimation of the αi effects, and later the prediction of the genomic values (a). Taking Δpi as an indicator of the relative magnitude of αi (the higher the Δpi, the higher the chance of the αi being higher), it is possible to define the proportion Δpipm (with their respective positive or negative signals); where Δpm is the average delta-p computed using the absolute value of Δpi multiplied by the average α (αm), which results in the αi effects. Thus, $αi=(Δpi/Δpm) αm$.

The αm value can be obtained using the theoretical expression of genetic gain. For a given locus i in sub-population one, $u1i=(p1i−q1i) αi=(2p1i−1) αi$, where $i=1,2,…,n$; and in subpopulation two, $u2i=(p2i−q2i) αi=(2p2i−1) αi$. Thus, $(u1i−u2i)=(2p1i−2p2i) αi=2Δpiαi$. Summing all loci, $(u1−u2)=2Σ(Δpiαi)$. The quantity ∆pm is defined by Falconer (1989)Falconer, D.S. 1989. Introduction to Quantitative Genetics. Longman Scientific and Technical, New York, NY, USA. as: $Δpm=sq(1−q)$, where $s=αk/σp$, where s is the selective coefficient of the locus, assuming an additive-dominance model. Therefore, $Δpm=αkpq/σP2$ where k is the standardized selection intensity and $σP2$ the phenotypic variance. The genetic gain is defined as $Gs=kσa2/σP2$. Thus, $k/σP2=Gs/σa2$. Replacing Δpm in the expression it is possible to obtain $Δpm=αpqGs/σa2$.

The genetic gain associated with the first differential selection is given by $Gs1=ha2(u1−u0)$, whereas the second one is $Gs2=ha2(u0−u2)$. These gains are symmetrical and must be summed to obtain the total gain 2Gs, given by: $2Gs=Gs1+Gs2=ha2(u1−u2)$. Thus, $Gs=(12)ha2(u1−u2)$. Replacing Gs into ∆pm, $Δpm=αpqGs/σa2$ and $Δpm=αpq(1/2)ha2(u1−u2)/σa2$; where $σa2=2pqα2$. Replacing and dividing by ∆pm, we have:

$1 = [ α p q ( 1 / 2 ) h a 2 ( u 1 − u 2 ) ] / ( 2 p q α 2 Δ p m ) ,$

thus, $1=(1/2)ha2(u1−u2)/(2αΔpm)$. Isolating α from the last expression, $α=(1/2)ha2(u1−u2)/(2Δpm)$) for n locus and the average a per locus is given by:

$α m = 0.5 h a 2 ( u 1 − u 2 ) / [ n m a r k e r s 2 Δ p m ] .$

Therefore, the additive markers effects (α) are non-parametrically estimated in the training population (individuals with known phenotypes and genotypes). Based on these effects, the genetic values (a) are estimated by the following expression: $a^=Wa^$, being W the incidence matrix of the additive genetic effects of the markers. This approach does not demand an iterative computational method and uses only genetic distance concepts (magnitude of ∆pi) and the genetic gains of two subpopulations. Additionally, it also does not use differential shrinkage based on allele frequencies, which benefits the loci with higher MAFs. In summary, the Delta-p method can be implemented by taking the followings steps: i) subdivision of the training population into two, according to the phenotype adjusted for systematic effects; ii) calculation of ∆pi and ∆pm; iii) calculation of αm; and iv) calculation of $αi=(Δpi/Δpm) αm$.

For the dominance effects, the quantity Δ(2pq)i (difference between frequencies of heterozygous genotypes of subpopulations one and two) is used instead of ∆pi. The better the phenotype and the higher the frequency of heterozygous genotypes in the population, the more favorable the value of the genotype corrected for the additive genetic effects (i.e., higher the of dominance deviations). In order to infer the effects of dominance, the interest is to find this association (frequency of heterozygous and values of phenotypes). For this, the phenotype must be returned as the frequency of heterozygous, corrected for additive effects. The incidence matrix S is inherent to the frequency of heterozygous among the three genotype classes. Then, gm is defined as the genotype value of the heterozygous, which is corrected for the additive effect (αi). It gives the dominance deviation of the locus i (δi), and the genomic dominance value (vector d) is defined as $d^=Sδ^$, where S is the previously mentioned incidence matrix.

## Delta-p/G-BLUP Index method

Based on the principle of combined selection (Lush, 1947Lush, J.L. 1947. Family merit and individual merit as basis for selection. America Naturalist 81: 241-261.; Falconer, 1989Falconer, D.S. 1989. Introduction to Quantitative Genetics. Longman Scientific and Technical, New York, NY, USA.), the G-BLUP can be jointly used with the Delta-p method. A selection index of the form $I=b1a^1+b2a^2$ is defined based on information from predicted additive genomic values from G-BLUP ($a^1$) and Delta-p ($a^2$), weighted by b1 and b2 coefficients, respectively. The index can be also represented by $I=b1f(y,W)+b2f((u1−u2),Δp)$; where y is the vector of individual phenotypic data, W the centralized allelic dosage matrix, (u1u2) the contrast of averages from subpopulations and ∆p the differential frequency between subpopulations. The weights of the index are given by b = P−1C as follows:

$C = [ σ a 1 2 r a ^ 1 a 2 σ a 2 2 r a ^ 2 a 2 2 ] , P = [ σ a 1 2 r a ^ 1 a 2 cov ( a 1 , a 2 ) r a ^ 1 a 2 r a ^ 2 a 2 cov ( a 1 , a 2 ) r a ^ 1 a 2 r a ^ 2 a 2 σ a 2 2 r a ^ 2 a 2 ] = [ σ a 1 2 r a ^ 1 a 2 σ Δ 2 r a ^ 1 a 2 r a ^ 2 a 2 σ a 2 2 r a ^ 1 a 2 r a ^ 2 a 2 σ a 2 2 r a ^ 2 a 2 ]$

where C is a matrix of variance and reliability (r2) of estimated genomic values through Delta-p and G-BLUP, P a covariance matrix between the estimated genomic values through Delta-p and G-BLUP and a the true additive genetic value of individuals.

Dividing P and C by $σa12$, and considering $Δ2=σa22/σa12$ it is possible to obtain:

$C = [ r a ^ 1 a 2 Δ 2 r a ^ 2 a 2 2 ] ; P = [ r a ^ 1 a 2 Δ 2 r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 r a ^ 2 a 2 ] .$

The weights (b1 and b2) and accuracy (rIa) between the index and the true additive genetic value a are given respectively by:

$b 1 = 1 − r a ^ 2 a 2 Δ 2 1 − r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 ; b 2 = 1 − r a ^ 1 a 2 1 − r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 ; r I a = [ r a ^ 1 a 2 + r a ^ 2 a 2 Δ 2 − 2 r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 1 − r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 ] 1 / 2 = [ 1 − ( 1 − r a ^ 1 a 2 ) ( 1 − r a ^ 2 a 2 Δ 2 ) 1 − r a ^ 1 a 2 r a ^ 2 a 2 Δ 2 ] 1 / 2 .$

Under this approach, the reliabilities of $a^1$($ra^1a2$, obtained by G-BLUP) and $a^2$ ($ra^2a2$ obtained by Delta-p) and the proportion of variances of $a^2$ and $a^1$ are required. If $Δ2=σa22/σa12$ provides a value higher than one, then ∆2 must be at least 1. Another option is to calculate cov(a1, a2) directly through $cov(a^1,a^2)$). Considering that the dominance effects of the methods are analogous to that described above only the genomic value due to the dominance should be used ($d^$) instead of the additive genomic value ($a^$) in derived equations.

## Computer resources

All computational implementations of the proposed methods were performed using the R software program (version 3.3.1). The G-BLUP was performed in the rrBLUP package through the mixed.solve function. The algorithm used for the development of the Delta-p method is available in http://www.novoscursos.ufv.br/departamentos/ufv/det/www/?page_id=1374.

## Comparisons of the proposed methods Simulated dataset

Each type of population was simulated 10 times under the same parameter settings, which preserved the same features and provided samples that were effectively from the same conceptual population. Nine replicates were used as training populations, and one replicate was used as a validation population. The estimations based on each of the nine replicates were validated by obtaining estimates of the efficiency measures. The efficiency measures of genomic prediction were given by the accuracies ($ra^a$ and $rd^d$), bias ($bya^$e$byd^$), recovered additive and dominance genomic heritabilities ($ha2$e$hd2$) and relative efficiency (RE).

The accuracy is defined as the correlation between the predicted values and the parametric genetic values; and the prediction bias is defined as the regression coefficient associated with this correlation. For regression coefficients below one (< 1) it is understood that predicted values had been overestimated, whereas for coefficients above one (> 1), it was concluded that predicted values had been underestimated. The recovered additive genomic heritability is given by

$h a 2 = σ a 2 σ a 2 + σ d 2 + σ e 2 ,$

where $σa2=∑i=1n2piqimi2$ is the additive genomic variance, $mi2$ the square of the ith marker effect with allelic frequencies equal to pi and qi. The genomic heritability due to dominance is given by

$h d 2 = σ d 2 σ a 2 + σ d 2 + σ e 2 ,$

where $σd2=∑i=1n(2piqidi)2$ and di is the genotypic value of the heterozygous. The relative efficiency is given by the ratio between the accuracies from the compared methods. All these measures were obtained for each replicate in each scenario and the general results were reported as average values.

## Real dataset

We obtained the measures of efficiency using five-fold cross validation (k = 5), in which in each fold the training population consisted of 296 individuals (80 % of total) and the validation population consisted of 74 individuals (20 % of total). The efficiency measures used were the predictive ability ($ra^y$), which consisted of the correlation between the estimated additive genomic values and the phenotypic values of the validation population, and prediction bias.

# Results and Discussion

The average and the respective estimated standard deviations for the accuracy, bias, heritabilities and relative efficiency were obtained by both Delta-p and G-BLUP and are shown in Table 2. We note an increase in accuracy with the use of the G-BLUP. The correlation between the genomic values predicted by the Delta-p and the G-BLUP method was high ($ra^pa^g$), ranging from 0.80 to 0.88. Although the G-BLUP is widely applied to GWS and recommended for the prediction of additive and dominant genomic values (Azevedo et al., 2015Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105.; De los Campos et al., 2013De los Campos, G.; Hickey, J.M.; Pong-Wong, R.; Daetwyler, H.D.; Callus, M.P.L. 2013. Whole genome regression and prediction methods applied to plant and animal breeding. Genetics 193: 327-345.; Gianola, 2013Gianola, D. 2013. Priors in whole-genome regression: the Bayesian alphabet returns. Genetics 194: 573-596.; Gianola et al., 2009Gianola, D.; De los Campos, G.; Hill, W.G.; Manfredi, E.; Fernando, R. 2009. Additive genetic variability and the Bayesian alphabet. Genetics 183: 347-363.), this study proposed combining the values derived from the Delta-p and G-BLUP methods to construct an index giving different weights to each of them which can be seen as an improvement in the G-BLUP additive and dominant genomic predictions. For the additive model, the increase in accuracy using the Delta-p/G-BLUP index (0.80 to 0.87) was 0.03 units on average compared to the accuracy of the G-BLUP method (0.77 to 0.85). For the additive-dominant model, this increase was more evident, being 0.08 units. The Delta-p/G-BLUP index, considering the additive model, provided a relative average efficiency of 103 % to 105 %, while in the additive-dominant model this efficiency increased to 109 % to 116 %. It is noted that the Delta-p/G-BLUP index improved the prediction of additive and dominant genomic values through G-BLUP, since the index provided greater accuracy. It is important to report that this advantage in terms of extra percentage points has no additional computational costs.

Table 2
Average and standard deviation of accuracy ($ra^a$), bias ($bya^$), additive heritability ($ha2$), correlation between additive genomic values estimated by Delta-p and Genomic Best Linear Unbiased Predictor - G-BLUP ($ra^1a^2$), proportion between the additive variances estimated (Δ2), weights (b1 and b2), accuracy (rIa) and bias (byIa) of the Delta-p/G-BLUP index method, relative efficiency (RE) between the index and G-BLUP considering simulated dataset.

As already mentioned, the Delta-p/G-BLUP index method is based on the combined selection index. According to Lush (1947)Lush, J.L. 1947. Family merit and individual merit as basis for selection. America Naturalist 81: 241-261. and Falconer (1989)Falconer, D.S. 1989. Introduction to Quantitative Genetics. Longman Scientific and Technical, New York, NY, USA., this process outperformed other selection methods, for example, individual or among or within family selection. Similarly, in the present study, the Delta-p/G-BLUP index method, which is based on the combined selection index, outperformed the G-BLUP and Delta-p methods. Gains of 5 % in accuracy are already significant in genetic improvement, often equivalent to the gain that is obtained in a complete cycle of improvement (Resende et al., 2015Resende, M.D.V.; Ramalho, M.A.P.; Guilherme, S.R.; Abreu, A.F.B. 2015. Multi generation index in the within progenies bulk method for breeding of self-pollinated plants. Crop Science 55: 1202-1211.).

The differences between allele frequencies were effective in estimating the marker effects. This result corroborates the principle of the Delta-p method. Weller et al. (2014)Weller, J.I.; Glick, G.; Ezra, E.; Seroussi, E.; Shemesh, M.; Zeron, Y.; Ron, M. 2014. Predictive ability of selected subsets of single nucleotide polymorphisms (SNPs) in a moderately sized dairy cattle population. Animal 8: 208-216. and Weller (2016)Weller, J.I. 2016. Genomic Selection in Animals. John Wiley, Hoboken, NJ, USA. reported that marker selection based on differences between allelic frequencies between groups of young and older bulls, i.e. between two contrasting subpopulations, can be effective in genomic selection in dairy cattle.

Considering the additive model, it was observed that for all scenarios, the G-BLUP method was the one that presented an estimated additive heritability close to the parametric value (Table 2). Considering the additive-dominant model, it was observed that for the scenarios whose traits are controlled by genes of small effects (scenarios 5 and 6), the Delta-p method presented the additive heritability estimates closest to the parametric values. On the other hand, for the scenarios whose traits are controlled by genes of both small and sizeable effects (scenarios 7 and 8), the G-BLUP method presented better performance.

Considering the additive model, the G-BLUP and Delta-p methods had non-biased predicted additive genomic values (Table 2), except for the Delta-p method in scenarios 1 and 3, where the method underestimated these values ($bya^>1$). Considering the additive-dominant model, the G-BLUP method provided non-biased additive genomic values, except for scenario 8, where the method overestimated these values ($bya^<1$). In contrast, the Delta-p method presented predicted additive genomic values which were overestimated for scenarios 5, 6 and 8; and underestimated for scenario 7. Already, the Delta-p/G-BLUP index overestimated the predicted additive genomic values for all scenarios considered (byIa < 1). However, the index had biases considerably closer to one in the additive model than in the additive-dominant model (Table 2).

Table 3 shows the average results obtained through a simulated dataset, the methods for the genomic values due to dominance considering the additive-dominant model and also the weights, accuracy and efficiency obtained through the Delta-p/G-BLUP index. In this case, the improvement in the prediction of genomic values with the use of the index (Table 3) was less evident. Only 0.01 units improved in prediction accuracy compared to the accuracy of G-BLUP. Bennewitz and Meuwissen (2010)Bennewitz, J.; Meuwissen, T.H.E. 2010. The distribution of QTL additive and dominance effects in porcine F2 crosses. Journal of Animal Breeding and Genetics 127: 171-179., Hill et al. (2008)Hill, W.G.; Goddard, M.E.; Visscher, P.M. 2008. Data and theory point to mainly additive genetic variance for complex traits. PLoS Genetics 4: e1000008., and Wellmann and Bennewitz (2012)Wellmann, R.; Bennewitz, J. 2012. Bayesian models with dominance effects for genomic evaluation of quantitative traits. Genetics Research 94: 21-37. discuss the relevance of the inclusion of dominance in Quantitative Genomics. However, according to the results reported by Azevedo et al. (2015)Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105., genomic values due to dominance are difficult to predict and are consequently associated with low accuracy values and higher bias. This corroborates the low accuracy level recorded for the Delta-p/G-BLUP index and G-BLUP when considering only dominance genomic values.

Table 3
Average and standard deviation of accuracy ($rd^d$), bias ($byd^$), additive heritability ($hd2$), correlation between dominant genomic values estimated by Delta-p and Genomic Best Linear Unbiased Predictor - G-BLUP ($rd^1d^2$), proportion between the dominant variances estimated (∆2), weights (b1 and b2), accuracy (rId) and bias (byId) of the Delta-p/G-BLUP index method, relative efficiency (RE) between the index and G-BLUP considering simulated dataset.

Although the dominance genomic values, considering the additive-dominant model, were inaccurately estimated, it can be seen that the index and the G-BLUP accuracy were similar. This result agrees with Azevedo et al. (2015)Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105., which used several methods applied to GWS and obtained results of similar dominance accuracy. Thus, it should be noted that the main criteria for the comparison of methods in these scenarios is the recovered heritabilities, due to dominance and prediction bias.

When analyzing the heritabilities due to dominance, it was observed that for the four scenarios, the closest estimates were provided by the G-BLUP method. However, in scenarios 5 and 8, the heritabilities, due to dominance, were slightly underestimated and in scenario 7 the heritabilities, due to dominance, were slightly overestimated. The Delta-p method underestimated the genomic values due to dominance in the scenarios whose traits were controlled by genes with small effects (scenarios 5 and 6) and overestimated the scenarios whose traits were controlled by genes with both small and sizeable effects. The G-BLUP overestimated the values in all scenarios, with bias being considerably lower in relation to the Delta-p method. In addition, the Delta-p/G-BLUP index also overestimated the genomic values due to dominance for all scenarios considered (byIa < 1) (Table 3).

Although we have previously analyzed the additive genomic values and because of the dominance, considering the additive-dominant model, it is also important to present the results referring to the total genomic value (a + d), since this is related to individual performance (Table 4). Analogous to Tables 2 and 3, the G-BLUP presented higher accuracy in the prediction of total genomic values, when compared to the Delta-p method. The increase in accuracy using the Delta-p/G-BLUP index (0.65 to 0.71) was 0.07 units, on average, compared to the accuracy of the G-BLUP method (0.56 to 0.65), showing again the efficiency of the index. The Delta-p/G-BLUP index provided a relevant relative average efficiency that was 108 % to 114 % (Table 4). It was observed that, under the four scenarios considered, the G-BLUP had a better performance in recovering the total genomic heritability than the Delta-p. However, these were still underestimated (Table 4). In addition, the three methods produced predicted total genomic values that were overestimated ($byg^$ < 1). Although the Delta-p/G-BLUP index presented biased total genomic values, this method has the great advantage of obtaining greater accuracy when compared to the G-BLUP method.

Table 4
Average and standard deviation of accuracy ($rg^g$), bias ($byg^$), additive heritability ($hg2$), correlation between total genomic values estimated by Delta-p and Genomic Best Linear Unbiased Predictor - G-BLUP ($rg^1g^2$), proportion between the total variances estimated (∆2), weights (b1 and b2), accuracy (rIg) and bias (byIg) of the Delta-p/G-BLUP index method, relative efficiency (RE) between the index and G-BLUP considering real dataset.

To elucidate the importance of the Delta-p/G-BLUP index for predicting genomic selection, the results for an applied breeding program of Asian rice are presented in Table 5. Analogous to previous results considering additive model (Table 2), G-BLUP showed more accurate genomic values than Delta-p. The Delta-p method presented values showing overestimated additive genomic values ($bya^<1$), except for the Panicle number per plant, while the G-BLUP presented values that were less biased. The Delta-p/G-BLUP index, considering the real dataset, provided a relative average efficiency of 101 % to 116 %, being on average 8 % more efficient than the G-BLUP. Eight % is highly expressive and under genomic selection performed annually, such as in rice cultivation, these gains are cumulative and grow rapidly (Resende et al., 2015Resende, M.D.V.; Ramalho, M.A.P.; Guilherme, S.R.; Abreu, A.F.B. 2015. Multi generation index in the within progenies bulk method for breeding of self-pollinated plants. Crop Science 55: 1202-1211.). Thus, the index proposed in this study can be satisfactorily implemented in rice breeding programs.

Table 5
Predictive ability ($ra^y$), bias ($bya^$), additive heritability ($ha2$), correlation between additive genomic values estimated by Delta-p and Genomic Best Linear Unbiased Predictor - G-BLUP ($ra^1a^2$), proportion between the additive variances estimated (∆2), weights (b1 and b2), predictive ability (rIa) and bias (byI) of the Delta-p/G-BLUP index method, relative efficiency (RE) between the index and G-BLUP considering real dataset.

Thus, the Delta-p/G-BLUP index is defined based on information from predicted genomic values from G-BLUP and Delta-p. This approach is rarely applied in the context of GWS. This method, as much for the simulated data as for the rice data, showed itself to be superior to the G-BLUP method. In addition, it led to genomic, additive and dominance values that are more accurate than the traditional G-BLUP and required reduced computational time. It can also be easily implemented in the context of GWS.

# Conclusions

The proposed method, the Delta-p/G-BLUP index, is easy to implement in the context of genomic selection. It demands reduced computational cost and leads to genomic, additive, dominance and total values that are more accurate than the G-BLUP method. On the other hand, the Delta-p method provided lower values of accuracy and dominance bias estimates than the G-BLUP method. However, for additive effects, under almost all scenarios evaluated, the estimates were more biased for the Delta-p model than for the G-BLUP. For the traits evaluated in rice it was concluded that the Delta-p/G-BLUP index was superior to the G-BLUP method, and increased predictive ability.

# Acknowledgments

The authors thank the following Brazilian funding organizations: the Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES) and Fundação Arthur Bernardes (FUNARBE).

# References

• Ammiraju, J.S.S.; Luo, M.; Goicoechea, J.L.; Wang, W.; Kudrna, D.; Mueller, C.; Talag, J.; Kim, H.; Sisneros N.B.; Blackmon, B.; Fang, E.; Tomkins, J.B.; Brar, D.; Mackill, D.; MacCouch, S.; Kurata, N.; Lambert, G.; Galbraith, D.W.; Arumuganathan, K.; Rao, K.; Walling, J.G.; Gill, N.; Yu, Y.; Sanmiguel, P.; Soderlund, C.; Jackson, S.; Wing, R.A. 2006. The Oryza bacterial artificial chromosome library resource: construction and analysis of 12 deep-coverage large-insert BAC libraries that represent the 10 genome types of the genus Oryza. Genome Research 16: 140-147.
• Azevedo, C.F.; Resende, M.D.V.; Silva, F.F.; Viana, J.M.S.; Valente, M.S.F.; Resende Junior, M.F.R.; Muñoz, P. 2015. Ridge, Lasso and Bayesian additive-dominance genomic models. BMC Genetics 16: 105.
• Bennewitz, J.; Meuwissen, T.H.E. 2010. The distribution of QTL additive and dominance effects in porcine F2 crosses. Journal of Animal Breeding and Genetics 127: 171-179.
• De los Campos, G.; Hickey, J.M.; Pong-Wong, R.; Daetwyler, H.D.; Callus, M.P.L. 2013. Whole genome regression and prediction methods applied to plant and animal breeding. Genetics 193: 327-345.
• Falconer, D.S. 1989. Introduction to Quantitative Genetics. Longman Scientific and Technical, New York, NY, USA.
• Gianola, D. 2013. Priors in whole-genome regression: the Bayesian alphabet returns. Genetics 194: 573-596.
• Gianola, D.; De los Campos, G.; Hill, W.G.; Manfredi, E.; Fernando, R. 2009. Additive genetic variability and the Bayesian alphabet. Genetics 183: 347-363.
• Goddard, M.E. 2009. Genomic selection: prediction of accuracy and maximization of long term response. Genetics 136: 345-357.
• Goddard, M.E.; Hayes, B.J.; Meuwissen, T.H.E. 2011. Using the genomic relationship matrix to predict the accuracy of genomic selection. Journal of Animal Breeding and Genetics 128: 409-421.
• Habier, D.; Fernando, R.L.; Dekkers, J.C.M. 2007. The impact of genetic relationship on genome-assisted breeding values. Genetics 117: 2389-2397.
• Hill, W.G.; Goddard, M.E.; Visscher, P.M. 2008. Data and theory point to mainly additive genetic variance for complex traits. PLoS Genetics 4: e1000008.
• Kempthorne, O. 1973. An Introduction to Genetic Statistics. Iowa State University Press, Ames, IA, USA.
• Lush, J.L. 1947. Family merit and individual merit as basis for selection. America Naturalist 81: 241-261.
• Meuwissen, T.H.E.; Hayes, B.J.; Goddard, M.E. 2001. Prediction of total genetic value using genome-wide dense marker maps. Genetics 157: 1819-1829.
• Resende, M.D.V. 2015. Quantitative and Population Genetics = Genética Quantitativa e de Populações. Editora Suprema, Visconde do Rio Branco, MG, Brazil (in Portuguese).
• Resende, M.D.V.; Ramalho, M.A.P.; Guilherme, S.R.; Abreu, A.F.B. 2015. Multi generation index in the within progenies bulk method for breeding of self-pollinated plants. Crop Science 55: 1202-1211.
• Subudhi, P.K.; Sasaki, T.; Khush, G.S. 2006. Rice. v.1, p. 1-9. In: Kole, C., ed. Genome mapping and molecular breeding in plants. Springer, Berlin, Germany.
• Van Raden, P.M. 2008. Efficient methods to compute genomic predictions. Journal of Dairy Science 91: 4414-4423.
• Viana, J.M.S. 2004. Quantitative genetics theory for non-inbred populations in linkage disequilibrium. Genetics and Molecular Biology 27: 594-601.
• Viana, J.M.S. 2011. Program for molecular and quantitative data analysis Real Breeding = Programa para análises de dados moleculares e quantitativos Real Breeding. Universidade Federal de Viçosa, Viçosa, MG, Brazil (in Portuguese).
• Viana, J.M.S.; Piepho, H.P.; Silva, F.F. 2016. Quantitative genetics theory for genomic selection and efficiency of breeding value prediction in open-pollinated populations. Scientia Agricola 73: 243-251.
• Vitezica, Z.G.; Varona, L.; Legarra, A. 2013. On the additive and dominance variance and covariance of individuals within the genomic selection scope. Genetics 195: 1223-1230.
• Weller, J.I. 2016. Genomic Selection in Animals. John Wiley, Hoboken, NJ, USA.
• Weller, J.I.; Glick, G.; Ezra, E.; Seroussi, E.; Shemesh, M.; Zeron, Y.; Ron, M. 2014. Predictive ability of selected subsets of single nucleotide polymorphisms (SNPs) in a moderately sized dairy cattle population. Animal 8: 208-216.
• Wellmann, R.; Bennewitz, J. 2012. Bayesian models with dominance effects for genomic evaluation of quantitative traits. Genetics Research 94: 21-37.
• Zhao, K.; Tung, C.W.; Eizenga, G.C.; Wright, M.H.; Ali, M.L.; Price, A.H., Norton, J. G.; Islam, A.R.; Reynolds, A.; Mezey, J.; McClung, A.M.; Bustamante, C.D.; McClung, A.M 2011. Genome-wide association mapping reveals a rich genetic architecture of complex traits in Oryza sativa. Nature Communications 2: 467-477.

### Edited by

Edited by: Roberto Fritsche Neto

# Publication Dates

• Publication in this collection
Jul-Aug 2019

# History

14 Oct 2017
• Accepted
07 Apr 2018
Escola Superior de Agricultura "Luiz de Queiroz" USP/ESALQ - Scientia Agricola, Av. Pádua Dias, 11, 13418-900 Piracicaba SP Brazil, Phone: +55 19 3429-4401 / 3429-4486 - Piracicaba - SP - Brazil
E-mail: scientia@usp.br