5

I use the R code:

dat<-data.frame(p1=c(0,1,1,0,0), GAMMA.1=c(1,2,3,4,3), VAR1=c(2,2,1,3,4), GAMMA.2=c(1,1,3,4,1))
form <- p1 ~ GAMMA.1:VAR1 + GAMMA.2:VAR1
mod <- glm(formula=form, data=dat, family=binomial)
(coef <- coefficients(mod))

# (Intercept) GAMMA.1:VAR1 VAR1:GAMMA.2 
#   1.7974974   -0.2563667   -0.2181079

As we can see the names of coef for the interaction GAMMA.2:VAR1 is not in the same order as in form (we have VAR1:GAMMA.2 instead). For several reasons, I need the output

# (Intercept) GAMMA.1:VAR1   GAMMA.2:VAR1
#   1.7974974   -0.2563667   -0.2181079

without changing the names of the coefficients afterwards. Specifically, I want the same names for the coefficients as I used in the form object (without switching as in the code above). Can I tell glm() not to switch the names of the interactions?

Giuseppe
  • 786
  • 1
  • 7
  • 18

2 Answers2

4

The answer is no, not without a lot of rewriting functions. The order of the label of interaction terms is determined by the terms.formula function, which itself is determined by the termsform function buried deep in the C code. There are no parameters that you can pass termsform that give you the behaviour that you want (although keep.order looked promising, it does not do what you want).

You would have to rewrite the terms.formula function to "swap back" the names after output from termsform, and then override the terms.formula function with your patched version, but are you sure you want that? It will be far easier to change the names of the coefficients afterwards.

You could also use terms.formula preemptively, and determine how your formula would be reordered, using and then create a mapping vector.

dat<-data.frame(p1=c(0,1,1,0,0), GAMMA.1=c(1,2,3,4,3), VAR1=c(2,2,1,3,4), GAMMA.2=c(1,1,3,4,1))
form <- p1 ~ GAMMA.1:VAR1 + GAMMA.2:VAR1
new.names<-labels(terms(form,data=dat,keep.order=TRUE))
names(new.names)<-as.character(form[[3]][-1])
new.names
# GAMMA.1:VAR1   GAMMA.2:VAR1 
# "GAMMA.1:VAR1" "VAR1:GAMMA.2"

You could use that vector to map names if you had the need later on.

nograpes
  • 18,623
  • 1
  • 44
  • 67
  • More specifically the culprit is `stats:::C_termsform` (see [model.c](https://github.com/wch/r-source/blob/trunk/src/library/stats/src/model.c), lines 970++). – jay.sf Apr 10 '23 at 07:32
1

I have two possible workarounds for you.

One workaround comes from the observation that terms within interaction labels are ordered according to their order of appearance in the formula. In your example, the order is GAMMA.1, VAR1, GAMMA.2. You may be able to rewrite your formula with a different ordering so that the formula and coefficient names match:

form <- p1 ~ VAR1:GAMMA.1 + VAR1:GAMMA.2
mod <- glm(formula=form, data=dat, family=binomial)
coefficients(mod)

# (Intercept) VAR1:GAMMA.1 VAR1:GAMMA.2 
#   1.7974974   -0.2563667   -0.2181079

Another workaround is to rename the coefficients according to the formula when you pull them out:

rhs_terms <- c("(Intercept)",as.character(form[[3]][2:length(form[[3]])]))
(coef <- setNames(coefficients(mod), rhs_terms))

# (Intercept) GAMMA.1:VAR1 GAMMA.2:VAR1 
#   1.7974974   -0.2563667   -0.2181079
pangia
  • 1,164
  • 10
  • 11