方差分析要求各比较组除了处理因素不同外，其他对观察指标有影响的因素应齐同或均衡，即要求控制对观察指标有影响的其他因素。但在实际工作中，有时有些因素无法加以控制或因为某些原因未加控制，例如要评价三种降压药物的降压效果，患者的初始血压水平对服药一段时间后血压下降量有影响，但是患者的初始血压是难以控制的，这时候如果不考虑初始血压的差异直接用方差分析是不恰当的。这时就需要用到协方差分析，初始血压对血压下降量的影响可以看作协变量。协方差分析就是扣除协变量影响后，比较不同组均数之间的差别。
协方差分析要满足两个条件：一是和方差分析一样，理论上要求观察变量服从正态分布且互相独立，各样本方差齐；二是各总体客观存在因变量对协变量的线性回归关系且斜率相等。
我们采用multcomp包中的litter数据集。某研究员将怀孕小鼠随机分成4组，分别接受不同剂量（0、5、50、500）的药物处理，称量各小鼠幼崽的体重数据如表1，变量分别为药物剂量（dose），幼崽体重（weight），怀孕时间（gesttime）。问各组怀孕小鼠产下的幼崽的体重是否有差异。
我们知道小鼠的怀孕时间对幼崽体重影响很大，但我们无法控制小鼠的怀孕时间，因此怀孕时间就是一个协变量，在分析药物对幼崽体重的影响时，我们要去除怀孕时间的影响。
表1 数据资料
| id | dose |
weight |
gesttime |
1 |
0 |
28.05 |
22.5 |
2 |
0 |
33.33 |
22.5 |
3 |
0 |
35.80 |
21.5 |
… |
… |
… |
… |
11 |
500 |
24.55 |
22.0 |
12 |
500 |
33.78 |
22.5 |
经正态性和方差齐性检验，数据符合正态性和方差齐性，更多正态性和方差齐性检验方法可见正态性检验和t检验。除此之外，协方差分析还要假定回归斜率相等。
在本例中，我们假定4个处理组通过怀孕时间来预测出生体重的回归斜率都相同。可用ANCOVA模型中gesttime*does的交互项对回归斜率的同质性进行检验。若交互项显著，则表示时间和幼崽出生体重的关系依赖于药物剂量的水平，代码和结果如下：
library(multcomp) #载人"multcomp包"
data(litter) #引用包中"litter"数据集
attach(litter) #将"litter"添加到R的搜索路径中
library(car) #载入"car"包
qqPlot(lm(weight~dose,data=litter),simulate=TRUE,main="Q-Q Plot") #正态性检验
bartlett.test(weight~dose,data=litter) #方差齐性检验
library(multcomp)
fit1<-aov(weight~gesttime*dose,data=litter) #回归斜率的同质性检验
summary(fit1)
图1 回归斜率同质性检验结果
由图1可以看出，怀孕时间和剂量的交互效应不显著，支持了斜率相等的假设，说明药物与怀孕时间对幼崽体重的影响是独立的。接下来就可以进行协方差分析了。
数据满足正态性、方差齐性和回归斜率同质性，接下来我们做协方差分析，代码和结果如下：
fit2<-aov(weight~gesttime+dose,data=litter) #协方差分析
summary(fit)
图2 协方差分析结果
协方差分析的F检验结果显示，怀孕时间与幼崽出生体重相关，控制怀孕时间，药物剂量与出生体重相关。
协方差分析比较的是修正均数，即去除了协变量效应后的组均值，可以使用effects包中的effect()函数计算控制了怀孕时间影响后，不同药物剂量组的出生体重的均值。代码和结果如下：
install.packages("effects")
library(effects)
effect("dose",fit2)
图3 调整组均值
图3显示了各组去除怀孕时间影响后的均值，我们可以用调整后的均值进行多重比较，分析不同药物剂量水平的差异，具体方法可见《组间差异的的多重比较》。
小结
协方差分析是将线性回归分析和方差分析结合起来的一种统计分析方法，需要注意的是，调整后的组均值实际上是假设协变量固定在其总体均数时的观察变量的均值，其差别与实际均数间的差别并不是一回事，只是用于合理比较。
参考来源
1. 孙振球,徐勇勇.医学统计学:第4版[M].北京:人民卫生出版社.2014.
2. [美]Robert I. Kabacoff著．R语言实战(第2版) [M].王小宁等译．北京:人民邮电出出版社.2016.
制作：周燕、蔡敏
初审：何冠豪、胡建雄
审核：肖建鹏、刘涛
指导：马文军
《数据分析和应用》致力于为全国各地公共卫生与医学工作者（机构）提供专业可靠的统计咨询、研究设计、数据分析、高通量测序数据和序列分析、调研报告等服务（详细可见公众号菜单栏），欢迎有需要的人员和机构与我们联系。
