grbf.m

竞争学习的matlab工具箱

% MATLAB implementation of global SOM-based RBF time series predictor.
%
% Based on the prediction method described in the following paper:
%
% Barreto, G.A. and Souza, L.G.M. (2006). 
% "Adaptive filtering with the self-organizing map: A performance
% comparison", vol. 19, pp. 785-798.
%
% Author: Guilherme A. Barreto
% Date: September 28st 2006

clear; clc; close all;

%--------------- Organize the training data ------------
load train_series_data.dat;    % Load the time series to be clustered

%% Building the input vectors from an univariate time series
p=5;       % Dimension of the input vector (length of the time window)
lap=p-1;      % Amount of overlapping between consecutive input vectors
Dw=buffer(train_series_data,p,lap);   % Build the data vectors
if lap>0,
     Dw=Dw(:,p:end)';  % Eliminate the first 'p-1' vectors with zeros)
else Dw=Dw';
end
Dw=fliplr(Dw);

Ynext=Dw(2:end,1);        % Target (one-step-ahead) values 
Dw=[Ynext Dw(1:end-1,:)]; % Input regressor vectors for Ynext

Dw=Dw+0.01*randn(size(Dw));  % Add some gaussian noise to the data
%-------------------------------------------------------
% SOM initialization and training
Mx = 16;             % Number of neurons in the X-dimension
My = 1;              % Number of neurons in the Y-dimension
msize = [Mx My];     % Size of 2-D SOM map
MASK=[0; ones(p,1)];   
sMap = som_randinit(Dw, 'msize', msize,'lattice','rect','shape','sheet');
sMap = som_seqtrain(sMap,Dw,'radius',[floor(0.5*max(msize)) 0],'sample_order','random','neigh','gaussian','trainlen',50);

%--------------- Organize the testing data ------------
load test_series_data.dat;    % Load testing time series
Dw=buffer(test_series_data,p,lap); % Build input data vectors

if lap>0,
     Dw=Dw(:,p:end)';  % get rid of the first 'p-1' vectors with zeros
else Dw=Dw';
end
Dw=fliplr(Dw);

Ynext=Dw(2:end,1);        % Target (one-step-ahead) values 
Dw=[Ynext Dw(1:end-1,:)]; % Input regressor vectors for Ynext

Dw=Dw+0.01*randn(size(Dw));  % Add some gaussian noise to the data
%------------------------------------------------------
[LEN_DATA DIM_INPUT]=size(Dw);  % Data matrix size (1 input vector per row)

G=[]; w=[];
s2=0.10;  % Spread of the basis functions
for t=1:LEN_DATA,
    for i=1:Mx*My,
        act = norm(Dw(t,2:end)-sMap.codebook(i,2:end));
        G(i) = exp(-act*act/(2*s2*s2));   % output of the i-th basis function
        w(i) = sMap.codebook(i,1);        % hidden-to-output layer weight
    end
    Yhat(t) = dot(w,G)/sum(G);  % predicted value
    error(t)=Ynext(t)-Yhat(t);  % prediction error
end


% Plot target and predicted time series
figure; plot(Ynext,'r-'); hold on; plot(Yhat,'b-'); hold off

% Normalized mean squared error
NMSE=var(error)/var(Ynext)