合理进行均值比较——单因素两水平和多水平设计定量资料多元方差分析.pdf

四川精神卫生 2023 年第 36 卷第 S1期http：/合理进行均值比较单因素两水平和多水平设计定量资料多元方差分析胡纯严1，胡良平1，2*（1.军事科学院研究生院，北京 100850；2.世界中医药学会联合会临床科研统计学专业委员会，北京 100029*通信作者：胡良平，E-mail：）【摘要】本文目的是介绍与单因素两水平和多水平设计定量资料多元方差分析有关的基本概念、计算方法、两个医学实例以及SAS实现。基本概念包括试验因素与属性因素、方差分析与协方差分析、自变量与协变量、矩阵与行列式；计算方法涉及一般统计量、4个多元方差分析检验统计量和多个均值向量的两两比较；两个医学实例分别为“2种脱机方式下患者脱机前有关定量指标的测定结果”和“3种不同类型地区的8岁男童生长发育调查数据”。借助SAS对两个实例中的定量资料分别进行多元方差分析和协方差分析，并对如何合理选取协变量进行阐释和讨论。【关键词】单因素多水平设计；协变量；多元方差分析；多元协方差分析；两两比较中图分类号：R195.1 文献标识码：A doi：10.11886/scjsws20230319002 Reasonably carry out mean value comparison：MANOVA of the quantitative data collected from single factor two-level and multi-level designHu Chunyan1，Hu Liangping1，2*（1.Graduate School，Academy of Military Sciences PLA China，Beijing 100850，China；2.Specialty Committee of Clinical Scientific Research Statistics of World Federation of Chinese Medicine Societies，Beijing 100029，China*Corresponding author:Hu Liangping，E-mail：）【Abstract】The purpose of the paper was to introduce the basic concepts，calculation methods，two medical examples and SAS implementation related to the multivariate analysis of variance（MANOVA）for the quantitative data with single factor two-level and multi-level design.Basic concepts included experimental factors and attribute factors，analysis of variance and covariance，independent variables and covariates，matrix and determinant.Calculation methods involved the general statistics，four MANOVA test statistics and pairwise comparison of multiple mean vectors.The two medical examples were measurement results of quantitative indicators related to patients before taking offline under two types of offline modes and investigation data on growth and development of 8-year-old boys in three different types of areas.With the help of SAS software，the multivariate analysis of variance and covariance were carried out on the quantitative data in the two cases，and how to reasonably select covariates was explained and discussed.【Keywords】Single factor multilevel design；Covariate；Multivariate analysis of variance；Multivariate analysis of covariance；Pairwise comparison单因素两水平和多水平设计是指试验中只涉及一个具有k（k2）水平的试验因素，全部受试对象被完全随机地分入该试验因素的k个水平组中去，或者从该试验因素的各个水平所对应的k个子总体中随机抽取受试对象的试验安排方法。本文将介绍与单因素两水平和多水平设计定量资料多元方差分析有关的基本概念、计算方法、两个医学实例及其SAS实现的方法。1 基本概念 1.1试验因素与属性因素在医学研究中，影响结果变量取值的因素非常多，这些因素被统称为“影响因素”。其中，有些因素是研究者希望着重考查的、可能还是研究者施加给受试对象的，这类因素被称为“试验因素”或“处理因素”；除试验因素外，在特定医学研究项目中出现的其他因素被统称为“非试验因素”。在非试验因素中，与受试对象自身条件或心理素质有关的因素，通常被称为“属性因素”，例如性别、年龄、职业、学历、家庭经济状况等。1.2方差分析与协方差分析方差分析最早由 Fisher1于 1918年提出，此法是基于对总离均差平方和的分解以及对总自由度302023 年第 36 卷第 S1期四川精神卫生http：/的分解，构造反映不同变异来源的均方（也称为方差），再基于不同变异来源的均方与误差的均方构造出检验统计量F，进而基于F分布推断某给定因素不同水平条件下的均值之间差异是否有统计学意义。Levene于1960年首先提出将方差分析用于方差齐性检验2，Brown和Forsythe于1974年对Levene提出的方法进行了扩展3。协方差分析是将线性回归分析与方差分析相结合的统计分析方法，先采用线性回归方法消除对结果变量有影响的某个或某些定量的自变量的影响，之后再进行方差分析，其主要作用是平衡定量的干扰因素对观测结果的影响。在协方差分析中，协变量可以只有一个，也可以有多个。根据协变量的数量和结果变量的数量，可以将单因素多水平设计定量资料的多元协方差分析称为带有k（k1）个协变量的单因素多水平设计m（m2）元定量资料的协方差分析。1.3自变量与协变量在进行差异性分析，常把影响结果取值的变量称为“影响因素”；在进行回归分析或协方差分析时，常把影响结果取值的变量称为“自变量”。在一个实际问题中，可能存在许多自变量，人们又习惯把定性自变量称为“分类变量”，而把定量自变量称为“协变量”，旨在消除或平衡定量自变量在分类自变量所形成的各组之间的影响，以便研究者着重考查各分类自变量对结果变量的影响。1.4矩阵与行列式1.4.1矩阵的定义由mn个数aij（i=1，2，m；j=1，2，n）排成的矩形阵列，称为 mn 矩阵，简记为A=(aij)m n。见式（1）。A=a11a12a1na21a22a2nam1am2amn（1）在式（1）的括号内，横排称为矩阵的行，竖排称为矩阵的列；aij称为矩阵A的第i行第j列元素。当m=n时，称此时的矩阵A为方阵。1.4.2行列式的定义n阶行列式是由 n2个数 aij（i，j=1，2，n）所确定的一个Dn，表达式见式（2）。Dn=|a11a12a1na21a22a2nan1an2ann=j1jn(-1)N()j1jna1j1a2j2anjn（2）式（2）中，“j1jn”表示对所有n个元素排列求和。Dn是n！项的代数和，每项是n个因子的乘积。n个因子取自不同的行与不同的列，每项前的符号由排列j1jn的逆序数N（j1jn）决定。【说明】逆序数的含义如下：在一个排列中，如果一对数的前后位置与大小顺序相反，则称它们之间有一个逆序。一个排列中，逆序的总数称为这个排列的逆序数。2 计算方法 2.1一般统计量单因素两水平和多水平设计定量资料多元方差分析的计算方法与一元方差分析类似，也是将定量结果的总变异分解为组间变异与组内变异，只是多元方差分析中参数和统计量是以向量或矩阵的形式呈现的。实质上多元方差分析是将一元方差分析中比较若干个一维均值的原理和方法拓展到比较若干个m维均值向量的过程。在单因素设计中，设试验因素有g个水平，也就是有g个试验组；第k个试验组的样本含量为nk，试验中的总例数为n=k=1gnk；需要观测的定量结果变量的个数为m。用k表示第k组的总体均值向量，Xk表示第k组的样本均值向量，X表示所有观测的样本均值向量，Xki表示第 k组第 i个观测向量。此时，第k组的离均差矩阵Lk和协方差矩阵Sk可表示如下，分别见式（3）和式（4）。Lk=i=1nk(Xki-Xk)(Xki-Xk)（3）Sk=Lknk-1（4）与一元方差分析相似，多元方差分析中全部观测的总变异由总离均差矩阵T表示，组间变异由组间离均差矩阵H表示，组内变异由组内离均差矩阵E 表示。当总变异被分解为组间变异与组内变异时，相应的总离均差矩阵被分解为组间离均差矩阵与组内离均差矩阵两部分，见式（5）。T=H+E（5）总离均差矩阵、组间离均差矩阵和组内离均差31四川精神卫生 2023 年第 36 卷第 S1期http：/矩阵可以被表示如下，分别见式（6）、式（7）、式（8）。T=k=1gi=1nk(Xki-X)(Xki-X)（6）H=i=1nknk(Xk-X)(Xk-X)（7）E=k=1gi=1nk(Xki-Xk)(Xki-Xk)=k=1gLk=k=1g(nk-1)Sk（8）2.2Wilks检验统计量根据组间离均差矩阵H与组内离均差矩阵E，可计算 Wilks 检验统计量。该统计量是 Wilks4于1932年提出的一种广义方差比，其计算公式见式（9）。Wilks=|E|H+E（9）式（9）中分子和分母都是行列式。当 Wilks检验统计量取值较小时，说明组间变异H大于随机误差的效应E，此时有理由怀疑原假设H0：1=2=g的正确性。Rao5于1951年给出由Wilks求F值的近似算法，H和E通过Wilks值可转变为F值，从而以简便方法实现多元方差分析。常见情况下由Wilks求F值的关系式见表1。2.3其他3个检验统计量在多元方差分析中，除Wilks统计量之外，还可计算其他3个检验统计量，分别为Pillai迹、Hotelling-Lawley迹以及Roy极大根6。与Wilks一样，这些统计量都是计算F统计量的过渡检验统计量。在大多数情况下，由不同检验统计量转换为F统计量后的结果相差无几。当响应变量组的第一个（或最大）特征根完全解释了响应变量组的共同变异时，4个检验统计量对应的 F 值相等，并且精确服从 F 分布。当响应变量组的第一个特征根不能完全解释响应变量组的共同变异时，4个检验统计量对应的F值不相等，并且它们不能精确地服从F分布，而是逼近F分布，此时，Roy的F值是上限值。为稳定起见，通常选用Wilks检验法或Hotelling-Lawley检验法。2.4多个均值向量的两两比较对单因素多水平设计定量资料进行多元方差分析之后，如果各组均值向量之间的差异有统计学意义，还需要进行各组均值向量之间的两两比较。多个均值向量两两比较所采用的方法，与两个均值向量比较的HotellingsT2检验在形式上完全相同，唯一不同的是，在此所用的是多元方差分析中的合并类内协方差矩阵，而不是由被比较的两组资料来计算误差的合并协方差矩阵。两两比较的检验统计量仍然是F统计量，其计算公式见式（10）。F=n-g-p+1(n-g)pT2ij，1=p，2=n-g-p+1（10）式（10）中T2ij的计算公式见（11）。T2ij=ninjni+nj(Xi-Xj)S-1(Xi-Xj)，i，j=1，2，g，i j，S=1n-gE（11）式（10）中，p为矩阵H+E的秩。式（11）中，S为多元方差分析中的合并类内协方差矩阵。3 实例与SAS实现 3.1问题与数据结构3.1.1两个医学问题及数据【例1】某研究欲比较智能监护（SmartCare）脱机模式和间断停用呼吸机两种方法对老年慢性阻塞性肺疾病（COPD）患者脱机成功率的影响。受试对象为某医院老年呼吸科收治的行机械通气的男性呼吸衰竭患者，共28例，机械通气时间均超过3天，年龄 7091 岁，随机分为智能监护脱机组（SC 组，表1常见情况下基于Wilks值求F值的公式Table 1Under common circumstances，the formula for calculating the value of F based on the value of Wilks 响应变量数m=1m=2m1m1组数g2g2g=2g=2转换关系式F=(n-gg-1)(1-)F=(n-g-1g-1)(1-)F=(n-m-1m)(1-)F=(n-m-2m)(1-)F统计量的自由度1=g-1，2=n-g1=2(g-1)，2=2(n-g-1)1=m，2=n-m-11=2m，2=2(n-m-2)322023 年第 36 卷第 S1期四川精神卫生http：/n=13）和间断停用呼吸机组（SBT 组，n=15）进行脱机治疗。部分指标记录结果见表27。假定资料满足参数检验的前提条件，试比较SBT组和SC组在脱机前机械通气时间、急性生理和慢性健康状态评分（APACHEII）、痰液pH值和二氧化碳分压（PaCO2）四项指标的均值向量差异是否有统计学意义。【例2】分别对大城市、县城和农村这三类地区的8岁男童进行随机抽样，调查8岁男童的身高、体重和胸围，样本含量均为30，三类地区8岁男童的生长发育调查数据见表37。假定资料满足参数检验的前提条件，试比较不同类型地区的男童身高、体重和胸围均值向量之间的差异是否有统计学意义。3.1.2对数据结构的分析例1中有一个二水平的试验因素，即脱机方式包括间断停用呼吸机（SBT）与智能监护脱机（SC），这两组分别包含 15 例和 13 例患者。研究者观测并收集了两组患者的年龄、脱机前机械通气时间、APACHEII、痰液pH值和PaCO2这5个定量变量的取值，其中，年龄为定量的协变量，已假定其他4个定量变量属于定量结果变量。故表2资料属于具有一个协变量的单因素两水平设计四元定量资料。值得注意的是，从临床基础知识角度考量，本例中的“脱机前机械通气时间、APACHEII、痰液pH值和PaCO2”这4个定量变量都应该被视为“脱机成功率”的影响因素，而不属于定量结果变量。这里，假设它们为定量结果变量，其目的是借用此资料来介绍单因素两水平设计四元方差分析方法。例2中有一个三水平的试验因素，即地区类型（大城市、县城和农村），这三组均分别包含30例8岁男童。研究者对三组都观测并收集了身高、体重和胸围这3个定量变量的数值。故表3资料属于单因素三水平设计三元定量资料。3.1.3创建SAS数据集分析例1资料所需要的SAS数据步程序如下：data a1；do group=SBT，SC；do i=1 to 15；input id age time apac ph paco；output；end；end；cards；1 78 33 10.2 7.5 36.42 83 10 10.7 7.5 45.63 90 19 11.0 7.5 39.8 27 85 16 8.5 7.4 45.928 77 4 4.2 7.4 27.5表2两组患者脱机前有关定量指标的测定结果Table 2Measurement results of relevant quantitative indicators of two groups of patients before offline组别SBT组SC组患者编号1215161728年龄（岁）788384858477脱机前机械通气时间（天）331024124APACHEII（分）10.210.710.225.723.84.2痰液pH值7.57.57.37.47.37.4PaCO2（mmHg）36.445.631.643.675.827.5注：APACHEII，急性生理和慢性健康状态评分；PaCO2，二氧化碳分压；患者编号115属于SBT组，1628属于SC组表3三类地区的8岁男童的生长发育调查数据Table 3Investigation data of growth and development of 8-year-old boys in three different types of areas编号122930身高（cm）大城市119.8121.7136.8141.2体重（kg）22.621.526.731.7胸围（cm）60.555.560.563.5身高（cm）县城125.1127.0109.9125.4体重（kg）23.021.517.823.5胸围（cm）62.059.056.460.8身高（cm）农村118.3121.3131.2122.6体重（kg）20.420.025.524.1胸围（cm）54.454.358.358.933四川精神卫生 2023 年第 36 卷第 S1期http：/ ；run；【变量说明】变量id、age、time、apac、ph、paco分别代表编号、年龄、脱机前机械通气时间、APACHEII、痰液pH值、PaCO2；第327行之间的数据已省略，第28行之后加了两行缺失数据，这是因为SC组的样本含量比SBT组少了2例。分析例2资料所需要的SAS数据步程序如下：data a2；length group$7；do i=1 to 30；do group=city，county，country；input height weight chest；output；end；end；cards；119.8 22.6 60.5 125.1 23.0 62.0 118.320.4 54.4121.7 21.5 55.5 127.0 21.5 59.0 121.320.0 54.3 136.8 26.7 60.5 109.9 17.8 56.4 131.2 25.5 58.3141.2 31.7 63.5 125.4 23.5 60.8 122.6 24.1 58.9；run；【变量说明】group 代表地区类型，city、county、country分别代表大城市、县城、农村，height、weight、chest分别代表身高、体重、胸围。3.2用SAS实现统计分析3.2.1分析例1中的资料不考虑协变量age对结果变量的影响，设所需要的SAS过程步程序如下：proc glm data=a1；class group；model time apac ph paco=group/nouni ss3；manova H=group；means group；run；quit；【SAS 程序说明】若在 model语句中，删除选项“nouni”，可先输出一元方差分析计算结果8，后输出多元方差分析计算结果；manova语句是实现多元方差分析的重要语句，通过关键词“H=”指定影响因素。【SAS输出结果及解释】T2=0.690，F=3.970，分子和分母的自由度分别为4和23，P=0.014，说明由4个定量指标组成的均值向量在 SBT组与 SC组之间 差 异 有 统 计 学 意 义。SBT 组 的 均 值 向 量 为15.800，10.087，7.447，39.533；SC组的均值向量为 15.769，20.038，7.377，42.746。考虑协变量age对结果变量的影响，设所需要的SAS过程步程序如下：proc glm data=a1；class group；model time apac ph paco=age group/nouni ss3；manova H=age group；means group；run；quit；【SAS输出结果及解释】对因素group而言，T2=0.574，F=3.160，分子和分母的自由度分别为 4 和22，P=0.034，说明由4个定量指标组成的均值向量在SBT组与SC组之间差异有统计学意义。对协变量age而言，T2=0.294，F=1.620，分子和分母的自由度分别为4和22，P=0.206，说明由4个定量指标组成的均值向量在各年龄组之间差异无统计学意义。【结论】由脱机前机械通气时间、APACHEII、痰液pH值、PaCO2这4个定量指标组成的均值向量在 SBT组与 SC 组之间差异有统计学意义；根据两组均值向量的各分量的数值 SBT 组为 15.800，10.087，7.447，39.533，SC 组为15.769，20.038，7.377，42.746 可初步判定，它们之间的差异主要由第二个分量（SBT组和SC组APACHEII的均值分别为10.087和20.038）的差异所致。严格地说，需要进行单变量差异性检验，才可得出明确结论。因篇幅所限，此处从略。另外，本例中协变量age对结果的影响非常小，当考虑或不考虑协变量时，4个定量指标组成的均值向量在SBT组与SC组之间差异均有统计学意义。3.2.2分析例2中的资料将身高、体重和胸围视为3个定量结果变量，设所需要的SAS过程步程序如下：proc sort data=a2；by group；342023 年第 36 卷第 S1期四川精神卫生http：/run；proc glm data=a2；class group；model height weight chest=group/nouni ss3；contrast city vs country group 1-1 0；contrast city vs county group 1 0-1；contrast country vs county group 0 1-1；manova H=group；means group；run；quit；【SAS输出结果及解释】Wilks=0.550，F=9.880，分子和分母的自由度分别为6和170，P0.001，说明由3个定量指标组成的均值向量在三组之间差异有统计学意义。均值向量在三组中任何两组之间的比较结果如下：city 组与 country 组比较，Wilks=0.692，F=12.590，分子和分母的自由度分别为3和85，P0.001，说明由3个定量指标组成的均值向量在city组与country组之间差异有统计学意义；city组与 county 组比较，Wilks=0.798，F=7.170，分子和分母的自由度分别为 3 和 85，P0.001，说明由3个定量指标组成的均值向量在city组与county组之间差异有统计学意义；country组与county组比较，Wilks=0.734，F=10.280，分子和分母的自由度分别为3和85，P0.001，说明由3个定量指标组成的均值向量在country组与county组之间差异有统计学意义。city 组、country 组、county 组的均值向量身 高 均 值，体 重 均 值，胸 围 均 值 分 别 为：126.520，23.517，60.020、120.917，22.527，57.423、120.517，20.767，57.883。将身高视为协变量，将体重和胸围视为2个定量结果变量，设所需要的SAS过程步程序如下：proc glm data=a2；class group；model weight chest=height group/nouni ss3；contrast city vs country group 1-1 0；contrast city vs county group 1 0-1；contrast country vs county group 0 1-1；manova H=height group；means group；run；quit；【SAS 输出结果及解释】对因素 group 而言，Wilks=0.704，F=8.150，分子和分母的自由度分别为4和170，P0.001，说明由2个定量指标组成的均值向量在三组之间差异有统计学意义；对协变量height 而言，Wilks=0.426，F=57.150，分子和分母的自由度分别为2和85，P0.001，说明由2个定量指标组成的均值向量在不同身高组之间差异有统计学意义。均值向量在三组中任何两组之间的比较结果如下：city组与country组比较，Wilks=0.828，F=8.830，分子和分母的自由度分别为 2 和 85，P0.001，说明由 3个定量指标组成的均值向量在city组与country组之间差异有统计学意义；city组与 county 组比较，Wilks=0.978，F=0.980，分子和分母的自由度分别为2和85，P=0.381，说明由3个定量指标组成的均值向量在city组与county组之间差异无统计学意义；country 组与 county 组比较，Wilks=0.734，F=15.370，分子和分母的自由度分别为2和85，P0.001，说明由3个定量指标组成的均值向量在 country组与 county组之间差异有统计学意义。city组、country组、county组的均值向量 体重均值，胸围均值 分别为：23.517，60.020、22.527，57.423、20.767，57.883。【结论】在以身高、体重、胸围3个定量指标作为结果变量时，三元方差分析结果表明，3个均值向量在三组之间总的差异有统计学意义；而且，三组中任何两组之间的差异均有统计学意义。在以身高为协变量，以体重和胸围这两个定量指标作为结果变量时，二元方差分析结果表明，3个均值向量在三组之间总的差异有统计学意义；3个均值向量在各年龄组之间总的差异也有统计学意义。但是，3个均值向量在三组中任何两组之间的差异并非均有统计学意义，其中，city（大城市）与county（县城）之间差异无统计学意义。4 讨论与小结 4.1讨论多元方差分析是从整体上来评价均值向量在组间的差异是否有统计学意义，即使得到差异有统计学意义的结果，并不意味着每个定量指标在组间的差异也都有统计学意义。另外，哪些定量变量适合作为协变量，取决于专业知识和研究目的。例如，在本文例1中，年龄与其他4个定量指标（脱机前机械通气时间、APACHEII、痰液 pH 值和 PaCO2）在专业上是有区别的，因为年龄很可能是影响4个定量指标的“试验因素”之一，而不是反映治疗效果的。又例如，在本文例2中，是将身高视为一个协变量合适，还是应该将身高与体重和胸围一起视为定35四川精神卫生 2023 年第 36 卷第 S1期http：/量结果变量合适，取决于研究目的。若研究目的是希望把“身高、体重和胸围”当作反映“生长发育”的3个重要指标，则将它们视为三元定量结果变量合适；若研究目的是希望把“体重和胸围”当作反映“生长发育”的2个重要指标，则将它们视为二元定量结果变量、将“身高”视为协变量合适。4.2小结本文介绍了与单因素两水平和多水平设计定量资料多元方差分析有关的基本概念、计算方法、两个医学实例及其SAS实现。基本概念包括试验因素与属性因素、方差分析与协方差分析、自变量与协变量、矩阵与行列式；计算方法包括一般统计量、4个多元方差分析检验统计量和多个均值向量的两两比较；两个医学实例分别涉及单因素两水平与单因素三水平设计多元定量资料。基于SAS实现了单因素两水平与单因素三水平设计定量资料的多元方差分析和协方差分析。参考文献1Fisher RA.The correlation between relatives on the supposition of Mendelian inheritanceJ.Earth and Environmental Science Transactions of The Royal Society of Edinburgh，1918，52（2）：399-433.2胡良平.医学统计学：运用三型理论分析定量与定性资料M.北京：人民军医出版社，2009：37-41.Hu LP.Medical statistics：applying the three types theory to analyze quantitative and qualitative dataM.Beijing：Peoples Military Medical Publishing House，2009：37-41.3Brown MB，Forsythe AB.Robust tests for equality of variancesJ.J Am Stat Assoc，1974，69（346）：364-367.4Wilks SS.Certain generalization in the analysis of varianceJ.Biometrika，1932，24（3-4）：471-494.5Rao CR.An asymptotic expansion of the distribution of Wilks criterionJ.Bull Inst Internat Statist，1951，33（Part）：177-180.6SAS Institute Inc.SAS/STAT15.1 users guide M.Cary，NC：SAS Institute Inc，2018：3957-4142.7胡良平.面向问题的统计学：（3）试验设计与多元统计分析M.北京：人民卫生出版社，2012：353-389.Hu LP.Problem-oriented statistics：（3）experimental design and multivariate statistical analysisM.Beijing：Peoples Medical Publishing House，2012：353-389.8胡纯严，胡良平.如何正确运用方差分析：单因素多水平设计定量资料一元方差分析 J.四川精神卫生，2022，35（1）：16-20.Hu CY，Hu LP.How to use analysis of variance correctly：an analysis of variance for the univariate quantitative data collected from the design of a single factor with multi-levelJ.Sichuan Mental Health，2022，35（1）：16-20.（收稿日期：2023-03-19）（本文编辑：陈霞）36