Genomic heritability from a relationship matrix

The restricted log-likelihood profiled over the total variance depends only on h^2, so the estimate is found by evaluating that profile on a deterministic grid and refining the grid maximum by golden-section search. sigma_e^2 = 0 forces h^2 = 1 exactly, which is the degenerate case the tests use.

Usage

Hertbg(y, K, grid = 41, refine = 40)

Arguments

y: Phenotypes.
K: Genomic relationship matrix, n x n.
grid: Number of grid points on (0, 1).
refine: Golden-section refinement steps.

Value

List with estimate, h2, var_g, var_e, total_var, loglik, grid_h2, grid_loglik, n, method.

Details

Formula: h^2 = sigma_g^2 / (sigma_g^2 + sigma_e^2); profile -0.5 [log|V| + log(1'V^{-1}1) + (n - 1) log(y' P y)] with V = h^2 K + (1 - h^2) I.

References

VanRaden (2008), Journal of Dairy Science 91(11):4414-4423. doi:10.3168/jds.2007-0980

Examples

set.seed(1)
X <- matrix(rnorm(20), 10, 2)
K <- X %*% t(X)
y <- rnorm(10)
Hertbg(y, K)
#> $estimate
#> [1] 0.3799601
#> 
#> $h2
#> [1] 0.3799601
#> 
#> $var_g
#> [1] 0.2548769
#> 
#> $var_e
#> [1] 0.4159222
#> 
#> $total_var
#> [1] 0.6707991
#> 
#> $loglik
#> [1] -8.772388
#> 
#> $grid_h2
#>  [1] 0.00000100 0.02500095 0.05000090 0.07500085 0.10000080 0.12500075
#>  [7] 0.15000070 0.17500065 0.20000060 0.22500055 0.25000050 0.27500045
#> [13] 0.30000040 0.32500035 0.35000030 0.37500025 0.40000020 0.42500015
#> [19] 0.45000010 0.47500005 0.50000000 0.52499995 0.54999990 0.57499985
#> [25] 0.59999980 0.62499975 0.64999970 0.67499965 0.69999960 0.72499955
#> [31] 0.74999950 0.77499945 0.79999940 0.82499935 0.84999930 0.87499925
#> [37] 0.89999920 0.92499915 0.94999910 0.97499905 0.99999900
#> 
#> $grid_loglik
#>  [1] -10.630110 -10.176966  -9.850849  -9.605833  -9.416498  -9.267474
#>  [7]  -9.148852  -9.053916  -8.977937  -8.917480  -8.869988  -8.833521
#> [13]  -8.806583  -8.788011  -8.776892  -8.772508  -8.774301  -8.781838
#> [19]  -8.794796  -8.812946  -8.836144  -8.864326  -8.897503  -8.935771
#> [25]  -8.979309  -9.028392  -9.083410  -9.144889  -9.213531  -9.290263
#> [31]  -9.376320  -9.473370  -9.583707  -9.710589  -9.858825 -10.035909
#> [37] -10.254441 -10.538090 -10.940057 -11.630191 -21.754331
#> 
#> $n
#> [1] 10
#> 
#> $method
#> [1] "REML profile in h^2 over V = h^2 K + (1 - h^2) I, golden-section refined; VanRaden (2008) G"
#> 
#> attr(,"class")
#> [1] "morie_rich_result" "list"