Analysis of Models Using Eigenvalues and Eigenvectors 

We will look at a few modesl similar to the ones presented in your text that use the eigenvalues and eigenvectors to study different discrete models.  

The first model will be a slight variation of the model presented in Example 3.12 on Deer Population.   

Example 1.  Deer Population 

We want to know how to maintain the deer population if no hunting is allowed.  We wish to know first what the population will by at the end of say 20 years and then predict the steady state if it exists.  Fawns and yearlings will be lumped into one category called young and the other will be the adult population.  Surveys have shown that 30% of the young change to adults each year and that 20% of the young die each year.  This tells us that 50% of the young in year k are still classified as young in year k+1.  Also, it is observed that adults reproduce at a rate of 20% of the adult population per year.  Also, it is observed that in the absence of hunting 12% of the adult population die naturally each year.  An extensive survey shows that there is a population of 500 young and 1000 adults currently in the ares.  Let y(k) represent the number of deer in the young category in kear k and let a(k) represent the number of deer in the adult category in year k.  The discrete model of the set of equation is then given by 

                           y(k+1) = .50y(k) + .20 a(k) 

                           a(k+1) = .30y(k) + .88 a(k) 

>	with(linalg):

 

>	with(LinearAlgebra):

 

>	A:=<<.50,.30>|<.20,.88>>;

 

>	EA:=Eigenvectors(A,output='list');

 

>	p1:=EA[1][3][1];

 

>	p2:=EA[2][3][1];

 

>	x0:=matrix([[500],[1000]]);

 

>	P:=<p1|p2>;

 

>	Lambda:=MatrixInverse(P).A.P;

 

We want to look at x(k+1) = A x(k) or x(k+1) = x(0).   we can then look at   and then we have x(k+1)= x(0).   To see what happens to the population at the end of 15 years then we set k+1 =15 and look at the right-hand side of the equation.  

>	x15:=P.Lambda^(15).MatrixInverse(P).x0;

 

The steady state is when we let times go to infinity which means that the eigenvalue that is less than one will disappear.   

>	xs:=P.<<1,0>|<0,0>>.MatrixInverse(P).x0;

 

Can we reacy any conclusion about allowing hunting?   

Lets see what happens if they allow a harvest of 5% of the adult population each year.  

This would reduce the adult population by another 5 %  hence our coefficient matrix would be 

>	A[2,2]:=.83;

 

>

A;

 

>	EAa:=Eigenvectors(A,output='list');

 

>	EAa[2][3][1];

 

>	Pa:=<p1|p2>;

 

>

 

>	Lambda:=MatrixInverse(Pa).A.Pa;

 

>	x15:=Pa.Lambda^(15).MatrixInverse(Pa).x0;

 

Example 2 

Weather model.  We will discuss Exercise 3.2,8 in your text.  I will give a brief discription of the  

model at this time.  The model will predict that there is precipitation, it is cloudy or it is sunny. 

It is determined that the probability that a day with precipitiation is followed by another day with precipitation is 1/4 and that the probability that a day with precipitation is followed by a cloudy day is 1/2.  This also tell us that the probability that a day with precipitation is follweed by a sunny day is 1/4.  It is also determined that the probability that a cloudy day is followed by a day with precipitation is 1/4 and the probability that a cloudy day is followed by a cloudy day is 1/4.  This then tells us that the probability that a cloudy dayis followed by a sunny day is 1/2.  Lastly, it is determined that the probability that a sunny day is followed by a day with precipitation is 1/8 and that the probability that a sunny day is followed by a cloudy day is 1/8.  This also, tells us that the probability that a sunny day is followed by a sunny day is 3/4.  Let p1(k), p2(k) and p3(k) denote the probability of precipition, cloudy and sunny respectively on the kth day.  Let p1(k+1), p2(k+1) and p3(k+1) denote the probability of preciptation, cloudy, or sunny respectively on the next day or the (k+1)st day.  This leads to the set of equations.  

    p1(k+1) = 1/4 p1(k) +1/4 p2(k) +1/8 p2(k) 

    p2(k+1) = 1/2 p1(k) + 1/4 p2(k) +1/8 p2(k) 

    p3(k+1) = 1/4 p1(k)  + 1/2 p2(k) +3/4 p2(k) 

We then have the model p(k+1) = W p(k) where W is the transition matrix.  Let us start with initial probabilities of p(0) = [1/4,1/4,1/2];  What is the prediction of the weather 5 days later.  

>	W:=<<1/4,1/2,1/4>|<1/4,1/4,1/2>|<1/8,1/8,3/4>>;

 

>	EW:=evalf(Eigenvectors(W,output='list'));

 

>	p1:=EW[1][3][1];

 

>	p2:=EW[2][3][1];

 

>	p3:=EW[3][3][1];

 

>	P:=<p1|p2|p3>;

 

>	Lambda:=simplify(MatrixInverse(P).W.P);

 

>	p0:=<1/4,1/4,1/2>;

 

>	p5:=P.(Lambda^5).MatrixInverse(P).p0;

 

>

 

A Few Properties of Matrices Related to Eigenvalues 

Def A real  matrix A is said to be normal it .   

If the matrix is conplex then we would have that the transpose is replaced by conjugate transpose.   

We see right away that every symmetric matrix satisfies the about property and hence is normal 

Can you think of any other normal matrices.   

Try the following  

                                         A =  

>

AT:=A^%T;

 

>

AT.A;

 

>

A.AT;

 

Can we get a little more general.  

>	Aa:=<<0,2>|<-2,0>>;

 

>	AaT:=Aa^%T;

 

>

AaT.Aa;

 

>

Aa.AaT;

 

What type of matrix are we now looking at?   

 so we wee that and so any ski-symmetric matrix is normal. 

Also how about orthogonal and Unitary matrices?   

Let us now look at matrix polynomials.  i.e. 

Given an matrix A we want to look at a matrix polynomial of the form 

  p(A) = + + ...
(this is called a matrix polynomial).  Then the eigenvalues of the matrix polynomial p(A) are  

, .., .   

>	A:=<<1,0,0>|<0,3,1>|<0,1,3>>;

 

>	EA:=Eigenvectors(A,output='list');

 

>	I3:=IdentityMatrix(3);

 

>	PA:=3*A^2 + 5*A + 3*I3;

 

>	lambda[1]:=EA[1][1];

 

>	lambda[2]:=EA[2][1];

 

>	lambda[3]:=EA[3][1];

 

>	Eigenvalues(PA);

 

>	Pl1:=3*lambda[3]^2 + 5*lambda[3] +3;

 

>	Pl2:=3*lambda[2]^2 + 5*lambda[2]+3;

 

>	Pl3:=3*lambda[1]^2 + 5*lambda[1] +3;

 

>

 

Recursion Formula and Matrices 

I tried to make a correction to this example and something happened when I tried to export the file as an html file.  Look in you text for this type of example.  

We want to study the following recursion formula using matrix theorm.  In partiuclar we want to use eigenvalues and eigenvectors to find the value for large values of n.   

Let = 0 and and for n=1,2,.. the following values of the sequence are given by the formula 

                                               
We would like to find the value of say without having to find all the values between and .    We will write the rcursion formula as a system of equation but this will require that we introduce another equation.  

=  

We can now write the system in matrix form 

=  and then 

=  

So we need to replace the matrix by a diagonal matrix if possible .  This means we need to find the eigenvalues and eigenvectors of the coefficient matrix.  

>	with(LinearAlgebra):

 

>	A:=<<0,3/4>|<1,1/4>>;

 

>	EA:=Eigenvectors(A,output='list');

 

>	p1:=EA[1][3][1];

 

>	p2:=EA[2][3][1];

 

>	P:=<p1|p2>;

 

>	PIn:=MatrixInverse(P);

 

>	Lambda:=PIn.A.P;

 

>	L50:=evalf(Lambda^50);

 

>	T50:=P.L50.PIn.<0,1>;

 

>

 

Exercises  

Exercises 3.1.13 (page 341), Exercise 3.2.14 (page 367), Exercise 3.3.2 (page) 379  Exercises are due  December 1, 2008 .   Each of the problems need to be worked by finding the eigenvalues and eigenvectors and then diagonalizing the coefficient matrix like the examples in the lecture.   Again they can be turned in at the end of class that day or e-mail a worksheet to me.  The Take Home part of the exam is already on the web sight.  You can start on it on December 3, 2008.  It should take a little more than one class period.  It is due on December 12, 2008.      

>