因子分析概述和应用简例

1、1因子分析概述和应用简例1 因子分析概述4 方差最大正交旋转5 因子得分6 对应分析5 应用简例2 R型因子分析3 主因子的解7 Q型因子分析21 因子分析概述一、因子分析的基本概念 1900年C.斯皮尔曼发表了对学生考试成绩分析的著名文章，可认为是因子分析的开始。1957年由Krumbein把该方法引入到地质学研究，后来，又由Imbrie对该方法在地质学中的应用和发展作了大量的基础工作。 因子分析是研究变量间相关关系、样品间相似关系、变量与样品间成因联系以及探索它们之间产生上述关系之内在原因的一些多元统计分析方法的总称。 根据它们的研究对象，因子分析大致可分为R型因子分析、Q型因子分析和对应

2、分析三种类型。在此仅介绍基于相关系数和相似系数统计量下的因子分析。3(一)变量的相关矩阵设n个样品m个变量观测值标准化后的矩阵形式如下（在因子分析后续的讨论中，若无特殊说明，数据均约定是标准化的）：对于任意两个变量xi(xi1, xi2, , xin)和xj(xj1, xj2, , xjn)，它们的相关系数为： 4m个变量的相关系数构成的mm的矩阵为： 且rij=rji，r11r22rmm1，并称R为变量的相关矩阵 5(二)样品的相似矩阵 任意两个样品xi(x1i, x2i, , xmi)和xj(x1j, x2j, , xmj)的相似系数为 n个样品的相似系数矩阵形式为 且q11q22qnn1

3、，qijqji，称Q为样品的相似矩阵 6(三) R型因子分析 R型因子分析是从研究变量的相关(相似)系数矩阵R的内部结构出发，找出m个控制变量相关性的主因子fi (i=1,2,m)，由此把 xi 表示为fi的线性组合，即 ： 在进行综合地质研究时，如果用前p(p m)个主因子就能解释原始数据8090%以上的信息，那么上式可改写为：7 用上式可极大地化简研究系统(降维分析问题，特别是当p=2时，可以在二维空间中绘图)进一步探索变量的成因联系及空间变化规律的控制因素。通常称上式为R型因子分析模型。 Q型因子分析是从研究样品相似(相关)系数矩阵R的内部结构出发，找出n个控制样品相似性的主因子，若用前

4、p( pn)个主因子化简样品研究系统，得到如下Q型因子分析模型:（四） Q型因子分析8f1, f2, f3,称为公因子, ei称为特殊因子。 因子分析模型中,各变量(或样品)中共同出现的因子fi (i=1,2, p)称公因子, 是相互独立的理论变量,可将其理解为p 维空间中相互垂直的p个坐标轴。 单一变量(或样品)中特有的ei(i=1,2,m或n)叫做特殊因子,它们之间以及它们与所有公因子之间都是相互独立的。 它是把上述两种因子分析方法结合起来，在同一个空间里研究样品和变量的分布规律，更有利于进行地质解释的一种多元统计分析方法。 （五）对应分析9二. 因子分析在地质研究中的作用 第二，指示成因

5、推理方向。因子分析能够把庞杂纷乱的原始数据按成因上的联系进行归纳、整理、精炼和分类，理出几条客观的成因线索，为地质人员提供逻辑推理方向，启发思考相应的成因结论。 第三，分解叠加的地质过程。现实观测到的地质现象往往是多种成因过程叠加的产物，因子分析提供了一个分解叠加过程进而识别每个单一地质过程的手段。 第一，压缩原始数据。因子分析能在数量上大大精简原始数据但又不损失数据中包含的成因信息。从而有利于地质人员进行综合分析。10三、 因子分析在地质研究中的应用 沉积盆地物源区的研究、沉积物的粒度分析、沉积相研究、地层分析、古环境与古生态的研究、石油与天然气的成因研究、油田水的化学研究、有机地球化学研究

6、及石油、天然气物化探资料分析等各个方面。11 假设有n个样品，每个样品有m个变量，它们的观测值记为数据矩阵：1 因子分析概述 一、主因子 当m=2时，在x1、x2坐标系内作出样品的散点图如图9-1所示。(9-1)12x1x2图9-1 主因子示意图 由图可见, 在x1、x2坐标系内, x1、x2的方差几乎没有差别。 若将坐标系旋转一个角度，得到新坐标系 f1、f2 , 这时f1、f2的方差具有明显差别, f1的方差占了总方差的绝大部分。 因此新变量f1能够反映原始数据所包含的绝大部分信息。由坐标旋转变换知新老变量间的关系为:或(9-2)f2f113 若f2的方差小于允许误差时,它就失去了存在的必

7、要,此时就可用f1代替两个原始变量,而式(9-2)改为:或(9-3)其中ei服从均值为0，方差为2的正态分布。 (9-4)(i=1,2,m) 如果每个样品有m个变量，仿照式(9-2)可以写出m个新变量： 可以证明，上式系数的平方和满足：14 把m个原始变量表示为p个主因子的线性组合，当p小于m，特别是 p=2时，可以在二维空间对变量作图，进而对变量的相关性及成因联系进行研究。 在此条件下，由原始变量经线性组合而得到的新变量fi叫做主因子或综合变量。 组合成新变量有什么用途呢？15二、因子分析 因子分析是研究变量间的相关性、样品间的相似性、两者成因联系、探索它们之间产生相关(相似)关系之内因的一

8、些多元统计分析方法的总称。 因子分析的任务之一是找出p( pm)个主因子，以此化简研究系统。 根据研究对象不同,将因子分析分为R型、Q型因子分析和对应分析三种方法。161. R型因子分析 (9-5) R型因子分析是从研究变量的相关(相似)系数矩阵R的内部结构出发，找出m个控制变量相关性的主因子fi (i=1,2,m)，由此把 xi 表示为fi的线性组合，即 ： 其中称因子载荷矩阵。17(9-6) 在进行综合地质研究时，如果用前p(p m)个主因子就能解释原始数据8090%以上的信息，那么式(9-5)可改写为： 用上式可极大地化简研究系统(降维分析问题，特别是当p=2时，可以在二维空间中绘图)进

9、一步探索变量的成因联系及空间变化规律的控制因素。通常称式(9-6)为R型因子分析模型。18 Q型因子分析是从研究样品相似(相关)系数矩阵R的内部结构出发，找出n个控制样品相似性的主因子，若用前p( pn)个主因子化简样品研究系统，得到如下Q型因子分析模型:3.公因子2. Q型因子分析 因子分析模型中,各变量(或样品)中共同出现的因子fi (i=1,2, p)称公因子, 是相互独立的理论变量,可将其理解为p 维空间中相互垂直的p个坐标轴。(9-7)f1, f2, f3,称为公因子, ei称为特殊因子。19 单一变量(或样品)中特有的ei(i=1,2,m或n)叫做特殊因子,它们之间以及它们与所有公

10、因子之间都是相互独立的。 它是把上述两种因子分析方法结合起来，在同一个空间里研究样品和变量的分布规律，更有利于进行地质解释的一种多元统计分析方法。 因子载荷矩阵中的aij是第i个变量(或样品)在第j 个因子轴上的负荷。如果把xi视为p维空间的一个向量，那么aij 则是xi在fj轴上的投影。 4.对应分析20 由上述可知，进行因子分析的关键是确定主因子载荷矩阵，在它的基础上选定前p个主因子则可得到因子分析模型。21四、基本定理及主因子解 因子分析的一个基本问题是要由变量之间的相关系数矩阵R去确定因子载荷矩阵A。由前所述有 1.主因子解 222324 R*叫做约相关矩阵。因子分析的基本问题之一就是

11、在已知约相关矩阵R*的条件下求解因子载荷矩阵A，使得R*AA。这也是因子分析的基本定理。 实际中a1,a2,am很难确定，R*也不知道。因此常常用R代替R*， 即由RAA求解因子载荷矩阵A。满足RAA的因子解很多,如果A是它的解，那么A与正交矩阵T的乘积矩阵AT也是它的解，这是因为 这就告诉我们有许许多多的因子解。25 按下面的方法求得的因子解称作主因子解。 其做法是:根据变量的相关选出第一个因子 F1，使其在各变量的公因子方差中所作的方差贡献最大，然后消去该因子影响，再从剩余的相关中选出与F1不相关的因子F2，使其在各变量的剩余公因子方差中方差贡献最大，如此继续挑选直至各变量公因子方差分解完

12、毕为止。比如，选取第一个主因子F1，要使它的方差贡献S1ai12(i=1,2,.,m)在R=AA之下达到最大，可用条件极值中的拉格朗日乘数法求出。这里不详细给出其求解过程，只给出它的结果如下： 26 设变量相关系数矩阵R的特征值12, , p p+1 , m 0 对应的特征向量为如果选定前p个主因子，那么可以证明前p个主因子载荷矩阵为:27相应的R型因子分析模型为：282.主因子个数 用Jacobi法可求出相关(相似)系数矩阵的全部特征值和单位特征向量。若m个特征值从大到小排列为1,2,m，计算特征值累计百分比：然后根据C确定主因子数p。293.例在冀北坳陷东部某烃源岩层采得烃源岩样品4个，每

13、个样品测得3个地化指标，共12个原始数据，将原始数据的各指标变量标准化，得到标准化数据阵为:要求对上述标准化数据进行主成分分析，考查各指标变量间关系，计算各主因子的因子计量和因子贡献。 30解：1、计算三个变量的相关系数矩阵2、计算主因子解，首先求得R的特征值为，方差贡献的累计百分数分别为、1。 313、确定公因子数p=2 前两个特征值所对应的特征向量为如下矩阵的两个行向量 得到主因子荷载矩阵为32二、Q型主因子载荷矩阵及主因子个数 设相似系数矩阵Q的特征值12p+1n对应的特征向量为V1V2Vp+1Vn。如果选定前p个主因子，那么可以得到与R型主因子载荷类似的Q型主因子载荷矩阵为:相应的Q型

14、因子分析模型为：1. 主因子载荷矩阵332.主因子个数 用Jacobi法可求出相似(相关)系数矩阵的全部特征值和单位特征向量。若n个特征值从大到小排列为1,2,n, 计算特征值累计百分比：然后根据C确定主因子数p。343 方差最大正交旋转与因子得分一、R型因子分析方差最大正交旋转 确定了p 个主因子后，为便于解释主因子的地质意义，还要对主因子轴进行旋转,使第j个公因子的代表性变量在fj轴上的系数等于或趋近于1，在其它公因子轴上的系数等于或趋近于0。1.方差最大正交旋转的实现方法 从数学上看,上述即对主因子矩阵A1实施正交变换。常用的方法是方差最大正交旋转,这种旋转方法使p个因子保持彼此正交,并

15、且使因子载荷矩阵中的各因子载荷平方后的方差达最大。对于R型因子载荷矩阵A1来说,对因子fj的简化效果,用因子载荷平方的方差进行描述:35当V 达到最大时，p个因子都得到了最优简化。(9-8) 其中bij是经过正交旋转后的因子载荷矩阵B的元素,为避免负值取bij2。若Vj为最大,即第j个因子得到了最优简化,此时第j个公因子的代表性变量在fj轴上的系数等于或趋近1,而在其他公因子轴上的系数等于或趋近于0。对p个因子的简化效果用因子载荷平方的方差之和进行度量：36 考虑到各个研究对象的公因子方差间的差异,可以用(bij2/hi2)代替bij2,并对式(9-8)两边乘以m，记为：(9-9) 实际是要求

16、V*达到最大。于是问题归结为求一个正交矩阵Tpp, 使B=A1T 满足V* 为最大的条件。 对于因子平面fg fq可以对A1作如下正交变换:1 . . . cos -sin . . . sin cos . . . 100Tgq=gqgq37B中的元素分别为： 对A1进行一次正交旋转变换,相当于对因子平面fgfq旋转了一个角度,每次旋转必须满足使式(9-9)的值为最大,得到矩阵: 对于选定的p个主因子,必须对A1中所有p列配对旋转,总共需要旋转p(p-1)/2次,全部旋转完毕得到载荷矩阵:(9-10)38 由B2计算出V2* 重复循环旋转，可以得到一个有界序列： V1*V2* 对于任意小正数，当

17、|Vk*-Vk+1*|时,结束计算,最后得到: 经过一个循环旋转得到B1 ,由此可按式(9-9)计算V1*。 在上个循环旋转的基础上，从B1出发完成下个旋转循旋转得到B2, 即:39 前面介绍了实现方差最大正交旋转的的基本思想，但在任何一次旋转变换中，都涉及旋转角度的问题。确定旋转角度的具体步骤如下：2.旋转角度(不讲) (1)把式(9-10)代入式(9-9)，得到函数V*();式中:; (2)求V*对的一阶导数并令其为0，可以解得:(9-11)40 (3)将函数V*()展开化简，剩下包含sin4、sin22 的项，使V*()成为以/2位周期的函数，因而4在/2范围内，通常在-/4 /4之间。

18、 (4)由V*对的二阶导数小于0，可得：由(9-12)可知,与E同号。 Q型因子分析方差最大正交旋转与R型中的旋转方法类似,不同之处是因子载荷矩阵A1=(aij)n p，因此,只要将上述公式中的m改为n即可。(9-12)41二、因子得分 因子分析是将变量(或样品)表示为公因子的线性组合。当然,也可以将公因子表示为变量(或样品)的线性组合,即:(9-13) 称fj为因子得分函数,表示对原始变量信息量的聚集量。如果取p=2,则可将每个样品的m个变量代入(9-13),计算出因子得分fi、fj ( ij , i,j=1,2, p),在二维空间作图对样品进行分类或者进行地质解释。1. 因子得分函数422

19、.因子得分函数的确定若公因子已提取全部信息的80%以上，则可采用简化因子模型计算：F=(AA)-1AX有些情况下(如公因子相对较少，提取的信息比较少)，直接计算因子得分的不完全适用，可以采用回归估计的方法可以采用。rij是原始变量的相关系数。43三、 因子分析的步骤以R型因子分析为例说明因子分析计算步骤(1)建原始数据矩阵。(2)对原始数据作数据标准化变换。(3)计算m个变量间的相关系数，建相关系数矩阵R。(4)求出R的特征值(12,m0并按大小排列)及相应于j的单位特征向量uj (j=1,2,m)。(5)确定公因子个数p。 (6)求出主因子载荷矩阵A=aij , (7)计算诸公因子方差hi2

20、。(8)将因子载荷矩阵A作方差最大正交旋转，求出旋转后的因子载荷矩阵，仍记为B。(9)计算因子得分。444 对应分析 R或Q型因子分析用最少的公因子提取研究对象的绝大部分信息，简化研究系统，从而在低维空间中研究样品的空间分布规律和变量成因联系，进行地质解释和推断。因此这两种因子分析在地质学中有着广泛的应用。但它毕竟存在以下不足： R型和Q型因子分析把研究样品空间分布规律和研究变量共生组合关系分隔开来。事实上，样品的特征经常由变量来揭示。如储层的渗透性问题，既要研究渗透性地质成因的空间分布规律，又要研究储层样品的特征，那么前者一、对应分析的概念不能在同一空间研究样品和变量45就要研究样品，而后者

21、就要利用变量来解释。这说明样品和变量的研究是不可分割的，也就是说应该找出一种把两种因子分析统一起来的研究方法。 地质研究中，样品数量通常远比变量数目多，计算样品相似系数时既让费机时，又占用大量的内存空间。 为了使地质数据在同一尺度下参与地质分析，往往是将变量进行标准化处理，然而对样品就不好作标准化了。 因此对变量和样品的非对等标准化处理，造成地质数据的尺度不同，影响对地质规律的认识。占机时和内存多地质数据的尺度不同46 鉴于上述原因，在两种因子分析的基础上又发展了一种新的多元统计分析方法对应分析。它采用一种数据处理方法,把两种因子分析结合起来，在同一空间内对变量和样品进行划类，由此既可研究样品

22、的分布规律，又可通过变量对样品进行地质解释。二、对应分析的数据变换 假设有n个样品,每个样品有m个变量,它们的观测值记为:47 其中xij0(i=1,2,m; j=1,2,n),并且每一行和每一列上至少有一个数不为0。1.原始数据的总和标准化矩阵X中各元素之和为 对矩阵X中各元素除以T,称数据总和标准化。结果如:48记P中第j列元素之和为第j列上的每个元素再除以p. j ，得:2.样品坐标的相对比例49 3. 样品点间的加权距离 上式表示m维空间中的n个样品点，每个样品点的坐标是样品中各变量的相对比例。 考虑到指标数量级的差异，采用加权距离表示样品点l与k之间的距离,以此表示它们之间的差异。(

23、第i行的和) 由上可知，为了计算样品点之间的加权距离，只要将坐标改为：50(9-14)4.变量坐标的相对比例矩阵P中第i行元素之和：第i行上的每个元素再除以pi.，得:51其中 上式表示n维空间中的m个变量点，每个变量点的坐标是样品中变量的相对比例。 52 5.变量点间的加权距离 上式表示n维空间中的m个变量点，任意两点l与k之间的加权距离为： 由上可知，为了计算变量点之间的距离，只要将坐标改为：(9-15)即可。53 为了对样品进行地质解释，利用式(9-14)和(9-15)进一步研究样品和变量的关系。(9-14)(9-15)54三、协方差矩阵1. 变量的协方差矩阵(1) 变量的均值 若将矩阵

24、P中的元素pij视为概率,那么pi. 、p.j 就是边缘概率,因此m维空间中样品点第i个变量的概率均值为： 55 (2)变量的协方差矩阵 第i与第j个变量的协方差为其中56那么m个变量的协方差矩阵为：其中Z=(zik)mn 2.样品的协方差矩阵 (1)样品的均值 在n维空间中，第k样品的概率均值为：57 (2)样品的协方差矩阵 第l和第k个样品的协方差为其中58那么n个样品的协方差矩阵为：其中 四、因子载荷矩阵 由线性代数可知，矩阵ZZ与ZZ有相同的非零特征值12p(pm)，并且每个j(1jp)，若对应的 uj是ZZ的单位特征向量，那么vj=Zuj 是ZZ 相应的单位特征向量；反之，若vj是Z

25、Z的特征值j对应的单位特征向量，那么uj=Zvj则是ZZ 所相应的单位特征向量。59 上述结果表明，当求得变量的协方差矩阵S的特征值j (j=1,2,., p)和与其对应的特征向量uj(j=1,2,., p)后，便可得到R型因子分析的因子载荷矩阵，再由vj=Zuj可直接求得Q型因子分析的因子载荷矩阵。 此外，S与S*有相同的特征值,这些特征值表示各公因子所提供的方差,因而变量空间中的p个因子与样品空间中的p个因子在总方差中所占的百分比完全相同，因此用同样的因子轴既可以表示变量，又可以表示样品，把R型和Q型因子分析统一起来。60 1. R型因子分析因子载荷矩阵 若取S 的前p个特征值j( j=1

26、,2,p)，与之对应的特征向量为uj( j=1,2,p)，那么因子载荷矩阵为： 2. Q型因子载荷矩阵 Q型因子载荷矩阵为：61其中： 五、对应分析计算步骤 1. 求Z矩阵 把原始数据矩阵X变换为Z矩阵62 2. R型因子分析 求变量的协方差矩阵S的特征值j(j=1,2,m) ，按其累计百分比取前p个主因子，计算因子载荷矩阵。 4. 制作对应分析图 在因子平面内以因子载荷为坐标所做的变量和样品散点图。其中： 根据vj=Zuj 计算Q型因子载荷矩阵。3. Q型因子分析63样品号石英含量%长石含量%岩屑含量%其他含量%鉴定命名123456789101112131415161718828386887

27、8847661605857554746585361508878759262529273027257171217995121211131315131615262935302533102230200000000020石英砂岩石英砂岩石英砂岩硬砂质石英砂岩硬砂质石英砂岩硬砂质石英砂岩硬砂质石英砂岩长石砂岩长石砂岩长石砂岩长石砂岩长石砂岩混合砂岩混合砂岩硬砂岩长石质硬砂岩长石质硬砂岩长石质硬砂岩5 应 用 简 例 例1 砂岩分类图9-2砂岩样品与碎屑对应分析平面图石英长石岩屑石英砂岩硬砂岩杂砂岩长石砂岩f2f164 例2 古潜山油气藏与成油地质参数分析 对已探明的一些古潜山油气藏进行划类，研究它们的资

28、源量与成油地质参数的关系，是评价未知古潜山油气资源量的基础。与其它类型的油气藏一样，油源、供油、储集、圈闭和保存条件等是形成古潜山油气藏的基本条件，将其拟定为8项成油地质参数： 由于存在同一个生油凹陷向不同古潜山圈闭提供油气的问题，因此本例未把油源条件列入成油地质参数。对已探明的26个古潜山油气藏及其的成油地质参数进行对应分析，结果如下：序号地质参数名称1古潜山到凹陷生油区中心的距离，km;2储集体的总孔隙度，%;3古潜山上伏盖层的厚度，km;4供油窗口，km2;5生油岩与不整合面的接触面积，km2;6圈闭闭合高度，km;7圈闭闭合面积，km2;8古潜山最浅埋藏深度，km.65图9-3 油气藏

29、与成油地质参数对应分析平面图距生油中心远，古潜山不整合面与生油岩接触面积小，断层供油，总孔隙度低(5百万吨)。距生油中心近、靠不整合面输导形成的油气藏(500百万吨)距生油中心较近，古潜山不整合面与生油岩接触面积较小，总孔隙度高 (10百万吨) f2f166 例3 云南某地钾盐矿床水化学资料对应分析1. 目的 2. 资料 20个三种类型的盐泉样品,每个盐泉有7项水化学分析指标,分别是:矿化度(g/l)、Br/CI(103)、K/盐(103)、K/CI(103)、Na/K、Mg/CI(102)、Na/CI。原始数据见下表。对应分析图如图9-4。 了解样品和变量的关系,合理评价盐泉及其成因联系。6

30、7盐泉序号矿化度/g/l(Br/Cl)103(K/盐)103（K/Cl）103Na/K(Mg/Cl)102Na/Cl123456789101112131415161718192011.85345.5963.5253.681 48.28 17.9567.370 4.2236.44216.234 10.585 23.5355.398 283.148316.604 307.310 322.515256.580304.092240.4460.4800.5260.0860.3700.386 0.2800.506 0.340 0.1900.3900.420 0.2300.1200.1480.3170.1730.3120.2970.2830.04214.36013.85024.40013.57014.5009.75013.6003.8004.7003.4002.4002.6002.8001.7631.4531.6271.3820.8990.7890