Neural Network Models

Single Layer Perceptron

Here we develop the training code for the single layer perceptron.

Setup

The sigmoid function we will be using, with a = 10 gives Figure 14.1:

g[p_] := 1/(1 + Exp[-a p])

a = 10;

Plot[g[p] /. a-> 10,{p,-1,1},PlotStyle->GrayLevel[0.5]]

The function TrainOnce takes a matrix A, representing the current "state" of the fitting model, and a pair of associated data, xin and yin, a steepest descent step size, eta, and the tolerance eps, and iterates using (14.25), minimizing the residual.  One could normalize the step by making deltaA a unit vector, that would give a better interpretation for eta.

TrainOnce[Amat_, xin_, yin_,eta_, eps_] := Module[{z, Aout, fval,deltaA,MM,NN},
NN = Length[xin];
MM = Length[yin];
Aout = Amat;
z = g[Aout.xin];
fval = (z - yin).(z - yin);
While[fval ≥ eps,
deltaA = Table[2 (z[[i]] - yin[[i]]) z[[i]] (1 - z[[i]]) a xin[[j]],{i,1,MM},{j,1,NN}];
Aout = Aout - eta deltaA;
z = g[Aout.xin];
fval = (z - yin).(z - yin);
];
Return[Aout];
]

As a test, here is a single x-y pairing, and starting matrix with random entries - we use a special value for a, the parameter appearing in the sigmoid, to probe the nonlinear region of that function:

xin = {.1,.2,.3,.4,.5};
yin = {.1,.2,.3};
a = 10/(Length[xin] Max[xin]);
SeedRandom[1234];
Amat = Table[Random[Real,{-1,1}],{i,1,Length[yin]},{j,1,Length[xin]}]

Run the steepest descent to generate a new matrix:

Atest = TrainOnce[Amat, xin, yin, 1.1, (1.0 10^(-10))^2]

We can run the neural network now that it is trained on a single pair:

g[Atest.xin]

Example - Figure 14.2, 14.3

Now we'll define a function that trains on sets of data - we give tables xdata and ydata with the implicit association that xdata[[i]] -> ydata[[i]] under our map, that's the target we'll use to set the matrix values.  Note that we start with a reasonable value of a, given the input data (so that the exponential doesn't get too small too fast), and a matrix A that is drawn with random values on -1 -> 1.  TrainOnData is an implementation of Algorithm 14.1:

TrainOnData[xdata_, ydata_, eta_, eps_] := Module[{Aout, MM, NN, index},
NN = Length[xdata[[1]]];
MM = Length[ydata[[1]]];
a = 10./(Length[xdata[[1]]] Max[Abs[Flatten[xdata]]]);
Aout = Table[Random[Real,{-1,1}],{i,1,MM},{j,1,NN}];
For[index = 1, index ≤ Length[xdata], index = index+1,
Aout = TrainOnce[Aout, xdata[[index]], ydata[[index]], eta, eps];
];
Return[Aout];
]

Here we'll set up a model problem (described in Section 14.3) that takes vectors of length 10 and returns the 10th entry.  Restricting our vectors to have values drawn in discrete steps from .1, .2, .3, to .9 and shifting the input data to live symmetrically on -.4 to .4 makes the encoding in the neural network easier.  We make 10 vectors and then randomly repeat them and their correct output 1000 times to get mixx (the input) and mixy (the correct output):

SeedRandom[11111];
trainsize = 10;
xdata = Table[Table[.1 * Random[Integer,{1,9}],{j,1,10}],{k,1,trainsize}];
ydata = Table[{xdata[[j,10]]},{j,1,Length[xdata]}];
xdata = xdata - .5;
itab = Table[Random[Integer, {1,trainsize}],{j,1,100 trainsize}];
mixx = Table[xdata[[itab[[j]]]],{j,1,Length[itab]}];
mixy = Table[ydata[[itab[[j]]]],{j,1,Length[itab]}];

The matrix A is trained on this data, using eta = .5, and eps = (10^(-12)):

Amat = TrainOnData[mixx, mixy, .5, 1.0 10^(-12)];

We get perfect agreement with our training data provided we interpret our results as living on {.1,.2,.3,.4,.5,.6,.7,.8,.9}, which is reasonable in this case

Table[Flatten[Abs[.1 Round[10 g[Amat.xdata[[j]]]]   - ydata[[j]]]],{j,1,Length[ydata]}]

Below, we plot the residual values without rounding to the nearest .1 - this is Figure 14.2

ListPlot[Flatten[Table[Flatten[g[Amat.xdata[[j]]]  - ydata[[j]]],{j,1,Length[ydata]}]],PlotStyle->GrayLevel[0.5]]

Now let's test on "new" data, we'll take 100 input vectors, and run the network - the resulting residual plot is Figure 14.3.

testingsize = 100;
testdata = Table[Table[.1 * Random[Integer,{1,9}],{j,1,10}],{k,1,testingsize}];
testans = Table[{testdata[[j,10]]},{j,1,Length[testdata]}];
testdata = testdata-.5;
ListPlot[Flatten[Table[.1 Round[10 g[Amat.testdata[[j]]]]  - testans[[j]],{j,1,Length[testdata]}]],PlotStyle->GrayLevel[0.5],PlotRange->{-.1,.1}]

What percentage of the time is the neural network model correct on the test data?

tdat = Flatten[Table[.1 Round[10 g[Amat.testdata[[j]]]]  - testans[[j]],{j,1,Length[testdata]}]];
N[Length[Select[tdat,#==0.&]]/Length[tdat]]

so ~80% for this set of test data.

The matrix in this case can be roughly interpreted as picking out the last component of the input vector (these values appear in (14.27)):

Amat

Multi-Layer Perceptron

The multi-layer perceptron has the same basic form - but now we need to specify a hidden layer size, k - TrainOnce again takes the matrices (now, A and B) and iterates the steepest descent to convergence (at a minimum) for a single x-y pair - TrainOnData does this for multiple associated inputs.  There are now two sigmoids, with different parameters a and ab:

g[p_] := 1/(1 + Exp[-a p])

gb[p_] := 1/(1 + Exp[-ab p])

Plot[gb[p] /. ab-> 10,{p,-1/2,1/2}]

and our output is produced via: gb[A g[B xin] - kappa] where kappa is used to shift the data for the hidden layer - we take it to be 1/2, for the most part.

TrainOnce[Amat_, Bmat_, xin_, yin_,eta_, eps_] := Module[{z, Aout, Bout,fval,deltaA,deltaB, MM,NN,KK , zt,sigma,Omega},
NN = Length[xin];
MM = Length[yin];
KK = Length[Bmat];
Aout = Amat;
Bout = Bmat;
zt = g[Bout.xin]- kappa;
z = gb[Aout.zt];
fval = (z - yin).(z - yin);
While[fval ≥ eps,
sigma = Table[2 ab (z[[ell]] - yin[[ell]]) z[[ell]] (1 - z[[ell]]) ,{ell,1,MM}];
deltaA = Table[sigma[[ell]]zt[[i]],{ell,1,MM},{i,1,KK}];
deltaB = Table[ Sum[sigma[[ell]] Aout[[ell,i]] a (zt[[i]] + kappa)(1-zt[[i]] - kappa),{ell,1,MM}] xin[[j]],{i,1,KK},{j,1,NN}];
Aout = Aout - eta deltaA;
Bout= Bout- eta deltaB;
zt = g[Bout.xin] - kappa;
z = gb[Aout.zt];
fval = (z - yin).(z - yin);
];
Return[{Aout,Bout}];
]

I’ve commented out the print line below, but it can be useful when testing different values of eta and eps.

TrainOnData[xdata_, ydata_, Ksize_, eta_, eps_] := Module[{Aout,Bout, MM, NN, KK, index},
NN = Length[xdata[[1]]];
MM = Length[ydata[[1]]];
KK = Ksize;
a = 10./(Length[xdata[[1]]] Max[Abs[Flatten[xdata]]]);
Bout = Table[Random[Real,{-1,1}],{i,1,KK},{j,1,NN}];
Aout = Table[Random[Real,{-1,1}],{i,1,MM},{j,1,KK}];
For[index = 1, index ≤ Length[xdata], index = index+1,
{Aout,Bout} = TrainOnce[Aout,Bout,  xdata[[index]], ydata[[index]], eta, eps];
(*If[Mod[index,10] == 0, Print[index];];*)
];
Return[{Aout,Bout}];
]

Example 1

We can try out our original test case in this expanded setting - here, the hidden layer size is 2.

SeedRandom[11111];
trainsize = 10;
xdata = Table[Table[.1 * Random[Integer,{1,9}],{j,1,10}],{k,1,trainsize}];
ydata = Table[{xdata[[j,10]]},{j,1,Length[xdata]}];
xdata = xdata - .5;
itab = Table[Random[Integer, {1,trainsize}],{j,1,500 trainsize}];
mixx = Table[xdata[[itab[[j]]]],{j,1,Length[itab]}];
mixy = Table[ydata[[itab[[j]]]],{j,1,Length[itab]}];

ab = 10.0;
kappa =.5;

{Amat,Bmat} = TrainOnData[mixx, mixy, 2, .5, 1.0 10^(-12)];

Table[Flatten[Abs[.1 Round[10 gb[Amat.(g[Bmat.xdata[[j]]]-kappa)]]   - ydata[[j]]]],{j,1,Length[ydata]}]

Once again, the results are “perfect” for the training data, provided we round appropriately - and even when we don’t, they are pretty good, as shown in the residual plot below:

ListPlot[Flatten[Table[Flatten[gb[Amat.(g[Bmat.xdata[[j]]]-kappa)]  - ydata[[j]]],{j,1,Length[ydata]}]]]

On new data, the network works fairly well:

testingsize = 100;
testdata = Table[Table[.1 * Random[Integer,{1,9}],{j,1,10}],{k,1,testingsize}];
testans = Table[{testdata[[j,10]]},{j,1,Length[testdata]}];
testdata = testdata-.5;
ListPlot[Flatten[Table[.1 Round[10 gb[Amat.(g[Bmat.testdata[[j]]]-kappa)]]  - testans[[j]],{j,1,Length[testdata]}]],PlotStyle->GrayLevel[0.5]]

Example 2 - Discrete Fourier Transform for Noisy Signal - Section 14.6, Figure 14.4

Here we can train a network to pick out frequencies from data - take a noisy signal with well-defined frequency, the input will be the noisy signal, the target output, the associated frequency (we limit the allowed frequencies to a set of four, as described in Section 14.6).

SeedRandom[1111];
testsize = 100;
S = 5;
xdata = {};
ydata = {};
ivals = {1,2,3,4};
For[index = 1, index ≤ testsize, index = index+1,
ii = Random[Integer,{1,4}];
dtab = N[Table[S  Random[Real,{-1,1}] + Cos[(2 Pi ) (ivals[[ii]]) t],{t,0,1,1/100}]];
xdata = Join[xdata,{dtab}];
ydata = Join[ydata,{{2.0 ivals[[ii]]/10}}];
];

Here’s an example of the original signal (solid curve) together with the signal plus noise (points) that we ultimately send in as training data:

G0 = ListPlot[xdata[[1]],PlotStyle->GrayLevel[0.5]];
G1 = ListLinePlot[Table[Cos[2 Pi  5 ydata[[1,1]]  t],{t,0,1,1/100}],PlotStyle->GrayLevel[0]];
Show[G0,G1]

Train the neural network on the 100 pairs of x-y data generated above, with hidden layer of size 2:

ab = 10./(2. .5);
{Amat,Bmat} = TrainOnData[xdata, ydata, 2, .5,10^(-12)];

We’ll check the residuals, rounded appropriately:

check = Flatten[Table[Round[gb[Amat.(g[Bmat.xdata[[jj]]]-kappa)],.2] - ydata[[jj]],{jj,1,Length[xdata]}]];
ListPlot[check]

This is the data that appears in (14.37) (percent perfect, percent off by 1, -1, 2, -2, 3, -3, and 4, -4 respectively):

Print[N[Length[Select[Chop[check],# == 0.&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.6&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.6&]]/Length[check]],
" ",N[Length[Select[Chop[check],# == 0.8&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.8&]]/Length[check]]]

Now for a hidden layer of size 4:

ab = 10./(4 .5);
{Amat,Bmat} = TrainOnData[xdata, ydata, 4, .5,10^(-12)];

A plot of the residual:

check = Flatten[Table[Round[gb[Amat.(g[Bmat.xdata[[jj]]]-kappa)],.2] - ydata[[jj]],{jj,1,Length[xdata]}]];
ListPlot[check]

and the data that makes up (14.38):

Print[N[Length[Select[Chop[check],# == 0.&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.6&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.6&]]/Length[check]],
" ",N[Length[Select[Chop[check],# == 0.8&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.8&]]/Length[check]]]

And the same tests for a hidden layer of size 8:

ab = 10./(8 .5);
{Amat,Bmat} = TrainOnData[xdata, ydata, 8, .5,10^(-12)];

check = Flatten[Table[Round[gb[Amat.(g[Bmat.xdata[[jj]]]-kappa)],.2] - ydata[[jj]],{jj,1,Length[xdata]}]];
ListPlot[check]

Below is the data for (14.39):

Print[N[Length[Select[Chop[check],# == 0.&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.6&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.6&]]/Length[check]],
" ",N[Length[Select[Chop[check],# == 0.8&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.8&]]/Length[check]]]

Moving on to new data, with our 8 hidden layer, trained network:

testsize = 100;
S = 5;
xdata = {};
ydata = {};
ivals = {1,2,3,4};
For[index = 1, index ≤ testsize, index = index+1,
ii = Random[Integer,{1,4}];
dtab = N[Table[S  Random[Real,{-1,1}] + Cos[(2 Pi ) (ivals[[ii]]) t],{t,0,1,1/100}]];
xdata = Join[xdata,{dtab}];
ydata = Join[ydata,{{2.0 ivals[[ii]]/10}}];
];

We generated 100 new x-y pairs, and the residuals are:

check = Flatten[Table[Round[gb[Amat.(g[Bmat.xdata[[jj]]]-kappa)],.2] - ydata[[jj]],{jj,1,Length[xdata]}]];
ListPlot[check]

with “statistics”, given in (14.40):

Print[N[Length[Select[Chop[check],# == 0.&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.2&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.4&]]/Length[check]]," ",N[Length[Select[Chop[check],# == 0.6&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.6&]]/Length[check]],
" ",N[Length[Select[Chop[check],# == 0.8&]]/Length[check]]," ",N[Length[Select[Chop[check],# == -0.8&]]/Length[check]]]