Introduction to Neural Networks

Multi-layer non-linear networks
Gradient descent: Learning by error back-propagation

Backprop simulation example: XOR

It is well-known that with appropriate weights, 2 weight layers with 3 hidden units can solve the XOR problem. But this is still a tough problem to learn, mainly because it requires that two very different inputs map to the same output. Let's try learning the weights.

Reading in packages

So far, we've only used standard packages that are part of the Mathematica  package. For these exercises, we will use a publicly available package written by James A. Freeman which can be downloaded from: http://library.wolfram.com/infocenter/Books/3485/

You can find out where your current default directory is by typing:

Directory[]

Then you can set your default directory (for example) to read custom packages, or to save your data in a particular place:
 SetDirectory["myfavoritedirectory"]

...but you'll have to determine where your Backpropagation.m file is, and set the directory appropriately. After you've set the directory, you can read in the package:

<<Backpropagation.m

That said, it is usually easier to bring up a window for file browsing using:

Import[SystemDialogInput["FileOpen"],"Package"];

and then find and select Backpropagation.m.

Standard backprop

We will first try a straightforward implementation of the algorithm described above , called bpnStandard[], which is in  Backpropagation.m.  You can open up the package and see how this and other functions are defined.  But we replicate it here to show that the basic operations of error backpropagation are rather straightforward. You don't have to execute the following function because if you've read it in, it is defined in Backpropagation.m.

The bpnStandard backprop function

Some sample inputs:

ioPairsXOR = { {{0.9,0.9},{0.1}}, {{0.1,0.1},{0.1}}, {{0.1,0.9},{0.9}},
{{0.9,0.1},{0.9}} };

bpnStandard[inNumber_,hidNumber_,outNumber_,ioPairs_,eta_,numIters_] :=
  Module[{errors,hidWts,outWts,ioP,inputs,outDesired,hidOuts,  outputs, outErrors,outDelta,hidDelta},
	hidWts = Table[Table[Random[Real,{-0.1,0.1}],{inNumber}],{hidNumber}];
	outWts = Table[Table[Random[Real,{-0.1,0.1}],{hidNumber}],{outNumber}];
 	errors = Table[
                        (* select ioPair *)
    ioP=ioPairs[[Random[Integer,{1,Length[ioPairs]}]]];
    inputs=ioP[[1]];
    outDesired=ioP[[2]];
                        (* forward pass *)
    hidOuts = sigmoid[hidWts.inputs];
    outputs = sigmoid[outWts.hidOuts];
                        (* determine errors and deltas *)
    outErrors = outDesired-outputs;

outDelta= outErrors (outputs (1-outputs));

Note that the identity that we derived earlier relating the derivative of the sigmoid to the sigmoid was used in the above Mathematica expression--(outputs (1-outputs).

hidDelta=(hidOuts (1-hidOuts)) Transpose[outWts].outDelta; (* update weights *)

outWts += eta Outer[Times,outDelta,hidOuts];
hidWts += eta Outer[Times,hidDelta,inputs];

(Note the C-style notation: x += a, is the same as: x = x + a)

(* add squared error to Table *)
       outErrors.outErrors,{numIters}];  (* end of Table *)
 Return[{hidWts,outWts,errors}];
 ];

Running the algorithm

We will first try a standard backpropagation network with 2 input units, 3 hidden layer units, and 1 output unit  (same network diagram as at the beginning of this notebook). The learning constant eta, we will set to 5. And let's try it for 1500 iterations.
bpnStandard[] expects a list with the input/output pairs set up as follows (for an XOR training set).

?bpnStandard

ioPairsXOR = { {{0.9,0.9},{0.1}}, {{0.1,0.1},{0.1}}, {{0.1,0.9},{0.9}},
{{0.9,0.1},{0.9}} };

Timing[outs=bpnStandard[2,3,1,ioPairsXOR,5,1500];]

Did the net converge? No. The errors wiggle all over and never settle down near zero.

ListPlot[outs[[3]],Joined->True]

We can see specifically where it is failing by calling bpnTest[]:

bpnTest[outs[[1]],outs[[2]],ioPairsXOR];

Improving the standard algorithm by preventing overlearning of certain patterns

The above network had a hard time learning the first pattern. Even with more iterations, you may discover the network to be stuck in a local minimum of the error function. 
You could try more units in the hidden layer. This probably won't help much.
One trick to improve the odds of convergence is avoid over-learning certain patterns. The function bpnMomentumSmart sets a maximum acceptable error for a given pattern to 0.1, and then if the error is less than that for a given iteration, the weights are not updated. The idea is to concentrate more on learning the associations where the net shows some obstinance. (This network also uses a momentum term, the weight of this term is determined by alpha, which below is set to 0.9. Momentum is discussed below)
But you have to be lucky. I tried the following several times, before the algorithm nailed it:

outs={0,0,0,0};
Timing[
outs=bpnMomentumSmart[2,3,1,ioPairsXOR,2.0,0.9,1300];]

ioPairsXOR

bpnTest[outs[[1]],outs[[2]],ioPairsXOR];

ListPlot[outs[[3]],Joined->True]

Momentum

A standard modification to backprop that typically has a significant effect on learning speed is a momentum term. The idea, as the name suggests, is to keep the weight changes moving in about the same direction that they have been going. For example, for standard backprop, the hidden units were updated as:

hidWts += eta Outer[Times,hidDelta,inputs];

With momentum, the weights are updated in the same way except that alpha times the previous weight update is added in:

hidLastDelta = 
eta Outer[Times,hidDelta,inputs]+ alpha hidLastDelta;
hidWts += hidLastDelta;

The momentum term was included in bpnMomentumSmart[].

Exercise: Can you find better learning parameters eta, and momentum term alpha?

Exercise: Try learning OR and AND.

ioPairsOR = { {{0.9,0.9},{0.9}}, {{0.1,0.1},{0.1}}, {{0.1,0.9},{0.9}},
{{0.9,0.1},{0.9}} };

Timing[outs=bpnStandard[2,3,1,ioPairsOR,5,1500];]

bpnTest[outs[[1]],outs[[2]],ioPairsOR];

ListPlot[outs[[3]],Joined->True]

outs[[1]]
outs[[2]]

Exercise: Can you find a learning sequence which improves the odds of standard backpropagation finding a solution to the XOR problem?

Appendix

Logistic function - a smooth, differentiable non-linearity

As we saw above, it can be useful to have a non-linearity which is smooth enough to be differentiable. D[] returns an expression for the derivative, so to define a function that is the derivative of another we write:

Df1[x_] := D[f[t],t] /.t->x
Df1[x_] := Evaluate[D[f[t],t]]

The more direct way is to use "operators" or functionals that take functions as inputs and return functions as outputs. Derivative[], or f'[x]  return functions:

(*Df[x_] := Derivative[f[x],x]*)
Df[x_] := f'[x]

Note that the derivative has a particularly simple expression in terms of f[x], which you can verify by comparing Simplify[Df[x]] with Simplify[hh[x]], or Simplify[Df[x]-hh[x]]:

hh[x_] := f[x](1-f[x]);

Here is a plot of the sigmoid and its derivative:

Plot[{Df[x],f[x]},{x,-3,3}]

References

Rumelhart, D. E., McClelland, J. L., & Group, and the PDP Research Group (1986). Parallel Distributed Processing. Explorations in the Microstructure of Cognition . Cambridge,  Massachusetts: MIT Press. See chapter 9 in Vol. I, page 318.
Freeman, J. A. (1994). Simulating Neural Networks with Mathematica . Reading, MA: Addison-Wesley Publishing Company.
http://www.mathsource.com/cgi-bin/MathSource/Publications/BookSupplements/Freeman-1993/0205-906