一、基本原理
在机器学习中,当维数较高时会出现数据样本稀疏、距离计算困难等问题,这类问题被称为维数灾难(curse of dimensionality)。缓解维数灾难的一个重要途径是降维(dimension reduction),即通过变换将高维空间转变为低维子空间。
主成分分析(Principal Component Analysis, PCA)是一种最常用的降维方法。若要使用超平面对高维空间进行表达,一般具有以下性质:最大重构性,即样本点到这个超平面距离都足够近;最大可分性,即样本点在超平面上投影尽可能分开。下面即从这两个方面说明:
最大重构性:
原样本点与投影重构的样本点距离为:
其中重构样本点,投影,是新坐标系,是标准正交基向量
考虑到是标准正交基向量,是协方差矩阵,故主成分分析的优化目标为:
最大可分性:
样本点在超平面上的投影是,若使投影点尽可能分开,即使投影后样本点方差最大化,即,次时优化目标可写为:
两者思想实现了统一,下图则显示了这种关系:
PCA算法描述如下:
PCA算法
二、算法实现
以下为PCA算法的具体实现过程,通过奇异值分解(singular value decomposition, SVD)对协方差矩阵做特征值分解,并选取其中主要分量进行数据降维。将二维数据降到一维。
import numpy as np
import matplotlib.pyplot as plt
import scipy.io as spio
def PCA_2D():
data_2d = spio.loadmat('data.mat')
X = data_2d['X']
m = X.shape[0]
X_copy = X.copy()
X_norm,mu,sigma = Normalize(X_copy)
cov_mat = np.cov(X_norm.T)
U,S,V = np.linalg.svd(cov_mat)
plt.plot(X[:,0],X[:,1],'bo')
plt.plot(np.array(mu[0],(mu+S[0]*U[:,0])[0]),np.array(mu[1],(mu+S[0]*U[:,0])[1]),'r-')
plt.show()
K = 1
Z = mapData(X_norm,U,K)
X_rec = recoverData(Z,U,K)
plt.plot(X_norm[:,0],X_norm[:,1],'bo')
plt.plot(X_rec[:,0],X_rec[:,1],'ro')
for i in range(m):
plt.plot([X_norm[i,:][0],X_rec[i,:][0]],[X_norm[i,:][1],X_rec[i,:][1]],'--')
plt.axis('square')
plt.show()
def Normalize(X):
n = X.shape[1]
mu = np.zeros((1,n))
sigma = np.zeros((1,n))
mu = np.mean(X,axis=0)
sigma = np.std(X,axis=0)
for i in range(n):
X[:,i] = (X[:,i]-mu[i])/sigma[i]
return X,mu,sigma
def mapData(X_norm,U,K):
Z = np.zeros((X_norm.shape[0],K))
U_reduce = U[:,0:K]
Z = np.dot(X_norm,U_reduce)
return Z
def recoverData(Z,U,K):
X_rec = np.zeros((Z.shape[0],U.shape[0]))
U_reduce = U[:,0:K]
X_rec = np.dot(Z,U_reduce.T)
return X_rec
if __name__=='__main__':
PCA_2D()
下图绘制了原始数据图和数据降维的映射关系,说明了降维过程及其主方向。
原始数据图 降维后的映射关系
PCA同样可用于图片的压缩,以下为通过sklearn的实现过程。
import numpy as np
import matplotlib.pyplot as plt
import scipy.io as spio
from sklearn.decomposition import pca
from sklearn.preprocessing import StandardScaler
def PCA_FACE():
image_data = spio.loadmat('data_faces.mat')
X = image_data['X']
draw_img(X[0:25,:])
scaler = StandardScaler()
scaler.fit(X)
x_train = scaler.transform(X)
K = 200
model = pca.PCA(n_components=K).fit(x_train)
Z = model.transform(x_train)
U_reduce = model.components_
draw_img(U_reduce[0:25,:])
X_rec = np.dot(Z,U_reduce)
draw_img(X_rec[0:25,:])
def draw_img(imgdata):
num = 0
m,n = imgdata.shape
width = int(round(np.sqrt(n)))
height = int(n/width)
rows = int(np.floor(np.sqrt(m)))
cols = int(np.ceil(m/rows))
pad = 1
draw_arr = -np.ones((pad+rows*(height+pad),pad+cols*(width+pad)))
for i in range(rows):
for j in range(cols):
max_val = max(np.abs(imgdata[num,:]))
draw_arr[pad+i*(height+pad):pad+i*(height+pad)+height,pad+j*(width+pad):pad+j*(width+pad)+width] = imgdata[num,:].reshape(height,width,order='F')/max_val
num += 1
plt.imshow(draw_arr)
plt.axis('off')
plt.show()
if __name__=='__main__':
PCA_FACE()
下图分别为原图、压缩后的图、恢复后的图,说明了PCA的在图片压缩的作用。
原图 压缩后的图 恢复后的图另外,可以通过mglearn模块直观说明PCA的过程。
import mglearn
mglearn.plots.plot_pca_illustration()
参考资料
[1] https://github.com/lawlite19/MachineLearning_Python
[2] 周志华 著. 机器学习. 北京:清华大学出版社,2016
[3] 李航 著. 统计学习方法. 北京:清华大学出版社,2012
[4] 史春奇等 著. 机器学习算法背后的理论与优化. 北京:清华大学出版社,2019
[5] Peter Harrington 著. 李锐等 译. 机器学习实战. 北京:人民邮电出版社,2013
何处是归程?长亭更短亭。——李白《菩萨蛮》
网友评论