[R-br] dose-resposta

[Cesar Rabak]

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.
> ​
>
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20161222/4588fc8e/attachment.html>