yang_rbf.m

一种最小二乘法的RBF神经网络

clear all
close all

global sp k w2 b2 unitcenters

alldatanum=400;
traindatanum=300;
testdatanum=100;
sp=8;
errorlimit=0.1;

u=rand(1,alldatanum+10);
y=zeros(1,alldatanum+10);
for i=4:alldatanum+3
    y(1,i)=1.5*y(1,i-1)*y(1,i-2)*y(1,i-3)/(3+y(1,i-1)+y(1,i-2)*y(1,i-2)+...
             y(1,i-3)*y(1,i-3)*y(1,i-3))+u(1,i)+u(1,i-1)+u(1,i-2);
end

alldatain=[];
alldataout=[];
for i=4:alldatanum+3
    newin=[y(1,i-1);y(i-2);y(1,i-3);u(1,i);u(1,i-1);u(1,i-2)];
    alldatain=[alldatain newin];
    alldataout=[alldataout y(1,i)];
end

traindatain=alldatain(:,1:traindatanum);
traindataout=alldataout(:,1:traindatanum);
testdatain=alldatain(:,traindatanum+1:alldatanum);
testdataout=alldataout(:,traindatanum+1:alldatanum);

figure
hold on
grid
j=traindatanum+1:alldatanum
plot(j,testdatain,'k-')
xlabel('input t');
ylabel('output u');

figure
hold on
grid
j=traindatanum+1:alldatanum
plot(j,testdataout,'b-')
xlabel('input t');
ylabel('output y');

[indim,maxunitnum]=size(traindatain);

distance=dist(traindatain',traindatain);
hiddenunitout=radbas(distance/sp);

posselected=[];
vectorsselected=[];
hiddenunitoutselected=[];
errorhistory=[];

vectorsselectfrom=hiddenunitout;

dd=sum((traindataout.*traindataout)')';
for k=1:maxunitnum
    pp=sum(vectorsselectfrom.*vectorsselectfrom)';
    denominator=dd*pp';
    [xxx,selectednum]=size(posselected);
    if selectednum>0,
        [lin,xxx]=size(denominator);
        denominator(:,posselected)=ones(lin,1);
    end
    angle=((traindataout*vectorsselectfrom).^2)./denominator;
    
    [value,pos]=max(angle);
    posselected=[posselected pos];
    
    hiddenunitoutselected=[hiddenunitoutselected;hiddenunitout(pos,:)];
    hiddenunitoutex=[hiddenunitoutselected;ones(1,traindatanum)];
    w2ex=traindataout*pinv(hiddenunitoutex);
    w2=w2ex(:,1:k);
    b2=w2ex(:,k+1);
    nnout=w2*hiddenunitoutselected+b2;
    SSE=sumsqr(traindataout-nnout)
    
    errorhistory=[errorhistory SSE];
    
    if SSE<errorlimit,break,end
    
   newvector=vectorsselectfrom(:,pos);
   projectionlen=newvector'*vectorsselectfrom/(newvector'*newvector);
   vectorsselectfrom=vectorsselectfrom-newvector*projectionlen;
end

unitcenters=traindatain(:,posselected);

testdistance=dist(unitcenters',testdatain);
testhiddenunitout=radbas(testdistance/sp);
testnnout=w2*testhiddenunitout+b2;
plot(j,testnnout,'r-')

k
unitcenters
w2
b2

figure
hold on
grid
[xx,num]=size(errorhistory);
plot(1:num,errorhistory,'r-');