Desdobramento de um fatorial triplo
Caros, Peço auxilio para um problema de um fatorial em que há um problema em atender aos pressupostos. Agradeço alguma sugestão. Maurício m_seca_raiz = c(0.000000000, 0.000000000, 0.000460000, 0.000000000, 0.000000000, 0.000140000, 0.000000000, 0.000300000, 0.000150000, 0.000000000, 0.000120000, 0.000500000, 0.000560000, 0.000300000, 0.000000000, 0.000160000, 0.000160000, 0.000600000, 0.000200000, 0.000200000, 0.000440000, 0.000000000, 0.000080000, 0.000240000, 0.000380000, 0.000080000, 0.000000000, 0.000250000, 0.000960000, 0.000800000, 0.000320000, 0.000760000, 0.000920000, 0.000366667, 0.000720000, 0.000540000, 0.000440000, 0.000620000, 0.000100000, 0.000600000, 0.000200000, 0.000280000, 0.000000000, 0.000520000, 0.000480000, 0.000000000, 0.000280000, 0.000280000, 0.000800000, 0.000280000, 0.000000000, 0.000680000, 0.000700000, 0.000380000, 0.000000000, 0.000360000, 0.000300000, 0.000340000, 0.000280000, 0.000360000, 0.000340000, 0.001780000, 0.000250000, 0.000260000, 0.000300000, 0.000280000, 0.000160000, 0.000120000, 0.000125000, 0.000080000, 0.000220000, 0.000320000) meio = factor((rep(c("m2","m3"),each=12, times=3))) sacarose = factor(rep(c("7,5","15,0","30,0"), each=24)) carvão = factor((rep(c("c1","c2"),each=6, times=6))) tratamento = meio:sacarose:carvão head(data.frame(m_seca_raiz, meio, sacarose, carvão, tratamento),18) hist(m_seca_raiz) boxplot(m_seca_raiz ~ tratamento,ylab="Massa seca da raiz", xlab="Tratamento", cex.axis=0.85 ) max(m_seca_raiz) modelo_m_seca_raiz = aov(m_seca_raiz ~ tratamento) summary(modelo_m_seca_raiz) hist(modelo_m_seca_raiz$residuals) shapiro.test(modelo_m_seca_raiz$residuals) require(car) leveneTest(modelo_m_seca_raiz$residuals ~ tratamento) anova(modelo_m_seca_raiz)$'Mean Sq'[2] res_pad_m_seca_raiz = modelo_m_seca_raiz$residuals/sqrt(anova(modelo_m_seca_raiz)$'Mean Sq'[2]) plot(modelo_m_seca_raiz$fitted.values, res_pad_m_seca_raiz, ylim = c(-5,5)) abline(h=c(-3,3),lty=3) #Elimando a observação 62 shapiro.test(modelo_m_seca_raiz$residuals[-62]) leveneTest(modelo_m_seca_raiz$residuals[-62] ~ tratamento[-62]) #Transformação Box-Cox require(forecast) #require(fpp) BoxCox.lambda(m_seca_raiz+0.5, method = "guerrero", lower=-3,upper = 3) BoxCox.lambda(m_seca_raiz+0.5, method = "loglik", lower=-3,upper = 3) #lambda = BoxCox.lambda(m_seca_raiz,lower=-3) #lambda m_seca_raizt = BoxCox(m_seca_raiz+0.5, lambda=-3) hist(m_seca_raizt) boxplot(m_seca_raizt ~ tratamento, ylab="m_seca_raiz (transformada)", xlab="Tratamentos") modelo_m_seca_raizt = aov(m_seca_raizt ~ tratamento) summary(modelo_m_seca_raizt) hist(modelo_m_seca_raizt$residuals) plot(modelo_m_seca_raizt$residuals) plot(modelo_m_seca_raizt$fitted.values, modelo_m_seca_raizt$residuals) shapiro.test(modelo_m_seca_raizt$residuals) leveneTest(modelo_m_seca_raizt$residuals ~tratamento) anova(modelo_m_seca_raizt)$'Mean Sq'[2] res_pad_m_seca_raizt = modelo_m_seca_raizt$residuals/sqrt(anova(modelo_m_seca_raizt)$'Mean Sq'[2]) plot(modelo_m_seca_raizt$fitted.values, res_pad_m_seca_raizt, ylim = c(-6,6)) abline(h=c(-3,3),lty=3)
Sem maiores esclarecimentos sobre os dados, análise, etc, o problema parece ser apenas uma observação com valor extremo se comparado com as demais. Ela salta aos olhos já na análise exploratória. Daí a importância de uma boa análise exploratória. Ainda que seja composta de gráficos simples, negar a sua utilidade e partir para a análise sem ver de fato os dados pode gerar preocupações e escrita de código desnecessários. Aproveitando o ensejo, eu tenho uma posição de reserva com relação a aplicação de sequencias de testes para a avaliação dos pressupostos. Eu tenho preferência pela análise visual porque ela dá um panorama amplo e detalhado da condição dos pressupostos com sugestões de curso de ação (é um outlier? é relação média variância? é assimetria?). Eu considero a análise de dados uma tarefa artesanal, por isso não se pode deixar de aplicar no olhar clínico do analista. Por outro lado, se fosse para automatizar a análise de dados como análise em fluxo contínuo de várias respostas, ou análise diária, grande volume de dados, pouco tempo disponível, extrema necessidade de eliminar subjetividade (como se estivesse diante de uma corte), etc etc, aí eu usaria tais testes como pré-teste para que não seja necessário correr a todo momento a análise diagnóstico pente fino. da <- data.frame(m_seca_raiz, meio, sacarose, carvão, tratamento) # Ocorrência dos pontos experimentais. ftable(xtabs(~meio + sacarose + carvão, data = da)) library(lattice) # Inspeção dos dados antes da análise. xyplot(m_seca_raiz ~ sacarose | carvão, groups = meio, auto.key = TRUE, data = da, type = c("p", "a")) # ATTENTION: Análise exploratória revela um ponto suspeito na cela # c1:30,0:m3. Existe conhecimento de algo que tenha provocado um valor # maior? E.g. parcela comprometida, erro de dosagem, erro de mensuração, # etc? # Valor extremo será eliminado. Isso não complica a análise. m0 <- lm(m_seca_raiz ~ sacarose * carvão * meio, data = da, subset = m_seca_raiz < 0.0015) par(mfrow = c(2, 2)) plot(m0) layout(1) # NOTE: Pressupostos OK! Apenas uma muito sutil indicação de relação # média variância mas que não é alarmante. # Verifica sugestão de transformação. range(m0$model[, 1]) MASS::boxcox(update(m0, . + 0.00001 ~ .)) abline(v = 0.5, col = "red") # Transformar por raíz deve melhorar quase nada, todavia, vamos # experimentar. m1 <- lm((m_seca_raiz + 0.00001)^0.5 ~ sacarose * carvão * meio, data = da, subset = m_seca_raiz < 0.0015) par(mfrow = c(2, 2)) plot(m1) layout(1) # NOTE: leve melhoria na relação média-variância. # Médias marginais (erros padrões diferentes pelo desbalanceamento). emm <- emmeans::emmeans(m1, specs = "sacarose") emm # Contrastes entre média marginais. emmeans::contrast(emm, method = "pairwise") À disposição. Walmes.
participantes (2)
-
Maurício Lordêlo -
Walmes Zeviani