PCA主成分分析去噪与降维

最新推荐文章于 2024-07-20 17:25:24 发布

NoTime4Emotion

最新推荐文章于 2024-07-20 17:25:24 发布

阅读量1.9w

点赞数 9

分类专栏：机器学习文章标签： PCA 主成分分析降维降噪

本文链接：https://blog.csdn.net/qq_42192910/article/details/85881257

版权

机器学习专栏收录该内容

2 篇文章 2 订阅

订阅专栏

PCA主成分分析去噪与降维

（项目学习笔记及代码）

一、目的

1、降低维数（主要目的）

主成分分析 (Principal Component Analysis，简称 PCA)是最常用的一种降维方法。（具体可见西瓜书）

2、去噪

PCA也可以用作噪音的过滤，因为任何一个成分的变化影响都远远大于随机噪声的影响的，所以针对于噪声，各种的成分相对不受影响。即可以用主成分来重构原始数据。

3、原理

假设 $\mathbb{R}^{d}$ 空间中有 $m$ 个样本点 $X=\{x_1,x_2,...,x_m\}$ 。这些点即可以理解为我们数据的特征，如图片的像素，每一个数据 $x_i$ 都有 $d$ 维， $x_i \in \mathbb{R}^{d}$ 。
$\begin{aligned} x_i={\left[ \begin{matrix} x_i^1 \\ x_i^2 \\ \cdot \\ \cdot \\ \cdot \\ x_i^d \\ \end{matrix} \right] }_{d \times 1} \tag{1} \end{aligned}$
于是有：
$\begin{aligned} X=(x_1,x_2,...,x_m)_{d \times m}^\mathrm{T}={\left[ \begin{matrix} x_1^1 & x_2^1 & \cdot & \cdot & \cdot & x_1^d \\ x_1^2 & x_2^2 & \cdot & \cdot & \cdot & x_2^d \\ \cdot & \cdot & \cdot & & & \cdot \\ \cdot & \cdot & & \cdot & & \cdot \\ \cdot & \cdot & & & \cdot & \cdot \\ x_m^1 & x_m^2 & \cdot & \cdot & \cdot & x_m^d \\ \end{matrix} \right] }_{m \times d} \tag{2} \end{aligned}$
其中 $x_i^j$ 表示是第 $i$ 个样本的第 $j$ 维的特征，单独的列向量 $x_i$ 即是 $d$ 维的。每一行代表一个数据样本。
首先对 $x_i$ 数据进行归一化（保证其均值为0）。
$\bar{x}= \frac{1}{m} \sum_{i=1}^m x_i={\left[ \begin{matrix} \bar{x}^1 \\ \bar{x}^2 \\ \cdot \\ \cdot \\ \cdot \\ \bar{x}^d \\ \end{matrix} \right] }_{d \times 1} \tag{3}$
对所有的数据样本点列向量 $x_i$ 求对应元素的均值，得到d维的均值列向量 $\bar{x}$ ，然后用原始样本数据减去均值即得到了归一化后的结果。

$\begin{aligned} x_i & =x_i-\bar{x}={\left[ \begin{matrix} x_i^1 - \bar{x}^1 \\ x_i^2 - \bar{x}^2 \\ \cdot \\ \cdot \\ \cdot \\ x_i^d - \bar{x}^d \\ \end{matrix} \right] }_{d \times 1} \tag{5} \end{aligned}$

对数据样本 $x_i\in \mathbb{R}^{d}$ 进行降维压缩得到结果 $z_i\in \mathbb{R}^{d'}$ ， $(d^{'} < d)$ 。由 $x_i$ 得到 $z_i$ 的过程可以如下表示，权重 $\mathbf{W}=\{w_1,w_2,...,w_{d'} \}$ 作用下的矩阵线性变换。这就是我们降维后的结果，也是PCA工作的目标：
$z_i^j=w_j^\mathrm{T} \cdot x_i={\left[ \begin{matrix} w_j^1 \\ w_j^2 \\ \cdot \\ \cdot \\ \cdot \\ w_j^d \\ \end{matrix} \right] }_{d \times 1}^\mathrm{T} \cdot {\left[ \begin{matrix} {x}^1 \\ {x}^2 \\ \cdot \\ \cdot \\ \cdot \\ {x}^d \\ \end{matrix} \right] }_{d \times 1} \tag{6}$
其中 $z_i^j$ 是降低维度后的原始样本数据 $x_i$ 对应的第 $j$ 维的元素值；而 $w_j$ 是 $d$ 维的列向量。
从最大可分（最大投影方差）的角度出发，同时当 $w_j \|_2^2=1$ 时， $w_j^\mathrm{T} \cdot x_i$ 即是样本数据点在该矩阵方向上得投影。因为 $\alpha$ 到 $\beta$ 的投影即是 $|\alpha|\cos \theta$ ，而向量点乘有 $\alpha \cdot \beta=|\alpha||\beta|\cos \theta$ ，当 $|\beta|=1$ 时，点乘就是投影了。
我们的目标就是找到恰当的 $\mathbf{W}$ 将这个投影 $w_j^\mathrm{T} \cdot x_i$ 的方差最大化，也就是让样本点在新的基准下划分得越分开越好。 于是有样本方差来度量：
（均值已经归一化为0了）
$\begin{aligned} S &= \frac{1}{m-1} \sum_{i=1}^m (w_j^\mathrm{T}x_i-w_j^\mathrm{T}\bar{x})^2\\ &= \frac{1}{m-1} \sum_{i=1}^m (w_j^\mathrm{T}(x_i-\bar{x}))^2\\ &= \frac{1}{m-1} \sum_{i=1}^m (w_j^\mathrm{T}(x_i-\bar{x})) \cdot (w_j^\mathrm{T}(x_i-\bar{x}))^\mathrm{T} \\ &= \frac{1}{m-1} \sum_{i=1}^m w_j^\mathrm{T}(x_i-\bar{x}) \cdot (x_i-\bar{x})^\mathrm{T}w_j \\ & \rightarrow w_j^\mathrm{T} \left(\sum_{i=1}^m (x_i-\bar{x}) \cdot (x_i-\bar{x})^\mathrm{T} \right)w_j \\ & \rightarrow w_j^\mathrm{T} \left(\text{Cov} (x) \right)w_j \\ & \rightarrow w_j^\mathrm{T}S w_j \\ \tag{7} \end{aligned}$
所以投影的方差最后转换成了样本矩阵的协方差。
于是得到优化目标：
$\argmax_{w_j} w_j^\mathrm{T}S w_j \\ \text{s.t. } w_j \cdot w_j^\mathrm{T}=1 \tag{8}$
利用拉格朗日数乘法：
$L(w_j,\lambda)=w_j^\mathrm{T}S w_j + \lambda (w_j \cdot w_j^\mathrm{T}-1) \tag{9}$
再对每一个 $w_j$ 求偏微分 $\frac{\partial L}{\partial w_j}=0$ ：
$Sw_j-\lambda w_j=0\\ Sw_j=\lambda w_j \\ \tag{10}$
于是乎，形如此形式有： $w_j$ 即是特征向量，而 $\lambda$ 是特征值。
再带回优化目标有：
$w_j^\mathrm{T}S w_j=w_j^\mathrm{T}\lambda w_j = \lambda \tag{11}$
这意味这最大化的结果值（方差最大化的结果值）就是数据样本协方差矩阵的特征值。 可以理解到特征值在某种意义上表达的是方差的意义。
依次类推，此后得到的每一个的 $w_j$ 还需要保证相互之间点乘为0是正交的。
于是一般即可以理解为找到数据样本 $X$ 的协方差矩阵 $S$ ，是一个实对称矩阵，那么其一定是可以相似对角化的，找到其最大的特征值 $\lambda$ 即是最大的投影方差，其对应的特征向量 $w_j$ 即是线性乘法因子。而最后选择最大的特征值即可。一般来说目前都直接奇异值分解法求解，因为奇异值分解数据矩阵 $X$ 的过程也是会构造出 $XX^\mathrm{T}$ ，也就是协方差矩阵。

二、代码

1、本文代码是通过加载原始的matlab 中的一个数组文件作为数据集的（数据集可以是任意的数据，但注意最后用plt库绘图时，注意下标的匹配）
2、其中有些地方代码已经注释的很清楚了。
3、若进行降噪即最后重构（restore）即可，若进行降维，那么经过转换后（transform）的即是新的维度下的结果。
4、代码引用了部分搜索到的，不过有些确实找不到了而没有引用！

1、sklearn库方式

# -*- coding: utf-8 -*-
"""
Created on Tue Jan  8 10:39:09 2019

@author: QSY
"""


import matplotlib.pyplot as plt
import scipy.io as sci
from sklearn.decomposition import PCA
import numpy as np

"""Read data from  *.mat including two kinds"""
data_tarin = sci.loadmat('data_train.mat')
# data_test = sci.loadmat('data_test.mat')
data_train_noisy = data_tarin['noiseSignal']
data_train_theory = data_tarin['clearSignal']

#选取数据进行操作
X=data_train_noisy[0:1,:]
#进行行列转换
X=list(map(list, zip(*X)))

#创建横坐标
p=[]
for i in range(434): 
    j=i
    p.append(j)


def ReadReductionofcomponents(X):
    
    pca=PCA(n_components=1)
    pca.fit(X)
    #进行操作转换后的结果
    x_reduction =pca.transform(X)
    
    #进行复原恢复后的结果
    x_restore=pca.inverse_transform(x_reduction)   
    x_restore=x_restore[:,0:1]  
  
    #X仅选择某一列的数据进行操作
    X=np.array(X)
    X=X[:,0:1]
    
    #进行线绘制   
    plt.plot(p, x_restore, color='r', linewidth=0.3, alpha=0.6,label="PCA")
    plt.plot(p, X, color='b', linewidth=0.2, alpha=0.6,label="init")
    plt.legend()#显示图例名字
    #pdf保存
    plt.savefig('1.png')
    plt.show()   
    
ReadReductionofcomponents(X)

2、数学方式

# -*- coding: utf-8 -*-
"""
Created on Tue Jan  8 16:24:57 2019

@author: QSY
"""


import numpy as np
import scipy.io as sci
import matplotlib.pyplot as plt


"""Read data from  *.mat including two kinds"""
data_tarin = sci.loadmat('data_train.mat')
other_data=sci.loadmat('pca.mat')
# data_test = sci.loadmat('data_test.mat')
data_train_noisy = data_tarin['noiseSignal']
data_train_theory = data_tarin['clearSignal']
perfect_data=other_data['U']
#选取数据进行操作
X=data_train_noisy[0:3,:]
X=list(map(list, zip(*X)))

X=perfect_data
#零均值化
def zeroMean(dataMat):      
    #按列求均值，即求各个特征的均值
    meanVal=np.mean(dataMat,axis=0)     
    newData=dataMat-meanVal
    return newData,meanVal
 
def pca(dataMat,n):
    
    newData,meanVal=zeroMean(dataMat)
    covMat=np.cov(newData,rowvar=0)    #求协方差矩阵
    
    eigVals,eigVects=np.linalg.eig(np.mat(covMat))#求特征值和特征向量,特征向量是按列放的，即一列代表一个特征向量
    eigValIndice=np.argsort(eigVals)            #对特征值从小到大排序
    n_eigValIndice=eigValIndice[-1:-(n+1):-1]   #最大的n个特征值的下标
    n_eigVect=eigVects[:,n_eigValIndice]        #最大的n个特征值对应的特征向量
    lowDDataMat=newData*n_eigVect               #低维特征空间的数据
    reconMat=(lowDDataMat*n_eigVect.T)+meanVal  #重构数据
    return lowDDataMat,reconMat

#newData,meanVal=zeroMean(X)
lowDDataMat,reconMat=pca(X,n=1)

p=[]
for i in range(24): 
    j=i#np.log10(i)
    p.append(j)

plt.plot(p, lowDDataMat, color='r', linewidth=0.2, alpha=0.6,label="after PCA")
#plt.plot(p, X, color='b', linewidth=0.2, alpha=0.6,label="init")
plt.legend()#显示图例名字
#pdf保存
plt.savefig('1.pdf')
plt.show()