# [S] All Contrasts for a Factor in a glm Fit

John Wallace (jrw@fish.washington.edu)
Fri, 29 May 1998 23:16:17 -0700 (PDT)

The recent thread on the reference category in factors inspired me to
debug a function I wrote last year that gives all the contrasts (without
redundancy) of a particular factor within a glm fit (factors in aov and lm
fits should also work).

The function, all.contr(), is attached below.

Keeping with Prof. Ripley's quine data example:

> library(MASS)
> options(contrasts=c("contr.treatment", "contr.poly"))
> attach(quine)

> all.contr(glm(Days ~ Sex + Age + Eth, family=poisson), Age)

Level F0 contrasted with level(s) F1, F2, F3,

Value Std. Error t value Pr(>|t|)
AgeF1 -0.2262100 0.2532980 -0.8930589 0.3733582
AgeF2 0.3490078 0.2259269 1.5447816 0.1246566
AgeF3 0.2968470 0.2381471 1.2464857 0.2146677

Level F1 contrasted with level(s) F2, F3,

Value Std. Error t value Pr(>|t|)
AgeF2 0.3490078 0.2259269 1.544782 0.1246566
AgeF3 0.2968470 0.2381471 1.246486 0.2146677

Level F2 contrasted with level(s) F3,

Value Std. Error t value Pr(>|t|)
AgeF3 0.296847 0.2381471 1.246486 0.2146677

If, instead of attach'ing, you used the data argument in the glm() function,
then you must use the data argument in the all.contr() function:

> all.contr(glm(Days ~ Sex + Age + Eth, family=poisson,
data=quine), Age, data = quine)

That's all there is to it. (Oh, and it is S v4 ready with some help at
the top of the function.)

```--
John Wallace
University of Washington                    ^    ^    ^
Fisheries Research Institute               / \  / \  / \   ^
Box 357980                                 / \  / \   |   / \
Seattle, WA 98195-7980                      |    |  o__~  / \
PHONE   (206) 543-1513                  @ @         /\/\   |
FAX     (206) 685-7471                   ~
E-MAIL        jw@u.washington.edu              o
WWW      http://www.fish.washington.edu/people/jrw/Wallace.html
o  _///_ //
<`)=  _<<
\\\  \\
--------------------------------------

"all.contr" <- function(fit, variable, data)
{
#   'fit' is an aov, lm, or glm object.
#   'variable' is a factor for which all the contrasts are wanted.
#    Use 'data' if the data argument was used to make 'fit'.
#
#   DATE WRITTEN:  March 1997 REVISED: May 1998
#   Author:  John R. Wallace (jrw@fish.washington.edu)
#
var.name <- deparse(substitute(variable))
if(!missing(data))
assign(var.name, data[, var.name], frame = 1)
var.lev.names <- levels(variable)
contrasts(variable) <- contr.treatment(var.lev.names)
dummy <- contr.treatment(var.lev.names, contrasts = F)
CALL <- summary(fit)\$call
fit.names <- dimnames(summary.lm(eval(CALL))\$coef)[[1]]
var.levels <- paste(var.name, var.lev.names, sep = "")
if(length(var.levels) == 2) {
fit.names[var.name == fit.names] <- var.levels[2]
}
fit.rows <- grep(var.levels, fit.names)
var.cols <- c(1, grep(fit.names, var.levels))
var.num <- length(var.cols)
for(i in 1:(var.num - 1)) {
cat("\n Level", var.lev.names[var.cols[i]],
"contrasted with level(s)", paste(
var.lev.names[var.cols[(i + 1):var.num]],
", ", sep = "", collapse = ""), "\n\n")
contrasts(variable) <- dummy[,  - var.cols[i]]
print(summary.lm(eval(CALL))\$coef[fit.rows[i:(
var.num - 1)],  , drop = F])
}
invisible()
}

------------------------------------
-----------------------------------------------------------------------
This message was distributed by s-news@wubios.wustl.edu.  To unsubscribe
send e-mail to s-news-request@wubios.wustl.edu with the BODY of the
message:  unsubscribe s-news
```