
Caros, Foi realizado um experimento no qual deseja-se saber qual composto extraído de uma determinada espécie vegetal é mais eficiente para matar uma determinado tipo de inseto. Foram utilizadas 32 placas cada uma contendo 10 insetos vivos. Em cada placa, foram realizadas 4 repetições das 8 diferentes doses do composto de cada espécie. O número de insetos mortos foi registrado após um determinado período de tempo. O objetivo é encontrar qual a dose, para cada um das três espécies, que mata 50% dos insetos. Alguma sugestão de como encontrar os modelos para as três espécies? dose = rep(c(0,0.15625,0.31250,0.62500,1.2500, 2.5000, 5.0000, 10.000), each=4) n_insetos = rep(10,32) m_especie1 = 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) m_especie2 = 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) m_especie3 = 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) pm_especie1 = m_especie1/n_insetos pm_especie2 = m_especie2/n_insetos pm_especie3 = m_especie3/n_insetos dose.fator = as.factor(dose) boxplot(pm_especie1~dose.fator, xlab= "Dose", ylab="Proporção de mortos espécie 1",cex.axis=.85) boxplot(pm_especie2~dose.fator, xlab= "Dose", ylab="Proporção de mortos espécie 2",cex.axis=.85) boxplot(pm_especie3~dose.fator, xlab= "Dose", ylab="Proporção de mortos espécie 3",cex.axis=.85) bartlett.test(pm_especie1~dose.fator) bartlett.test(pm_especie2~dose.fator) bartlett.test(pm_especie3~dose.fator) modelo1 = glm(cbind(good=m_especie1, bad=n_insetos-m_especie1)~dose,family="binomial", data=dados) par(mfrow=c(2,2)) plot(modelo1) shapiro.test(modelo1$residuals) modelo2 = glm(cbind(good=m_especie2, bad=n_insetos-m_especie2)~dose,family="binomial", data=dados) par(mfrow=c(2,2)) plot(modelo2) shapiro.test(modelo2$residuals) modelo3 = glm(cbind(good=m_especie3, bad=n_insetos-m_especie3)~dose,family="binomial", data=dados) par(mfrow=c(2,2)) plot(modelo2) shapiro.test(modelo3$residuals)

Veja se na Task View de farmacocinética algum dos pacotes indicados não lhe serve. HTH -- Cesar Rabak 2016-12-20 22:26 GMT-02:00 Maurício Lordêlo via R-br < r-br@listas.c3sl.ufpr.br>:
Caros, Foi realizado um experimento no qual deseja-se saber qual composto extraído de uma determinada espécie vegetal é mais eficiente para matar uma determinado tipo de inseto. Foram utilizadas 32 placas cada uma contendo 10 insetos vivos. Em cada placa, foram realizadas 4 repetições das 8 diferentes doses do composto de cada espécie. O número de insetos mortos foi registrado após um determinado período de tempo. O objetivo é encontrar qual a dose, para cada um das três espécies, que mata 50% dos insetos. Alguma sugestão de como encontrar os modelos para as três espécies?
dose = rep(c(0,0.15625,0.31250,0.62500,1.2500, 2.5000, 5.0000, 10.000), each=4) n_insetos = rep(10,32) m_especie1 = 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) m_especie2 = 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) m_especie3 = 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) pm_especie1 = m_especie1/n_insetos pm_especie2 = m_especie2/n_insetos pm_especie3 = m_especie3/n_insetos
dose.fator = as.factor(dose)
boxplot(pm_especie1~dose.fator, xlab= "Dose", ylab="Proporção de mortos espécie 1",cex.axis=.85) boxplot(pm_especie2~dose.fator, xlab= "Dose", ylab="Proporção de mortos espécie 2",cex.axis=.85) boxplot(pm_especie3~dose.fator, xlab= "Dose", ylab="Proporção de mortos espécie 3",cex.axis=.85)
bartlett.test(pm_especie1~dose.fator) bartlett.test(pm_especie2~dose.fator) bartlett.test(pm_especie3~dose.fator)
modelo1 = glm(cbind(good=m_especie1, bad=n_insetos-m_especie1)~dose,family="binomial", data=dados) par(mfrow=c(2,2)) plot(modelo1) shapiro.test(modelo1$residuals)
modelo2 = glm(cbind(good=m_especie2, bad=n_insetos-m_especie2)~dose,family="binomial", data=dados) par(mfrow=c(2,2)) plot(modelo2) shapiro.test(modelo2$residuals)
modelo3 = glm(cbind(good=m_especie3, bad=n_insetos-m_especie3)~dose,family="binomial", data=dados) par(mfrow=c(2,2)) plot(modelo2) shapiro.test(modelo3$residuals)
_______________________________________________ R-br mailing list R-br@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.

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.

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@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.

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@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@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

Cesar, Na análise existem duas variáveis: `Dose` (fator) e `dose` (numérica). Os gráficos foram feitos com a primeira (Dose) mantendo o que o Maurício fez mas o ajuste foi com a `dose` numérica pois só assim é possível determinar a DL_50. Se permite uma observação, posso estar errado inclusive e não me recordo de uma definição formal, mas acredito que na ANCOVA é quando se utiliza uma variável contínua (fazendo uma regressão) mas ela entra como "bloco", ou seja, para explicar variação nos dados porém sem haver interesse do pesquisador em explicar esse efeito. Esse é o caso do peso inicial de animais e idade em experimentos de nutrição (aves, suínos, bovinos, etc) e o caso de plantas por parcela (stand) em experimentos agronômicos no campo. No caso de usar um fator contínuo, como é o caso desse experimento do Maurício, não sei se é aplicável o nome ANCOVA porque a dose é um fator de interesse. Caso eu não tenha entendido a sua observação, deixe-me saber. À disposição. Walmes.

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@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@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@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

Muito obrigado Walmes pelo script enviado. Aos demais também pelos comentários e sugestões. Gostaria de fazer somente mais alguns questionamentos: 1. No gráfico gerado pela função "bwplot", como modificar na barra superior onde aparece o nome das espécies (m1, m2 e m3) para que estes nomes fiquem em itálico? bwplot( 100*(mort/n) ~ Dose | esp, data = da, pch = "|", as.table = TRUE,layout = c(NA, 1), ylab = "percentual de mortos", xlab=expression(paste("Concentração","(",~mu,"l","/ml"," ",")")) ) 2. Os pesquisadores da área de Agronomia e afins tem uma certa "idolatria" pelo R2. Acredito que não há necessidade deste coeficiente no caso do ajuste de um MLG. Como explicar isso para eles? 3. No caso de um ajuste usando a função "aov" ou "lm", a tabela da ANOVA é exposta da maneira como os pesquisadores da área de Agronomia conseguem entender. No caso do MLG há uma forma de saída destes resultados semelhante a isso? Maurício Em 28 de dezembro de 2016 17:53, Cesar Rabak via R-br < r-br@listas.c3sl.ufpr.br> escreveu:
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@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@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@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
_______________________________________________ R-br mailing list R-br@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.

Para o texto em itálico, considere o CMR abaixo. library(lattice) # Nomes das faixas em itálico. xyplot(Sepal.Length ~ Petal.Length | Species, data = iris, strip = strip.custom(par.strip.text = list(font = 3))) Sobre R², se realmente você quer apresentar um, pode considerar o quadrado da correlação entre observado e predito que justamente é o R² em modelos lineares. Esse não é único R² usado para GLM, tem uns como porcentagem de deviance explicada, enfim. # Correlação entre observado e predito ao quadrado. cor(iris$Sepal.Length, fitted(m0))^2 summary(m0)$r.squared À disposição. Walmes.

Obrigado Walmes! Problema resolvido. Maurício Em 10 de janeiro de 2017 11:00, Walmes Zeviani <walmeszeviani@gmail.com> escreveu:
Para o texto em itálico, considere o CMR abaixo.
library(lattice)
# Nomes das faixas em itálico. xyplot(Sepal.Length ~ Petal.Length | Species, data = iris, strip = strip.custom(par.strip.text = list(font = 3)))
Sobre R², se realmente você quer apresentar um, pode considerar o quadrado da correlação entre observado e predito que justamente é o R² em modelos lineares. Esse não é único R² usado para GLM, tem uns como porcentagem de deviance explicada, enfim.
# Correlação entre observado e predito ao quadrado. cor(iris$Sepal.Length, fitted(m0))^2 summary(m0)$r.squared
À disposição. Walmes.
participantes (4)
-
Cesar Rabak
-
Marcelo Laia
-
Maurício Lordêlo
-
Walmes Zeviani