# NOT RUN {
avec <- c(5, 10) # Alter these values parametrically
ivec <- c(3, 15) # Inflate these values
tvec <- c(6, 7) # Truncate these values
set.seed(1); pobs.mix <- pstr.mix <- 0.1
gdata <- data.frame(x2 = runif(nn <- 1000))
gdata <- transform(gdata, shape.p = loglink(0.1+0.1*x2, inverse = TRUE))
gdata <- transform(gdata,
y1 = rgaitzeta(nn, shape.p, alt.mix = avec, pobs.mix = pobs.mix,
inf.mix = ivec, pstr.mix = pstr.mix, truncate = tvec))
gaitzeta(alt.mix = avec, inf.mix = ivec, max.support = max.support)
with(gdata, table(y1))
spikeplot(with(gdata, y1), col = "pink2", log = "x")
gaitzxfit <- vglm(y1 ~ x2, crit = "coef", trace = TRUE, data = gdata,
gaitzeta(inf.mix = ivec, truncate = tvec,
alt.mix = avec, eq.ap = TRUE, eq.ip = TRUE))
head(fitted(gaitzxfit, type.fitted = "Pstr.mix"))
head(predict(gaitzxfit))
t(coef(gaitzxfit, matrix = TRUE)) # Easier to see with t()
summary(gaitzxfit, HDEtest = FALSE) # summary(gaitzxfit) is better
# }
Run the code above in your browser using DataLab