Caso 1
dados <- c(4.2,4.5,5.2,3.8,3.7,3.5,4.2,4,3.9,
4.6,4.7,5,4.4,3.5,3.1,4.2,3.8,3.9,
4.5,4.3,6.8,4.8,3.1,3.4,5.2,3.7,3.7,
4.4,4.7,5.8,3.9,3.7,3.3,5.1,4.1,4)
cult<- factor(rep(letters[1:3], 3*4))
prepar<- factor(rep(1:3, each = 3, 4))
bloco<- factor(rep(1:4, each = 9))
data.frame(prepar, dados)
## prepar dados
## 1 1 4.2
## 2 1 4.5
## 3 1 5.2
## 4 2 3.8
## 5 2 3.7
## 6 2 3.5
## 7 3 4.2
## 8 3 4.0
## 9 3 3.9
## 10 1 4.6
## 11 1 4.7
## 12 1 5.0
## 13 2 4.4
## 14 2 3.5
## 15 2 3.1
## 16 3 4.2
## 17 3 3.8
## 18 3 3.9
## 19 1 4.5
## 20 1 4.3
## 21 1 6.8
## 22 2 4.8
## 23 2 3.1
## 24 2 3.4
## 25 3 5.2
## 26 3 3.7
## 27 3 3.7
## 28 1 4.4
## 29 1 4.7
## 30 1 5.8
## 31 2 3.9
## 32 2 3.7
## 33 2 3.3
## 34 3 5.1
## 35 3 4.1
## 36 3 4.0
# modelo
mod <- aov(dados ~ prepar*cult+ Error(bloco/prepar))
summary(mod)
##
## Error: bloco
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 3 0.5297 0.1766
##
## Error: bloco:prepar
## Df Sum Sq Mean Sq F value Pr(>F)
## prepar 2 8.912 4.456 69.24 7.16e-05 ***
## Residuals 6 0.386 0.064
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## cult 2 1.322 0.6608 3.636 0.047181 *
## prepar:cult 4 6.107 1.5267 8.399 0.000525 ***
## Residuals 18 3.272 0.1818
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Soma de quadrados
medparcela <- tapply(dados, paste(bloco, prepar), mean)
sqparcela <- 3*sum((medparcela - mean(dados))^2)
medbloco <- tapply(dados, bloco, mean)
sqbloco <- 9*sum((medbloco - mean(dados))^2)
medprepar <- tapply(dados, prepar, mean)
sqprepar <- 12*sum((medprepar - mean(dados))^2)
sqresa <- sqparcela - sqbloco -sqprepar
sqresa
## [1] 0.3861111
sqbloco
## [1] 0.5297222
sqprepar
## [1] 8.911667
summary(mod)
##
## Error: bloco
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 3 0.5297 0.1766
##
## Error: bloco:prepar
## Df Sum Sq Mean Sq F value Pr(>F)
## prepar 2 8.912 4.456 69.24 7.16e-05 ***
## Residuals 6 0.386 0.064
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## cult 2 1.322 0.6608 3.636 0.047181 *
## prepar:cult 4 6.107 1.5267 8.399 0.000525 ***
## Residuals 18 3.272 0.1818
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
medcult <- tapply(dados, cult, mean)
sqcult <- 12*sum((medcult - mean(dados))^2)
medint <- tapply(dados, paste(prepar, cult), mean)
sqint <- 4*sum((medint - mean(dados))^2)- sqcult - sqprepar
sqtot <- sum((dados- mean(dados))^2)
sqresb <- sqtot- sqparcela - sqint- sqcult
sqcult
## [1] 1.321667
sqint
## [1] 6.106667
sqresb
## [1] 3.271667
summary(mod)
##
## Error: bloco
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 3 0.5297 0.1766
##
## Error: bloco:prepar
## Df Sum Sq Mean Sq F value Pr(>F)
## prepar 2 8.912 4.456 69.24 7.16e-05 ***
## Residuals 6 0.386 0.064
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## cult 2 1.322 0.6608 3.636 0.047181 *
## prepar:cult 4 6.107 1.5267 8.399 0.000525 ***
## Residuals 18 3.272 0.1818
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
df <- data.frame(dados, cult, bloco, prepar)
# Prepar -> cult -> valores
# Within, resb é suficiente para fazer os testes (teste para as sub parcela)
interaction.plot(prepar, cult, dados)
## testes se tem diferença sig. em cada prepar
sqprep1 = summary(aov(dados~cult ,data = df[df$prepar==1,]))[[1]][1,2]
sqprep2 = summary(aov(dados~cult ,data = df[df$prepar==2,]))[[1]][1,2]
sqprep3 = summary(aov(dados~cult ,data = df[df$prepar==3,]))[[1]][1,2]
sqprep1+sqprep2+sqprep3
## [1] 7.428333
sqint +sqcult
## [1] 7.428333
## Teste se existe alguma diferença dentro de cada prepar
1-pf(c(sqprep1,sqprep2,sqprep3)/(sqresb/18),2,18)
## [1] 1.580133e-05 1.187933e-03 1.853098e-03
## Testes 2 a 2 para cada nivel de prepar
# 4 é quantidade de blocos por preparo/cult
qtukey(.95, 3, 18)*sqrt((sqresb/18)/4)
## [1] 0.7693814
### nivel de prepar = 1
so <- df[df$prepar== 1, ]
dist(tapply(so$dados, so$cult, mean))
## a b
## b 0.125
## c 1.275 1.150
# cult -> prepar -> valores
# terar que criar um residuo composto
interaction.plot(cult, prepar, dados)
## Testes 2 a 2 para o cultivo, se os preparos sao diferentes
k = 3 # quantidade de sub parcelas
sqcult1 = summary(aov(dados~prepar ,data = df[df$cult=="a",]))[[1]][1,2]
sqcult2 = summary(aov(dados~prepar ,data = df[df$cult=="b",]))[[1]][1,2]
sqcult3 = summary(aov(dados~prepar ,data = df[df$cult=="c",]))[[1]][1,2]
sqcult1+sqcult2+sqcult3
## [1] 15.01833
sqint + sqprepar
## [1] 15.01833
qmresa <- sqresa/6
qmresb <- sqresb/18
qmrescomp <- 1/3*(qmresa+ 2*qmresb)
glcomp <- (qmresa+ 2*qmresb)^2/(qmresa^2/6+ (2*qmresb)^2/18)
# Testes
# Teste se para cada nivel de cult tem diferenca no prep
1-pf(c(sqcult1,sqcult2,sqcult3)/qmrescomp, 3, round(glcomp))
## [1] 5.954688e-02 8.741906e-06 1.102674e-12
## Testes 2 a 2 para cada nivel de cult
# 4 é quantidade de blocos por preparo/cult
# Tukey
qtukey(.95, 3, round(glcomp))*sqrt(qmrescomp/4)
## [1] 0.6687644
# fisher
qt(0.975, round(glcomp))*sqrt(2*qmrescomp/4)
## [1] 0.5524196
so <- df[df$cult== "c", ]
dist(tapply(so$dados, so$prepar, mean))
## 1 2
## 2 2.375
## 3 1.825 0.550