如何解释二维列联表？

机器算法验证 r 分类数据 解释 列联表 对数线性
2022-03-23 21:26:29

我试图了解如何解释通过泊松 GLM 拟合的列联表的对数线性模型。

考虑 CAR 中的这个例子（Fox 和 Weisberg，2011，第 252 页）。

require(car)
data(AMSsurvey)
(tab.sex.citizen <- xtabs(count ~ sex + citizen, data=AMSsurvey))

产量：

        citizen
sex      Non-US  US
  Female    260 202
  Male      501 467

然后我们拟合（相互）独立的模型：

AMS2 <- as.data.frame(tab.sex.citizen)
(phd.mod.indep <- glm(Freq ~ sex + citizen, family=poisson, data=AMS2))
pchisq(2.57, df=1, lower.tail=FALSE)

输出：

> (phd.mod.indep <- glm(Freq ~ sex + citizen, family=poisson, data=AMS2))

Call:  glm(formula = Freq ~ sex + citizen, family = poisson, data = AMS2)

Coefficients:
(Intercept)      sexMale    citizenUS  
     5.5048       0.7397      -0.1288  

Degrees of Freedom: 3 Total (i.e. Null);  1 Residual
Null Deviance:      191.5 
Residual Deviance: 2.572    AIC: 39.16
> pchisq(2.57, df=1, lower.tail=FALSE)
[1] 0.1089077

p 值接近 0.1，表明拒绝独立性的证据不足。但是，让我们假设我们有足够的证据拒绝 NULL（即，对于我们的目的，0.10 p 值表示两个变量之间的关联）。

问题：那么，我们如何解释这个对数线性模型？

（我们是否适合饱和模型（即update(phd.mod.indep, . ~ . + sex:citizen)）？我们是否解释了估计的回归系数？在 CAR 中，由于拒绝 NULL 的证据不足，他们停在这一点上，但我有兴趣了解解释这一点的机制简单的对数线性模型，好像“交互”很重要......）

2个回答

一般来说，2 路列联表没有太多内容，但是您正试图以如此详细的级别解开它，以至于出现了一些混淆。通常，对于像这样的简单列联表，人们只想知道变量 (sex和citizen) 是否独立。要评估这一点，您可以运行卡方检验： 

chisq.test(tab.sex.citizen)
#  Pearson's Chi-squared test with Yates' continuity correction
# 
# data:  tab.sex.citizen
# X-squared = 2.389, df = 1, p-value = 0.1222

您可以通过执行针对完整（饱和）模型运行的 Poisson GLM 的嵌套模型测试来执行此测试的似然比版本（上面的卡方测试是分数测试）： 

phd.mod.indep <- glm(Freq ~ sex + citizen, family=poisson, data=AMS2)
phd.mod.sat   <- glm(Freq ~ sex * citizen, family=poisson, data=AMS2)
anova(phd.mod.indep, phd.mod.sat, test="LRT")
# Analysis of Deviance Table
# 
# Model 1: Freq ~ sex + citizen
# Model 2: Freq ~ sex * citizen
#   Resid. Df Resid. Dev Df Deviance Pr(>Chi)
# 1         1     2.5721                     
# 2         0     0.0000  1   2.5721   0.1088
deviance(phd.mod.indep)                                        # [1] 2.572123
deviance(phd.mod.sat)                                          # [1] 3.308465e-14
1-pchisq(deviance(phd.mod.indep)-deviance(phd.mod.sat), df=1)  # [1] 0.1087617

还可以得到交互项的 Wald 检验（即独立性检验）： 

summary(phd.mod.sat)
# ...
# Coefficients:
#                   Estimate Std. Error z value Pr(>|z|)    
# (Intercept)        5.56068    0.06202  89.663  < 2e-16 ***
# sexMale            0.65592    0.07643   8.582  < 2e-16 ***
# citizenUS         -0.25241    0.09379  -2.691  0.00712 ** 
# sexMale:citizenUS  0.18214    0.11373   1.602  0.10926    
# ...
#     Null deviance: 1.9148e+02  on 3  degrees of freedom
# Residual deviance: 3.3085e-14  on 0  degrees of freedom
# AIC: 38.586
# ...

（要了解分数与 Wald 与似然比检验的相关信息，请在此处查看我的答案：为什么我的 p 值在逻辑回归输出、卡方检验和 OR 的置信区间之间存在差异？）

但是请注意，您进行测试的phd.mod.indep方式不正确（请参阅此处：使用 null 和模型偏差测试 GLM 模型）。该模型对空模型的检验是检验所有细胞在总体中是否具有相同的概率。它将按如下方式实施： 

1-pchisq(191.5-2.57, df=3-1)  # [1] 0

撇开所有细胞概率是否相等的测试（我怀疑这对您很感兴趣），如果您认为sex和之间存在关联citizen，您将不会解释模型phd.mod.indep。那将是一个错误指定的模型。

除了@gung 的回答，我发现马赛克图在评估与分类数据的相对关系方面非常有用。一旦我们建立（实际上，我们在这里假设）背离独立性，就可以使用包vcd以图形方式显示此类列联表。

马赛克图：

mosaic(tab.sex.citizen, shade=TRUE, legend=TRUE)

在此处输入图像描述

关联图：

assoc(tab.sex.citizen, shade=TRUE, legend=TRUE)

在此处输入图像描述

从第二个情节可以看出，非公民的女性比例明显高于美国公民。同样，女性的非公民比例高于男性。

当脱离独立性明确时，图表会更清晰，所以这里有一个关于泰坦尼克号数据的示例（请参阅本教程）。

data(Titanic)
sex.survived <- margin.table(Titanic, c(2,4))
ftable(sex.survived)

产量：

> ftable(sex.survived)
       Survived   No  Yes
Sex                      
Male            1364  367
Female           126  344

然后是情节：

mosaic(sex.survived, shade=TRUE, legend=TRUE)
assoc(sex.survived, shade=TRUE, legend=TRUE)

在此处输入图像描述

在此处输入图像描述

从皮尔逊残差（观察到的频率与预期频率的偏差的度量，或数据中未由对数线性模型解释的位），很明显，对于男性来说，有一个不成比例的高数字没有幸存下来，并且幸存下来的人数极少。同样，在幸存者中，女性的数量出乎意料地高。

其它你可能感兴趣的问题
上一篇傅里叶变换和多元法线 下一篇使用 Hedonic 回归对价格进行建模