earg = list(cutpoint = 2)
polf("p", earg = earg, short = FALSE)
polf("p", earg = earg, tag = TRUE)
p = seq(0.01, 0.99, by=0.01)
y = polf(p, earg = earg)
y. = polf(p, earg = earg, deriv = 1)
max(abs(polf(y, earg = earg, inv = TRUE) - p)) # Should be 0
par(mfrow=c(2,1), las = 1)
plot(p, y, type = "l", col = "blue", main = "polf()")
abline(h=0, v=0.5, col = "red", lty = "dashed")
plot(p, y., type = "l", col = "blue",
main = "(Reciprocal of) first POLF derivative")
# Rutherford and Geiger data
ruge = data.frame(yy = rep(0:14,
times=c(57,203,383,525,532,408,273,139,45,27,10,4,0,1,1)))
with(ruge, length(yy)) # 2608 1/8-minute intervals
cutpoint = 5
ruge = transform(ruge, yy01 = ifelse(yy <= cutpoint, 0, 1))
earg = list(cutpoint=cutpoint)
fit = vglm(yy01 ~ 1, binomialff(link = "polf", earg = earg), ruge)
coef(fit, matrix = TRUE)
exp(coef(fit))
# Another example
pdat = data.frame(x2 = sort(runif(nn <- 1000)))
pdat = transform(pdat, x3 = runif(nn))
pdat = transform(pdat, mymu = exp( 3 + 1 * x2 - 2 * x3))
pdat = transform(pdat, y1 = rpois(nn, lambda=mymu))
cutpoints = c(-Inf, 10, 20, Inf)
pdat = transform(pdat, cuty = Cut(y1, breaks=cutpoints))
with(pdat, plot(x2, x3, col=cuty, pch=as.character(cuty)))
with(pdat, table(cuty) / sum(table(cuty)))
fit = vglm(cuty ~ x2 + x3, fam = cumulative(link = "polf",
reverse = TRUE, parallel = TRUE, intercept.apply = TRUE,
mv = TRUE, earg = list(cutpoint=cutpoints[2:3])),
pdat, trace = TRUE)
head(fit@y)
head(fitted(fit))
head(predict(fit))
coef(fit)
coef(fit, matrix = TRUE)
constraints(fit)
fit@misc$earg
Run the code above in your browser using DataLab