standard error for heritability estimates

[kidd...@gmail.com]

Hello Karl,

I was wondering if it possible to work out standard error for heritability from est_herit or is there some other way to estimate how confident we can be about the heritability estimate.

Thanks

Senthil

[Karl Broman]

I’ve not looked at this. Seems like you could use equation 2 from Visscher and Goddard (2015),  which just depends on eigenvalues of kinship matrix, sample size, and heritability. Or use a “parametric bootstrap” (simulate multivariate normal data with the estimated covariance matrix). 

karl

[Karl Broman]

I took a stab at it here: https://github.com/rqtl/qtl2/issues/193#issuecomment-790711138

# estimated SE of estimated heritability
# Visscher & Goddard (2015), https://doi.org/10.1534/genetics.114.171017
#
# lambda_i = eigenvalues of kinship matrix
# a = sum[ (lam_i-1)^2 / (1 + h^2(lam_i - 1))^2 ]
# b = sum[ (lam_i - 1) / (1 + h^2(lam_i - 1)) ]

# N = sample size
#

# var(h_hat^2) = 2 / (a - b^2/N)

se_herit <-
    function(pheno, kinship)
{

    h2 <- as.numeric(est_herit(pheno, kinship))

    d <- decomp_kinship(kinship)

    lam <- d$values
    a <- sapply(h2, function(hsq) sum( (lam-1)^2 / (1 + hsq*(lam-1))^2 ) )
    b <- sapply(h2, function(hsq) sum( (lam-1) / (1 + hsq*(lam-1)) ) )
    n <- length(lam)

    setNames( sqrt(2 / (a - b^2/n)), colnames(pheno) )
}

[kidd...@gmail.com]

Hello Karl,

Thanks a lot for the quick answer and putting together code.

Will run it on some data and let you know how it goes.

Thanks

Senthil