前言
通过模型算法,熟练对Matlab和python的应用。
学习视频链接:
https://www.bilibili.com/video/BV1EK41187QF?p=16&vd_source=67471d3a1b4f517b7a7964093e62f7e6
一、主成分分析法
1.主成分分析法简介
在分析研究多变量的课题时,变量太多就会增加课题的复杂性。人们自然希望变量个数较少而得到的信息较多。在很多情形,变量之间是有一定的相关关系的,可以解释为这两个变量反映此课题的信息有一定的重叠。主成分分析是对于原先提出的所有变量,将重复的变量(关系紧密的变量)删去多余,建立尽可能少的新变量,使得这些新变量是两两不相关的,且这些新变量在反映课题的信息方面尽可能保持原有的信息。设法将原来变量重新组合成一组新的互相无关的几个综合变量,同时根据实际需要从中可以取出几个较少的综合变量尽可能多地反映原来变量的信息的统计方法叫做主成分分析或称主分量分析,也是数学上用来降维的一种方法。
2.主成分分析法原理
- PCA的主要目标是将特征维度变小,同时尽量减少信息损失。就是对一个样本矩阵,一是换特征,找一组新的特征来重新标识;二是减少特征,新特征的数目要远小于原特征的数目。
- 通过PCA将 n 维原始特征映射到 k 维(k<n)上,称这 k 维特征为主成分,需要强调的是,不是简单地从 n 维特征中去除其余(n-k )维特征,而是重新构造出全新的 k 维正交特征,且新生成的 k 维数据尽可能多地包含原来 n 维数据的信息。例如,使用PCA将20个相关的特征转化为5个无关的新特征,并且尽可能保留原始数据集的信息。
- PCA 是一种线性降维方法,即通过某个投影矩阵将高维空间中的原始样本点线性投影到低维空间,以达到降维的目的,线性投影就是通过矩阵变换的方式把数据映射到最合适的方向。
- 降维的几何意义可以理解为旋转坐标系,取前 k 个轴作为新特征。
- 降维的代数意义可以理解为 m × n 阶的原始样本 X,与 n × k 阶矩阵W做矩阵乘法运算 X × W ,即得到 m × k 阶低维矩阵Y,这里的 n × k 阶矩阵 W 就是投影矩阵。
3.主成分分析法思想
- 1、假设有 n 个样本,p 个指标,则可构成大小为 n × p 的样本矩阵 x :
x = [ x 11 x 12 ⋯ x 1 p x 21 x 22 ⋯ x 2 p ⋮ ⋮ ⋱ ⋮ x n 1 x n 2 ⋯ x n p ] = ( x 1 , x 2 , ⋯ , x p ) x=\begin{bmatrix}x_{11}&x_{12}&\cdots&x_{1p}\\x_{21}&x_{22}&\cdots&x_{2p}\\\vdots&\vdots&\ddots&\vdots\\x_{n1}&x_{n2}&\cdots&x_{np}\end{bmatrix}=\begin{pmatrix}x_1,x_2,\cdots,x_p\end{pmatrix} x= x11x21⋮xn1x12x22⋮xn2⋯⋯⋱⋯x1px2p⋮xnp =(x1,x2,⋯,xp) - 2、假设我们想找到新的一组变量
z
1
,
z
2
,
⋅
⋅
⋅
,
z
m
(
m
≤
p
)
z_1,z_2,\cdotp\cdotp\cdotp,z_m(m\leq p)
z1,z2,⋅⋅⋅,zm(m≤p) 且它们满足:
{ z 1 = l 11 x 1 + l 12 x 2 + ⋯ + l 1 p x p z 2 = l 21 x 1 + l 22 x 2 + ⋯ + l 2 p x p ⋮ z m = l m 1 x 1 + l n 2 x 2 + ⋯ + l n p x p \begin{cases}z_1=l_{11}x_1+l_{12}x_2+\cdots+l_{1p}x_p\\z_2=l_{21}x_1+l_{22}x_2+\cdots+l_{2p}x_p\\\vdots\\z_m=l_{m1}x_1+l_{n2}x_2+\cdots+l_{np}x_p\end{cases} ⎩ ⎨ ⎧z1=l11x1+l12x2+⋯+l1pxpz2=l21x1+l22x2+⋯+l2pxp⋮zm=lm1x1+ln2x2+⋯+lnpxp
系数 l i j l_{ij} lij 的确定原则:
(1) z i z_i zi与 z j ( i ≠ j z_j( i\neq j zj(i=j ; i , j = 1 , 2 ; i, j= 1, 2 ;i,j=1,2 , m , m ,m)相互无关
(2) z 1 z_1 z1是与 x 1 , x 2 , . . . , x p x_1,x_2,...,x_p x1,x2,...,xp 的一切线性组合中方差最大者
(3) z 2 与 z 1 z_2与z_1 z2与z1不相关的 x 1 , x 2 , . . . , x p x_1,x_2,...,x_p x1,x2,...,xp 的一切线性组合中方差最大者
(4)以此类推, z m z_m zm是与 z 1 , z 2 , . . . , z m − 1 z_1,z_2,...,z_{m-1} z1,z2,...,zm−1不相关的 x 1 , x 2 , . . . , x p x_1,x_2,...,x_p x1,x2,...,xp 的所有线性组合中方差最大者
(5)新变量指标 z 1 , z 2 , . . . , z m z_1,z_2,...,z_m z1,z2,...,zm分别称为原变量指标 x 1 , x 2 , . . . , x p x_1,x_2,...,x_p x1,x2,...,xp的第一,第二,…,第 m m m主成分
4.PCA的计算步骤
-
1、假设有 n 个样本,p 个指标,则可构成大小为 n × p 的样本矩阵 x :
x = [ x 11 x 12 ⋯ x 1 p x 21 x 22 ⋯ x 2 p ⋮ ⋮ ⋱ ⋮ x n 1 x n 2 ⋯ x n p ] = ( x 1 , x 2 , ⋯ , x p ) x=\begin{bmatrix}x_{11}&x_{12}&\cdots&x_{1p}\\x_{21}&x_{22}&\cdots&x_{2p}\\\vdots&\vdots&\ddots&\vdots\\x_{n1}&x_{n2}&\cdots&x_{np}\end{bmatrix}=\begin{pmatrix}x_1,x_2,\cdots,x_p\end{pmatrix} x= x11x21⋮xn1x12x22⋮xn2⋯⋯⋱⋯x1px2p⋮xnp =(x1,x2,⋯,xp)
(1)首先对其进行标准化处理:
按列计算均值 x ˉ j = 1 n ∑ i = 1 n x i j \bar{x}_{j}=\frac{1}{n}\sum_{i=1}^{n}x_{ij} xˉj=n1∑i=1nxij 和标准差 S j = ∑ i = 1 n ( x i j − x j ‾ ) 2 n − 1 S_{j}=\sqrt{\frac{\sum_{i=1}^{n}(x_{ij}-\overline{x_{j}})^{2}}{n-1}} Sj=n−1∑i=1n(xij−xj)2 ,计算得标准化数据 X i j = x i j − x j ‾ s j X_{ij}=\frac{x_{ij}-\overline{x_{j}}}{s_{j}} Xij=sjxij−xj,原始样本矩阵经过标准化变为: X = [ X 11 X 12 ⋯ X 1 , p X 21 X 22 ⋯ X 2 , p ⋮ ⋮ ⋱ ⋮ X n 1 X n 2 ⋯ X n p ] = ( X 1 , X 2 , ⋯ , X p ) X=\begin{bmatrix}X_{11}&X_{12}&\cdots&X_{1,p}\\X_{21}&X_{22}&\cdots&X_{2,p}\\\vdots&\vdots&\ddots&\vdots\\X_{n1}&X_{n2}&\cdots&X_{np}\end{bmatrix}=\left(X_1,X_2,\cdots,X_p\right) X= X11X21⋮Xn1X12X22⋮Xn2⋯⋯⋱⋯X1,pX2,p⋮Xnp =(X1,X2,⋯,Xp)
(2)计算标准化样本的协方差矩阵/样本相关系数矩阵
R = [ r 11 r 12 ⋯ r 1 p r 21 r 22 ⋯ r 2 p ⋮ ⋮ ⋱ ⋮ r n 1 r n 2 ⋯ r n p ] {R}=\begin{bmatrix}r_{11}&r_{12}&\cdots&r_{1p}\\r_{21}&r_{22}&\cdots&r_{2p}\\\vdots&\vdots&\ddots&\vdots\\r_{n1}&r_{n2}&\cdots&r_{np}\end{bmatrix} R= r11r21⋮rn1r12r22⋮rn2⋯⋯⋱⋯r1pr2p⋮rnp
其中 r i j = 1 n − 1 ∑ k = 1 n ( X k i − X i ‾ ) ( X k j − X j ‾ ) {r_{ij}}=\frac{1}{n-1}\sum_{k=1}^{n}\Big(X_{ki}-\overline{X_{i}}\Big)\Big(X_{kj}-\overline{X_{j}}\Big) rij=n−11∑k=1n(Xki−Xi)(Xkj−Xj),即 R = 1 n − 1 X T X R=\frac1{n-1}X^TX R=n−11XTX
(3)计算 R 的特征值和特征向量- 特征值: λ 1 ≥ λ 2 ≥ . . . ≥ λ p ≥ 0 \lambda_{1}\geq\lambda_{2}\geq...\geq\lambda_{p}\geq0 λ1≥λ2≥...≥λp≥0
- 特征向量: a 1 = [ a 11 a 21 ⋮ a p 1 ] , a 2 = [ a 12 a 22 ⋮ a p 2 ] , ⋯ , a p = [ a 1 p a 2 p ⋮ a p p ] a_1=\begin{bmatrix}a_{11}\\a_{21}\\\vdots\\a_{p1}\end{bmatrix},a_2=\begin{bmatrix}a_{12}\\a_{22}\\\vdots\\a_{p2}\end{bmatrix},\cdots, a_p=\begin{bmatrix}a_{1p}\\a_{2p}\\\vdots\\a_{pp}\end{bmatrix} a1= a11a21⋮ap1 ,a2= a12a22⋮ap2 ,⋯,ap= a1pa2p⋮app
(4)计算主成分贡献率以及累计贡献率
- 贡献率: α i = λ i ∑ k = 1 p λ k ( i = 1 , 2 , … , p ) \alpha_{i}=\frac{\lambda_{i}}{\sum_{k=1}^{p}\lambda_{k}}\left(i=1 ,2 ,\dots ,p\right) αi=∑k=1pλkλi(i=1,2,…,p)
- 累计贡献率: ∑ G = ∑ k − 1 i λ k ∑ k = 1 p λ k ( i = 1 , 2 , … , p ) \sum G{=}\frac{\sum_{k-1}^{i}\lambda_{k}}{\sum_{k=1}^{p}\lambda_{k}}\left(i=1 ,2 ,\ldots ,p\right) ∑G=∑k=1pλk∑k−1iλk(i=1,2,…,p)
(5)写出主成分
一般取累计贡献率超过80%的特征值所对应的第一、第二、…、第 m(m≤p) 个主成分。
第 i 个主成分: F i = a 1 i X 1 + a 2 i X 2 + . . . + a p i X p ( i = 1 , 2 , . . . , m ) F_{i}=a_{1i}X_{1}+a_{2i}X_{2}+...+a_{pi}X_{p}\left(i=1 ,2 ,... ,m\right) Fi=a1iX1+a2iX2+...+apiXp(i=1,2,...,m)(6)根据系数分析主成分代表的意义
对于某个主成分而言,指标前面的系数越大,代表该指标对于该主成分的影响越大。(7)利用主成分的结果进行后续的分析
- 主成分得分
- 聚类分析
- 回归分析
-
探究棉花单产和五个指标之间的关系:
将五个指标:种子费、化肥费、农药费、机械费、灌溉费作为指标矩阵 x ,按照主成分分析的步骤对指标矩阵进行降维。
二、代码实现----Matlab
clear;clc;
load data1.mat
[n,m] = size(x);
% 标准化处理
X = zscore(x);
% 计算协方差矩阵
R = cov(X);
% 计算R的特征值和特征向量
[V,D] = eig(R);
eigenvalues = diag(D);
[eigenvalues_sorted, sort_index] = sort(eigenvalues, 'descend');
V_sorted = V(:, sort_index);
D_sorted = eigenvalues_sorted;
% 计算主成分贡献率和累计贡献率
alpha = D_sorted ./ sum(D_sorted);
alpha_cumsum = cumsum(D_sorted) / sum(D_sorted); % cumsum()是求累加值的函数
三、代码实现----python
import numpy as np
import pandas as pd
df = pd.read_excel("棉花产量论文作业的数据.xlsx","Sheet1",engine="openpyxl")
x = df.values[:,2:] # 提取表格数据构造矩阵
n, m = x.shape
# 标准化矩阵
x_mean = np.mean(x,0)
S = np.std(x,0)
X = (x - np.tile(x_mean,(n,1))) / np.tile(S,(n,1))
# 计算协方差
R = np.cov(X,rowvar=False)
# 计算R的特征值,特征向量
D,V = np.linalg.eig(R)
sort_index = np.argsort(D)
descending_indices = sort_index[::-1]
sort_D = D[descending_indices]
sort_V = V[:,descending_indices]
# 计算主成分贡献率
alpha = sort_D / sum(sort_D)
# 累计贡献率
add_alpha = np.cumsum(sort_D) / sum(sort_D)
总结
本文介绍了主成分分析法对指标矩阵降维,并分别使用Matlab和python进行代码编写。