O modelo para estes dados seria
\[y_{ij}=\mu+\tau_i+\varepsilon_{i,j},\; i=1,...,5\] assim temos que pela análise de variância temos como hipótese:
valores <- c(7,7,15,11,9,
12,17,12,18,18,
14,19,19,18,18,
19,25,22,19,23,
7,10,11,15,11)
trat<- factor(rep(1:5, each = 5))
modelo <- aov(valores~trat)
# data.frame(valores, trat)
summary(modelo)
## Df Sum Sq Mean Sq F value Pr(>F)
## trat 4 475.8 118.94 14.76 9.13e-06 ***
## Residuals 20 161.2 8.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Logo, pelo teste, temos evidência com 5% de significância, que existe algum \(\tau_i\) diferente (média do tratamento).
Os costrantes são do tipo que agrega os valores em subgrupos, no qual é possivel com 5 contranstes obter 5-1 constrastes:
Para construir contrastes ortogonais, temos que construido 4 constrantes, com o intuito, de baixo, médio, alto sobre o Cotton weight percent:
\[\begin{array}{c} A=\sum a_i \bar{y_i} \Rightarrow \sum a_i=0\\ B=\sum b_i \bar{y_i} \Rightarrow \sum b_i=0\\ C=\sum c_i \bar{y_i} \Rightarrow \sum c_i=0\\ D=\sum d_i \bar{y_i} \Rightarrow \sum d_i=0\\ \sum a_i.c_i = 0\\ \sum b_i.c_i = 0\\ \sum b_i.d_i = 0\\ \sum c_i.d_i = 0 \end{array}\]
Ou seja, temos como hipótese:
Onde se analisarmos A, \(\sum a_i \bar{y_i}\sim N(\sum a_i \mu_i, \sum a_i^2\sigma^2/n)\)
no qual os quadrados médios são dados por exemplo no contranste \(A\),
\[E(A)=E(\sum a_i \bar{y}_{i.})=\sum a_i E(\bar{y}_{i.})=\sum a_i \mu_i\] \[Var(A)=Var(\sum a_i \bar{y}_{i.})=\sum Var( a_i \bar{y}_{i.})=\sum a_i^2Var( \bar{y}_{i.})=\sum a_i^2\sigma^2/n\] Assim como estatistica de teste temos como hipotese de não haver efeito na comparação dos grupos ortogonais, ou seja, \(\sum a_i \bar{y_i}=0\)
\[{\sum a_i \bar{y_i} - 0\over \sqrt{\sum a_i^2\sigma^2/n}}\sim N(0, 1)\] No qual, para o estimador de \(\hat\sigma^2=QMRes\), temos que
\[{\sum a_i \bar{y_i} - 0\over \sqrt{\sum a_i^2QMRes/n}}\sim t(an-a)\] ou ainda
\[{(\sum a_i \bar{y_i} - 0)^2\over \sum a_i^2QMRes/n}\sim F(1,an-a)\]
medias <- apply(matrix(valores, ncol = 5, byrow=T), 1, mean) # y_barra
a = 5
n = 5
A = c(0, 1, 1, -1, -1)
B = c(0, 1, -1, 0, 0)
C = c(0, 0, 0, 1, -1)
D = c(4, -1, -1, -1, -1)
sqA <- sum(A*medias)^2/(sum(A^2)/n)
sqB <- sum(B*medias)^2/(sum(B^2)/n)
sqC <- sum(C*medias)^2/(sum(C^2)/n)
sqD <- sum(D*medias)^2/(sum(D^2)/n)
QMres <- summary(modelo)[[1]][2,3]
sqA
## [1] 0.45
sqB
## [1] 12.1
sqC
## [1] 291.6
sqD
## [1] 171.61
Logo, realizando os testes, temos que, como é uma comparação 2 a 2, onde possui normalidade e não se sabe a variância, efetuamos o teste T, mas ao elevar ao quadrado temos que a estatistica do teste tem distribuição F.
round(1-pf(c(sqA,sqB, sqC, sqD)/QMres, 1, a*n - a),4)
## [1] 0.8156 0.2347 0.0000 0.0002
Assim, pelo teste F, com 95% de confiança, rejeitamos as hipóteses dos contrastes, B e C, e não rejeitamos A e C.
sqA+sqB+sqC+sqD
## [1] 475.76
summary(modelo)[[1]][1,2]
## [1] 475.76
dist(medias)
## 1 2 3 4
## 2 5.6
## 3 7.8 2.2
## 4 11.8 6.2 4.0
## 5 1.0 4.6 6.8 10.8
O teste de Turkey, temos com hipótese o nível exato de significância pensando no global (onde no teste de fisher, isso não acontece, havendo um nível de significâcia menor que o adotado), para a comparação 2 a 2 das médias dos tratamentos.
# Valor critico, valores maiores que isso são rejeitados com nivel
# de significancia igual a 5%
qtukey(.95, nmeans = 5, df = a*n - a)
## [1] 4.231857
# P-valor
round((1-ptukey(dist(medias)/sqrt(QMres/n),nmeans = a, df = a*n-a)),4)
## 1 2 3 4
## 2 0.0385
## 3 0.0026 0.7372
## 4 0.0000 0.0189 0.2101
## 5 0.9798 0.1163 0.0091 0.0001
TukeyHSD(modelo)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = valores ~ trat)
##
## $trat
## diff lwr upr p adj
## 2-1 5.6 0.2270417 10.9729583 0.0385024
## 3-1 7.8 2.4270417 13.1729583 0.0025948
## 4-1 11.8 6.4270417 17.1729583 0.0000190
## 5-1 1.0 -4.3729583 6.3729583 0.9797709
## 3-2 2.2 -3.1729583 7.5729583 0.7372438
## 4-2 6.2 0.8270417 11.5729583 0.0188936
## 5-2 -4.6 -9.9729583 0.7729583 0.1162970
## 4-3 4.0 -1.3729583 9.3729583 0.2101089
## 5-3 -6.8 -12.1729583 -1.4270417 0.0090646
## 5-4 -10.8 -16.1729583 -5.4270417 0.0000624
plot(TukeyHSD(modelo))
O teste de fisher faz a comparação multipla dos valores 2 a 2.
# Valor critico, valores maiores que isso são rejeitados com nivel
# de significancia igual a 5%, mas há um incremento no nível de significacia
# global
qt(.975, df = a*n - a)*sqrt(QMres*2/n)
## [1] 3.745452
# P-valor
round(2*(1-pt(dist(medias)/sqrt(QMres*2/n), df = a*n-a)),4)
## 1 2 3 4
## 2 0.0054
## 3 0.0003 0.2347
## 4 0.0000 0.0025 0.0375
## 5 0.5838 0.0186 0.0012 0.0000
Como o modelo segue \(y_{i,j}=\mu+\tau_i+\varepsilon_{i,j}\), a variável aleatória é \(\varepsilon\sim N(0,\sigma^2)\), logo isolando esta variável, podemos verificar se os resíduos possui normalidade, ou seja, através de uma exploração empírica de um gráfico ou um teste de normalidade como shapiro-wilk:
# approximadamente 0
sum(modelo$residuals)
## [1] -2.664535e-15
hist(modelo$residuals, breaks = 10)
shapiro.test(modelo$residuals)
##
## Shapiro-Wilk normality test
##
## data: modelo$residuals
## W = 0.94387, p-value = 0.1818
Assim, com o teste de shapiro-wilk, não rejeitamos a normalidade dos dados, ao nível de siginificâcia de 5%.