Trunc(c(1, 8), 2)
set.seed(1) # The following example is based on the normal
mymean <- 20; m.truth <- 3 # approximation to the Poisson.
gdata <- data.frame(y1 = round(rnorm((nn <- 1000), mymean,
sd = sqrt(mymean / m.truth))))
org1 <- with(gdata, range(y1)) # Original range of the raw data
m.max <- 5 # Try multipliers 1:m.max
logliks <- numeric(m.max)
names(logliks) <- as.character(1:m.max)
for (i in 1:m.max) {
logliks[i] <- logLik(vglm(i * y1 ~ offset(rep(log(i), nn)),
gaitdpoisson(truncate = Trunc(org1, i)), data = gdata))
}
sort(logliks, decreasing = TRUE) # Best to worst
if (FALSE) par(mfrow = c(1, 2))
plot(with(gdata, table(y1))) # Underdispersed wrt Poisson
plot(logliks, col = "blue", type = "b", xlab = "Multiplier")
Run the code above in your browser using DataLab