Potential Bug / different defaults for Intercept / Reference Levels when using `Formula.get_model_matrix()` with categoricals

[s3alfisc]

Hi @matthewwardrop, I am wondering if the following behavior is a bug or not:

model_matrix() and Formula.get_model_matrix() handle reference levels for categoricals differently: model_matrix includes an intercept and drops a reference level by default, while get_model_matrix() does not.

from formulaic import Formula, model_matrix
import pandas as pd

df = pd.DataFrame({"a": [1, 2, 3], "b": [1, 2, 3], "c": [1, 2, 3]})

_, rhs = model_matrix("a ~ C(b)", data = df)
rhs
# Intercept C(b)[T.2]   C(b)[T.3]
# 0 1.0 0   0
# 1 1.0 1   0
# 2 1.0 0   1

Formula(fml1 = "a ~ C(b)").get_model_matrix(df)
#    .rhs:
#           C(b)[T.1]  C(b)[T.2]  C(b)[T.3]
#        0          1          0          0
#        1          0          1          0
#        2          0          0          1

If intended behavior, is it possible to mimic model_matrix() behavior for Formula.get_model_matrix()?

[s3alfisc]

Thought about this for a bit and concluded it is likely intended behavior =)

[matthewwardrop]

Hi @s3alfisc ,

Just wanted to clarify a few things here:

• model_matrix("a ~ C(b)", df) and Formula("a ~ C(b)").get_model_matrix(df, context=locals()) will give identical output (unless you are using globals too).
• Nested formula are treated as not requiring the intercept by default; that is Formula(nested="a ~ C(b)").get_model_matrix(df, context=locals()) will have a nested model matrix at .nested which does not have the intercept. This is because in most cases where this is used, it is unlikely you want duplicate intercepts. Note that multipart formula like Formula("y ~ x1 | x2 | x3") do not have this behaviour (the rhs formula will all have intercepts).
• You can customize the behaviour of this using the _parser and _nested_parser arguments to Formula; e.g. Formula(nested="a ~ C(b)", _nested_parser=DefaultFormulaParser(include_intercept=True). You can read more in the Formula docstring.
• You can also bypass this logic by nesting formulae explicitly: Formula(nested=Formula("y~x"))

Hope that helps.