求解这段指令的意思，越详细越好，谢谢啦
 

function [Kp,T2]=KPCA(ax,ay)
[Nx]=size(ax);
mean_X =mean(ax);
axb=ax;
std_X=std(ax);
ax=ax-mean_X(ones(Nx,1),:);
std_X(find(std_X==0))=1;%数据预处理
ax=ax./std_X(ones(Nx,1),:);
c=10000;
% gama=0.05;
% ni=1;
% F1=ax(1,:);
% F=F1';
% for ib=2:Nx
% for i=1:ni
% for j=1:ni
%
% batar1(ib).block(i,j)=exp(-norm(ax(i,:)-ax(j,:))^2/c);
% end
% batar2(i,1)=exp(-norm(ax(i,:)-ax(ib,:))^2/c);
% batar3(1,i)=exp(-norm(ax(ib,:)-ax(i,:))^2/c);
% end
% s1=exp(-norm(ax(ib,:)-ax(ib,:))^2/c);
% batar= batar3(1,i)*inv(batar1(ib).block)* batar2(i,1);
% s=(s1- batar)/s1;
% if s>sqrt(gama)
% ni=ni+1;
% F=[F ax(ib,:)'];
% end
%
%
% end
% AX=F'%训练集基的提取结束
[N]=size(ax,1);
for i=1:N
for j=1:N
K(i,j)=exp(-norm(ax(i,:)-ax(j,:))^2/c);%求核矩阵
end
end
n1=ones(N,N);
N1=1/N*n1;
Kp=K-N1*K-K*N1+N1*K*N1;
[u,s,v]=svd(Kp);
lamda=s*(1/(N-1));
P=v;
lamda=diag(lamda);
B=length(find(lamda>1e-10)); %求非零的特征值个数
%求主元个数
npc=1;
while sum(lamda(1:npc))/sum(lamda(1:B))<0.9
npc=npc+1;
end
%求特征空间有效维数
DimFS=1;
while sum(lamda(1:DimFS))/sum(lamda(1:B))<=0.99
DimFS=DimFS+1;
end
lamda=diag(lamda);
for i=1:B

% P(:,i)=P(:,i)/norm(P(:,i)*s(i,i));
% P(:,i)=P(:,i)/(norm(P(:,i))*sqrt(N*lamda(i,i)));
P(:,i)=P(:,i)/(sqrt((N-1)*lamda(i,i)));
end
[Ny]=size(ay,1);
mean_X =mean(axb);
std_X = std(axb);
[num_sample] = Ny;
ay = ay-mean_X(ones(num_sample,1),:);
ay = ay./std_X(ones(num_sample,1),:);
% mean_y = mean(ay);
% std_y=std(ay);
% ay = ay-mean_y(ones(Ny,1),:);
% std_y(find(std_y==0))=1;%数据处理
% ay = ay./std_y(ones(Ny,1),:);

for i=1:Ny
for j=1:N
Ky(i,j)=exp(-norm(ay(i,:)-ax(j,:))^2/c);
end
end
t1=ones(1,N);
t11=1/N*t1;
for i=1:Ny
kp1(i,:)= Ky(i,:)-t11*K- Ky(i,:)*N1+t11*K*N1;
end
for i=1:Ny
for k=1:B
t(i,k)=P(:,k)'*kp1(i,:)';
end
end
% 求T2,SPE
% covtyb=inv(t'*t);
for i=1:Ny
T2(i)=t(i,1:npc)*inv(lamda(1:npc,1:npc))*t(i,1:npc)'; %也可以
% SPE(i)=t(i,1:npc)*t(i,1:npc)';
% T2(1,i)=t(i,1:npc)*(covtyb(1:npc,1:npc))*t(i,1:npc)';

% SPE(i)=t(i,(npc+1):B)*t(i,(npc+1):B)';
SPE(i)=1-(2/N)*ones(N,1)'*Ky(i,:)'+(1/(N^2))*ones(N,1)'*K*ones(N,1)-t(i,1:npc)*t(i,1:npc)';
end
%T2,SPE控制线
t2cl=npc*(N-1)*(N+1)*icdf('f',0.99,npc,N-npc)/(N*(N-npc));
for i=1:3
theta(i)=trace((lamda(npc+1:DimFS,npc+1:DimFS))^i);
end
h0=1-2*theta(1)*theta(3)/(3*theta(2)^2);
ca=icdf('norm',0.99,0,1);
s_cl=theta(1)*(ca*sqrt(2*theta(2)*h0^2)/theta(1)+1+theta(2)*h0*(h0-1)/theta(1)^2)^(1/h0)*4;
for i=1:Ny
t2cl1(i)=t2cl;
end
for i=1:Ny
s_cl1(i)=s_cl;
end
figure
plot(1:Ny,T2(1:Ny),'k');
hold on;
plot(1:Ny,t2cl1(1:Ny),'r');
title('T2');
hold off;
figure
plot(1:Ny,SPE(1:Ny),'k')
hold on;
plot(1:Ny, s_cl1(1:Ny),'r');
title('SPE');
hold off;