da_algorithm.m

主成分分析和偏最小二乘SquaresPrincipal成分分析（ PCA ）和偏最小二乘（ PLS ）

function da=da_algorithm(modelinfo);
% Main computational algorithm for discr. analysis

% last modified 4.02.05

da=[];
win1=modelinfo.win1;
win2=modelinfo.win2;
coeffs=modelinfo.coeffs;
funcs=modelinfo.funcs;
X=modelinfo.X;
dTime=X.time;
crd_number=length(X.name);
for i=1:crd_number
    assignin('caller',sprintf('x%d',i),X.data(:,i));
end
R=zeros(length(X.data(:,1)),length(funcs));
for i=1:length(funcs)
    R(:,i)=evalin('caller',funcs{i});
end
for i=1:crd_number
    evalin('caller',['clear ' sprintf('x%d',i)]);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if ~strcmp(modelinfo.mode,'auto')
	handles.fig=figure('name','Processing...','numbertitle','off','color',[1 1 1]);
	handles.crd_ax=subplot(311);
     plot(X.time,X.data(:,1));
	handles.fisher_ax=subplot(312);
	title('Fisher criterion');
	grid on;
	handles.df_ax=subplot(313);
	title('Discriminant function  (-c0)');
	grid on;
	set([handles.crd_ax handles.fisher_ax handles.df_ax],'xlim',[min(dTime) max(dTime)]);
	set([handles.crd_ax handles.fisher_ax handles.df_ax],'nextplot','add');
	set([handles.crd_ax handles.fisher_ax handles.df_ax],'fontsize',8);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
data_len=length(X.data(:,1));
Fisher_length=data_len-win2-win1+1;
Fisher=zeros(Fisher_length,1);
Discr_func=Fisher(:,1);
Coeff_data=zeros(Fisher_length,length(funcs));
fisher_max=1;
df_max=1;
df_min=0;

k=1;
if strcmp(modelinfo.mode,'manual')
	for i=win1:data_len-win2
        draw_windows(win1,win2,i,dTime,handles);
        drawnow;
        W1=R(i-win1+1:i,:);
        W2=R(i+1:i+win2,:);
        [h,resh]=da2(W1,W2);
        Fisher(k)=h;
        Coeff_data(k,:)=resh;
        Discr_func(k)=sum(R(i,:).*resh);
        if k>1
            set(0,'currentfigure',handles.fig);
            set(handles.fig,'currentaxes',handles.fisher_ax);
		    if max(h)>fisher_max
		        ylim([0 1.1*max(h)]);
		        fisher_max=1.1*max(h);
		    end
            plot(dTime(i-1:i),Fisher(k-1:k,1),'erasemode','none','LineWidth',2);
            set(0,'currentfigure',handles.fig);
            set(handles.fig,'currentaxes',handles.df_ax);
		    if Discr_func(k)>df_max
		        ylim([df_min Discr_func(k)+0.1*abs(Discr_func(k))]);
		        df_max=Discr_func(k)+0.1*abs(Discr_func(k));
		    end
		    if Discr_func(k)<df_min
		        ylim([Discr_func(k)-0.1*abs(Discr_func(k)) df_max]);
		        df_min=Discr_func(k)-0.1*abs(Discr_func(k));
		    end
            plot(dTime(i-1:i),Discr_func(k-1:k),'erasemode','none');
        end
        k=k+1;
	end
else
	for i=win1:data_len-win2
        W1=R(i-win1+1:i,:);
        W2=R(i+1:i+win2,:);
        [h,resh]=da2(W1,W2);
        Fisher(k,:)=h;
        Coeff_data(k,:)=resh;
        Discr_func(k)=sum(R(i,:).*resh);
        k=k+1;
	end
end    
time=dTime(win1:data_len-win2);
IX=[];
for i=1:Fisher_length
    if (~all(isfinite(Coeff_data(i,:))))|(~isfinite(Discr_func(i)))
        IX=[IX i];
        fprintf('infinity solution removed at time %d\n',time(i));
    end
end

if ~isempty(IX)
    Coeff_data(IX,:)=[];
    Fisher(IX,:)=[];
    time(IX)=[];
    Discr_func(IX)=[];
end
da.df_txt=modelinfo.df;
da.win1=win1;
da.win2=win2;
da.functions_txt=modelinfo.funcs;
da.coefficients_txt=modelinfo.coeffs;
da.time=time;
da.fisher=Fisher; clear Fisher;
da.coeff_data=Coeff_data; clear Coeff_data;
da.discriminant_func=Discr_func; clear Discr_func;
da.x=X;

if ~strcmp(modelinfo.mode,'auto')
	set(handles.fig,'name','Processing...Done');
	h = findobj(handles.crd_ax,'Type','line','color','red');
	delete(h);
end



function draw_windows(win1,win2,common_point,X,handles);
axes(handles.crd_ax);
h = findobj(handles.crd_ax,'Type','line','color','red');
delete(h);

line1_x=[X(common_point-win1+1) X(common_point)];
line2_x=[X(common_point+1) X(common_point+win2)];
Ylim=get(handles.crd_ax,'ylim');
line(line1_x,Ylim(1)*ones(2,1),'Color','r','LineWidth',2,'erasemode','xor');
line(line1_x,Ylim(2)*ones(2,1),'Color','r','LineWidth',2,'erasemode','xor');
line([line1_x(1) line1_x(1)],Ylim,'Color','r','LineWidth',2,'erasemode','xor');
line([line1_x(2) line1_x(2)],Ylim,'Color','r','LineWidth',2,'erasemode','xor');
line(line2_x,Ylim(1)*ones(2,1),'Color','r','LineWidth',2,'erasemode','xor');
line(line2_x,Ylim(2)*ones(2,1),'Color','r','LineWidth',2,'erasemode','xor');
line([line2_x(1) line2_x(1)],Ylim,'Color','r','LineWidth',2,'erasemode','xor');
line([line2_x(2) line2_x(2)],Ylim,'Color','r','LineWidth',2,'erasemode','xor');


function [h,alf] = da2(dat1,dat2);

m1=mean(dat1);
m2=mean(dat2);
m3=length(dat1(:,1))/(length(dat1(:,1))+length(dat2(:,1)))*m1+length(dat2(:,1))/(length(dat1(:,1))+length(dat2(:,1)))*m2;
W=length(dat1(:,1))/(length(dat1(:,1))+length(dat2(:,1)))*cov(dat1)+length(dat2(:,1))/(length(dat1(:,1))+length(dat2(:,1)))*cov(dat2);
B=(m1-m3)'*(m1-m3)+(m2-m3)'*(m2-m3);

[V,D]=eig(B,W,'chol');
[h,I] = max(diag(D));
alf = V(:,I)';