In variani/biglmmz: Low Rank Linear Mixed Models (LMMs) Powered by 'bigstatsr'

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

About

library(biglmmz)

Estimating the heritability

The linear mixed model (LMM) is expected to recover the true heritability of 0.8 given the sample size of 1,500 individuals and 200 causal genetic variants.

N <- 1500; M <- 200; h2 <- 0.8

Zg <- sapply(1:M, function(i) rbinom(N, 2, 0.5)) # allele freq. = 0.5

col_means <- colMeans(Zg, na.rm = TRUE)
col_freq <- col_means / 2  # col_means = 2 * col_freq
col_sd <- sqrt(2 * col_freq * (1 - col_freq))

Z <- sweep(Zg, 2, col_means, "-")
Z <- sweep(Z, 2, col_sd , "/")

b <- rnorm(M, 0, sqrt(h2/M))
y <- Z %*% b + rnorm(N, 0, sqrt(1 - h2))

# biglmmg accepts genotypes in the FBM format
G <- as_FBM(Zg)
mod <- biglmmg(y, G = G)
mod$gamma

Calculating the effective sample size (ESS)

# varianace components = s2 * c(h2, 1 - h2)
h2 <- mod$gamma
s2 <- mod$s2

ess(G, h2 = h2, s2 = s2)

variani/biglmmz documentation built on Dec. 15, 2020, 7:58 a.m.