%功能:演示回归算法算法在计算机视觉中的应用
%实现如何利用偏最小二乘回归模型实现数据拟合;
%环境:Win7,Matlab2012b
%Modi: NUDT-VAP30
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all;
close all;
clc;
%多重相关性计算
load('data.mat'); %读入数据文件data.mat
cr=corrcoef(data) %计算变量之间的相关系数
%建立偏最小二乘回归模型
%提取可能的主成分
X=data(:,1:3);
Y=data(:,3:6);
E0=stand(X)
F0=stand(Y)
A=rank(E0)
[W,C,T,U,P,R]=plspcr(E0,F0)
%主成分解释能力分析
%计算主成分累计复测定系数
RA=plsra(T,R,F0,A)
%计算主成分的信息解释能力
[Rdx,RdX,RdXt,Rdy,RdY,RdYt]=plsrd(E0,F0,T,A)
%计算第一主成分间的相关性
%通过t1/u1图像直观的考查第一主成分间的相关性
cr=plsutcor(U,T)
%计算PLS回归方程的系数
%标准化因变量关于主成分t1的经验回归系数
TCOEFF=R(:,1)
%标准化因变量关于标准化自变量的经验回归系数
SCOEFF=pls(1,5,W,P,R)%1为建模的主成分个数,5为自变量个数
%计算原始因变量关于原始自变量的经验回归系数
[COEFF,INTERCEP]=plsiscoeff(X,Y,SCOEFF)%对标准化的回归系数进行逆标准化处理,输出原始自变量对因变量的回归系数及常数项
%变量投影重要性分析与模型的改进
result=plsresult(W,RdY,RdYt,1)%result表示第j个自变量对因变量的解释能力模
************************************************************************************************************************************
function SCOEFF = pls(h,p,W,P,R)
%PLS求偏最小二乘法回归方程的系数
%SCOEFF = pls(h,p,W,P,R)
%h-用于建模的主成分个数
%p-自变量个数
%W-模型效应权重p×rankE0矩阵
%P-模型效应载荷量p×rankE0矩阵
%R-因变量载荷量q×rankE0矩阵
%SCOEFF--偏最小二乘法回归方程的系数p×q矩阵
for byk=1:h
if byk==1
WX(:,byk)=W(:,byk);
SCOEFF=WX(:,byk)*R(:,byk)';
else
I=eye(p);
ww=eye(p);
for bykbyk=1:byk-1
ww=ww*(1-W(:,bykbyk)*P(:,bykbyk)');
end
WX(:,byk)=ww*W(:,byk);
end
SCOEEF=WX(:,byk)*R(:,byk)';
end
**************************************************************************************、
function [COEFF,INTERCEP] = plsiscoeff(X,Y,B)
%PLESTSCOEFF标准化回归系数逆标准化处理,输出原始自变量对因变量的回归系数及常数项
% [COEFF,INTERCEP]=plsiscoeff(X,Y,B)
%X-原始自变量数据
%Y-原始因变量数据
%B-标准化变量回归方程的系数
%COEFF-原始变量回归方程的系数
%INTERCEP-原始变量回归方程的常数项
[xrow,xcol]=size(X);
[yrow,ycol]=size(Y);
for i=1:ycol
bykCOEFF(:,i)=B(:,i).*std(Y(:,i));
end
for j=1:ycol
for i=1:xcol
COEFF(i,j)=bykCOEFF(i,j)/std(X(:,i));
end
end
INTERCEP=mean(Y)-(mean(X)*COEFF);
******************************************************************************
function [W,C,T,U,P,R]=plspcr(E0,F0)
%PLSPCR提取PLS建模过程所有可能的主成分
%[W,C,T,U,P,R]=plspcr(E0,F0)
%E0-自变量标准化的样本数据n×p矩阵;
%F0-因变量标准化的样本数据n×p矩阵;
%W -模拟效应权重p×rankE0矩阵
%C -因变量权重q×rankE0矩阵
%T -自变量系统主成分得分n×rankE0矩阵
%U -因变量系统主成分得分n×rankE0矩阵
%P -模型效应载荷量p×rankE0矩阵
%R -因变量载荷量q×rankE0矩阵
A=rank(E0);
W=[];
C=[];
T=[];
U=[];
P=[];
R=[];
for byk=1:A
%提取主轴与主成分
EFFE=E0'*F0*F0'*E0;
FEEF=F0'*E0*E0'*F0;
[w,LAMBDA]=eigs(EFFE,1,'lm')
[c,LAMBDA]=eigs(FEEF,1,'lm')
t1=E0*w;
u1=F0*c;
W=[W,w];
C=[C,c];
T=[T,t1];
U=[U,u1];
%计算残差
p1=(E0'*t1)/norm(t1)^2;
E1=E0-t1*p1';
E0=E1;
r1=(F0'*t1)/norm(t1)^2;
F1=F0-t1*r1';
F0=F1;
P=[P,p1];
R=[R,r1];
end
****************************************************************************
function RA = plsra(T,R,F0,rankE0)
%PLSRA求出主成分的累积复测定系数
%RA=plsra(T,R,F0,rankE0)
%T—自变量系统主成分得分N×rankE0矩阵
%R—因变量载荷量q×rankE0矩阵
%F0—因变量标准化的样本数据n×q矩阵
%rankE0-plspcr提取的主成分个数
RAAA=[];
for byk=1:rankE0
RAbyk=sum(norm(T(:,byk)).^2*norm(R(:,byk)).^2)./(norm(F0))^2;
RAAA=[RAAA,RAbyk];
end
RA=cumsum(RAAA);
***********************************************************************************
function [Rdx,RdX,RdXt,Rdy,RdY,RdYt] = plsrd(E0,F0,T,h)
%E0-标准化后的自变量数据
%F0-标准化后的因变量数据
%T-自变量系统主成分得分n×rankE0矩阵
%h-用于建模或希望进行解释能力分析的主成分个数
%Rdx-各主成分对于某自变量的解释能力
%RdX-各主成分对于自变量组的解释能力
%RdXt-全部主成分对自变量组的累积解释能力
%Rdy-各主成分对于某因变量的解释能力
%RdY-各主成分对因变量组的解释能力
%RdYt-全部主成分对因变量组的累积解释能力
%--成分对自变量解释能力分析--
[nr,nx]=size(E0);
[nr,ny]=size(F0);
Rdx=zeros(nx,h);
t1=zeros(nr,1);
x1=zeros(nr,1);
for xj=1:nx
for ti=1:h
t1=T(:,ti);
x1=E0(:,xj);
cc=(corrcoef(t1,x1)).^2;
Rdx(xj,ti)=cc(1,2);
end
end
Rdx;
RdX=sum(Rdx)./nx;
RdXt=sum(RdX);
%--成分对因变量解释能力的分析--
Rdy=zeros(ny,h);
y1=zeros(nr,1);
for yj=1:ny
for ti=1:h
t1=T(:,ti);
y1=F0(:,yj);
rr=(corrcoef(t1,y1)).^2;
Rdy(yj,ti)=rr(1,2);
end
end
Rdy;
RdY=sum(Rdy)./ny;
RdYt=sum(RdY);
******************************************************************************************
function result = plsresult( W,RdY,RdYt,h )
[nx,wk]=size(W);
result=zeros(1,nx);
for j=1:nx
for hh=1:h
Whj=W(j,hh);
tvip(hh)=RdY(hh)*Whj.^2;
end
S_tvip=sum(tvip);
result(j)=sqrt((nx/RdYt)*S_tvip);
end
bar(result,'c');
title('变量投影重要性result图')
*****************************************************************
function cr=plsutcor(U,T)
%PLSUTCOR绘制t1/u1图,并给出二者的相关系数
%cr=plsutcor(U,T)
%U-因变量提取的成分
%T-自变量提取的成分
%cr-自变量与因变量的相关系数
u1=U(:,1);
t1=T(:,1);
ut=[u1,t1];
cr=corrcoef(ut)
plot(t1,u1,'o')
lsline
title('t1/u1图')
xlabel('t1')
ylabel('u1')
***********************************************************************************
function X0 = stand(X)
% STAND将数据矩阵X逐列进行标准化处理,输出标准化数据X0
% X是原始数据矩阵,X0是标准化后的数据矩阵
zeros(size(X));
[nr,nx]=size(X);
for mk=1:nr
X0(mk,:)=(X(mk,:)-mean(X))./std(X);
end
**************************************************************