主成分分析pca和白化

最近在学习dnn，autoEncoding有点问题，先把pca做了，这个不太难

作业见http://deeplearning.stanford.edu/wiki/index.php/Exercise:PCA_in_2D

原始数据集

-6.7644914e-01  -6.3089308e-01  -4.8915202e-01  -4.8005424e-01  -3.7842021e-01  -3.3788391e-01  -3.2023528e-01  -3.1108837e-01  -2.3145555e-01  -1.9623727e-01  -1.5678926e-01  -1.4900779e-01  -1.0861557e-01  -1.0506308e-01  -8.0899829e-02  -7.1157518e-02  -6.3251073e-02  -2.6007219e-02  -2.2553443e-02  -5.8489047e-03  -4.3935323e-03  -1.7309716e-03   7.8223728e-03   7.5386969e-02   8.6608396e-02   9.6406046e-02   1.0331683e-01   1.0531131e-01   1.1493296e-01   1.3052813e-01   1.6626253e-01   1.7901863e-01   1.9267343e-01   1.9414427e-01   1.9770003e-01   2.3043613e-01   3.2715844e-01   3.2737163e-01   3.2922364e-01   3.4869293e-01   3.7500704e-01   4.2830153e-01   4.5432503e-01   5.4422436e-01   6.6539963e-01
-4.4722050e-01  -7.4778067e-01  -3.9074344e-01  -5.6036362e-01  -3.4291940e-01  -1.3832158e-01   1.2360939e-01  -3.3934986e-01  -8.2868433e-02  -2.4759514e-01  -1.0914760e-01   4.2243921e-01  -5.2329327e-02  -2.0126541e-01   1.3016657e-01   1.2293321e-01  -3.4787750e-01  -1.4584897e-01  -1.0559656e-01  -5.4200847e-02   1.6915422e-02  -1.1069762e-01   9.0859816e-02   1.5269096e-01  -9.4416463e-02   1.5116385e-01  -1.3540126e-01   2.4592698e-01   5.1087447e-02   2.4583340e-01  -5.9535372e-02   2.9704742e-01   1.0168115e-01   1.4258649e-01   1.0662592e-01   3.1698532e-01   6.1577841e-01   4.3911172e-01   2.7156501e-01   1.3572389e-01   3.1918066e-01   1.5122962e-01   3.4979047e-01   6.2316971e-01   5.2018811e-01

第一步：获取数据

# -*- coding: utf-8 -*
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
import math
# %%================================================================
# %% Step 0: Load data
# %  We have provided the code to load data from pcaData.txt into x.
# %  x is a 2 * 45 matrix, where the kth column x(:,k) corresponds to
# %  the kth data point.Here we provide the code to load natural image data into x.

def normal(x):
return 1.0/math.sqrt(x)

normal_ufunc = np.vectorize(normal,otypes=[np.float])

x = np.loadtxt('pcaData.txt')
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter (x[0,:] ,x[1,:])
plt.show()

第二步：平移使数据的均值为0

# make mean to zero
xmean = np.mean(x,axis=1).reshape((-1,1))
y = x -xmean

第三步：计算协方差矩阵

这里需要注意，cov函数需要传递参数bias=1，n个数据的协方差需要除以n，若不传递bias=1，将除以(n-1)否则

# 协方差矩阵
xx = np.cov(y,bias=1)
print (u'原始协方差矩阵',xx)

第四步：计算特征值特征向量

因为协方差矩阵是对称矩阵，因此它的特征值等于奇异值，故采用奇异值分解（svd）来计算

# 特征值或奇异值分解
u,sigma,v = np.linalg.svd(xx)
print ('u sigma v=',u,sigma,v)

第五步：旋转数据

U是特征向量构成的正交矩阵，用u左乘原始数据矩阵相当于求原始数据在以u的列向量为标准正交基下的坐标

因这里是2维数据，直接将其压缩为1维

# 旋转
xrot = np.dot(np.transpose (u),y)
xrotcov = np.cov(xrot,bias=1)
print ('xrot=',xrot)
print (u'旋转后协方差矩阵',xrotcov)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(xrot[0,:],xrot[1,:]-xrot[1,:])
plt.show()

关于pca的压缩就到这里。

下面是pca白化

白化的目标，特征之间相关性最低，所有特征有相同的方差。通常使协方差矩阵为单位矩阵

svd中分解出来的sigma为对角矩阵，这里sigma=array([ 0.16071768, 0.01528232])

为了使每个输入特征具有单位方差，我们可以直接使用

作为缩放因子来缩放每个特征

。具体地，我们定义白化后的数据

如下：



实践中需要实现PCA白化或ZCA白化时，有时一些特征值

在数值上接近于0，这样在缩放步骤时我们除以

将导致除以一个接近0的值；这可能使数据上溢
(赋为大数值)或造成数值不稳定。因而在实践中，我们使用少量的正则化实现这个缩放过程，即在取平方根和倒数之前给特征值加上一个很小的常数

：



当

在区间

上时,
一般取值为

。

# pca白化
epsilon = 1e-5
scal = np.diag (normal_ufunc(sigma+epsilon))
xpcawhite = np.dot(scal,xrot)
xpcawhitecov = np.cov(xpcawhite,bias=1)
print ('xpcawhite',xpcawhite)
print('xpcawhitecov',xpcawhitecov)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(xpcawhite[0,:],xpcawhite[1,:])
plt.show()

绘制出

，我们得到:



注意pca白化会缩放特征范围

zca白化

用任意正交矩阵左乘pca白化后的数据，得到的任然是单位协方差矩阵数据，这就是zca白化

# zcawhite
xzcawhite = np.dot(u,xpcawhite)
xzcawhitecov = np.cov(xzcawhite,bias=1)
print('xzcawhite',xzcawhite)
print ('xzcawhitecov',xzcawhitecov)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(xzcawhite[0,:],xzcawhite[1,:])
plt.show()