第8章 Stata主成分分析与因子分析

     

目录

8.1主成分分析

8.2因子分析

1.主成因子法

2.主因子法

3.迭代公因子方差的主因子法

4.最大似然因子法

---

     在进行数据统计分析时，还往往会遇见变量特别多的情况，而且很多时候这些变量之间还存在着很强的相关关系或者说变量之间存在着很强的信息重叠，如果我们直接对数据进行分析，一方面会带来工作量的无畏的加大，另一方面还会出现一些模型应用的错误，于是主成分分析与因子分析应运而生。这两种分析方法的基本思想都是在不损失大量信息的前提下，利用较少的独立变量来替代原来的变量进行进一步的分析。

8.1主成分分析

        基本思想是将众多的初始变量整合成少数几个互相无关的主成分变量，而这些新变量尽可能地包含了初始变量的全部信息，然后利用这些新的变量来代替以前的变量。

        数据（案例8.1）一共V1-V19 19个变量，分别代表了年份、全国人口（万人）、农林牧渔业总产值（亿元）、工业总产值（亿元）、国内生产总值（亿元）、全社会投资总额（亿元）、货物周转量（亿吨千米）、社会消费品零售总额（亿元）、进出口贸易总额（亿元）、原煤（亿元）、发电量（亿千瓦时）、原油（万吨）、钢（万吨）、汽车（万辆）、布（亿米）、糖（万吨）、粮食（万吨）、棉花（万吨）、油料（万吨）。使用主成分分析对这些指标提取主成分并写出提取主成分与这些指标之间的表达式。

correlate V2-V19  #对18个变量进行相关性分析

        结果如图所示。这里得到的是所有变量之间的方差-协方差矩阵。我们可以发现很多变脸之间的相关关系是非常强的甚至有的超过了90%，这说明变量之间存在着相当数量的信息重叠。我们进行主成分分析把众多的初始变量整合成少数几个相互之间无关的主成分变量是非常必要的。

pca V2-V19 #对18个变量进行主成分分析

        结果如图所示，（Component）表示的是系统提取的主成分名称，可以发现，我们stata总共提取了18个主成分。Eigenvalue列表示的是系统提取的主成分特征值，特征值的大小意味着该主成分的解释能力，特征值越大解释能力越强，可以发现Stata提取的18个主成分中只有前7个是有效的。Proportion列表示的是系统提取的主成分的方差贡献率，方差贡献率也便是主成分的解释能力，可以发现第一个主成分的方差贡献率是0.8023，表示该主成分解释了所有变量的80.23%的信息。

        这个图站试的是主成分特征向量矩阵，以表明各个主成分在各个变量上的载荷，从而可以的出主成分的表达式。在表达式中各个变量已经不是原始变量而是标准变量。其中第一个 特征值的主成分表达式是：

        comp1 = 0.1377*V2+0.2605*V3+···+0.1913*V19

        在第一主成分中，除了粮食变量（V17）以外的变量系数比较大，可以看成是反映哪些变量的综合指标，在第二主成分中，粮食变量系数比较大，可以看作是反应粮食的综合指标。因为主成分分析只不过是一种矩阵变换，所以各个主成分并不一定具有实际意义，本例中各个主成分的内在含义就不是很明确。

案例延伸：

1.只保留特征值大于1的主成分

pca V2-V19,mineigen(1)

        结果如图所示，值得说明的是，上图最后一列Unexplained表示的是该变量未被系统提取的两个主成分的解释的信息比例，例如V2未被解释的信息比例就是72.55%这种信息丢失的情况i昂是我们舍弃其他主成分必然付出的代价。

2.限定提取的主成分个数

pca V2-V19,components(1)

        图例与结果不再过多阐释。

8.2因子分析

        因子分析在一定程度上可被视做主成分分析的深化和拓展，它对相关问的研究更为深入透彻。因子分析的基本原理是将具有一定相关关系的多个变量综合为数量较少的几个因子，从而研究一组错综复杂关系的实测指标是如何受少数几个内在的独立因子所支配的，所以它属于多元分析中处理降维问题的一种常用统计方法。

        因子分析法有很多种，在这里我们介绍四种，分别是主成因子法、主因子法、迭代公因子方差的主因子法、最大似然因子法。

        数据（案例8.2）本例中有7个变量，分别是年份、工业总产值、国内生产总值、货物周转量、原煤、发电量和原油。我们把这些变量分别设为V1-V7。试用因子分析法对这些指标提取公因子并写出提取公因子与这些指标之间的表达式。

1.主成因子法

factor V2-V7 , pcf ##本命令是对这6个变量进行因子分析

        结果如图所示，我们可以看到一共有9（Number of obs = 9）个样本参与了分析，提取保留的因子共有两个（Retained factors = 2），模型LR检验的卡方值为100.47，P值为0.0000（ LR test: independent vs. saturated:  chi2(15) =  100.47 Prob>chi2 = 0.0000），可以看出模型是非常显著的。图上半部分最左列（Factor）说明的是因子名称，可以看出模型共提取了6个因子 Eigenvalue列表示的是提取因子的特征情况，只有前面两个因子的特征值是大于1 的，Proportion列表表示的是提取因子的方差贡献率。Cumluative列表示的是提取因子的累计方差贡献率。

        图下半部分说明的是模型的因子载荷矩阵以及比变量未被解释部分。其中，Variable列表示的是变量名称，Factor1、Factor2两列分别说明的是提取的前两个主因子（特征值大于1的）对各个变量的解释程度。本例中Factor1主要解释的是V2\V3\V4\V6\V7这5个变量的信息，Factor主要解释的是V5变量的信息。Uniqueness列表示变量未被提取的前两个主因子解释的部分，可以发现在舍弃其他主因子的情况下，信息的损失量是很小的。

rotate #本命令的含义是对因子结构进行旋转

        结果如图所示，上图展示的是对因子结构进行旋转的结果。经过学者们的眼界表明，旋转操作有助于进一步简化因子结构。Stata支持的旋转方式有两种；一种是最大方差正交旋转，一般适用于相互独立的因子或者成分，也是系统默认的情况；另一种是promax斜交旋转，允许因子或者成分子间存在相关关系。此处我们选择系统默认的方式，当然我们后面的操作也证明了这样做的恰当性。

        上图包括三部分内容，第一部分说明的是因子旋转模型的一般情况，从图中我们可以看出共有9个样本参与了分析，提取保留的因子共有2，模型LR检验的卡方值为100.47，P值为0.0000，模型非常显著。最左列（Factor）说明的是因子名称，可以看出模型旋转后共提取了两个因子。Proportion列表示的是提取因子的方差贡献率，

        第二部分说明的是模型的因子载荷矩阵以及比变量未被解释部分。其中，Variable列表示的是变量名称，Factor1、Factor2两列分别说明的是旋转提取的两个主因子对各个变量的解释程度。本例中Factor1主要解释的是V2\V3\V4\V6\V7这5个变量的信息，Factor主要解释的是V5变量的信息。Uniqueness列表示变量未被提取的前两个主因子解释的部分，可以发现在舍弃其他主因子的情况下，信息的损失量是很小的。

        第三部分是因子旋转矩阵的一般情况，提取的两个因子不存在相关关系。

loadingplot,factor(2)yline(0)xline(0) #本命令的含义是绘制因子旋转后的因子载荷图。

        我们可以看出factor1主要解释的是V2\V3\V4\V6\V7,factor2主要解释的是V5

predict f1 f2 #本命令的含义是显示因子得分系数矩阵

        这样的话我们就可以写出各公因子的表达式，值得一提的是在表达式中各个变量已经不是原始变量而是标准化变量。

        F1=0.19062*V2+0.21609*V3+···+0.21185*V7

        F2=0.31358*V2+···

list V1 f1 f2 #本命令的含义是估计因子分析后各个样本的因子得分情况，请看一下数据浏览就知道了，这个操作只是显示。 真正的估计命令是上面的 predict f1 f2

        结果如图所示，展示的是因子得分数据。

correlate f1 f2 #本命令的含义是查看两个主因子之间的相关系数矩阵

        我们可以看到这两个主因子之间几乎没有任何关系，这也说明了我们在前面对因子进行旋转的操作环节中采用最大方差正交旋转方式是明智的。值得说明的是途中f1与f2的相关系是-0.0000，并非不正确，这是因为stata只保留了四位小数所导致的，如果真是数据是-0.0001那么结果显示的是-0.0000。

scoreplot,mlabel(V1)yline(0)xline(0) #本命令的含义是展示每个样本的因子得分示意图

        从图中可以看出所有样本被分到了4个象限，其中第一象限包括了96、97、98年这三年的两个因子得分比较高，第二象限包括94、95年其中这两个年在因子2上得分较高，在因子1分上得分较低。后面的不再阐述。

estat kmo #本命令展示的是本例因子分析的KMO检验结果。

        KMO检验是为了判单数据是否适合进行因子分析，其取值范围是0~1。其中0.9~1表示极好，0.8~0.9表示可奖励，0.7~0.8表示还好，0.6~0.7表示中等，0.5~0.6表示糟糕，0~0.5表示不可接受。在本例总体（Overall）KMO的取值为0.6566，表示可以进行因子分析。各个变量的KMO值也大多在0.6以上，所以本例是比较适合因子分析的，模型的构建是有意义的。

screeplot #本命令的含义是绘制因子分析的碎石图

        通过碎石图我们可以非常直观的观测出提取因子的特征值的大小情况。图的横轴表示的是系统提取因子的名称，并且已经按照特征值大小进行降序排列，纵轴表示因子特征值的大小情况。

2.主因子法

factor V2-V7,pf #本命令的含义是使用主因子法对V2-V7变量进行因子分析
 
rotate #本命令的含义是对因子结构进行旋转
 
loadingplot,factor(2)yline(0)xline(0) #本命令的含义是绘制因子旋转后的因子载荷图
 
perdict f1 f2 f3 f4 #本命令的含义是显示因子得分系数矩阵
 
list V1 f1 f2 f3 f4 #显示因子得分系数矩阵
 
correlate f1 f2 f3 f4 #展示提取的主因子的相关系数矩阵
 
scoreplot ，mlabel(V1) yline(0) xline(0) #展示每个样本的得分示意图
 
estat kmo #本命令的含义是显示KMO检验结果
 
screeplot #绘制因子分析的碎石图

3.迭代公因子方差的主因子法

factor V2-V7,ipf #本命令的含义是使用主因子法对V2-V7变量进行因子分析
 
rotate #本命令的含义是对因子结构进行旋转
 
loadingplot,factor(2)yline(0)xline(0) #本命令的含义是绘制因子旋转后的因子载荷图
 
perdict f1 f2 f3 f4 f5 #本命令的含义是显示因子得分系数矩阵
 
list V1 f1 f2 f3 f4 f5  #显示因子得分系数矩阵
 
correlate f1 f2 f3 f4 f5 #展示提取的主因子的相关系数矩阵
 
scoreplot ，mlabel(V1) yline(0) xline(0) #展示每个样本的得分示意图
 
estat kmo #本命令的含义是显示KMO检验结果
 
screeplot #绘制因子分析的碎石图

4.最大似然因子法

factor V2-V7,pf #本命令的含义是使用主因子法对V2-V7变量进行因子分析
 
rotate #本命令的含义是对因子结构进行旋转
 
loadingplot,factor(2)yline(0)xline(0) #本命令的含义是绘制因子旋转后的因子载荷图
 
perdict f1 f2 f3  #本命令的含义是显示因子得分系数矩阵
 
list V1 f1 f2 f3  #显示因子得分系数矩阵
 
correlate f1 f2 f3 #展示提取的主因子的相关系数矩阵
 
scoreplot ，mlabel(V1) yline(0) xline(0) #展示每个样本的得分示意图
 
estat kmo #本命令的含义是显示KMO检验结果
 
screeplot #绘制因子分析的碎石图

         最大似然因子法、主因子法、迭代公因子方差的主因子法和主成分因子法的结果解释几乎一样大家只需要稍微变通即可。留给自己思考和记忆的空间。不过最大似然因子法有几个需要提醒的地方：

        该检验有助于确定合适的因子数目。第一部分说明的是因子分析经过迭代计算后在第4次（Iteration 4:   log likelihood = -3.4269988）达到饱和，此时系统提取的主因子个数是3个。然后就是第二部分 BIC（Schwarz's BIC  =  39.8124）的值为39.8124，AIC（(Akaike's) AIC   =     36.854）为35.854，模型LR（LR test: independent vs. saturated:  chi2(15) =  100.47）检验的卡方值为100.47，P值为0.0000，意为模型是非常显著的