p <- seq(0.01, 0.99, by = 0.01)
cloglog(p)
max(abs(cloglog(cloglog(p), inverse = TRUE) - p)) # Should be 0
p <- c(seq(-0.02, 0.02, by = 0.01), seq(0.97, 1.02, by = 0.01))
cloglog(p) # Has NAs
cloglog(p, bvalue = .Machine$double.eps) # Has no NAs
p <- seq(0.01, 0.99, by = 0.01)
plot(p, logit(p), type = "l", col = "limegreen", ylab = "transformation",
lwd = 2, las = 1, main = "Some probability link functions")
lines(p, probit(p), col = "purple", lwd = 2)
lines(p, cloglog(p), col = "chocolate", lwd = 2)
lines(p, cauchit(p), col = "tan", lwd = 2)
abline(v = 0.5, h = 0, lty = "dashed")
legend(0.1, 4, c("logit", "probit", "cloglog", "cauchit"),
col = c("limegreen", "purple", "chocolate", "tan"), lwd = 2)
# This example shows that a cloglog link is preferred over the logit
n <- 500; p <- 5; S <- 3; Rank <- 1 # Species packing model:
mydata <- rcqo(n, p, S, eq.tol = TRUE, es.opt = TRUE, eq.max = TRUE,
family = "binomial", hi.abundance = 5, seed = 123, Rank = Rank)
fitc <- cqo(attr(mydata, "formula"), ITol = TRUE, data = mydata,
fam = binomialff(mv = TRUE, link = "cloglog"), Rank = Rank)
fitl <- cqo(attr(mydata, "formula"), ITol = TRUE, data = mydata,
fam = binomialff(mv = TRUE, link = "logit"), Rank = Rank)
# Compare the fitted models (cols 1 and 3) with the truth (col 2)
cbind(concoef(fitc), attr(mydata, "ccoefficients"), concoef(fitl))
Run the code above in your browser using DataLab