PCA降维算法

目前图像特征的提取主要有两种方法：传统图像特征提取方法 和 深度学习方法。

1. 传统的特征提取方法：基于图像本身的特征进行提取（PCA）；
2. 深度学习方法：基于样本自动训练出区分图像的特征分类器；

特征选择（feature selection）和特征提取（Feature extraction）都属于降维（Dimension reduction）

特征提取：

特征提取的主要方法：主要目的是为了排除信息量小的特征，减少计算量等：

主成分分析（PCA）；

PCA主成分分析是传统特征提取方法当中最常用的算法，没有之一，传统的特征提取方法没有深度学习那么严谨，智能，准确度也没那么高，但它在做法上没有深度学习那么复杂，因为深度学习还需要训练参数等等，比较耗时。很多情况下并不需要用到深度学习的知识，简单点说，就是杀鸡焉用牛刀。很多人知道PCA是用来降维的，但是却不了解PCA算法的原理，于是深入了解PCA算法原理就成了一个分水岭。基于此，本文主要介绍PCA算法的原理。

PCA算法是如何实现的

简单来说，就是将数据从原始的空间中转换到新的特征空间中，例如原始的空间是三维的(x,y,z)，x、y、z分别是原始空间的三个基，我们可以通过某种方法，用新的坐标系(a,b,c)来表示原始的数据，那么a、b、c就是新的基，它们组成新的特征空间。在新的特征空间中，可能所有的数据在c上的投影都接近于0，即可以忽略，那么我们就可以直接用(a,b)来表示数据，这样数据就从三维的(x,y,z)降到了二维的(a,b)。

问题的关键在于如何求新的基（a，b，c）？

步骤：

1. 对原始数据零均值化（中心化），
2. 求协方差矩阵，
3. 对协方差矩阵求特征向量和特征值，这些特征向量组成了新的特征空间。

1. PCA——零均值化（中心化）：

把所有的变量都去减掉他们的均值，使得减完之后得到的数据集均值为0，其实本质就是一个平移的过程，平移后得到的所有数据的中心是（0，0）

公式推导
设均值为E
(1/n) * (x1 + x2 + … +xn) = E
所有的变量都减去均值
(1/n) * (x1 - E + x2 - E + … + xn - E)
= (1/n) * ((x1 + x2 + … +xn) -nE)
= (1/n) * (x1 + x2 + … +xn) - E
= 0

只有“中心化”数据之后，计算得到的方向才能比较好的“概括”原来的数据。
下图描述了“中心化”的几何意义，即将样本集的中心平移到坐标系的原点0上

对于上左图，蓝色的圈代表数据域，凭肉眼观察，按理来说数据的发展方向是从左上到右下，正好是跟它现在的均值是垂直的。而我们一般而言一组数据的均值是能够代表这组数据的趋势的，但显然现在的均值并不能代表该数据的发展趋势。所以当我们做了零均值化过后，就和我们的预期达到一致了。所以我们做零均值化的原因就在于优化，使得我们在计算的时候更加准确，如果不做这一步，有可能会导致误差变大。

2. PCA——PCA降维的几何意义

对于一组数据，如果它在某一坐标轴上的方差越大，说明坐标点越分散，该属性就能够比较好地反映源数据。

这里举个例子便于理解，假设现在有3个人的人脸特征，他们的特征都有三个维度（眼睛，鼻子，嘴巴），假设对于鼻子这一维度的数据3个人都是一样的，都是大鼻子，那么这个时候鼻子这个维度的方差就为0，这个维度的离散程度就很小。那么如果我们要做分类任务，我们从鼻子这个维度所获得的信息量就基本为0，也就是说我们无法从这个维度得知一张人脸属于谁。所以我要想通过鼻子这一维度来区分出这3个人的话，我需要一个大的方差，比如说他们的鼻子分别是大，中，小鼻子，此时的方差就比较大了，这时就可以通过鼻子这一维度来获得比较大的信息量，通过这一维度来进行人脸识别。这样一来鼻子这一维度就是一个有效的特征。基于此。引出PCA降维算法的目标：

1. 降维后同一维度的方差最大
2. 不同维度之间的相关性为0

因为我们的目的在于用尽量少维度来表示物体特征，那么我一定是希望每个维度的特征越明显越好，所以要求方差最大。并且我不能有a和b能够推出c的情况，因为这样的话就形成了冗余特征，所以不同维度之间不能有相关性

由于不同维度之间相关性应为0，由此便引出协方差概念

3. PCA——协方差

协方差是用来度量两个随机变量关系的统计量
同一元素的协方差就表示该元素的方差，不同元素之间的协方差就表示它们的相关性
协方差>0 , 表示X和Y正相关
协方差=0 , 表示X和Y不相关（我们要找的就是这个点）
协方差<0 , 表示X和Y负相关

在理解协方差时应将其与方差结合起来，因为方差描述的是变量和均值之间的关系，那对于协方差而言，如果两个变量之中有一个是均值，那此时协方差就变成了方差，所以方差实际上是协方差的一个特例，当X和Y为同一特征时，此时的协方差就是方差。所以说同一元素的协方差就表示该元素的方差，不同元素之间的协方差就表示它们的相关性。所以说我们最后得到的几个特征，两两之间的协方差都应该为0。
这里有一个问题在于协方差只能处理二维的数据，但是我们的整个数据集当中，特征维度肯定不止2维，可能有100维，那我们两两之间都要做协方差，所以这时就引入了矩阵来组织多维数据，所以我们就用到了协方差矩阵。

4. PCA——协方差矩阵

定义：

这里其实就是把两两之间的协方差都列举出来了
协方差矩阵的特点：

• 协方差矩阵计算的是不同维度之间的协方差，而不是不同样本的
• 样本矩阵的每行是一个样本，每列为一个维度，所以我们要按列计算均值
• 协方差矩阵的对角线就是各个维度上的方差

这里需要区分样本矩阵和协方差矩阵

特别的，如果做了中心化，则协方差矩阵为（中心化矩阵的协方差矩阵公式，m为样本个数）：

5. PCA——对协方差矩阵求特征值，特征矩阵

A为n阶矩阵（在应用中，A为协方差矩阵），若数λ和n维非0列向量x满足Ax=λx，那么数λ称为A的特征值，x称为A的对应于特征值λ 的特征向量。
式Ax=λx也可写成( A-λE)x=0，E是单位矩阵，并且|A-λE|叫做A 的特征多项式。当特征多项式等于0的 时候，称为A的特征方程，特征方程是一个齐次线性方程组，求解特征值的过程其实就是求解特征方程的解。
对于协方差矩阵A，其特征值（ 可能有多个）计算方法为：

对数字图像矩阵做特征值分解，其实是在提取这个图像中的特征，这些提取出来的特征是一个个的向量，即 对应着特征向量。而这些特征在图像中到底有多重要，这个重要性则通过特征值来表示。
比如一个100x100的图像矩阵A分解之后，会得到一个100x100的特征向量组成的矩阵Q，以及一个100x100 的只有对角线上的元素不为0的矩阵E，这个矩阵E对角线上的元素就是特征值，而且还是按照从大到小排列 的（取模，对于单个数来说，其实就是取绝对值），也就是说这个图像A提取出来了100个特征，这100个特 征的重要性由100个数字来表示，这100个数字存放在对角矩阵E中。
所以归根结底，特征向量其实反应的是矩阵A本身固有的一些特征，而此时A表示的是协方差矩阵，所以此时的特征向量反映的是协方差矩阵的特征。本来一个矩阵就是一个线性变换， 当把这个矩阵作用于一个向量的时候，通常情况绝大部分向量都会被这个矩阵A变换得“面目全非”， 但是偏偏刚好存在这么一些向量，被矩阵A变换之后居然还能保持原来的样子，于是这些向量就可以 作为矩阵的核心代表了。 于是我们可以说：一个变换（即一个矩阵）可以由其特征值和特征向量完全表述，这是因为从数学上 看，这个矩阵所有的特征向量组成了这个向量空间的一组基底。而矩阵作为变换的本质其实就是把一个基底下的东西变换到另一个基底表示的空间中。

6. PCA——对特征值进行排序

7. PCA——评价模型的好坏，k值的确定

通过特征值的计算我们可以得到主成分所占的百分比，用来衡量模型的好坏。
对于前k个特征值所保留下的信息量计算方法如下：

8. PCA——代码实现

1. 手动实现

# -*- coding: utf-8 -*-
"""
使用PCA求样本矩阵X的K阶降维矩阵Z
"""import numpy as npclass CPCA(object):'''用PCA求样本矩阵X的K阶降维矩阵ZNote:请保证输入的样本矩阵X shape=(m, n)，m行样例，n个特征'''def __init__(self, X, K):''':param X,训练样本矩阵X:param K,X的降维矩阵的阶数，即X要特征降维成k阶'''self.X = X       #样本矩阵Xself.K = K       #K阶降维矩阵的K值self.centrX = [] #矩阵X的中心化self.C = []      #样本集的协方差矩阵Cself.U = []      #样本矩阵X的降维转换矩阵self.Z = []      #样本矩阵X的降维矩阵Zself.centrX = self._centralized()self.C = self._cov()self.U = self._U()self.Z = self._Z() #Z=XU求得def _centralized(self):'''矩阵X的中心化'''print('样本矩阵X:\n', self.X)centrX = []mean = np.array([np.mean(attr) for attr in self.X.T]) #样本集的特征均值print('样本集的特征均值:\n',mean)centrX = self.X - mean ##样本集的中心化print('样本矩阵X的中心化centrX:\n', centrX)return centrXdef _cov(self):'''求样本矩阵X的协方差矩阵C'''#样本集的样例总数ns = np.shape(self.centrX)[0]#样本矩阵的协方差矩阵CC = np.dot(self.centrX.T, self.centrX)/(ns - 1)print('样本矩阵X的协方差矩阵C:\n', C)return Cdef _U(self):'''求X的降维转换矩阵U, shape=(n,k), n是X的特征维度总数，k是降维矩阵的特征维度'''#先求X的协方差矩阵C的特征值和特征向量a,b = np.linalg.eig(self.C) #特征值赋值给a，对应特征向量赋值给b。函数doc：https://docs.scipy.org/doc/numpy-1.10.0/reference/generated/numpy.linalg.eig.html print('样本集的协方差矩阵C的特征值:\n', a)print('样本集的协方差矩阵C的特征向量:\n', b)#给出特征值降序的topK的索引序列ind = np.argsort(-1*a)#构建K阶降维的降维转换矩阵UUT = [b[:,ind[i]] for i in range(self.K)]U = np.transpose(UT)print('%d阶降维转换矩阵U:\n'%self.K, U)return Udef _Z(self):'''按照Z=XU求降维矩阵Z, shape=(m,k), n是样本总数，k是降维矩阵中特征维度总数'''Z = np.dot(self.X, self.U)print('X shape:', np.shape(self.X))print('U shape:', np.shape(self.U))print('Z shape:', np.shape(Z))print('样本矩阵X的降维矩阵Z:\n', Z)return Zif __name__=='__main__':'10样本3特征的样本集, 行为样例，列为特征维度'X = np.array([[10, 15, 29],[15, 46, 13],[23, 21, 30],[11, 9,  35],[42, 45, 11],[9,  48, 5],[11, 21, 14],[8,  5,  15],[11, 12, 21],[21, 20, 25]])K = np.shape(X)[1] - 1print('样本集(10行3列，10个样例，每个样例3个特征):\n', X)pca = CPCA(X,K)

2. 利用numpy库封装好的函数简化代码

#coding=utf-8import numpy as np
class PCA():def __init__(self,n_components):self.n_components = n_componentsdef fit_transform(self,X):self.n_features_ = X.shape[1]# 求协方差矩阵X = X - X.mean(axis=0)self.covariance = np.dot(X.T,X)/X.shape[0]# 求协方差矩阵的特征值和特征向量eig_vals,eig_vectors = np.linalg.eig(self.covariance)# 获得降序排列特征值的序号idx = np.argsort(-eig_vals)# 降维矩阵self.components_ = eig_vectors[:,idx[:self.n_components]]# 对X进行降维return np.dot(X,self.components_)# 调用
pca = PCA(n_components=2)
X = np.array([[-1,2,66,-1], [-2,6,58,-1], [-3,8,45,-2], [1,9,36,1], [2,10,62,1], [3,5,83,2]])  #导入数据，维度为4
newX=pca.fit_transform(X)
print(newX)                  #输出降维后的数据

3. 利用sklearn库简化

sklearn是一个第三方库，里面提供了一个PCA函数，调用它可以直接进行降维

#coding=utf-8import numpy as np
from sklearn.decomposition import PCA
X = np.array([[-1,2,66,-1], [-2,6,58,-1], [-3,8,45,-2], [1,9,36,1], [2,10,62,1], [3,5,83,2]])  #导入数据，维度为4
pca = PCA(n_components=2)   #降到2维
pca.fit(X)                  #训练
newX=pca.fit_transform(X)   #降维后的数据
# PCA(copy=True, n_components=2, whiten=False)
print(pca.explained_variance_ratio_)  #输出贡献率
print(newX)                  #输出降维后的数据

9. PCA——鸢尾花实例

我们通过Python的sklearn库来实现鸢尾花数据 进行降维，数据本身是4维的，降维后变成2维。
其中样本总数为150，鸢尾花的类别有三种。

#!/usr/bin/env python
# encoding=gbkimport matplotlib.pyplot as plt
import sklearn.decomposition as dp
from sklearn.datasets.base import load_irisx,y=load_iris(return_X_y=True) # 加载数据，x表示数据集中的属性数据，y表示数据标签
pca=dp.PCA(n_components=2) # 加载PCA算法，设置降维后主成分数目为2
reduced_x=pca.fit_transform(x) # 对原始数据进行降维，保存在reduced_x中
red_x,red_y=[],[]
blue_x,blue_y=[],[]
green_x,green_y=[],[]
for i in range(len(reduced_x)): # 按鸢尾花的类别将降维后的数据点保存在不同的表中if y[i]==0:red_x.append(reduced_x[i][0])red_y.append(reduced_x[i][1])elif y[i]==1:blue_x.append(reduced_x[i][0])blue_y.append(reduced_x[i][1])else:green_x.append(reduced_x[i][0])green_y.append(reduced_x[i][1])
plt.scatter(red_x,red_y,c='r',marker='x')
plt.scatter(blue_x,blue_y,c='b',marker='D')
plt.scatter(green_x,green_y,c='g',marker='.')
plt.show()

结果：

10. PCA——算法优缺点

优点：

1. 完全无参数限制的。在PCA的计算过程中完全不需要人为的设定参数或是根据任何经验模型对计算 进行干预，最后的结果只与数据相关，与用户是独立的。
2. 用PCA技术可以对数据进行降维，同时对新求出的“主元”向量的重要性进行排序，根据需要取前面 最重要的部分，将后面的维数省去，可以达到降维从而简化模型或是对数据进行压缩的效果。同时 最大程度的保持了原有数据的信息。
3. 计算方法简单，易于在计算机上实现。

缺点：

1. 如果用户对观测对象有一定的先验知识，掌握了数据的一些特征，却无法通过参数化等方法对处理 过程进行干预，可能会得不到预期的效果，效率也不高。