[R-br] Comparações múltiplas usando MMC (mean–mean multiple comparisons) plots

ASANTOS alexandresantosbr em yahoo.com.br
Quinta Maio 10 10:22:43 BRT 2012


*Bom dia pessoal,

          Estou tendo problemas em visualizar os resultados de uma 
comparação múltipla usando GLM, com o uso MMC plots do pacote HH. Bom, 
tenho quatro grupos funcionais de cupins que são humívoros (humi), 
ceifadores(ceifa), intermediários (inter) e xilófagos (xilo) que são 
meus tratamentos e como variável resposta o número de ocorrência de 
espécies de cada grupo (acorr). Para coleta da variável resposta, 
fizemos 6 transectos que são as repetições. Para ver se existe diferença 
entre os grupos funcionais fiz:*

  #
mod2<-glm(ocorr~grupo, family="quasipoisson", data=trat1_90antes)## Um 
plot(mod2) prévio mostrou que quasipoisson foi a distribuição mais 
adequada aos dados que poisson
summary(mod2)
glm(formula = ocorr ~ grupo, family = "quasipoisson", data = trat1_90antes)
Deviance Residuals:
     Min       1Q   Median       3Q      Max
-1.9148  -0.4946   0.1213   0.4178   1.3813
Coefficients:
             Estimate Std. Error t value Pr(>|t|)
(Intercept)   1.5041     0.1624   9.260 1.13e-08 ***
grupohumi    -1.0986     0.3249  -3.382  0.00296 **
grupointer   -0.8979     0.3019  -2.974  0.00750 **
grupoxilo    -0.6568     0.2780  -2.363  0.02838 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for quasipoisson family taken to be 0.7123377)
     Null deviance: 28.663  on 23  degrees of freedom
Residual deviance: 16.889  on 20  degrees of freedom
AIC: NA
Number of Fisher Scoring iterations: 5
anova(mod2, test="Chi") ## Anova do modelo
Analysis of Deviance Table
Model: quasipoisson, link: log
Response: ocorr
Terms added sequentially (first to last)
       Df Deviance Resid. Df Resid. Dev  Pr(>Chi)
NULL                     23     28.663
grupo  3   11.774        20     16.889 0.0008834 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
plot(mod2)

##Teste de comparação múltipla
summary(glht(mod2, linfct=mcp(grupo="Tukey")))

          Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Fit: glm(formula = ocorr ~ grupo, family = "quasipoisson", data = 
trat1_90antes)

Linear Hypotheses:
                    Estimate Std. Error z value Pr(>|z|)
humi - ceifa == 0   -1.0986     0.3249  -3.382  0.00405 **
inter - ceifa == 0  -0.8979     0.3019  -2.974  0.01491 *
xilo - ceifa == 0   -0.6568     0.2780  -2.363  0.08255 .
inter - humi == 0    0.2007     0.3794   0.529  0.95122
xilo - humi == 0     0.4418     0.3606   1.225  0.60632
xilo - inter == 0    0.2412     0.3401   0.709  0.89176
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)

 > ##Resultados
 > tapply(trat1_90antes$ocorr, trat1_90antes$grupo,mean, na.rm = 
TRUE)##Media dos tratamentos
    ceifa     humi    inter     xilo
4.500000 1.500000 1.833333 2.333333
 > 
tapply(trat1_90antes$ocorr,trat1_90antes$grupo,sd,na.rm=TRUE)/sqrt(tapply(trat1_90antes$ocorr,trat1_90antes$grupo,length))###Erros 
padrão da media
     ceifa      humi     inter      xilo
0.6708204 0.5000000 0.5426274 0.3333333

*Aqui tornou-se difícil saber a diferença entre os grupos funcionais em 
termos de atribui a, b, c, as médias, então resolvi fazer uma análise 
gráfica usando o pacote HH:*

 > ##Visualizando as comparações
 > contrast.1<-
+ if.R(r=glht(mod2, linfct=mcp(grupo="Tukey")),
+      s=multicomp(mod2, plot=FALSE))
 > plot(contrast.1)


*Que me deu um gráfico pouco informativo, porém seguindo o help do mmc 
{HH}, encontrei que o gráfico de comparações deveria ser feito usando as 
funções glht.mmc e multicomp.mmc, que requerem o uso de aov(), então fiz:*

 >
 > ## Contraste por mmc
 > contrast.aov <- aov(ocorr ~ grupo, data=trat1_90antes)
 > contrast.mmc <-
+ if.R(r=glht.mmc(contrast.aov, linfct=mcp(grupo="Tukey")),
+      s=multicomp.mmc(contrast.aov, plot=FALSE))
 > plot(contrast.mmc)
 > #

*Pois bem, agora tenho a comparação entre os grupos, mas não utilizando 
a distribuição de erros quasipoisson e não sei o que fazer, pois o 
pacote HH não me permite utilizar os objetos criados usando glm(), para 
se fazer as comparações múltiplas. Alguém do grupo teria alguma sugestão,

Obrigado,

Alexandre*

-- 
Alexandre dos Santos
Engenheiro Florestal, Dr.
Universidade Federal de Lavras
Departamento de Entomologia
Laboratório de Entomologia Florestal
Caixa Postal 3037
37200-000 - Lavras/MG
Fone: +55 (35) 9223-0304

-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20120510/b63175ce/attachment.html>


Mais detalhes sobre a lista de discussão R-br