% x(N,M)  - input- M input feature vectors with size Nx1
% alpha   - input- step size
%
clear all
close all
clc;
M=100;  % # of equations constraint
N=5;    % # of unknowns (w0, w1, wN-1)
x=zeros(N,M);

for i=1:M
    x(:,i)=round(rand(N,1));
end
% -- to balance the number of 0 examples and 1 exmaples --
is0_token=0;
for i=1:M
    if sum(x(:,i))==0
        is0_token=is0_token+1;
    end
end
x=[zeros(N,M-is0_token),x];
M=size(x,2);
rank = randperm(M);
x=x(:,rank);
% --------------------------------------------------------------------------------
t=zeros(M,1);tp=t;
for i=1:M
    if sum(x(:,i))==0
        t(i)=1;
    end
end

x=[ones(1,M);x];  N=N+1;        % add 1st row of x with 1's
%% ----------nn node-------------------------
% Define the nodes in each layer (do not include bias)
% N is the feature number of input, 1 is the label number (0 or 1)
layer_temp=input('please input the number of nodes in each hidden layer (not include the input feathuires and output labels) (e.g:[15,10,5]): \n ');
layer_size=[N,layer_temp,1]; %layer_size=[N,1] for one node, layer_size=[N,10,1] for two node
layer_num=numel(layer_size); %layer number include the input and output
for ilayer=1:layer_num-1
    ww{ilayer}=rand(layer_size(ilayer+1),layer_size(ilayer));
end
alpha=0.01;         % step size
lim=.00001;         % stopping criterion
nit=5000;res=zeros(nit,1);
%if set more layers and more node, seems fall into local minima, need to
%increase nit 2~4 order of magnitude
kk=0;
for k=2:nit     % Looping over iterations
    [grad,res(k)]=gradientnn(M,x,t',ww,layer_size);
    for ilayer=1:layer_num-1
        ww{ilayer}=ww{ilayer}-alpha*grad{ilayer};
    end
    kk=kk+1
end
% if not converge, do more iteration
while (res(nit)/res(2))>0.05
    res=[res;zeros(nit,1)];
     nit=nit*2;
   
    for k=nit/2:nit     % Looping over iterations
        [grad,res(k)]=gradientnn(M,x,t',ww,layer_size);
        for ilayer=1:layer_num-1
            ww{ilayer}=ww{ilayer}-alpha*grad{ilayer};
        end
        kk=kk+1
    end
end
figure(3);Displaynn(M,res,nit,ww,layer_size,x,t);pause(0.5)

%% ------ single node------------------------------------
layer_size=[N,1]; 
alpha=.01;         % step size
w=rand(N,1);        % starting model
lim=.00001;         % stopping criterion
nit=5000;res=zeros(nit,1);
kk=0; beta=1;ioption=0;
for k=2:nit     % Looping over iterations
    [grad,res]=gradient(M,N,x,t,w,beta,res,k,ioption);
    w=w-alpha*grad;  % Gradient descent update
    if res(k)<=lim;k=nit;end
    kk=kk+1
end
figure(1);Display(M,res,nit,w,layer_size,x,t,beta);pause(0.5)
%% ----------two node-------------------------
alpha=.01;         % step size
N1=10;
layer_size=[N,10,1]; 
w1=rand(N,N1);
w2=rand(N1,1);        % starting model
lim=.00001;         % stopping criterion
nit=5000;res=zeros(nit,1);
kk=0; beta=1;ioption=0;
for k=2:nit     % Looping over iterations
    [grad1,grad2,res]=gradient2(M,x,t,w1,w2,k,res,ioption);
    w1=w1-alpha*grad1;  % Gradient descent update
    w2=w2-alpha*grad2;  % Gradient descent update
    if res(k)<=lim;k=nit;end
    kk=kk+1
end
figure(2);Display2(M,res,nit,w1,w2,layer_size,x,t);pause(0.5)