quasibinomialff()
quasibinomialff(link = "probit")
shunua = hunua[sort.list(with(hunua, altitude)), ] # Sort by altitude
fit = vglm(agaaus ~ poly(altitude, 2), binomialff(link = cloglog), shunua)
plot(agaaus ~ jitter(altitude), shunua, col = "blue", ylab = "P(Agaaus = 1)",
main = "Presence/absence of Agathis australis", las = 1)
with(shunua, lines(altitude, fitted(fit), col = "orange", lwd = 2))
# Fit two species simultaneously
fit2 = vgam(cbind(agaaus, kniexc) ~ s(altitude), binomialff(mv = TRUE), shunua)
with(shunua, matplot(altitude, fitted(fit2), type = "l",
main = "Two species response curves", las = 1))
# Shows that Fisher scoring can sometime fail. See Ridout (1990).
ridout = data.frame(v = c(1000, 100, 10), r = c(4, 3, 3), n = c(5, 5, 5))
(ridout = transform(ridout, logv = log(v)))
# The iterations oscillates between two local solutions:
glm.fail = glm(r/n ~ offset(logv) + 1, weight = n,
binomial(link = cloglog), ridout, trace = TRUE)
coef(glm.fail)
# vglm()'s half-stepping ensures the MLE of -5.4007 is obtained:
vglm.ok = vglm(cbind(r, n-r) ~ offset(logv) + 1,
binomialff(link = cloglog), ridout, trace = TRUE)
coef(vglm.ok)
Run the code above in your browser using DataLab