set.seed(123)
# Variable: Pesos
peso = c(runif(40, 130, 150),
runif(40, 150, 180),
runif(40, 80, 130))
# Variable: Compuestos Fenolicos
cf = c(rnorm(40, 20, 1),
rnorm(40, 12, 1.2),
rnorm(40, 16, 1.15))
# Factor 1: Tiempo de Coccion
tiemp = gl(3, 40, 120, c(0, 10, 12))
# Factor 2: Vierdades
varie = gl(2, 20, 120, c('var1', 'var2'))
df = data.frame(varie, tiemp, cf, peso)
head(df)
## varie tiemp cf peso
## 1 var1 0 20.37964 135.7516
## 2 var1 0 19.49768 145.7661
## 3 var1 0 19.66679 138.1795
## 4 var1 0 18.98142 147.6603
## 5 var1 0 18.92821 148.8093
## 6 var1 0 20.30353 130.9111
View(df)
boxplot(df$cf ~ df$varie)

boxplot(df$cf ~ df$tiemp)
lattice::bwplot(df$cf ~ df$tiemp | df$varie)


library(collapsibleTree)
collapsibleTreeSummary(df,
hierarchy = c('varie', 'tiemp'),
collapsed = FALSE)
mod1 = aov(cf ~ tiemp*varie, df)
summary(mod1)
## Df Sum Sq Mean Sq F value Pr(>F)
## tiemp 2 1399.9 700.0 602.599 <2e-16 ***
## varie 1 5.4 5.4 4.639 0.0334 *
## tiemp:varie 2 0.1 0.0 0.033 0.9678
## Residuals 114 132.4 1.2
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
tapply(df$cf, df$varie, mean)
## var1 var2
## 15.74983 16.17367
tapply(df$cf, df$tiemp, mean)
## 0 10 12
## 20.12759 11.76136 15.99629
modt = aov(cf+peso ~ tiemp*varie, df)
summary(modt)
## Df Sum Sq Mean Sq F value Pr(>F)
## tiemp 2 59372 29686 271.832 <2e-16 ***
## varie 1 271 271 2.483 0.118
## tiemp:varie 2 67 34 0.307 0.736
## Residuals 114 12450 109
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(mod1, 'tiemp')
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = cf ~ tiemp * varie, data = df)
##
## $tiemp
## diff lwr upr p adj
## 10-0 -8.366227 -8.938526 -7.793928 0
## 12-0 -4.131296 -4.703594 -3.558997 0
## 12-10 4.234931 3.662632 4.807230 0
TukeyHSD(modt, 'tiemp')
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = cf + peso ~ tiemp * varie, data = df)
##
## $tiemp
## diff lwr upr p adj
## 10-0 13.53163 7.982541 19.08072 2e-07
## 12-0 -38.94101 -44.490101 -33.39192 0e+00
## 12-10 -52.47264 -58.021731 -46.92355 0e+00
res1 = mod1$residuals
# Prueba de normalidad en los residuales
shapiro.test(res1)
##
## Shapiro-Wilk normality test
##
## data: res1
## W = 0.98475, p-value = 0.1949
trat = interaction(df$tiemp, df$varie)
bartlett.test(res1, trat)
##
## Bartlett test of homogeneity of variances
##
## data: res1 and trat
## Bartlett's K-squared = 7.1214, df = 5, p-value = 0.2118
plot(res1)

plot(mod1$fitted.values, res1)

rest = modt$residuals
# Prueba de normalidad en los residuales
shapiro.test(rest)
##
## Shapiro-Wilk normality test
##
## data: rest
## W = 0.99105, p-value = 0.6311
bartlett.test(rest, trat)
##
## Bartlett test of homogeneity of variances
##
## data: rest and trat
## Bartlett's K-squared = 32.322, df = 5, p-value = 5.129e-06
plot(rest)

plot(mod1$fitted.values, rest)

library(ggplot2)
mdf = aggregate(list(cf = df$cf),
list(varie = df$varie,
tiemp = df$tiemp),
mean)
ggplot(mdf)+
aes(x = tiemp, y = cf, group = varie, color = varie)+
# geom_boxplot()+
geom_line()

ggplot(mdf)+
aes(x = varie, y = cf, group = tiemp, color = tiemp)+
# geom_boxplot()+
geom_line(linetype = 'dashed')

mod1b = aov(cf ~ tiemp + varie, df)
summary(mod1b)
## Df Sum Sq Mean Sq F value Pr(>F)
## tiemp 2 1399.9 700.0 612.819 <2e-16 ***
## varie 1 5.4 5.4 4.718 0.0319 *
## Residuals 116 132.5 1.1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
shapiro.test(mod1b$residuals)
##
## Shapiro-Wilk normality test
##
## data: mod1b$residuals
## W = 0.98569, p-value = 0.2366
bartlett.test(mod1b$residuals, trat)
##
## Bartlett test of homogeneity of variances
##
## data: mod1b$residuals and trat
## Bartlett's K-squared = 7.1214, df = 5, p-value = 0.2118
df$bloq = gl(2, 10, 120, c('biol','alim'))
head(df)
## varie tiemp cf peso bloq
## 1 var1 0 20.37964 135.7516 biol
## 2 var1 0 19.49768 145.7661 biol
## 3 var1 0 19.66679 138.1795 biol
## 4 var1 0 18.98142 147.6603 biol
## 5 var1 0 18.92821 148.8093 biol
## 6 var1 0 20.30353 130.9111 biol
collapsibleTreeSummary(df,
hierarchy = c('bloq','varie', 'tiemp'),
collapsed = FALSE)
mod2 = aov(cf ~ bloq * tiemp * varie, df)
summary(mod2)
## Df Sum Sq Mean Sq F value Pr(>F)
## bloq 1 0.1 0.1 0.066 0.7981
## tiemp 2 1399.9 700.0 589.774 <2e-16 ***
## varie 1 5.4 5.4 4.541 0.0354 *
## bloq:tiemp 2 0.8 0.4 0.351 0.7044
## bloq:varie 1 0.1 0.1 0.055 0.8148
## tiemp:varie 2 0.1 0.0 0.032 0.9684
## bloq:tiemp:varie 2 3.3 1.6 1.375 0.2573
## Residuals 108 128.2 1.2
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
mod2b = aov(cf ~ bloq + tiemp * varie, df)
summary(mod2b)
## Df Sum Sq Mean Sq F value Pr(>F)
## bloq 1 0.1 0.1 0.067 0.7968
## tiemp 2 1399.9 700.0 597.666 <2e-16 ***
## varie 1 5.4 5.4 4.601 0.0341 *
## tiemp:varie 2 0.1 0.0 0.033 0.9680
## Residuals 113 132.3 1.2
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
lattice::bwplot(cf ~ tiemp | varie + bloq, df)
