Atualiza conteúdo da inst/doc.

parent 4599baf6
Pipeline #6379 passed with stage
in 33 minutes and 37 seconds
## ----setup, include=FALSE-----------------------------------------
## ----setup, include=FALSE------------------------------------------------
source("_setup.R")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
library(MRDCr)
## -----------------------------------------------------------------
# help(soja)
ls("package:MRDCr")
......@@ -28,7 +26,7 @@ xyplot(nvag + ngra ~ K, groups = umid, outer = TRUE,
soja <- soja[-74, ]
## ---- message=FALSE-----------------------------------------------
## ---- message=FALSE------------------------------------------------------
# Considerar K categórico.
soja <- transform(soja, K = factor(K))
......@@ -96,7 +94,7 @@ lapply(K,
test = adjusted(type = "fdr"))
})
## ---- message=FALSE-----------------------------------------------
## ---- message=FALSE------------------------------------------------------
m1 <- glm(ngra ~ bloc + umid * K, data = soja, family = poisson)
par(mfrow = c(2, 2))
......@@ -154,7 +152,7 @@ lapply(K,
test = adjusted(type = "fdr"))
})
## ---- echo=FALSE, results="hold"----------------------------------
## ---- echo=FALSE, results="hold"-----------------------------------------
cat(format(Sys.time(),
format = "Atualizado em %d de %B de %Y.\n\n"))
sessionInfo()
......
This diff is collapsed.
## ----setup, include=FALSE-----------------------------------------
## ----setup, include=FALSE------------------------------------------------
source("_setup.R")
## ---- echo=FALSE, include=FALSE-----------------------------------
devtools::load_all()
## ---- results="hide", message=FALSE-------------------------------
## ---- results="hide", message=FALSE--------------------------------------
# Pacotes requeridos.
library("lmtest")
library("boot")
......@@ -14,16 +11,13 @@ library("RColorBrewer")
library("sandwich")
library("hnp")
library("knitr")
library("MRDCr")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
data(postura)
str(postura)
summary(postura)
# names(postura) <- c("npost", "trat", "linh")
# postura <- postura[, c(2, 3, 1)]
# use_data(postura, overwrite = TRUE)
# tab <- xtabs(~trat + npost, data = postura)
bwplot(npost ~ linh | trat,
data = postura,
main = "Mudanças de postura vs tratamento e linhagem",
......@@ -46,7 +40,7 @@ mdp <- aggregate(npost ~ trat + linh,
})
mdp
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Ajuste do modelo Poisson.
mP <- glm(npost ~ trat + linh,
data = postura,
......@@ -68,12 +62,12 @@ X2 <- sum(resid(mP, type = "pearson")^2)
phichap <- X2/df.residual(mP)
phichap
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Diagnóstico do modelo - gráficos padrão do R.
par(mfrow = c(2, 2))
plot(mP)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
envelope <- function(modelo) {
dados <- na.omit(modelo$data)
nsim <- 100
......@@ -104,12 +98,12 @@ envelope <- function(modelo) {
points(quantis, r1, pch = 16, cex = 0.75)
}
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Gráfico quantil-quantil com envelopes simulados.
layout(1)
envelope(mP)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Modelo quasi poisson (V(mu) = mu).
mQP0 <- glm(npost ~ trat + linh,
data = postura,
......@@ -142,7 +136,7 @@ qqnorm(resid(mQP1, type = "deviance"),
pch = 20, main = "QL", las = 1)
qqline(resid(mQP1, type = "deviance"))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Estimador sanduíche.
estrb <- coeftest(mP, vcov = sandwich)
estrb
......@@ -153,7 +147,7 @@ mB <- Boot(mP)
# Resultados obtidos via bootstrap.
summary(mB)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
erroz <- rbind(summary(mP)$coefficients[2,2:3],
summary(mQP0)$coefficients[2,2:3],
summary(mQP1)$coefficients[2,2:3],
......@@ -179,7 +173,7 @@ kable(quadres,
format = "markdown",
caption = "Comparativo dos modelos ajustados")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
postura.ex <- postura[-c(8, 18, 28), ]
# Ajustando o modelo Poisson sem as três observações.
......@@ -218,7 +212,7 @@ exp(coef(mPe)[2])
# O efeito de intervenção aumenta, e torna-se mais significativo,
# mediante exclusão dos outliers.
## ---- echo=FALSE, results="hold"----------------------------------
## ---- echo=FALSE, results="hold"-----------------------------------------
cat(format(Sys.time(),
format = "Atualizado em %d de %B de %Y.\n\n"))
sessionInfo()
......
This diff is collapsed.
## ----setup, include=FALSE-----------------------------------------
## ----setup, include=FALSE------------------------------------------------
source("_setup.R")
## ---- results = "hide", message = FALSE---------------------------
## ---- results = "hide", message = FALSE----------------------------------
# Pacotes necessários.
library(lattice)
library(MASS)
library(effects)
library(knitr)
library(MRDCr)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Dez primeiras linhas da base.
head(seguros, 10)
str(seguros)
seguros$lexpo <- log(seguros$expos)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Distribuição do números de sinistros.
tb <- table(seguros$nsinist)
tb
......@@ -58,7 +57,7 @@ tabidsex <- cbind(tabidsex, taxa)
kable(tabidsex, align = "c",
caption = "**Taxa de sinistros segundo sexo e idade**")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
seguros <- na.omit(seguros)
mP <- glm(nsinist ~ sexo + idade + I(idade^2) + valor +
offset(log(expos)),
......@@ -72,7 +71,7 @@ phichap <- X2/mP$df.residual
# Indicador de superdispersão.
phichap
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
envelope <- function(modelo) {
dados <- na.omit(modelo$data)
nsim <- 100
......@@ -103,7 +102,7 @@ envelope <- function(modelo) {
points(quantis, r1, pch = 16, cex = 0.75)
}
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Diagnóstico do modelo - gráficos.
par(mfrow = c(2, 2))
plot(mP)
......@@ -111,24 +110,24 @@ plot(mP)
par(mfrow = c(1, 1))
envelope(mP)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
mPo <- glm(nsinist ~ sexo + idade + I(idade^2) + valor +
log(expos),
data = seguros, family = poisson)
summary(mPo)
anova(mP, mPo, test = "Chisq")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
mBNo <- glm.nb(nsinist ~ sexo + idade + I(idade^2) + valor +
log(expos), data = seguros)
summary(mBNo)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Diagnóstico do modelo - gráficos.
par(mfrow = c(2, 2))
plot(mBNo)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
dadosnb3 <-
seguros[, c("sexo", "idade", "valor", "expos", "nsinist")]
dadosnb3$lexpo <- log(seguros$expos)
......@@ -167,11 +166,11 @@ envelope <- function(modelo) {
points(quantis, r1, pch = 16, cex = 0.75)
}
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
par(mfrow = c(1, 1))
envelope(mBNo)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
intervalos <- confint(mBNo)
estimat <- cbind(mBNo$coefficients, intervalos)
colnames(estimat)[1] <- "Estimativa pontual"
......@@ -181,9 +180,8 @@ kable(round(estimat, 5), align = "c",
caption = paste("**Estimativas pontuais e intervalos de",
"confiança - Modelo Binomial Negativo**"))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
efeitos <- allEffects(mBNo, given.values = c(lexpo = 0))
trellis.par.set(list(axis.text = list(cex = 1.2)))
plot(efeitos[[2]],
type = "response",
......@@ -203,7 +201,9 @@ plot(efeitos[[4]],
xlab = "Valor (x1000 reais)",
ylab = "Taxa de sinistros")
## -----------------------------------------------------------------
dev.off()
## ------------------------------------------------------------------------
# Poisson sem ajuste de covariáveis.
n <- nrow(seguros)
......@@ -252,13 +252,16 @@ rownames(matfreq) <- c("Amostra",
"BN não ajustada por covariáveis",
"BN ajustada por covariáveis")
## ---- results="markup"--------------------------------------------
## ---- results="markup"---------------------------------------------------
kable(matfreq, format = "markdown",
caption = paste("Frequências amostrais e frequências",
"ajustadas para o número de sinistros"))
## ---- echo=FALSE, results="hold"----------------------------------
## ---- echo=FALSE, results="hold"-----------------------------------------
cat(format(Sys.time(),
format = "Atualizado em %d de %B de %Y.\n\n"))
sessionInfo()
## ---- include=FALSE------------------------------------------------------
detach("package:effects")
This diff is collapsed.
## ----setup, include=FALSE-----------------------------------------
## ----setup, include=FALSE------------------------------------------------
source("_setup.R")
## ---- message=FALSE, error=FALSE, warning=FALSE-------------------
## ---- message=FALSE, error=FALSE, warning=FALSE--------------------------
# Definições da sessão.
# devtools::load_all("../")
library(lattice)
library(latticeExtra)
library(bbmle)
......@@ -14,7 +13,7 @@ library(doBy)
library(multcomp)
library(MRDCr)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
grid <- expand.grid(lambda = c(2, 8, 15),
alpha = c(-0.05, 0, 0.05))
y <- 0:35
......@@ -60,7 +59,7 @@ for (a in seq_along(alpha)) {
add = TRUE, xname = "lambda", col = col[a], lwd = 2)
}
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Gráfico do espaço paramétrico de lambda x alpha.
......@@ -87,7 +86,7 @@ levelplot(sum ~ lambda + alpha,
# Já que lambda * alpha > -1, então alpha = -1/lambda dá a fronteira.
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Função de log-Verossimilhança da Poisson Generalizada na
# parametrização de modelo de regressão.
MRDCr::llpgnz
......@@ -128,7 +127,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Carregando e explorando os dados.
......@@ -203,7 +202,7 @@ dev.off()
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Testes de hipótese.
......@@ -273,7 +272,7 @@ xyplot(fit ~ K | umid, data = pred,
desloc = 8 * scale(as.integer(pred$modelo), scale = FALSE),
panel = panel.cbH)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
xyplot(ngra ~ K | umid, data = soja, layout = c(NA, 1),
......@@ -338,7 +337,7 @@ dev.off()
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Teste de Wald para a interação.
a <- c(0, attr(model.matrix(m0), "assign"))
ai <- a == max(a)
......@@ -405,7 +404,7 @@ xyplot(fit ~ K | umid, data = pred,
desloc = 8 * scale(as.integer(pred$modelo), scale = FALSE),
panel = panel.cbH)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Número de grãos por vagem (offset).
......@@ -482,7 +481,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -556,7 +555,7 @@ xyplot(fit ~ K | umid, data = pred,
desloc = 8 * scale(as.integer(pred$modelo), scale = FALSE),
panel = panel.cbH)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Número de capulhos em função do nível de desfolha artificial e fase
# fenológica do algodoeiro.
......@@ -622,7 +621,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -689,12 +688,12 @@ p2 <- xyplot(fit ~ des | est, data = pred, groups = modelo,
panel = panel.superpose)
# p2
## ---- fig.width=7, fig.height=3.5---------------------------------
## ---- fig.width=7, fig.height=3.5----------------------------------------
update(p1, type = "p", layout = c(NA, 1),
key = key, spread = 0.07) +
as.layer(p2, under = TRUE)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
data(nematoide, package = "MRDCr")
......@@ -763,7 +762,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -884,7 +883,7 @@ qqmath(~values | ind, data = rp,
})
## ---- echo=FALSE, results="hold"----------------------------------
## ---- echo=FALSE, results="hold"-----------------------------------------
cat(format(Sys.time(),
format = "Atualizado em %d de %B de %Y.\n\n"))
sessionInfo()
......
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
## ----setup, include=FALSE-----------------------------------------
## ----setup, include=FALSE------------------------------------------------
source("_setup.R")
## ---- message=FALSE, error=FALSE, warning=FALSE-------------------
## ---- message=FALSE, error=FALSE, warning=FALSE--------------------------
# Definições da sessão.
library(lattice)
library(latticeExtra)
......@@ -14,7 +14,7 @@ library(doBy)
library(multcomp)
library(MRDCr)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Função densidade na parametrização original.
dgcnt
......@@ -48,7 +48,7 @@ useOuterStrips(xyplot(py ~ y | factor(lambda) + factor(alpha),
var.name = expression(alpha == ""),
sep = ""))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# A média da Gamma-Count.
......@@ -121,7 +121,7 @@ xyplot(m ~ lambda, groups = alpha, type = "l",
layer(panel.abline(a = -0.5, b = 1, lty = 2)) +
layer(panel.abline(a = 1, b = 1, lty = 2))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
h <- function(...) {
dgamma(...)/(1 - pgamma(...))
}
......@@ -156,7 +156,7 @@ legend("topright", legend = sprintf("%0.3f", shape),
col = col, lty = 1, lwd = 2, bty = "n",
title = expression(alpha))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Função de log-Verossimilhança da Poisson Generalizada na
# parametrização de modelo de regressão.
......@@ -199,7 +199,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Carregando e explorando os dados.
......@@ -276,7 +276,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -348,7 +348,7 @@ xyplot(fit ~ K | umid, data = pred,
desloc = 8 * scale(as.integer(pred$modelo), scale = FALSE),
panel = panel.cbH)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
xyplot(ngra ~ K | umid, data = soja, layout = c(NA, 1),
......@@ -430,7 +430,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -504,7 +504,7 @@ xyplot(fit ~ K | umid, data = pred,
desloc = 8 * scale(as.integer(pred$modelo), scale = FALSE),
panel = panel.cbH)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Número de grãos por vagem (offset).
......@@ -574,7 +574,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -658,7 +658,7 @@ xyplot(ngra/nvag ~ K | umid, data = soja, layout = c(NA, 1),
scale = FALSE),
panel = panel.cbH))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Número de capulhos em função do nível de desfolha artificial e fase
# fenológica do algodoeiro.
......@@ -729,7 +729,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -801,12 +801,12 @@ p2 <- xyplot(fit ~ des | est, data = pred, groups = modelo,
panel = panel.superpose)
# p2
## ---- fig.width=7, fig.height=3.5---------------------------------
## ---- fig.width=7, fig.height=3.5----------------------------------------
update(p1, type = "p", layout = c(NA, 1),
key = key, spread = 0.07) +
as.layer(p2, under = TRUE)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
data(nematoide, package = "MRDCr")
......@@ -878,7 +878,7 @@ corrplot.mixed(V, lower = "number", upper = "ellipse",
tl.cex = 0.8, col = brewer.pal(9, "Greys")[-(1:3)])
dev.off()
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
# Tamanho das covariâncias com \alpha.
each(sum, mean, max)(abs(V[1, -1]))
......@@ -933,7 +933,7 @@ xyplot(nema/off ~ cult, data = nematoide,
scale = FALSE),
panel = panel.cbH))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
#-----------------------------------------------------------------------
# Log-verossimilhança para o número de gols dos times em uma partida.
......@@ -1177,7 +1177,7 @@ colnames(dwl) <- sprintf(paste(levels(db$home)[i], collapse = " %s "),
c("=", ">", "<"))
t(dwl)
## ---- echo=FALSE, results="hold"----------------------------------
## ---- echo=FALSE, results="hold"-----------------------------------------
cat(format(Sys.time(),
format = "Atualizado em %d de %B de %Y.\n\n"))
sessionInfo()
......
This diff is collapsed.
## ----setup, include=FALSE-----------------------------------------
## ----setup, include=FALSE------------------------------------------------
source("_setup.R")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
library(MRDCr)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
data(peixe)
peixe$campista <- as.factor(peixe$campista)
levels(peixe$campista) <- c("Não", "Sim")
......@@ -13,7 +13,7 @@ str(peixe)
## help(peixe)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Estudo observacional
ftable(with(peixe, table(npessoas, ncriancas, campista)))
......@@ -55,7 +55,7 @@ xyplot(lnpeixes ~ ncriancas + npessoas,
type = c("p", "g", "spline"))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
##======================================================================
## Ajuste de modelos hurdle
......@@ -91,7 +91,7 @@ m1ZBN <- zeroinfl(f1, data = peixe, dist = "negbin")
m2ZBN <- zeroinfl(f2, data = peixe, dist = "negbin")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Via log-verossimilhança
cbind("Poisson" = sapply(list(m1P, m2P), logLik),
......@@ -103,7 +103,7 @@ cbind("Poisson" = sapply(list(m1P, m2P), logLik),
)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Testes de razão de verossimilhanças para o efeito de interação
anova(m1BN, m2BN)
......@@ -111,7 +111,7 @@ lmtest::lrtest(m1HBN, m2HBN)
lmtest::lrtest(m1ZBN, m2ZBN)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Teste de Vuong para diferença entre os modelos BN e HUBN
vuong(m1BN, m1HBN)
......@@ -120,7 +120,7 @@ vuong(m1BN, m1HBN)
vuong(m1ZBN, m1HBN)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Estimativas dos parâmetros e testes de Wald
summary(m1BN)
......@@ -128,7 +128,7 @@ summary(m1HBN)
summary(m1ZBN)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Ajuste de preditor do modelo proposto
## Retira o efeito de campista no preditor para as contagens não nula
......@@ -141,7 +141,7 @@ lmtest::lrtest(m1HBN, m3HBN)
vuong(m1BN, m3HBN)
## ----diag, cache = TRUE, fig.height = 4---------------------------
## ----diag, cache = TRUE, fig.height = 4----------------------------------
## Uma pequena avaliação dos resíduos
p1 <- xyplot(residuals(m3HBN) ~ fitted(m3HBN),
......@@ -167,13 +167,13 @@ print(p2, split = c(2, 1, 3, 1), more = TRUE)
print(p3, split = c(3, 1, 3, 1), more = FALSE)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Estimativas dos parâmetros e testes de Wald
summary(m3HBN)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
##----------------------------------------------------------------------
......@@ -213,7 +213,7 @@ barchart(py ~ y | factor(ncriancas), data = da,
))
## ----boot, cache = TRUE-------------------------------------------
## ----boot, cache = TRUE--------------------------------------------------
## Intervalos de confiança para predição
......@@ -250,7 +250,7 @@ xyplot(fit ~ npessoas | campista,
))
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Modelo Hurdle (binomial e binomial negativa)
m3HBN$formula
......@@ -279,17 +279,17 @@ cbind("glm_binomial" = coef(mzero),
# "zeroinfl" = m3HBN$theta)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
data(seguros)
str(seguros)
## help(seguros)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
summary(seguros)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
library(pscl)
......@@ -330,7 +330,7 @@ m1ZBN <- zeroinfl(f1, data = seguros, dist = "negbin")
m2ZBN <- zeroinfl(f2, data = seguros, dist = "negbin")
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Via log-verossimilhança
cbind("Poisson" = sapply(list(m0P, m1P, m2P), logLik),
......@@ -342,7 +342,7 @@ cbind("Poisson" = sapply(list(m0P, m1P, m2P), logLik),
)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Testes de razão de verossimilhanças
lmtest::lrtest(m0HP, m1HP, m2HP)
......@@ -350,14 +350,14 @@ lmtest::lrtest(m0HP, m1HP, m2HP)
lmtest::lrtest(m0HBN, m1HBN, m2HBN)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Para a componente de contagens não negativas é necessário um modelo
## Binomial Negativa, considerando superdispersão dos dados?
vuong(m1HP, m1HBN)
## -----------------------------------------------------------------
## ------------------------------------------------------------------------
## Estimativas do modelo proposto