gui_ssa_forecasting.m

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

function gui_ssa_forecasting(groups_No);

% last modified 11.02.05

global SSA GSD_GLOBALS

grN=length(SSA.groups);
str={};
ylabel_str={};
curr_g=1;
for i=1:grN
    if ~isempty(SSA.groups{i})
        str{curr_g}=['group ' num2str(curr_g) ': '];
        ylabel_str{curr_g}=['group ' num2str(curr_g)];
        for j=1:length(SSA.groups{i})
            str{curr_g}=[str{curr_g} ' ' num2str(SSA.groups{i}(j))];
        end
        curr_g=curr_g+1;
    end
end

if nargin<1
	[selection,ok] = listdlg('PromptString','Select group:',...
                    'SelectionMode','multiple','ListString',str);
else
    ok=1;
    selection=groups_No;
end


if ok
	full_base_data=SSA.base_data;
    
    forec_comp=[];
    base_data=0;
    for i=1:length(selection)
        for j=1:length(SSA.groups{selection(i)})
            forec_comp=[forec_comp SSA.groups{selection(i)}(j)];
            base_data=base_data+ssa_get_elementary_vector(SSA.singular_numbers(SSA.groups{selection(i)}(j)),SSA.U(:,SSA.groups{selection(i)}(j)),SSA.V(:,SSA.groups{selection(i)}(j)),SSA.delay);
        end
    end
    
	prompt = {'Enter number of steps'};
	dlg_title = 'Input for SSA forecasting';
	num_lines= 1;
	def     = {'20'};
	answer  = inputdlg(prompt,dlg_title,num_lines,def);
	if isempty(answer)
        return
	end
	step_No=str2num(answer{1});
    f=full_base_data;
	for i=1:step_No
        f=ssa_my_forecasting(f,forec_comp);
	end
    
	tempHandleF=figure('Units','characters','Name',['SSA forecasting. Significant components: ' num2str(forec_comp)],...
                'NumberTitle','off','color',[1 1 1]);
    new_time=1:length([base_data; f(end-step_No+1:end)]);
    
	plot(new_time, [base_data; f(end-step_No+1:end)]);
	hold on
	plot(1:length(full_base_data), full_base_data,'r');
    xlim([new_time(1) new_time(end)]);
	GSD_GLOBALS.fig_handles=[GSD_GLOBALS.fig_handles tempHandleF];
end



function forecasted=ssa_my_forecasting(f,r);

global SSA

L=length(SSA.singular_numbers);
delay=SSA.delay;
N=length(f);
K=N-delay*(L-1);  %% K must be positive!!!

X=zeros(L,K);
for i=1:L
    X(i,:)=f(1+delay*(i-1):delay*(i-1)+K)';
end

C=X*X'/K;

[V,D] = eig(C);
D=diag(D);
[D,IX]=sort(-D);
V=V(:,IX);

v=V(end,r);
V_star=V(1:end-1,r);

Q=f(end-L+2:end);
Q=Q(:);

novoe=v*inv(V_star'*V_star)*V_star'*Q;
forecasted=[f; novoe];