devtools::load_all()
library(gammmbest)
library(mbest)

set.seed(1)

n_groups <- 1000
n_obs <- 4

laplace_ids <- sample(c(-1, 1), n_groups, replace = TRUE) * rexp(n_groups) / sqrt(2)
gaussian_ids <- rnorm(n_groups)

group <- rep(seq_len(n_groups), each = n_obs)
noise <- rnorm(n_obs * n_groups)
y <- laplace_ids[group] + noise

mbest_model <- mhglm(y ~ 1 + (1 | group), data = data.frame(group = group))

gammmbest_model_alpha0 <- gammmbest_fit(x = matrix(1, nrow = n_obs * n_groups), 
    y = y, group = group, family = gaussian(), fixef_index = 1, ranef_index = 1,
    cov_supp = matrix(1), alpha = 0)

gammmbest_model_alpha1 <- gammmbest_fit(x = matrix(1, nrow = n_obs * n_groups), 
    y = y, group = group, family = gaussian(), fixef_index = 1, ranef_index = 1,
    cov_supp = matrix(1), alpha = 1)

gammmbest_model_alpha0.5 <- gammmbest_fit(x = matrix(1, nrow = n_obs * n_groups), 
    y = y, group = group, family = gaussian(), fixef_index = 1, ranef_index = 1,
    cov_supp = matrix(1), alpha = 0.5)



kschmaus/gammmbest documentation built on May 7, 2019, 9 p.m.