资源描述
吴 海 龙 湖南大学化学生物传感与计量学国家重点实验室 E-mail: hlwu,主成分分析 Principal Component Analysis,Chemometrics ？,The chemical discipline that uses mathematical, statistical and other methods employing formal logic (i) to design or select optimal measurement procedures and experiments, and (ii) to provide maximum relevant chemical information by analyzing chemical data.,一次测量得到一个值 例如： 温度、压力、pH、 单波长的吸光度 等,单变量数据,分析仪器的高性能化 例如： UV-可见分光吸收光谱 IR、 NIR、荧光光谱 GC、LC、MS NMR、等 样品-浓度数据 样品-变量-时间 等等,多变量数据,分析过程,测定-得到数据-数据解析-信息-结论 = 单变量（矢量型数据）： 平均值、标准差、数据检验、。 多变量（矩阵型数据）：主成分分析（抽象因子分析，数据简约） 、 聚类分析、方差分析、。 变量组合： （矢量-矢量、矢量-矩阵、矩阵-矩阵） 多元线性回归、典型相关分析、主成分回归、。 偏最小二乘回归,分析化学中常用方法,定性分析：分类、判别 定量分析：工作曲线法、多元校正、QSAR、。 建模： Ys= Xs B+ E, 如： Rs= Cs S + Es. 已知： Rs、 Cs，求： S. 预测： Yu= Xu B+ E, 如： Ru= Cu S + Es. 已知： Ru、 S，求： Cu.,1 Principal Component Analysis (PCA),主成分分析 (PCA),对一个矩阵，利用其变量之间的共线性， 对数据进行简约。这样， 可直观显示（图示） 可提取抽象因子 有效克服因严重共线性引起的不稳定算法带来的计算误差放大，即病态,PCA 例子,# B T Bmc Tmc 1 48 26 13 12 2 44 20 9 6 3 40 24 5 10 4 38 18 3 4 5 32 9 -3 -5 6 28 6 -7 -8 7 26 5 -9 -9 8 24 4 -11 -10 mean 35 14 0 0 B: 苯， T: 二甲苯 ； Bmc和Tmc为减去平均值后的值,利用GC得到的8个样品中苯和二甲苯的含量,Tmc vs Bmc,X矩阵中含有8个样品和两个变量 方差协方差阵为： 77.71 76.29 76.29 80.86 也即，它使用平均化的X矩阵XTX(列数-1) 特征矢量 # 特征值 B的系数 T的系数 方差百分数 1 155.59 0.6998 0.7144 98.1 2.98 0.7144 -0.6998 1.9 方差= /（样品数-1）,PCA 方差协方差矩阵,夹角余弦 cos(ij)= ( xik xjk)/( xik2 xjk2) 相关系数 cos(ij)=(xik-mi)(xjk-mj)/( (xik-mi)2 (xjk-mj)2) 其中mi和mj分别表示第i和第j个样本的均值.,PCA,计算特征矢量 t=x1p1+x2p2 (p12+p22)=1 p=(p1,p2)t 约束条件极大极小化，使用Lagrange乘子法, I为2x2的单位阵,PCA,得到特征值： ( 77.71-155.59 )p1 + 76.29p2 =0 76.29p1 + ( 80.86-155.59 )p2 =0 p1=0.6998 p2=0.7144 p1=0.7144 p2=-0.6998,Lamda: 1 2 155.5876 2.9838 Eigenvector: 0.6998 0.7144 0.7144 -0.6998,特征矢量1和2,特征值与特征矢量,特征值与特征矢量相关连，表示方差的程度 特征矢量是根据方差最大化原理进行原始数据的变量的线性组合 得到的特征矢量的方差比前一个特征矢量的更小，也就是依次递减 此时，特征矢量相互正交，也就是无相关性,PCA 的术语,特征矢量又可以叫做 载荷轴、主成分、潜变量、抽象因子 所谓得分就是原始数据在主成分轴上的投影，也就是下一个图中第一主成分上的圈,PC1和PC2的表示,PCA的优点,它能找到表现原始数据阵最重要的变量的组合 通过表示最大的方差，能有效的直观反映样本之间的关系 能从最初的几个主成分的得分来近似反映原始的数据阵,用PC1对X阵近似,# Bapprox Tapprox Bmc Tmc Score 1. 12.36 12.62 13 12 17.67 2. 7.41 7.56 9 6 10.58 3. 7.45 7.60 5 10 10.64 4. 3.47 3.54 3 4 4.96 5. -3.97 -4.05 -3 -5 -5.67 6. -7.43 -7.58 -7 -8 -10.61 7. -8.91 -9.09 -9 -9 -12.73 8. -10.39 -10.60 -11 -10 -14.84,用PC1对X阵近似,( 77.71-155.59 )p1 + 76.29p2 =0 76.29p1 + ( 80.86-155.59 )p2 =0 p1=0.6998 p2=0.7144 p1=0.7144 p2=-0.6998 计算得分的方法 例：13 x 0.6998+12 x 0.7144=17.67 利用PC1对X做近似 例：17.67 x 0.6998=12.36 X残差 例：12.36-13=-0.64,X 残差,# Bresid Tresid Bmc Tmc 1. -0.64 0.62 13 12 2. -1.59 1.56 9 6 3 2.45 -2.40 5 10 4. 0.47 -0.46 3 4 5. -0.97 0.95 -3 -5 6. -0.43 0.42 -7 -8 7. 0.09 0.09 -9 -9 8. 0.61 -0.60 -11 10,PC1对X阵的贡献,T (SPC1X0.7144=) 12.6236 7.5616 7.6034 3.5413 -4.0516 -7.5825 -9.0927 -10.6030,B (SPC1X0.6998=) 12.3656 7.4071 7.4480 3.4689 -3.9688 -7.4275 -8.9069 -10.3863,PC1 的得分Scores 17.6702 10.5846 10.6430 4.9570 -5.6714 -10.6138 -12.7278 -14.8418,PC2对X阵的贡献,T (SPC2*（-0.6998）=) -0.6225 -1.5611 2.3975 0.4591 -0.9488 -0.4182 0.0920 0.6021,B (SPC2*0.7144=) 0.6355 1.5937 -2.4475 -0.4686 0.9686 0.4269 -0.0939 -0.6147,PC2 的得分Scores 0.8896 2.2308 -3.4260 -0.6560 1.3558 0.5976 -0.1314 -0.8604,得分的平方和,(n-1) x Lamda1= PC1得分的平方和 (n-1) x Lamda2= PC2得分的平方和 Lamda反映的是相应主成分的方差大小,常用的 PCA,在该情况下，特征值就是得分的平方和 例： (17.67)2+(10.58)2+(10.64)2+(4.96)2+ (-5.67)2+(-10.61)2+(-12.73)2+(-14.84)2=1089 得分与载荷与上例相同,PCA,如果对X不进行预处理，则： 得分与载荷与前例不同,PCA,在实际计算中，PCA的计算常采用NIPALS (Nonlinear Iterative Partial Least Squares)方法 NIPALS方法并不是计算所有的因子，仅仅计算最初的 k个主成分,以上为8x2的矩阵，可变为2x2的矩阵,= 在实际应用中， 对于一个矩阵，Xmxn 其每一维的变量都大于2,即 m2,n2 怎么办？,Single value decomposition,主成分分析投影的数学意义,因为 X =USVt 即 XV=US 亦即 XV= X v1, v2,., vA = US 可见矩阵US=T (亦称非标准化的得分矩阵) 的每一个元素实际是每一个样本向量xit(i=1,2,.,n)对荷载矩阵V中的每一相互正交的荷载矢量上的投影坐标（内积本质上就是投影），它反映了样本与样本之间的相互关系; 同理可得，荷载矩阵的每一个元素实际是每一个变量向量xj(j=1,2,.,d)对得分矩阵中的每一相互正交的得分矢量上的投影坐标，它反映了变量与变量之间的相互关系。,主成分分析数学几何意义,=,投影,主成分分析的数学与几何意义示意图,2 Principal Component Regression (PCR),Lambert-Beer Law,Emxn,Ymxn,PCR 方法是采用多元统计中的成分分析方法，先对混合物量测矩阵 Y 矩阵直接进行分解，然后只取其中的主成分来进行回归分析，故有主成分回归之称。,方法模型,C= PY Y=Y0+E Y0 =U*S*Vt* Y0+=V*(S*)-1Ut* P=CY0+=C V*(S*)-1Ut* C未知=PY未知,C，Y分别是浓度矩阵和混合物测量矩阵 Y0表示只含混合物的量测值和一部分植入误差矩阵 Y0+是Y0 的广义逆 P 是回归系数矩阵 U，Vt 分别为标准列正交和行正交矩阵，即Scores 和 Loadings E为误差矩阵 U*，S*， Vt*分别为U，S，Vt的前n个特征值和特征矢量作为主成分，其余作为误差丢弃,PCR Advantages,Does not require wavelength selection. Any number can be used; usually the whole spectrum, or large regions. Larger number of wavelengths gives averaging effect, making model less susceptible to spectral noise. PCA data compression allows using inverse regression to calculate model coefficients; can calibrate only for constituents of interest. Can be used for very complex mixtures since only knowledge of constituents of interest is required. Can sometimes be used to predict samples with constituents (contaminants) not present in the original calibration mixtures.,PCR Disadvantages,Calculations are slower than most Classical methods. Optimization requires some knowledge of PCA; models are more complex to understand and interpret. No guarantee PCA vectors directly correspond to constituents of interest. Generally, a large number of samples are required for accurate calibration. Collecting calibration samples can be difficult; must avoid collinear constituent concentrations.,3 Partial Least Squares Regression (PLSR),PLSR 方法不仅直接对混合物量测矩阵 Y 进行正交分解，而且同时对浓度矩阵 C 进行正交分解。,方法模型,Y=USVt=U*S*Vt*+ Er =T*V*t+ Er C=PGQt=P*G*Qt*+Ec =R*Qt*+Ec,C，Y分别是浓度矩阵和混合物测量矩阵 P 是回归系数矩阵 U，Vt (P,Qt)分别为标准列正交和行正交矩阵，即Scores 和 Loadings Er，Ec为误差矩阵 U*，S*， Vt*分别为U，S，Vt的前n个特征值和特征矢量作为主成分，其余作为误差丢弃,NIPLS算法步骤,1. 初始迭代矢量： 取C某一列 r 2. r代替T*的t vt: vt=rtY/(rtr) 3. 归一化v： vt新= vt旧 / | vt旧| 4. vt t: t=Yv/(vtv) 5. t代替R*的r qt: qt=ttC/(ttt) 6. 归一化q： qt新= qt旧 / | qt旧| 7. Qt r: r =Cq/(qtq) 8.判断t是否收敛|t前一轮-t后一轮|是否小于给定常数。否，则重复2-7；是，则继续9。,9. 进行回归运算：b1=rt1/(tt1t1) 10. EY,1=Y-t1vt1 ; EC,1=C-r1qt1=C-b1t1qt1 11. EY,1用代替Y， EC,1代替C，返回主因子迭代过程9-10，算下一个主因子：r2, qt2, t2, vt2 得到EY,2 ，EC,2，直到得到所有的Y和C的主因子(可用交互校验法和方差比较来决定) 12. 由储存的Vt* Y未知(Y=T*Vt*) t未知 r未知 ( b=rt/(ttt) ) 由储存的Qt* C未知(C=R*Qt*),The geometric representation of PLSR. The X-matrix can be represented as N points in the K dimensional space where each column of X (x_k) defines one coordinate axis. The PLSR model defines an A-dimensional hyper-plane, which in turn, is defined by one line, one direction, per component. The direction coefficients of these lines are p_ak. The coordinates of each object, i, when its ak data (row i in X) are projected down on this plane are t_ia. These positions are related to the values of Y.,PLS Disadvantages,Calculations are slower that most Classical methods, especially PLS-1. Models are more abstract, thus more difficult to understand and interpret. Generally, a large number of samples are required for accurate calibration. Collecting calibration samples can be difficult; must avoid collinear constituent concentrations.,4 Pattern Recognition (PR),Feature extraction methods,Pretreatment methods for pattern recognition,变换法 xij,new= xij,old xij,new=lg(xij,old) 以达到改变数据标度的目的。 变量组合法 将原来的变量按一定方式，如变量相加、变量相减或变量相比等进行组合以产生新的变量。,Pattern recognition by latent projections,Principal component analysis(PCA) Partial least squares (PLS) SIMCA (classification based on principal component analysis),主成分分析的数学与几何意义示意图,SIMCA计算过程框图,从此框图我们可以看出，SIMCA实际上是在循环地使用主成分分析方法，它先是用主成分分析方法来得到整个样本的分类，然后，在此基础上，分别建立各类样本的主成分模型，继用它们来检验未知样本，判别未知样本的类别。由于整个SIMCA计算过程可在投影图上直接进行，所以，只要计算机程序的人机对话界面建立得好，此方法可完全在图上直接操作完成，使用起来非常方便。将SIMCA编入商业软件的例子很多，比较著名的有Unscrumber和Serius的软件。,令样本数据矩阵表示为： X= =USVt xjt(j=1,2,.,n)称为样本向量，则一个样本就由样本矩阵中的一个行向量来表示。 xjt=xct+k tikvkt+eit,主成分分析建模算法可由下述步骤所组成：,1) 定义样本到主成分模型的距离： si = eitei /(d-A) 这是因为样本矢量表达式的前一部分，即 tikvkt 正是主成分荷载矢量的线性组合，属于可由主成分模型表出的部分，在样本矢量表达式中没被主成分模型所包含的就是误差矢量，所以，样本到主成分模型的距离就是误差矢量的点积，上式中的分母是样本矢量在建立了主成分模型，即降维后的自由度，由原样本空间的维数减去主成分数而得。si在SIMCA方法中常被称为第i个样本的残余标准偏差 (RSD) 。,主成分分析建模算法可由下述步骤所组成：,2) 定义整个类的均值标准偏差 sc： 如令s= s1, s2, ., si,., snt，则整个类的均值标准偏差 sc 可由下式给出 sc= sts / (n-A-1) 上式以(n-A-1)为分母的理由是为了整个类的均值标准偏差不受该类内样本数多少的影响，而且，还考虑了由于引入A个主成分所失去的自由度。,主成分分析建模算法可由下述步骤所组成：,3) 为了判别某一样本与某一类主成分模型的差别或相似程度，SIMCA方法的提出者Wold引入了F-检验统计量。从上面的讨论可知，某一样本的残余标准偏差的平方si2和整个类的均值标准偏差的平方sc2分别具有自由度(d-A)和(d-A)(n-A-1)，据此，我们可以采用自由度为(d-A)和(d-A)(n-A-1) 的F-检验来判别某一样本与某一类主成分模型的差别或相似程度。F-检验或记为Fcrit可用来求得一个判定某样本是否属于某类的残余标准偏差平方的上界值： smax2=sc2 Fcrit 用于F-检验的置信水平一般取0.05或0.01即可。,置信水平取0.05或0.01的由一个主成分建立的类模型,由极端得分值表述的由一个主成分构成的类模型,对于任意未知样本xit，都可给出它的主成分模型: xit=xct+k tikvkt+eit 然后，计算它的误差项eit，继算出它到所定义的主成分模型的距离： si = eitei /(d-A) 将它平方值与smax2进行比较，即用前述F-检验即可判别它是否属于所定义的主成分模型。从此可以看出，只要求得了该未知样本矢量的误差项，判别就很容易了。所以，判别的过程，就只要是误差项的计算过程。,未知样本误差项的计算可由下列步骤组成：,1)令 ei(0)= xit-xct 在此xct为所求类的样本均值向量，为已知值。 2) 对a=1, 2, .,A个主成分进行下述循环计算，以求得该样本在所求类的得分tia(a=1,2,.,A)： tia = ei(a-1)va ei(a) = ei(a-1)- tiava 最终所得的ei(A)即是未知样本误差项。继用下式算出它与所求主成分模型的距离： si = eitei /(d-A) 注意到上述求未知样本误差项的过程始终都用到了主成分分析所得的荷载矢量的正交性质，即 vitvj=0 (ij)。,
展开阅读全文