解釋一下核主成分分析(Kernel Principal Component Analysis, KPCA)的公式推導過程~
KPCA，中文名稱”核主成分分析“，是對PCA演算法的非線性擴充套件，言外之意，PCA是線性的，其對於非線性資料往往顯得無能為力，例如，不同人之間的人臉影象，肯定存在非線性關係，自己做的基於ORL資料集的實驗，PCA能夠達到的識別率只有88%，而同樣是無監督學習的KPCA演算法，能夠輕鬆的達到93%左右的識別率（雖然這二者的主要目的是降維，而不是分類，但也可以用於分類），這其中很大一部分原因是，KPCA能夠挖掘到資料集中蘊含的非線性資訊。
今天突然心血來潮，想重新推導一下KPCA的公式，期間遇到了幾個小問題，上部落格查閱，發現目前並沒有一個專注於KPCA公式推導的文章，於是決定寫一篇這樣的部落格（轉載請註明：http://blog.csdn.net/wsj998689aa/article/details/40398777）。
1. 理論部分
KPCA的公式推導和PCA十分相似，只是存在兩點創新：
1. 為了更好地處理非線性資料，引入非線性對映函式
2. 引入了一個定理：空間中的任一向量（哪怕是基向量），都可以由該空間中的所有樣本線性表示，這點對KPCA很重要，我想大概當時那個大牛想出KPCA的時候，這點就是它最大的靈感吧。話說這和”稀疏“的思想比較像。
假設中心化後的樣本集合X（d*N，N個樣本，維數d維，樣本”按列排列“），現將X對映到高維空間，得到
於是乎！假設D（D >> d）維向量
和PCA太像了吧？這個時候，在利用剛才的定理，將特徵向量
然後，在把
進一步，等式兩邊同時左乘
你可能會問，這個有啥用？
這樣做的目的是，構造兩個
第二個等號，是源於核函式的性質，核函式比較多，有如下幾種：
於是，公式進一步變為如下形式：
兩邊同時去除K，得到了PCA相似度極高的求解公式：
求解公式的含義就是求K最大的幾個特徵值所對應的特徵向量，由於K為對稱矩陣，所得的解向量彼此之間肯定是正交的。
但是，請注意，這裡的
這時有的朋友可能會問，這個時候，如果給定一個測試樣本
是這樣的，既然我們可以得到高維空間的一組基
於是呼~就可以對
2. 實驗部分
做了一些模擬實驗，分別比較了PCA與KPCA之間的效果，KPCA基於不同核函式的效果，二者對於原始資料的要求，以及效果隨著引數變化的規律。
1）下面展示的是“無重疊的”非線性可分資料下，PCA與KPCA（基於高斯核）的區別，注意，原始資料是二維資料，投影之後也是二維資料
2）下面展示的是“部分重疊的”非線性可分資料下，PCA與KPCA的區別
3）下面展示的是“無高斯擾動的”非線性可分資料下，PCA與KPCA的區別
4）下面展示的是上述三類資料下，基於多項式核函式的KPCA效果
5）下面展示的是在“部分重疊的”非線性可分資料下，基於多項式核函式的KPCA在不同多項式引數下的效果圖
3. 實驗結論
1. 從2.1中我們可以看出，PCA與KPCA對於非線性資料各自的處理能力，仔細觀察PCA其實只對原始資料進行了旋轉操作，這是由於其尋找的是資料的“主要分佈方向”。KPCA可以將原始資料投影至線性可分情況，其原因就是第一部分所說的內容。 2. 至於為何將資料分為“無重疊”，“部分重疊”，“無高斯擾動”，是自己在試驗中發現，對於部分重疊的資料，KPCA不能將資料投影至完全線性可分的程度（2.3第三幅圖中，不同類別資料仍舊存在重疊現象），這說明KPCA只是個無監督的降維演算法，它不管樣本的類別屬性，只是降維而已。 3. 這裡提供了高斯核與多項式核的效果，我們很容易發現，二者的效果有很大不同，這直觀地說明不同核函式具有不同的特質。並且，針對於無高斯擾動資料，始終沒有找到引數p，有可能針對這類資料，多項式核函式無能為力。 4. 2.5中展示了多項式核的引數影響，我們可以發現，往往p值是偶數時，資料可以做到近似線性可分，p是奇數時，資料分佈的形態也屬於另外一種固定模式，但是不再是線性可分。4. 程式碼
前面給出了自己對KPCA的理論解釋，以及做的一些基礎實驗，不給出實現程式碼，就不厚道了，程式碼如下所示，一部分是KPCA演算法程式碼，另一部分是實驗程式碼。function [eigenvalue, eigenvectors, project_invectors] = kpca(x, sigma, cls, target_dim)
% kpca進行資料提取的函式
psize=size(x);
m=psize(1); % 樣本數
n=psize(2); % 樣本維數
% 計算核矩陣k
l=ones(m,m);
for i=1:m
for j=1:m
k(i,j)=kernel(x(i,:),x(j,:),cls,sigma);
end
end
% 計算中心化後的核矩陣
kl=k-l*k/m-k*l/m+l*k*l/(m*m);
% 計算特徵值與特徵向量
[v,e] = eig(kl);
e = diag(e);
% 篩選特徵值與特徵向量
[dump, index] = sort(e, 'descend');
e = e(index);
v = v(:, index);
rank = 0;
for i = 1 : size(v, 2)
if e(i) < 1e-6
break;
else
v(:, i) = v(:, i) ./ sqrt(e(i));
end
rank = rank + 1;
end
eigenvectors = v(:, 1 : target_dim);
eigenvalue = e(1 : target_dim);
% 投影
project_invectors = kl*eigenvectors; %計算在特徵空間向量上的投影
endfunction compare
clear all;
close all;
clc;
% 生成非線性可分的三類資料
if exist('X1.mat')
load 'X1.mat'
load 'X2.mat'
load 'X3.mat'
figure(1)
plot(X1(1, :),X1(2, :) ,'ro')
hold on;
plot(X2(1, :),X2(2, :),'g*')
hold on;
plot(X3(1, :),X3(2, :),'b.')
hold on;
title('原始資料');
xlabel('第一維');
ylabel('第二維');
saveas(gcf, '原始資料圖.jpg')
else
[X1, X2, X3] = generate_data();
save 'X1.mat' X1
save 'X2.mat' X2
save 'X3.mat' X3
end
X = [X1 X2 X3];
[nFea, nSmps] = size(X);
nClsSmps = nSmps / 3;
% PCA
[vec_pca, Y_pca, value_pca] = princomp(X');
Y_pca = Y_pca';
figure(2);
plot(Y_pca(1, 1 : nClsSmps),Y_pca(2, 1 : nClsSmps), 'ro');
hold on;
plot(Y_pca(1, nClsSmps + 1 : 2 * nClsSmps),Y_pca(2, nClsSmps + 1 : 2 * nClsSmps), 'g*');
hold on;
plot(Y_pca(1, 2 * nClsSmps + 1 : end),Y_pca(2, 2 * nClsSmps + 1 : end), 'b.');
hold on;
title('PCA');
xlabel('第一維');
ylabel('第二維');
saveas(gcf, 'PCA投影圖.jpg')
% KPCA
percent = 1;
var = 2; % 1 代表高斯核，2代表多項式核，3代表線性核
sigma = 6; % 核引數
[vec_KPCA, value_KPCA, Y_pca] = kpca(X', sigma, var, 2);
Y_pca = Y_pca';
figure(3);
plot(Y_pca(1, 1 : nClsSmps),Y_pca(2, 1 : nClsSmps), 'ro');
hold on;
plot(Y_pca(1, nClsSmps + 1 : 2 * nClsSmps),Y_pca(2, nClsSmps + 1 : 2 * nClsSmps), 'g*');
hold on;
plot(Y_pca(1, 2 * nClsSmps + 1 : end),Y_pca(2, 2 * nClsSmps + 1 : end), 'b.');
hold on;
str = strcat('KPCA', '(p =', num2str(sigma), ')');
title(str);
xlabel('第一維');
ylabel('第二維');
str = strcat(str, '.jpg')
saveas(gcf, str)
end