Modelos Lineales Generalizados
Modelo de regresión logístico
Pulso.data <- read.table(file = "http://tarwi.lamolina.edu.pe/~clopez/Categoricos/Pulso.txt", header = T)
head(Pulso.data)
## Pulso Peso
## 1 0 140
## 2 0 145
## 3 0 160
## 4 0 190
## 5 0 155
## 6 0 165
attach(Pulso.data)
Pulso.m1 <- glm(Pulso ~ Peso, family = binomial(link = logit))
summary(Pulso.m1)
##
## Call:
## glm(formula = Pulso ~ Peso, family = binomial(link = logit))
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.0072 -0.7881 -0.6682 -0.4573 1.9954
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.45452 1.59763 0.910 0.363
## Peso -0.01832 0.01127 -1.626 0.104
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 101.214 on 91 degrees of freedom
## Residual deviance: 98.377 on 90 degrees of freedom
## AIC: 102.38
##
## Number of Fisher Scoring iterations: 4
Pulso.m2 <- glm(Pulso ~ Peso, family = binomial(link = probit))
summary(Pulso.m2)
##
## Call:
## glm(formula = Pulso ~ Peso, family = binomial(link = probit))
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.0111 -0.7915 -0.6662 -0.4324 2.0186
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.887258 0.930572 0.953 0.3404
## Peso -0.011177 0.006485 -1.724 0.0848 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 101.214 on 91 degrees of freedom
## Residual deviance: 98.247 on 90 degrees of freedom
## AIC: 102.25
##
## Number of Fisher Scoring iterations: 4
detach(Pulso.data)
Modelo loglineal de Poisson
Cangrejo.data <- read.table(file = "http://tarwi.lamolina.edu.pe/~clopez/Categoricos/Cangrejo.txt", header = T)
head(Cangrejo.data)
## Color Cond Ancho Sat Peso
## 1 3 3 28.3 8 3050
## 2 4 3 22.5 0 1550
## 3 2 1 26.0 9 2300
## 4 4 3 24.8 0 2100
## 5 4 3 26.0 4 2600
## 6 3 3 23.8 0 2100
attach(Cangrejo.data)
Cangrejo.m1 <- glm(Sat ~ Ancho, family = poisson(link = log))
summary(Cangrejo.m1)
##
## Call:
## glm(formula = Sat ~ Ancho, family = poisson(link = log))
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.8526 -1.9884 -0.4933 1.0970 4.9221
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.30476 0.54224 -6.095 1.1e-09 ***
## Ancho 0.16405 0.01997 8.216 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 632.79 on 172 degrees of freedom
## Residual deviance: 567.88 on 171 degrees of freedom
## AIC: 927.18
##
## Number of Fisher Scoring iterations: 6
Modelo binomial negativo
library(MASS)
Cangrejo.m2 <- glm(Sat ~ Ancho, family = negative.binomial(theta = 1.0,link = "identity"), start = coef(Cangrejo.m1))
summary(Cangrejo.m2)
##
## Call:
## glm(formula = Sat ~ Ancho, family = negative.binomial(theta = 1,
## link = "identity"), start = coef(Cangrejo.m1))
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.8227 -1.4281 -0.2631 0.4644 2.1628
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -11.63290 0.98973 -11.75 <2e-16 ***
## Ancho 0.55395 0.04713 11.75 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(1) family taken to be 0.8997994)
##
## Null deviance: 224.93 on 172 degrees of freedom
## Residual deviance: 202.89 on 171 degrees of freedom
## AIC: 752.08
##
## Number of Fisher Scoring iterations: 8
Mínimos cuadrados ponderados iterativos
#Función de enlace identidad
Cangrejo.m3 <- glm(Sat ~ Ancho, family = poisson(link = identity), start = c(10, 0.5))
summary(Cangrejo.m3)
##
## Call:
## glm(formula = Sat ~ Ancho, family = poisson(link = identity),
## start = c(10, 0.5))
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.9113 -1.9598 -0.5405 1.0406 4.7988
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -11.52567 0.67700 -17.02 <2e-16 ***
## Ancho 0.54925 0.02966 18.52 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 632.79 on 172 degrees of freedom
## Residual deviance: 557.71 on 171 degrees of freedom
## AIC: 917.01
##
## Number of Fisher Scoring iterations: 23
y <- Sat
x <- Ancho
z <- y
X <- cbind(1, c(x))
beta <- c(10, 0.5)
i = 1
n = 23
beta0 <- rep(0, n)
beta1 <- rep(0, n)
while (i <= n)
{
eta <- beta%*%t(X)
w <- 1/eta
W <- diag(c(w))
b <- solve(t(X)%*%W%*%X)%*%t(X)%*%W%*%z
beta0[i] <- b[1]
beta1[i] <- b[2]
beta <- t(b)
i <- i + 1
}
data.frame(cbind(beta0, beta1))
## beta0 beta1
## 1 -10.53632 0.5116343
## 2 -11.22481 0.5378140
## 3 -11.36222 0.5430388
## 4 -11.41842 0.5451760
## 5 -11.44867 0.5463262
## 6 -11.46746 0.5470408
## 7 -11.48022 0.5475257
## 8 -11.48941 0.5478750
## 9 -11.49631 0.5481377
## 10 -11.50168 0.5483417
## 11 -11.50595 0.5485042
## 12 -11.50942 0.5486361
## 13 -11.51228 0.5487450
## 14 -11.51468 0.5488362
## 15 -11.51671 0.5489133
## 16 -11.51844 0.5489791
## 17 -11.51993 0.5490358
## 18 -11.52122 0.5490849
## 19 -11.52235 0.5491278
## 20 -11.52334 0.5491654
## 21 -11.52421 0.5491985
## 22 -11.52498 0.5492278
## 23 -11.52567 0.5492538
#Función de enlace log
Cangrejo.m4 <- glm(Sat ~ Ancho, family = poisson(link = log), start = c(0,0))
summary(Cangrejo.m4)
##
## Call:
## glm(formula = Sat ~ Ancho, family = poisson(link = log), start = c(0,
## 0))
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.8526 -1.9884 -0.4933 1.0970 4.9221
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.30476 0.54223 -6.095 1.1e-09 ***
## Ancho 0.16405 0.01996 8.217 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 632.79 on 172 degrees of freedom
## Residual deviance: 567.88 on 171 degrees of freedom
## AIC: 927.18
##
## Number of Fisher Scoring iterations: 7
beta <- c(0, 0)
i = 1
n = 7
beta0 <- rep(0, n)
beta1 <- rep(0, n)
while (i <= n)
{
eta <- beta%*%t(X)
z <- eta + y/exp(eta) - 1
w <- exp(eta)
W <- diag(c(w))
b <- solve(t(X)%*%W%*%X)%*%t(X)%*%W%*%t(z)
beta0[i] <- b[1]
beta1[i] <- b[2]
beta <- t(b)
i <- i + 1
}
data.frame(cbind(beta0, beta1))
## beta0 beta1
## 1 -11.424351 0.5073769
## 2 -10.185336 0.4379312
## 3 -7.281620 0.3181471
## 4 -4.366641 0.2054990
## 5 -3.385293 0.1672320
## 6 -3.305312 0.1640671
## 7 -3.304757 0.1640451