[R-br] Heterocedasticidade e teste de média

[Walmes Zeviani]

Pessoal,

Estes dados apresentam um padrão peculiar onde a solução transformação
Box-Cox não resolve. O Wagner mostrou um ajuste que modela a variância com
uma função. Aqui estou repetindo o que ele fez, usando funções do pacote
nlme, para fazer teste de razão de verossimilhança. Como a relação média
variância nos dados não é exatamente uma função monótona, o ajuste com
função Exp para a variância não resultou em melhoria. Então eu modelei a
variância em função de grupos. Pro caso de comparar os tratamentos com a
testemunha, reordene os níveis do fator para a testemunha ser o primeiro
nível. Assim, o summary vai exibir estes contrastes.

Fenois <- c(337.311, 344.874, 342.353, 325.546, 333.950, 330.588, 328.067,
            328.067, 318.824, 331.429, 333.950, 334.790, 336.471, 338.151,
            342.353, 259.160, 252.437, 268.403, 265.882, 266.723, 287.731,
            88.571, 88.571,  90.252,  41.513,  52.437,  49.076,  88.571,
            88.571,  90.252,  64.202,  60.000,  61.681)
Cor <- factor(c(rep("ambar",6), rep("ambar_claro",3), rep("ambar",6),
                rep("ambar_claro",6),rep("branco",6),
                rep("extra_ambar_claro",3),rep("branco",3)))

dados <- data.frame(Fenois, Cor)

xtabs(~Cor, data = dados)

library(lattice)
library(latticeExtra)

xyplot(Fenois ~ Cor, data = dados)

# ATENÇÃO. O padrão a relação média variância não é uma função
# monótona, ou seja, pelo gráfico, a variância não aumenta com a média
# pois o trat 3 tem mesma média do 4 mas com variância bem maior o trat
# 1 que é o de maior média tem variância baixa.

mv <- aggregate(Fenois ~ Cor,
                data = dados,
                FUN = function(y) {
                    cbind(m = mean(y), v = var(y))
                })
str(mv)

xyplot(Fenois[, 2] ~ Fenois[, 1],
       data = mv,
       xlab = "Média", ylab = "Variância",
       type = c("p", "r")) +
    layer(panel.text(x = x, y = y, labels = mv$Cor, pos = 3))

m0 <- lm(Fenois ~ Cor, data = dados)

par(mfrow = c(2, 2))
plot(m0)
layout(1)

MASS::boxcox(m0)

# A Box-Cox sugere transformação quando existe relação média variância
# por função monótona ou problema de simetria. Esse padrão ela corrige.

library(nlme)

m1 <- gls(Fenois ~ Cor, data = dados, weights = varExp(), method = "ML")

logLik(m0)
logLik(m1)

anova(m1, m0)

# Razão entre modelos aponta que a função de variância exponencial (que
# é monótona) não difere do modelo nulo para a variância.

anova(m1)

# Criar um fator que represente tratamento com a mesma variância. Nesse
# caso serão formados dois grupos.

dput(levels(dados$Cor))
dados$grupo <- factor(ifelse(dados$Cor %in% c("ambar_claro", "branco"),
                             "alta", "baixa"))

xyplot(Fenois ~ Cor, data = dados, groups = grupo, auto.key = TRUE)

m2 <- gls(Fenois ~ Cor, data = dados,
          weights = varIdent(form = ~1 | grupo),
          method = "ML")

anova(m2, m0)
anova(m2)

summary(m2)

​À
​ disposição
Walmes.​
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20161105/84913c2b/attachment.html>