Exercícios de Delineamento Inteiramente Casualizado

Prof. Fernando de Souza Bastos

16 de outubro de 2018

Exercício 4.1

Para comparar a produtividade de quatro variedades de milho, um agrônomo tomou vinte parcelas similares e distribuiu, inteiramente ao acaso, cada uma das 4 variedades em 5 parcelas experimentais. A partir dos dados experimentais fornecidos abaixo, é possível concluir que existe diferença significativa entre as variedades com relação a produtividade, utilizando o nível de significância de 5%?

\[\begin{array}{ccccc} \hline &&&Variedades&&\\ \hline & A & B & C & D \\ \hline & 25 & 31 & 22 &33 \\ & 26 & 25 & 26 &29 \\ & 20 & 28 & 28 &31 \\ & 23 & 27 & 25 &34 \\ & 21 & 24 & 29 &28 \\ \hline Totais& 115 & 135 & 130 &155 \\ \hline Médias& 23 & 27 & 26 &31 \\ \hline \end{array}\] 
library(readr)
if(!is.null(dev.list())) dev.off()
## null device 
##           1
da <- read_delim("https://raw.githubusercontent.com/maf261/maf261.github.io/master/Exercicios/Exercicios_noR/Exer_4_1.txt", 
                 "\t")
da
## # A tibble: 20 x 2
##    Var    Prod
##    <chr> <int>
##  1 A        25
##  2 A        26
##  3 A        20
##  4 A        23
##  5 A        21
##  6 B        31
##  7 B        25
##  8 B        28
##  9 B        27
## 10 B        24
## 11 C        22
## 12 C        26
## 13 C        28
## 14 C        25
## 15 C        29
## 16 D        33
## 17 D        29
## 18 D        31
## 19 D        34
## 20 D        28
str(da)
## Classes 'tbl_df', 'tbl' and 'data.frame':    20 obs. of  2 variables:
##  $ Var : chr  "A" "A" "A" "A" ...
##  $ Prod: int  25 26 20 23 21 31 25 28 27 24 ...
##  - attr(*, "spec")=List of 2
##   ..$ cols   :List of 2
##   .. ..$ Var : list()
##   .. .. ..- attr(*, "class")= chr  "collector_character" "collector"
##   .. ..$ Prod: list()
##   .. .. ..- attr(*, "class")= chr  "collector_integer" "collector"
##   ..$ default: list()
##   .. ..- attr(*, "class")= chr  "collector_guess" "collector"
##   ..- attr(*, "class")= chr "col_spec"
da$Var <- as.factor(da$Var)
is.factor(da$Var)
## [1] TRUE
#Resumo total dos dados
summary(da)
##  Var        Prod      
##  A:5   Min.   :20.00  
##  B:5   1st Qu.:24.75  
##  C:5   Median :26.50  
##  D:5   Mean   :26.75  
##        3rd Qu.:29.00  
##        Max.   :34.00
#Média dos dados por variedade
res <- aggregate(Prod ~ Var, data = da, FUN = mean);
kable(res,format = "markdown",digits = 2, padding = 3,align = "c")
Var Prod
A 23
B 27
C 26
D 31

Sempre é interessante fazer uma análise descritiva dos dados para avaliar linearidade, independencia, normalidade e homocedásticidade.

O pacote lattice de Sarkar (2017) do R é valioso para tais análises:

library(lattice)

xyplot(Prod~Var, data=da, jitter.x=TRUE)

xyplot(Prod~reorder(Var,Prod), data=da, jitter.x=TRUE)

dotplot(Prod~Var, data=da, jitter.x=TRUE)

dotplot(Prod~reorder(Var,Prod), data=da, jitter.x=TRUE)

bwplot(Prod~Var, data=da, pch="|")

bwplot(Prod~reorder(Var,Prod), data=da, pch="|")

#De outra forma

if (!require("pacman")) install.packages("pacman")
pacman::p_load(readr, dplyr, tibble, ggplot2, ggthemes, car, agricolae)
ggplot(da, aes(Var, Prod)) +
  geom_boxplot() +
  theme_bw() +
  theme_few()

O gráfico dos resíduos versus valores ajustados (valores preditos) é uma das principais técnicas utilizadas para verificar as suposições dos resíduos. Além da detecção de heteroscedasticidade, esse gráfico pode indicar que não existe uma relação linear entra as variáveis explicativas com a variável resposta por meio de alguma tendência nos pontos. Por exemplo, se os pontos do gráfico formam uma parábola, é indicativo que termos de segundo grau sejam necessários.

Para o diagnóstico de heteroscedasticidade, tentamos encontrar alguma tendência no gráfico. Por isso, se os pontos estão aleatoriamente distribuídos em torno do 0, sem nenhum comportamento ou tendência, temos indícios de que a variância dos resíduos é homoscedástica. Já a presença de “funil” é um indicativo da presença de heteroscedasticidade.

Para obtenção dos resíduos padronizados fazemos uma mudança de escala dos resíduos, tend uma linha horizontal temos que a variabilidade é constante

Pontos ao redor da linha significam normalidade.

Mesma avaliação do plot 1, mais utilizado quando utilizamos regressão linear

m0 <- lm(Prod~Var, data=da)
par(mfrow=c(2,2))
plot(m0)

layout(1)

\[ \begin{align} H_{0}&:\mu_{1}=\mu_{2}=\mu_{3}=\mu_{4}\\ H_{1}&: \textrm{Não}\ H_{0}\ (\textrm{Pelo menos uma média difere das demais}) \end{align} \]

m0 <- aov(Prod~Var, data=da)
anova(m0)
## Analysis of Variance Table
## 
## Response: Prod
##           Df Sum Sq Mean Sq F value   Pr(>F)   
## Var        3 163.75  54.583  7.7976 0.001976 **
## Residuals 16 112.00   7.000                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
tukey <- HSD.test(m0, "Var",alpha=0.05)
tukey
## $statistics
##   MSerror Df  Mean       CV      MSD
##         7 16 26.75 9.890659 4.787402
## 
## $parameters
##    test name.t ntr StudentizedRange alpha
##   Tukey    Var   4         4.046093  0.05
## 
## $means
##   Prod      std r Min Max Q25 Q50 Q75
## A   23 2.549510 5  20  26  21  23  25
## B   27 2.738613 5  24  31  25  27  28
## C   26 2.738613 5  22  29  25  26  28
## D   31 2.549510 5  28  34  29  31  33
## 
## $comparison
## NULL
## 
## $groups
##   Prod groups
## D   31      a
## B   27     ab
## C   26      b
## A   23      b
## 
## attr(,"class")
## [1] "group"

Para deixar mais visual ainda, podemos construir um gráfico de barras com a média de cada tratamento e adicionar a sua letra correspondente ao teste de Tukey:

tukey$groups %>% 
  rownames_to_column(var = "trt") %>% 
  ggplot(aes(reorder(trt, Prod, function(x) -mean(x)), Prod)) +
    geom_bar(stat = "identity") +
    geom_text(aes(label = groups), vjust = 1.8, size = 9, color = "white") +
    labs(x = "Var", y = "Médias") +
    theme_few()

Duncan<- duncan.test(m0, "Var", alpha=0.05)
Duncan
## $statistics
##   MSerror Df  Mean       CV
##         7 16 26.75 9.890659
## 
## $parameters
##     test name.t ntr alpha
##   Duncan    Var   4  0.05
## 
## $duncan
##      Table CriticalRange
## 2 2.997999      3.547280
## 3 3.143802      3.719797
## 4 3.234945      3.827638
## 
## $means
##   Prod      std r Min Max Q25 Q50 Q75
## A   23 2.549510 5  20  26  21  23  25
## B   27 2.738613 5  24  31  25  27  28
## C   26 2.738613 5  22  29  25  26  28
## D   31 2.549510 5  28  34  29  31  33
## 
## $comparison
## NULL
## 
## $groups
##   Prod groups
## D   31      a
## B   27      b
## C   26     bc
## A   23      c
## 
## attr(,"class")
## [1] "group"

Para deixar mais visual ainda, podemos construir um gráfico de barras com a média de cada tratamento e adicionar a sua letra correspondente ao teste de Tukey:

Duncan$groups %>% 
  rownames_to_column(var = "trt") %>% 
  ggplot(aes(reorder(trt, Prod, function(x) -mean(x)), Prod)) +
    geom_bar(stat = "identity") +
    geom_text(aes(label = groups), vjust = 1.8, size = 9, color = "white") +
    labs(x = "Var", y = "Médias") +
    theme_few()

library(ExpDes.pt)
dic <- dic(da$Var,
  da$Prod,
  quali = TRUE,
  mcomp = "tukey",
  nl = FALSE,
  hvar = 'bartlett',
  sigT = 0.05,
  sigF = 0.05)
## ------------------------------------------------------------------------
## Quadro da analise de variancia
## ------------------------------------------------------------------------
##            GL     SQ     QM     Fc     Pr>Fc
## Tratamento  3 163.75 54.583 7.7976 0.0019756
## Residuo    16 112.00  7.000                 
## Total      19 275.75                        
## ------------------------------------------------------------------------
## CV = 9.89 %
## 
## ------------------------------------------------------------------------
## Teste de normalidade dos residuos 
## Valor-p:  0.2358736 
## De acordo com o teste de Shapiro-Wilk a 5% de significancia, os residuos podem ser considerados normais.
## ------------------------------------------------------------------------
## 
## ------------------------------------------------------------------------
## Teste de homogeneidade de variancia 
## valor-p:  0.9981235 
## De acordo com o teste de bartlett a 5% de significancia, as variancias podem ser consideradas homogeneas.
## ------------------------------------------------------------------------
## 
## Teste de Tukey
## ------------------------------------------------------------------------
## Grupos Tratamentos Medias
## a     D   31 
## ab    B   27 
##  b    C   26 
##  b    A   23 
## ------------------------------------------------------------------------
#Gráficos de Resíduos do Pacote ExpDes.pt
plotres(dic)

#Teste Duncan

dic <- dic(da$Var,
  da$Prod,
  quali = TRUE,
  mcomp = "duncan",
  nl = FALSE,
  hvar = 'bartlett',
  sigT = 0.05,
  sigF = 0.05)
## ------------------------------------------------------------------------
## Quadro da analise de variancia
## ------------------------------------------------------------------------
##            GL     SQ     QM     Fc     Pr>Fc
## Tratamento  3 163.75 54.583 7.7976 0.0019756
## Residuo    16 112.00  7.000                 
## Total      19 275.75                        
## ------------------------------------------------------------------------
## CV = 9.89 %
## 
## ------------------------------------------------------------------------
## Teste de normalidade dos residuos 
## Valor-p:  0.2358736 
## De acordo com o teste de Shapiro-Wilk a 5% de significancia, os residuos podem ser considerados normais.
## ------------------------------------------------------------------------
## 
## ------------------------------------------------------------------------
## Teste de homogeneidade de variancia 
## valor-p:  0.9981235 
## De acordo com o teste de bartlett a 5% de significancia, as variancias podem ser consideradas homogeneas.
## ------------------------------------------------------------------------
## 
## Teste de Duncan 
## ------------------------------------------------------------------------
## Grupos  Tratamentos  Medias
## a     D       31 
##  b    B       27 
##  bc   C       26 
##   c   A       23 
## ------------------------------------------------------------------------

Exercício 4.2

Um treinador de corrida rústica, objetivando melhorar o desempenho de seus atletas, testou três novas técnicas de preparação. Para tanto trabalhou com um grupo de 15 atletas completamente homogêneos para as características essenciais. A designação das técnicas de preparação aos atletas foi feita totalmente ao acaso e de tal forma que o número de atletas avaliados em cada uma das técnicas fosse o mesmo. Os resultados obtidos, após um determinado período de tempo de aprendizado da técnica pelos atletas, foram os seguintes (minutos / 25 Km):

\[\begin{array}{ccccc} \hline &&\textrm{Técnicas de Preparação}&&\\ \hline \textrm{Repetições} & 1 & 2 & 3 \\ \hline 1 & 130 & 125 & 135 \\ 2 & 129 & 131 & 129 \\ 3 & 128 & 130 & 131 \\ 4 & 126 & 129 & 128 \\ 5 & 130 & 127 & 130 \\ \hline Totais& 643 & 642 & 653 \\ \hline \textrm{Médias}& 128.6 & 128.4 & 130.6 \\ \hline \end{array}\] 
library(readr)
da <- read_delim("https://raw.githubusercontent.com/maf261/maf261.github.io/master/Exercicios/Exercicios_noR/Exer_4_2.txt", 
                 "\t")
da
## # A tibble: 15 x 2
##    Tec   Tempo
##    <chr> <int>
##  1 X1      130
##  2 X1      129
##  3 X1      128
##  4 X1      126
##  5 X1      130
##  6 X2      125
##  7 X2      131
##  8 X2      130
##  9 X2      129
## 10 X2      127
## 11 X3      135
## 12 X3      129
## 13 X3      131
## 14 X3      128
## 15 X3      130
str(da)
## Classes 'tbl_df', 'tbl' and 'data.frame':    15 obs. of  2 variables:
##  $ Tec  : chr  "X1" "X1" "X1" "X1" ...
##  $ Tempo: int  130 129 128 126 130 125 131 130 129 127 ...
##  - attr(*, "spec")=List of 2
##   ..$ cols   :List of 2
##   .. ..$ Tec  : list()
##   .. .. ..- attr(*, "class")= chr  "collector_character" "collector"
##   .. ..$ Tempo: list()
##   .. .. ..- attr(*, "class")= chr  "collector_integer" "collector"
##   ..$ default: list()
##   .. ..- attr(*, "class")= chr  "collector_guess" "collector"
##   ..- attr(*, "class")= chr "col_spec"
da$Tec <- as.factor(da$Tec)
is.factor(da$Tec)
## [1] TRUE
#Resumo total dos dados
summary(da)
##  Tec        Tempo      
##  X1:5   Min.   :125.0  
##  X2:5   1st Qu.:128.0  
##  X3:5   Median :129.0  
##         Mean   :129.2  
##         3rd Qu.:130.0  
##         Max.   :135.0
#Média dos dados por variedade
res <- aggregate(Tempo ~ Tec, data = da, FUN = mean);
kable(res,format = "markdown",digits = 2, padding = 3,align = "c")
Tec Tempo
X1 128.6
X2 128.4
X3 130.6

Usando o pacote lattice de Sarkar (2017):

library(lattice)

xyplot(Tempo~Tec, data=da, jitter.x=TRUE)

xyplot(Tempo~reorder(Tec,Tempo), data=da, jitter.x=TRUE)

dotplot(Tempo~Tec, data=da, jitter.x=TRUE)

dotplot(Tempo~reorder(Tec,Tempo), data=da, jitter.x=TRUE)

bwplot(Tempo~Tec, data=da, pch="|")

bwplot(Tempo~reorder(Tec,Tempo), data=da, pch="|")

m0 <- lm(Tempo~Tec, data=da)
par(mfrow=c(2,2))
plot(m0)

layout(1)

\[ \begin{align} H_{0}&:\mu_{1}=\mu_{2}=\mu_{3}\\ H_{1}&: \textrm{Não}\ H_{0}\ (\textrm{Pelo menos uma média difere das demais}) \end{align} \] - Usando a função interna do R

m0 <- lm(Tempo~Tec, data=da)
anova(m0)
## Analysis of Variance Table
## 
## Response: Tempo
##           Df Sum Sq Mean Sq F value Pr(>F)
## Tec        2   14.8     7.4  1.3962  0.285
## Residuals 12   63.6     5.3
library(ExpDes.pt)
dic <- dic(da$Tec,
  da$Tempo,
  quali = TRUE,
  mcomp = "tukey",
  nl = FALSE,
  hvar = 'bartlett',
  sigT = 0.01,
  sigF = 0.01)
## ------------------------------------------------------------------------
## Quadro da analise de variancia
## ------------------------------------------------------------------------
##            GL   SQ  QM     Fc Pr>Fc
## Tratamento  2 14.8 7.4 1.3962 0.285
## Residuo    12 63.6 5.3             
## Total      14 78.4                 
## ------------------------------------------------------------------------
## CV = 1.78 %
## 
## ------------------------------------------------------------------------
## Teste de normalidade dos residuos 
## Valor-p:  0.9262885 
## De acordo com o teste de Shapiro-Wilk a 5% de significancia, os residuos podem ser considerados normais.
## ------------------------------------------------------------------------
## 
## ------------------------------------------------------------------------
## Teste de homogeneidade de variancia 
## valor-p:  0.663661 
## De acordo com o teste de bartlett a 5% de significancia, as variancias podem ser consideradas homogeneas.
## ------------------------------------------------------------------------
## 
## De acordo com o teste F, as medias nao podem ser consideradas diferentes.
## ------------------------------------------------------------------------
##   Niveis Medias
## 1     X1  128.6
## 2     X2  128.4
## 3     X3  130.6
## ------------------------------------------------------------------------
plotres(dic)

Exercício 4.4

da = data.frame(Trat=rep(c("A","B","C"),each=8),
                Resp=c(96,95,100,108,120,110.5,97  ,92.5,
                  90,93,89 ,88 ,87 ,92.5 ,87.5,85  ,
                  86,85,105,105,90 ,100  ,95  ,95))
head(da)
##   Trat  Resp
## 1    A  96.0
## 2    A  95.0
## 3    A 100.0
## 4    A 108.0
## 5    A 120.0
## 6    A 110.5
str(da)
## 'data.frame':    24 obs. of  2 variables:
##  $ Trat: Factor w/ 3 levels "A","B","C": 1 1 1 1 1 1 1 1 2 2 ...
##  $ Resp: num  96 95 100 108 120 ...
m0 <- lm(Resp~Trat,data=da)
par(mfrow=c(2,2))
plot(m0)

layout(1)

Escolha um número que seja próximo, o ideal é que seja um número do tipo \(1/a,\) onde \(a \in \mathbb{Z}.\)

MASS::boxcox(m0)

MASS::boxcox(lm(Resp~Trat,data=da),lambda=seq(-10,1,by=.1))

MASS::boxcox(lm(Resp~Trat,data=da),lambda=seq(-10,1,by=.1))
abline(v=-4,col="red")
abline(v=-4.25,col="blue")

m1 <- lm(Resp^(-4)~Trat,data = da)
plot(m1)

library(ExpDes.pt)
dic <- dic(da$Trat,
  da$Resp^(-4),
  quali = TRUE,
  mcomp = "tukey",
  nl = FALSE,
  hvar = 'bartlett',
  sigT = 0.05,
  sigF = 0.05)
## ------------------------------------------------------------------------
## Quadro da analise de variancia
## ------------------------------------------------------------------------
##            GL         SQ         QM   Fc     Pr>Fc
## Tratamento  2 1.6047e-16 8.0234e-17 7.64 0.0032122
## Residuo    21 2.2054e-16 1.0502e-17               
## Total      23 3.8101e-16                          
## ------------------------------------------------------------------------
## CV = 25.08 %
## 
## ------------------------------------------------------------------------
## Teste de normalidade dos residuos 
## Valor-p:  0.6995647 
## De acordo com o teste de Shapiro-Wilk a 5% de significancia, os residuos podem ser considerados normais.
## ------------------------------------------------------------------------
## 
## ------------------------------------------------------------------------
## Teste de homogeneidade de variancia 
## valor-p:  0.1587246 
## De acordo com o teste de bartlett a 5% de significancia, as variancias podem ser consideradas homogeneas.
## ------------------------------------------------------------------------
## 
## Teste de Tukey
## ------------------------------------------------------------------------
## Grupos Tratamentos Medias
## a     B   1.606925e-08 
## ab    C   1.296106e-08 
##  b    A   9.735807e-09 
## ------------------------------------------------------------------------
#Avalie os resultados
plotres(dic)

Referências

Ferreira, Eric Batista, Portya Piscitelli Cavalcanti, and Denismar Alves Nogueira. 2018. ExpDes.pt: Pacote Experimental Designs (Portuguese). https://CRAN.R-project.org/package=ExpDes.pt.

Sarkar, Deepayan. 2017. “Lattice: Trellis Graphics for R.” março. https://cran.r-project.org/web/packages/lattice/index.html.