El documento describe un experimento con dos factores fijos y aleatorios. Se analizan diferentes modelos lineales mixtos para el conjunto de datos Machines, incluido un modelo con interacciones aleatorias entre máquina y trabajador. Los resultados muestran que este último modelo es el más apropiado según las pruebas de significancia y los gráficos de verificación.
0 calificaciones0% encontró este documento útil (0 votos)
14 vistas1 página
El documento describe un experimento con dos factores fijos y aleatorios. Se analizan diferentes modelos lineales mixtos para el conjunto de datos Machines, incluido un modelo con interacciones aleatorias entre máquina y trabajador. Los resultados muestran que este último modelo es el más apropiado según las pruebas de significancia y los gráficos de verificación.
El documento describe un experimento con dos factores fijos y aleatorios. Se analizan diferentes modelos lineales mixtos para el conjunto de datos Machines, incluido un modelo con interacciones aleatorias entre máquina y trabajador. Los resultados muestran que este último modelo es el más apropiado según las pruebas de significancia y los gráficos de verificación.
El documento describe un experimento con dos factores fijos y aleatorios. Se analizan diferentes modelos lineales mixtos para el conjunto de datos Machines, incluido un modelo con interacciones aleatorias entre máquina y trabajador. Los resultados muestran que este último modelo es el más apropiado según las pruebas de significancia y los gráficos de verificación.
Descargue como PDF, TXT o lea en línea desde Scribd
Descargar como pdf o txt
Está en la página 1/ 1
Práctica 4
Maximiliano Vaca Montejano
364897 maximiliano.vaca@uabc.edu.mx
Profesor: Andrés García Medina
andres.garcia.medina@uabc.edu.mx Fecha de entrega: jueves 11 de abril, 2024 (12pm).
Ejercicio 1 (20 pts)
Resuelva los siguientes problemas conceptuales:
:D
Ejercicio 2 (40 pts)
Analice los datos Machines del ejemplo guia 3 (ec. 4) a través de la función lme
del ejemplo guía:
Consideremos ahora un experimento en el que cada observación se agrupa según dos factores.En la figura se muestra un diagrama esquemático de dicho diseño.
No description has been provided for this image
Supongamos que un factor se va a modelar como un efecto fijo y el otro como un efecto aleatorio. Entonces, un modelo para la k-ésima observación en el nivel i de efecto fijo A y el nivel j de efecto aleatorio B es
yijk = μ + αi + bj + (ab)ij + eijk
(ec. 4)
donde
2 bj ∼ N (0, α ) b
2 (αb)ij ∼ N (0, σ ) αθ
2 ϵijk ∼ N (0, σ )
y todas las v.a. son mutuamente independientes.
Ademas, μ es la media poblacional global, αi son los I efectos fijos para el factor A , y bj representan los J efectos aleatorios para el factor B.
Por otro lado, (α ∗ b)ij son las interacciones.
(a) Intente encontrar el modelo más apropiado, teniendo cuidado de examinar los gráficos de verificación del modelo apropiado
In [3]: library(nlme) data(Machines)
head(Machines, 5)
A nffGroupedData: 5 × 3 Worker Machine score
<ord> <fct> <dbl>
1 1 A 52.0
2 1 A 52.8
3 1 A 53.1
4 2 A 51.8
5 2 A 52.8
verificamos que hay interacciones
In [4]: attach( Machines ) # make variables in Machines available by name
interaction.plot( Machine, Worker, score, las = 1) # Figure 1.10 detach() # undo the effect of ‘attach( Machines )’
partimos del modelo mas simple, sin interacciones
con un efecto fijo para cada maquina y un efecto aleatorio para cada trabajador
fm1Machine lme.formula(fixed = score ~ Machine, data = Machines, random = ~1 | Worker) 1 5 296.8782 306.5373 -143.4391 NA NA
fm2Machine lme.formula(fixed = score ~ Machine, data = Machines, random = ~1 | Worker/Machine) 2 6 227.6876 239.2785 -107.8438 1 vs 2 71.19063 3.24324e-17
nuestro modelo fm2, el de la ec4, es apropiado en comparacion con el mas simple
(b) Asegúrese de probar si la interacción de la ec. 4 es apropiada.
In [9]: m1 <- lm(score ~ Worker*Machine, Machines)
In [14]: plot(m1, which = 1)
si, la interaccion de la ec. 4 es apropiada y es a la que habiamos llegado desde lme con fm2
(c) De manera similar, pruebe si sería apropiada una estructura de efectos aleatorios más compleja: específicamente una en la que la interacción máquina-trabajador esté correlacionada con la del trabajador.
In [16]: fm3Machine <- update(fm1Machine, random = ~ Machine | Worker)
fm2Machine lme.formula(fixed = score ~ Machine, data = Machines, random = ~1 | Worker/Machine) 1 6 227.6876 239.2785 -107.8438 NA NA
fm3Machine lme.formula(fixed = score ~ Machine, data = Machines, random = ~Machine | Worker) 2 10 228.3112 247.6295 -104.1556 1 vs 2 7.37635 0.117287
el modelo fm3, donde la interacción máquina-trabajador está correlacionada con la del trabajador parece ser apropiada en comparacion con el modelo mas simple y se muestra bastante competente respecto al segundo modelo fm2, aunque este se siga prefiriendo.
(d) Si algún dato parece particularmente problemático en los gráficos de verificación, repita el análisis y vea si las conclusiones cambian.
In [20]: #nuevo conjunto de datos que excluye los puntos 18, 17 y 26
Machines_new <- Machines[-c(18, 17, 26), ]
Machines_new
In [22]: attach( Machines_new ) # make variables in Machines available by name
interaction.plot( Machine, Worker, score, las = 1) # Figure 1.10 detach() # undo the effect of ‘attach( Machines )’
yijk = βj + bi + ϵijk
In [23]: fm1Machine_new <- lme( score ~ Machine, data = Machines_new, random = ~ 1 | Worker )
el modelo mejoró al deshacernos del ruido en las graficas
no me parece que lleve a un sobreajuste ya que eliminamos pocos datos
Nos aseguramos nuevamente de probar si la interacción de la ec. 4 es apropiada.
In [30]: m1_new <- lm(score ~ Worker*Machine, Machines_new)
In [32]: plot(m1_new, which = 1)
las mismas conclusiones
de manera similar, pruebe si sería apropiada una estructura de efectos aleatorios más compleja: específicamente una en la que la interacción máquina-trabajador esté correlacionada con la del trabajador.
In [33]: fm3Machine_new <- update(fm1Machine_new, random = ~ Machine | Worker)
fm2Machine_new lme.formula(fixed = score ~ Machine, data = Machines_new, random = ~1 | Worker/Machine) 1 6 196.6804 207.9076 -92.34019 NA NA
fm3Machine_new lme.formula(fixed = score ~ Machine, data = Machines_new, random = ~Machine | Worker) 2 10 196.3857 215.0977 -88.19285 1 vs 2 8.294684 0.08136082
mismas conclusiones
Ejercicio 3 (40 pts)
Repitamos el ejemplo guiado de los rieles usando un enfoque de máxima verosimilitud. En este caso: Σθ Σ_φ = Iσ^2, porloquelosparámetros(θ,φ) = (log σ,log σ_b). Lasiguientef uncióntomacomoentrada :(θ,φ),X,Z,y 2 = Iσ , b
(a)Explique cada una de las lineas del código de la función con base en la teoría revisada en el curso
In [39]: llm <- function(parameters, X, Z, y) {
sigma.b <- exp(parameters[1]) #regresa el 1er parámetro a la escala original (desde logaritmo) sigma <- exp(parameters[2]) #regresa el 2do parámetro a la escala original (desde logaritmo)
n <- length(y) # # de observaciones
pr <- ncol(Z) # # de efectos aleatorios pf <- ncol(X) # # de efectos fijos
X1 <- cbind(X, Z) #combina las matrices de efectos fijos y aleatorios
ipsi <- c(rep(0, pf), rep(1/sigma.b^2, pr)) #vector de inversas de las varianzas de los efectos fijos y aleatorios
b1 <- solve(crossprod(X1)/sigma^2 + diag(ipsi), t(X1) %*% y/sigma^2) #estimación de los coeficientes de los efectos fijos y aleatorios
ldet <- sum(log(diag(chol(crossprod(Z)/sigma^2 + diag(ipsi[-(1:pf)]))))) #logaritmo del determinante de la matriz de covarianza de los efectos aleatorios
attr(l, "b") <- as.numeric(b1) #añade los coeficientes estimados como un atributo de la log-verosimilitud
return(-l) #devuelve la log-verosimilitud negativa
}
(b) Determine X,Z ,y de manera explicita.
lo dejé para hoy y ya no alcanzo (11:40)
(c) Utilice la función optim para optimizar la función de verosimilitud (llm) fijando como valores iniciales: parameters = c(0, 0)
In [40]: parameters <- c(0, 0)
res_optim <- optim(parameters, llm, X = X, Z = Z, y = y)
Error in fn(par, ...): object 'y' not found
Traceback:
1. optim(parameters, llm, X = X, Z = Z, y = y) 2. (function (par) . fn(par, ...))(c(0, 0)) 3. fn(par, ...)
(d) Determine los valores ́ optimos de (σ2b , σ2 ), y compare con los obtenidos en el ejemplo guía correspondiente. ¿A que se deben las posible diferencias numéricas?
In [41]: sigma.b.opt <- exp(res_optim$par[1])
sigma.opt <- exp(res_optim$par[2])
Error in eval(expr, envir, enclos): object 'res_optim' not found
