[R-br] dose-resposta

Cesar Rabak cesar.rabak em gmail.com
Quarta Dezembro 28 18:53:47 BRST 2016


Marcelo,

Primeiro a sua observação mais "preocupante": os pacotes da Task View terem
sido criticados pelos membros da banca. . .

Como não tenho ideia de que pacote nem qual crítica, fico apenas ressabiado
que essas críticas não sejam mais divulgadas, nem que seja para a gente não
recair no(s) erro(s) que eles possam ter descoberto...

Como o Walmes explicou (inclusive em post posterior a minha observação) o
código teria que ser outro, apropriado para a ANCOVA, que do ponto de vista
de sintaxe é similar ao da ANOVA desde que a variável contínua seja mantida
dessa forma no dataframe.

HTH
--
Cesar Rabak


2016-12-22 18:23 GMT-02:00 Marcelo Laia <marcelolaia em gmail.com>:

> Prezado Cesar,
>
> Neste caso, ANCOVA, como ficaria o código que o Walmes gentilmente
> publicizou?
> Tenho trabalhado com dose resposta e os pacotes do Task View foram
> criticados
> por membros de bancas. Achei a solução do Walmes muito interessante. Por
> isso,
> gostaria de testar a sua abordagem, também!
>
> Marcelo
>
> On 22/12/16 at 01:21, Cesar Rabak via R-br wrote:
> > Walmes,
> >
> > Apenas um observação: não seria mais interessante aproveitar mais a
> > informação da dose se a tratássemos como variável contínua fazendo uma
> > ANCOVA?
> >
> > Minha interpretação apressada do seu código me faz crer que a dose virou
> > fator também. . .
> >
> > []s
> > --
> > Cesar Rabak
> >
> >
> > On Thu, Dec 22, 2016 at 9:46 AM, Walmes Zeviani <walmeszeviani em gmail.com
> >
> > wrote:
> >
> > > Ajuste um modelo com as três espécies como se fosse um fatorial. Para
> > > obter as DLs é mais fácil usar a parametrização de "estimativas
> separadas
> > > por nível". Veja o código abaixo.
> > >
> > > rm(list = ls())
> > > da <- data.frame(
> > >     dose = rep(c(0, 0.15625, 0.3125, 0.625, 1.25, 2.5, 5, 10), each =
> 4),
> > >     n = rep(10, 32),
> > >     m1 = c(1, 3, 4, 0, 5, 2, 5, 5, 4, 4, 2, 2, 0, 3, 3, 5, 10,
> > >            10, 7, 0, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10,
> > >            10),
> > >     m2 = c(0, 3, 4, 4, 1, 0, 4, 0, 0, 1, 0, 0, 3, 2, 0, 3, 3, 2,
> > >            0, 3, 3, 3, 3, 6, 9, 4, 7, 2, 10, 10, 10, 10),
> > >     m3 = c(4, 2, 2, 3, 4, 6, 8, 6, 4, 6, 6, 5, 5, 9, 5, 8, 10,
> > >             8, 5, 5, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10,
> > >             10))
> > > da$Dose <- as.factor(da$dose)
> > >
> > > da <- reshape2::melt(da,
> > >                      id.vars = c("dose", "Dose", "n"),
> > >                      variable.name = "esp",
> > >                      value.name = "mort")
> > > str(da)
> > >
> > > library(lattice)
> > >
> > > bwplot(mort/n ~ Dose | esp,
> > >        data = da,
> > >        pch = "|",
> > >        as.table = TRUE,
> > >        layout = c(NA, 1))
> > >
> > > xyplot(mort/n ~ Dose | esp,
> > >        data = da,
> > >        type = c("p", "a"),
> > >        as.table = TRUE,
> > >        layout = c(NA, 1))
> > >
> > > m0 <- glm(cbind(mort, n - mort) ~ esp * dose,
> > >           data = da,
> > >           family = quasibinomial)
> > >
> > > par(mfrow = c(2, 2))
> > > plot(m0)
> > > layout(1)
> > >
> > > summary(m0)
> > >
> > > anova(m0, test = "F")
> > >
> > > # Em um modelo binomial do tipo ~ \beta_0 + \beta_1 * x, a DL_{50} é
> > > # -\beta_0/\beta_1.
> > >
> > > # Ajuste do modelo com estimativas separadas para cada nível.
> > > m1 <- glm(cbind(mort, n - mort) ~ 0 + esp/dose,
> > >           data = da,
> > >           family = quasibinomial)
> > >
> > > # Deviances iguais porque são modelos iguais.
> > > deviance(m0)
> > > deviance(m1)
> > >
> > > summary(m1)
> > >
> > > beta <- matrix(coef(m1), ncol = 2)
> > > dl <- -beta[, 1]/beta[, 2]
> > >
> > > pred <- with(da,
> > >              expand.grid(esp = levels(esp),
> > >                          dose = seq(min(dose), max(dose),
> > >                                     length.out = 30)))
> > >
> > > pred$p <- predict(m1, newdata = pred, type = "response")
> > >
> > > library(latticeExtra)
> > >
> > > xyplot(mort/n ~ dose | esp,
> > >        data = da,
> > >        as.table = TRUE,
> > >        layout = c(NA, 1)) +
> > >     as.layer(xyplot(p ~ dose | esp,
> > >                     data = pred,
> > >                     type = "l",
> > >                     col = 2,
> > >                     lwd = 2)) +
> > >     layer(panel.abline(v = dl[which.packet()],
> > >                        h = 0.5,
> > >                        lty = 2))
> > >
> > > À disposição.
> > > Walmes.
> > > ​
> > >
>
> > _______________________________________________
> > R-br mailing list
> > R-br em listas.c3sl.ufpr.br
> > https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
> > Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forne?a
> c?digo m?nimo reproduz?vel.
>
>
>
> --
> Marcelo
>
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20161228/52388926/attachment.html>


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