% Example 8.3
clear;
x=[0 1]; y=[-1 1];
alpha=0.2;   %alpha=0.5;
wx=1; wh=1; bh=0; bz=0;
for i=1:50
 h=wx*x+bh;
 z=wh*h+bz;
 MSE(i)=sum((y-z).^2)/2;
 if MSE(i) < 1E-6, break, end
 dwh = -(y(1)-z(1))*h(1)    - (y(2)-z(2))*h(2);
 dbz = -(y(1)-z(1))         - (y(2)-z(2));
 wh = wh -alpha*dwh;
 bz = bz -alpha*dbz;
 dwx = -(y(1)-z(1))*wh*x(1) - (y(2)-z(2))*wh*x(2);
 dbh = -(y(1)-z(1))*wh      - (y(2)-z(2))*wh;
 wx = wx -alpha*dwx;
 bh = bh -alpha*dbh;
end
epoch=1:i; plot(epoch,MSE)
%-----------------------------------------------------------------------------
%
% Example 8.4
x=[0    5   10   15   20];           %[1 x 5] training input data
y=[1 0.99 0.99 0.94 0.95];          %[1 x 5] training output data
nh=1;                                       % num. of hidden node
net=feedforwardnet(nh);          % create two-layer network model
net=configure(net,x,y);                        % this is optional
[netModel,trainRecor]=train(net,x,y);     % train the model ‘net’
xNew=-5:0.1:25;                                % new input points
z=sim(netModel,xNew);                        % simulation results
bestE=trainRecor.best_epoch;
bestP=trainRecor.best_vperf; % or trainRecor.vperf(bestE+1)
hold on; plot(bestE,bestP,'o');
%-----------------------------------------------------------------------------
%
% Example 8.5
%Plotting the actual function at grid points
clear; n=51;
x=linspace(0.5,1.5,n); y=linspace(0.5,1.5,n);
[X,Y]=meshgrid(x,y);
Z=zeros(n);
for i=1:n, for j=1:n
 Z(i,j)=Rosenbrock([y(j) x(i)]);
end, end
figure(1); cs=surface(X,Y,Z); hold on;
%
%LHS samples
rng default;
xx=0.5+lhsdesign(20, 2,'iterations',5000,'smooth','off');
yy=Rosenbrock(xx);
plot3(xx(:,1), xx(:,2), yy,'or'); hold off
%
% FFNN
nh=5;                                       % num. of hidden node
net=feedforwardnet(nh);          % create two-layer network model
net=configure(net,xx',yy');                    % this is optional
[netModel,trainRecor]=train(net,xx',yy'); % train the model 'net'
ZNN=zeros(n);
for i=1:n, for j=1:n
 ZNN(i,j)=sim(netModel,[y(j) x(i)]');
end, end
figure(2); cs=surface(X,Y,ZNN); hold on;
plot3(xx(:,1), xx(:,2), yy,'or'); hold off
%
RMSE=sqrt(sum(sum(abs(Z-ZNN)))/(n*n))
%-----------------------------------------------------------------------------
%
% Figure 8.15
clear; close all; rng default
a=0; b=4*pi; n=51; ns=15; nh=5; 
x=linspace(a,b,n);                             % Test grid points
z=0.2*x+sin(x);                                   % True function
%
xx=linspace(a,b,ns);                     % Equally spaced samples
yy=0.2*xx+sin(xx);
%
figure(1); hold on
for i=1:20                                   % Repetition of FFNN
 net=feedforwardnet(nh);         % create two-layer network model
 net=configure(net,xx,yy);
 if i==1, wb=getwb(net); end             
 net = setwb(net,wb);                    % Fix initial parameters
 %
 net.divideFcn = 'divideind';                    % Fixed datasets
 net.divideParam.trainInd=[1 2 3 4 6 7 9 11 12 14 15];
 net.divideParam.valInd=[5 10];
 net.divideParam.testInd=[8 13];
 %
 [netModel,trainRecor]=train(net,xx,yy);  % train the model 'net'
 zNN=sim(netModel,x);                    % FFNN model predictions
 plot(x,zNN,'k');
 %
 RMSE(i)=sqrt((z-zNN)*(z-zNN)'/n);            % RMSE at test grid
end
plot(x,z,'r',xx, yy,'or'); hold off
m=mean(RMSE), s=std(RMSE)   % Mean and standard deviation of RMSE
%-----------------------------------------------------------------------------
%
% Example 8.7
clear; close all; rng default
a=0; b=4*pi; n=51; ns=15; nh=10; 
x=linspace(a,b,n);                             % Test grid points
z=0.2*x+sin(x);                                   % True function
%
xx=linspace(a,b,ns);                     % Equally spaced samples
v=randn(1,ns); v=0.3*(v-mean(v))/std(v); 
yy=0.2*xx+sin(xx)+v;                           % Add random noise
%
net=feedforwardnet(nh);          % create two-layer network model
net=configure(net,xx,yy);
all=1:ns; valInd=[]; testInd=[];
net.divideFcn = 'divideind';                    % Fixed datasets
net.divideParam.valInd=valInd;
net.divideParam.testInd=testInd;
net.divideParam.trainInd=setdiff(setdiff(all,valInd),testInd);
[netModel,trainRecor]=train(net,xx,yy);  % train the model 'net'
zNN=sim(netModel,x);                    % FFNN model predictions
%
plot(x,zNN,'k',x,z,'r',xx, yy,'or'); hold off
RMSE=sqrt((z-zNN)*(z-zNN)'/n)                % RMSE at test grid
%-----------------------------------------------------------------------------
%
% Example 8.8
clear; close all; rng default
a=0; b=4*pi; n=51; ns=15; nh=10; 
x=linspace(a,b,n);                             % Test grid points
z=0.2*x+sin(x);                                   % True function
%
xx=linspace(a,b,ns);                     % Equally spaced samples
yy=0.2*xx+sin(xx);
%
figure(1); hold on
for i=1:20                                   % Repetition of FFNN
 v=randn(1,ns); v=0.3*(v-mean(v))/std(v); 
 yy=0.2*xx+sin(xx)+v;                          % Add random noise
 %
 net=feedforwardnet(nh,'traincgf');% create two-layer network model
 net=configure(net,xx,yy);
 net.performParam.regularization=0.01;
 %
 all=1:ns; valInd=[]; testInd=[];
 net.divideFcn = 'divideind';                    % Fixed datasets
 net.divideParam.valInd=valInd;
 net.divideParam.testInd=testInd;
 net.divideParam.trainInd=setdiff(setdiff(all,valInd),testInd);
 %
 [netModel,trainRecor]=train(net,xx,yy);  % train the model 'net'
 zNN=sim(netModel,x);                    % FFNN model predictions
 plot(x,zNN,'k');
 %
 RMSE(i)=sqrt((z-zNN)*(z-zNN)'/n);            % RMSE at test grid
end
plot(x,z,'r',xx, yy,'or'); hold off
m=mean(RMSE), s=std(RMSE)   % Mean and standard deviation of RMSE
%-----------------------------------------------------------------------------
%
% Example 8.10
k=1.0; %k=0.5;
d=0:0.1:10;
for i=1:5
 z(i,:)=-(2*d-3*i).^2+2*i;
end
tmax=max(0,z);
smax=log(1+sum(exp(k*z)))/k;
plot(d,z(1,:),'k--',d,z(2,:),'k--',d,z(3,:),'k--',d,z(4,:),'k--'...
    ,d,z(5,:),'k--',d,tmax,'b',d,smax,'r')
axis([0 10 -1 10])
%-----------------------------------------------------------------------------