Broad-sense heritability in plant breeding

Maria Belen Kistner & Flavio Lozano-Isla

Broad-sense heritability (\(H^2\))

Broad-sense heritability (\(H^2\)) is defined as the proportion of phenotypic variance that is attributable to an overall genetic variance for the genotype (Schmidt et al., 2019b). There are usually additional interpretations associated with \(H^2\): (i) It is equivalent to the coefficient of determination of a linear regression of the unobservable genotypic value on the observed phenotype; (ii) It is also the squared correlation between predicted phenotypic value and genotypic value; and (iii) It represents the proportion of the selection differential (\(S\)) that can be realized as the response to selection (\(R\)) (Falconer and Mackay, 2005).

There are two main reasons why heritability on an entry-mean basis is of interest in plant breeding (Schmidt et al., 2019a):

  1. It is plugged into the breeder’s Equation to predict the response to selection.
  2. It is a descriptive measure used to assess the usefulness and precision of results from cultivar evaluation trials.

Breeder´s equation

\[\Delta G=H^2S\] Where:

Usual Problems

In practice, most trials are conducted in a multienvironment trial (MET) presente unbalanced data as not all cultivars are tested at each environment or simply when plot data is lost or when the number of replicates at each location varies between genotypes (Schmidt et al., 2019b). However, the standard method for estimating heritability implicitly assumes balanced data, independent genotype effects, and homogeneous variances.

How calculate the Heritability?

According Schmidt et al. (2019a), the variance components could be calculated in two ways:

1) Two stages approach

For the two stage approach, in the first stage each experiment is analyzed individually according their experiment design (Lattice, CRBD, etc) (Zystro et al., 2018). And for the second stage environments are denotes a year-by-location interaction. This approach assumes a single variance for genotype-by-environment interactions (GxE), even when multiple locations were tested across multiple years.

Model

\[y_{ikt}=\mu\ +\ G_i+E_t+GxE_{it}+\varepsilon_{ikt}\]

Phenotypic variance

\[\sigma_p^2=\sigma_g^2+\frac{\sigma_{g\cdot e}^2}{n_e}+\frac{\sigma_{\varepsilon}^2}{n_e\cdot n_r}\]

2) One stage approach

For the one stage approach only one model is used for the MET analysis. The environmental effects are included via separate year, and location main interaction effects.

\[y_{ikt}=\mu+G_i+Y_m+E_q+YxE_{mq}+GxY_{im}+GxE_{iq}+GxYxE_{imq}+\varepsilon_{ikmq}\]

Phenotypic variance

\[\sigma_p^2=\sigma_g^2+\frac{\sigma_{g\cdot e}^2}{n_e}+\frac{\sigma_{g\cdot y}^2}{n_y}+\frac{\sigma_{g\cdot y\cdot e}^2}{n_y\cdot n_e}+\ \frac{\sigma_{\epsilon}^2}{n_e\cdot n_y\cdot n_r}\]

Differentes heritability calculations

Differentes heritability calculation
Standart Cullis Piepho
\(H^2=\frac{\sigma_g^2}{\sigma_p^2}=\frac{\Delta G}{S}\) \(H_{Cullis}^2=1-\frac{\overline{V}_{\Delta..}^{^{BLUP}}}{2\cdot\sigma_g^2}\) \(H_{Piepho}^2=\frac{\sigma_g^2}{\sigma_g^2+\frac{\overline{V}_{\Delta..}^{BLUE}}{2}}\)

Heritability function in the package

For calculate the standard heritability in MET experiments the number of location and replication should be include manually in the function H2cal(). In the case of difference number of replication in each experiments, take the maximum value (often done in practice) (Schmidt et al., 2019b).

For remove the outliers the function implemented is the Method 4 used for Bernal-Vasquez et al. (2016): Bonferroni-Holm using re-scaled MAD for standardizing residuals (BH-MADR).

Load packages

library(inti)
library(agridat)

H2cal function

 dt <- john.alpha
 hr <- H2cal(data = dt
            , trait = "yield"
            , gen.name = "gen"
            , rep.n = 3
            , fix.model = "rep + (1|rep:block) + gen"
            , ran.model = "rep + (1|rep:block) + (1|gen)"
            , emmeans = TRUE
            , plot_diag = TRUE
            , plot_dots = "rep"
            , outliers.rm = FALSE
            )

Model information

hr$model %>% summary()
#> Linear mixed model fit by REML ['lmerMod']
#> Formula: yield ~ rep + (1 | rep:block) + (1 | gen)
#>    Data: dt.rm
#> Weights: weights
#> 
#> REML criterion at convergence: 93.2
#> 
#> Scaled residuals: 
#>      Min       1Q   Median       3Q      Max 
#> -1.78956 -0.37299  0.02882  0.37245  2.65803 
#> 
#> Random effects:
#>  Groups    Name        Variance Std.Dev.
#>  gen       (Intercept) 0.14290  0.3780  
#>  rep:block (Intercept) 0.07022  0.2650  
#>  Residual              0.08162  0.2857  
#> Number of obs: 72, groups:  gen, 24; rep:block, 18
#> 
#> Fixed effects:
#>             Estimate Std. Error t value
#> (Intercept)   4.5182     0.1451  31.136
#> repR2         0.2978     0.1738   1.714
#> repR3        -0.4140     0.1738  -2.382
#> 
#> Correlation of Fixed Effects:
#>       (Intr) repR2 
#> repR2 -0.599       
#> repR3 -0.599  0.500

Variance component

hr$tabsmr %>% kable(caption = "Variance component table")
Variance component table
variable rep geno env year mean std min max V.g V.gxl V.gxy V.e h2.s h2.c h2.p
yield 3 24 1 1 4.479517 0.4153925 3.4992 5.107699 0.1429021 0 0 0.0816171 0.8400678 0.7979968 0.7966375

Best Linear Unbiased Predictors (BLUPs)

hr$blups %>% kable(caption = "BLUPs")
BLUPs
gen yield
G01 5.019434
G02 4.523213
G03 3.733687
G04 4.524376
G05 4.993200
G06 4.562890
G07 4.209302
G08 4.580480
G09 3.708318
G10 4.428877
G11 4.321815
G12 4.744008
G13 4.749915
G14 4.761650
G15 4.942950
G16 4.719215
G17 4.596328
G18 4.408069
G19 4.807826
G20 4.179281
G21 4.774382
G22 4.542339
G23 4.341252
G24 4.265192

Best Linear Unbiased Estimators (BLUEs)

hr$blues %>% kable(caption = "BLUEs")
BLUEs
gen yield SE df lower.CL upper.CL
G01 5.107699 0.1990122 44.32700 4.706700 5.508699
G02 4.478532 0.1990122 44.32700 4.077533 4.879531
G03 3.499200 0.1990122 44.32700 3.098200 3.900199
G04 4.490094 0.1990122 44.32700 4.089095 4.891094
G05 5.037210 0.1989050 44.28463 4.636417 5.438004
G06 4.536662 0.1989050 44.28463 4.135868 4.937456
G07 4.111136 0.1990122 44.32700 3.710137 4.512135
G08 4.527634 0.1990122 44.32700 4.126634 4.928633
G09 3.502181 0.1989050 44.28463 3.101387 3.902975
G10 4.373200 0.1990122 44.32700 3.972200 4.774199
G11 4.283264 0.1990122 44.32700 3.882265 4.684263
G12 4.755276 0.1990122 44.32700 4.354277 5.156276
G13 4.757913 0.1989050 44.28463 4.357119 5.158707
G14 4.775662 0.1989050 44.28463 4.374868 5.176456
G15 4.969112 0.1989050 44.28463 4.568318 5.369905
G16 4.730131 0.1990122 44.32700 4.329132 5.131130
G17 4.602612 0.1989050 44.28463 4.201819 5.003406
G18 4.361692 0.1990122 44.32700 3.960693 4.762692
G19 4.840328 0.1989050 44.28463 4.439534 5.241122
G20 4.039985 0.1989050 44.28463 3.639191 4.440779
G21 4.795007 0.1989050 44.28463 4.394214 5.195801
G22 4.527545 0.1989050 44.28463 4.126751 4.928339
G23 4.252449 0.1989050 44.28463 3.851655 4.653243
G24 4.153874 0.1990122 44.32700 3.752875 4.554873

References

Bernal-Vasquez, A.-M., H.-F. Utz, and H.-P. Piepho. 2016. Outlier detection methods for generalized lattices: A case study on the transition from ANOVA to REML. Theoretical and Applied Genetics 129(4): 787–804. doi: 10.1007/s00122-016-2666-6.

Falconer, D.S., and T.F. Mackay. 2005. Introduction to quantitative genetics (Pearson Prentice Hall, editor). Fourth.

Schmidt, P., J. Hartung, J. Bennewitz, and H.-P. Piepho. 2019a. Heritability in Plant Breeding on a Genotype-Difference Basis. Genetics 212(4): 991–1008. doi: 10.1534/genetics.119.302134.

Schmidt, P., J. Hartung, J. Rath, and H.-P. Piepho. 2019b. Estimating Broad-Sense Heritability with Unbalanced Data from Agricultural Cultivar Trials. Crop Science 59(2): 525–536. doi: 10.2135/cropsci2018.06.0376.

Zystro, J., M. Colley, and J. Dawson. 2018. Alternative Experimental Designs for Plant Breeding. Plant Breeding Reviews. John Wiley & Sons, Ltd. pp. 87–117