R 中的 GLM 子模型测试：为什么在删除一个连续协变量后所有统计数据仍然相同？
GLM submodel testing in R: why all statistics are still the same after remove one continuous covariate?
我正在进行子模型测试。较小的模型嵌套在较大的模型中。与较小的模型相比，较大的模型有一个连续变量。我使用似然比检验。结果很奇怪。两个模型具有相同的统计信息，例如残差和 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
您以两种不同的方式指定了同一个模型。为了展示这一点，我将首先解释一下引擎盖下发生的事情，然后我将遍历你的两个模型的系数以表明它们是相同的，最后我将以一些更高层次的直觉和解释结束。
创建不同交互项的不同公式
首先，您的两个模型之间的唯一区别是第一个模型包含 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_age 作为与 veh_value 的交互包含时，所有 4 个级别都包含在内，但 veh_value 未包含在模型中。
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_age = 1时veh_age条件下veh_value的系数是-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.
的值来调节它
我正在进行子模型测试。较小的模型嵌套在较大的模型中。与较小的模型相比，较大的模型有一个连续变量。我使用似然比检验。结果很奇怪。两个模型具有相同的统计信息，例如残差和 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
您以两种不同的方式指定了同一个模型。为了展示这一点，我将首先解释一下引擎盖下发生的事情，然后我将遍历你的两个模型的系数以表明它们是相同的，最后我将以一些更高层次的直觉和解释结束。
创建不同交互项的不同公式
首先，您的两个模型之间的唯一区别是第一个模型包含 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_age 作为与 veh_value 的交互包含时，所有 4 个级别都包含在内，但 veh_value 未包含在模型中。
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_age = 1时veh_age条件下veh_value的系数是-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.