主成分分析法是一个非监督的机器学习算法,主要用于数据的降维。通过降维,可以发现更便于人类理解的特征。
求解目标
主成分分析的步骤:
- 对样本进行demean处理(使所有样本的均值为0)
-
取一个轴的方向 w = (w1,w2...,wn),使我们的样本,映射到w之后,使下式最大
由于均值为0,则只需要使下式最大
映射的过程可以如下示意
w为单位向量,则有
则最终目标为求取一个w,使得下式最大
上式中为向量相乘,假设数据有n个维度,展开来是这个样子
那么,这就变成了一个目标函数的最优化问题,可以使用梯度上升法解决
这个过程看起来跟线性回归很像,其实是不同的,需要注意
梯度上升
梯度上升的过程与梯度下降是类似地,需要先求导
像前面一样,可以化为矩阵运算的形式。设X为这样的矩阵
则求导可以写成这样的形式。这里就不推导了
整理一下就是这个样子
求取n个主成分
上面的操作中求取w的是第一个主成分,称为第一主成分。如果要求取第二主成分,则需要将数据在第一个主成分上的分量去掉,得到的数据再求取主成分,就得到了第二主成分。
比如说有这么一个数据集
找到第一主成分的方向w1。将数据在第一主成分上的分量去掉,得到的数据如下。再求取一次主成分,就得到了第二主成分的方向。
假设我们找到了k个主成分的方向,如果想将数据从n维映射到k维(n>k),则可以如下进行。令Wk为计算出的k个主成分的方向。
则可以进行矩阵运算进行降维。X有m个样本n个方向,映射后有m个样本k个方向
同样,可以将映射到低维的数据重映射到高维,但是会损失一些信息,结果跟原来是不一样的
编程实现
"""
Created by 杨帮杰 on 11/4/2018
Right to use this code in any way you want without
warranty, support or any guarantee of it working
E-mail: yangbangjie1998@qq.com
Association: SCAU 华南农业大学
"""
import numpy as np
class PCA:
def __init__(self, n_components):
"""初始化PCA"""
assert n_components >= 1, "n_components must be valid"
self.n_components = n_components
self.components_ = None
def fit(self, X, eta=0.01, n_iters=1e4):
"""获得数据集的前n个主成分"""
assert self.n_components <= X.shape[1], \
"n_components must not be greater than the feature number of X"
def demean(X):
return X - np.mean(X, axis=0)
def f(w, X):
return np.sum((X.dot(w) ** 2)) / len(X)
def direction(w):
return w / np.linalg.norm(w)
def first_components(X, initial_w, eta=0.01, n_iters=1e4, epsilon=1e-8):
w = direction(initial_w)
cur_iter = 0
while cur_iter < n_iters:
gradient = df(w, X)
last_w = w
w = w + eta*gradient
w = direction(w)
if(abs(f(w, X) - f(last_w, X)) < epsilon):
break
cur_iter += 1
return w
X_pca = demean(X)
self.components_ = np.empty(shape=(self.n_components, X.shape[1]))
for i in range(self.n_components):
initial_w = np.random.random(X_pca.shape[1])
w = first_components(X_pca, initial_w, eta, n_iters)
self.components_[i,:] = w
X_pca = X_pca - X_pca.dot(w).reshape(-1,1) * w
return self
def transform(self, X):
"""将给定的X,映射到各个主成分分量中"""
assert X.shape[1] == self.components_.shape[1]
return X.dot(self.components_.T)
def inverse_transform(self, X):
"""将给定的X,反向映射回原来的特征空间"""
assert X.shape[1] == self.components_.shape[0]
return X.dot(self.components_)
def __repr__(self):
return "PCA(n_components = %d" % self.n_components