【问题标题】：GLM submodel testing in R: why all statistics are still the same after remove one continuous covariate?R中的GLM子模型测试：为什么在删除一个连续协变量后所有统计数据仍然相同？
【发布时间】：2017-12-14 13:25:46
【问题描述】：

我正在进行子模型测试。较小的模型嵌套在较大的模型中。与较小的模型相比，较大的模型具有一个连续变量。我使用似然比检验。结果相当奇怪。两种模型具有相同的统计数据，例如残差偏差和 df。我还发现两个模型具有相同的估计系数是 std.errors。事实怎么可能？

summary(m2221)

Call:
glm(formula = clm ~ veh_age + veh_body + agecat + veh_value:veh_age + 
veh_value:area, family = "binomial", data = Car)

Deviance Residuals: 
Min       1Q   Median       3Q      Max  
-0.9245  -0.3939  -0.3683  -0.3437   2.7095  

Coefficients:
                    Estimate Std. Error z value Pr(>|z|)    
(Intercept)        -1.294118   0.382755  -3.381 0.000722 ***
veh_age2            0.051790   0.098463   0.526 0.598897    
veh_age3           -0.166801   0.094789  -1.760 0.078457 .  
veh_age4           -0.239862   0.096154  -2.495 0.012611 *  
veh_bodyCONVT      -2.184124   0.707884  -3.085 0.002033 ** 
veh_bodyCOUPE      -0.850675   0.393625  -2.161 0.030685 *  
veh_bodyHBACK      -1.105087   0.374134  -2.954 0.003140 ** 
veh_bodyHDTOP      -0.973472   0.383404  -2.539 0.011116 *  
veh_bodyMCARA      -0.649036   0.469407  -1.383 0.166765    
veh_bodyMIBUS      -1.295135   0.404691  -3.200 0.001373 ** 
veh_bodyPANVN      -0.903032   0.395295  -2.284 0.022345 *  
veh_bodyRDSTR      -1.108488   0.826541  -1.341 0.179883    
veh_bodySEDAN      -1.097931   0.373578  -2.939 0.003293 ** 
veh_bodySTNWG      -1.129122   0.373713  -3.021 0.002516 ** 
veh_bodyTRUCK      -1.156099   0.384088  -3.010 0.002613 ** 
veh_bodyUTE        -1.343958   0.377653  -3.559 0.000373 ***
agecat2            -0.198002   0.058382  -3.391 0.000695 ***
agecat3            -0.224492   0.056905  -3.945 7.98e-05 ***
agecat4            -0.253377   0.056774  -4.463 8.09e-06 ***
agecat5            -0.441906   0.063227  -6.989 2.76e-12 ***
agecat6            -0.447231   0.072292  -6.186 6.15e-10 ***
veh_age1:veh_value -0.000637   0.026387  -0.024 0.980740    
veh_age2:veh_value  0.035386   0.031465   1.125 0.260753    
veh_age3:veh_value  0.114485   0.036690   3.120 0.001806 ** 
veh_age4:veh_value  0.189866   0.057573   3.298 0.000974 ***
veh_value:areaB     0.044099   0.021550   2.046 0.040722 *  
veh_value:areaC     0.021892   0.019189   1.141 0.253931    
veh_value:areaD    -0.023616   0.024939  -0.947 0.343658    
veh_value:areaE    -0.013506   0.026886  -0.502 0.615415    
veh_value:areaF     0.057780   0.026602   2.172 0.029850 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 33767  on 67855  degrees of freedom
Residual deviance: 33592  on 67826  degrees of freedom
AIC: 33652

Number of Fisher Scoring iterations: 5

 summary(m222)

Call:
glm(formula = clm ~ veh_value + veh_age + veh_body + agecat + 
    veh_value:veh_age + veh_value:area, family = "binomial", 
    data = Car)

  Deviance Residuals: 
  Min       1Q   Median       3Q      Max  
-0.9245  -0.3939  -0.3683  -0.3437   2.7095  

 Coefficients:
                    Estimate Std. Error z value Pr(>|z|)    
(Intercept)        -1.294118   0.382755  -3.381 0.000722 ***
veh_value          -0.000637   0.026387  -0.024 0.980740    
veh_age2            0.051790   0.098463   0.526 0.598897    
veh_age3           -0.166801   0.094789  -1.760 0.078457 .  
veh_age4           -0.239862   0.096154  -2.495 0.012611 *  
veh_bodyCONVT      -2.184124   0.707884  -3.085 0.002033 ** 
veh_bodyCOUPE      -0.850675   0.393625  -2.161 0.030685 *  
veh_bodyHBACK      -1.105087   0.374134  -2.954 0.003140 ** 
veh_bodyHDTOP      -0.973472   0.383404  -2.539 0.011116 *  
veh_bodyMCARA      -0.649036   0.469407  -1.383 0.166765    
veh_bodyMIBUS      -1.295135   0.404691  -3.200 0.001373 ** 
veh_bodyPANVN      -0.903032   0.395295  -2.284 0.022345 *  
veh_bodyRDSTR      -1.108488   0.826541  -1.341 0.179883    
veh_bodySEDAN      -1.097931   0.373578  -2.939 0.003293 ** 
veh_bodySTNWG      -1.129122   0.373713  -3.021 0.002516 ** 
veh_bodyTRUCK      -1.156099   0.384088  -3.010 0.002613 ** 
veh_bodyUTE        -1.343958   0.377653  -3.559 0.000373 ***
agecat2            -0.198002   0.058382  -3.391 0.000695 ***
agecat3            -0.224492   0.056905  -3.945 7.98e-05 ***
agecat4            -0.253377   0.056774  -4.463 8.09e-06 ***
agecat5            -0.441906   0.063227  -6.989 2.76e-12 ***
agecat6            -0.447231   0.072292  -6.186 6.15e-10 ***
veh_value:veh_age2  0.036023   0.034997   1.029 0.303331    
veh_value:veh_age3  0.115122   0.039476   2.916 0.003543 ** 
veh_value:veh_age4  0.190503   0.058691   3.246 0.001171 ** 
veh_value:areaB     0.044099   0.021550   2.046 0.040722 *  
veh_value:areaC     0.021892   0.019189   1.141 0.253931    
veh_value:areaD    -0.023616   0.024939  -0.947 0.343658    
veh_value:areaE    -0.013506   0.026886  -0.502 0.615415    
veh_value:areaF     0.057780   0.026602   2.172 0.029850 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 33767  on 67855  degrees of freedom
Residual deviance: 33592  on 67826  degrees of freedom
AIC: 33652

anova(m2221,m222,  test ="LRT")###
Analysis of Deviance Table

Model 1: clm ~ veh_age + veh_body + agecat + veh_value:veh_age + 
veh_value:area
Model 2: clm ~ veh_value + veh_age + veh_body + agecat + veh_value:veh_age + 
veh_value:area
Resid. Df Resid. Dev Df Deviance Pr(>Chi)
1     67826      33592                     
2     67826      33592  0        0

【问题讨论】：

您能否更新您的帖子以包含summary(m221) 和summary(m22) 的输出？

标签： r glm

【解决方案1】：

您以两种不同的方式指定了同一个模型。为了说明这一点，我将首先解释一下幕后发生的事情，然后我将遍历您的两个模型的系数以表明它们是相同的，最后我将以一些更高层次的直觉和解释结束。

不同的公式创建不同的交互项

首先，您的两个模型之间的唯一区别是第一个模型包含 veh_value 作为预测变量，而第二个模型不包含。但是，veh_value 与两个模型中的其他预测器交互。

让我们考虑一个简单的可重现的例子，看看当我们这样做时 R 做了什么。我将使用model.matrix 来简单地输入两个不同的公式，并查看 R 作为结果创建的连续变量和虚拟变量。

colnames(model.matrix(am ~ mpg + factor(cyl):mpg, mtcars))
#[1] "(Intercept)"      "mpg"              "mpg:factor(cyl)6" "mpg:factor(cyl)8"

colnames(model.matrix(am ~ factor(cyl):mpg, mtcars))
#"(Intercept)"      "factor(cyl)4:mpg" "factor(cyl)6:mpg" "factor(cyl)8:mpg"

请注意，在第一次调用中我包含了连续预测器 mpg，而我没有在第二次调用中包含它（您正在做的一个更简单的示例）。

现在请注意，第二个模型矩阵在第一个模型中包含一个额外的交互变量 (factor(cyl)4:mpg)。换句话说，因为我们没有直接在模型中包含mpg，所以cyl的所有级别都包含在交互中。

你的模型是一样的

您的模型与上面的简单示例基本上做同样的事情，我们可以证明一天结束时的系数在加在一起时实际上是相同的。

在您的第一个模型中，当 veh_value 被包含为与 veh_value 的交互但 veh_value 不包含在模型中时，所有 4 个级别的 veh_age 都包含在内。

veh_age1:veh_value -0.000637   0.026387  -0.024 0.980740    
veh_age2:veh_value  0.035386   0.031465   1.125 0.260753    
veh_age3:veh_value  0.114485   0.036690   3.120 0.001806 ** 
veh_age4:veh_value  0.189866   0.057573   3.298 0.000974 ***

在您的第二个模型中，当与 veh_value 交互时，仅包含 3 个级别的 veh_age，因为模型中包含 veh_value。

veh_value:veh_age2  0.036023   0.034997   1.029 0.303331    
veh_value:veh_age3  0.115122   0.039476   2.916 0.003543 ** 
veh_value:veh_age4  0.190503   0.058691   3.246 0.001171 **

现在，看到模型实际上是相同的关键部分。只需遍历veh_age 的所有级别即可轻松展示。

首先考虑veh_age = 1

对于这两个模型，veh_value 上的系数以 veh_age 为条件，而 veh_age = 1 为 -0.000637

# For first model
veh_age1:veh_value -0.000637   0.026387  -0.024 0.980740    

# For second model
veh_value          -0.000637   0.026387  -0.024 0.980740

现在考虑veh_age = 2

对于这两种模型，veh_value 上的系数以 veh_age 为条件，而 veh_age = 2 为 0.035386

# For first model
veh_age2:veh_value  0.035386   0.031465   1.125 0.260753    

# For second model - note sum of the two below is 0.035386
veh_value          -0.000637   0.026387  -0.024 0.980740
veh_value:veh_age2  0.036023   0.034997   1.029 0.303331

直觉

当您包含交互 veh_value:veh_age 时，您实际上是在说您希望 veh_value 一个连续变量的系数以 veh_age 一个分类变量为条件。包括交互veh_value:veh_age 和veh_value 作为预测变量是在说同样的事情。您想知道veh_value 的系数，但您想根据veh_age 的值对其进行调节。

【讨论】：