# NOT RUN {
set.seed(1)
x <- matrix(rnorm(1000 * 25), ncol = 25)
eta <- 0.1 + 0.25 * x[,1] - 0.25 * x[,3] + 0.75 * x[,5] -0.35 * x[,6] #0.25 * x[,1] - 0.25 * x[,3]
y <- 1 * (eta > rnorm(1000))
yp <- rpois(1000, eta ^ 2)
yg <- rgamma(1000, exp(eta) * 1.75, 1.75)
# binomial
system.time(gl1 <- glm.fit(x, y, family = binomial()))
system.time(gf1 <- fastglmPure(x, y, family = binomial(), tol = 1e-8))
system.time(gf2 <- fastglmPure(x, y, family = binomial(), method = 1, tol = 1e-8))
system.time(gf3 <- fastglmPure(x, y, family = binomial(), method = 2, tol = 1e-8))
system.time(gf4 <- fastglmPure(x, y, family = binomial(), method = 3, tol = 1e-8))
max(abs(coef(gl1) - gf1$coef))
max(abs(coef(gl1) - gf2$coef))
max(abs(coef(gl1) - gf3$coef))
max(abs(coef(gl1) - gf4$coef))
# poisson
system.time(gl1 <- glm.fit(x, yp, family = poisson(link = "log")))
system.time(gf1 <- fastglmPure(x, yp, family = poisson(link = "log"), tol = 1e-8))
system.time(gf2 <- fastglmPure(x, yp, family = poisson(link = "log"), method = 1, tol = 1e-8))
system.time(gf3 <- fastglmPure(x, yp, family = poisson(link = "log"), method = 2, tol = 1e-8))
system.time(gf4 <- fastglmPure(x, yp, family = poisson(link = "log"), method = 3, tol = 1e-8))
max(abs(coef(gl1) - gf1$coef))
max(abs(coef(gl1) - gf2$coef))
max(abs(coef(gl1) - gf3$coef))
max(abs(coef(gl1) - gf4$coef))
# gamma
system.time(gl1 <- glm.fit(x, yg, family = Gamma(link = "log")))
system.time(gf1 <- fastglmPure(x, yg, family = Gamma(link = "log"), tol = 1e-8))
system.time(gf2 <- fastglmPure(x, yg, family = Gamma(link = "log"), method = 1, tol = 1e-8))
system.time(gf3 <- fastglmPure(x, yg, family = Gamma(link = "log"), method = 2, tol = 1e-8))
system.time(gf4 <- fastglmPure(x, yg, family = Gamma(link = "log"), method = 3, tol = 1e-8))
max(abs(coef(gl1) - gf1$coef))
max(abs(coef(gl1) - gf2$coef))
max(abs(coef(gl1) - gf3$coef))
max(abs(coef(gl1) - gf4$coef))
# }
Run the code above in your browser using DataCamp Workspace