1.1

O modelo para esse tipo de análise é:

\[y_{ij}=\mu+\tau_{i}+\varepsilon_{ij}\]

onde \(\varepsilon_{ij}\sim N(0,\sigma^2)\)

No qual testamos se os valores dos tratamentos \(\tau_i\) são iguais ou seja:

valores <- c(99.1, 98.1, 99.8, 99.6, 101,
            117.4, 117.4, 113, 113, 111,
            117.5, 114.9, 113.2, 121.3, 110.4,
            97, 94.4, 94, 94, 94.8)
  
  
trat<- factor(rep(1:4, each = 5))


modelo <- aov(valores~trat)


summary(modelo)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## trat         3 1629.5   543.2   76.57 1.05e-09 ***
## Residuals   16  113.5     7.1                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Conforme a tabela acima, com um nível de significância de 5%, temos que existe algum valor \(\tau_i\) que seja diferente os demais, ou seja, existe alguma média do grupo diferente entre os demais.

1.2

Para avaliar a normalidade, utilizamos o resíduo como modo geral de todos os dados (em vez de testar normalidade nos grupos de tratamento). Assim, temos como no gráfico abaixo:

hist(modelo$residuals)

notamos que é dificil definir se possui normalidade, assim, testamos a hipótese de normalidade de shapiro-wilk. No qual

shapiro.test(modelo$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$residuals
## W = 0.96651, p-value = 0.6801

Assim com o resultado acima, com o nível de significância de 0.05, não há evidência de rejeitar a normalidade.

1.3

Para verificar a hipótese de homogêniedade de variância, temos o teste de Bartlet

bartlett.test(valores, trat)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  valores and trat
## Bartlett's K-squared = 8.3324, df = 3, p-value = 0.03962

no qual com há indicios com, 5% de significância, de rejeitar a homogêniedade de variância.

Teste F

y_t <- matrix(valores, ncol = 5, byrow = T)

var_y <- apply(y_t, 1, var)
  
1-pf(max(var_y)/min(var_y),4,4)  
## [1] 0.0105634

Com o mesmo nível de significância, rejeitamos a homogêniedade de variância.

Teste levene

library(lawstat)
## Warning: pacote 'lawstat' foi compilado no R versão 4.4.2
levene.test(valores, trat, location = 'mean')
## 
##  Classical Levene's test based on the absolute deviations from the mean
##  ( none not applied because the location is not set to median )
## 
## data:  valores
## Test Statistic = 4.0607, p-value = 0.02532

Com o mesmo nível de significância, rejeitamos a homogêniedade de variância.

para a análise gráfica, temos

boxplot(valores~trat)

Assim, temos que há evidência gráfica de diferença nas variâncias.

1.4

Para testar a igualdade de variância quando não há normalidade nos resíduos temos que pelo teste de kruskal-wallis:

Onde as hipóteses são:

kruskal.test(valores, trat)
## 
##  Kruskal-Wallis rank sum test
## 
## data:  valores and trat
## Kruskal-Wallis chi-squared = 16.248, df = 3, p-value = 0.001009

Assim, rejeitamos a hipótese de igualdade de distribuição acumulada, ou seja, existe um trataménto que é diferente entre os outros.

1.5

require(MASS)
## Carregando pacotes exigidos: MASS
## Warning: pacote 'MASS' foi compilado no R versão 4.4.2
resultados <- boxcox(modelo, lambda = seq(-10, 0, 0.1))

lambda <-  resultados$x[which.max(resultados$y)]

Como pelo intervalo de confiança não contem o \(\lambda =1\), então é aconselhável fazer a transformação…

new_model <-aov(valores^(lambda)~trat)

summary(new_model)
##             Df    Sum Sq   Mean Sq F value   Pr(>F)    
## trat         3 6.791e-22 2.264e-22   150.5 6.19e-12 ***
## Residuals   16 2.410e-23 1.500e-24                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Teste de variância
bartlett.test(valores^(lambda)~trat)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  valores^(lambda) by trat
## Bartlett's K-squared = 0.9857, df = 3, p-value = 0.8047
# Teste de normalidade
shapiro.test(new_model$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  new_model$residuals
## W = 0.96641, p-value = 0.6781

Ou seja, depois da transformação com o \(\hat{\lambda}= -5.4\), temos como evidência em m nível de confiança de 95%, que as variâncias são iguais e que há normalidade nos dados, assim, validando o modelo ANOVA, no qual há eviência de que existe algum \(\tau_i\) diferente entre os demais.

SIMULAÇÃO

1.1

Temos dois tipos de erro

  • Erro tipo I: \(P(\text{Rejeitar}H_0|H_0\text{ Verdadeira})\)
  • Erro tipo II: \(P(\text{Não rejeitar }H_0|H_0\text{ Falso})\)
set.seed(10)

simulacoes = 1000
mu = 100
sigma2 = 16
n=5

t = c(5, -5, 0, 0)


trat = factor(rep(1:length(t), each = n))

taus = rep(t, each = n)

p_valores = c()


for (simu in 1:simulacoes) {
  e = rnorm(n*length(t), 0, sqrt(sigma2))

  y = mu + taus + e
  
  p_valores[simu] = summary(aov(y~trat))[[1]][1,5]
}


mean(p_valores > 0.05)
## [1] 0.137

Sabendo que a distribuição F não tem padronização, a estatistica de teste quando mede o erro tipo II, temos que dispõe do numerador seguindo uma não centralidade, ou seja

\[{MS_{treat}\over MS_{res}} = {\sum_{i=1}^a (\bar{y}_{i.}-\bar{y}_{..})^2 /(a-1) \over \sum_{j=1}^n\sum_{i=1}^a (\bar{y}_{ij}-\bar{y}_{i.})^2 /(an-a)}\sim F(a-1,an-1, ncp)\;,ncp={n\sum_{i=1}^a \tau_i^2\over\sigma^2}\]

Assim, para calcular analíticamente o Erro tipo II, temos que, com nível de significância de 0.05:

  • Achar valor crítico, na distribuição central \(f_{crit}\sim F(a-1,an-1, 0)\) com \(\alpha=0.05\)

  • Achar valor a probabilidade acumulada até o ponto crítico sob a hipótese H1, ou seja, que a distribuição não é central com a acumulada até o valor crítico da central, assim temos o valor do erro tipo II.

n=5;a=4
taus  = c(5, -5, 0, 0)
sigma2 = 16

alpha =  0.05

# nao centralidade
ncp = n*sum(taus^2)/sigma2

f_critc = qf(1-alpha, a-1, a*n-a)

# ERRO  TIPO II
pf(f_critc, a-1, a*n-a, ncp)
## [1] 0.1534387

Assim o valor análiticamente do erro tipo II é 0.15, que representa o valor acumulado até a reta vertical da linha vermelha.

xx <- seq(0,30,0.1)


plot(xx, df(xx, a-1, a*n-a), type='l')
abline(v=f_critc)

lines(xx, df(xx, a-1, a*n-a, ncp), col=2)