#' Generacion de datos simulados
#'
#' @description generate data
#'
#' @param n.fixed.effects number of fixed effects
#' @param n.random.effects n random effects
#' @param n.sample size of sample
#' @importFrom fastDummies dummy_columns
#' @importFrom stats rnorm
#' @importFrom MASS mvrnorm
#' @importFrom magic adiag
#' @export
data.generation <- function(n.fixed.effects, n.random.effects, n.sample) {
"
generacion de datos artificiales con los que probaremos la estimacion del modelo
jerarquico
"
n <- n.sample
p <- n.fixed.effects
m <- n.random.effects
sigma2 <- 25
tau00_2 <- 20
tau11_2 <- 10
# simulacion construccion de la data de entrada
group <- sort(sample(1:m, n, replace = T))
# agregar intercepto
dummy <- as.matrix(dummy_columns(group)[, -1])
# variable longitudinal
l <- matrix(rnorm(n = n, mean = 0, sd = 1), nrow = n, ncol = 1)
d <- diag(n)
diag(d) <- l
U <- combine_matrix(dummy, d %*% dummy)
# punto de validacion de la matriz
stopifnot(sum(sum(U)) == n + sum(l))
# covariables aleatorias para testear
A <- matrix(rnorm(n = n * p, mean = 10, sd = 5), n, p)
A <- d %*% A
# agregar intercepto
X <- cbind(matrix(rep(1, n), ncol = 1), A)
# matriz Q que representa la covarianza entre el factor aleatorio del intercepto con la pendiente, en
# este caso, independientes.
Q <- matrix(nrow = 2, ncol = 2)
Q[1, 1] <- tau00_2
Q[2, 2] <- tau11_2
Q[2, 1] <- 0
Q[1, 2] <- 0
# matriz de varianza covarianza de los factores aleatorios
G <- do.call(adiag, rep(list(Q), ncol(U)/2))
# generar los factores aleatorios
eta <- matrix(mvrnorm(mu = rep(0, ncol(U)), Sigma = G), nrow = ncol(U), ncol = 1)
# generar los factores fijos, para testear generamos solo unos
beta <- matrix(rep(1, p + 1), nrow = p + 1, ncol = 1)
# error
eps <- matrix(rnorm(n = n, mean = 0, sd = sqrt(sigma2)), ncol = 1)
# generar una variable de respuesta a partir de los datos
Y <- X %*% beta + U %*% eta + eps
list(X = X, U = U, Y = Y, beta = beta, eta = eta, tau00_2 = tau00_2, tau11_2 = tau11_2, sigma2 = sigma2,
group = group)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.