Introduction to Neural Networks
U. Minn. Psy 5038

Daniel Kersten

Introduction

Last time

Modeled aspects of spatial lateral inhibition in the visual system:
	Simple linear network: intensity->spatial filter->response
	intensity and response were represented as vectors
	the spatial filter was represented as a matrix

If the spatial filter is a neural network, then the elements of the matrix can be interperted as strengths of synaptic connections. The idea is to represent the synaptic weights for the ith output neuron by the values in the ith row of the matrix. The matrix of weights is sometimes called a connection matrix just because non-zero weights are between connected neurons, and neurons that aren't connected get fixed zero weights.
This model is what we called a linear feedforward network, and used the generic neuron model. 
We generalized the linear discrete model to a linear continous one with feedback. Although more complex, its steady-state solution is identical to the linear feedforward model. 
We studied the continous system by approximating it as a discrete time system with ϵ = Δt as a free parameter.

Today

There is a large body of mathematical results on linear algebra and matrices, and it is worth our while to spend some time going over some of the basics of matrix manipulation. 
We will first review  matrix arithmetic (addition and multiplication). We will review the analog of division, namely finding the inverse of a matrix--something that was used in Lecture 5 to show how the steady-state solution of the feedback model of lateral inhibition was equivalent to a feedforward model with the appropriate weights.
You may wonder at times how all this is used in neural modeling. But as this course goes on, we will see how otherwise obscure notions of things like an "outer product" between two vectors, or the "eigenvectors" of a matrix are meaningful for neural networks. For example, the "outer product" between two vectors can be used in modeling learning, and the eigenvectors corresponding to a "matrix of memories" can represent stored prototypes, and in dimensionality reduction.

Basic matrix arithmetic

Definition of a matrix: a list of scalar lists

Traditional and Standard form output format defaults.

Defining arrays or lists of function outputs using indices

As we've already seen, Table[] can be used to generate a matrix, or list of lists. E.g.

H = Table[i^2 + j^2, {i, 1, 3}, {j, 1, 3}]

And we can view it in the traditional MatrixForm either by: Table[i^2+j^2,{i,1,3},{j,1,3}]//MatrixForm, or

MatrixForm[H]

Before going on, use the Cell menu to set: Cell->Default Output FormatType->TraditionalForm. This will cause vectors and matrices to be printed out in traditional mathematics format. Then we don't have to explicitly type MatrixForm every time to see the output in traditional matrix form. But if you want to see the output in standard Mathematica form, use //StandardForm  or StandardForm[ ].

Defining arrays or lists of symbols

Although most of the time, we'll be working with numerical matrices, Mathematica also allows one to specify arrays or lists of variables.
An mxn matrix has m rows, and n columns. Here is a 3x4 matrix of symbolic elements w[i,j]:

Array[w, {3, 4}]

You can produce the same array with: Table[w[i,j],{i,1,3},{j,1,4}]]. You can also have Mathematica display a subscripted list

L=Table[(Subscript[w, i, j]), {i, 1, 3}, {j, 1, 4}]

Mathematica lists can have elements of almost any type:

(Subscript[w, 2, 3]) = 3.0;
(Subscript[w, 3, 4]) = Pi;
(Subscript[w, 2, 4]) = "hi";
(Subscript[w, 1, 1]) = 3;

L

Adding, subtracting and multiplying by a scalar

As with vectors, matrices are added, subtracted, and multiplied by a scalar component by component:

A = {{a,b},{c,d}};
B = {{x,y}, {u,v}};

Add A to B:

A+B

Subtract B from A:

A-B

Multiply A by 3:

3 A

Multiplying two matrices

We have already seen how to multiply a vector by a matrix: we replace the ith row of the output vector by the inner product of the ith row of the matrix with the vector.
In order to multiply a matrix A, by another matrix B to get C = AB, we calculate the ijth component of the output matrix by taking the inner product of the ith row of A with the jth column of B:

A.B
%//MatrixForm

Note that AB is not necessarily equal to BA:

B.A

Laws of commutation, association and distribution

In particular, look at the element in the upper left of the matrix BA above--there is no reason, in general, for ax+bu to equal ax + cy. That is, matrix multiplication does not commute.

Apart from commutation for matrix multiplication, the usual laws of commutation, association, and distribution that hold for scalars hold for matrices. Matrix addition and subtraction do commute. Matrix multiplication is associative, so (AB)C = A(BC). The distributive law works too: A(B+C) = AB + AC

Non-square matrices

It is not necessary for A and B to be square matrices (i.e. have the same number of rows as columns) to multiply them. But if A is an mxn matrix, then B has to be an nxp matrix in order for AB to make sense. For example, here F is a 3x2 matrix, and G is a 2x4 matrix.

F = {{a,b},{c,d},{e,f}}
G = {{p,q,r,s},{t,u,v,w}}

Check the dimensions with

Dimensions[F]
Dimensions[G]

Because F has 2 columns, and G has 2 rows, it makes sense to multiply G by F:

Print[F, G, " =",  F . G]

However, because the number of columns of G (4)  do not match the number of rows of F (3), G.F is not well-defined:

Try this:

Print[F, G, " =",  G . F]

Inverse of a Matrix

Dividing a matrix by a matrix: the identity matrix & matrix inverses

The matrix corresponding to 1 or unity is the identity matrix. Like 1, the identity matrix is fundamental enough, that Mathematica  provides a special function to generate n-dimensional identity matrices. Here is a 2x2:

IdentityMatrix[2]

It is easy to show that the identity matrix plays the role for matrix arithmetic that the scalar 1 plays for scalar arithmetic.

How can one divide one matrix, say B, by another, say A? We can divide numbers, x by y, by multiplying x times the inverse of y, i.e. 1/y. So to do the equivalent of dividing B by A, we need to find a matrix Q such that when A is multiplied by Q, we get the matrix equivalent of unity, i.e. the identity matrix. Then "B/A" can be achieved by calculating the matrix product: B.Q.

A = {{a,b},{c,d}};

Mathematica  provides a built-in function to compute matrix inverses:

Q = Inverse[A]

We can test to see whether the product of a A and Q is the identity matrix, but Mathematica  won't go through the work of simplifying the algebra in this case, unless we specifically ask it to.

Q.A

Simplify[Q.A]

Here is a simple numerical example. Verify that R is the inverse of B by seeing whether the product is the Identity matrix.

B = {{1,-1},{3,2}};
R = Inverse[B]

Singular and badly conditioned matrices

What if one row is a scaled version of another row?

B1 = {{1.5,1},{3,2.0}}

Then the rows of the matrix are not linearly independent. In this case, the inverse is not defined.

Ask for the inverse of B1

Inverse[B1]

Sometimes the rows are almost, but not quite, linearly dependent (because the elements are represented as approximate floating point approximations to the actual values). Mathematica may warn you that the matrix is badly conditioned. Mathematica  may try to find a solution to the inverse , but you should be suspicious of the solution. In general, one has to be careful of badly conditioned matrices.

Let's try finding the inverse of the following matrix:

B2 = {{-2,-1},{4.00000000000001,2.0}};
Inverse[B2]

Were we lucky or not?

B2 . Inverse[B2]

Determinant of a matrix

There is a scalar function of a matrix called the determinant. If a matrix has an inverse, then its determinant is non-zero. Does B2 have an inverse?

Det[B2]

Det[B]

Matrix transpose

Transpose definition

We will use the transpose operation quite a bit in this course.  It interchanges the rows and columns of a matrix:

Y  =  Table[(Subscript[y, i, j]), {i, 1, 3}, {j, 1, 4}]

In StandardForm on the input line, the transpose is written:

Transpose[Y]

The output default is TradtionalForm. On the input line, in TraditionalForm, the transpose is written:

(Y^T) // StandardForm

What do the ouputs, X and X\(TraditionalForm\`\^T\) look like in StandardForm? (Make sure your input form is appropriate)

Getting parts of a matrix

We can pull out the ith row of a matrix in Mathematica by simply writing W[[i]]. For example, the 2nd row of Y is:

Y[[2]]

What about ith column of a matrix? There is no equally simple way of getting the column of a matrix in Mathematica, but we can use the transpose operation to do it. Transpose[W] [[i]] produces the ith column of matrix W. For example, the 3nd column of Y is:

Transpose[Y][[3]]

Try putting //ColumnForm after the above expression.

Pull out the element (Subscript[x, 2, 3])from Y

Symmetric matrices

For a square matrix, the diagonal elements remain the same under transpose.
If the transpose of a matrix equals itself, HT = H, H is said to be a symmetric matrix. Symmetric matrices occur quite often in physical systems (e.g. the force on particle i by particle j is equal to the force of j on i). This means that the elements of a symmetric matrix H look like they are reflected about the diagonal. We constructed a symmetric matrix in Lecture 5 when we modeled the exponential drop-off in lateral inhibition:

H = Table[N[Exp[-Abs[i-j]],1],{i,5},{j,5}]

The notion of a Hermitian matrix is a generalization of symmetric, where H* = H, and H* is the transpose of the complex conjugate of H.

Neural networks and symmetric connections

Do neural systems have symmetric connections? Most real neural networks probably do not. Although, when the nature of the processing would not be expected to favor a particular asymmetry, we might expect that there should be symmetric connections on average. We made this assumption when setting up our lateral inhibition weight matrix, as seen above for H. Symmetric matrices have so many nice properties, that they are appealing to neural modelers on theoretical grounds. We'll see this later when we study the Hopfield networks. Lack of symmetry can have profound effects on the dynamics of non-linear networks and can produce chaotic trajectories of the state vector.

Outer product of two vectors

We've already seen that the inner product of two vectors produces a scalar. In the next lecture, we will begin our discussion of Hebbian learning in neural networks. In this context, the outer product of an input and output vector will be used to model synaptic modification. Consider two vectors:

fv = Table[(Subscript[f, i]), {i, 1, 3}]
gv = Table[(Subscript[g, i]), {i, 1, 3}]

The outer product is just all of the pairwise products of the elements of f and g arranged in a nice (and special) order:

h = Outer[Times,fv,gv]

To get a preview, the idea in Hebbian learning is that an adjustment to the strength of  a connection between say the ith and jth neurons is proportional to the product of the strength of pre- and post-synaptic activities, represented by (Subscript[g, i])and (Subscript[f, j]).

The outer product is also written in traditional row and column format as: f gT

(fg^T) = (({{(Subscript[f, 1])

                (Subscript[f, 2])}, {(Subscript[f, 3])}}))((Subscript[g, 1]) (Subscript[g, 2]) 
            (Subscript[g, 3]))

We'll also use outer product to calculate autocorrelation matrices.

Eigenvectors and eigenvalues

Basic idea

Eigenvalues and eigenvectors crop up in many ways. For example, the solution to a set of linear first-order (i.e. no derivative orders bigger than 1) differential equations  can be found in terms of eigenvalues and eigenvectors. For a simplified version of the limulus equations

((df)/dt)=A.f

the solution is determined by the eigenvalues and eigenvectors of A and the initial conditions. 
Eigenvectors also crop up in statistics and related neural network algorithms. For example, principal components analysis (PCA) is used in dimensionality reduction--also important in neural networks and "natural computation".

So what are eigenvalues, eigenvectors?An eigenvector, x,  of a matrix, A,  is vector that when you multiply it by A, you get an output vector that points in the same (or opposite, depending on the sign of λ) direction as x:
	Ax = λx
where λ is an eigenvalue--a scalar that adjusts the length change of x. The eigenvalues are found by solving: 
	Det[A-λId]==0, for λ where Id is the identity matrix.

Let's find the eigenvalues of:

A = {{1,2},{3,4}};

Id = IdentityMatrix[2]

Solve[Det[A - λ*Id] == 0, λ]

Now pick one of the solutions for λ, find the vector x={(Subscript[x, 1]),(Subscript[x, 2])} that satisfies the basic definition, A.x=λx

Solve[{A . {(Subscript[x, 1]), (Subscript[x, 2])} == 
          (1/2) (5 - (Sqrt[33]))*{{1, 0}, {0, 1}} . {(Subscript[x, 1]), (Subscript[x, 2])}}, {(Subscript[x, 1]), 
        (Subscript[x, 2])}]

So any multiple of

x = {1,  -(1/6) ((Sqrt[33]) + 3)}

is a valid eigenvector.

Because the Det[ ]==0 is a polynomial of order n, there can't be any more than n distinct eigenvectors for an nxn matrix--and there may be less.

We can similarly find the second distinct eigenvector for λ -> (1/2) (5 + (Sqrt[33])), but Mathematica provides a simpler way of doing things.

Built-in functions for finding eigenvalues and eigenvectors

Eigenvalues[ ]

The eigenvalues are given by:

Eigenvalues[A]

Eigenvalues and eigenvector elements do not have to be real numbers. They can be complex, that is an element can be the sum of a real and imaginary number. Imaginary numbers are represented by multiples (or fractions) of I, the square root of -1:

Sqrt[-1]
Sqrt[-1]//OutputForm

B = {{1,2},{-3,4}};
Eigenvalues[B]

Eigenvectors[ ]

The built-in Mathematica function Eigenvectors[A] returns the eigenvectors of matrix A as the rows of a matrix, lets call eig:

eig = Eigenvectors[A]

Verify that eig[[1]] and A.eig[[1]] lie along the line (i.e. in the same or opposite directions) by taking the dot product of the unit vectors pointing in the directions of each. Use this function:

normalize[x_] := N[x/Sqrt[x.x]];

Preview of eigenvectors and dimensionality reduction

A measure of the structure of an ensemble of signals (e.g. of vectors) is the degree of correlation between their elements. Signals, such as sounds and visual images, have correlational structure that is taken advantage of in sound and image compression algorithms. One simple kind of statistical structure is characterized by the degree to which one can predict one element of a vector from a nearby element. For images, the color of a pixel at location i is a pretty good predictor of the color at location j = i+1. As j gets far from i, however, the prediction gets worse. It is possible to characterize the degree of predictability by a matrix whose ijth element is big if the ith pixel of an image is a good predictor of the jth. One measure of predictability is the autocorrelation matrix where, for example, the i,jth element represents the degree of correlation between the ith and jth pixels, i.e. that are |j = i| apart. An autocorrelation matrix has real elements and is symmetric. The eigenvectors of the matrix capture the subspace with most of the variance, i.e. where the "action" is. The eigenvalues correspond to the amount of "action"--i.e. the variance in the eigenvector directions.

Verify that the eigenvectors of the symmetric matrix above (H) are orthogonal.

Next time

Linear systems analysis & introduction to learning/memory

In the next lecture, we'll apply what we've learned about matrices and vectors to two problems. We'll first take an in-depth look at the properties of the linear neural network to see how it relates to general linear systems. We'll also show the advantages of using the eigenvectors of the linear transformation matrix as a useful coordinate system to represent input and output pattern vectors.
We will introduce the topic of learning and memory, and see how the outer product rule can be used to model associative learning.

Preview of linear systems analysis

Linear systems are an important, and tractable class of input/output models. Matrix operations on vectors provide the prototype linear system.
Consider the generic 2-layer network. It consists of a weighted average of the inputs (stage 1), followed by a point-nonlinearity (the squash function of stage 2), and added noise (stage 3). Although later we will see how the non-linearity enables computations that are not possible without it, useful functions can be realized with just the linear or stage 1 part of the network. We've seen one application already with the model of the limulus eye. In the next lecture, we will see how linear networks can be used to model associative memory. But first, let us take what we've learned so far about modeling linear networks and look at in the general context of linear systems theory. Our 2-layer net is a matrix of weights that operates on a vector of input activities ty computing a weighted sum. One property of such a system is that it satisfies the fundamental definition of a "linear system", which we will define shortly.
The world of input/output systems can be divided up into linear and non-linear systems. Linear systems are nice because the mathematics that describes them is not only well-known, but also has a mature elegance. On the other hand, it is a fair statement to say that most real-world systems are not linear, and thus hard to analyze...but fascinating if for that reason alone. Scientists were lucky with the limulus eye. That nature is usually non-linear doesn't mean one shouldn't familiarize oneself with the basics of linear system theory. Many times non-linear systems can be approximated by a linear one over some restricted range of parameter values. 
So what is a "linear system"?
The notion of a "linear system" is a generalization of the input/output properties of a straight line passing through zero. The matrix equation W.x == y is a linear system. This means that if W is  a matrix, x1 and x2 are vectors, and a and b are scalars:
                     W.(a x1 + b x2) = a W.x1 + b W.x2

This is a consequence of the laws of matrix algebra.The idea of a linear system has been generalized beyond matrix algebra.  Imagine we have a box that takes inputs such as f, and outputs g = T[f].

The abstract definition of a linear system is that it satsifies:
	                     T[a f + b g] = a T[f] + b T[g]
where T is the transformation that takes the sum of scaled inputs f, g (which can be functions or vectors) to the sum of the scaled transformation of f and g. The property, that the output of a sum is the sum of the outputs, is sometimes known as the superposition principle for linear systems. The fact that linear systems show superposition is good for doing theory, but as we will see later, it limits the kind of computations that can be done with linear systems, and thus with linear neural network models.

References

Linear Algebra and Its Applications by Gilbert Strang, 3rd Edition , Hardcover, 505 pages, Published by Hbj College & School Div. Publication date: February 1988