合理进行多元分析——主成分分析.pdf

四川精神卫生 2023 年第 36 卷第 S1期http：/合理进行多元分析主成分分析胡纯严1，胡良平1，2*（1.军事科学院研究生院，北京 100850；2.世界中医药学会联合会临床科研统计学专业委员会，北京 100029*通信作者：胡良平，E-mail：）【摘要】本文目的是介绍与主成分分析有关的基本概念、计算方法、两个实例以及SAS实现。基本概念包括相关矩阵、特征值与特征向量、主成分变量、主成分表达式和主成分的性质；计算方法涉及特征值与特征向量的求法、主成分分析的计算原理以及系数估计和个数的确定；两个实例中的资料分别为“20例肝病患者的4项肝功能指标的测定结果”和“23种肿瘤类期刊的文献计量学指标的调查结果”；借助SAS对两个实例中的定量资料进行了主成分分析，并基于主成分的计算结果分别实现了样品聚类和样品排序，并对输出结果作出了解释。【关键词】特征值；特征向量；主成分分析；样品聚类；样品排序中图分类号：R195.1 文献标识码：A doi：10.11886/scjsws20230605001Reasonably carry out multivariate analysis：principal component analysisHu 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 this article was to introduce the basic concepts，calculation methods，two examples and SAS implementation related to the principal component analysis.Basic concepts included correlation matrix，eigenvalues and eigenvectors，principal component variables，principal component expressions and principal component properties.The calculation method involved the calculation of eigenvalues and eigenvectors，the calculation principle and the coefficient estimation and the number determination of the principal component.The data in the two examples were measurement results of 4 liver function indicators in 20 patients with liver disease and survey results of literature metrology indicators in 23 tumor journals.With the help of SAS software，the principal component analysis was carried out on the quantitative data in the two cases，and based on the calculation results of the principal components，the sample clustering and sample sorting were respectively realized，and a reasonable explanation was given for the output results.【Keywords】Eigenvalue；Eigenvector；Principal component analysis；Sample clustering；Sample sorting在生物医学和临床研究中，研究者经常收集到单组设计多元定量资料。如何选择合适的统计分析方法处理这种定量资料，是研究者经常面临的一个棘手的统计问题，因为可用于处理这种多元定量资料的多元统计分析方法约有十几种。本文将介绍一种最简单、最基础的多元统计分析方法，即主成分分析。1 基本概念 1.1相关矩阵设具有同质性的n个个体，测量其m个定量指标（记为x1，x2，xm）的取值，采用Pearson相关分析公式 见式（1）计算出任意两个定量变量之间的相关系数（记为rij，i，j=1，2，m）1，将它们按一定顺序排列成一个 mm 方阵，见式（2），此方阵被称为相关矩阵。rij=k=1n()xik-x ik()xjk-x jkk=1n()xik-x ik2k=1n()xjk-x jk2（1）R=r11r12r1mr21r22r2mrm1rm2rmm（2）式（1）中，k代表观测或个体的编号。式（2）中，rij=rji，i，j=1，2，m，即R是一个对称的矩阵。1.2特征值与特征向量对于一个 m阶矩阵 见式（3），如果存在一个数0和非零向量X0使式（4）成立，则称0为矩阵A的482023 年第 36 卷第 S1期四川精神卫生http：/特征值，称 X0为矩阵 A 对应于特征值 0的特征向量2。A=a11a12a1ma21a22a2mam1am2amm（3）AX0=0X0（4）1.3主成分变量主成分变量（简称主成分）是一个不能被直接观测的隐变量。在一个单组设计m元定量资料中，有m个主成分，其中，每一个主成分都是由m个原变量线性组合而成的，但它们彼此互相独立；通常，它们所包含的信息量是不相等的，且满足如下关系，见式（5）。Pm Pm-1 P2 P1（5）式（5）中，Pi(i=1，2，m)代表第i个主成分。主成分由原始变量的线性组合而成，故主成分也被称为线性主成分。当对数据进行非线性变换后再进行主成分分析时，就称为非线性主成分分析3。1.4主成分表达式主成分的常见表达式有以下 4 种，见式（6）、式（7）、式（8）、式（9）。Z1=a11()X1-X1+a12()X2-X2+a1m()Xm-XmZ2=a21()X1-X1+a22()X2-X2+a2m()Xm-XmZm=am1()X1-X1+am2()X2-X2+amm()Xm-Xm（6）Z1=b11x1+b12x2+b1mxmZ2=b21x1+b22x2+b2mxmZm=bm1x1+bm2x2+bmmxm（7）Z1=c11x1+c12x2+c1mxmZ2=c21x1+c22x2+c2mxmZm=cm1x1+cm2x2+cmmxm（8）x1=c11Z1+c21Z2+cm1Zmx2=c12Z1+c22Z2+cm2Zmxm=c1mZ1+c2mZ2+cmmZm（9）式（6）、式（7）、式（8）、式（9）中，Xi为第i个原始变量，xi为Xi的标准化变量，见式（10）。xi=Xi-XiSi（10）式（10）中，Si为第 i 个变量的样本标准差，i=1，2，m。式（6）、式（7）、式（8）、式（9）中，各aij、bij、cij都是线性组合的系数，称为因子负荷量。其中，aij、bij分别是基于原始数据与标准化数据计算所得到的第i个特征向量中的第j个元素。cij的计算见式（11）。cij=ibij（11）式（11）中，i为第i个特征值；bij为式（7）中的系数；i，j=1，2，m。1.5主成分的性质各主成分之间互不相关，即任何两个主成分之间的相关系数为0，若原变量服从正态分布，则各主成分之间互相独立；全部主成分所反映的信息等于全部原变量的总信息；各主成分的作用大小不等：第一主成分的作用大于等于第二主成分，第二主成分的作用大于等于第三主成分，以此类推，最后一个主成分的作用最小。因篇幅所限，主成分的其他性质从略。2 计算方法 2.1特征值与特征向量的求法第一步，构造一个特征矩阵，见式（12）：E-A=-a11-a12-a1m-a21-a22-a2m-am1-am2-amm（12）式（12）中，是一个未知量；E是一个对角线上元素全为 1 的 mm 对角矩阵；矩阵 A 的定义见前文式（3）。第 二 步，构 造 一 个 的 n 齐 次 多 项 式，见式（13）。f()=|E-A|=|-a11-a12-a1m-a21-a22-a2m-am1-am2-amm（13）式（13）中的“|*|”代表由“*”确定的行列式，它是一个具体的数值。第三步，构造一个矩阵 A 的特征方程组，见式（14）。f()=|E-A|=0（14）第四步，求出特征方程组的全部根，即特征值k（k=1，2，m）。第五步，把特征值逐个代入齐次线性特征方程组，见式（15）。(kE-A)X=0（15）49四川精神卫生 2023 年第 36 卷第 S1期http：/式（15）中，X是由全部定量变量构成的向量，见式（16）。X=x1x2xm（16）求出方程组的基础解系，则这个基础解系的非零线性组合就是矩阵A对应于特征值k（k=1，2，m）的全部特征向量。因篇幅所限，具体求解过程见文献 2，4。2.2主成分分析的计算原理设有n例儿童的身高（X1）和体重（X2）两个观测指标，显然它们之间有较强的相关性。以 X1为横轴，以X2为纵轴，绘制散点图如图1所示。可见这n个点的分布在一条直线的近旁，呈现出直线化趋势。它们沿X1轴和X2轴方向都具有较大的变异度，个体在某个方向上的变异度可用该方向上相应观测变量的方差来表示。以此直线作为新的横轴Z1，再作一条垂直于Z1的直线作为纵轴Z2，在平面上，这n点的变异主要集中在Z1方向上，在Z2方向上变异很小。所以，研究这n个对象的变异，可以只考虑Z1值的大小，忽略Z2值的大小。也就是说，若取Z1作为第一主成分，则Z1就反映了原始指标X1和X2所包含的大部分信息。若将X1和X2标准化后的指标记为Y1、Y2，则Z1、Z2与Y1、Y2有以下关系，见式（17）和式（18）。Z1=l11Y1+l12Y2（17）Z2=l21Y1+l22Y2（18）Z1、Z2是 Y1、Y2的线性函数，显然也是 X1、X2的线性函数，且Z1、Z2不相关。称Z1为第一主成分，Z2为第二主成分，并称这种分析方法为主成分分析法。2.3主成分表达式中系数的估计为了计算简便，主成分的计算一般从求相关矩阵出发2，5。以表 1为例，由 m个变量 X1，X2，Xm的n个样本观测值求出主成分。计算步骤如下。第一步，对各原始指标数据进行标准化，见式（19）。然后用标准化的数据来计算主成分。记X为标准化后的数据矩阵，见式（20）。xij=Xij-XSj，j=1，2，m（19）X=x11x21xn1x12x22xn2x1mx2mxnm（20）第二步，求出X的相关矩阵R 标准化后，X的相关矩阵即为协方差矩阵Cov（X），见式（21）。R=Cov(X)=r11r21rm1r12r22rm2r1mr2mrmm=1r21rm1r121rm2r1mr2m1（21）第三步，求出相关矩阵的特征值和特征值所对应的特征向量。求主成分问题，实际上就是要求出X的协方差矩阵Cov（X）（这里即为X的相关矩阵R）的特征值和特征向量。由于R为半正定矩阵，故可由R的特征方程|R-I|=0求得m个非负特征值，这些特征值按从大到小排序为：12m。特征方程见式（22）。()R-iI ai=0i=1，2，maiai=1（22）求出式（22）的解，得每一个特征值对应的单位特征向量ai=（ai1，ai2，aim），可写出主成分的表达式，见式（23）。Zi=aiX=ai1x1+ai2x2+aimxm，i=1，2，m（23）2.4主成分个数的确定主成分分析的目的是采用降维的思想，降低数据的维数2，5。故一般只取前k（km）个主成分。通常采用累计贡献率（c）的方法确定k值，即根据实际问题给出一个c值，一般c在70%85%之间。若前k1个主成分的累计贡献率小于c，而前k个主成分的累计贡献率c，则取k个主成分。图1主成分分析示意图Figure 1Schematic diagram of principal component analysis502023 年第 36 卷第 S1期四川精神卫生http：/3 实例与SAS实现 3.1问题与数据结构3.1.1两个实际问题及数据【例1】某医学院的研究者测得20例肝病患者的4 项肝功能指标的具体值，各指标分别为：转氨酶（X1）、肝大指数（X2）、硫酸锌浊度（X3）和胎甲球（X4）。试采用某种多元统计分析方法处理资料，并依据一定的理由对20例肝病患者进行分类。数据见表16。【例2】某研究者提供了某年我国23种肿瘤类期刊文献计量学方面8个定量指标的调查数据，这些指标分别为：载文量（X1）、基金论文比（X2）、总被引频次（X3）、影响因子（X4）、5年影响因子（X5）、即年指标（X6）、被引半衰期（X7）和 Web 即年下载率（X8）。试采用某种多元统计分析方法处理资料，并依据一定的理由对23种肿瘤类期刊进行排序。各指标的具体数据见表27。3.1.2对数据结构的分析在例1中，20例患者均为肝病患者，研究者获取了这些患者肝功能状态的4项定量指标的值，故这是一个单组设计4元定量资料。在例2中，23种期刊均为肿瘤类期刊，研究者收集了关于这些期刊文献计量学方面的8项定量指标的数值，故这是一个单组设计8元定量资料。3.1.3创建SAS数据集分析例 1 资料，设所需要的 SAS 数据步程序如下：data a1；input X1-X4；obs=_n_；cards；表120例肝病患者的4项肝功能指标的测定结果Table 1Measurement results of four liver function indicators in 20 patients with liver disease患者编号121920转氨酶（X1）401020120肝大指数（X2）2.01.51.02.0硫酸锌浊度（X3）551220胎甲球（X4）2030600表223种肿瘤类期刊的文献计量学指标的调查结果Table 2Survey results of literature metrology indicators of 23 tumor journals刊名中华肿瘤杂志癌症中国肿瘤临床中华放射肿瘤学杂志肿瘤中国肿瘤中国肿瘤临床与康复肿瘤防治研究实用肿瘤杂志实用癌症杂志中国癌症杂志肿瘤防治杂志中国肺癌杂志癌变 畸变 突变实用肿瘤学杂志临床肿瘤学杂志中国肿瘤生物治疗杂志肿瘤研究与临床现代肿瘤医学白血病 淋巴瘤河南肿瘤学杂志肿瘤学杂志四川肿瘤防治载文量（X1）23431650710219124325530219825118850917212023332582256336200274188143基金论文比（X2）0.350.490.330.170.400.130.050.250.170.150.140.180.240.450.090.060.700.070.090.110.040.130.04总被引频次（X3）2 7051 9351 710942702660595585566546526476412341302298296246243231230207110影响因子（X4）1.4150.7420.4201.0110.4700.3580.2060.2800.3260.2960.3550.2300.6030.4060.1370.3180.3870.1630.2590.1590.1000.2330.1105年影响因子（X5）1.3940.8790.6731.2900.5250.3670.2280.3320.3320.2940.4190.2440.6430.4520.2200.2620.4590.1630.1970.2020.0970.1860.132即年指标（X6）0.1200.1040.0260.0290.0210.0580.0200.0230.0350.0120.0320.0240.0580.1670.0090.0370.0370.0080.0630.0250.0000.0050.000被引半衰期（X7）6.04.35.35.34.73.74.14.75.33.93.83.02.95.35.52.84.43.92.54.44.63.54.4Web即年下载率（X8）38.90035.20022.3007.50025.00015.5007.00024.70024.10019.30025.40015.80021.70031.2008.60015.10011.20018.00012.40015.80010.20015.40012.30051四川精神卫生 2023 年第 36 卷第 S1期http：/40 2.0 5 2010 1.5 5 30120 3.0 13 50250 4.5 18 0120 3.5 9 5010 1.5 12 5040 1.0 19 40270 4.0 13 60280 3.5 11 60170 3.0 9 60180 3.5 14 40130 2.0 30 50220 1.5 17 20160 1.5 35 60220 2.5 14 30140 2.0 20 20220 2.0 14 1040 1.0 10 020 1.0 12 60120 2.0 20 0；run；分析例 2 资料，设所需要的 SAS 数据步程序如下：data a1；input X1-X8；Jour=_n_；cards；234 0.35 2705 1.415 1.394 0.120 6.0 38.900316 0.49 1935 0.742 0.879 0.104 4.3 35.200 188 0.13 207 0.233 0.186 0.005 3.5 15.400143 0.04 110 0.110 0.132 0.000 4.4 12.300；run；3.2用SAS实现统计分析3.2.1分析例1中的资料设所需要的SAS程序如下8：proc princomp out=aaa prefix=z；var X1-X4；run；data a2；set aaa；maxz=max（of z1-z4）；if maxz=z1 then do；zz=z1；c1=obs；end；if maxz=z2 then do；zz=z2；c2=obs；end；if maxz=z3 then do；zz=z3；c3=obs；end；if maxz=z4 then do；zz=z4；c4=obs；end；maxz=round（maxz，0.001）；proc print data=a2；var maxz c1-c4；run；【SAS输出结果及解释】由4个变量两两之间的相关系数按一定规律排列出来，所形成的方阵称为相关矩阵。见表3。X1与X2之间的相关性最高，r值为0.695；X1与X3之间的相关性次之，r值为0.220。相关矩阵的特征值见表4。若希望累计贡献率在85%以上，至少需要前3个主成分。与 4 个特征值对应的特征向量的计算结果见表5。由表5可知，4列上的数值代表4个特征向量的元素，按前文的式（7）可写出4个主成分表达式，下面仅呈现第1个主成分的表达式，见式（24）。z1=0.700 x1+0.699x2+0.083x3+0.163x4（24）值得一提的是：SAS输出结果中的英文字母都是大写的，而且是正体；但表5中的系数是基于标准表3相关矩阵Table 3Correlation matrix矩阵X1X2X3X4X11.0000.6950.2200.025X21.000-0.1480.135X31.0000.071X41.000注：X1，转氨酶；X2，肝大指数；X3，硫酸锌浊度；X4，胎甲球表4相关矩阵的特征值Table 4Eigenvalues of the correlation matrix特征值编号1234特征值1.7181.0940.9810.207相邻特征值之差0.6250.1120.774-比例0.4300.2730.2450.052累计贡献率43.0%70.3%94.8%100.0%表5与特征值对应的特征向量Table 5Eigenvectors corresponding to the eigenvalues变量名称X1X2X3X4特征变量z10.7000.6990.0830.163特征变量z20.095-0.2840.9040.305特征变量z3-0.2400.058-0.2700.931特征变量z4-0.6660.6640.319-0.121522023 年第 36 卷第 S1期四川精神卫生http：/化变量计算的结果，故呈现在表达式中，应将英文字母写成小写形式。基于临床专业知识，结合各主成分表达式中系数的绝对值和正负号，解读各主成分所代表的含义：第一个主成分表达式中，变量 x1和 x2的系数最大，说明第一主成分受控于转氨酶和肝大指数，这两项指标的数值大，意味着肝病患者处于急性炎症状态。第二个主成分表达式中，变量x3的系数最大，说明第二主成分受控于硫酸锌浊度，这项指标的数值大，意味着肝病患者处于慢性炎症状态。第三个主成分表达式中，变量x4的系数最大，说明第三主成分受控于胎甲球，这项指标的数值大，意味着肝病患者可能处于肝癌可疑状态。第四个主成分表达式中，变量 x1和 x2的系数的绝对值最大，但符号相反，由于第四个主成分的贡献率很小（5.170%），仅供参考，临床上认为可能指向急性肝萎缩。患者分类的主要输出结果见表6。由1类4类列可知，第一类包含编号为4、8、9、11、15、17这6位肝病患者，他们属于急性炎症肝病患者；第二类包含编号为7、12、13、14、16这5位肝病患者，他们属于慢性炎症肝病患者；第三类包含编号为2、3、5、6、10、19这6位肝病患者，他们属于疑似肝癌患者；第四类包含编号为1、18、20这3位肝病患者，他们属于疑似肝萎缩肝病患者。结论：在本例中，求出了4个主成分，从第一到第四个主成分的贡献率依次为 43.0%、27.3%、24.5%和5.2%。将每位肝病患者在4项定量指标上的取值分别代入所求出的4个主成分表达式中进行计算（注意：对原变量需要做标准化变换），并按绝对值最大作为分类的“标准”，可得到4类分类结果。3.2.2分析例2中的资料设所需要的SAS过程步程序如下：proc princomp out=aaa prefix=z；var X1-X8；run；data a2；set aaa；zt=0.531*z1+0.157*z2+0.121*z3+0.086*z4；zt=round（zt，0.001）；proc rank descending data=a2 out=bbb；var zt；ranks order；run；proc print data=bbb noobs；var jour zt order；run；【SAS程序说明】由proc princomp过程计算的结果可知，前4个主成分的累计贡献率为89.45%，大于85.00%，故取前4个主成分，并以它们各自的贡献率为权重，进行加权平均，求出与每种期刊对应的一个综合评价指标zt的值。实现的语句为：“zt=0.531*z1+0.157*z2+0.121*z3+0.086*z4；”，该语句中的系数分别代表前4个主成分的贡献率。【SAS 输出结果及解释】因篇幅所限，由 proc princomp过程输出的结果从略。下面仅输出基于主成分计算结果进行综合评价的结果，见表7。两列“总分”代表各期刊的综合评分值；两列“排序”代表期刊的排序位次，数值越小，表明综合评分值越大。结论：在本例中，基于8项文献计量学指标的数值实施了主成分分析，获得了累计贡献率达到89.45%的前4个主成分；再基于这4个主成分以及它们各自的贡献率，采取加权算法获得综合评价指标zt；进而，基于zt的表达式，算出每种期刊的综合得分值，并对综合得分值排序，最终获得了全部23种肿瘤类期刊的排列顺序。这是将“无序样品（本文指 期刊）”转化成“有序样品”的过程，属于传统综合评价方法的一种扩展9。表620例肝病患者分成4类的结果Table 6Results of 20 liver disease patients divided into 4 categories患者编号12345678910最大值0.1850.4300.7762.0760.9491.0260.8022.2932.0221.2121类-4-89-2类-7-3类-23-56-104类1-患者编号11121314151617181920最大值1.1212.1090.3373.0240.7070.4830.232-0.1181.3970.3411类11-15-17-2类-121314-16-3类-19-4类-18-20注：每位患者可基于4个主成分计算出4个得分值，取其中最大值作为该患者的最终得分值53四川精神卫生 2023 年第 36 卷第 S1期http：/4 讨论与小结 4.1讨论主成分分析方法是众多多元统计分析方法的基础，虽然计算过程和方法比较简单，但从它的计算结果中可以引申出多种非常有价值的应用。由本文的例1可知，主成分分析可以被用来实现“样品聚类”；由本文的例2可知，主成分分析可以被用来实现“传统综合评价”，使无序样品变成有序样品，从而实现样品排序。事实上，主成分分析还有许多其他应用，因篇幅所限，不再赘述。4.2小结本文介绍了与主成分分析有关的 5 个基本概念、计算方法和两个实例以及SAS实现。5个基本概念包括相关矩阵、特征值与特征向量、主成分变量、主成分表达式和主成分的性质；计算方法涉及特征值与特征向量的求法、主成分分析的计算原理、系数估计和个数确定；两个实例涉及的资料分别是“20例肝病患者的4项肝功能指标的测定结果”和“23 种肿瘤类期刊的文献计量学指标的调查结果”。基于两个实例，分别实现了样品聚类和样品排序。参考文献1Rice JA.Mathematical statistics and data analysis M.2 版.北京：机械工业出版社，2003：129-134.Rice JA.Mathematical statistics and data analysisM.2nd edition.Beijing：China Machine Press，2003：129-134.2Johnson DE.应用多元统计分析方法 M.北京：高等教育出版社，2005：77-92，93-146.Johnson DE.Applied multivariate methods for data analysis M.Beijing：Higher Education Press，2005：77-92，93-146.3余锦华，杨维权.多元统计分析与应用 M.广州：中山大学出版社，2005：189-209.Yu JH，Yang WQ.Multivariate statistical analysis and applicationM.Guangzhou：Sun Yat-sen University Press，2005：189-209.4郭大钧.大学数学手册 M.济南：山东科学技术出版社，1985：300-302.Guo DJ.Handbook of college mathematics M.Jinan：Shandong Science and Technology Press，1985：300-302.5Johnson RA，Wichern DW.实用多元统计分析 M.6 版.北京：清华大学出版社，2008：430-480.Johnson RA，Wichern DW.Applied multivariate statistical analysisM.6th edition.Beijing：Tsinghua University Press，2008：430-480.6胡良平.现代统计学与SAS应用 M.北京：军事医学科学出版社，1996：316-323.Hu LP.Modern statistics and SAS applicationsM.Beijing：Military Medical Science Press，1996：316-323.7胡良平.面向问题的统计学：（3）试验设计与多元统计分析M.北京：人民卫生出版社，2012：19-39.Hu LP.Problem-oriented statistics：（3）experimental design and multivariate statistical analysisM.Beijing：Peoples Mental Publishing House，2012：19-39.8SAS Institute Inc.SAS/STAT15.1 users guide M.Cary，NC：SAS Institute Inc，2018：7891-7932.9胡良平，黄国平.医学科研设计方法与关键技术 M.成都：四川大学出版社，2017：349-360.Hu LP，Huang GP.Medical research design methods and key technologiesM.Chengdu：Sichuan University Press，2017：349-360.（收稿日期：2023-06-05）（本文编辑：吴俊林）表7对23种肿瘤期刊进行综合评价的结果Table 7Comprehensive evaluation results of 23 cancer journals期刊编号123456789101112总分3.5462.0441.3930.8790.398-0.273-0.7670.1100.092-0.391-0.184-0.399排序12346121878131114期刊编号1314151617181920212223-总分-0.0740.749-0.655-0.756-0.182-0.813-0.892-0.738-0.964-0.969-1.154-排序95151710192016212223-54