[R-es] Comparaciones múltiples en ANOVA anidadp
Luciano Selzer
luciano.selzer en gmail.com
Dom Feb 5 04:14:06 CET 2012
Hola José, puedo estár equivocado porque son las 12 de la noche ...
Pero creo que el modelo equivalente en lmer es
lmer(VR ~ Tratamiento + (1|Individuo), data = Medidas)
El otro modelo implica que tus individuos son únicos para cada tratamiento.
Un saludo
Luciano
El día 4 de febrero de 2012 16:45, José Trujillo Carmona
<trujillo en unex.es> escribió:
> Explico un poco más el problema con lmer:
>
> Si utilizo lmer, todo parece funcionar, pero el análisis de la varianza no
> es consistente con el de aov:
>
>> AnovaModel.4 <- lmer(VR ~ Tratamiento + (1|Tratamiento/Individuo),
>> data=Medidas)
>
>
>> Pares <- glht(AnovaModel.3, linfct = mcp(Tratamiento = "Tukey"))
>> summary(Pares)
>
> Simultaneous Tests for General Linear Hypotheses
>
> Multiple Comparisons of Means: Tukey Contrasts
>
>
> Fit: lmer(formula = VR ~ Tratamiento + (1 | Tratamiento/Individuo), data =
> Medidas)
>
> Linear Hypotheses:
> Estimate Std. Error z value Pr(>|z|)
> 2 - 1 == 0 -0.14933 0.93547 -0.160 1.000
> 3 - 1 == 0 -0.13078 0.93547 -0.140 1.000
> 4 - 1 == 0 0.24593 0.93547 0.263 0.999
> 5 - 1 == 0 -1.08025 0.93547 -1.155 0.777
> 3 - 2 == 0 0.01855 0.93547 0.020 1.000
> 4 - 2 == 0 0.39526 0.93547 0.423 0.993
> 5 - 2 == 0 -0.93091 0.93547 -0.995 0.858
> 4 - 3 == 0 0.37671 0.93547 0.403 0.994
> 5 - 3 == 0 -0.94947 0.93547 -1.015 0.849
> 5 - 4 == 0 -1.32618 0.93547 -1.418 0.616
> (Adjusted p values reported -- single-step method)
>
>> anova(AnovaModel.4)
> Analysis of Variance Table
>
> Df Sum Sq Mean Sq F value
> Tratamiento 4 2.4993 0.62482 0.5819
>
>
>> AnovaModel.2 <- aov(VR ~ Tratamiento + Error(Individuo),data=Medidas)
>
>> summary(AnovaModel.2)
>
> Error: Individuo
> Df Sum Sq Mean Sq F value Pr(>F)
> Tratamiento 4 9.166 2.2915 3.473 0.0502 .
> Residuals 10 6.599 0.6599
> ---
> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
>
> Error: Within
> Df Sum Sq Mean Sq F value Pr(>F)
> Residuals 30 36.35 1.212
>
> No comprendo como habría de escribir la función en lmer para que fuese
> equivalente al de aov que creo congruente con los textos:
>
> Bioestaticas Analysis de Zar
> Bimetry de Sokal y Rohlf
> Bioestadística de Steel y Torrie
>
>
> Reitero mi agradecimiento por aguantar
>
>
>
> El 04/02/12 20:05, José Trujillo Carmona escribió:
>
>> Dispongo de un experimento en el que cinco tratamientos ha sido aplicados
>> a cinco grupos de voluntarios.
>>
>> En cada grupo había tres personas y a cada persona se le tomaron 3
>> medidas,
>>
>> En total dispongo de 45 medidas, pero evidentemente no son independientes
>> entre sí. Si no tomo en cuenta que las medidas de la misma persona son más
>> parecidas entre sí (bloques anidados) estaría incurriendo en
>> pseudoreplicación.
>>
>> Los datos y el análisis se podrían hacer de la siguiente forma:
>>
>> Generación de datos para ejemplo:
>> > Medidas <- as.data.frame(matrix(rnorm(45*1, mean=10, sd=1), ncol=1))
>> > colnames(Medidas) <- "VR"
>> > Medidas$Tratamiento <- as.factor(rep(1:5,rep(9,5)))
>> > Medidas$Individuo <- as.factor(rep(1:15,rep(3,15)))
>>
>> El análisis de la varianza sugerido podría ser:
>>
>> > AnovaModel.1 <-aov(VR ~ Tratamiento + Tratamiento/Individuo,
>> > data=Medidas)
>> > summary(AnovaModel.1)
>> Df Sum Sq Mean Sq F value Pr(>F)
>> Tratamiento 4 1.15 0.2871 0.25 0.907
>> Tratamiento:Individuo 10 12.39 1.2395 1.08 0.407
>> Residuals 30 34.44 1.1479
>>
>>
>> Pero este análisis incurre en evidente pseudoreplicación: El valor F es el
>> cociente entre el Cuadrado Medio del Tratamiento y los Residuos como si las
>> 45 mediciones fuesen independientes entre sí.
>>
>> Para tener en cuenta que la variabilidad inducida en la respuesta por los
>> Tratamientos debe ser medida entre individuos y no entre mediciones el
>> análisis procedente podría ser:
>>
>> > AnovaModel.2 <- aov(VR ~ Tratamiento + Error(Individuo),data=Medidas)
>> > summary(AnovaModel.2)
>>
>> Error: Individuo
>> Df Sum Sq Mean Sq F value Pr(>F)
>> Tratamiento 4 1.148 0.2871 0.232 0.914
>> Residuals 10 12.395 1.2395
>>
>> Error: Within
>> Df Sum Sq Mean Sq F value Pr(>F)
>> Residuals 30 34.44 1.148
>>
>> Pero el objeto AnovaModel.2 no es un objeto aov sino aovlist:
>> > attr(AnovaModel.2,"class")
>> [1] "aovlist" "listof"
>>
>> Como consecuencia no puedo utilizarlo ni para la función glht ni para la
>> función TukeyHSD.
>>
>> Puedo "extraer" un aov:
>>
>> > ls.str(pat="AnovaModel.2")
>> AnovaModel.2 : List of 3
>> $ (Intercept):List of 9
>> $ Individuo :List of 9
>> $ Within :List of 6
>> > AnovaModel.3<-AnovaModel.2$Individuo
>> > attr(AnovaModel.3,"class")
>> [1] "aov" "lm"
>>
>> A este modelo ya puedo pedirle límites de confianza adecuados para los
>> parámetros del modelo (diferencias entre Tratamientos):
>>
>> > confint(AnovaModel.3)
>> 2.5 % 97.5 %
>> Tratamiento[T.2] -1.485182 0.8535804
>> Tratamiento[T.3] -1.072891 1.2658714
>> Tratamiento[T.4] -1.258890 1.0798723
>> Tratamiento[T.5] -1.453857 0.8849054
>>
>> Pero el nuevo modelo sigue siendo inválido para glht o para TuleyHSD:
>>
>> > Pares <- glht(AnovaModel.3, linfct = mcp(Tratamiento = "Tukey"))
>> [1] ERROR: no 'model.matrix' method for 'model' found!
>>
>> > TukeyHSD(AnovaModel.3, "Tratamiento")
>> [2] ERROR: this fit does not inherit from "lm"
>>
>> Para no alargarme en mis indagaciones. Para los análisis de la varianza
>> podría haber utilizado las funciones lme del paquete nlme o lmer del paquete
>> lme4, pero tampoco proporcionan objetos válidos para glht o TukeyHSD.
>>
>> ¿Alguien sabe como abordar las comparaciones múltiples con medidas
>> repetidas o anova anidado?
>>
>> Gracias de antemano.
>>
>
> --
> _____---^---_____
>
> Univ. de Extremadura
> Dept. Matemáticas.
> Despacho B29
> Tf: + 34 924 289 300
> Ext. 86823
>
> _______________________________________________
> R-help-es mailing list
> R-help-es en r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-help-es
Más información sobre la lista de distribución R-help-es