Efeito da densidade de nematóides em cultivares de aveia
Santino Aleandro da Silva
Walmes Marques Zeviani
Experimento de responsabilidade de Patricia Meiriele Marini, aluna do curso de mestrado do IAPAR, que avaliou o dano causado pelo nematóide Meloidogyne incognita em cultivares de aveia (Avena sativa), a tendência de comportamento populacional do nematóide (fator de reprodução e número de nematóides por grama de raiz) e o nível populacional a partir do qual este nematóide causa dano acima de 10% ao desenvolvimento da planta em relação ao nível 0 (testemunha). O experimento foi conduzido em casa de vegetação no Instituto Agranômico do Paraná (IAPAR), sob orientação de Andressa C. Z. Machado. O delineamento foi inteiramente casualizado em arranjo experimental do tipo fatorial duplo (3x10) com 4 repetições por cela. Foram utilizados vasos com 3000 cm³ de solo, sendo conduzida uma planta de aveia por vaso. O significado das variáveis está descrito abaixo:
cultivar: fator de 3 níveis, cultivares de aveia;nivel: numérico de 10 níveis, número de espécimes de nematóide/cm³de solo;mfr: resposta numérica, massa fresca de raiz (g);mfpa: resposta numérica, massa fresca de parte aérea (g);mspa: resposta numérica, massa seca de parte aérea (g);fr: resposta numérica, fator de reprodução (população final/população inicial, exceto no nível 0 - sem nematóide);nema: resposta numérica, nematóide/grama de raiz (exceto no nível 0 - sem nematóide).
##-----------------------------------------------------------------------------
## Definições da sessão.
## Lista de pacotes a serem usados na sessão.
pkg <- c("lattice", "latticeExtra", "gridExtra", "nlme",
"doBy", "multcomp", "plyr", "reshape", "splines",
"wzRfun")
sapply(pkg, require, character.only=TRUE)## lattice latticeExtra gridExtra nlme doBy multcomp
## TRUE TRUE TRUE TRUE TRUE TRUE
## plyr reshape splines wzRfun
## TRUE TRUE TRUE TRUE## Instruções para instalação do pacote wzRfun em:
## https://github.com/walmes/wzRfun
sessionInfo()## R version 3.1.1 (2014-07-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=pt_BR.UTF-8
## [4] LC_COLLATE=en_US.UTF-8 LC_MONETARY=pt_BR.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=pt_BR.UTF-8 LC_NAME=C LC_ADDRESS=C
## [10] LC_TELEPHONE=C LC_MEASUREMENT=pt_BR.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] splines grid stats graphics grDevices utils datasets base
##
## other attached packages:
## [1] wzRfun_0.2 reshape_0.8.5 plyr_1.8.1 multcomp_1.3-6
## [5] TH.data_1.0-3 mvtnorm_0.9-9997 doBy_4.5-10 MASS_7.3-34
## [9] survival_2.37-7 nlme_3.1-117 gridExtra_0.9.1 latticeExtra_0.6-26
## [13] RColorBrewer_1.0-5 lattice_0.20-29 rmarkdown_0.2.68 knitr_1.6
##
## loaded via a namespace (and not attached):
## [1] digest_0.6.4 evaluate_0.5.5 formatR_1.0 htmltools_0.2.6 lme4_1.1-7
## [6] Matrix_1.1-4 methods_3.1.1 minqa_1.2.3 nloptr_1.0.4 Rcpp_0.11.0
## [11] sandwich_2.3-0 stringr_0.6.2 tools_3.1.1 yaml_2.1.13 zoo_1.7-10## Sem cores nos gráficos da lattice.
trellis.device(color=FALSE)
##-----------------------------------------------------------------------------
## Ler.
## list.files(pattern="*.txt")
da <- read.table("http://www.leg.ufpr.br/~walmes/data/aveia_nematoide.txt",
header=TRUE, sep="\t")
str(da)## 'data.frame': 120 obs. of 7 variables:
## $ cultivar: Factor w/ 3 levels "Afrodite","Slava",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ nivel : num 0 0 0 0 0.0625 0.0625 0.0625 0.0625 0.125 0.125 ...
## $ mfr : num 24.7 32 18.5 28.8 18.8 ...
## $ mfpa : num 42.8 51.5 36.4 60.1 24.3 ...
## $ mspa : num 5.24 6.08 4.3 7.7 3.94 4.64 4.61 3.52 4.29 4.87 ...
## $ fr : num NA NA NA NA 1.2 0.8 0.86 0.43 0 0 ...
## $ nema : int NA NA NA NA 12 8 7 4 0 0 ...##-----------------------------------------------------------------------------
## Tabela de frequência dos registros.
xtabs(~nivel+cultivar, da)## cultivar
## nivel Afrodite Slava Torena
## 0 4 4 4
## 0.0625 4 4 4
## 0.125 4 4 4
## 0.25 4 4 4
## 0.5 4 4 4
## 1 4 4 4
## 2 4 4 4
## 4 4 4 4
## 8 4 4 4
## 16 4 4 4##-----------------------------------------------------------------------------
## Ver.
db <- melt(data=da, id.vars=c("cultivar","nivel"))
str(db)## 'data.frame': 600 obs. of 4 variables:
## $ cultivar: Factor w/ 3 levels "Afrodite","Slava",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ nivel : num 0 0 0 0 0.0625 0.0625 0.0625 0.0625 0.125 0.125 ...
## $ variable: Factor w/ 5 levels "mfr","mfpa","mspa",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ value : num 24.7 32 18.5 28.8 18.8 ...## Transformação potência não elimina os zeros.
xyplot(value~nivel^0.25|variable, groups=cultivar, data=db,
type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
xlab="Densidade", auto.key=TRUE)
xyplot(value~nivel|variable, groups=cultivar, data=db,
type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
xlab="Densidade", auto.key=TRUE)
## ** Com o log() os zeros desaparacem.
xyplot(value~log2(nivel)|variable, groups=cultivar, data=db,
type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
xlab="Densidade", auto.key=TRUE)
##-----------------------------------------------------------------------------
## Criar uma indicadora da testemunha que é o nível zero.
## Pelos gráficos tem-se a impressão de que na dose zero foge à
## tendência de média considerando os demais níveis. Para isso, será
## criada uma nova variável (dummy) indicadora de ser nível 0.
da$zero <- as.integer(da$nivel==0)
## Criar uma versão fator do nível/densidade de nematóides.
da$niv <- factor(da$nivel)
levels(da$niv)## [1] "0" "0.0625" "0.125" "0.25" "0.5" "1" "2" "4" "8"
## [10] "16"## Representar nível na potência 0.25 e no log2().
da$np <- da$nivel^0.25
da$nl <- log2(da$nivel)
da$nl[!is.finite(da$nl)] <- -5Análise da MFR
##-----------------------------------------------------------------------------
## Ajuste.
## m0 <- lm(mfr~cultivar*niv, data=subset(da, nivel>0))
m0 <- lm(mfr~cultivar*niv, data=da)
par(mfrow=c(2,2)); plot(m0); layout(1)
MASS::boxcox(m0)
## Transformação Box-Cox aponta que, mesmo os resíduos não apresentando
## padrões fortes de fuga aos pressupostos, uma transformação raíz
## quadrada da resposta pode ser considerada.
## xyplot(sqrt(mfr)~np, groups=cultivar, data=da,
## type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
## xlab="Densidade", auto.key=TRUE)
## xyplot(sqrt(mfr)~nl, groups=cultivar, data=da,
## type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
## xlab="Densidade", auto.key=TRUE)
m0 <- update(m0, sqrt(mfr)~.)
par(mfrow=c(2,2)); plot(m0); layout(1)
## Estabilizou/anulou a relação média-variância.
##-----------------------------------------------------------------------------
## Quadro de anova.
anova(m0)## Analysis of Variance Table
##
## Response: sqrt(mfr)
## Df Sum Sq Mean Sq F value Pr(>F)
## cultivar 2 32.7 16.36 57.76 < 2e-16 ***
## niv 9 30.2 3.36 11.86 4.1e-12 ***
## cultivar:niv 18 10.5 0.58 2.06 0.014 *
## Residuals 90 25.5 0.28
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##-----------------------------------------------------------------------------
## Ajuste de um modelo reduzido considerando efeito de nível (^0.25) com
## um spline.
##
## m1 <- update(m0, .~cultivar*ns(x=nl, df=5))
## anova(m1, m0)
## str(da)
##
## pred <- with(da,
## expand.grid(
## cultivar=levels(cultivar),
## nl=seq(min(nl, na.rm=TRUE), max(nl, na.rm=TRUE),
## length.out=50))
## )
## pred <- cbind(pred,
## predict(m1, newdata=pred, interval="confidence"))
##
## xyplot(sqrt(mfr)~nl, groups=cultivar, data=da,
## as.table=TRUE, scales=list(y="free"),
## xlab="Densidade", auto.key=TRUE)+
## as.layer(
## xyplot(fit~nl, groups=cultivar, data=pred, type="l")
## )
## Não se obteve um ajuste consideravel usando polinômio ou
## splines. Qualquer uma dessas abordagens exigiu um alto grau da função
## para que não apresenta-se falta de ajuste. Sendo assim, optou-se por
## representar o efeito de nível por categorias.
## Causas para o efeito de nivel ser de alto grau:
## 1) de fato ser de alto grau, embora isso seja pouco sustentável.
## 2) existir uma causa de variação somada ao efeito em cada nível de
## densidade de nematóides.
##-----------------------------------------------------------------------------
## Estimativa das médias ajustadas.
X <- LSmatrix(m0, effect=c("cultivar","niv"))
grid <- equallevels(attr(X, "grid"), da)
Xm <- by(X, INDICES=grid$cultivar, FUN=as.matrix)
Xm <- lapply(Xm, "rownames<-", levels(da$niv))
L <- lapply(Xm, apmc, model=m0, focus="niv", test="fdr")
L## $Afrodite
## niv estimate lwr upr cld
## 1 0 5.072 4.544 5.601 a
## 2 0.0625 4.497 3.968 5.026 ab
## 3 0.125 3.989 3.460 4.517 ab
## 4 0.25 4.046 3.517 4.575 ab
## 5 0.5 3.840 3.312 4.369 b
## 6 1 4.193 3.664 4.721 ab
## 7 2 4.710 4.181 5.238 ab
## 8 4 4.426 3.898 4.955 ab
## 9 8 3.679 3.150 4.207 b
## 10 16 3.963 3.435 4.492 ab
##
## $Slava
## niv estimate lwr upr cld
## 1 0 3.813 3.284 4.342 ab
## 2 0.0625 3.484 2.956 4.013 b
## 3 0.125 3.836 3.308 4.365 ab
## 4 0.25 4.056 3.527 4.584 ab
## 5 0.5 3.486 2.958 4.015 b
## 6 1 3.975 3.446 4.503 ab
## 7 2 4.176 3.648 4.705 ab
## 8 4 4.593 4.064 5.122 a
## 9 8 1.773 1.244 2.302 c
## 10 16 3.645 3.116 4.174 ab
##
## $Torena
## niv estimate lwr upr cld
## 1 0 3.573 3.045 4.102 ab
## 2 0.0625 2.443 1.915 2.972 cd
## 3 0.125 3.100 2.571 3.628 bc
## 4 0.25 3.750 3.221 4.278 b
## 5 0.5 2.801 2.273 3.330 bc
## 6 1 3.241 2.712 3.770 bc
## 7 2 3.361 2.832 3.889 bc
## 8 4 3.046 2.517 3.574 bc
## 9 8 1.603 1.075 2.132 d
## 10 16 2.743 2.214 3.272 acgrid <- ldply(L); names(grid)[1] <- "cultivar"
grid <- equallevels(grid, da)
## Passando o inverso da tranformação feita nos dados.
grid[,c("estimate","lwr","upr")] <-
sapply(grid[,c("estimate","lwr","upr")], function(x) x^2)
##-----------------------------------------------------------------------------
## Gráfico com resultados.
segplot(niv~lwr+upr|cultivar, centers=estimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab="Massa fresca de raízes (g)")
Análise da MFPA
##-----------------------------------------------------------------------------
## Ajuste.
m0 <- lm(mfpa~cultivar*niv, data=da)
par(mfrow=c(2,2)); plot(m0); layout(1)
MASS::boxcox(m0)
##-----------------------------------------------------------------------------
## Quadro de anova.
anova(m0)## Analysis of Variance Table
##
## Response: mfpa
## Df Sum Sq Mean Sq F value Pr(>F)
## cultivar 2 904 452 4.78 0.011 *
## niv 9 10496 1166 12.34 1.6e-12 ***
## cultivar:niv 18 1941 108 1.14 0.328
## Residuals 90 8504 94
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1m0 <- update(m0, .~cultivar+niv)
##-----------------------------------------------------------------------------
## Estimativa das médias ajustadas.
X <- LSmatrix(m0, effect=c("cultivar","niv"))
grid <- equallevels(attr(X, "grid"), da)
Xm <- by(X, INDICES=grid$cultivar, FUN=as.matrix)
Xm <- lapply(Xm, "rownames<-", levels(da$niv))
L <- lapply(Xm, apmc, model=m0, focus="niv", test="fdr")
grid <- ldply(L); names(grid)[1] <- "cultivar"
grid <- equallevels(grid, da)
##-----------------------------------------------------------------------------
## Gráfico com resultados.
segplot(niv~lwr+upr|cultivar, centers=estimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab="Massa fresca de parte aérea (g)")
##-----------------------------------------------------------------------------
## Teste para os efeitos principais.
Xm <- LSmatrix(m0, effect="cultivar")
rownames(Xm) <- levels(da$cultivar)
apmc(Xm, model=m0, focus="cultivar")## cultivar estimate lwr upr cld
## 1 Afrodite 35.58 32.50 38.66 a
## 2 Slava 32.51 29.43 35.59 ab
## 3 Torena 28.86 25.78 31.95 bXm <- LSmatrix(m0, effect="niv")
rownames(Xm) <- levels(da$niv)
apmc(Xm, model=m0, focus="niv", test="fdr")## niv estimate lwr upr cld
## 1 0 38.51 32.884 44.14 ac
## 2 0.0625 27.00 21.374 32.63 de
## 3 0.125 29.39 23.765 35.02 de
## 4 0.25 35.81 30.181 41.44 bcd
## 5 0.5 34.26 28.631 39.89 cd
## 6 1 43.94 38.314 49.57 ab
## 7 2 45.53 39.904 51.16 a
## 8 4 32.93 27.303 38.56 cd
## 9 8 12.53 6.901 18.16 f
## 10 16 23.27 17.647 28.90 eAnálise da MSPA
##-----------------------------------------------------------------------------
## Ajuste.
m0 <- lm(mspa~cultivar*niv, data=da)
par(mfrow=c(2,2)); plot(m0); layout(1)
MASS::boxcox(m0)
##-----------------------------------------------------------------------------
## Quadro de anova.
anova(m0)## Analysis of Variance Table
##
## Response: mspa
## Df Sum Sq Mean Sq F value Pr(>F)
## cultivar 2 24.6 12.28 10.5 7.9e-05 ***
## niv 9 157.4 17.49 15.0 1.4e-14 ***
## cultivar:niv 18 27.4 1.52 1.3 0.21
## Residuals 90 105.2 1.17
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1m0 <- update(m0, .~cultivar+niv)
##-----------------------------------------------------------------------------
## Estimativa das médias ajustadas.
X <- LSmatrix(m0, effect=c("cultivar","niv"))
grid <- equallevels(attr(X, "grid"), da)
Xm <- by(X, INDICES=grid$cultivar, FUN=as.matrix)
Xm <- lapply(Xm, "rownames<-", levels(da$niv))
L <- lapply(Xm, apmc, model=m0, focus="niv", test="fdr")
grid <- ldply(L); names(grid)[1] <- "cultivar"
grid <- equallevels(grid, da)
##-----------------------------------------------------------------------------
## Gráfico com resultados.
segplot(niv~lwr+upr|cultivar, centers=estimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab="Massa seca de parte aérea (g)")
##-----------------------------------------------------------------------------
## Teste para os efeitos principais.
Xm <- LSmatrix(m0, effect="cultivar")
rownames(Xm) <- levels(da$cultivar)
apmc(Xm, model=m0, focus="cultivar")## cultivar estimate lwr upr cld
## 1 Afrodite 4.884 4.536 5.231 a
## 2 Slava 4.381 4.033 4.728 a
## 3 Torena 3.777 3.430 4.124 bXm <- LSmatrix(m0, effect="niv")
rownames(Xm) <- levels(da$niv)
apmc(Xm, model=m0, focus="niv", test="fdr")## niv estimate lwr upr cld
## 1 0 5.024 4.390 5.658 bc
## 2 0.0625 3.874 3.240 4.508 d
## 3 0.125 3.550 2.916 4.184 d
## 4 0.25 5.085 4.451 5.719 ac
## 5 0.5 4.123 3.488 4.757 cd
## 6 1 5.248 4.613 5.882 ab
## 7 2 6.047 5.413 6.682 a
## 8 4 5.194 4.560 5.828 ab
## 9 8 1.916 1.282 2.550 e
## 10 16 3.411 2.777 4.045 dAnálise da FR
##-----------------------------------------------------------------------------
## Ajuste.
## xyplot((fr^0.3)~log2(nivel), groups=cultivar, data=da,
## type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
## xlab="Densidade", auto.key=TRUE)
## m0 <- lm(c(fr+0.01)~cultivar*niv, data=da)
## par(mfrow=c(2,2)); plot(m0); layout(1)
## MASS::boxcox(m0)
m0 <- lm(fr~cultivar*niv, data=da)
par(mfrow=c(2,2)); plot(m0); layout(1)
## Transformação logtrans: y_trans = log(y+alpha).
MASS::logtrans(m0, alpha=seq(0.001, 0.1, by=0.01))
abline(v=0.001)
m0 <- update(m0, log(fr+0.01)~.)
par(mfrow=c(2,2)); plot(m0); layout(1)
##-----------------------------------------------------------------------------
## Quadro de anova.
anova(m0)## Analysis of Variance Table
##
## Response: log(fr + 0.01)
## Df Sum Sq Mean Sq F value Pr(>F)
## cultivar 2 377 188.5 285.03 < 2e-16 ***
## niv 8 142 17.7 26.82 < 2e-16 ***
## cultivar:niv 16 54 3.4 5.11 3.5e-07 ***
## Residuals 81 54 0.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##-----------------------------------------------------------------------------
## Estimativa das médias ajustadas.
X <- LSmatrix(m0, effect=c("cultivar","niv"))
grid <- equallevels(attr(X, "grid"), da)
Xm <- by(X, INDICES=grid$cultivar, FUN=as.matrix)
Xm <- lapply(Xm, "rownames<-", levels(da$niv)[-1])
L <- lapply(Xm, apmc, model=m0, focus="niv", test="fdr")
L## $Afrodite
## niv estimate lwr upr cld
## 1 0.0625 -0.2451 -1.054 0.5640 a
## 2 0.125 -4.6052 -5.414 -3.7961 f
## 3 0.25 -1.7184 -2.527 -0.9093 b
## 4 0.5 -2.7217 -3.531 -1.9126 bd
## 5 1 -2.2378 -3.047 -1.4287 bc
## 6 2 -3.3427 -4.152 -2.5336 cdf
## 7 4 -2.9776 -3.787 -2.1685 bde
## 8 8 -3.7093 -4.518 -2.9002 df
## 9 16 -4.1572 -4.966 -3.3481 ef
##
## $Slava
## niv estimate lwr upr cld
## 1 0.0625 2.00727 1.19816 2.8164 ab
## 2 0.125 1.90931 1.10021 2.7184 ab
## 3 0.25 3.15637 2.34727 3.9655 a
## 4 0.5 2.06144 1.25234 2.8705 ab
## 5 1 2.23848 1.42938 3.0476 ab
## 6 2 0.90513 0.09602 1.7142 bc
## 7 4 1.19884 0.38974 2.0079 bc
## 8 8 -0.09123 -0.90033 0.7179 cd
## 9 16 -0.84710 -1.65620 -0.0380 d
##
## $Torena
## niv estimate lwr upr cld
## 1 0.0625 0.7950 -0.01406 1.6041 bc
## 2 0.125 1.9256 1.11651 2.7347 ab
## 3 0.25 2.4412 1.63211 3.2503 a
## 4 0.5 2.7138 1.90473 3.5229 a
## 5 1 0.8533 0.04419 1.6624 bc
## 6 2 1.0988 0.28970 1.9079 bc
## 7 4 0.4880 -0.32111 1.2971 c
## 8 8 -1.2590 -2.06813 -0.4499 d
## 9 16 -2.4115 -3.22062 -1.6024 dgrid <- ldply(L); names(grid)[1] <- "cultivar"
grid <- equallevels(grid, da)
## Passando o inverso da tranformação feita nos dados.
ci <- sapply(grid[,c("estimate","lwr","upr")], function(x) exp(x)-0.01)
colnames(ci) <- paste0("p", colnames(ci))
grid <- cbind(grid, ci)
##-----------------------------------------------------------------------------
## Gráfico com resultados.
segplot(niv~plwr+pupr|cultivar, centers=pestimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab=expression(Fator~de~reprodução))
## Na escala em que foi realizada a análise.
segplot(niv~lwr+upr|cultivar, centers=estimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab=expression(Fator~de~reprodução~(y[t]==log(y+0.01))))
Análise da NEMA
##-----------------------------------------------------------------------------
## Ajuste.
## xyplot((nema^0.3)~log2(nivel), groups=cultivar, data=da,
## type=c("p", "smooth"), as.table=TRUE, scales=list(y="free"),
## xlab="Densidade", auto.key=TRUE)
## m0 <- lm(c(fr+0.01)~cultivar*niv, data=da)
## par(mfrow=c(2,2)); plot(m0); layout(1)
## MASS::boxcox(m0)
m0 <- lm(nema~cultivar*niv, data=da)
par(mfrow=c(2,2)); plot(m0); layout(1)
## Transformação logtrans: y_trans = log(y+alpha).
MASS::logtrans(m0, alpha=seq(0.001, 0.1, by=0.01))
abline(v=0.05)
m0 <- update(m0, log(nema+0.05)~.)
par(mfrow=c(2,2)); plot(m0); layout(1)
##-----------------------------------------------------------------------------
## Quadro de anova.
anova(m0)## Analysis of Variance Table
##
## Response: log(nema + 0.05)
## Df Sum Sq Mean Sq F value Pr(>F)
## cultivar 2 600 300.0 292.56 < 2e-16 ***
## niv 8 148 18.5 18.00 3.7e-15 ***
## cultivar:niv 16 69 4.3 4.21 7.5e-06 ***
## Residuals 81 83 1.0
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##-----------------------------------------------------------------------------
## Estimativa das médias ajustadas.
X <- LSmatrix(m0, effect=c("cultivar","niv"))
grid <- equallevels(attr(X, "grid"), da)
Xm <- by(X, INDICES=grid$cultivar, FUN=as.matrix)
Xm <- lapply(Xm, "rownames<-", levels(da$niv)[-1])
L <- lapply(Xm, apmc, model=m0, focus="niv", test="fdr")
L## $Afrodite
## niv estimate lwr upr cld
## 1 0.0625 1.982 0.97418 2.989 ab
## 2 0.125 -2.996 -4.00318 -1.988 c
## 3 0.25 1.980 0.97219 2.987 ab
## 4 0.5 1.066 0.05859 2.073 a
## 5 1 2.789 1.78163 3.797 ab
## 6 2 1.705 0.69790 2.713 ab
## 7 4 3.124 2.11669 4.132 b
## 8 8 3.183 2.17604 4.191 b
## 9 16 3.055 2.04750 4.062 b
##
## $Slava
## niv estimate lwr upr cld
## 1 0.0625 4.757 3.750 5.765 c
## 2 0.125 5.151 4.144 6.159 cd
## 3 0.25 7.003 5.996 8.010 b
## 4 0.5 6.907 5.899 7.914 b
## 5 1 7.491 6.483 8.498 ab
## 6 2 6.744 5.737 7.752 bd
## 7 4 7.548 6.541 8.556 ab
## 8 8 8.891 7.883 9.898 a
## 9 16 7.314 6.306 8.321 ab
##
## $Torena
## niv estimate lwr upr cld
## 1 0.0625 4.275 3.267 5.282 c
## 2 0.125 5.594 4.586 6.601 cd
## 3 0.25 6.420 5.413 7.428 bd
## 4 0.5 7.975 6.968 8.983 b
## 5 1 6.520 5.512 7.527 bd
## 6 2 7.372 6.364 8.379 ab
## 7 4 7.683 6.675 8.690 ab
## 8 8 7.858 6.851 8.866 ab
## 9 16 6.218 5.210 7.225 adgrid <- ldply(L); names(grid)[1] <- "cultivar"
grid <- equallevels(grid, da)
## Passando o inverso da tranformação feita nos dados.
ci <- sapply(grid[,c("estimate","lwr","upr")], function(x) exp(x)-0.01)
colnames(ci) <- paste0("p", colnames(ci))
grid <- cbind(grid, ci)
##-----------------------------------------------------------------------------
## Gráfico com resultados.
## Na escala da resposta.
segplot(niv~plwr+pupr|cultivar, centers=pestimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab=expression(Nematóides~por~cm^{-3}~de~solo))
## Na escala em que foi realizada a análise.
segplot(niv~lwr+upr|cultivar, centers=estimate,
data=grid, draw=FALSE, horizontal=FALSE,
scales=list(x=list(rot=90)), layout=c(3,1),
xlab=expression(Nematóides~por~cm^{-3}~de~solo),
ylab=expression(Nematóides~por~cm^{-3}~de~solo~(y[t]==log(y+0.05))))
##-----------------------------------------------------------------------------
## Análise considerando um glm Poisson.
## m0 <- glm(nema~cultivar*niv, data=da, family=quasipoisson)
## par(mfrow=c(2,2)); plot(m0); layout(1)
## summary(m0)
## anova(m0, test="F")
## Possivelmente essa contagem não seja de fato uma contagem mas um
## número obtido pela multiplicação por um fator de escala de uma
## contagem baseada em uma fração de cm³. Isso é comum na determinação
## de concentração de bactérias e demais organismos microscópicos.