两个因子变量之间具有相等约束的线性模型
此问题与https://stats.stackexchange.com/questions/3143/linear-model-with-constraints有关，但情况略有不同。两个因子变量之间具有相等约束的线性模型
我有一个简单的双因子线性模型，连续结果Y。 factor1具有〜350个分类值，并且factor2具有相同的〜350个分类。我想限制每个级别上的系数总和为零跨越这两个因素。
（这样做的原因是，factor1和factor2每个级别中的任何训练例如正或负进入，但从来没有在相同的例子中出现两次。）
这里是示出这种情况的示例数据集，那里有四个级别每个因素：
Y factor1 factor2
1 -1.2470416 A B
2 4.3368592 C D
3 1.0005147 D A
4 -2.8309146 A C
5 1.7501315 B D
6 -0.8372193 B A
7 3.3542627 C A
8 4.3319422 D C
9 1.4937895 D B
10 2.0951559 A D
11 -2.6610207 C D
12 -4.9917367 D B
13 2.2424169 D A
14 1.0205409 C A
15 -3.4584576 C B
统计模型我想估计是： $$Ý_ {（I，J）} = \ alpha_i- \ beta_j + \ varepsilon _ {（I，J） } $$ 其中$（i，j）$是取决于配对的结果。 factor1标记$ i $和factor2标记$ j $。如果组A出现在factor2中，则A上的参数应该等于factor1中出现的负数。因此，我想为所有$ i $和$ j $设置$ \ alpha $等于$ \ beta $。
我可以在lm()如下估计这个模型的（无意义的）版本很容易：
Y <- c(-1.2470416, 4.3368592 , 1.0005147 , -2.8309146 , 1.7501315 , -0.8372193 , 3.3542627 , 4.3319422 , 1.4937895 , 2.0951559 , -2.6610 207 , -4.9917367 , 2.2424169 , 1.0205409 , -3.4584576)
factor1 <- c("A" , "C" , "D" , "A" , "B" , "B" , "C" , "D" , "D" , "A" , "C" , "D" , "D" , "C" , "C")
factor2 <- c("B", "D", "A", "C", "D", "A", "A", "C", "B", "D", "D", "B", "A", "A", "B")
DF <- data.frame(Y,factor1,factor2)
lm(Y~factor1+factor2,data=DF)
，我得到下面的输出：
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.5363 2.5856 0.207 0.841
factor1B -0.4579 3.1121 -0.147 0.887
factor1C 0.4047 2.4925 0.162 0.875
factor1D 1.8737 2.4098 0.778 0.459
factor2B -3.6252 2.2050 -1.644 0.139
factor2C -0.7226 2.8903 -0.250 0.809
factor2D 0.7561 2.2094 0.342 0.741
需要注意的是，从理论上说，factor1C应等于-factor2C由我的模型决定。在简单的lm()输出中情况并非如此，因为我没有施加任何限制。
所以我想要做的是什么估计
Y ~ factor1 + factor2 [subject to factor1+factor2=0 for each level of factor1, factor2]
用简单的英语，这将是像
model2 <- lm(Y~factor1-factor2, data=DF)
但是，这当然不是R如何解释表达（因为在model语句中放一个减号告诉R从模型中排除该变量）。
我读过对比，但我不认为有办法做到这一点。我也阅读了glmc，但没有看到一个简单的方法将它纳入这个多层次的因素。此外，我不清楚生成新的factor3 = factor1-factor2是针对此特定场景的明确操作。最后，我尝试运行model3 <- lm(Y+factor2 ~ factor1, data=DF)但收到错误。
我的感觉是，我需要通过循环每个变量的级别来创建一个约束矩阵。我对R很新，我不确定这是如何完成的。任何帮助，将不胜感激。
请注意，这是很容易做到这一点在Stata，如下：
input ID y factor1 factor2
1 -1.2470416 1 2
2 4.3368592 3 4
3 1.0005147 4 1
4 -2.8309146 1 3
5 1.7501315 2 4
6 -0.8372193 2 1
7 3.3542627 3 1
8 4.3319422 4 3
9 1.4937895 4 2
10 2.0951559 1 4
11 -2.6610207 3 4
12 -4.9917367 4 2
13 2.2424169 4 1
14 1.0205409 3 1
15 -3.4584576 3 2
end
constraint 1 2.factor1 = -2.factor2
constraint 2 3.factor1 = -3.factor2
constraint 3 4.factor1 = -4.factor2
cnsreg y i.factor1 i.factor2, constraints(1/3)
这给下面的输出：
Constrained linear regression Number of obs = 15
F( 3, 11) = 0.73
Prob > F = 0.5554
Root MSE = 2.9875
(1) 2.factor1 + 2.factor2 = 0
(2) 3.factor1 + 3.factor2 = 0
(3) 4.factor1 + 4.factor2 = 0
------------------------------------------------------------------------------
y | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
factor1 |
B | 2.104393 1.439085 1.46 0.172 -1.063011 5.271798
C | .5222649 1.377463 0.38 0.712 -2.509511 3.55404
D | .6589209 1.266188 0.52 0.613 -2.127941 3.445783
|
factor2 |
B | -2.104393 1.439085 -1.46 0.172 -5.271798 1.063011
C | -.5222649 1.377463 -0.38 0.712 -3.55404 2.509511
D | -.6589209 1.266188 -0.52 0.613 -3.445783 2.127941
|
_cons | .5054862 .829675 0.61 0.555 -1.320616 2.331589
------------------------------------------------------------------------------
怎样才能做到上面R中？
正如https://stats.stackexchange.com/questions/3143/linear-model-with-constraints中最受欢迎（但未被接受）的答案所指出的，通过创建一个新变量可以很容易地解决这个问题，这个变量是“单热”编码因子的差异。
在Stata，可以做到这一点如下：
* one-hot encode each of the factors
qui tab factor1, gen(f1dum)
qui tab factor2, gen(f2dum)
* generate difference in one-hot vectors
forv x=1/4{
gen fdiffdum`x' = f1dum`x'-f2dum`x'
}
* regress y on differenced one-hot vectors
reg y fdiffdum2 fdiffdum3 fdiffdum4
其中给出以下输出：
Source | SS df MS Number of obs = 15
-------------+---------------------------------- F(3, 11) = 0.73
Model | 19.5429062 3 6.51430205 Prob > F = 0.5554
Residual | 98.1766922 11 8.92515383 R-squared = 0.1660
-------------+---------------------------------- Adj R-squared = -0.0614
Total | 117.719598 14 8.40854274 Root MSE = 2.9875
------------------------------------------------------------------------------
y | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
fdiffdum2 | 2.104393 1.439085 1.46 0.172 -1.063011 5.271798
fdiffdum3 | .5222648 1.377463 0.38 0.712 -2.509511 3.55404
fdiffdum4 | .6589209 1.266188 0.52 0.613 -2.127941 3.445783
_cons | .5054862 .829675 0.61 0.555 -1.320616 2.331589
------------------------------------------------------------------------------
在R，一个能做到这一点如下：
factor1mat <- model.matrix(~factor1, DF)
factor2mat <- model.matrix(~factor2, DF)
factordiffmat <- factor1mat - factor2mat
summary(lm(Y~factordiffmat, data=DF))
Coefficients: (1 not defined because of singularities)
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.5055 0.8297 0.609 0.555
factordiffmat(Intercept) NA NA NA NA
factordiffmatfactor1B 2.1044 1.4391 1.462 0.172
factordiffmatfactor1C 0.5223 1.3775 0.379 0.712
factordiffmatfactor1D 0.6589 1.2662 0.520 0.613
我不清楚你是否想要限制factor1和factor2上的_coefficients_等于零，或者_values_是否被约束为0 ... – MichaelChirico
我对这个问题的理解是'factor1'和'factor2'完全是多重共线性的。所以你只能包括一个或另一个...... – MichaelChirico
我不认为这是一个关于R/R代码的问题，因此我不认为这个问题在这方面是脱离主题的。 OTOH，我并不真正关注你的情况，或者它是如何激发你怀疑的解决方案的。对于这个问题，我不清楚你的建议解决方案*是什么（例如，我分享@ MichaelChirico的困惑）。它可能有助于开发一个简单的例子，只需几个级别和一个示例数据集就可以了，然后添加一些额外的解释。 – gung