在前面文章中介绍了单因素协方差分析(One-Way ANCOVA)的假设检验理论，本篇文章将实例演示在Stata软件中实现单因素协方差分析的操作步骤。
关键词：Stata; 协方差分析; 单因素协方差分析; 平行性检验; 平行线检验
一、案例介绍
为研究A、B、C三种饲料对增加小鼠体重的影响，将初始体重相近的45只小鼠随机分成三组，分别喂养A、B、C三种饲料，但在实验设计时未对小鼠的进食量加以限制，现测得三组小鼠的进食量(Food)和所增体重(Weight)，请推断A、B、C三种饲料对小鼠的增重效果是否有差别？部分数据见图1。本文案例可从“附件下载”处下载。
二、问题分析
本案例的分析目的是比较A、B、C三种饲料对增加小鼠体重的影响，即A、B、C三组体重增量是否存在差异。但显然，每个老鼠的进食量会对体重增量产生影响，因此针对这种情况可将进食量作为体重增量的影响因素进行单因素协方差分析(One-Way ANCOVA)。但需要满足9个条件：
条件1：观察变量为连续变量。本研究中观察变量为体重增量，为连续变量，该条件满足。
条件2：自变量为二分类或多分类变量。本研究中自变量为A、B、C三组，该条件满足。
条件3：协变量是连续变量。本研究中协变量为进食量，为连续变量，该条件满足。
条件4：各研究对象之间具有独立的观测值。本研究中各个研究对象均为独立样本，不存在互相干扰的情况，该条件满足。
条件5：观察变量不存在显著异常值，该条件需要通过软件分析后判断。
条件6：各组内因变量的残差服从正态(或近似正态)分布。
条件7：组间因变量的残差方差齐，该条件需要通过软件分析后判断。
条件8：各组内协变量和因变量之间存在线性关系，该条件需要通过软件分析后判断。
条件9：各组内协变量和因变量的回归直线平行，即通过平行性检验，该条件需要通过软件分析后判断。
三、软件操作及结果解读
(一) 适用条件判断
条件5—条件7需要使用残差进行判断，首先使用单因素协方差分析生成因变量残差。
1. 生成因变量残差
(1) 软件操作
anova weight group c.food predict e, residual
(2) 结果解读
通过“Data Editor”视图可以看到生成了一列新的变量“e (残差) ”(图2)。
2. 条件5判断(异常值判断)
(1) 软件操作
绘制残差箱线图。
graph box e, over(group)
(2) 结果解读
图3残差的箱线图未提示任何异常值和极端值，满足条件5。
3. 条件6判断(正态性检验)
(1) 软件操作
单因素协方差分析需要检验因变量残差的整体正态性，结果如图4所示。
swilk e
(2) 结果解读
图4显示因变量残差正态性检验的P值为0.682，>0.1，可认为残差服从正态分布。此外，本案例也可以绘制残差的正态分位图，结果也提示残差满足正态性要求(请读者自行操作)。综上，本案例满足条件6。
4. 条件7判断(方差齐性检验)
(1) 软件操作
robvar e, by(group)
(2) 结果解读
由图5可知，F=2.063，P=0.140，提示本案例组间因变量残差的方差齐，满足条件7。
5. 条件8判断(线性关系检验)
(1) 软件操作
①绘制拟合图，结果如图6所示。
scatter weight food || lfit weight food ||, by(group,total row(1))
②计算各组内weight和food的Pearson相关系数，结果如图7所示。
bysort group: pwcorr weight food, sig
(2) 结果解读
由图6可知，1、2、3三组内food与weight均呈现线性关系；图7给出了相关性检验结果，三组的相关系数分别为0.844、0.836、0.849，P值均<0.001。综上，本案例各组内协变量和因变量之间均存在线性关系，满足条件8。
6. 条件9判断(平行性检验)
平行性检验通过判断group与food的交互项是否有统计学意义决定。
(1) 软件操作
anova weight group c.food group#c.food
(2) 结果解读
在图8所示结果中，模型整体的FModel=81.72，P <0.001，表示整体模型检验具有统计学意义。“group#food”的交互项检验结果为，Fgroup#food =0.19，P=0.831，表示交互项无统计学意义，即满足平行性检验。因此，本案例满足条件9。
(二) 统计描述
1. 软件操作
tabstat weight food,by(group) stat(count mean median sd min max) col(stat) long
2. 结果解读
图9中列出了1、2、3三组体重增量分别为36.633±20.369、46.427±21.389、119.027±31.944 g，三组的进食量分别为273.220±46.247、274.147±41.661、493.587±73.797 g。
(三) 单因素协方差分析
1. 软件操作
①根据上面的分析结果可知，“group#food”的交互项无统计学意义，则不应该放在模型中。结果如图10所示。
anova weight group c.food
②计算估计边际均值，并绘制边际均值图。结果如图11、图12所示。
margins group, plot
2. 结果解读
图10结果显示，模型(Model)的F =141.70，P <0.001，表示整体模型检验具有统计学意义。协变量food的P<0.001，说明进食量对小鼠体重的增加具有影响；而自变量group的P=0.149，提示三组小鼠体重的增重效果差异无统计学意义。
图11的结果显示了三组的估算边际均值、标准误及均值的95% 置信区间，1、2、3三组体重增量的估算边际均值分别为64.649 (95%CI: 55.552~73.746)、74.091 (95%CI: 65.037~83.144)、63.347 (95%CI: 50.054~76.640) g，可见三组体重增量的估算边际均值与三组体重增量的均值(图9)存在较大差异。图12的估算边际均值图也绘制了三组体重增量的情况，可见组2最高、组1其次、组3最低。但这种差异是否有统计学意义，还需要依据统计学推断结果做出判断。
(四) 事后检验(两两比较)
本案例中由于三组间整体比较差异无统计学意义，因此无须再进行事后两两比较。假设本案例三组间整体比较差异有统计学意义，其分析步骤如下。
1. 软件操作
pwcompare group, effects pwcompare group, effects bonferroni
2. 结果解读
图13-1和图13-2分别显示了未校正和用“Bonferroni (邦弗罗尼)”校正法的两两比较的结果。可知，三组中任意两组比较差异均无统计学意义(P>0.05)。
四、结论
本研究采用单因素协方差分析，判断在调整进食量后组别为1、2、3的三种饲料对小鼠的增重效果是否有差别。通过对模型残差绘制箱线图提示，数据不存在异常值；通过Shapiro-Wilk检验，提示各组因变量残差满足正态性要求；通过Levene’s检验，提示各组因变量残差满足方差齐性要求；通过绘制散点图和相关性分析，提示各组内协变量和因变量之间存在线性关系；通过平行性检验，发现“group#food”的交互项无统计学意义(F=0.19，P=0.831)，提示满足平行性检验要求。
1、2、3三组小鼠体重增量分别为36.633±20.369、46.427±21.389、119.027±31.944 g。单因素协方差分析显示在调整进食量后1、2、3三组体重增量的估算边际均值分别为64.649 (95%CI: 55.552~73.746)、74.091 (95%CI: 65.037~83.144)、63.347 (95%CI: 50.054~76.640) g，差异无统计学意义(Fgroup=2.00，Pgroup=0.149)。本研究结果提示1、2、3三种饲料对小白鼠体重的增重效果无差异。
五、知识小贴士
- 协方差分析是针对在实验设计阶段难以控制其取值水平，或者无法严格控制的因素，在统计分析阶段对其进行统计控制的一种分析方法，实质为线性回归分析和方差分析的结合。适用于完全随机设计、随机区组(配伍)设计、拉丁方设计、析因设计等类型的方差分析。
- 协方差分析一般要求协变量在组间的观察范围相差不宜太大(分析前最好先对协变量均数间的差别作假设检验)，否则修正后的边际均值的差值可能会落在回归线的延长线上，回归线外推后，是否仍然满足平行线和线性关系尚不可知，其协方差分析的结论可能不一定准确。
六、分析小技巧
- 正态性检验：对协方差分析的正态性条件考察应检验因变量残差的正态性，此时检验多组整体残差的正态性即可，无须检验每组残差的正态性。在Stata中使用“swilk”命令可进行“Shapiro-Wilk (夏皮罗-威尔克正态性)”检验，也可使用“qnorm”命令绘制正态分位图检验正态性。关于正态性检验的注意事项详见文章(医学统计学核心概念及重要假设检验的软件实现(2/4)——正态性假设检验的SPSS实现)。
- 方差齐性检验：对协方差分析的方差齐性条件考察应检验组间因变量残差是否相等。在Stata中使用“robvar”命令可进行“Levene’s检验”，其结果中“W0”即“Levene’s F statistic(F值)”，“Pr>F”即“P值”。关于方差齐性检验的更多内容请阅读(医学统计学核心概念及重要假设检验的软件实现(4/4)——方差齐性检验及SPSS实现)。