model.matrix()：为什么在这种情况下我会失去对对比度的控制
model.matrix(): why do I lose control of contrast in this case
假设我们有一个玩具数据框：
x <- data.frame(x1 = gl(3, 2, labels = letters[1:3]),
x2 = gl(3, 2, labels = LETTERS[1:3]))
我要构建模型矩阵
# x1b x1c x2B x2C
# 1 0 0 0 0
# 2 0 0 0 0
# 3 1 0 1 0
# 4 1 0 1 0
# 5 0 1 0 1
# 6 0 1 0 1
作者：
model.matrix(~ x1 + x2 - 1, data = x,
contrasts.arg = list(x1 = contr.treatment(letters[1:3]),
x2 = contr.treatment(LETTERS[1:3])))
但实际上我得到：
# x1a x1b x1c x2B x2C
# 1 1 0 0 0 0
# 2 1 0 0 0 0
# 3 0 1 0 1 0
# 4 0 1 0 1 0
# 5 0 0 1 0 1
# 6 0 0 1 0 1
# attr(,"assign")
# [1] 1 1 1 2 2
# attr(,"contrasts")
# attr(,"contrasts")$x1
# b c
# a 0 0
# b 1 0
# c 0 1
# attr(,"contrasts")$x2
# B C
# A 0 0
# B 1 0
# C 0 1
我在这里有点困惑：
- 我已经传入显式对比矩阵以降低第一个因子水平；
- 我已要求放弃拦截。
那为什么我得到一个有 5 列的模型矩阵？我怎样才能得到我想要的模型矩阵？
每当我们在 R 级别失去对某些东西的控制时，在 C 级别肯定会有一些默认的、不可更改的行为。 model.matrix.default() 的 C 代码可以在 R 中找到源包位于：
R-<release_number>/src/library/stats/src/model.c
我们可以在这里找到解释：
/* If there is no intercept we look through the factor pattern */
/* matrix and adjust the code for the first factor found so that */
/* it will be coded by dummy variables rather than contrasts. */
让我们用一个数据框对此做一个小测试
x <- data.frame(x1 = gl(2, 2, labels = letters[1:2]), x2 = sin(1:4))
如果我们在RHS上只有x2，我们可以成功拦截：
model.matrix(~ x2 - 1, data = x)
# x2
#1 0.8414710
#2 0.9092974
#3 0.1411200
#4 -0.7568025
如果 RHS 上只有 x1，则不应用对比度：
model.matrix(~ x1 - 1, data = x)
# x1a x1b
#1 1 0
#2 1 0
#3 0 1
#4 0 1
当我们同时拥有 x1 和 x2 时，不应用对比度：
model.matrix(~ x1 + x2 - 1, data = x)
# x1a x1b x2
#1 1 0 0.8414710
#2 1 0 0.9092974
#3 0 1 0.1411200
#4 0 1 -0.7568025
这意味着虽然两者之间存在差异：
lm(y ~ x2, data = x)
lm(y ~ x2 - 1, data = x)
没有区别
lm(y ~ x1, data = x)
lm(y ~ x1 - 1, data = x)
或
lm(y ~ x1 + x2, data = x)
lm(y ~ x1 + x2 - 1, data = x)
之所以有这样的行为，并不是为了保证数值的稳定性，而是为了保证估计/预测的敏感性。如果我们在对 x1 应用对比时真的放弃截距，我们最终会得到一个模型矩阵：
# x1b
#1 0
#2 0
#3 1
#4 1
效果是我们将级别 a 的估计限制为 0。
在这个post中：，我们有一个数据集：
# Y X1 X2
#1 1.8376852 TRUE TRUE
#2 -2.1173739 TRUE FALSE
#3 1.3054450 FALSE TRUE
#4 -0.3476706 TRUE FALSE
#5 1.3219099 FALSE TRUE
#6 0.6781750 FALSE TRUE
此数据集中没有联合存在(X1 = FALSE, X2 = FALSE)。但从广义上讲，model.matrix() 必须做一些安全和明智的事情。假设训练数据集中没有两个因子水平的联合存在意味着它们不需要预测是有偏见的。如果我们在应用对比时真的放弃拦截，这种联合存在被限制在 0。但是，post 的 OP 故意想要这种非标准行为（出于某种原因），在这种情况下，给出了一个可能的解决方法在我的回答中。
假设我们有一个玩具数据框：
x <- data.frame(x1 = gl(3, 2, labels = letters[1:3]),
x2 = gl(3, 2, labels = LETTERS[1:3]))
我要构建模型矩阵
# x1b x1c x2B x2C
# 1 0 0 0 0
# 2 0 0 0 0
# 3 1 0 1 0
# 4 1 0 1 0
# 5 0 1 0 1
# 6 0 1 0 1
作者：
model.matrix(~ x1 + x2 - 1, data = x,
contrasts.arg = list(x1 = contr.treatment(letters[1:3]),
x2 = contr.treatment(LETTERS[1:3])))
但实际上我得到：
# x1a x1b x1c x2B x2C
# 1 1 0 0 0 0
# 2 1 0 0 0 0
# 3 0 1 0 1 0
# 4 0 1 0 1 0
# 5 0 0 1 0 1
# 6 0 0 1 0 1
# attr(,"assign")
# [1] 1 1 1 2 2
# attr(,"contrasts")
# attr(,"contrasts")$x1
# b c
# a 0 0
# b 1 0
# c 0 1
# attr(,"contrasts")$x2
# B C
# A 0 0
# B 1 0
# C 0 1
我在这里有点困惑：
- 我已经传入显式对比矩阵以降低第一个因子水平；
- 我已要求放弃拦截。
那为什么我得到一个有 5 列的模型矩阵？我怎样才能得到我想要的模型矩阵？
每当我们在 R 级别失去对某些东西的控制时，在 C 级别肯定会有一些默认的、不可更改的行为。 model.matrix.default() 的 C 代码可以在 R 中找到源包位于：
R-<release_number>/src/library/stats/src/model.c
我们可以在这里找到解释：
/* If there is no intercept we look through the factor pattern */
/* matrix and adjust the code for the first factor found so that */
/* it will be coded by dummy variables rather than contrasts. */
让我们用一个数据框对此做一个小测试
x <- data.frame(x1 = gl(2, 2, labels = letters[1:2]), x2 = sin(1:4))
如果我们在RHS上只有
x2，我们可以成功拦截：model.matrix(~ x2 - 1, data = x) # x2 #1 0.8414710 #2 0.9092974 #3 0.1411200 #4 -0.7568025如果 RHS 上只有
x1，则不应用对比度：model.matrix(~ x1 - 1, data = x) # x1a x1b #1 1 0 #2 1 0 #3 0 1 #4 0 1当我们同时拥有
x1和x2时，不应用对比度：model.matrix(~ x1 + x2 - 1, data = x) # x1a x1b x2 #1 1 0 0.8414710 #2 1 0 0.9092974 #3 0 1 0.1411200 #4 0 1 -0.7568025
这意味着虽然两者之间存在差异：
lm(y ~ x2, data = x)
lm(y ~ x2 - 1, data = x)
没有区别
lm(y ~ x1, data = x)
lm(y ~ x1 - 1, data = x)
或
lm(y ~ x1 + x2, data = x)
lm(y ~ x1 + x2 - 1, data = x)
之所以有这样的行为，并不是为了保证数值的稳定性，而是为了保证估计/预测的敏感性。如果我们在对 x1 应用对比时真的放弃截距，我们最终会得到一个模型矩阵：
# x1b
#1 0
#2 0
#3 1
#4 1
效果是我们将级别 a 的估计限制为 0。
在这个post中：
# Y X1 X2
#1 1.8376852 TRUE TRUE
#2 -2.1173739 TRUE FALSE
#3 1.3054450 FALSE TRUE
#4 -0.3476706 TRUE FALSE
#5 1.3219099 FALSE TRUE
#6 0.6781750 FALSE TRUE
此数据集中没有联合存在(X1 = FALSE, X2 = FALSE)。但从广义上讲，model.matrix() 必须做一些安全和明智的事情。假设训练数据集中没有两个因子水平的联合存在意味着它们不需要预测是有偏见的。如果我们在应用对比时真的放弃拦截，这种联合存在被限制在 0。但是，post 的 OP 故意想要这种非标准行为（出于某种原因），在这种情况下，给出了一个可能的解决方法在我的回答中。