PCA的数学原理和python实现

最近学习了一下PCA，具体原理网址如下：

CodingLabs - PCA的数学原理

主要原理是通过线性变换将原始数据变换为一组各维度线性无关的表示，其中将方差最大的方向作为主要特征。提取数据的主要特征分量，可用于高维数据的降维

主要算法步骤是

PCA算法

设有m条n维数据。

1）将原始数据按列组成n行m列矩阵X

2）将X的每一行（代表一个属性字段）进行零均值化，即减去这一行的均值

3）求出协方差矩阵

4）求出协方差矩阵的特征值及对应的特征向量

5）将特征向量按对应特征值大小从上到下按行排列成矩阵，取前k行组成矩阵P

6）Y=PX即为降维到k维后的数据

参考了网上的代码，用python实现了一遍，对图像进行处理，用列作为特征，行作为样本。

如下图451*300的图像，300个样本，每个样本451个特征。

按照pca算法，提取前k个主要特征合并出的图像如下

k=5

 

 k=15

k=25

 

k=35

 

 python代码如下

import numpy as np
from skimage import io,data

img=data.chelsea()
print(img.shape)
io.imsave("cat.bmp",img)
img=img/255
io.imshow(img)
io.show()

def pca_m(data,k):
    n_samples, n_features = data.shape  # 每列是一个特征，每一行是一个样本
    # 求解一个特征的均值，即每一列的平均值，shape(n_samples,1)
    mean = (data.sum(axis=1) / n_samples).reshape(-1, 1)
    # 去中心化,shape(n_samples,1)
    normal_data = data - mean
    # 得到协方差矩阵
    matrix_ = np.dot(np.transpose(normal_data), normal_data)
    eig_val, eig_vec = np.linalg.eig(matrix_)
    # 得到最主要的k个特征
    eigIndex = np.argsort(eig_val)
    eigVecIndex = eigIndex[:-(k + 1):-1]
    feature = eig_vec[:, eigVecIndex]

    new_data = np.dot(normal_data, feature)
    # 将降维后的数据映射回原空间
    rec_data = np.dot(new_data, np.transpose(feature)) + mean
    return rec_data


R,G,B=img[:,:,0],img[:,:,1],img[:,:,2]

for k in range(5,80,10):
    R_new,G_new,B_new=pca_m(R,k),pca_m(G,k),pca_m(B,k)
    img_new = np.zeros(img.shape)
    img_new[:, :, 0] = R_new
    img_new[:, :, 1] = G_new
    img_new[:, :, 2] = B_new

    io.imshow(img_new)
    io.show()
    imgimg = np.maximum(img_new,0)
    imgsave=np.uint8(imgimg*255)
    io.imsave("cat_k"+str(k)+".bmp",imgsave)