parametric H2 bigger than 1 upon selection

[marcio-resende]

Hi all,

I noticed that when I calculate the heritability of the target population in each different cycle of the breeding, sometimes the parametric broad-sense h2, calculated as (varG/varP), is bigger than 1. This happens particularly when we apply selection of individuals and within families. See a code below that reproduces this issue in an F4 population.

Any suggestions of what may be happening or if I am doing anything wrong? Best regards,

Marcio

> library(AlphaSimR)
> nQtl = 300
> nSnp = 3000
> nParents = 200
> founderPop <- runMacs(nInd = nParents, nChr = 10, segSites = nQtl+nSnp, inbred = TRUE,
+                       species = "MAIZE", split = 30, ploidy = 2L,
+                       manualCommand = NULL, manualGenLen = NULL)
> SP <- SimParam$new(founderPop)
> SP$restrSegSites(nQtl,nSnp)
> SP$addTraitADG(nQtl,mean=0,var=30,
+                varGxE=20,meanDD=0,varDD=0.3)
> SP$setVarE(h2=0.65)
> Parents_public = newPop(founderPop[1:100])
> F1_public=randCross(Parents_public, 50)
> 
> F2_public=self(F1_public, nProgeny = 25, parents = NULL, keepParents = FALSE)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.4690354
> F2_public = selectWithinFam(F2_public, nInd = 6, selectTop = TRUE)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.8819334
> 
> F3_public=self(F2_public, nProgeny = 8, parents = NULL, keepParents = FALSE) 
> F3_public = selectWithinFam(F3_public, nInd = 2, selectTop = TRUE)
> 
> F4_public=self(F3_public, nProgeny = 10, parents = NULL, keepParents = FALSE) 
> V1 = varG(F4_public)/varP(F4_public)
> F4_public = setPheno(F4_public,reps=2)
> V2 = varG(F4_public)/varP(F4_public)
> F4_public = selectWithinFam(F4_public, nInd = 1, selectTop = TRUE)
> V3 = varG(F4_public)/varP(F4_public)
> F4_public = selectInd(F4_public, nInd = length(F4_public@id) * 0.2, selectTop = TRUE)
> V4 = varG(F4_public)/varP(F4_public)
> V1;V2;V3;V4
          [,1]
[1,] 0.5204006
          [,1]
[1,] 0.6761066
          [,1]
[1,] 0.8925379
         [,1]
[1,] 1.359076

[gregorgorjanc]

@marcio-resende thanks for this report! I think what you are seeing here is that selection on phenotype is reducing phenotypic variance, but since selection is not with 100% accuracy, you still have plenty of genetic variance, which in fact can be larger than the phenotypic variance in the end. I get the below results when I select on "gv" or on "pheno" - I think plots like the ones shown below can be informative what is happening

> Parents_public = newPop(founderPop[1:100])
> F1_public = randCross(Parents_public, 50)
> 
> use = "gv"
> 
> F2_public = self(F1_public, nProgeny = 25, parents = NULL, keepParents = FALSE)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.3833625
> 
> F2_public = selectWithinFam(F2_public, nInd = 6, selectTop = TRUE, use = use)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.3163724
> 
> F3_public = self(F2_public, nProgeny = 8, parents = NULL, keepParents = FALSE) 
> F3_public = selectWithinFam(F3_public, nInd = 2, selectTop = TRUE, use = use)
> 
> F4_public = self(F3_public, nProgeny = 10, parents = NULL, keepParents = FALSE) 
> V1 = varG(F4_public)/varP(F4_public)
> 
> F4_public = setPheno(F4_public,reps=2)
> V2 = varG(F4_public)/varP(F4_public)
> 
> F4_public = selectWithinFam(F4_public, nInd = 1, selectTop = TRUE, use = use)
> V3 = varG(F4_public)/varP(F4_public)
> 
> F4_public = selectInd(F4_public, nInd = length(F4_public@id) * 0.2, selectTop = TRUE, use = use)
> V4 = varG(F4_public)/varP(F4_public)
> 
> V1;V2;V3;V4
          [,1]
[1,] 0.5841276
          [,1]
[1,] 0.7205692
          [,1]
[1,] 0.6752042
          [,1]
[1,] 0.1888792
> 
> gv = gv(F4_public)
> pheno = pheno(F4_public)
> range = range(c(gv, pheno))
> plot(gv ~ pheno, xlim = range, ylim = range)
> 
> use = "pheno"
> 
> F2_public = self(F1_public, nProgeny = 25, parents = NULL, keepParents = FALSE)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.5584471
> 
> F2_public = selectWithinFam(F2_public, nInd = 6, selectTop = TRUE, use = use)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.7877201
> 
> F3_public = self(F2_public, nProgeny = 8, parents = NULL, keepParents = FALSE) 
> F3_public = selectWithinFam(F3_public, nInd = 2, selectTop = TRUE, use = use)
> 
> F4_public = self(F3_public, nProgeny = 10, parents = NULL, keepParents = FALSE) 
> V1 = varG(F4_public)/varP(F4_public)
> 
> F4_public = setPheno(F4_public,reps=2)
> V2 = varG(F4_public)/varP(F4_public)
> 
> F4_public = selectWithinFam(F4_public, nInd = 1, selectTop = TRUE, use = use)
> V3 = varG(F4_public)/varP(F4_public)
> 
> F4_public = selectInd(F4_public, nInd = length(F4_public@id) * 0.2, selectTop = TRUE, use = use)
> V4 = varG(F4_public)/varP(F4_public)
> 
> V1;V2;V3;V4
          [,1]
[1,] 0.5854095
          [,1]
[1,] 0.7404762
          [,1]
[1,] 0.9117285
         [,1]
[1,] 1.306791
> gv = gv(F4_public)
> pheno = pheno(F4_public)
> range = range(c(gv, pheno))
> plot(gv ~ pheno, xlim = range, ylim = range)

Variation of genetic and phenotype values when selecting on true genetic values

Variation of genetic and phenotype values when selecting on phenotype values

[marcio-resende]

Thank you Gregor,

This is interesting. Would these results imply a negative covariance between gv and the residuals?

[gregorgorjanc]

@marcio-resende yes - selection on a phenotype induces negative correlation between components of the phenotype (genetic value and environmental effect) - this is the Bulmer effect, though Bulmer was thinking only (I think!) on negative correlation between QTL. See this plot on relationship between gv and pheno (left) and gv and residual (pheno - gv) before selection (top row) and after selection (bottom row) - code after the plot

> F2_public = self(F1_public, nProgeny = 25, parents = NULL, keepParents = FALSE)
> varG(F2_public)/varP(F2_public)
          [,1]
[1,] 0.4730515
> 
> par(mfrow = c(2, 2))
> gv = gv(F2_public)
> pheno = pheno(F2_public)
> res = pheno - gv
> range = range(c(gv, pheno))
> range2 = range(c(gv, res))
> plot(gv ~ pheno, xlim = range, ylim = range)
> plot(gv ~ res, xlim = range2, ylim = range2)
> abline(coef(lm(gv ~ res)), col="red")
> cor(gv, res)
           [,1]
[1,] -0.0168661
> 
> F2_public = selectWithinFam(F2_public, nInd = 6, selectTop = TRUE, use = use)
> varG(F2_public)/varP(F2_public)
       [,1]
[1,] 0.8127
> 
> gv = gv(F2_public)
> pheno = pheno(F2_public)
> res = pheno - gv
> plot(gv ~ pheno, xlim = range, ylim = range)
> plot(gv ~ res, xlim = range2, ylim = range2)
> abline(coef(lm(gv ~ res)), col="red")
> cor(gv, res)
          [,1]
[1,] -0.290182