function demo_mlp_learning()
syms c   % adaptable parameters in an MLP network
M=20;d=2;N=800
% initialization 
c0=rand(1,M*(d+1)+M+1)*2-1;
size(c0)
% preparation of training data
% uniform sampling
% Substitute to the target function g
z=rand(N,2)*4*pi-2*pi; 
y=g_sin(z(:,1),z(:,2))+rand(N,1)*0.02-0.01;  
size(y)
plot3(z(:,1),z(:,2),y,'.');
h = g_hat2(z(:,1),z(:,2),c0,M);
% Levenberg-Marquardt method
options = optimoptions('fsolve','Algorithm', 'levenberg-marquardt')
% specification of a nonlinear system for MLP_learning
f = @(c)learning_mlp(c,z(:,1),z(:,2),y,M);
% apply fsolve
% Verification of c_zero
c_zero = fsolve(f,c0,options);
mean((learning_mlp(c_zero,z(:,1),z(:,2),y,M)).^2)
mean((g_hat2(z(:,1),z(:,2),c_zero,M)-y).^2)
% testing phase
z_test=rand(100,2)*2*pi-pi; 
y_test=g_sin(z_test(:,1),z_test(:,2));  
mean((g_hat2(z_test(:,1),z_test(:,2),c_zero,M)-y_test).^2)
end

% a nonlinear system for MLP_learning

function F = learning_mlp(c,x1,x2,y,M)
A = [x1 x2 ones(length(x1),1)];
d=2;W=[];
%W=reshape(c(1:M*(d+1)),d+1,M);
for j=1:M
   W=[W c((j-1)*(d+1)+1:j*(d+1))'];
end
v =tanh(A*W);
F = v*c(M*(d+1)+1:end-1)'+ones(length(x1),1)*c(end);
F =F-y;
end


function h = g_sin(x1,x2)
C1=[1 1/2 -1/2]'; %weight
A = [x1 x2 ones(length(x1),1)];
h = sin(A*C1); %activation function
end

function h = g_hat2(x1,x2,c,M)
A = [x1 x2 ones(length(x1),1)];
d=2;W=[];
%W=reshape(c(1:M*(d+1)),d+1,M);
for j=1:M
   W=[W c((j-1)*(d+1)+1:j*(d+1))'];
end
v =tanh(A*W);
h = v*c(M*(d+1)+1:end-1)';
h=h+ones(length(x1),1)*c(end);
end

