library(MASS)
library(ggplot2)
library(gridExtra)

Tamano de muestra y variables independientes

La estructura de la población y la edad de la menarquia se simulan de la misma forma que en el ejemplo anterior.

Variables dependientes

Estatura

A diferencia de la primera simulación, en ésta si se toma en cuenta la variación materna, es decir, este modelo toma en cuenta que existen diferencias entre familias en la relación entre la menarquia y la estatura (i.e. una \(\alpha\) y \(\beta\) diferente por familia). Utilizando los mismos parámetros que anteriormente (Onland-Mored et al., 2005), el modelo del cual se simula la estatura es el siguiente:

\[estatura_{ij} = \alpha_j + \beta_j*menarquia_{ij} + e_{ij}\\ \begin{bmatrix} \alpha_j\\ \beta_j \end{bmatrix} \sim \mathcal{MVN}\left( \begin{bmatrix} \alpha_{0}\\ \beta_{0} \end{bmatrix}, S\right)\]

# Parametros
alfa.0 <- 161.5  # cm
beta.0 <- 0.31
sigma_e <- 0.1
e <- rnorm(sum(n_off),0,sigma_e)  # error aleatorio
rho <- matrix(c(1, -0.3, -0.3, 1), nrow = 2)
s_alpha <- 0.2
s_beta <- 0.4 
Sigmas <- c(s_alpha, s_beta)
factor_var <- 0.1 # representa la variabilidad ambiental
S <- factor_var * diag(Sigmas) %*% rho %*% diag(Sigmas)

# Simulacion de los parametros para cada madre
alfa_beta <- mvrnorm(N_mom, mu=c(alfa.0,beta.0),Sigma=S) 

# Simulacion de la estatura (cm) para las hijas de cada madre
estatura <- vector(mode="numeric",length = sum(n_off))
estatura[1:n_off[1]] <- alfa_beta[1,1] +
  alfa_beta[1,2]*menarquia[1:n_off[1]] +
  e[1:n_off[1]]  # error aleatorio

cum_n_off <- cumsum(n_off)
for(j in 2:N_mom){
  estatura[(cum_n_off[j-1]+1):cum_n_off[j]] <- alfa_beta[j,1] +
    alfa_beta[j,2]*menarquia[(cum_n_off[j-1]+1):cum_n_off[j]] +
    e[(cum_n_off[j-1]+1):cum_n_off[j]]  # error aleatorio
}

# Base de datos de las variables creadas
d <- data.frame(estatura = estatura,
                menarquia = menarquia,
                mother = as.factor(mother))

Visualización