R/estimate_ranef.R

estimate_ranef <- function(xz_svd, y, beta, prior_cov, family,
        alpha) {
    ranef_all_groups <- Map(estimate_ranef_single_group, xz_svd, y,
                            list(beta), list(prior_cov), list(family),
                            list(alpha))
    ranef_all_groups <- Reduce(rbind, ranef_all_groups)
    tibble::as_tibble(ranef_all_groups)
}

estimate_ranef_single_group <- function(xz_svd, y, beta, prior_cov, family,
		alpha) {
	U <- xz_svd$U
	D <- xz_svd$D
	V1 <- xz_svd$V1
	V2 <- xz_svd$V2

	UD <- U %*% D
	x <- UD %*% t(V1)
	z <- UD %*% t(V2)
    gamma <- chol(solve(prior_cov))
    z <- z %*% solve(gamma)

	precision <- (1 - alpha) * diag(ncol(z))
	offset <- x %*% beta

	fit <- fit_lbfgs(z, y, precision, family, offset = offset, alpha = alpha)
	solve(gamma, fit$par)
}
kschmaus/gammmbest documentation built on May 7, 2019, 9 p.m.