主成分分析PCA

数据降维的意义

1.维数灾难

高维数据空间给机器学习算法带来的计算量是十分巨大的。事实上，在高维情形下出现的数据样本稀疏、距离计算困难等问题，是所有机器学习方法共同面临的严重障碍，被称为“维数灾难”。

缓解维数灾难的一个重要途径就是降维，即通过某种数学变换将原始高维属性空间转变为一个低维“子空间”，这个子空间内样本密度大幅提高，距离计算也变得更为容易。

为什么能进行降维呢？这是因为在很多时候，人们观测或收集到的数据样本虽是高维的，但与学习任务密切相关的也许仅是某个低维分布，即高维空间的一个低维“嵌入”，在这个低维空间学习起来效率更高，如下图所示。

2.数据可视化

当数据处于高维度时，我们比较难直观地发现其特性和各数据维度的关系。但当通过降维处理使得数据只有两维或者三维时可以通过画图寻找其中存在的规律。

例如下面这个反应各个国家经济发展状态的例子，不同国家对应二维图上的一个点。

    

PCA流程

PCA其实是在寻找一个超平面，这个超平面具有如下性质：

1.最近重构性：样本点到这个超平面的距离都足够近

2.最大可分性：样本点在这个超平面上的投影尽可能分开

得到某一数据样本，根据以上原则建立代价函数，求代价函数的极小值。（数学过程较复杂）

         PCA算法流程，较其数学过程简单太多了：

怎么选择低维数k?

一般选择k的方式是根据我们所要求保留的数据差异性。例如原数据为X（m*n），投影为Xapprox（m*k，另n-k维为0）。保留数据差异性定义如下（99%保留）：

在PCA过程2中的协方差矩阵记作COV，过程对COV做特征值分解，一般使用SVD（奇异值分解，更加稳定），即：

[U,S,V] = SVD（COV）；

其中U=[w1,w2,..wn]为n个特征向量，S为特征向量对应的特征值矩阵（主对角为特征值，其他位置为0），如下图。

k为降维数，则取前k个特征向量P = [w1,w2...wk]即为投影矩阵，Xapprox = X*P;保留差异性等于（1-保留k个特征向量特征值的和比上所有特征值的和），即：

由此可确定K。

如何恢复数据维度？

X' (m*n) = Xapprox(m*k)*P'(k*n)即可，但是注意此时恢复的数据与原数据不再相同，差异性由投影时舍去的信息决定。  

PCA和线性回归的差异

直观上，PCA和线性回归很像。以二维数据降到一维为例，PCA和线性回归都是找一条直线，尽量使得各个数据与直线距离最近，那么差异在哪呢？可以看下图

差异就在于线性回归最小化的是坐标点到直线的垂直距离(即||δy||)，而PCA最小化的是欧式距离（即||δx,δy||）。如果不懂的话再好好想想回归函数和PCA的定义。

Pca的Matlab实现：

function [mappedX, mapping] = pca(X, no_dims)%PCA Perform the PCA algorithm%%   [mappedX, mapping] = pca(X, no_dims)%% The function runs PCA on a set of datapoints X. The variable% no_dims sets the number of dimensions of the feature points in the % embedded feature space (no_dims >= 1, default = 2). % For no_dims, you can also specify a number between 0 and 1, determining % the amount of variance you want to retain in the PCA step.% The function returns the locations of the embedded trainingdata in % mappedX. Furthermore, it returns information on the mapping in mapping.%%% This file is part of the Matlab Toolbox for Dimensionality Reduction.% The toolbox can be obtained from http://homepage.tudelft.nl/19j49% You are free to use, change, or redistribute this code in any way you% want for non-commercial purposes. However, it is appreciated if you % maintain the name of the original author.%% (C) Laurens van der Maaten, Delft University of Technology    if ~exist('no_dims', 'var')        no_dims = 2;    end% Make sure data is zero mean    mapping.mean = mean(X, 1);X = bsxfun(@minus, X, mapping.mean);% Compute covariance matrix    if size(X, 2) < size(X, 1)        C = cov(X);          else        C = (1 / size(X, 1)) * (X * X');        % if N>D, we better use this matrix for the eigendecomposition    end% Perform eigendecomposition of CC(isnan(C)) = 0;C(isinf(C)) = 0;    [M, lambda] = eig(C);        % Sort eigenvectors in descending order    [lambda, ind] = sort(diag(lambda), 'descend');    if no_dims < 1        no_dims = find(cumsum(lambda ./ sum(lambda)) >= no_dims, 1, 'first');        disp(['Embedding into ' num2str(no_dims) ' dimensions.']);    end    if no_dims > size(M, 2)        no_dims = size(M, 2);        warning(['Target dimensionality reduced to ' num2str(no_dims) '.']);    endM = M(:,ind(1:no_dims));    lambda = lambda(1:no_dims);% Apply mapping on the data    if ~(size(X, 2) < size(X, 1))        M = bsxfun(@times, X' * M, (1 ./ sqrt(size(X, 1) .* lambda))');     % normalize in order to get eigenvectors of covariance matrix    end    mappedX = X * M;        % Store information for out-of-sample extension    mapping.M = M;mapping.lambda = lambda;