本文转载自统计之都，文章翻译自   Jonas Kristoffer Lindeløv   的  Common statistical tests are linear models (or: how to teach stats) ，翻译工作已获得原作授权。本翻译工作首发于统计之都网站和微信公众号上。

7 比率：卡方检验是对数线性模型

回想一下，使用对数转换可以简单地处理 比率（proportion） ，举例来说， x  每一个单位的增加，都会带来  y  的一定数量的百分比的增加。借助这个特点，可以最简单有效地理解计数数据和列联表。

7.1 拟合优度检验

7.1.1  理论：作为对数线性模型

模型：一个单独截距来预测 log（y）。
我提议你参考阅读 列联表一节 ，基本上它就是 双因素拟合优度检验 。

7.1.2  示例数据

本例子中，需要一些宽格式的计数数据： 
# Data in long formatD <- data.frame(
  mood = c("happy", "sad", "meh"),
  counts = c(60, 90, 70)
)# Dummy coding for the linear model
D$mood_happy <- ifelse(D$mood == "happy", 1, 0)
D$mood_sad <- ifelse(D$mood == "sad", 1, 0)

7.1.3  R代码：拟合优度检验

现在，让我们展示一下，拟合优度检验其实是单因素方差分析的对数线性变换版本。设置 family = poisson() ，它默认使用 log 作为 联系函数（link function） （family = poisson(link='log')）。 
# Built-in test
a <- chisq.test(D$counts)

# As log-linear model, comparing to an intercept-only model
full <- glm(counts ~ 1 + mood_happy + mood_sad, data = D, family = poisson())
null <- glm(counts ~ 1, data = D, family = poisson())
b <- anova(null, full, test = "Rao")

# Note: glm can also do the dummy coding for you:
c <- glm(counts ~ mood, data = D, family = poisson())
来看看结果：

注意，其中  anova(..., test = 'Rao')  表示用 Rao 得分 检验，又称为拉格朗日乘子检验（Lagrange Multiplier test，LM test）来计算 p 值。当然也可以使用  test='Chisq'  或  test='LRT' ，它们计算近似的 p 值。你可能会认为我们在这里作弊了，偷偷地对卡方模型进行后续处理，实际上， anova  仅仅指定了 p 值的计算方式，内部对数线性模型仍然发生在  glm  中。
顺带一说，如果只有两个计数变量，而且样本量较大（ N > 100 ），模型会有效地近似于 二项检验（binomial test） binom.test 。这个样本量比通常情况要更大，所以我认为这不是经验准则，也不会在此进一步探索。

7.2 列联表

7.2.1  理论：作为对数线性模型

这里的理论会变得更令人费解，我会简单地写一下，从而你可以感受到它其实就是对数线性 双因素方差分析模型 。来，开始探索……
对于双因素列联表，计数变量 yy 的模型使用了列联表的边缘比率来建模。为什么这是可行的呢？答案比较高深，我们不会在这里详解，但读者可以通过查阅 Christoph Scheepers 的相关幻灯片来获得精彩的解答。这个模型包含了很多计数变量和回归系数 Ai 和 B i ：
yi=N⋅xi(Ai/N)⋅zj(Bj/N)⋅xij/((Aixi)/(Bjzj)/N)
多复杂呀！！！这里， i 是行标号， j 列标号，Xsomething是相应行和或列和的值：N=∑y。请记得， y 是计数变量，所以  N 是总计数值。
我们可以通过定义 比率 来简化以上记号：αi=xi(Ai/N)， βi =xj(Bi/N) ， αiβj =xij/(Aixi)/(Bjzj)/N 。重写模型如下：
yi=N⋅αi⋅βj⋅αiβj
嗯，好多了。然而，这里依然有很多乘法项，使得我们很难从直观上理解变量之间是如何交互的。如果还记得 log(A⋅B)=log(A)+log(B)， 那么两边取对数就清晰易懂了，可得：
log(yi)=log(N)+log(αi)+log(βj)+log(αiβj)

太爽了！现在我们可以直观地理解回归系数（都是比率）是怎样独立地影响到  y 的。这就是为什么对数变换对比率数据如此有效。注意到，这其实就是双因素方差分析模型加上一些对数变换，这就回到了所熟悉的线性模型———只是对系数的解释发生了变化而已！此外，我们不能继续使用 R 里的  lm  函数了。

7.2.2  示例数据

我们需要一些“长”格式的数据，并且需要保存为表格格式，才能作为 chisq.test 的输入： 

# Contingency data in long format for linear model
D <- data.frame(
  mood = c("happy", "happy", "meh", "meh", "sad", "sad"),
  sex = c("male", "female", "male", "female", "male", "female"),
  Freq = c(100, 70, 30, 32, 110, 120)
)

# ... and as table for chisq.test
D_table <- D %>%
  tidyr::spread(key = mood, value = Freq) %>% # Mood to columns
  select(-sex) %>% # Remove sex column
  as.matrix()

# Dummy coding of D for linear model (skipping mood=="sad" and gender=="female")
# We could also use model.matrix(D$Freq~D$mood*D$sex)
D$mood_happy <- ifelse(D$mood == "happy", 1, 0)
D$mood_meh <- ifelse(D$mood == "meh", 1, 0)
D$sex_male <- ifelse(D$sex == "male", 1, 0)

7.2.3  R代码：卡方检验

接下来看看卡方检验和对数线性模型之间的等价性。这个过程和上述 双因素方差分析 过程非常相似： 
# Built-in chi-square. It requires matrix format.
a <- chisq.test(D_table)
# Using glm to do a log-linear model, we get identical results when testing the interaction term:
full <- glm(Freq ~ 1 + mood_happy + mood_meh + sex_male + 
              mood_happy * sex_male + mood_meh * sex_male, data = D, family = poisson())
null <- glm(Freq ~ 1 + mood_happy + mood_meh + sex_male, data = D, family = poisson())
b <- anova(null, full, test = "Rao") # Could also use test='LRT' or test='Chisq'

# Note: let glm do the dummy coding for you
full <- glm(Freq ~ mood * sex, family = poisson(), data = D)
c <- anova(full, test = "Rao")

# Note: even simpler syntax using MASS:loglm ("log-linear model")
d <- MASS::loglm(Freq ~ mood + sex, D)
代码里用了 summary(full)，你可以取消以上代码的折叠，查看回归模型拟合的系数的原始值。作为对数线性模型，这些系数表示了：如果跳转到一个类别的话，y 将会获得多少 比例上的提升 。

8 资料来源和更多的等价模型

下面是本文内容的部分资料来源，还包含了很多本文没有提到的等价性模型：

  • Cross Validated 网站上， 我的原始想法 。

    https://stats.stackexchange.com/questions/303269/common-statistical-tests-as-linear-models

  • 对于“非参”检验， 我之前提出的疑问 和有用的答案。

    https://stats.stackexchange.com/questions/210529/are-parametric-tests-on-rank-transformed-data-equivalent-to-non-parametric-test

  • StackOverflow 网站上，关于 t 检验和方差分析的 问题和回答 。

    https://stats.stackexchange.com/questions/59047/how-are-regression-the-t-test-and-the-anova-all-versions-of-the-general-linear

  • Christoph Scheepers 的幻灯片 ，介绍了卡方检验如何被理解为对数线性模型。

    https://uni-tuebingen.de/fileadmin/Uni_Tuebingen/SFB/SFB_833/A_Bereich/A1/Christoph_Scheepers_-_Statistikworkshop.pdf

  • Philip M. Alday 的笔记，里面包括了卡方、二项、多项、泊松分布作为对数线性模型和 logistic 模型的理解。文中介绍的“等价性”没有我在本文展示的那么精确，因此我没有在本文详细介绍。然而，它们对理解这些检验是有帮助的！

    https://rpubs.com/palday/glm-test

  • Kristoffer Magnusson  的文章使用  lme4::lmer  混合模型（ mixed model ）介绍了  RM-ANOVA  和增长模型（ growth model ）。

    https://rpsychologist.com/r-guide-longitudinal-lme-lmer

  • Thom Baguley 的文章介绍了 Friedman 检验。这篇文章实际上启发了我开始思考“非参”检验的线性模型等价形式，而且最终推动我写下了本文章。

    https://seriousstats.wordpress.com/2012/02/14/friedman/

9 教材和课程大纲

大部分高等统计书籍（和一些入门书籍）也都同意“所有模型都是 GLMM（广义线性混合效应模型） 的观点”。然而，线性模型部分通常都是概念上提了一下，而没有清晰地指出细节。我想通过简练的方式把线性模型当作工具。幸运地，大部分对初学者友好的教材后来都合并了：

  • Russ Poldrack 的开源书籍 "Statistical Thinking for the 21st century"（从关于建模的第 5 章开始）

    http://statsthinking21.org/fitting-models-to-data.html

  • Jeff Rouder 的课程笔记，介绍了仅使用 R^2 和 BIC 来对比模型。它避开了所有关于 p 值、F 值等等的繁琐问题。完整的材料和幻灯片。

    https://jeffrouder.blogspot.com/2019/03/teaching-undergrad-stats-without-p-f-or.html

    https://drive.google.com/drive/folders/1CiJK--bAuO0F-ug3B5I3FvmsCdpPGZ03

我说一下对我所做的事情的看法。我已使用了本文的一部分进行教学，并获得了巨大的成功，但是这并不是完整的教学过程，因为我并没有分派到教授整个课程。
我会花费 50% 的时间在数据的线性模型上，因为它包含了学生所需知道的 70%（以下的第 1 点）。剩下来的课程则是关于当你有一个组、两个组等等数据的时候会发生什么事情。
注意，主流统计课程的开始部分都是关于采样和假设检验的理解；我这里把这部分移动到后面，这样，学生可以基于之前学习的知识来进行理解，而不是一上来就面对各种陌生的概念。

  1. 回归的基础知识

    1. 回想高中的知识： 然后获得对斜率和截距的非常好的直觉。理解到这条式子能用所有的变量名称来重写：如 money = profit * time + starting_money或   或去除系数之后可写成 y ~ x + 1。如果听众接受程度高的话，可以探索这些模型是如何解微分方程的，并指出 y 是如何随着 x 的变化而变化的。

    2. 扩展到多元回归模型。记得这时候要带有非常多的生活例子和练习，从而使这些概念变得直觉上非常容易理解。让听众感叹于这些简洁的模型都可以用来描述非常大的数据集。

    3. 介绍对于非数值型数据如何进行秩转换，并进行各种尝试。

    4. 教授三种前提假设：数据点的独立性，残差分布的正态性和方差齐性 （homoscedasticity）。

    5. 参数的置信（confidence）/可信（credible）区间。指出极大似然估计（Maximum-Likelihood estimate）很难计算，因此区间估计更为重要。

    6. 对以上简单的回归模型，简要地介绍 R^2。顺便提及一下，这就是 Pearson 和 Spearman 相关系数。

  2. 特殊情况 #1：一个或两个均值（t 检验、Wilcoxon 检验、Mann-Whitney 检验）：

    1. 单均值： 当只有一个 x 值的时候，回归模型简化成了 y = b。如果 y 不是数值型的，你可以进行秩转换。应用模型假设（只有一个 x，因此方差齐性不适用于这里）。顺便提及一下，这些仅有截距的模型也分别可称为单样本 t 检验和 Wilcoxon 符号秩检验。

    2. 双均值： 如果我们把两个变量一起放在 x 轴，两者均值之差就是斜率。很好！这就能用我们称为瑞士军刀的线性模型来解决。应用模型的假设条件，检查两个组的方差是否相等，相等即方差齐性。这模型称为独立 t 检验。构造一些例子，做一些练习，也许还能加上 Welch 检验，再加上秩转换 ---- 变成所谓的 Mann-Whitney U 检验的版本。

    3. 配对样本： 违反了独立性假设。通过计算配对组的差值，这就转化成了 2.1（单截距）的等价形式，尽管这种情况有另外的名称：配对 t 检验和 Wilcoxon 配对组检验。

  3. 特殊情况 #2：三个或多个均值（方差分析（ANOVA））

    1. 对类别转化后的 示性变量 ： 类别的每一个取值范围对应的回归系数，是如何通过乘以一个二元（binary）示性变量，来对每个取值范围对应的截距来进行建模的。 （ How one regression coefficient for each level of a factor models an intercept for each level when multiplied by a binary indicator.） 这只是我们为了使数据能用线性模型建模，而扩展了在 2.1 所做的事情而已。

    2. 一个变量的均值： 单因素方差分析（one-way ANOVA） .

    3. 两个变量的均值： 双因素方差分析（two-way ANOVA） .

  4. 特殊情况 #3：三个或多个比率（卡方检验）

    1. 对数变换： 通过对数变换，把“多元乘法”模型转化成线性模型，从而可以对比率进行建模。关于对数线性模型和对比率的卡方检验的等价性，可以查阅这个非常优秀的介绍。此外，还需要介绍 (log-) odds ratio（一般翻译为“比值比”或“优势比”）。当“多元乘法”模型使用对数变换转化为“加法”模型之后，我们仅加上来自 3.1 的示性变量技巧，就会在接下来发现模型等价于 3.2 和 3.3 的方差分析----除了系数的解释发生了变化。

    2. 单变量的比率： 拟合优度检验.

    3. 双变量的比率： 列联表.

  5. 假设检验：

    1. 视为模型比较的假设检验： 假设检验用于全模型和某个参数固定了（通常为 0，也即从模型中去除）的模型进行比较，而不是对模型进行估计。比如说，在 t 检验  把两个均值之一固定为零之后，我们探究单独一个均值（单样本 t 检验）对两个组的数据的解释程度。如果解释程度比较好，那么我们更倾向于这个单均值模型，而不是双均值模型，因为前者更为简单。假设检验其实是比较多个线性模型，来获得更多的定量描述。单参数的检验，假设检验包含的信息更少。但是，同时对多个参数（如方差分析的类别变量）进行检验的话，模型比较就会变得没有价值了。

    2. 似然比： 似然比是一把瑞士军刀，它适用于单样本 t 检验到 GLMM 等情况。BIC 对模型复杂度进行惩罚。还有，加上先验（prior）的话，你就能得到贝叶斯因子（Bayes Factor）。一个工具，就能解决所有问题。我在上文方差分析中使用了似然比检验。

10 不足之处

一些需要澄清的简化前提：

  1. 我没在这里覆盖到前提假设的内容。这会在另一篇文章揭晓！但是所有检验都很可能有三个预定假设：a)  数据点的独立性， b)  残差的正态性， c)  同方差性（homoscedasticity）。

  2. 我假定所有的零假设是缺失了效应的情况，但是所有原理都和非 0 的零假设所一致的。

  3. 我没有讨论推断内容。因为大家都会关心 p 值，因此我在比较中提到了 p 值，从而简短地展示了背后的模型等价性。参数的估计值也会展示出相同的等价性。如何进行推断则是另一个话题了。我个人是贝叶斯学派的，但是展示贝叶斯学派内容的话，会减少这篇文章的受众。此外，构造 稳健模型 是更好的选择，但是它无法揭示模型的等价性。

  4. 本文列表依然缺失了很多其它著名的检验，有可能在以后添加进来。比如说符号检验（sign test）（要求很大的 N 从而可以有效地使用线性模型来近似），Friedman 检验 -- 即在 rank(y) 上的 RM-ANOVA，McNemar 检验，和二项（Binomial）/多项（Multinomial）检验。如果你认为它们需要在本文提及到，欢迎在本文档的 Github 仓库  提交对应说明！

Github：https://github.com/lindeloev/tests-as-linear/

11 附录：翻译稿评论

11.1  译者评论

译者：相对于统计检验来说，线性回归实际上是更符合直觉的。想当年某检验实在让笔者百思不得其解，某师姐指点迷津：“你实在搞不懂可以看成是线性回归对系数的检验，我们如此这般建造一个 X ……”让笔者茅塞顿开。故听朋友推荐本文之后，笔者毛遂自荐承接了翻译任务。希望各位读者能从本文感受统计的威力和它带来的喜悦。如各位读者有指正或建议之处，热烈欢迎于主站或微信文下留言评论。

11.2  审稿人讨论

审稿人：我（黄湘云）看完这篇文章的感受是怀疑自己读了个“假”大学，开个玩笑哈！感觉这篇文章是继《心理学的危机》后又一篇需要找个地方安安静静读几天的，文中很多检验的细节都略过了，更加数学严格的检验介绍估计得去看《数理统计引论》陈希孺著的这本书才能明白。这篇文章的覆盖面起码是一个学期的课，如果把本文没有详述的其他检验补充进来，特别是加上检验的数学推导和一些实际应用案例后，估计能成一本书。我是学线性模型的（此线性模型非大多数人了解的彼线性模型），看完之后有点汗颜和如梦初醒，惊奇于作者独具一格的视角。不足之处是有些地方还不够通俗，比如列联表作为对数线性模型来理解，一点也不直接，作者也略去了！这里的列联表其实是指我们通常教科书上的独立性检验。拟合优度检验和独立性检验的检验统计量的极限分布都是卡方分布，故而都归纳在卡方检验下。
文章介绍了那么多的检验问题，实际上都可以归为统计学三大检验 --- 似然比检验、 Wald 检验、（Rao）Score（得分）检验 --- 在线性模型下的特殊情况。数理统计的教材往往是利用似然比这把瑞士军刀展开介绍的。似然比在假设检验中地位相当于极大似然估计在参数估计中的地位，相当于正态分布在抽样分布中的地位。抽样分布、参数估计和假设检验合称统计推断。学数理统计的人往往不愿去记那么多的检验名称，比如 t 检验、F 检验和卡方检验，特别是诸多名人检验，因为本质上那只是似然比统计量在不同的条件下呈现的极限分布不同而已。三大检验的渐近等价性可参考

  • Dennis D. Boos and L. A. Stefanski.  Essential Statistical Inference: Theory and Methods. 2013. Springer, New York. 

    <10> Chapter 3. Likelihood-Based Tests and Confidence Regions. 10>

  • Vito M. R. Muggeo and Gianfranco Lovison (2014) The “Three Plus One” Likelihood-Based Test Statistics: Unified Geometrical and Graphical Interpretations,  The American Statistician,  68:4, 302-306,  DOI: 10.1080/00031305.2014.955212

Rao 在 1948 年给出得分检验的渐近性，文中提及 Rao 得分检验这个名称只是强调 Rao 在得分检验中的贡献，有些书籍中提及 Rao 检验基本等同于得分检验，而拉格朗日乘子检验由 Aitchison 和 Silvey 于 1958 年在经济学领域独立提出来的，所以在统计学文献中见最多的是得分检验，经济学文献中多描述为拉格朗日乘子检验。这二者是殊途同归，都对检验统计量的得分函数做泰勒展开，其极限分布都为卡方分布。
列联表是分类数据的一种方便的组织形式，与之相关的检验和前面带连续变量的线性回归模型的检验是有本质区别的，列联表是与多项分布联系起来的，这里没有残差，线性回归模型往往对残差做了独立同正态分布的假设。
译者简介
黄俊文 ，本科中山大学，硕士 University of California, San Diego。现于业界从事数据分析工作。
作者： Jonas Kristoffer Lindeløv
审校： 蔡占锐，黄湘云，谢益辉
编辑： 孙腾飞，任焱

统计之都官网： cosx.org

统计之都公众号： CapStat

相关文章