Kernel PCA

参考：http://blog.csdn.net/wolenski/article/details/7961703

主要思想：

1.PCA中的内积运算可以通过核函数完成

2. 假设向量X1,X2通过某种映射F得到高维空间中的Y1,Y2，利用核技巧，可以使高维空间中的内积运算<Y1,Y2>通过核函数Q(X1,X2)即可得到结果。

假设存在训练样本B ，那么根据PCA的原理就可以得到公式

C v=r *v

其中C为样本B的协方差矩阵，r,v分别为特征值和特征向量

根据特征值可由样本线性张成（不进行具体推导），最终得到

借助核技巧 既可以化简成

                Ka=λa

计算核函数空间中的特征值和特征向量，通过排序特征值，既可以得到原始样本在高维空间中的主成份。其中K为核函数映射得到的协方差矩阵

以上只是证明了核函数的可操作性。

实际使用核函数应当：

显然在推导传统PCA时需要先对训练样本数据进行中心化再求协方差，因此样本数据中心化后映射到核函数空间F时仍需要进行中心化得到Kc，然而映射方法的显示结果我们并不知道，因此可以通过一些列推导，我们可以通过未中心化的映射数据得到矩阵K，再对矩阵按照一个公式进行中心化，绕过协方差矩阵在运算前的中心化，直接得到中心化结果Kc

中心化方法网上有很多：Kc=K-1(N)K-K1(N)+1(N)K1(N)

对Kc求特征值和特征空间，排序得到核函数空间主成份的特征向量空间Vc

当得到新测试数据 t 时，按照传统PCA，并根据线性张成，得到 u*t= ∑α(x’ t)  因此很明显(x’ t) 内积可以使用核函数运算，同样得到另一映射后的新矩阵Q，该矩阵同样需要中心化得到Qc,将该矩阵Qc乘以训练样本核函数空间特征向量空间Vc，既可以得到测试样本在高维空间中的主成份。这样既可以解决一些线性难分问题。

杭州趋数网络科技有限公司 数据挖掘/数据分析师花京华

参考Matlab代码

http://blog.sina.com.cn/s/blog_74e58b1001019uej.html

function [train_kpca,test_kpca] = kpcaFordata(train,test,threshold,rbf_var)
%% Data kpca processing
%% 源地址：http://blog.sina.com.cn/lssvm
%% 函数默认设置
if nargin <4
rbf_var=10000;%？
end
if nargin <3
threshold = 90;
end
%% 数据处理
patterns=zscore(train); %训练数据标准化
test_patterns=zscore(test); %测试数据标准化
train_num=size(patterns,1); %train_num是训练样本的个数
test_num=size(test_patterns,1);%test_num是测试样本的个数
cov_size = train_num; %cov_size是训练样本的个数
%% 计算核矩阵
for i=1:cov_size,
for j=i:cov_size,
K(i,j) = exp(-norm(patterns(i,:)-patterns(j,:))^2/rbf_var); %核函数 rbf_var ？？
K(j,i) = K(i,j);
end
end
unit = ones(cov_size, cov_size)/cov_size;%cov_size是样本的个数
%% 中心化核矩阵
K_n = K - unit*K - K*unit + unit*K*unit;% 中心化核矩阵
%% 特征值分解
[evectors_1,evaltures_1] = eig(K_n/cov_size);
[x,index]=sort(real(diag(evaltures_1))); %sort每行按从小到大排序，x为排序后结果，index为索引
evals=flipud(x) ;% flipud函数实现矩阵的上下翻转
index=flipud(index);
%% 将特征向量按特征值的大小顺序排序
evectors=evectors_1(:,index);
%% 单位化特征向量
% for i=1:cov_size
% evecs(:,i) = evectors(:,i)/(sqrt(evectors(:,i)));
% end
train_eigval = 100*cumsum(evals)./sum(evals);
index = find(train_eigval >threshold);

train_kpca = zeros(train_num, index(1)); %train_num是训练样本的个数
%% evecs单位化后的特征矩阵,K_n训练数据的中心化核矩阵
train_kpca=[K_n * evectors(:,1:index(1))];
%% 重建测试数据
unit_test = ones(test_num,cov_size)/cov_size;%cov_size是训练样本的个数
K_test = zeros(test_num,cov_size); %test_num是测试样本的个数， cov_size是训练样本的个数
for i=1:test_num, %test_num是测试样本的个数
for j=1:cov_size,%cov_size是训练样本的个数
K_test(i,j) = exp(-norm(test_patterns(i,:)-patterns(j,:))^2/rbf_var);
end
end
K_test_n = K_test - unit_test*K - K_test*unit + unit_test*K*unit;
test_kpca = zeros(test_num, index(1));%test_num是测试样本的个数
test_kpca = [K_test_n * evectors(:,1:index(1))];