Best way to extract link category from glm model?

I am writing a function that takes an object fulland reduced glmto summarize the results of the interaction of the variable of interest varofintand the interaction variable interaction_var(by executing lrtestand using the svycontrastobject fullto get the results varofintfor each level interaction_var). Sample data:

x <- data.frame(outcome=rbinom(100,1,.3),varofint=rnorm(100), interaction_var=sample(letters[1:3],100,replace=TRUE))

reduced <- glm(outcome~varofint+interaction_var,data=x)
full <- glm(outcome~varofint*interaction_var,data=x)

I would like to know the best way to extract the reference category for the specified ( full) glm model. I could do something like

levels(full$data$interaction_var)[1]

but will it be a "safe" method for extracting a reference data category into arguments contrasts? It seems that, given the choice of SAS contrast, this method may lead to a level interactionv_varthat is not used as a reference category in the model. Will the following be safer?

mf <- model.frame(full)
setdiff(rownames(contrasts(mf[, "interaction_var"])), colnames(contrasts(mf[, "interaction_var"])))

or similar

names(which(apply(contrasts(mf[, "interaction_var"]),1,function(.v){all(.v==0)})))

Am I missing an easier way to retrieve a link category?

+4
source share
1 answer

Here is the function for this task:

refCat <- function(model, var) {
  cs <- attr(model.matrix(model), "contrasts")[[var]]
  if (is.character(cs)) {
    if (cs == "contr.treatment")
      ref <- 1
    else stop("No treatment contrast")
  }  
  else {
    zeroes <- !cs
    ones <- cs == 1
    stopifnot(all(zeroes | ones))
    cos <- colSums(ones)
    stopifnot(all(cos == 1))
    ros <- rowSums(ones)
    stopifnot(sum(!ros) == 1 && sum(ros) != ncol(cs))
    ref <- which(!ros)
  }
  return(levels(model$data[[var]])[ref])  
}    

The function stops if the variable is varnot presented as treatment contrasts.

Examples:

refCat(reduced, "interaction_var")
# [1] "a"
refCat(full, "interaction_var")
# [1] "a"
+1
source

All Articles