x = runif(n <- 2000)
p0 = logit(-1 + 2*x, inverse=TRUE)
y1 = rposnegbin(n, munb=exp(0+2*x), k=exp(1)) # With covariates
y2 = rposnegbin(n, munb=exp(1+2*x), k=exp(1)) # With covariates
y1 = ifelse(runif(n) < p0, 0, y1)
y2 = ifelse(runif(n) < p0, 0, y2)
table(y1)
table(y2)
fit = vglm(cbind(y1,y2) ~ x, zanegbinomial, trace=TRUE)
coef(fit, matrix=TRUE)
fitted(fit)[1:9,]
predict(fit)[1:9,]
Run the code above in your browser using DataLab