# Nonparametric proportional odds model
pneumo <- transform(pneumo, let = log(exposure.time))
vgam(cbind(normal, mild, severe) ~ s(let),
cumulative(parallel = TRUE), data = pneumo, trace = TRUE)
# Nonparametric logistic regression
fit <- vgam(agaaus ~ s(altitude, df = 2), binomialff, data = hunua)
plot(fit, se = TRUE)
pfit <- predict(fit, type = "terms", raw = TRUE, se = TRUE)
names(pfit)
head(pfit$fitted)
head(pfit$se.fit)
pfit$df
pfit$sigma
# Fit two species simultaneously
fit2 <- vgam(cbind(agaaus, kniexc) ~ s(altitude, df = c(2, 3)),
binomialff(mv = TRUE), data = hunua)
coef(fit2, matrix = TRUE) # Not really interpretable
plot(fit2, se = TRUE, overlay = TRUE, lcol = 1:2, scol = 1:2)
ooo <- with(hunua, order(altitude))
with(hunua, matplot(altitude[ooo], fitted(fit2)[ooo,], ylim = c(0, 0.8),
xlab = "Altitude (m)", ylab = "Probability of presence", las = 1,
main = "Two plant species' response curves", type = "l", lwd = 2))
with(hunua, rug(altitude))
# The subset argument does not work here. Use subset() instead.
set.seed(1)
zdata <- data.frame(x2 = runif(nn <- 100))
zdata <- transform(zdata, y = rbinom(nn, 1, 0.5))
zdata <- transform(zdata, subS = runif(nn) < 0.7)
sub.zdata <- subset(zdata, subS) # Use this instead
if (FALSE)
fit4a <- vgam(cbind(y, y) ~ s(x2, df = 2), binomialff(mv = TRUE),
data = zdata, subset = subS) # This fails!!!
fit4b <- vgam(cbind(y, y) ~ s(x2, df = 2), binomialff(mv = TRUE),
data = sub.zdata) # This succeeds!!!
Run the code above in your browser using DataLab