# fast generalized linear model fitting ```r fastglmPure( x, y, family = gaussian(), weights = rep(1, NROW(y)), offset = rep(0, NROW(y)), start = NULL, etastart = NULL, mustart = NULL, method = 0L, tol = 1e-07, maxit = 100L ) ``` ## Arguments - `x`: input model matrix. Must be a matrix object - `y`: numeric response vector of length nobs. - `family`: a description of the error distribution and link function to be used in the model. For `fastglmPure` this can only be the result of a call to a family function. (See `family` for details of family functions.) - `weights`: an optional vector of 'prior weights' to be used in the fitting process. Should be a numeric vector. - `offset`: this can be used to specify an a priori known component to be included in the linear predictor during fitting. This should be a numeric vector of length equal to the number of cases - `start`: starting values for the parameters in the linear predictor. - `etastart`: starting values for the linear predictor. - `mustart`: values for the vector of means. - `method`: an integer scalar with value 0 for the column-pivoted QR decomposition, 1 for the unpivoted QR decomposition, 2 for the LLT Cholesky, 3 for the LDLT Cholesky, 4 for the full pivoted QR decomposition, 5 for the Bidiagonal Divide and Conquer SVD - `tol`: threshold tolerance for convergence. Should be a positive real number - `maxit`: maximum number of IRLS iterations. Should be an integer ## Returns A list with the elements - **coefficients**: a vector of coefficients - **se**: a vector of the standard errors of the coefficient estimates - **rank**: a scalar denoting the computed rank of the model matrix - **df.residual**: a scalar denoting the degrees of freedom in the model - **residuals**: the vector of residuals - **s**: a numeric scalar - the root mean square for residuals - **fitted.values**: the vector of fitted values ## Examples ```r 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)) ```

fastglmPure function