PCA(Principal Component Analysis)

PCA的本质就是找一些投影方向，使得数据在这些投影方向上的方差最大，而且这些投影方向是相互正交的。

这其实就是找新的正交基的过程，计算原始数据在这些正交基上投影的方差，方差越大，就说明在对应正交基上包含了更多的信息量。

而原始数据协方差矩阵的特征值越大，对应的方差越大，在对应的特征向量上投影的信息量就越大。反之，如果特征值较小，

则说明数据在这些特征向量上投影的信息量很小，可以将小特征值对应方向的数据删除，从而达到了降维的目的。

算法思想

将一组P维向量降为K维(K

各字段两两间协方差为0，而字段的方差则尽可能大——这不就是对角矩阵吗？

Y = XP，其中X为原始数据(每个维度按列组成)，Y为X对P做基变换后的数据，P的基行列组成。

设X、Y的协方差矩阵(对称矩阵)分别为C、D，则：D=PCP’。

我们希望D为对角矩阵，则D为C对P做相似对角化所得。

由于协方差矩阵为对称矩阵，则构成相似变换矩阵P的特征向量是正交的，实现逆变换。

实验步骤

(1)分别求每个维度的平均值，然后对于所有的样例，都减去对应维度的均值，得到去中心化的数据；

(2)求协方差矩阵C：用去中心化的数据矩阵乘上它的转置，然后除以(N-1)即可，N为样本数量；

(3)求协方差的特征值和特征向量；

(4)将特征值按照从大到小排序，选择前k个，然后将其对应的k个特征向量分别作为列向量组成特征向量矩阵；

(5)将样本点从原来维度投影到选取的k个特征向量，得到低维数据；

(6)通过逆变换，重构低维数据，进行复原。

import numpy as np

创建 eigValPct(eigVals, percentage)

通过方差的百分比来计算将数据降到多少维。

函数传入的参数是特征值eigVals和百分比percentage，返回需要降到的维度数num

def eigValPct(eigVals, percentage):

sortArray=np.sort(eigVals)[::-1] # 特征值从大到小排序

pct = np.sum(sortArray)*percentage

tmp = 0

num = 0

for eigVal in sortArray:

tmp += eigVal

num += 1

if tmp>=pct:

return num

创建 im_PCA(dataMat, percentage=0.9)

函数有两个参数，其中dataMat是已经转换成矩阵matrix形式的数据集，每列表示一个维度；

其中的percentage表示取前多少个特征需要达到的方差占比，默认为0.9。

注意：

np.cov(dataMat, rowvar=False)，按照rowvar的默认值，会把一行当成一个特征，一列当成一个样本。

def im_PCA(dataMat, percentage=0.9):

meanVals = np.mean(dataMat, axis=0)

meanRemoved = dataMat - meanVals

# 这里不管是对去中心化数据 or 原始数据计算协方差矩阵，结果都一样，特征值大小会变，但相对大小不会改变

# covMat = np.cov(dataMat, rowvar=False)

# 标准的计算需要除以(dataMat.shape[0]-1)，不算也不会影响结果，理由同上

covMat = np.dot(np.transpose(meanRemoved), meanRemoved)

eigVals, eigVects = np.linalg.eig(np.mat(covMat))

k = eigValPct(eigVals,percentage) #要达到方差的百分比percentage，需要前k个向量

print('K =', k)

eigValInd = np.argsort(eigVals)[::-1] #对特征值eigVals从大到小排序

eigValInd = eigValInd[:k]

redEigVects = eigVects[:,eigValInd] #主成分

lowDDataMat = meanRemoved*redEigVects #将原始数据投影到主成分上得到新的低维数据lowDDataMat

reconMat = (lowDDataMat*redEigVects.T)+meanVals #得到重构数据reconMat

return lowDDataMat, reconMat

注意：

图像Matrix格式必须转换为uint8格式，否则cv2.imshow时图像不能正常显示！np.array(reconMat, dtype=‘uint8’)。

但是，强制类型转化后会丢失信息，比如将 6. 变成 5 ，因为强制类型转化是直接采用截断二进制位的方式。。。

import cv2

img = cv2.imread('ava.jpg')

blue = img[:,:,0]

dataMat = np.mat(blue)

lowDDataMat, reconMat = im_PCA(dataMat, 1)

print('原始数据', blue.shape, '降维数据', lowDDataMat.shape)

# print(dataMat)

# print(reconMat)

# 格式必须转换为uint8格式，这里丢失了很多信息！！！

# reconMat = np.array(reconMat, dtype='uint8')

cv2.imshow('blue', blue)

cv2.imshow('reconMat', np.array(reconMat, dtype='uint8'))

cv2.waitKey(0)

K = 426

原始数据 (640, 426) 降维数据 (640, 426)

计算降维误差率

根据PCA的误差计算方式：信息丢失率 =

注意：np.dot(A, B)

(1)如果处理的是一维数组，则得到的是两数组的内积；

(2)如果是二维数组/矩阵之间的运算，得到的是矩阵积；

def PrintError(data, recdata):

sum1 = 0

sum2 = 0

D_value = data - recdata # 计算两幅图像之间的差值矩阵

# 计算两幅图像之间的误差率，即信息丢失率

for i in range(data.shape[0]):

sum1 += np.dot(data[i],data[i])

sum2 += np.dot(D_value[i], D_value[i])

print('丢失信息量：', sum2)

print('原始信息量：', sum1)

print('信息丢失率：', sum2/sum1)

print(dataMat)

print(reconMat)

PrintError(np.array(blue, dtype='double'), np.array(reconMat, dtype='double'))

[[ 6 6 6 ... 11 13 13]

[ 6 6 6 ... 12 15 15]

[ 7 7 7 ... 13 17 17]

...

[16 16 16 ... 29 29 29]

[16 16 16 ... 29 29 29]

[16 16 16 ... 29 29 29]]

[[ 6. 6. 6. ... 11. 13. 13.]

[ 6. 6. 6. ... 12. 15. 15.]

[ 7. 7. 7. ... 13. 17. 17.]

...

[16. 16. 16. ... 29. 29. 29.]

[16. 16. 16. ... 29. 29. 29.]

[16. 16. 16. ... 29. 29. 29.]]

丢失信息量： 7.140138346734739e-16

原始信息量： 784219850.0

信息丢失率： 9.104766152928595e-25

使用sklearn的PCA方法

参数说明：

sklearn.decomposition.PCA(n_components=None, copy=True, whiten=False);

n_components: PCA算法中所要保留的主成分个数n，也即保留下来的特征个数n。赋值为int，比如n_components=1，将把原始数据降到一个维度。赋值为string，比如n_components=‘mle’，将自动选取特征个数n，使得满足所要求的方差百分比。

copy: 是否在原始数据的副本上进行运算。

whiten: 白化，使得每个特征具有相同的方差。

from sklearn.decomposition import PCA

pca = PCA(n_components=426).fit(blue)

# 降维

x_new = pca.transform(blue)

# 还原降维后的数据到原空间

recdata = pca.inverse_transform(x_new)

print(recdata)

# 计算误差

PrintError(np.array(blue, dtype='double'), np.array(reconMat, dtype='double'))

cv2.imshow('sklearn-recdata', np.array(recdata, dtype='uint8'))

cv2.waitKey(0)

[[ 6. 6. 6. ... 11. 13. 13.]

[ 6. 6. 6. ... 12. 15. 15.]

[ 7. 7. 7. ... 13. 17. 17.]

...

[16. 16. 16. ... 29. 29. 29.]

[16. 16. 16. ... 29. 29. 29.]

[16. 16. 16. ... 29. 29. 29.]]

丢失信息量： 7.140138346734739e-16

原始信息量： 784219850.0

信息丢失率： 9.104766152928595e-25

Reference

计算过程

https://www.cnblogs.com/lzllovesyl/p/5235137.html

PCA详解，原理及步骤，包括数学基础

https://my.oschina.net/gujianhan/blog/225241#OSC_h2_1

PCA，图像压缩还原

https://blog.csdn.net/Vincent_zbt/article/details/88648739

scikit-learn中PCA的使用方法

https://blog.csdn.net/u012162613/article/details/42192293/