# Capítulo 6. Introducción a los modelos
# lineales generalizados (GLM)

#  6.1. Estructura de los GLM
#  6.1.2. La función de vínculo

(y <- 5:15)
log(y)
(p <- seq(0.01, 0.99, 0.05))
log(p/(1-p))

#  6.3. Estimación de parámetros por máxima verosi-
#         militud
#  Cuadro 6.2. (continuación)

dist <- seq(0, 200, 20)
n <- length(dist)
Xmat <- model.matrix(~dist)
beta.vec <- c(4, -0.02)
lin.pred <- Xmat[,] %*% beta.vec
mu <- exp(lin.pred)
set.seed(1)
y <- rpois(n, mu)
plot(y ~ dist, ylab = "Riqueza de especies",
xlab = "Distancia al margen (m)")

range.intercept <- seq(0, 5, 0.1)
range.slope <- seq(-3, 0, 0.01)


MLE <- NULL
for (i in 1:length(range.intercept)) {
for (j in 1:length(range.slope)) {
pred.i <- exp(range.intercept[i] + range.slope[j] * dist)
logL.i <- sum(dpois(y, pred.i, log = T))
res.i <- c(range.intercept[i], range.slope[j],logL.i)
MLE <- rbind(MLE, res.i)
}
}

MLE <- as.data.frame(MLE)
names(MLE) <- c("intercept", "slope", "logL")
MLE[which.max(MLE$logL),]
beta.vec
iter <- length(range.intercept) * length(range.slope)
iter