Diseño Factorial Completo en arreglo Completamente al Azar (DFCCA)
set.seed(123)
# Respuesta: Volumen de raices
vol = c(rnorm(40, 6, 1),
rnorm(20, 8, 1.3))
# Factor 1: variedad de papa (2 niveles)
variedad = gl(2, 30, 60, c('Colombia', 'Pastusa'))
# Factor 2: ferilizante foliar (3 niveles)
ferti = gl(3, 10, 60, c('0cc','5cc','10cc'))
df = data.frame(vol, variedad, ferti)
head(df)
## vol variedad ferti
## 1 5.439524 Colombia 0cc
## 2 5.769823 Colombia 0cc
## 3 7.558708 Colombia 0cc
## 4 6.070508 Colombia 0cc
## 5 6.129288 Colombia 0cc
## 6 7.715065 Colombia 0cc
medias = tapply(vol, list(variedad, ferti), mean); medias
## 0cc 5cc 10cc
## Colombia 6.074626 6.208622 5.575441
## Pastusa 6.322045 7.988670 8.288192
desv = tapply(vol, list(variedad, ferti), sd); desv
## 0cc 5cc 10cc
## Colombia 0.9537841 1.038073 0.9308092
## Pastusa 0.5273024 1.407274 1.1133786
cv = 100*desv/medias; cv
## 0cc 5cc 10cc
## Colombia 15.701117 16.71987 16.69481
## Pastusa 8.340694 17.61587 13.43331
library(lattice)
bwplot(vol~ ferti | variedad)
xyplot(vol~ ferti | variedad)
Modelo:
Hipotesis: \[H_O:(tau*beta)_ij=0\]
Nota: Las dos siguientes solo son interpretables si no ahy interaccion
mod1 = aov(vol ~ variedad*ferti, df)
summary(mod1)
## Df Sum Sq Mean Sq F value Pr(>F)
## variedad 1 37.45 37.45 35.375 2.07e-07 ***
## ferti 2 9.18 4.59 4.334 0.01797 *
## variedad:ferti 2 15.49 7.75 7.318 0.00154 **
## Residuals 54 57.17 1.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ifelse(mod1$p.value<0.05,"Hay interaccion","No hay interaccion")
## logical(0)
Como el p-value de la interaccion es menor a 0.05, entonces se rechaza la hipotesis de No interaccion
interaction.plot(df$ferti,
df$variedad,
df$vol)
df_mean = aggregate(list(mvol = df$vol),
list(variedad = df$variedad,
ferti = df$ferti),
mean)
library(ggplot2)
ggplot(df_mean)+
aes(ferti, mvol, color = variedad)+
geom_line(aes(group = variedad), size = 1)+
geom_point(size = 2)
ggplot(df_mean)+
aes(variedad, mvol, color = ferti)+
geom_line(aes(group = ferti), size = 1)+
geom_point(size = 2)
¿Es esta interacción remediable?
df2 = df[df$ferti != '0cc',]
head(df2)
## vol variedad ferti
## 11 7.224082 Colombia 5cc
## 12 6.359814 Colombia 5cc
## 13 6.400771 Colombia 5cc
## 14 6.110683 Colombia 5cc
## 15 5.444159 Colombia 5cc
## 16 7.786913 Colombia 5cc
mod2 = aov(vol ~ variedad*ferti, df2)
summary(mod2)
## Df Sum Sq Mean Sq F value Pr(>F)
## variedad 1 50.46 50.46 39.088 3.2e-07 ***
## ferti 1 0.28 0.28 0.216 0.645
## variedad:ferti 1 2.17 2.17 1.685 0.203
## Residuals 36 46.48 1.29
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
df_mean2 = aggregate(list(mvol = df2$vol),
list(variedad = df2$variedad,
ferti = df2$ferti),
mean)
ggplot(df_mean2)+
aes(ferti, mvol, color = variedad)+
geom_line(aes(group = variedad), size = 1)+
geom_point(size = 2)
Al no tener interacción se puede hacer comparación de medias
TukeyHSD(mod2, 'variedad')
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = vol ~ variedad * ferti, data = df2)
##
## $variedad
## diff lwr upr p adj
## Pastusa-Colombia 2.246399 1.517693 2.975106 3e-07
res = residuals(mod2)
trt = interaction(df2$variedad, df2$ferti)
shapiro.test(res)
##
## Shapiro-Wilk normality test
##
## data: res
## W = 0.97349, p-value = 0.4608
bartlett.test(df2$vol, trt)
##
## Bartlett test of homogeneity of variances
##
## data: df2$vol and trt
## Bartlett's K-squared = 1.6537, df = 3, p-value = 0.6473
df_mean2
## variedad ferti mvol
## 1 Colombia 5cc 6.208622
## 2 Pastusa 5cc 7.988670
## 3 Colombia 10cc 5.575441
## 4 Pastusa 10cc 8.288192
# aumento porcentual del rendimiento entre dosis fertilizante para variedad Pastusa
100 * (df_mean2$mvol[4]-df_mean2$mvol[2])/df_mean2$mvol[2]
## [1] 3.749335
Factorial Completo en Bloques Completos al Azar
set.seed(123)
# Respuesta: Volumen de raices
vol = c(rnorm(20, 6, 1),
rnorm(10, 8, 1.3))
# Factor 1: variedad de papa (2 niveles)
variedad = gl(2, 15, 30, c('Colombia', 'Pastusa'))
# Factor 2: ferilizante foliar (3 niveles)
ferti = gl(3, 5, 30, c('0cc','5cc','10cc'))
# Bloque: Nivel de alguna caracteristica del suelo
bloq = gl(5, 1, 30, c('--','-','+-','+','++'))
df = data.frame(vol, variedad, ferti, bloq)
head(df)
## vol variedad ferti bloq
## 1 5.439524 Colombia 0cc --
## 2 5.769823 Colombia 0cc -
## 3 7.558708 Colombia 0cc +-
## 4 6.070508 Colombia 0cc +
## 5 6.129288 Colombia 0cc ++
## 6 7.715065 Colombia 5cc --
library(collapsibleTree)
collapsibleTree(df, c('bloq', 'variedad','ferti'))
table(df$bloq, df$ferti, df$variedad)
## , , = Colombia
##
##
## 0cc 5cc 10cc
## -- 1 1 1
## - 1 1 1
## +- 1 1 1
## + 1 1 1
## ++ 1 1 1
##
## , , = Pastusa
##
##
## 0cc 5cc 10cc
## -- 1 1 1
## - 1 1 1
## +- 1 1 1
## + 1 1 1
## ++ 1 1 1
ggplot(df)+
aes(ferti, vol, color=variedad)+
geom_point(size = 2)+
geom_line(aes(group=variedad))+
facet_wrap(.~bloq, ncol = 5)
mod3 = aov(vol ~ bloq + variedad*ferti, df)
summary(mod3)
## Df Sum Sq Mean Sq F value Pr(>F)
## bloq 4 2.256 0.564 0.418 0.7934
## variedad 1 5.412 5.412 4.015 0.0588 .
## ferti 2 4.384 2.192 1.626 0.2216
## variedad:ferti 2 3.522 1.761 1.307 0.2929
## Residuals 20 26.958 1.348
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(mod3, 'variedad')
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = vol ~ bloq + variedad * ferti, data = df)
##
## $variedad
## diff lwr upr p adj
## Pastusa-Colombia 0.8494454 -0.03487277 1.733764 0.0588366
med_variedad = tapply(df$vol, df$variedad, mean)
med_variedad
## Colombia Pastusa
## 6.152384 7.001830
ggplot(df)+
aes(vol, fill = variedad)+
geom_density(alpha = 0.3)