Resumo

Nesse trabalho estuda-se a análise de variância, conhecida como ANOVA. Análise de variância (Anova) é uma técnica utilizada com a finalidade de comparar os efeitos de diferentes tratamentos. Para que tenhamos um conteúdo adequado sobre o assunto precisamos falar de:

Para essa análise, precisamos certificar que os erros são variáveis aleatórias independentes com distribuição normal de média zero e variância constante, o qual é chamado depressuposto de normalidade.

Esse método estátistico trata-se de investigar igualdade de médias em experimentos propostos. Usado também na regressão estátistica, o método é usado na aplicação de modelos estatístico onde particionamos a variância amostral em fatores. Na ANOVA aplicamos também o teste de hipótese, onde testamos a hipótese das igualdades das médias contra a dierença das mesmas. Na aplicação do método em modelos estatísticos usamos alguns testes alternativos na ANOVA que são testes de normalidade e teste de homoscedasticidade. A ANOVA tem uma forte importância na tomada de decisão, ajuda a verificar diferenças estatísticas na hora da comparação. Em resumo, a ANOVA trata-se de um método estatístico que permite realizar comparações simultâneas entre duas ou mais médias, ou seja, permite testar hipóteses sobre médias de distintas populações, tarta de um teste sobre a igualdade (ou não) dos valores esperados (médias) de uma determinada variável de interesse nas k subpopulações de interesse.

1 Introdução

Análise da Variância (ANOVA) é um método para testar a igualdade de duas ou mais médias populacionais, baseado na análise das variâncias amostrais. Os dados amostrais são separados em grupos seguindo uma característica (fator). Fator (ou tratamento): é uma característica que permite distinguir diferentes populações umas das outras. Cada fator contém dois ou mais grupos (classificações).

Podemos também descrever a ANOVA como uma coleção de modelos estatísticos no qual a variância amostral é particionada em diversos componentes devido a diferentes fatores (variáveis), que nas aplicações estão associados a um processo, produto ou serviço. A análise de variância consiste em decompor a variação total das observações do experimento em partes que podem ser atribuídas a causas conhecidas (tratamentos, etc) e em partes atribuídas a causas não controladas ou não controláveis (erro ou resíduo). Ou seja, a variação é vista de modo que:

Variação Total = Variação Controlada + Variação não Controlada.

2 Referencial Teórico

Para a construção da ANOVA vamos apresentar algumas regras a seguir, mais conhecidas como pressupostos.

  1. Os efeitos principáis devem ser aditivos: Nos experimentos, cada observação segue um modelo linear aditivo \(y_{ij}\) = \(\mu\) + \(\tau_{i}\) + \(\epsilon_{ij}\).

  2. Os erros de observação são independentes, cov(\(\epsilon_{ij}\)\(\epsilon_{ij}'\)) = 0.

  3. Os erros são homoscedásticos: Cada tratamento deve ter aproximadamente a mesma variância.

  4. Os erros são normalmente distribuídos: Para que a ANOVA seja válida, os erros devem ser originários da mesma população.

Em resumo, temos que \(\epsilon_{ij}\) ~ N(0,\(\sigma^{2}\)).

Além dos pressupostos que a ANOVA precisa obedecer, precisamos verificar as condições da hipótese e seu objetivo. O objetivo da hipótese é, em geral, verificar se existe diferença significativa entre pelo menos duas médias de tratamentos. Portanto, veremos a seguir as hipóteses testadas.

\(H_{0}\): \(\mu_{1}\) = \(\mu_{2}\) = … = \(\mu_{a}\)

\(H_{1}\): \(\mu_{i} \neq \mu_{i}'\) para pelo menos um \(i \neq i'\).

Uma forma equivalente de usar essa hipótese é aplicando \(\tau\), onde iremos usa-lo nas demais descrições da ANOVA. Logo,

\(H_{0}\): \(\tau_{1}\) = \(\tau_{2}\) = … = $_{a} $

\(H_{1}\): \(\tau_{i} \neq 0\) para pelo menos um \(i\).

2.1 Estimação dos parâmetros - Métodos dos minímos quadrados.

Nessa seção vamos mostrar que essas são prioridades para o ajuste por quadrados minímos decorrente do sistema de equações normais. Esse método é muito conhecido quando se trata da anova pois é não viesado e eficiente na hora de calcular o erro na reta da regressão, ela calcula a distãncia dos pontos (observações ou erros) da reta regressora. Esse método é eficiente pois como mostraremos na construção a seguir, o calculo proposto consiste em encontrar o melhor ajuste para um conjunto de dados tentando minimizar a soma dos quadrados das diferenças entre o valor estimado e os dados observados. O termo “ao quadrado” é importante pois aplicando essa função faz com que os dados obtenham o mesmo comportamento tanto para os dados positivos quanto para os negativos. Em resumo, queremos minimizar a distância aos valores dados, isto é importante em termos de aplicações, já que podemos ter valores obtidos, experimentalmente, com uma certa incerteza. Assim, a construção da equação é de forma,

\(L(\mu, \tau_{1},...,\tau_{a})\) = \(\sum_{i=1}^{a} \sum_{j=1}^{n} \epsilon_{ij}^2\) = \(\sum_{i=1}^{a} \sum_{j=1}^{n}(y_{ij}-\mu - \tau_{i})^2\)

Derivando-se \(L\) em relação a cada um dos parâmetros, tem-se

\(\frac{\partial L(\mu,\tau_{1},...,\tau_{a})}{\partial \mu}\) = -2\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij}-\mu - \tau_{i})\)

\(\frac{\partial L(\mu,\tau_{1},...,\tau_{a})}{\partial \tau_{i}}\) = -2\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij}-\mu - \tau_{i})\)

Igualando os resultados a zero, temos:

-2\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij}-\hat\mu - \hat\tau_{i})\) = 0 e -2\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij}-(\hat\mu) - \hat\tau_{i})\) = 0

Aplicando os somatórios, obtem-se o sistema de equações normais:

\(\sum_{i=1}^{a} \sum_{j=1}^{n} y_{ij}\) = \(an\hat\mu\) + \(n\sum_{i=1}^{a} \hat\tau_{i}\) \(\sum_{ij}^{n}\) = \(n\hat\mu + n\hat\tau_{i}\).

Dessa formas, vemos que o conjunto de equações normais não é linearmente independente. Não existe solução única para os parâmetro a serem estimados.

Dessa forma, uma solução para isso é, impor a restrição de que \(\sum_{i=1}^{a} \hat\tau_{i}\) = 0. Portanto, devido a isso temos que \(\hat\mu = \overline{y}..\) e \(\hat\tau_{i} = \overline{y}_{i.} - \overline{y}_{..}\), com \(i=1,...,a\).

Assim, decompondo as somas dos quadrados temos que a variabilidade total dos dados pode ser reescrita em função das médias de cada tratamento, assim

\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij} - \overline{y}_{..})^2\) = \(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij} - \overline{y}_{.} + \overline{y}_{.} - \overline{y}_{..})^2\).

Fazendo toda a manipulação com os somatórios, chegamos em

\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij} - \overline{y}_{..})^2\) = \(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij} - \overline{y}_{.})^2\) +

\(\sum_{i=1}^{a} \sum_{j=1}^{n} (\overline{y}_{.} - \overline{y}_{..})^2\).

Logo, desenvolvendo os quadrados, obtem-se que

\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij} - \overline{y}_{..})^2\) = SQT, que é a sooma de quadrados total, representando a variação total de y em torno de sua média.

\(\sum_{i=1}^{a} \sum_{j=1}^{n} (y_{ij} - \overline{y}_{.})^2\) = SQE, o que significa que obtivemos a soma de quadrado dos erros representando o y em torno da reta em um modelo.

\(\sum_{i=1}^{a} \sum_{j=1}^{n} (\overline{y}_{.} - \overline{y}_{..})^2\) = SQReg, representando a variação das esperanças especificas de y, dado x, em torno de sua média.

Assim, podemos expressar SQT = SQE + SQReg.

Para as observações, a variância é calculada como \(\frac{SQ}{gl}\), onde SQ representa a soma de quadrados dos desvios com respeito a média, e gl são os graus de liberdades.

Na ANOVA têm-se quantidades semelhantes a variâncias, chamadas de quadrados médios (QM), calculados como:

  • QMReg = \(\frac{SQReg}{glreg}\)

  • QME = \(\frac{SQE}{gl}\)

Observamos acima que não se trata exatamente de médias pois o denominador não é n. Das observações, vemos que:

Existe N = an observações totais, logo a SQT tem N - 1 graus de liberdade. Os tratamentos são representados por “a” observações, assim SQReg tem a-1 graus de liberddades. As replicações são representadas por “n”, assim proporcionam n-1 graus de liberdade para estimar o erro experimental. Dado que existem “a” tratamentos, SQE tem a(n-1) graus de liberdade (ou N-a gl).

Daí, QMReg = \(\frac{SQReg}{a-1}\) e QME = \(\frac{SQE}{N-a}\).

3 Tabela da ANOVA.

Onde SQTr representa SQReg e QMTr reepresenta QMReg.

4 Forma para testar a igualdade de várias médias

4.1 Modelo estatístico

O modelo para descrever os dados de um experimento pode ser escrito como \(y_{ij} = \mu + \epsilon_{ij}\), sendo \(i = 1,...,a\), \(j = 1, ..., n\) \((1)\)

onde,

  • \(y_{ij}\) representa a \(ij - ésima\) observação
  • \(\mu_{i}\) é a média do \(i-ésima\) nível do fator ou tratamento
  • \(\epsilon_{ij}\) é o erro aleatório que incorpora todas as outras fontes de variabilidade no experimento (medição, fatores não controlados, diferenças entre unidades experimentais, ruídos, etc). Assume-se que \(E(\epsilon_{ij})\) = 0 de modo que \(E(\epsilon_{ij})\) = \(\mu_{i}\)

Observação: O modelo na Equação (1) é chamado Modelo de médias.

4.2 Forma alternativa do modelo

Seja \(\mu_{1}\) = \(\mu + \tau\) \(i= 1,2,...,a\)

A Equação (1) torna-se o Modelo de efeitos, escrito como \(y_{ij} = \mu + \tau + \epsilon_{ij}\), \(\\ i = 1,...,a\), \(\\ j = 1,...,n\), \((2)\),

onde

  • \(y_{ij}\) é o valor observado na unidade j que recebeu o tratamento i

  • \(\mu\) é um parâmetro constante, comum a todos os tratamentos, chamado média geral (quando os dados são balanceados)

  • \(\tau\) é um parâmetro único que representa o efeito do i-ésimo tratamento

  • \(\epsilon\) é um componente do erro aleatório, associado à j-ésima repetição do i-ésimo tratamento.

  • Se os tratamentos foram selecionados especificamente pelo experimentador, o interesse é testar hipóteses sobre as médias dos tratamentos e as conclusões aplicam apenas aos níveis do fator considerados na análise. Interessa também estimar os parâmetros do modelo \((\mu, \tau_{i}\), \(\sigma^{2})\) \(=>\) Modelo de efeitos fixos.

  • Se os a tratamentos foram selecionados como uma amostra aleatória de uma população maior de tratamentos, o interesse é estender as conclusões a todos os tratamentos na população. Nesse caso, os \(\tau_{i}\) são variáveis aleatórias e os testes de hipóteses recaem sobre a variabilidade de \(\tau_{i}\) tentando-se estimar essa variabilidade \(=>\) Modelo de efeitos aleatórios ou Modelo de componentes de variância.

4.3 Modelo de efeitos aleátorios

\(y_{ij}\) = \(\mu\) + \(\tau_{i}\) + \(\epsilon_{ij}\), i = 1,…,a, j = 1,…,n, onde,

  • \(y_{ij}\) é a abservação do i-ésimo tratamento na j-ésima unidade experimental (variável resposta).
  • \(\mu\) é a média global.
  • \(\tau_{i}\) é o efeito do i-ésimo tratamento.
  • \(\epsilon_{ij}\) é o erro associado ao i-ésimo tratamento na j-ésima unidade experimental. Imcorpora todas as fontes de variabilidade não controladas no experimento.

Assim, temos que,

  • \(\tau_{i}\) e \(\epsilon_{ij}\) são variáveis aleatórias independentes entre si, com \(\tau_{i}\) ~ N(0,\(\sigma_{\tau}^{2}\)) e \(\epsilon_{ij}\) ~ N(0,\(\sigma^{2}\)).
  • A condição E(\(\tau_{i}\)) = 0 é similar á condição \(\sum_{i=1}^{a} \tau_{i}\) = 0. Estabelece que o efeito esperado do i-ésimo nível como um desvio de \(\mu\) é zero.
  • Se a variância dos efeitos dos tratamentos \(\tau_{i}\) e \(\sigma_{\tau}^{2}\), a variância da resposta é dada por var(\(y_{ij}\)) = \(\sigma_{\tau}^{2}\) + \(\sigma^{2}\). As variâncias \(\sigma_{\tau}^{2}\) e \(\sigma^{2}\) são chamadas componentes de variância.

Assim, nesse modelo os efeitos têm interpretação diferente. As hipóteses testadas são:

\(H_0\): \(\sigma_{\tau}^{2}\) = 0

\(H_1\): \(\sigma_{\tau}^{2}\) \(>\) 0

Logo, se \(\sigma_{\tau}^{2}\) = 0, todos os tratamentos são idênticos, mas se \(\sigma_{\tau}\) > 0, há variabilidade entre os tratamentos.

4.4 Modelos de efeitos fixos

\(y_{ij}\) = \(\mu + \tau_{i} + \epsilon_{ij}\), com i = 1,…,a, j = 1,…,n, onde

A seguir apresentaremos um experimento com um único fator.

5 Experimento com um único fator

Exempo 1, aplicação.

Assume-se que há tratamentos ou diferentes níveis de um único fator a serem comparados.

Tratamento(Níveis) Observações Total Média
1 \(y_{11} y_{12} · · · y_{1n}\) \(y_{1.}\) \(\overline{y}_{1}\)
2 \(y_{21} y_{22} · · · y_{2n}\) \(y_{2.}\) \(\overline{y}_{2}\)
. . . . . . .
. . . . . . .
. . . . . . .
a \(y_{a1} y_{a2} . . . y_{an}\) \(y_{a.}\) \(\overline{y}_{a}\)
\(y_{..}\) \(\overline{y}_{..}\)
  • \(y_{ij}\) representa a j-ésima observação do nível ou do tratamento \(i\).
  • \(y_{i}\)= \(\sum_{j=1}^{n}\), \(y_{i}\) = \(\sum_{j=1}^{n}\) \(\frac{y_{ij}}{n}\) \(\overline{y}_{..}\) = \(\sum ^{a}_{i=1} \sum _{j=1}^{n}\) \(\frac{y_{ij}}{kn}\)

Considere um experimento cujo objetivo é verificar se a inclusão de raízes e tubérculos, como suplementação de inverno na alimentação de vacas em lactação, aumenta a produção de leite. Consideram-se 24 animais, três tipos de suplementos e uma testemunha (placebo), que são:

  • Sem suplemento (S)
  • Mandioca (M)
  • Araruta (A)
  • Batata doce (B).

Para definir o tipo de suplemento que será dado a cada animal, realiza-se um sorteio aleatório enumerando cada um dos 24 animais (parcelas) que participarão do estudo (1 a 24) e, em seguida, colocam-se os tratamentos em uma sequência, como a dada a seguir:

  • S1 S2 S3 S4 S5 S6
  • M1 M2 M3 M4 M5 M6
  • A1 A2 A3 A4 A5 A6
  • B1 B2 B3 B4 B5 B6

Utilizando um gerador de números aleatórios, aloca-se o tipo de suplemento a cada animal. Suponha que a sequência de números aleatórios sorteada, tenha sido

  • 24 23 22 14 1 13
  • 6 20 8 7 9 4
  • 21 15 17 16 19 2
  • 11 5 10 3 18 12.

Assim, tem-se a configuração do experimento:

Vaca Trat Vaca Trat Vaca Trat Vaca Trat
1 \(S_{5}\) 7 \(M_{4}\) 13 \(S_{6}\) 19 \(S_{5}\)
2 \(A_{6}\) 8 \(M_{3}\) 14 \(S_{4}\) 20 \(S_{5}\)
3 \(B_{4}\) 9 \(M_{5}\) 15 \(A_{2}\) 21 \(S_{5}\)
4 \(M_{6}\) 10 \(B_{3}\) 16 \(A_{4}\) 22 \(S_{5}\)
5 \(B_{2}\) 11 \(B_{1}\) 17 \(A_{3}\) 23 \(S_{5}\)
6 \(M_{1}\) 12 \(B_{6}\) 18 \(B_{5}\) 24 \(S_{5}\)

Considerem-se as seguintes produções médias diárias (kg) de leite a 4% de gordura das vacas submetidas a administração de raízes e tubérculos, como suplementação de inverno na alimentação de vacas em lactação.

Id. Prod. Id. Prod. Id. Prod. Id. Prod.
1 22,81 7 25,12 13 23,54 19 35,04
2 35,19 8 24,36 14 25,42 20 22,37
3 20,37 9 22,94 15 32,47 21 35,42
4 24,80 10 26,54 16 34,48 22 23,43
5 24,37 11 22,15 17 35,04 23 21,07
6 23,40 12 24,06 18 19,54 24 19,58
Trat. Observações Média
S 19,58 21,07 23,43 25,42 22,81 23,54 22,64
M 23,40 22,37 24,36 25,12 22,94 21,56 23,29
A 35,42 32,47 34,48 33,79 35,04 35,19 34,39
B 22,15 24,37 26,54 20,37 19,54 24,06 22,83

Seja \(y_{ij}\) o valor da produção de leite da j-ésima vaca que recebeu o j-ésima tratamento. Os valores das produções (kg) de leite a 4% de gordura das vacas que participaram do estudo podem ser resumidos na forma:

Objetivo: Testar se há diferença na produção média de leite de acordo com o tipo de suplementação.

6 Análise descritiva

Resultados

suple<-factor(rep(c('A','B','M','S'), each=6))
x<-c(19.58,21.07,23.43,25.42,22.81,23.54,23.40,22.37,24.36,25.12,22.94,
     21.56,35.42,32.47,34.48,33.79,35.04,35.19,22.15,24.37,26.54,20.37,
     19.54,24.06)
dados<-data.frame(x,suple)

media<-tapply(x,suple,mean)
dp<-tapply(x,suple,sd)
var<-tapply(x,suple,var)
summary(x)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   19.54   22.32   23.80   25.79   28.02   35.42
media
##        A        B        M        S 
## 22.64167 23.29167 34.39833 22.83833
dp
##        A        B        M        S 
## 2.050360 1.301360 1.111529 2.645233
var
##        A        B        M        S 
## 4.203977 1.693537 1.235497 6.997257

Box-Plot

Medida descritiva

Tratamento Média Desvio padrão
S 22,64 2,05
M 23,29 1,30
A 34,40 1,11
B 22,84 2,65

Análise do exemplo 1

Pelo box-plot conseguimos ver que batata doce tem uma maior variabilidade, tem maior concentração de dados, mas não representa o melhor suplemento tubérculo. Deste modo vemos que por mais que araruta tenha uma menor concentração, menor variabilidade, ela representa a melhor raíz como suplemento de inverno na alimentação das vacas.

Complementando a análise, vamos aplicar o test t para amostras independentes e análisar todos os pares de médias resultantes. Para isso, iremos definir nossas hipóteses do teste. Assim,

\(H_{0}\) : \(\mu_{1}\) - \(\mu_{2}\) = 0 vs. \(H_{1}\) : \(\mu_{1}\) - \(\mu_ {2}\) \(\neq 0\).

Logo, aplicaremos a estatística de teste que é ela:

\(t_{0}\) = \(\frac{ (\overline{y1} - \overline{y2} ) - (\mu_{1} - \mu_{2})}{Sp\sqrt[]{\frac{1}{n1}+\frac{1}{n2}}}\),

onde

\(S{p}^{2}\) = \(\sqrt{\frac{(n_{1}-1)S_{}1{p}^{2} + (n_{2}-1)S_{2}{p}^{2}}{n_{1} +n_{2} - 2 }}\).

Feito isso, podemos então dar continuidade a análise onde mostraremos os resultados das hipóteses ditas acima.

Resultados

Seja a matriz de confusão dada a seguir

S M A B
S \(-\) ND \(**\) ND
M \(-\) \(-\) \(**\) ND
A \(-\) \(-\) \(-\) \(**\)
B \(-\) \(-\) \(-\) \(-\)

  • Onde ND significa SEM DIFERENÇA SIGNIFICATIVA

  • Onde ** significa DIFERENÇA SIGNIFICATIVA AO NÍVEL DE 5%.

Dado essas informações, podemos dizer que a solução proposta nos leva a distorção do erro tipo I. Ou seja, rejeitamos a hipotese nula dado que ela é verdadeira. O correto seria aceitar a hipótese nula. A probabilidade de cometer um erro do tipo I é \(\alpha\), que é o nível de significância que você definiu para seu teste de hipóteses. Um \(\alpha\) de 0,05 indica que você quer aceitar uma chance de 5% de que está errado ao rejeitar a hipótese nula

Vejamos a explicação para isso. Suponha que seja testada a igualdade das quatro médias usando comparações pareadas. Há 6 pares possíveis e, se a probabilidade de aceitar corretamente a hipótese nula para cada par testado é de (1 - \(\alpha\)) = 0,95, então a probabilidade de aceitar corretamente a hipótese nula para todos os 6 pares é \((0,95)^{6}\) = \(0,7359\), se os testes forem independentes.

Aplicação do modelo do exemplo 1

modelo<-aov(x~suple)
names(modelo)
##  [1] "coefficients"  "residuals"     "effects"       "rank"         
##  [5] "fitted.values" "assign"        "qr"            "df.residual"  
##  [9] "contrasts"     "xlevels"       "call"          "terms"        
## [13] "model"
coef(modelo)    
## (Intercept)      supleB      supleM      supleS 
##  22.6416667   0.6500000  11.7566667   0.1966667

Assim, usando como casela de referência o suplemento A, obtivemos um intercepto de 22,64, ou seja, o valor da inicialização do modelo, mais conhecida como \(\alpha\) ou \(\beta_{0}\). Tendo os outros coefcientes dados por 0.650, 11.756, 0.196 que são usualmente representados por \(\beta_{1}\), \(\beta_{2}\), \(\beta_{3}\) respectivamente.

Como dito nas seções acima, para a adequação do modelo e da ANOVA, devemos verificar a normalidade dos dados, então, aplicaremos 3 testes de normalidade nos dados.

require(nortest)
## Loading required package: nortest
lillie.test(modelo$res)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo$res
## D = 0.10192, p-value = 0.7484
shapiro.test(modelo$res)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$res
## W = 0.97917, p-value = 0.88
ad.test(modelo$res) 
## 
##  Anderson-Darling normality test
## 
## data:  modelo$res
## A = 0.23818, p-value = 0.7559

Assim, conforme mostram os testes de normalidade apicados no modelo, em todos eles não rejeitamos a hipótese nula que é dada por,

\(H_{0}\): p-valor > 0,05 vs \(H_{1}\): p-valor < 0.05.

Ou até mesmo \(H_{0}\): A amostra provém de uma população normal vs \(H_{1}\): A amostra não provém de uma população normal

Tendo conhecimento então das hipóteses e analisando os resultados dos testes propostos, não devemos rejeitar a hipótese nula, temos normalidade na amostra.

Outro pressuposto que devemos veridicar é a homoscedasticidade das variâncias. Uma maneira de testar é usando um teste específico, portanto, usaremos o teste de Bartlett.

require(car)
## Loading required package: car
## Loading required package: carData
bartlett.test(x~suple)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  x by suple
## Bartlett's K-squared = 4.2843, df = 3, p-value = 0.2324
leveneTest(modelo)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  3  1.9172 0.1593
##       20

Notem o p-valor do teste, ele resume o resultado de um teste de hipótese. Quando p-valor for menor que o nível de significância escolhido para o teste (em geral, convencionado em 0,05), rejeita-se a hipótese nula em favorecimento da hipótese alternativa - neste caso assume-se a desigualdade das variâncias; Caso ele seja maior que o nível de significância, falha-se em rejeitar H0 por falta de evidências, sendo assim, conclui-se que as variâncias são constantes ou homogêneas. No caso do exemplo, o p-valor do teste de Bartlett e de Levene foi maior que 0,05, assim há homogeneidade das variâncias.

Após apresentar todos os passos da análise de variância vamos por último apresentar a tabela da ANOVA que pode ser obtidar usando “summary(modelo)” que teremos toda a significância do modelo construido ou de forma a escrever “anova(modelo)” que apresenta a tabela da anova. Portanto, a tabela a seguir apresenta a aplicação da ANOVA no modelo proposto.

anova(modelo)
## Analysis of Variance Table
## 
## Response: x
##           Df Sum Sq Mean Sq F value    Pr(>F)    
## suple      3 593.82 197.939  56.033 6.495e-10 ***
## Residuals 20  70.65   3.533                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Nas saídas da análise, pode-se notar a existência de uma diferença significativa entre os padrões de eventos dos quatro tratamentos. Notem que o p-valor estimado foi igual à aproximadamente 0 (Pr(>F)). Isto nos leva a rejeitar a hipótese nula em favorecimento à hipótese alternativa, ou seja, pelo menos uma dos tratamentos apresenta um padrão médio distinto das demais.

É importante verificar o pressuposto de independência das observações. Este pressuposto pode ser avaliado através da análise dos resíduos do modelo ANOVA. Aplicando a função genérica plot() do R aos resíduos do modelo, que podem ser acessados pela função residuals() ao objeto que armazena os resultados do modelo ANOVA pode-se visualizar a distribuição dos resíduos e sua independência à cada observação. Outra forma de visualizar esta independência dos dados sobre a análise dos resíduos de um modelo está centrada na distribuição dos mesmos. Espera-se que os resíduos de um modelo bem ajustado possua uma distribuição normal, centrada em 0 e variância constante.

Logo vemos que os dados não correspondem ao esperado, os resíduos não estão concentrados em torno da média 0 e muita variabilidade nos residuos, O histograma não apresenta muita normalidade nos resíduos. Vemos que mesmo tendo passado no teste de normalidade, os dados apresentam problemas na variabilidade e na normalidade.

Por fim, resta identificar qual ou quais áreas possuem média de ventos distintas entre sí. Para isto, utilizamos a função TukeyHSD() do R para contrastar as diferenças.

post.hoc <- TukeyHSD(modelo)
post.hoc
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = x ~ suple)
## 
## $suple
##            diff        lwr       upr     p adj
## B-A   0.6500000  -2.387229  3.687229 0.9311689
## M-A  11.7566667   8.719438 14.793896 0.0000000
## S-A   0.1966667  -2.840562  3.233896 0.9978139
## M-B  11.1066667   8.069438 14.143896 0.0000000
## S-B  -0.4533333  -3.490562  2.583896 0.9747902
## S-M -11.5600000 -14.597229 -8.522771 0.0000000

Com base no resultado do teste de Tukey, pode-se concluir nem que todos os tratamentos são diferentes uns dos outros pois não tivemos para todas as diferenças o valor de p-adj menor que 0,05, mas vemos que o que obtivemos os resultados esperados são para o suplemento A e B, todas as diferenças com os dois suplementos são significativas. Assim, o melhor suplemento é o Araruta (A).

Exemplo 2, aplicação.

Nesse exemplo iremos testar a absorção do cimento. O intuíto é testar o quao bom é o cimento em relação as misturas utilizadas. Em cada um deles foi usado um tipo de areia e queremos testar em minutos o tempo que cada tipo leva pra secar. O experimento foi separado em 6 blocos para cada tipo de cimento.

trat<-factor(rep(c('T1','T2','T3','T4', 'T5'), each=6))

absorcao<-c(2370, 1678, 2592, 2283, 2910, 3020, 1282, 1527, 871, 1025, 825, 920, 562, 321, 636, 317, 485, 842, 173, 127, 132, 150, 129, 227, 193, 71, 82, 62, 96, 44)

dados<-data.frame(trat,absorcao)

media<-tapply(absorcao,trat,mean)

dp<-tapply(absorcao,trat,sd)

var<-tapply(absorcao,trat,var)

media
##         T1         T2         T3         T4         T5 
## 2475.50000 1075.00000  527.16667  156.33333   91.33333
dp
##        T1        T2        T3        T4        T5 
## 486.42276 274.87961 200.31517  38.75908  52.83812
var
##         T1         T2         T3         T4         T5 
## 236607.100  75558.800  40126.167   1502.267   2791.867

Box-Plot

boxplot(absorcao~trat,ylab="Absorção do cimento", xlab="Misturas", names=c("Mistura 1", "Mistura 2", "Mistura 3", "Mistura 4", "Mistura 5"), col="pink")
points(media,pch="t",col="red")

Aplicação do modelo do exemplo 2

modelo<-aov(absorcao~trat)
names(modelo)
##  [1] "coefficients"  "residuals"     "effects"       "rank"         
##  [5] "fitted.values" "assign"        "qr"            "df.residual"  
##  [9] "contrasts"     "xlevels"       "call"          "terms"        
## [13] "model"
coef(modelo)
## (Intercept)      tratT2      tratT3      tratT4      tratT5 
##    2475.500   -1400.500   -1948.333   -2319.167   -2384.167

tendo então como casela de referência o tratamento T1, obtivemos um intercepto de 2475.500 como sendo \(\beta_{0}\), tendo os outros coefcientes dados por -1400.500, -1948.333, -2319.167, -2384.167 que são usualmente representados por \(\beta_{1}\), \(\beta_{2}\), \(\beta_{3}\) e \(\beta_{4}\) respectivamente.

A seguir, verificaremos a normalidade dos dados como feito no exemplo anterior, uma vez que sabemos a importância para a análise.

require(nortest)
lillie.test(modelo$res)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo$res
## D = 0.15343, p-value = 0.06941
shapiro.test(modelo$res)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$res
## W = 0.89447, p-value = 0.006168
ad.test(modelo$res) 
## 
##  Anderson-Darling normality test
## 
## data:  modelo$res
## A = 1.1824, p-value = 0.003643

Desse modo, em nenhum dos testes aplicados nos dados apresentaram normalidade, o que pode acontecer muitas das vezes em um conjunto de dados. Tendo em vista disso, usaremos uma trasformação nos dados que é bastante usual quando não encontramos normalidade, afim de corresponder ao pressuposto de normalidade. A transformação que iremos usar é a de Box-Cox, bastante usada quando não temos normalidade nos dados.

Transformação BOX-COX

Em estatística, uma transformação de potência é uma família de funções que são aplicadas para criar a transformação monotônica de dados usando funções de potência. Esta é uma técnica de transformação de dados útil usada para estabilizar a variância, tornar os dados mais semelhantes à distribuição normal, melhorar a validade das medidas de associação (como a correlação de Pearson entre as variáveis) e para outros procedimentos de estabilização de dados

Tanto a forma linear quanto a logarítmica são dois casos particulares de uma família mais extensa de transformações não-lineares. A transformação de potência é definida como uma função de variação contínua, em relação ao parâmetro de potência \(\lambda\) (lambda), ou seja, xλ. Uma classe geral de transformação que pode ser utilizada é a de Box-Cox, definida por:

\(f_{\lambda}\) = \(\frac{x^{\lambda} - 1}{\lambda}\) para \(\lambda\) \(\neq\) 0.

\(f_{0}\) = \(log(x)\) para \(\lambda\) = 0

Assim, usaremos a trasformação Box-Cox para tentar corrigir os dados e transforma-los em normais.

Aplicação Box-Cox

require(MASS)
## Loading required package: MASS
boxcox(modelo) 

dados<-data.frame(trat,"absorca"=log(absorcao))
attach(dados)
## The following object is masked _by_ .GlobalEnv:
## 
##     trat

Medidas descritivas após aplicado a transformação

media<-tapply(absorca,trat,mean)

dp<-tapply(absorca,trat,sd)

var<-tapply(absorca,trat,var)

media
##       T1       T2       T3       T4       T5 
## 7.796392 6.954847 6.206162 5.029280 4.401294
dp
##        T1        T2        T3        T4        T5 
## 0.2123848 0.2413457 0.3867126 0.2265903 0.4987614
var
##         T1         T2         T3         T4         T5 
## 0.04510732 0.05824776 0.14954661 0.05134315 0.24876293

Agora faremos o Box-Plot dos dados após a transformação para que passamos verificar a diferença dos dados e como as classes sofrem alteração.

boxplot(absorca~trat,ylab="Absorção do cimento", xlab="Misturas", names=c("Mistura 1", "Mistura 2", "Mistura 3", "Mistura 4", "Mistura 5"), col="pink")
points(media,pch="t",col="red")

Assim vemos que a variabilidade das misturas sofrem alteração, a mistura 1 apresentava a maior variabilidade das misturas e agora vemos que a maior variabilidade está concentrada na mistura 3. As misturas que tinham os dados bem abaixos agora apresentam uma nitidez, não estão sendo puxados tornando a variabilidade bem baixa e concentrados no 0, após essa transformação, os dados tiveram alteração nas classes tornando mais visível e mais fácil para a interpretação.

Desse modo, ajustaremos um novo modelo com a trasformação nos dados.

Aplicação do modelo do exemplo 2

modelo<-aov(absorca~trat)
names(modelo)
##  [1] "coefficients"  "residuals"     "effects"       "rank"         
##  [5] "fitted.values" "assign"        "qr"            "df.residual"  
##  [9] "contrasts"     "xlevels"       "call"          "terms"        
## [13] "model"
coef(modelo)
## (Intercept)      tratT2      tratT3      tratT4      tratT5 
##   7.7963924  -0.8415451  -1.5902303  -2.7671127  -3.3950988

A casela de referência continuou sendo o tratamento T1, o intercepto passa a ter um valor de 7.796 e, \(\beta_{1}\), \(\beta_{2}\), \(\beta_{3}\) e \(\beta_{4}\) tendo valores de -0.8415, -1,5902, -2,7671, -3.3950 respectivamente.

Para verifircarmos se essa transformação foi efiiente, faremos os mesmos testes de normalidade usados nas aplicações anteriores, uma vez que usamos a trasnformação para tornar os dados normais. Desse modo,

require(nortest)
lillie.test(modelo$res)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo$res
## D = 0.09299, p-value = 0.7302
shapiro.test(modelo$res) 
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$res
## W = 0.97645, p-value = 0.7255
ad.test(modelo$res)
## 
##  Anderson-Darling normality test
## 
## data:  modelo$res
## A = 0.25301, p-value = 0.7115

Assim, vemos que a trasformação usual nos dados foi eficiente e agora respondemos ao pressuposto da ANOVA de normalidade, em todos os testes não se deve rejeitar a hipótese nula, os dados provém de uma população normal.

Portanto, daremos continuidade na análise e verificaremos a homoscedastícidade das variâncias usando o teste de Bartlett e de Levene.

require(car)
bartlett.test(absorca~trat) 
## 
##  Bartlett test of homogeneity of variances
## 
## data:  absorca by trat
## Bartlett's K-squared = 5.5244, df = 4, p-value = 0.2376
leveneTest(modelo)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  4  0.9975 0.4274
##       25

Logo o p-valor é maior que 0.05, assumimos então que as variâncias são constantes, não devemos rejeitar \(H_{0}\). Isso significa que as variâncias são iguais pelo teste de Levene.

Uma vez tendo sucesso em todos os pressupostos para a ANOVA, podemos então apresentar a tabela da ANOVA para análise.

anova(modelo)
## Analysis of Variance Table
## 
## Response: absorca
##           Df Sum Sq Mean Sq F value    Pr(>F)    
## trat       4 45.896 11.4740  103.74 3.454e-15 ***
## Residuals 25  2.765  0.1106                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Na saída da análise podemos encontrar também que existe a evidência de uma diferença significativa entre os padrões dos eventos de tratamento dos cimentos. Nota-se que o p-valor estimado está aproximadamente a 0, rejeitamos então a hipótese nula, assim, pelo menos um dos tratamentos apresenta um padrão médio distinto dos demais.

Para os resíduos analizaremos os gráficos conforme no exemplo 1.

Logo vemos que os dados correspondem ao esperado, os resíduos estão concentrados em torno da média 0 e se colocarmos um instervalo em -0.5 e 0.5 encontramos apenas 2 pontos fora do intervalo. O histograma mostra que os resíduos temdem a uma distribuição normal. Dessa forma vemos que mesmo que encontremos problemas nos dados na análise proposta, existem maneiras estatística de corrigi-la e adequa-la a análise.

E para finaliza, faremos o teste de Tukey nos dados para assim cocncluirmos qual é o tipo de cimento que tem o menor tempo de secagem.

post.hoc <- TukeyHSD(modelo)
post.hoc
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = absorca ~ trat)
## 
## $trat
##             diff       lwr        upr     p adj
## T2-T1 -0.8415451 -1.405449 -0.2776414 0.0015853
## T3-T1 -1.5902303 -2.154134 -1.0263267 0.0000001
## T4-T1 -2.7671127 -3.331016 -2.2032090 0.0000000
## T5-T1 -3.3950988 -3.959002 -2.8311952 0.0000000
## T3-T2 -0.7486853 -1.312589 -0.1847816 0.0052888
## T4-T2 -1.9255676 -2.489471 -1.3616640 0.0000000
## T5-T2 -2.5535537 -3.117457 -1.9896501 0.0000000
## T4-T3 -1.1768823 -1.740786 -0.6129787 0.0000193
## T5-T3 -1.8048684 -2.368772 -1.2409648 0.0000000
## T5-T4 -0.6279861 -1.191890 -0.0640825 0.0237089

Vemos então que em modo de significânci dois a dois, todos eles são diferentes umm dos outros. O melhor tipo de cimento é o T5. Podemos concluir também que o cimento T4 também seria uns dos melhores para utilizarmos.

7 Conclusões

Nesse trabalho observamos que para a construção da ANOVA não basta somente mostrar a tabela e seus resultados, vemos que segue toda uma aplicação nos dados respeitando os pressupostos e aplicações de cada teste. Nesse trabalho verificamos que mesmo que os pressupostos não sejam aceitos, temos outras maneiras estatísticas para transformar os dados e torna-los eficientes e que respeitem os pressupostos da ANOVA. Concluimos também que a anova é muito importante qando precisamos e desejamos comparar médias, verificar se são iguais ou diferentes pois sabemos que isso causa efeitos diferentes na análise, assim como na inferência de uma análise. Assim, foi apresentado passos que devem ser usados em toda a construção antes de tomar a decisão.

Referências

[1] - http://lite.acad.univali.br/rcurso/anova/

[2] - http://www.portalaction.com.br/anova

[3] - Damodar N. Gujarati Dawn C. Porter. Econometria Básica. 5ªe, 2011.

[4] - Portal action. ANOVA. http://www.portalaction.com.br/anova.

[5] - Sonia Vieira. Análise de Variância ( Anova ).

[6] - Reinaldo Charnet, Clarice Azevedo de Luna Freire, Eugênia M. Reginato Charnet, Heloísa Bonvino. Análise de Modelos de Regressão Linear com Aplicações. 2ªe (1 de janeiro de 2008).