function?qn=kpca（dtrainkernelqFaCa）;
%?KPCA?-核主成分分析
%
%?使用?[trainFeatbj]=kpca（datakernelp1p2）
%?輸入?data?-?原始數據文件名
%??????kernel?-?核函數
%??????p1?p2?-?核函數參數
%??????q?-?投影后的維數
%?輸出????trainFeat?-?映射到特征空間的樣本

%訓練樣本
%?global?pp1FaCa
global?pp1
[m?n]=size（dtrain）;
k=zeros（mm）;
%?計算核矩陣
?for?i=1:m
???????for?j=1:m
??????????k（ij）?=?svkernel（kerneldtrain（i:）dtrain（j:））;
???????end
?end
D=mean（k1）;%列均????運行時注意查看一下D和C是否與原來的J和J’相等
%?C=mean（k2）;%行均
E=mean（D）;
J=repmat（Dm1）;
%?J‘
%?J1=repmat（D‘1m）;
%?J2=repmat（C‘m1）;%利用這種方法得到的結果和利用J和J‘得到的結果是一樣的
k=k-J-J‘+repmat（Emm）;?%將核矩陣中心化
[pcsevalue]=eigs（kq‘LM‘）;?%核矩陣的特征向量和特征值;這里的q只是粗選q個最大的特征值，這樣應該可以減小后面的計算量
%注意：使用eigs函數取的最大的q個特征值后，matlab將自動將這些特征值當做是一個q階矩陣產生的
%的特征值。為此，在下面計算貢獻率時，它將q個特征值的貢獻率當做為100%，將和實際的樣本貢獻率發生很大的差別
%?evoctor=pcs;%%%%%%%%%
%EIGS（AKSIGMA）?and?EIGS（ABKSIGMA）?return?K?eigenvalues?based?on?SIGMA:
???????%‘LM‘?or?‘SM‘?-?Largest?or?Smallest?Magnitude
%特征向量規范化
sqrtL=1./sqrt（diag（evalue））;
invsqrtL=diag（sqrtL）;
clear?sqrtL
%sqrtL=diag（sqrt（evalue））;
%invsqrtL=diag（1./diag（sqrtL））;%注意，這里最后得到的invsqrtL變量是以ai為元素的對角矩陣！！！！！！！！！！！！
%?pcs=invsqrtL*pcs‘;%這里是是根據公式將a的模值等于特征值的根號分之一
pcs=invsqrtL*pcs‘;


value=sum（evalue）;
eev=value;
Tnew=pcs*k‘;%映射到主成分空間；取對應主元個數個ai系數
Tnew=Tnew‘;%%測試樣本經過投影映射后得到的主元矩陣%Tnew是m1*pnum維的
es=inv（diag（eev））;
for?i=1:m
????T2（i）=Tnew（i:）*es*Tnew（i:）‘;
end
c1=ones（m1）;
ep=（1/（m^2））*sum（sum（k））;
epp=repmat（epm1）;
SPE=c1-（2/m）.*sum（k2）+epp-diag（Tnew*Tnew‘）;%對于每一個新樣本的SPE值是按照《故障診斷與容錯控制》PPT上140頁的公式2.64計算的，
tts=length（SPE）;
ttt=1:1:tts;
figure（2）
plot（tttSPE‘g-o‘）;
num=ceil（0.08*m）;%去掉總數的10%，將其當作誤差
sortts1=sort（T2‘descend‘）;
sortts2=sort（SPE‘descend‘）;
ddt=sortts2（1:num）;
tt=length（ddt）;
qn=[];
for?i=1:tt
????qn=[qnfind（ddt（i）==SPE）];
end

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----
?????文件????????2201??2009-03-21?10:05??kKKKKKPCA.m