[6423] in s-news-athena
Summary regarding Null Deviance in GLM in Splus 3.0.
daemon@ATHENA.MIT.EDU (YONG CHANG )
Fri Feb 10 09:13:40 1995
From: ychang@welchlink.welch.jhu.edu (YONG CHANG )
To: s-news@utstat.toronto.edu (S Newsgroup)
Date: Fri, 10 Feb 1995 08:49:08 -0500 (EST)
This message is a followup on the message I posted two days ago regarding
some "unbelievable" null deviance I had encountered running Poisson regression
in Splus 3.0.
Prof. Brian Ripley and Trevor Hastie have offered help to me in a very timely
manner, and both pointed out that the error could be due to omission of offset
in the Null model in the old version 3.0 of Splus, which turned out to be
the case. I was using 3.0 version of Splus, instead of the latest 3.2(3.3?).
The following is the corrected script of "glm" function given by Prof. Ripley.
--------------------------------------------------------------
Here's the current (corrected) glm():
glm
function(formula = formula(data), family = gaussian, data = sys.parent(),
weights, subset, na.action, start = eta, control = glm.control(...),
method = "glm.fit", model = F, x = F, y = T, contrasts = NULL, ...)
{
call <- match.call()
m <- match.call(expand = F)
m$family <- m$method <- m$model <- m$x <- m$y <- m$control <- m$
contrasts <- m$... <- NULL
m[[1]] <- as.name("model.frame")
m <- eval(m, sys.parent())
Terms <- attr(m, "terms")
if(method == "model.frame")
return(m)
a <- attributes(m)
Y <- model.extract(m, response)
X <- model.matrix(Terms, m, contrasts)
w <- model.extract(m, weights)
if(!length(w))
w <- rep(1, nrow(m))
else if(any(w < 0))
stop("negative weights not allowed")
start <- model.extract(m, start)
offset <- model.extract(m, offset)
family <- as.family(family)
if(missing(method))
method <- attr(family, "method")
if(!is.null(method)) {
if(!exists(method, mode = "function"))
stop(paste("unimplemented method:", method))
}
else method <- "glm.fit"
glm.fitter <- get(method)
fit <- glm.fitter(x = X, y = Y, w = w, start = start, offset = offset,
family = family, maxit = control$maxit, epsilon = control$
# If an offset and intercept is present, iterations are needed to
# compute the Null deviance; these are done here, unless the model
# is NULL, in which case the computations have been done already
#
if(any(offset) && attr(Terms, "intercept")) {
null.deviance <- if(length(Terms)) glm.fitter(X[, "(Intercept)",
drop = F], Y, w, offset = offset, family =
family, maxit = control$maxit, epsilon =
control$epsilon, null.dev = NULL)$deviance
else fit$deviance
fit$null.deviance <- null.deviance
}
attr(fit, "class") <- c("glm", "lm")
fit$terms <- Terms
fit$formula <- as.vector(attr(Terms, "formula"))
fit$call <- call
if(model)
fit$model <- m
if(x)
fit$x <- X
if(!y)
fit$y <- NULL
fit
}