knitr::opts_chunk$set(echo = TRUE)

ANÁLISE DE VARIÂNCIA EM EXPERIMENTOS INTEIRAMENTE AO ACASO

Lembre-se que ao realizar uma análise de variência você está testando duas hipóteses:

A partir dos conceitos teóricos abordados em sala de aula, vamos ver na prática como realizamos com o software R esse tipo de análise. Inicialmente precisamos carregar nossos dados que serão submetidos a análise de variância. Verifique com seu professor como construir a tabela de dados a serem analisados.

library(readxl)
dados <- read_excel("~/Documentos/ESTATÍSTICA INDUSTRIAL/intacaso.xls")
head(dados)
## # A tibble: 6 x 2
##   VARIEDADE PRODUCAO
##   <chr>        <dbl>
## 1 A               25
## 2 A               26
## 3 A               20
## 4 A               23
## 5 A               21
## 6 B               31

Realizando a análise de variância

Precisamos verificar a linguagem R apropriada para a realização desta análise.

analise<-aov(PRODUCAO~VARIEDADE, data = dados)
summary(analise)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## VARIEDADE    3  163.8   54.58   7.798 0.00198 **
## Residuals   16  112.0    7.00                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Precisamos verificar agora o resultado obtido por essa análise. Observamos que o valor do Fcalculado foi de aproximadamente 7,80. O valor de p foi de 0,001, o que verifica-se pela legenda. Mas o que esses resultados nos informam afinal????

Quer dizer que as médias obtidas dependem de cada variedade da planta.

Análise de Variância em experimentos em blocos ao acaso

Como discutido anteriormente, quando tem-se uma variável de bloco, não é o objetivo principal do pesquisador investigá-la, mas ao mesmo tempo desconfia-se que a mesma pode estar influenciando o resultado.

library(readxl)
bloco <- read_excel("~/Documentos/ESTATÍSTICA INDUSTRIAL/bloco.xls", col_types = c("text", "text", "numeric"))
head(bloco)
## # A tibble: 6 x 3
##   VARIEDADE BLOCO PRODUCAO
##   <chr>     <chr>    <dbl>
## 1 A         I           34
## 2 A         II          26
## 3 A         III         33
## 4 A         IV          36
## 5 A         V           31
## 6 B         I           26

Agora que carregamos nossos dados vamos realizar a análise:

analise1<-aov(PRODUCAO~VARIEDADE+BLOCO, data = bloco)
summary(analise1)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## VARIEDADE    3    490   163.3   6.533 0.00722 **
## BLOCO        4    160    40.0   1.600 0.23753   
## Residuals   12    300    25.0                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

E agora como interpretamos os resultados obtidos?????

A lógica é a mesma, porém temos mais uma variável chamada bloco. Observe que somente a variedade apresentou valor significativo nas médias obtidas. A variável bloco apresentou valor p>0,05 e por isso não influencia a resposta. É importante destacar que essa análise é somente válida quando não há fenômeno de interação. Caso isso ocorra, o procedimento de análise deve ser readequado.

Comparação das médias

Agora que já sabemos realizar a análise de variância é necessário saber quem é diferente de quem? No primeiro exemplo em que verificamos que as médias obtidas de cada variedade são diferentes precisamos saber agora quais são diferentes umas das outras.

medias<-TukeyHSD(analise, "VARIEDADE")
medias
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = PRODUCAO ~ VARIEDADE, data = dados)
## 
## $VARIEDADE
##     diff        lwr       upr     p adj
## B-A    4 -0.7874018  8.787402 0.1192178
## C-A    3 -1.7874018  7.787402 0.3123298
## D-A    8  3.2125982 12.787402 0.0010547
## C-B   -1 -5.7874018  3.787402 0.9313122
## D-B    4 -0.7874018  8.787402 0.1192178
## D-C    5  0.2125982  9.787402 0.0391175

E agora??? Ao rodarmos o script gerou uma tabela com todas as comparações possíveis. Mas como saber quem é diferente de quem? Para isso você deve observar o valor indicado na coluna “padj”. Se o valor de cada comparação indicado nessa coluna for inferior a 0,05, quer dizer que as médias em comparação são diferentes.

Outra forma mais visual de verificar isso é através de gráfico. Observe o gráfico abaixo:

plot(medias)

O que você observou no gráfico? Pode-se se perceber que os valores que não interceptam o zero são as comparações que apresentam diferença entre suas médias.

Vamos tentar utilizar algum pacote que torne nossa vida um pouco mais fácil?

Pacote agricolae

library(agricolae)
pacote<-aov(PRODUCAO~VARIEDADE, data = dados)
summary(pacote)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## VARIEDADE    3  163.8   54.58   7.798 0.00198 **
## Residuals   16  112.0    7.00                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Realizamos pelo pacote verifique que o resultado obtido foi exatamente igual. Vamos ao teste de comparação de médias para visualizar.

teste.HSD = with(dados, HSD.test(pacote, 'VARIEDADE', main='PRODUCAO', console=T))
## 
## Study: PRODUCAO
## 
## HSD Test for PRODUCAO 
## 
## Mean Square Error:  7 
## 
## VARIEDADE,  means
## 
##   PRODUCAO      std r Min Max
## A       23 2.549510 5  20  26
## B       27 2.738613 5  24  31
## C       26 2.738613 5  22  29
## D       31 2.549510 5  28  34
## 
## Alpha: 0.05 ; DF Error: 16 
## Critical Value of Studentized Range: 4.046093 
## 
## Minimun Significant Difference: 4.787402 
## 
## Treatments with the same letter are not significantly different.
## 
##   PRODUCAO groups
## D       31      a
## B       27     ab
## C       26      b
## A       23      b
plot(teste.HSD)

Vamos tentar outros gráficos??

bar.group(teste.HSD$groups, ylim=c(0,40))

Melhorando ainda mais o gráfico:

bar.group(teste.HSD$groups, density=10, border="blue", ylim=c(0,40), 
          las=1, angle=45, col='red', main='Teste de Tukey',
          xlab='VARIEDADE', ylab='PRODUÇÃO (ton/ha)')
abline(h=0, col='black', lwd=1.9)