CHAPTER 13
DIFFERENCE EQUATIONS

There are many cases in which some economic variable is related to its lagged values by discrete time intervals.  Equations that express such a relationship between the present value of a variable and its values in the past are called difference equations.  For example, yt = byt-1 is a first-order, homogenous linear difference equation.  The problem here is to find the value of Y at time t, given that the value of y at some time in the past is known.
The order of a difference equation is the highest time difference in the equation.  For example, yt+1 + 3yt +2 = 0 and yt + 2yt+1 + yt+2 = 0 are first-order and second-order difference equations, respectively.  The first two sections of this chapter are concerned with finding solutions of the first-order linear difference equations, followed by the stability analysis of the first-order linear difference equations.  The following sections concentrate on the solutions of the second-order linear difference equations and the stability analysis of the second-order difference equations.  Finally, we study the solution and stability of the system of linear difference equations.

TOOL BOX
Mathematica  Commands				Function
<<DiscreteMath`RSolve`				Loads Discreate Math Package
RSolve[{y[t] - ay[t-1] ==b, y[0]==yo}, y[t], t]		Solves the difference equation yt-ayt-1=b, Y0= Yo, for Yt.
H[{y1_, y2_}] = {2y1+3y2,  y1- 2y2}			Defines any values of y1 and y2 as a set of functions 								defined in H
f[y1_, y2_, t_] = NestList[H, {y1,y2}, n]		Nestlists y1 and y2 of H n times
f[a, b, n]						Finds the first n values visited around y1= a, y2 = b. 

I.  First-Order, Homogenous Linear Difference Equations
Suppose that the value of a variable at time t is a constant porportion b of its value at time t-1.  This assumption gives rise to a first-order, homogenous linear difference equation, the general form of which is  yt = byt-1.  This equation yields the following sequence of values for y:  y1 = byo,  y2 = by1 = (b^2)yo,  y3 = by2 = (b^3)yo, . . . .,  yt = yt-1 = (b^t)yo.   It is obvious that, yt = (b^t)yo satisfies equation yt = byt-1, and is a solution of that equation.  If the initial value of y at time t = 0 is also given, then the initial-value problem yt = byt-1, yo = a, has a solution yt = a(b^t).   This equation satisfies both the initial condition and the difference equation.

To solve the difference equations using Mathematica, first load the package "DiscreteMath`RSolve' , then use the command "RSolve".   The form of the equations given in the "RSolve" command is similar to the form used in the "DSolve" command of the previous chapter.  You can use a single equation, equations with initial conditions or several equations. The solution is returned in the form of a list of replacement rules. 

Important Notice:  The package "DiscreteMath" should be loaded before any program is run in this chapter.  To load DiscreteMath, type <<DiscreteMath`RSolve` and press Shift Enter.  If you run a program before running <<DiscreteMath`RSolve`, the program will not run and attempts to load DiscreteMath will also fail.  In such cases, close the notebook, exit Mathematica,  and run <<DiscreteMath`RSolve` before running any program. 

Example:  Solve the homogeneous difference equation yt = byt-1, y(0) = yo.
Solution:  Following the solution rule of the linear homogeneous difference equations, y = yo (b^t).  Using Mathematica, the problem can be solved as, 
<<DiscreteMath`RSolve`
Clear[y]
RSolve[{y[t]==b y[t-1], y[0] == yo}, y[t], t] 
{{y[t] -> (b^t) yo}}

Example:  Solve the homogeneous difference equation yt = .2yt-1, yo = 7.
Solution:  Following the solution rule of the homogeneous difference equations,  y = 7 (2^t).  Using Mathematica, the problem can be solved as,
Clear[y]
RSolve[{y[t]==.2y[t-1], y[0]==7}, y[t], t]
\(\(\({{Null\ \((y[t] -> 
                7.`\ 0.2`\^t)\)}}\)\(
\)\)\)
Example:  Suppose that an initial amount of A is saved at an interest rate of i.  What is the balance in the account at every point at time?
Solution:  The problem here is the economic problem of compounding interest, where St = (1 + i)St-1, So = A.  Using the solution rule of the initial-value, first-order difference equation, the solution is St = A(((1 + i))^t).

Example:  What will be the value of $10000 invested at an interest rate of 6% compounded annually, after 5, 10, and 15 years?
Solution: Solve the difference equation yt = (1 + .06t)yt-1,  y(0) = 10000 and evaluate yt for t = 5, t = 10, and t = 15.  Plot the time path of the investment growth.
Clear[ y, rose, path1]
rose = RSolve[{y[t]==1.06y[t-1], y[0]==10000}, y[t], t];
y[t_]=y[t]/.rose[[1]]
\(10000.`\ 1.06`\^t\)
y[5]
13382.3
y[10]
17908.5
y[15]
23965.6
to plot the time path of  y over t, list the values of the y(t) for a certain range of t and list plot the table of values.  (Note:  If the coefficient of Yt-1 is positive, y(t) can be graphed directly using "Plot" command.  However, if it is negative, the discrete values of y(t) should be evaluated first using "Table" command, then the values can be graphed.)
path1=Table[y[t], {t, 0, 10}];
ListPlot[path1, PlotJoined->True]

Null (⁃Graphics⁃)

Do More Practice:
1- Change the coefficient of yt-1 to 2.  Re-run the program and observe the time path of y(t).  
2- Change the coefficient of yt-1 to -2.  Re-run the program and observe the time path of y(t).  
3- Change the coefficient of yt-1 to -.2.  Re-run the program and observe the time path of y(t).  

Example:  Contrary to Domar's model,  which is a continouous time growth model of the economy (see chapter 12), Harrod's model treats time as being discrete ( R. F. Harrod, "An Essay In Dynamic theory,"  Economic Journal, March, 1939; and Toward a Dynamic Economics, New York: Macmillan, 1952).   The other assumptions of the model, however, are similar to Domar's model.  That is, saving (S) is proportinal to income (y), St = syt, investment ( I ) is proportional to change in income, It = c(yt - yt-1), and the equilibrium condition, St = It.  The problem is to find the time path of income.  
Solving the three equations together yields yt = ((cyt-1)/(c - s)), which is a difference equation.  Therefore, given the initial income level yo,  the time path of y is determined as yt = yo(b^t),  where b = (c/(c - s)).

Example:   The actual saving in an economy, St, is a constant proportion .5 of the economy's income, (yt), and the economy's desired investment ( It ) is a constant proportion 8 of the change in income, DYt = Yt - Yt-1.  Assuming that the desired investment is realized find the time path of income, if Y(0) = 6000.
Solution:  The problem can be formulated as, St = .5Yt, It = 8(Yt - Yt -1), St = It.   Substituting for St and It from the correspoding equations, yields .5Yt = 8(Yt - Yt -1), or Yt = (8/7.5)Yt-1, which is a diffrence equation with an initial value of Y(0) = 6000.  To solve the equation and to graph the time path of Yt, use Mathematica as,
Clear[ y, sage, path2]
sage = RSolve[{y[t]==(8/7.5)y[t-1], y[0]==6000}, y[t], t];
y[t_] = y[t]/.sage[[1]]
\(6000.`\ 1.0666666666666667`\^t\ Null\^2\)
path2=Table[y[t], {t, 0, 10}]
ListPlot[path2, PlotJoined->True] 
{6000.,6400.,6826.67,7281.78,7767.23,8285.04,8837.38,9426.54,10055.,10725.3,11440.3}

⁃Graphics⁃

Problems:
1.   Find the order of the following difference equations.
	a)  yt = 3yt-1 + yt-2 + 1
	b)  yt = 2yt-2 + 5
	c)  yt = 3yt+1 + (e^t)
	d)  yt+1 = yt+2 + 2		
2.   Solve the following initial-value difference equations.
	a)  4yt + 3yt-1 = 0,	yo = 2
	b)  yt = yt-1,		yo = 10
	c)  yt+1 + 2yt = 0,	yo = 2
	d)  2yt+1 = yt,		yo = 10
	e)  yt = -.5yt-1,		yo = 10
		
3.   What are the compunded values of $20000 deposited at i = 8%, one year, two years, five years, and ten years from now.	
4.  Suppose that the actual saving in an economy, St, is a constant proportion .2 of the economy's income, Yt, and the economy's desired investment ( It ) is a constant proportion 7 of the difference between this year's income, Yt, and the last year's income, Yt-1.  Assuming that the desired investment is realized, find the time path of income.

 Stability Analysis of the Solutions of First Order, Homogenous Linear Difference Equations:  The solution of the difference equation Yt = bYt-1 is Yt = Yo(b^t),  where Yo and b are known constants.  The time path of Y, therefore, depends on the values of  Yo and b.    The initial value of Y at t = 0, or Yo decides the position of Y at the initial period.  The time path of Y, however, depends on the absolute value and the sign of b.  Cases that may arise are as follows:
 
 Case 1:   b > 1.    This gives rise to an divergent nonoscillatory time path.  
 
 Example:  Graph the time path of Yt = 2Yt-1, Yo = 7.
  Solution:
Clear[y, okra, path3] 
okra = RSolve[{y[t]==2y[t-1], y[0]==7}, y[t], t];
 y[t_]=y[t]/.okra[[1]]
7 (2^t)
 path3=Table[y[t], {t, 0, 10}];
 ListPlot[path3, PlotJoined->True]

Null (⁃Graphics⁃)
 
Do More Practice:
1- Change Yo to Yo = -7.  Re-run the program and compare the results.
2- Change Yo to Yo = 1.  Re-run the program and compare the results.
3- Change Yo yo Yo = -1.  Re-run the program and compare the results. 
  
Case 2:   b = 1.   This gives rise to a constant time path for Y .  
 
 Example:  Graph the time path of Yt = Yt-1, Yo = 7.
  Solution:
 Clear[y, iris, path4]
 iris = RSolve[{y[t]==y[t-1], y[0]==7}, y[t], t]
{{Null (y[t]->7)}}
 y[t_]=y[t]/.iris[[1]];
 path4=Table[y[t], {t, 0, 10}];
 ListPlot[path4, PlotJoined->True]
 

(Null^2) (⁃  Graphics  ⁃)

Do More Practice:
1- Change Yo to Yo = -7.  Re-run the program and compare the results.
2- Change Yo to Yo = 1.  Re-run the program and compare the results.
3- Change Yo yo Yo = -1.  Re-run the program and compare the results. 
                                             

 Case 3:   0 < b < 1.   This gives rise to a convergent nonoscillatory time path, where the value of y approaches zero as t -> ∞.

 Example:  Graph the time path of Yt = .5Yt-1, Yo = 7.
 Solution:
ClearAll[y, fern, path5]
fern=RSolve[{y[t]==.5y[t-1], y[0]==7}, y[t], t]
\({{Null\ \((y[t] -> 7.`\ 0.5`\^t)\)}}\)
y[t_]=y[t]/.fern[[1]];
path5=Table[y[t], {t, 0, 10}];
ListPlot[path5, PlotJoined->True]

(Null^2) (⁃  Graphics  ⁃)

Do More Practice:
1- Change Yo to Yo = -7.  Re-run the program and compare the results.
2- Change Yo to Yo = 1.  Re-run the program and compare the results.
3- Change Yo yo Yo = -1.  Re-run the program and compare the results. 

Case 4:    -1 < b < 0.    This gives rise to a convergent oscillatory time path.  
 
 Example:  Graph the time path of Yt = -.5Yt-1, Yo = 7.
  Solution:
 ClearAll[y, tulip, path6]
 tulip= RSolve[{y[t]==-.5y[t-1], y[0]==7}, y[t], t]
{{Null (y[t] -> 7.` (((-0.5`))^t))}}
 y[t_]=y[t]/.tulip[[1]];
 path6=Table[y[t], {t, 0, 10}];
 ListPlot[path6, PlotJoined->True]

⁃Graphics⁃
 
Do More Practice:
1- Change Yo to Yo = -7.  Re-run the program and compare the results.
2- Change Yo to Yo = 1.  Re-run the program and compare the results.
3- Change Yo yo Yo = -1.  Re-run the program and compare the results.  

Case 5:   b < -1.   This gives rise to a divergent oscillatory time path .
 
 Example:  Graph the time path of Yt = -2Yt-1, Yo = 7.
  Solution:
 ClearAll[y, cactus, path7]
 cactus = RSolve[{y[t]==-2y[t-1], y[0]==7}, y[t], t]
{{Null (y[t] -> 7 (((-2))^t))}}
 y[t_]=y[t]/.cactus[[1]];
 path7=Table[y[t], {t, 0, 10}];
 ListPlot[path7, PlotJoined->True]

⁃Graphics⁃
 
Do More Practice:
1- Change Yo to Yo = -7.  Re-run the program and compare the results.
2- Change Yo to Yo = 1.  Re-run the program and compare the results.
3- Change Yo yo Yo = -1.  Re-run the program and compare the results. 

Problems:
1.  Plot the values of Yt for the following difference equations   Discuss the convergence or divergence of the solutions as t ->∞.
	a)  4yt + 3yt-1 = 0,	yo = 2
	b)  yt = yt-1,		yo = 10
	c)  yt+1 + 2yt = 0,	yo = 2
	d)  2yt+1 = yt,		yo = 10
	e)  yt = -.5yt-1,		yo = 10
	
II.  First-Order Nonhomogenous Linear Difference Equations
One major property of every difference equation is that the variable Y of the equation is a function of time.  For every first-order nonhomogeneous linear difference equation such as Yt = bYt-1 + c, there is a steady-state or stationary solution for Y such that the dynamic system becomes stationary and Y has no tendency for change.  The steady-state value of Y is derived by setting the change is Y equal to zero, ΔYt = 0, or by setting Yt+1 = Yt .  That is, for the difference equation Yt = bYt-1 + c, the steady-state value of Y is Y = (c/(1  -  b)),  b ≠ 1.  

Starting from any value of Y other than the steady-state value, Y will have a time path that may be convergent to the steady state value or divergent from it.  One interesting feature of the difference equations is finding the solution of the system and analyzing the time path of Y.  

Consider the following nonhomogenous difference equation Yt = bYt-1 + c, where c is a constant.  The solution here follows the solution of difference equations of the previous section.  Equation Yt = bYt-1 + c gives rise to the following sequence of values for y,
	Y1 = bYo + c
	Y2 = bY1 +c = b(bYo + c) + c = (b^2)yo + bc + c
	Y3 = bY2 + c = b((b^2)Yo + bc + c) + c = (b^3)Yo + (b^2) c  +  bc  +  c
	          . . . . . . . . . . . . . . . . . . . . . . . . . .
	    
	Yt = bYt-1 + c = (b^t)Yo + (b^(t - 1))c + (b^(t - 2))c +  . .. . + c
	Yt = (b^t)Yo + c((b^(t - 1)) + (b^(t - 2)) +  . . . . .+ 1)
	
Since the series ((b^(t - 1)) + (b^(t - 2)) +  . . . . .+ 1) adds up to (((1 - (b^t)))/(1 - b)), the solution of the nonhomogeneous linear difference equation Yt = bYt-1 + c is Yt = (b^t)Yo + ((c (1 - (b^t)))/(1 - b)),  if b ≠ 1.  If b = 1, using the same argument, it can be proved that, Yt = Yo + ct.  
	     	
Example:  Application to Markets  
The Cobweb model of markets, which is generally true for the agricultural markets, assumes that quantity demanded at time t depends on the market price at t, Qdt = a + bPt, but the quantity supplied at time t depends on the price that prevailed in the previous period.  That is, Qst = c + dPt-1.  Finally, price adjusts in such away as to clear the market,  Qdt = Qst.  Find the steady-state equilibrium price and quantity and analyze the stability of the market under these conditions.
Solution:  Solving the three equations for the quilibrium condition yields, Pt = ((c - a)/b) + (d/b) Pt-1.  However, this equation does not reveal what the equilibrium price is.  The steady-state equilibrium condition requires that the period-by-period price remains the same, i.e., Pt = Pt-1 = Pe = ((c - a)/(b - d)).  To solve Pt = ((c - a)/b) + (d/b) Pt-1, which is a first-order non-homogeneous linear difference equation, apply the solution rule studied above.  The solution is Pt = ((((d/b)))^t) Po +  ((c - a)/b). ((1 - ((((d/b)))^t))/(1  -  (d/b))) = ((((d/b)))^t) Po + Pe(1 - ((((d/b)))^t)), or Pt = Pe + (Po - Pe)((((d/b)))^t).  
 In a conventional demand ad supply, d > 0, b < 0, and (d/b) is negative.  Therefore, the market will exhibit an oscillatory behavior.  That is If Po > Pe, price will be above the equilibrium point when t is an even number and below the equilibrium when t is an odd number.  The opposite will occur if Po < Pe.  it is also clear that the condition for market stability and convergence is -1 < (d/b) < 1.  That is, stability of the market requires that the absolute value of the slope of the demand curve exceeds the absolute value of the slope of the supply curve.
 
 Example:  Application to Markets  
 The market demand and supply of a product are given as Qdt = 100 -2Pt and Qst = -8 + 4Pt-1.  Given the initial market price of P(0) = 20, find the time path of price adjustment to its steady-state equilibrium.
 Solution:  The steady-state equilibrium is derived by setting Qdt = Qst, when Pt = Pt-1= P.  To find P using Mathematica,
 Clear[p, sumac, path8]
 Solve[100-2p==-8+4p, p]
{{Null (p->18)}}
 the steady state equilibrium price is P = 18.  To find the time path of price adjustment from P = 20 to P = 18, solve for Pt in terms of Pt-1 by setting Qdt = Qst.  This yields Pt = 54 - 2Pt-1, which is a first-order nonhomogeneous difference equation, with initial-value y(0) = 20.  To solve this equation using Mathematica, 
sumac = RSolve[{p[t]==54-2p[t-1], p[0]==20}, p[t], t]
{{p[
              t] -> -2 (-9 - 
                  (2^t) (e^(i π t)))}}
p[t_]=p[t]/.sumac[[1]];
path8=Table[p[t], {t, 0, 10}];
ListPlot[path8, PlotJoined->True]

(Null^2) (⁃  Graphics  ⁃)
 
Example:  Application to IS-LM Model  
Suppose that a dynamic macro model of an economy is expressed as Ct = bYt-1 + Co,  It = Io, and  Yt = Ct + It, where C is consumption, Y is income, and I is investment.  Analyze the time path of income.
Solution:  Solving for Yt from the equilibrium equation Yt = Ct + It yields, Yt = bYt-1 + Co + Io.  The steady-state equilibrium condition requires that the period-by-period income remain the same. i.e., Yt = Yt-1 = Ye = ((Co  +  Io)/(1  -  b)).  However, the time path of Y is determined by the solution of the difference equation Yt = bYt-1 + Co + Io, which is Yt = Yo(b^t) + (((Co  +  Io) (1  -  (b^t)))/(1  -  b))= Yo(b^t)+ ye(1 - (b^t)) = (Yo - Ye)(b^t) + Ye.  Assuming a conventional Keynesian consumption function with b < 1, the solution implies that if there is a displacement in income, the system will approach the equilibrium point gradually as t -> ∞.

Example:  Application to Macroeconomic Model  
A simple macro model of an economy is expressed as Ct = .75Yt-1 + 500, It = 200, and Yt = Ct + It.  Given that the initial income is Y(0) = 3000, graph the time path of income.   
Solution:  The equilibrium condition yields the initial-value diifference equation Yt = .75Yt-1 + 700, Y(0) = 3000.  Using Mathematica, the solution is,
ClearAll[y, catkin, path9]
catkin= RSolve[{y[t]==.75y[t-1]+700, y[0]==3000}, y[t], t]
\({{Null\ \((y[t] -> 
                1.3333333333333333`\ \((3000.`\ \((\(-2.250000000000013`\)\ 0.7500000000000006`\^t + 3.0000000000000124`\ 0.9999999999999993`\^t)\) - 
                      2300.`\ If[
                          t ≥ 
                            1, \(-3.000000000000015`\)\ 0.7500000000000006`\^t + 3.000000000000014`\ 0.9999999999999993`\^t, 
                          0])\))\)}}\)
y[t_]=y[t]/.catkin[[1]];
path9 = Table[y[t], {t, 0, 15}];
ListPlot[path9, PlotJoined->True]

⁃Graphics⁃
The time path of Y(t) shows that it approaches to its steady-state equilibrium of Y = 2800 in about 15 periods.

Do More Practice:
1- Change the MPC (the coefficient of Yt -1 in the consumption function) from .75 to .9.  Keep the other assumptions the same.  Re-run the program and compare the results.
2- Change the MPC from .75 to 1.2.  Keep the other assumptions the same.   Re-run the program and compare the results.
3- Change the initial value from Y(0) = 3000 to Y(0) = 2800.  Keep the other assumptions the same.  Re-run the program and compare the results.
4- Delete the autonomous spending of 700 in the consumption function.  Keep other assumptions the same.   Re-run the program and compare the results.


Problems:
1.  Solve the following difference equations and analyze the stability of each equation.
 	a)	Yt = Yt-1 + 5,			yo = 2
 	b)	5Yt+1 - Yt = 2,			yo = (1/2)
 	c)	3Yt + Yt-1 = 8,			yo = 10
 	d)	Yt = 2Yt-1 + 5,			yo = 5
 	e)	Yt = -3Yt-1 -10			yo = 2 		
 2.  Determine the price behavior of the following market systems and graphically show the price movements.
 	a)	Qdt = 45 - 5Pt
 		Qst = 10 + 2Pt-1
 		Po = 7 			
 	b)	Qdt = 200 - 3Pt
 		Qst = -10 + 4Pt-1
 		Po = 28 			
 	c)	Qdt = 150 -2Pt
 		Qst = 10 + 5Pt-1
 		Po = 15 	
 	d)	Qdt = 80 - 4Pt
 		Qst = -40 + 2Pt-1
 		Po =15 			
 3.  Suppose that the price adjustment in shock market follows the rule Pt+1 = Pt - .3(Qst  - Qdt).   Given the following market models find the time path of Pt and determine the stability of the models.
 	a)	Qdt = 21 - 2Pt
 		Qst = -3 + 6Pt 		         
 	b)	Qdt = 100  - 5Pt
 		 Qst = -20 + 3Pt 		          
 4.  A macro model of an economy is expressed as:  Ct = .75Yt-1 + 1000,  It = 100, and Yt = Ct + It, where C is consumption, I is investment, and Y is income.   Find the time path of income and graphically show the adjustment of income from Yo = 5000 to equilibrium. 		

III.  Second-Order, Homogenous Linear Difference Equations
The general form of a second order, homogenous linear difference equation is Yt + aYt-1 + bYt-2 = 0 where a and b are arbitrary constants.  The solution of the first order, homogenous difference equation of the previous section suggests that a solution of the form, Yt = (k^t) may be applicable here also.   Substitute Yt = (k^t) and its lagged equivalent into the difference equation Yt + aYt-1 + bYt-2 = 0 to obtain, 	(k^t)+ a(k^(t - 1)) + b(k^(t - 2)) = 0.  This equation is called the root equation or auxiliary equation.  Dividing both sides of the root equation by (k^(t - 2)) reduces the root equation to a quadratic equation (k^2) + ak + b = 0.  This equation has the roots k1 = ((-a  +  (Sqrt[(a^2)  -  4  b]))/2) and k2 = ((-a  +  (Sqrt[(a^2)  -  4  b]))/2).   Through direct substitution, it is easy to show that y1 = (k1^t),  y2 = (k2^t), or any linear combination of y1 and y2, such as Yt = c1Y1 + c2Y2 satisfies the difference equation Yt + aYt-1 + bYt-2 = 0 and is a solution of that equation.  However, depending upon the value and the sign of (a^2) - 4b in the solution of the root equation, three distinct cases can occur:

Case1:  (a^2) - 4b - 4b > 0:  Under this assumption, k1 and k2 are unequal real roots and y1 = (k1^t) or y2 = (k2^t), or any linear combination of y1 and y2 such as Yt = c1Y1 + c2Y2 = c1(k1^t)+ c2(k2^t) are solutions to the equation.

Example:  Solve the second order, homogenous difference equation Yt + Yt-1 - 6Yt-2 = 0.
Solution:  Form the auxiliary equation, (k^2) + k - 6 = 0.   This equation has two unequal real roots, k1 = 2 and k2 = -3.  Therefore, y1 = (2^t) or y2 = -(3^t), or any linear combinatin of y1 and y2 will be a solution of the equation.
Using Mathematica, the solution can be found as,
Clear[y]
RSolve[y[t]+y[t-1]-6y[t-2]==0, y[t], t]
{{Null (y[t] -> 
                (1/5) (((-3))^(2 + t)) (-C[1] + 2 C[2]) + 
                  (1/5) (2^(2 + t)) (C[1] + 3 C[2]))}}
If the initial values of  Yt at times t = 0 and t = 1 are also given, then the initial-value problem Yt + Yt-1 - 6Yt-2 = 0,  Yo = 6,  Y1 = 2 has a solution, Yt = c1(k1^t) + c2(k2^t) = c1(((2))^t) + c2(((-3))^t).   Substituting for the initial values yields the two equations, c1 + c2 = 6 and 2c1 -3c2 =2, the solution of which yields, c1 = 4 and c2 =2.  Therefore, Yt = 4(((2))^t) + 2(((-3))^t) is a solution to the initial-value problem.  It is easy to show that this solution satisfies both the initial value and the difference equation.
To solve the same initial-value problem using Mathematica, use the RSolve command as,
Clear[y]
RSolve[{y[t]+y[t-1]-6y[t-2]==0, y[0]==6, y[1]==2}, y[t], t]
{{Null (y[t] -> 
                2 ((((-3))^t) + (2^(1 + t))))}}

Case 2: (a^2) - 4b = 0:  Under this assumption, k1 = k2 = k = -(a/2) and Y1 = (k^t) is a solution to the difference equation.  However, to solve for the initial-value problem, we need to search for another solution.  It is easy to show that Y2 = t(k^t) is also a solution of the difference equation (this can easily be proved by substituting t(k^t) into the difference equation and showing that it satisfies the equation).   Therefore, any linear combination of y1 and y2, such as Yt = c1Y1 +c2Y2 = c1(k^t) + c2t(k^t) is also a solution  of the problem.

Example:  Solve the second order, homogenous difference equation Yt + 4Yt-1 + 4Yt-2 = 0 for the initial values Yo = 2, Y1 = 0.
Solution: Form the auxiliary equation, (k^2) + 4k + 4 = 0.  This equation has only one root, k = -2.   Therefore, the solution to the difference equation is y1 = (((-2))^t).    However, to solve for the initial values, we need one more solution.  We try the solution Y2 = t(((-2))^t).  Since this solution also satisfies the difference equation, a linear combination of Y1 and Y2,  Yt = c1Y1 + c2Y2 = c1(((-2))^t)(-2) + c2t(((-2))^t).
will also satisfy the difference equation.  Solving for the initial values yields, c1 = 2 and c2 = -2.  Therefore, the solution of the initial-value problem is Yt  = 2(((-2))^t) - 2t(((-2))^t).    By direct substitution, it can be shown that this solution satisfies both the inital values and the difference equation.
To solve the same problem using Mathematica,
Clear[y]
RSolve[{y[t]+4y[t-1]+4y[t-2]==0, y[0]==2, y[1]==0}, y[t], t]
\({{Null\ \((y[
                  t] -> \(-2\^\(1 + t\)\)\ e\^\(i\ π\ t\)\ \((\(-1\) + t)\))\)}}\)

Case 3:   (a^2) - 4b < 0:   Under this assumption, k1 and k2 are two distinct complex conjugate roots.   Solutions of difference equations with complex conjugate roots are discussed in chapter 14.

Example:  Application to Macroeconomics
The acceleration principle states that the net investment in the economy during the time period t (It) depends on the rate of change of aggregate income, It = v(Yt-1 - Yt-2), where v is a constant coefficient.  This principle, together with the consumption model Ct = bYt-1can be merged into a simple multiplier-accelerator model of the economy.  Solving the two equations above for the equilibrium condition Yt = Ct + It, yields Yt = (b+v)Yt-1 -vYt-2.  This equation is a second-order, homogeneous difference equation that can be solved to find the time path of Yt.  The solution and the stability of  the model follows the three cases discussed above.

Example:  Application to Macroeconomics
Suppose that in the example above investment adjusts to change in income at a constant proportion 8, (v = 8), and that marginal propensity to consume,  b = .9.   Given the initial values Yo = 6000 and Y1 = 6200, find the time path of Y. 
Solution:  The assumptions above yields the second-order difference equation Yt = 8.9Yt-1 - 8Yt-2 with the initial values Yo = 6000 and Y1 = 6200.  Using Mathematica, the solution and the time path of Yt can be found as:
ClearAll[y, lily, graph1]
lily = RSolve[{y[t]==8.9y[t-1]-8y[t-2], y[0]==6000, y[1]==6200}, y[t], t]
\({{Null\ \((y[t] -> 
                5983.574322551905`\ 1.014523322739623`\^t + 
                  16.42567744809575`\ 7.885476677260377`\^t)\)}}\)
y[t_]=y[t]/.lily[[1]];
graph1 = Table[y[t], {t, 0, 5}];
ListPlot[graph1, PlotJoined->True]

(Null^2) (⁃  Graphics  ⁃)

Do More Practice:
1- Change b to .75.  Keep other assumptions the same.  Re-run the program and compare the results.
2- Change v to 2.  Keep other assumptions the same.  Re-run the program and compare the results.
3- Change v to .7 and b to .9.  Keep other assumptions the same.  Re-run the program and compare the results.
4- Change v to .9 and b to .7.  Keep other assumptions the same.  Re-run the program and compare the results.

Problems:
1.   Solve the following difference equations and graph the time path of Yt..
	a)	Yt + 3Yt-1 + 2Yt-2 = 0,			Yo =2,	Y1 =3
	b)	Yt + Yt-1 - 2Yt-2 = 0,			Yo = 1,	Y1 = -1
	c)	Yt - 4Yt-2 = 0,				Yo = 1,	Y1 = 1
	d)	6Yt + 5Yt-1 + Yt-2 = 0,			Yo = 0,	Y1 = 1
	e)	9Yt - 12Yt-1 + 4yT-2 = 0,		Yo =1,	Y1 = 0
	f)	Yt + 8Yt-1 + 16Yt-2 = 0,		Yo = 2,	Y1 = -4		
2.   A simple Keynesian macro-model of the economy is described as:  It = 4.07(Yt-1 - Yt-2), Ct = .73Yt-1, and Yt = Ct + It where I is investment, C is consumption and Y is income.   Given that Yo = 4000 and y1 = 4150, find the time path of income and graph it.

IV.  Second Order, Nonhomogenous Linear Difference Equations
The solution of a second-order, nonhomgenous linear difference equation of the form Yt + aYt-1 + bYt-2 + c = 0 consists of two parts:  The solution of the complementary homogenous part, Yc, and the particular solution, Yp, of the nonhomogenous difference equation.  The general solution of the difference equation is the sum of the solution to the complementary homogenous equatiion and the particular solution.   That is, Yt = Yc + Yp.  To obtain the particular solution of the nonhomogeneous part, try a function f(t) = z such that z satisfies the difference equation for all t values.

Example:  Solve the second-order, nonhomogenous linear difference equation Yt - 2Yt-1 - 3Yt-2 + 8 = 0 for the initial values Yo = 3 and Y1 = 2. 
Solution:  The complementary homogenous equaion Yt - 2Yt-1 - 3Yt-2 = 0 has the auxiliary equation (k^2) - 2k - 3 =0, with the roots k1 = 3 and k2 = -1.  Therefore, the solution of the complementary homogenous part is Yc = c1(((3))^t) + c2(((-1))^t).  To obtain a particular solution, we try f(t) = z, such that z satisfies the difference equation, for all t.  That is, z - 2z - 3z + 8 =  0,   or  z = 2.  Therefore, the particular solution is Yp = 2.  The general solution of the difference eqaution is Yt = Yc + Yp = c1(((3))^t) + c2(((-1))^t)+ 2.   Given the initial values, c1 and c2 can be found from the general solution of the difference equation.  Substitute for Yo and y1 to get  3 = c1 + c2 + 2 and 2 = 3c1 - c2 + 2, respectively.  Solving these two equations together yields, c1 = (1/4)and c2 = (3/4).   The solution of the initial-value, nonhomogeneous difference equation is Yt = (1/4) (3^t) + (3/4)(((-1))^t) + 2.  It can be shown, through direct substitution, that this solution satisfies both the initial conditions and the difference equation.
To solve the same problem using Mathematica,
Clear[y]
RSolve[{y[t]-2y[t-1]-3y[t-2]+8==0, y[0]==3, y[1]==2}, y[t], t]
{{Null (y[t] -> 
                (1/4) (8 + 3 (((-1))^t) + (3^t)))}}

Example:  Application to Macroeconomics
The multiplier-accelerator model of the previous section can be extended to a second order nonhomogenous difference equation by including autonomous spending and/or by including govenrment spending.   Suppose that the investment at time t ( It ) is defined as a constant 4 proportion of the change in income, It = 4(Yt-1 - Yt-2); Consumption at time t is defined as Ct = .80Yt-1 + 400, and the government spending is Gt = 200.  Given the equilibrium condition Yt = Ct + It + Gt and the initial values Yo = 5000 and Y1 = 5100 find Yt and graph the time path of Yt. 
Solution:  Solving the equations of the model for the equilibrium condition yields Yt - 4.8Yt-1+ 4Yt-2 - 600 = 0 with the initial values Yo = 5000 and Y1 = 5100.  This is an initial-value, second-order, nonhomogeneous difference equation that can be solved for Yt.  Using Mathematica the solution is, 
ClearAll[y, fuchsia, graph2]
fuchsia = RSolve[{y[t]-4.8y[t-1]+4y[t-2]-600==0, y[0]==5000, y[1]==5100}, y[t], t]
\({{Null\ \((y[
                  t] -> \(-0.25`\)\ \((5000.`\ \((\(-19.99999999999723`\)\ 0.9999999999999956`\^t + 
                            23.67858991804094`\ 1.0733500838578451`\^t - 
                            7.678589918043711`\ 3.726649916142159`\^t)\) + 
                      If[t ≥ 1, 
                        87999.9999999865`\ 0.9999999999999956`\^t - 
                          126463.35274200485`\ 1.0733500838578451`\^t + 
                          38463.35274201838`\ 3.726649916142159`\^t, 
                        0])\))\)}}\)
y[t_]=y[t]/.fuchsia[[1]];
graph2=Table[y[t], {t, 0, 5}];
ListPlot[graph2, PlotJoined->True]

⁃Graphics⁃

Do More Practice:
1- Change b to .75.  Keep other assumptions the same.  Re-run the program and compare the results.
2- Change v to 2.  Keep other assumptions the same.  Re-run the program and compare the results.
3- Change v to .7 and b to .9.  Keep other assumptions the same.  Re-run the program and compare the results.
4- Change v to .9 and b to .7.  Keep other assumptions the same.  Re-run the program and compare the results.

Problems:
1.   Solve the following difference equations and graph the time path of Yt.
	a)  Yt - Yt-1 - 2Yt-2 + 4 = 0,		Yo =0, Y1 = 1
	b)  Yt - 4Yt-2 + 6 = 0,			Yo = 10, Y1 = 2
	c)  Yt+2 - 7Yt+1 + 10Yt + 12 = 0,	Yo = 10,  Y1 = 5
	d)  Yt+2 - 9Yt + 8 = 0			Yo = 2, Y1 = 3
	e)  Yt + 4Yt-1 + Yt-2 =18,		Yo = 1, Y1 = 2
	 f)  Yt - 4Yt-1 + 4Yt-2 -3 = 0,		Yo = 0, Y1 = 2		
2.  A simple Keynesian macro-model of an economy is described as It = 5.08(Yt-1 - Yt-2), Ct = .8Yt-1 + 200, Yt = Ct + It + Gt.  Given that Yo = 3000 and Y1 = 3200, find the time path of income when Go = 150.
		
Stability Analysis of the solutions of Second-Order, Linear Difference Equations:  The general solution of a second-order, nonhomogenous difference equation is Yt = c1(k1^t) + c2(k2^t) + c, where c1, c2 and c are constants and k1 and k2 are the roots of auxiliary equation corresponding to the homogenous part of the difference equation.  This solution generates a time path for Yt as t varies from 0 to ∞.  In this solution, the only values that change as t changes are (k1^t) and (k2^t).  Therefore, changes in Yt depend on the sign and absolute values of k1 and k2.

A common way to study the time path of Yt is to transform Yt = c1(k1^t) + c2(k2^t) + c to Yt = (k1^t)(c1 + c2((((k1/k2)))^t)) + c by factoring out the root that has the greater absolute value.   Since it is assumed that |k1| >  |k2|, ((((k1/k2)))^t) approaches zero as t approaches ∞.   Therefore, the time path of xYt can be determined by concentrating on the behavior of (k1^t), or the greater root only.

If k1 and k2 are real roots, (k1≠ k2), Yt will have a number of different time paths depending upon the sign and magnitude of k1 and k2.  In general, if the absolute value of the largest root is smaller than one, as t approaches ∞,  Yt will approach its long run equilibrium value.  Otherwise, Yt will move further away from equilibrium.  Some of those cases are studied in the following examples.

Example:  The case where k1 and k2 > 1.
Solve the difference equation Yt - 5Yt-1 + 6Yt-2 - 10 = 0, Yo = 8 and Y1 = 13.
Solution:  The solution of the difference equation is, Yt = (2^(t ))+ 2(((3))^t) + 5.  Factor out (3^t) to transform the solution to (3^t)(((((2/3)))^t) + 2) + 5.  As t approaches ∞, ((((2/3)))^t) approaches zero, while (3^t) approaches ∞.  Therefore, Yt will explode gradually in a positive direction.  The solution and the graph using Mathematica illustrate the time apth of Yt.
ClearAll[y, corymb, graph3]
corymb=RSolve[{y[t]-5y[t-1]+6y[t-2]-10==0, y[0]==8, y[1]==13}, y[t], t]
{{Null (y[t] -> 5 + (2^t) + 2 (3^t))}}
y[t_]=y[t]/.corymb[[1]];
graph3=Table[y[t], {t, 0, 10}];
ListPlot[graph3, PlotJoined->True]

⁃Graphics⁃

Example:  The case where 0 < k1 &  K2 < 1.
Solve the difference equaion 6Yt - 5Yt-1 + Yt-2 - 4 = 0, Yo = 10, Y1 = 5 and analyze the time path of Yt.
Solution:  The solution of the difference equation is Yt = 2((((1/2)))^(t ))+ 6((((1/3)))^t) + 2.   Factor out ((((1/2)))^t) to transform the solution to Yt = ((((1/2)))^t)(2 + 6((((2/3)))^t)) + 2.   As t approaches ∞, both ((((1/2)))^t) and ((((2/3)))^t) approach zero.   Therefore, Yt gradually approaches its long run equilibrium level, 2.  The solution and the graph using Mathematica illustrate the time path of Yt.
ClearAll[y, dafodil, graph4]
dafodil=RSolve[{6y[t]-5y[t-1]+y[t-2]-4==0, y[0]==10, y[1]==5}, y[t], t]
{{Null (y[t] -> 
                (2^(1 - t)) (3^(-t)) (3 (2^t) + (3^t) + (6^t)))}}
y[t_]=y[t]/.dafodil[[1]];
graph4=Table[y[t], {t, 0, 5}];
ListPlot[graph4, PlotJoined->True]

⁃Graphics⁃

Example:  The case where k1 & k2 < -1.
Solve the difference equation Yt + 5Yt-1 + 6Yt-2 - 24 = 0, Yo = 5, Y1 = -5.
Solution:  The solution of the difference equaion is Yt = 2(((-2))^t) + (((-3))^t) + 2.   The solution can be written as Yt = (((-3))^t)(1 + ((((2/3)))^t)) + 2.   As t approaches ∞,  ((((2/3)))^t) approaches zero and (((-3))^t) approaches ∞.  Therefore, Yt will explode with an oscillatory path.  The solution and the graph using Mathematica illustrate the time path of Yt.
ClearAll[y, mint, graph5]
mint=RSolve[{y[t]+5y[t-1]+6y[t-2]-24==0, y[0]==5, y[1]==-5}, y[t], t]
{{Null (y[t] -> 
                2 + (((-3))^t) + (((-1))^t) (2^(1 + t)))}}
y[t_]=y[t]/.mint[[1]];
graph5=Table[y[t], {t, 0, 10}];
ListPlot[graph5, PlotJoined->True]

⁃Graphics⁃

Example:  The case where -1  <  k1 &  K2  < 0.
Solve the difference equation 6Yt + 7Yt-1 + 2Yt-2 + 30 = 0, Yo = 3, Y1 = -5, and analyze the time path of Yt.
Solution:  The solution of the difference equatiion is Yt = 3(((((-2)/3)))^t) + 2(((((-1)/2)))^t) - 2.  Factor out (((((-2)/3)))^t) to get Yt = (((((-2)/3)))^t)(3 + 2((((3/4)))^t))  - 2.   As t approaches ∞, both (((((-2)/3)))^t) and ((((3/4)))^t) approach zero and Yt approaches its long run equilibrifum, -2.  However, since the roots are negative numbers, Yt approaches -2 with an oscillatory path around the equilibrium value.  The solution and the graph using Mathematica illustrate the time path of Yt.
ClearAll[y, poppy, graph6]
poppy=RSolve[{6y[t]+7y[t-1]+2y[t-2]+30==0, y[0]==3, y[1]==-5}, y[t], t]
\({{Null\ \((y[t] -> 
                3\ \((\(-\(2\/3\)\) + \((\((2\/3)\)\^t + 
                            2\^\(1 - t\)\/3)\)\ e\^\(i\ π\ t\))\))\)}}\)
y[t_]=y[t]/.poppy[[1]];
graph6=Table[y[t], {t, 0, 10}];
ListPlot[graph6, PlotJoined->True]

⁃Graphics⁃

Example:  The case where k1 = k2  > 1.
solve the difference equation Yt - 4Yt-1 + 4Yt-2 -2 = 0, Yo = 3, Y1 = 10 and plot the time path of Yt.
Solution:  The solution of the difference equation is Yt = (2^t) + 3t(2^t) + 2.  The solution can be written as Yt = (2^t)(1 + 3t) + 2.   As t approaches ∞, Yt will gradually explode in the positive direction.  The solution and the graph using Mathematica illustrate the time path of Yt.
ClearAll[y, basil, graph7]
basil=RSolve[{y[t]-4y[t-1]+4y[t-2]-2==0, y[0]==3, y[1]==10}, y[t], t]
{{Null (y[t] -> 2 + (2^t) + 3 (2^t) t)}}
y[t_]=y[t]/.basil[[1]];
graph7=Table[y[t], {t, 0, 10}];
ListPlot[graph7, PlotJoined->True]

⁃Graphics⁃

Example:  The case where 0  <  k1= k2  < 1.
Solve the difference equaiion 4Yt - 4Yt-1 + Yt-2 + 2 = 0, Yo = 4, Y1 = 5, and graph the time path of Yt..
Solution: The solution of the equation is Yt = 2(((((1/2)))^t) + 4t((((1/2)))^t)) + 2.  The solution can be written as Yt = ((((1/2)))^t)(2 + 4t) + 2.   As t approaches ∞,  Yt will approach its long run equilibrium, -2.  The solution and the graph using Mathematica illustrate the time path of Yt.
ClearAll[y, bay, graph8]
bay=RSolve[{4y[t]-4y[t-1]+y[t-2]+ 2==0, y[0]==4, y[1]==5}, y[t], t]
{{Null (y[t] -> (2^(1 - t)) (3 - (2^t) + 4 t))}}
y[t_]=y[t]/.bay[[1]];
graph8=Table[y[t], {t, 0, 10}];
ListPlot[graph8, PlotJoined->True]

⁃Graphics⁃

Example:  The case where  -1 < k1=k2 < 0.
Solve the difference equation 4Yt + 4Yt-1 + Yt-2 - 9 = 0, Yo = 3, Y1 = -2, and graph the time path of Yt.
Solution:  The solution of the difference equation is Yt = 2(((((-1)/2)))^t) + 4t(((((-1)/2)))^t) + 1.   The solution can be written as Yt = (((((-1)/2)))^t)(2 + 4t) + 1.    As t approaches ∞, (((((-1)/2)))^t) approaches zero and Yt approaches its long run equilibrium, 1.  The solution and the graph using Mathematica illustrate the time path of Yt.
ClearAll[y,thyme, graph9]
thyme = RSolve[{4y[t]+4y[t-1]+y[t-2]-9==0, y[0]==3, y[1]==-2}, y[t], t]
{{y[t] -> 
            (2^(-t)) (2 (((-1))^t) + (2^t) + 
                  4 (((-1))^t) t)}}
y[t_]=y[t]/.thyme[[1]];
graph9=Table[y[t], {t, 0, 10}];
ListPlot[graph9, PlotJoined->True]

⁃Graphics⁃

Example:  Application to Markets
The demand and supply for a product are expressed as Qdt = 20 - Pt and Qst = -4 + 1.1Pet, respectively.  Pt is the market price of the product and Pet is the expected price.  Given that expectations of future prices by the producers follow the simple adaptive expectations rule, Pet = (9/11)Pt-1 + (2/11)Pt-2, solve the market model and analyze the time path of Pt, if the initial prices are Po = 10 and P1 = 8.
Solution:  The market equilibrium requires that Pt = 24 - .6Pet.  Substituting for Pet yields the second-order difference equation Pt = 24 - .9Pt-1 - .2Pt-2.  Given the initial values Po = 10 and P1 = 8, this equation can be solved to trace its time path.  Using Mathematica,
ClearAll[p,fennel, tplot1]
fennel = RSolve[{p[t]==24-.9p[t-1]-.2p[t-2], p[0]==10, p[1]==8}, p[t],t]
\({{Null\ \((p[
                  t] -> \(-5.`\)\ \((1.3333333333333326`\ \((\(-0.5000000000000001`\))\)\^t - 1.7142857142857135`\ \((\(-0.4`\))\)\^t - 
                      1.6190476190476193`\ 1.`\^t + 
                      7.`\ If[
                          t ≥ 
                            2, \(-1.3333333333333324`\)\ \((\(-0.5000000000000001`\))\)\^t + 1.4285714285714277`\ \((\(-0.4`\))\)\^t - 
                            0.09523809523809523`\ 1.`\^t, 0])\))\)}}\)
p[t_]=p[t]/.fennel[[1]];
tplot1=Table[p[t], {t, 0, 10}];
ListPlot[tplot1, PlotJoined->True]

⁃Graphics⁃

Do More Practice:
1-  Suppose that the expectations of the future prices follows the rule Pet = .5Pt + .3Pt-1 + .2Pt-2.  R-run the program and compare the results.
2-  Suppose that demand and supply are given as Qdt = 20 - 1.2Pt and Qst = -4 + Pet, respectively.  The other assumptions stay the same as in example.  Re-run the program and compare the rsults. 

Problems:
1.   Solve the following second-order linear difference equations and graph the values of 	Yt for t = 0 to t = 10.  Analyze the time path of Yt as t approaches ∞.   
	a)	4Yt - Yt-2 -9 = 0,			Yo = 11,Y1 = 1
	b)	Yt + Yt-1 - 6Yt-2 + 4 = 0,		Yo =6,	Y1 = 1
	c)	Yt + 4Yt-1 + 4Yt-2 -18 = 0,		Yo = 3,	Y1 = -4
	d)	10Yt - 23Yt-1 + 12Yt-2 + 3 = 0,	Yo = 6,	Y1 = 6.1
	e)	2Yt - Yt-1 - 6Yt-2 + 4 = 0,		Yo = 6,	Y1 = 4
	f)	2Yt + 7Yt-1 + 6Yt-2 + 30 = 0,	Yo = .5,	Y1 = -6
	g)	6Yt + Yt-1 -Yt-2 + 12 = 0,		Yo = 3,	Y1 = -(17/6)
	h)	4Yt - 12Yt-1 + 9Yt-2 + 2 = 0,		Yo = 0,	Y1 = 4		
2.   The demand and supply equations of a product are given as, Qdt = 10 - Pt and Qst = -8 + .8Pet  where, Pt is the market price and Pet is the expected price.  Given that expectations of future prices by the producers follow the simple adaptive expectations model, Pet = (7/8)Pt-1 + (1/8)Pt-2, solve the market model and analyze the time path of Pt, if the initial prices are Po = 15 and P1 = 8.1.		
3.  The demand and supply equations for a product are given as Qdt = 30 - 2Pt and Qst = -20 + 8Pet where, Pet is the expected market price.  If the expectations of future prices by the producers follow the simple adaptive expectations rule, Pet = (5/8)Pt-1 + (3/8)Pt-2, solve the market model. given the initial prices Po = 6.5 and P1 = 3.   Analyze the time path of prices.  Will the market price approach its long run equilibrium?
4.  Consider a simple macro model of economy where Yt = Ct + It, Ct = 100 + .8Yt-1, and It = 10(Yt-1 - Yt-2).  The first equation is the equilibrium condition and expresses income at time t (Yt) equal to consumption at time t (Ct) and investment at time t {It).  The second equation expresses consumption as a function of lagged income and the third equation expresses investment as a function of the change in lagged income.  Given that the initial income levels are Yo =  450 and Y1 = 480, find the time path of Yt and analyze the stability of the system.

V.  System of First-Order Linear Difference Equations
Economic variables may be related to each other through a system of several discrete time dynamic inter-related variables and equations.  For example, a discrete time dynamic IS-LM model of the economy may be expressed by the system of equations ΔYt  = α[D(Y, r) - Y] and Δrt = β[L(Y,r) - M], where Y is the the income or output, r is the interest rate, M is the money supply, and α and β are the coefficients of adjustment.  The first equation represents equilibrium in the output market, IS.  It expresses changes in output, ΔYt = Yt+1 - Yt, as α proportion of the gap between aggregate demand, D(Y, r), and aggregate output, Y.  The second equation,  (dr/dt),represents the equilibrium in the money market, LM.  It expresses changes in the interest rate as a β propostion of the gap between the demand for money, L(Y,r), and the money supply, M.

The laws that govern the behavior of variables in this model and the system can be estimated by a system of linear difference equations.  That is, a model which assumes that the change of each variable is a linear combination of other variables.  

System of First-Order Linear Homogeneous Difference Equations:  The general form of a system of first-order linear homogeneous difference equation is ΔYt = AYt where, ΔY(t) is a vector with elements (({{ΔY1t}, {ΔY2t}})), A is a matrix of constant coefficients, and Y(t) is the vector with elements (({{Y1t}, {Y2t}})).

Equation ΔYt = AYt may be written as (({{Y1t + 1  -  Y1t}, {Y2t + 1  -  Y2t}}))= A(({{Y1t}, {Y2t}})) or (({{Y1t+1 }, {Y2t + 1}}))= A(({{Y1t}, {Y2t}}))+ I2(({{Y1t}, {Y2T}})) = (A + I2)(({{Y1t}, {Y2t}})), where I2 is a 2x2 identity matrix.   Representing (({{Y1t+1 }, {Y2t + 1}}))by Yt+1 and (({{Y1t}, {Y2t}}))by Yt, the model can be written in matrix form as Yt+1 =  (A + I2)Yt.   The resemblance of this equation to Yt+1 = kYt suggets that a trial solution of the form Yt+1 = ((Yok)^t) may be applicable to the system of difference equations.  Let us study the possibility by searching for a constant vector u and a number k such that, Yt+1 = u(k^t) satisfies Yt+1 =  (A + I2)Yt.  Since Yt+1 = u(k^t), for Y to be a solution it is necessary and sufficient that u(k^t)= (A + I2)u(k^(t - 1)).   This equation implies that, (A + I2)u - ku = ((A + I2) - kI)u = 0.  Therefore, k is the eigenvalue and u is the eigenvector of the matrix (A + I2) (see, Chapter 5, Special Topics in Matrix Algebra).
	
Example:  Solve the following system of linear homogeneous difference equations, ΔY1t = 2Y1t + Y2t and ΔY2t = 3Y1t, Y1o  = 1 Y2o = 1.
Solution:  The system in the matrix form is, (({{Δy1t}, {Δy2t}}))= (({{2, 1}, {3, 0}}))(({{y1t}, {y2t}})).  Find the eigenvalues and eigenvectors of A + I2= (({{3, 1}, {3, 1}})),
Eigensystem[(({{3, 1}, {3, 1}}))]
{{0,4},{{-1,3},{1,1}}}
the system has eigenvalues 0 and 4, with the corresponding eigenvectors (({{-1}, {3}})) and (({{1}, {1}})), respectively .   Therefore, the general solution of the system is (({{Y1t + 1}, {Y2t + 1}})) = c1(({{-1}, {3}}))(0^t) + c2(({{1}, {1}}))(4^t) , \(where\ c1\ and\ c2\ are\ two\ constants\ to\ be\ decided\ by\ the\ initial\ conditions . \ \ Substitution\ for\ the\ initial\ values\ Y1o\  = \ \(Y2o\  = \ \(1\ results\ in\ Y1t + 1\  = \ \(4\^t\ and\ Y2t + 1\  = \ 
                  4\^t . \ \ To\ solve\ the\ same\ problem\ using\ Mathematica\ type\)\)\)\),
{{0,4},{{-1,3},{1,1}}}
Clear[y1, y2]
RSolve[{y1[t+1]==3y1[t]+y2[t], y2[t+1]==y2[t]+3y1[t], y1[0]==1, y2[0]==1}, {y1[t], y2[t]}, t]
{{y1[t] -> (4^t), y2[t] -> (4^t)}}

Example:  Solve the system of difference equations ΔYt = (({{1, 3}, {0, 1}}))Yt for the initial values Y0 = (({{1}, {1}})).  Plot the time path of Y1, Y2,  and the relationship between y1 and y2 over time.
Solution:  A  + I2 matrix of the system is (({{2, 3}, {0, 2}})).  The system can be written as (({{Y1t + 1}, {Y2t + 1}})) = (({{2, 3}, {0, 2}}))(({{Y1t}, {Y2t}})).  Using the Mathematica command RSolve the problem can be solved as,
Clear[y1,  y2,  parsley,  oregano,  tplot2,  tplot3]
A=(({{1, 3}, {0, 1}}));
H[{y1_  ,  y2_}] = (A + IdentityMatrix[2]) . {y1,  
                y2}
{2 y1+3 y2,2 y2}
oregano= RSolve[{y1[t+1]==2y1[t]+3y2[t], y2[t+1]==2y2[t], y1[0]==y2[0]==1}, {y1[t], y2[t]}, t] 
{{y1[t] -> (2^(-1 + t)) (2 + 3 t), 
          y2[t] -> (2^t)}}
to find the time path of Y1 and Y2,
y1[t_]=y1[t]/.oregano[[1]];
tplot2=Table[y1[t], {t, 0, 6}];
ListPlot[tplot2, PlotJoined->True]

⁃Graphics⁃
to find the time path of y2, 
y2[t_]=y2[t]/.oregano[[1]];
tplot3=Table[y2[t], {t, 0, 6}];
ListPlot[tplot3, PlotJoined->True]

⁃Graphics⁃
The two graphs above represent the time path of y1 and y2 individually.  To graph the relationship between y1 and y2 in the y1-y2 space, use the "NestList" command of  Mathematica, as shown below.  This graph shows the dynamic relationship between y1 and y2 over time. 
f[y1_, y2_, t_]:=NestList[H, {y1, y2}, 10];
parsley= f[1, 1, 10];
ListPlot[parsley, PlotStyle->PointSize[.02]]

⁃Graphics⁃

 Do More Practice:  
1- Change the coefficients of the A matrix in the preceding example to (({{1, 0}, {1, 4}})), find the eigenvalues of the A + I2 matrix, and re-run the program.  Observe the changes in the eigenvalues and their relations to the solution.
2- Change the coefficients of the A matrix in the preceding example to (({{0, 1}, {3, 4}})), find the eigenvalues of the A + I2 matrix, and re-run the program.  Observe the changes in the eigenvalues and their relations to the solution.
3- Change the coefficients of the A matrix in the preceding example to (({{1, -1}, {4, 4}})), find the eigenvalues of the A + I2 matrix, and re-run the program.  Observe the changes in the eigenvalues and their relations to the solution.

Problems:
1.  For the  following matrices, find a solution to ΔYt = AYt.
a)  A =  (({{1, -3}, {0, 1}}))		b) A = (({{2, 0}, {1, 0}}))

c)  A = (({{1, 2}, {2, 1}}))			d)  A = (({{2, 3}, {1, 0}}))

e)  A = (({{0, 1}, {(5/3), (2/3)}}))		f)  A = (({{1, 4, 1}, {0, 2, 3}, {0, 0, 3}}))
		                 
g)  A = (({{1, 0, 1}, {0, 2, 1}, {4, 0, 1}}))		h)  A = (({{1, 1, 0}, {0, 0, 1}, {0, 1, 1}}))
2.  Solve problems a, b, c, and d of the problem 1 for the initial value, y(0) = (({{1}, {2}})). Plot the time path of Y1t and Y2t.
3.  Solve problems a and b of  the problem 1 for the inital value,  y(0) = (({{0}, {3}})).  Plot the time path of Y1t and Y2t.

System of First-Order Nonhomogeneous Linear Difference Equations:  The general form of a system of first-order nonhomogeneous linear difference equations is ΔYt = AYt + B where, ΔYt = AYt is the homogenous part of the system and B is a vector of constants.  If the A matrix is non-singular matrix, then the steady-state equilibrium for the system exists and can be derived by setting ΔYt = 0 and solving for Y1t and Y2t.  That is, (({{Y1t}, {Y2t}})) = (A^(-1)). B.  The solution of the system, however, follows the solution of the complementary homogeneous part, discussed in the previous section, plus the particular solution due to nonhomogenous part.

Example:   Solve the following system of nonhomogenous linear difference equations, ΔY1t = -.4Y1t + .2Y2t + 2.8 and ΔY2t = .8Y1t - .7Y2t - 3.8, Y1o  = Y2o = 8.  Plot the time path of Y1t, Y2t, and the relationship between Y1t and Y2t over time.
Solution:  The system in the matrix form is, (({{Δy1t}, {Δy2t}}))= (({{- .4, .2}, {.8, - .7}}))(({{y1t}, {y2t}}))+ (({{2.8}, {-3.8}})).   The steady-state equilibrium condition can be derived by setting (({{Δy1t}, {Δy2t}}))= (({{0}, {0}})), and solving for (({{y1t}, {y2t}})).  Using Mathematica,
Clear[y1,  y2, paprika,  rice,  tplot4,  tplot5,  
          tplot6]
A = (({{- .4, .2}, {.8, - .7}}));    B = (({{2.8}, {-3.8}}));   y = (({{y1t}, {y2t}}));
Solve[A . y + B == 0,  {y1t,  y2t}]
{{y1t->10.,y2t->6.}}
the steady-state equilibrium of the system is at the Y1 and Y2 values of  10 and 6, respectively.  To solve the difference equation, find A + I2 matrix, and solve for (({{Y1t + 1}, {Y2t + 1}}))in terms of (({{Y1t}, {Y2t}})).  That is, (({{Y1t + 1}, {Y2t + 1}})) = ( A + I2)(({{Y1t}, {Y2t}}))+ B.   Using Mathematica, 
H[{y1_  , y2_}]=(A+IdentityMatrix[2]).{y1, y2}+B
{{2.8+0.6 y1+0.2 y2},{-3.8+0.8 y1+0.3 y2}}
use the resulting equtions to write the "RSolve" command as,
paprika =RSolve[{y1[t+1]==2.8 + .6y1[t] + .2y2[t] , y2[t+1]== -3.8  + .8y1[t]+.3y2[t], y1[0]==y2[0]==8}, {y1[t], y2[t]}, t]
\({{y1[
              t] -> \(-8.`\)\ \((\(-1.0000000000000002`\)\ If[t == 0, 1, 
                        0] + If[
                      t ≥ 
                        1, \(-0.0006226121454243346`\)\ 0.0227998127341234`\^t + 7.333955945478759`\ 0.8772001872658766`\^t - 
                        8.333333333333336`\ 1.`\^t, 0])\) + 
              5.999999999999999`\ \((\(-1.0000000000000002`\)\ If[t == 1, 1, 
                        0] + If[
                      t ≥ 
                        2, \(-0.027307774528011825`\)\ 0.0227998127341234`\^t + 8.360641107861348`\ 0.8772001872658766`\^t - 8.333333333333336`\ 1.`\^t, 
                      0])\) + 
              0.8000000000000003`\ \((\(-1.0000000000000002`\)\ If[t == 2, 1, 
                        0] + If[
                      t ≥ 
                        3, \(-1.1977192464893176`\)\ 0.0227998127341234`\^t + 
                        9.531052579822653`\ 0.8772001872658766`\^t - 
                        8.333333333333336`\ 1.`\^t, 0])\), 
          y2[t] -> 
            50.00000000000009`\ \((\(-0.005471168060446508`\)\ 0.022799812734123404`\^t + 0.5321378347271122`\ 0.8772001872658766`\^t - 
                  0.36666666666666603`\ 1.`\^t - 
                  2.92`\ If[
                      t ≥ 
                        2, \(-0.023954384929786306`\)\ 0.022799812734123404`\^t + 0.19062105159645268`\ 0.8772001872658766`\^t - 
                        0.1666666666666664`\ 1.`\^t, 0])\)}}\)
to graph the time path of Y1,
y1[t_]=y1[t]/.paprika[[1]]
tplot4=Table[y1[t], {t, 0, 20}];
ListPlot[tplot4, PlotJoined->True] 
\(\(-8.`\)\ \((\(-1.0000000000000002`\)\ If[t == 0, 1, 0] + 
              If[t ≥ 
                  1, \(-0.0006226121454243346`\)\ 0.0227998127341234`\^t + 
                  7.333955945478759`\ 0.8772001872658766`\^t - 
                  8.333333333333336`\ 1.`\^t, 0])\) + 
        5.999999999999999`\ \((\(-1.0000000000000002`\)\ If[t == 1, 1, 0] + 
              If[t ≥ 
                  2, \(-0.027307774528011825`\)\ 0.0227998127341234`\^t + 
                  8.360641107861348`\ 0.8772001872658766`\^t - 
                  8.333333333333336`\ 1.`\^t, 0])\) + 
        0.8000000000000003`\ \((\(-1.0000000000000002`\)\ If[t == 2, 1, 0] + 
              If[t ≥ 
                  3, \(-1.1977192464893176`\)\ 0.0227998127341234`\^t + 
                  9.531052579822653`\ 0.8772001872658766`\^t - 
                  8.333333333333336`\ 1.`\^t, 0])\)\)

⁃Graphics⁃
it is obvious that as t -> ∞, Y1 converges to its steady-state equilibrium 10.  To graph the time path of Y2t,
y2[t_]=y2[t]/.paprika[[1]];
tplot5=Table[y2[t], {t, 0,20}];
ListPlot[tplot5, PlotJoined->True] 

⁃Graphics⁃
it is obvious that as t -> ∞, Y2 approaches its steady-state equilibrium 6.  To graph the first 100 points visited over time in adjustment towards the equilibrium point, use the following commands,
H[{y1_, y2_}]={2.8 +.6y1+.2y2, -3.6 + .8y1  + .3y2};
f[y1_, y2_, t_]:=NestList[H, {y1, y2}, 100];
rice = f[1, 1, 100];
tplot6 =ListPlot[rice, PlotStyle->PointSize[.02]]

⁃Graphics⁃
It is also possible to graph the two equations of the system in the y1-y2 space and study the dynamics of the adjustment over time from one point to another or the dynamics of the adjustment to steady-state equilibrium.  To do this, find y2 in terms of y1 for the two equations of the system by setting ΔY1t = ΔY2t = 0 and plot them together with the previous graph, as shown below,
Clear[y1, y2, eq1, eq2, tplot7]
eq1=Solve[-.4y1 + .2y2 +2.8==0, y2]
{{y2->-5. (2.8-0.4 y1)}}
eq2=Solve[.8y1 - .7y2 - 3.8==0, y2]
{{y2->1.42857 (-3.8+0.8 y1)}}
tplot7=Plot[{-5(2.8 - .4y1), 1.42857(-3.8 + .8y1)}, {y1, 0, 20}]
Show[tplot6, tplot7] 


⁃Graphics⁃

⁃Graphics⁃

Example:  Application to IS-LM Model 
A simple macro model of an economy is represented by the following equations: C = 400 + .8Y,  I = 100 - 8r, Y = C + I,  Md = .25Y - .9r, and Md = Ms = 500.  The first equation expresses consumption (C) as a function of income (y); the second equation expresses investment ( I ) as a function of interest rate (r); the third equation is the equilibrium condition in the goods market; the fourth equation expresses demand for money (Md) as a function of income and interest rate; and the last equation is the equilibrium condition in the money market.  Suppose that the change in income in the goods market and the change in interest rate in the money market are proportional to excess demand in the corresponding markets, with a proportionality factor of .4 for the goods market and a proportionality factor of .5 for the money markets.  Given that the initial values of y and r are Y(0) = 2200 and r(0) = 10, find the time path of Y and r and analyze the dynamic relatioship between Y and r over time.  
Solution:  The change in income (ΔYt) can be expressed in terms of the excess demand in the goods market as ΔYt = .4[ 500 + .8Yt  - 8rt - Yt] or ΔYt = -.08Yt - 3.2rt + 200 and the change in interest rate (Δrt) can be expressed in terms of the excess demand in the money market as Δrt = .5(.25Yt - .9rt - 500) or Δrt = .125Yt - .45rt - 250.   The system can be written in matrix form as (({{ΔYt}, {Δrt}}))= (({{- .08, -3.2}, {.125, - .45}}))(({{Yt}, {rt}}))+(({{200}, {-250}})).  By setting (({{ΔYt}, {Δrt}}))= 0, the long run equilirium can be derived as,
Clear[y, r, yucca, tplot8, tplot9]
Solve[{-.08y -3.2 r +200==0, .125y -.45r - 250==0}, {y, r}]
{{y->2041.28,r->11.4679}}
the steady-state equilibrium income and interest rate are 2,041.28 and 11.47, respectively.  To study the dynamic adjustment of income and interest rate from the initial levels to the equilibrium levels, solve the system of difference equations (({{ΔYt}, {Δrt}}))= (({{- .08, -3.2}, {.125, - .45}}))(({{Yt}, {rt}}))+ (({{200}, {-250}})).  To solve the system, find A + I2 matrix, and solve for (({{Yt + 1}, {rt + 1}}))in terms of (({{Yt}, {rt}})).  That is, (({{Yt + 1}, {rt + 1}})) = ( A + I2)(({{Yt}, {rt}}))+ B.   Using Mathematica, 
A =(({{- .08, -3.2}, {.125, - .45}}));   B = (({{200}, {-250}}));
H[{y_  ,  
                r_}] = (A + IdentityMatrix[2]) . {y,  r} + B
{{200-3.2 r+0.92 y},{-250+0.55 r+0.125 y}}
use the resulting equtions to write the "RSolve" command, 
yucca=RSolve[{y[t+1]==200 -3.2r[t] + .92y[t],r[t+1]==-250 +.55r[t]+.125y[t], y[0]==2100, r[0]==10}, {y[t], r[t]}, t]
\({{y[
              t] -> \(-2100.`\)\ \((\(-1.`\)\ If[t == 0, 1, 0] + 
                    If[t ≥ 
                        1, \((\(-0.08821100917431163`\) - 
                              1.076112880629833`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) - 0.604793353138078`\ i)\)\^\(\(-1\) + t\) - \((\(\(0.08821100917431163`\)\(\)\) - 
                              1.076112880629833`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) + 0.604793353138078`\ i)\)\^\(\(-1\) + t\) - 
                        2.293577981651377`\ 1.`\^\(\(-1\) + t\), \(\(0.`\)\(\)\) + 0.`\ i])\) + 
              3087.0000000000005`\ \((\(-1.`\)\ If[t == 1, 1, 0] + 
                    If[t ≥ 
                        2, \((\(-0.08821100917431163`\) - 
                              1.076112880629833`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) - 0.604793353138078`\ i)\)\^\(\(-2\) + t\) - \((\(\(0.08821100917431163`\)\(\)\) - 
                              1.076112880629833`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) + 0.604793353138078`\ i)\)\^\(\(-2\) + t\) - 
                        2.293577981651377`\ 1.`\^\(\(-2\) + t\), \(\(0.`\)\(\)\) + 0.`\ i])\) - 
              1877.0000000000002`\ \((\(-1.`\)\ If[t == 2, 1, 0] + 
                    If[t ≥ 
                        3, \((\(-0.08821100917431163`\) - 
                              1.076112880629833`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) - 0.604793353138078`\ i)\)\^\(\(-3\) + t\) - \((\(\(0.08821100917431163`\)\(\)\) - 
                              1.076112880629833`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) + 0.604793353138078`\ i)\)\^\(\(-3\) + t\) - 
                        2.293577981651377`\ 1.`\^\(\(-3\) + t\), \(\(0.`\)\(\)\) + 0.`\ i])\), 
          r[t] -> 
            1.1037527593818983`\ \((\(-10.`\)\ \((\((\(\(0.5859908256880737`\)\(\)\) - 
                              0.8442923992851531`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) - 0.604793353138078`\ i)\)\^t + \((\(\(0.5859908256880737`\)\(\)\) + 
                              0.8442923992851531`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) + 0.604793353138078`\ i)\)\^t - 2.0779816513761475`\ 1.`\^t)\) + 
                  6.699999999999989`\ If[
                      t ≥ 
                        1, \((\(\(0.5859908256880737`\)\(\)\) - 0.8442923992851531`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) - 0.604793353138078`\ i)\)\^\(\(-1\) + t\) + \((\(\(0.5859908256880737`\)\(\)\) + 
                              0.8442923992851531`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) + 0.604793353138078`\ i)\)\^\(\(-1\) + t\) - 
                        2.0779816513761475`\ 1.`\^\(\(-1\) + t\), \(\(0.`\)\(\)\) + 0.`\ i] - 
                  1.700000000000017`\ If[
                      t ≥ 
                        2, \((\(\(0.5859908256880737`\)\(\)\) - 0.8442923992851531`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) - 0.604793353138078`\ i)\)\^\(\(-2\) + t\) + \((\(\(0.5859908256880737`\)\(\)\) + 
                              0.8442923992851531`\ i)\)\ \((\(\(0.7350000000000001`\)\(\)\) + 0.604793353138078`\ i)\)\^\(\(-2\) + t\) - 
                        2.0779816513761475`\ 1.`\^\(\(-2\) + t\), \(\(0.`\)\(\)\) + 0.`\ i])\)}}\)
to graph the time path of Y,
y[t_]=y[t]/.yucca[[1]];
tplot8=Table[y[t], {t, 0, 100}];
ListPlot[tplot8, PlotJoined->True] 

⁃Graphics⁃
it is obvious that as t->∞, Y approaches its steady-state equilibrium value 2041.  To graph the time path of r,
r[t_]=r[t]/.yucca[[1]];
tplot9=Table[r[t], {t, 0,100}];
ListPlot[tplot9, PlotJoined->True] 

⁃Graphics⁃
its is obvious that as t->∞, r approaches its steady-state equilibrium value 11.47.  To graph the first 100 points visited in the dynamic adjustment from one point to the equilibrium point use the following commands,
H[{y_, r_}]:={200-3.2r+.92y, -250 +.55r +.125y}
f[y_, r_, t_]:=NestList[H, {y, r}, 100]
vanilla = f[2500, 8, 100];
pl1 =ListPlot[vanilla, PlotStyle->PointSize[.02], PlotRange->All]

⁃Graphics⁃
It is also possible to graph the IS and LM curves in y-r space and study the dynamics of the adjustment over time from one point to another or to the equilibrium point.  To do this, find r in terms of y for IS and LM curves and plot them.  Then, show the lisplot above and the IS-LM graphs together as shown below,
IS=Solve[500+.8y -8r==y, r]
LM=Solve[.25y-.9r==500, r]
{{r->0.125 (500.-0.2 y)}}
{{r->1.11111 (-500.+0.25 y)}}
pl2=Plot[{.125(500-.2y), 1.111(-500+.25y)}, {y, 2000, 2200}, PlotRange->All]

⁃Graphics⁃
Show[pl1, pl2] 


⁃Graphics⁃
it is obvious from this graph that the equilibrium point is asymptotically stable and as time passes the points converge to the steady-state equilibrium point.  The stability of the system of difference equations is discussed in the following section.

Do More Practice:  
1- Change the speed of the adjustments in the goods and money markets from .4 and .5 to .2 and .9, respectively.   Re-run the program and observe changes in the solution and the adjustment path.
2- Change the speed of the adjustments in the goods and money markets from .4 and .5 to .9 and .2, respectively.   Re-run the program and observe changes in the solution and the adjustment path.
3- Suppose in the IS-LM model of the example investment is independent of the interest rate.  Re-run the program under this assumption and compare the results.
4- Change the model as C = .7Y + 400, I = 100 - r, Y = C + I, Md = .3Y - 1.2r, and Md = Ms = 480.  Keep the speed of adjustments in the goods and money markets the same as in the example.   Re-run the program and observe changes in the solution and the adjustment path.

Stability Analysis of the System of Linear Difference Equations:   The eigenvalues of the A + I matrix of the system of linear difference equations reveal much about the time path of the variables of the system and overall stability of the equilibrium point.  The equilibrium point of a system of linear difference equations will asymptotically be stable if and only if the absolute value of the real parts of the eigenvalues of the A + I matrix are less than one.  

Example:  Analyze the stability of the equilibrium point of the IS-LM model of the preceding example.
Solution:  The A matrix of the system of the difference equations of the preceding example is A =(({{- .08, -3.2}, {.125, - .45}})).   Find the eigenvalues of the A + I2 matrix,
Clear[A]
A =(({{- .08, -3.2}, {.125, - .45}}));
Eigenvalues[A   +  IdentityMatrix[2]]
{0.735+0.604793 i,0.735-0.604793 i}
the real parts of the eigenvalues of the A + I2 matrix is .735 and the system is asymptotically stable.

Example:   Analyze the stability of the system of difference equations ΔY1t = .3Y1t - .1Y2t - 10, ΔY2t = -Y1t +1.2Y2t + 12, Y1(0) = 50, Y2(0) = 25.  Solve the system and show the time path of each variables and the relationship between y1 and y2 over time.
Solution:   The system can be written in the matrix form as (({{ΔY1t}, {ΔY2t}}))= (({{.3, - .1}, {-1, 1.2}}))(({{Y1t}, {Y2t}}))+ (({{-10}, {12}})).  To analyze the stability of the system, find the eigenvalues of the A + I2 matrix, 
Clear[A,  B]
A =(({{.3, - .1}, {-1, 1.2}}));   B = (({{-10}, {12}}));
Eigenvalues[A + IdentityMatrix[2]]
{2.3,1.2}
the eigenvalues of the A + I2 matrix are 2.3 and 1.2.  Since the eigenvalues are greater than one, the equilibrium point is not stable.   To solve the system, find Y1t+1 and Y2t+1 in terms of Y1t and Y2t and solve the system of difference equations.  Using Mathematica,  
Clear[y1, y2, tendrils,sessile,   pl3, pl4, pl5, pl6]
H[{y1_  , y2_}]=(A+IdentityMatrix[2]).{y1, y2}+B
{{-10+1.3 y1-0.1 y2},{12-y1+2.2 y2}}
use the resulting equtions to write the "RSolve" command, 
tendrils =RSolve[{y1[t+1]==-10 +1.3y1[t]-.1y2[t],y2[t+1]==12-y1[t]+2.2y2[t],y1[0]==50,y2[0]==25}, {y1[t], y2[t]}, t]
\({{y1[
              t] -> \(-50.`\)\ \((\(-0.9999999999999998`\)\ If[t == 0, 
                        1, 0] + 
                    If[t ≥ 
                        1, \(-3.8461538461539724`\)\ 1.0000000000000024`\^t + 
                        6.545454545454665`\ 1.1999999999999968`\^t - 
                        3.699300699300693`\ 2.3000000000000007`\^t, 0])\) + 
              172.5`\ \((\(-0.9999999999999998`\)\ If[t == 1, 1, 0] + 
                    If[t ≥ 
                        2, \(-3.8461538461539626`\)\ 1.0000000000000024`\^t + 
                        5.454545454545569`\ 1.1999999999999968`\^t - 
                        1.608391608391605`\ 2.3000000000000007`\^t, 0])\) - 
              133.30000000000004`\ \((\(-0.9999999999999998`\)\ If[t == 2, 1, 
                        0] + If[
                      t ≥ 
                        3, \(-3.8461538461539537`\)\ 1.0000000000000024`\^t + 
                        4.545454545454652`\ 1.1999999999999968`\^t - 
                        0.6993006993006976`\ 2.3000000000000007`\^t, 0])\), 
          y2[t] -> 
            0.3623188405797101`\ \((\(-748.3846153846375`\)\ 1.0000000000000024`\^t + 986.0727272727493`\ 1.1999999999999968`\^t - 
                  168.68811188811148`\ 2.3000000000000007`\^t - 
                  76.9`\ If[
                      t ≥ 
                        2, \(-10.615384615384915`\)\ 1.0000000000000024`\^t + 
                        12.545454545454843`\ 1.1999999999999968`\^t - 
                        1.9300699300699256`\ 2.3000000000000007`\^t, 
                      0])\)}}\)
to graph the time path of Y1,
y1[t_]=y1[t]/.tendrils[[1]];
pl3=Table[y1[t], {t, 0, 10}];
ListPlot[pl3, PlotJoined->True] 

⁃Graphics⁃
to graph the time path of Y2,
y2[t_]=y2[t]/.tendrils[[1]];
pl4=Table[y2[t], {t, 0,10}];
ListPlot[pl4, PlotJoined->True] 

⁃Graphics⁃
To graph Y1t against Y2t over time use the following commands,
H[{y1_, y2_}]:={-10+1.3y1-.1y2, 12 -y1 +2.2y2}
f[y1_, y2_, t_]:=NestList[H, {y1, y2}, 5]
sessile = f[40, 20,5];
pl5 =ListPlot[sessile, PlotStyle->PointSize[.02], PlotRange->All]

⁃Graphics⁃
It is also possible to graph the equations of the system for their steady-state equilibrium by settimg (({{ΔY1t}, {ΔY2t}}))= (({{0}, {0}})) and study the adjustment over time from any point on the Y1-Y2 space towards the equilibrium point.
Clear[IS, LM]
IS=Solve[-10+.3y1 -.1y2==0, y2]
LM=Solve[12 -y1+1.2y2==0, y2]
{{y2->10. (-10.+0.3 y1)}}
{{y2->-0.833333 (12.-1. y1)}}
pl6=Plot[{10(-10 + .3y1), -.83(12-y1)}, {y1, 10, 60}, PlotRange->All]
Show[pl5, pl6] 


⁃Graphics⁃

⁃Graphics⁃
the graphs of Y1, Y2, and the dynamic relations between Y1 and Y2 indicate that the variables and the model are all divergent.

Problems:
1.  For the following system of difference equations 

a-  (({{ΔY1t}, {ΔY2t}})) = (({{- .3, 0}, {1.2, .8}}))(({{Y1t}, {Y2t}})) + (({{9}, {-40}})),  (({{Y1 (0)}, {Y2 (0)}}))= (({{5}, {30}}))

b-  (({{ΔY1t}, {ΔY2t}})) = (({{- .3, -1}, {1.5, .2}}))(({{Y1t}, {Y2t}})) + (({{40}, {-30}})),  (({{Y1 (0)}, {Y2 (0)}}))= (({{12}, {35}}))

c-  (({{ΔY1t}, {ΔY2t}})) = (({{-1.2, .8}, {.5, .2}}))(({{Y1t}, {Y2t}})) + (({{-5}, {-10}})),  (({{Y1 (0)}, {Y2 (0)}}))= (({{10}, {20}}))

d-  (({{ΔY1t}, {ΔY2t}})) = (({{1, -1}, {.5, 2}}))(({{Y1t}, {Y2t}})) + (({{-24}, {-20}})),  (({{Y1 (0)}, {Y2 (0)}}))= (({{16}, {5}}))
a.  find the steady-state equilibrium point.
b.  solve the system of difference equations and graph the time path of Y1 and Y2.
c.  graph Y1 and Y2 on the Y1-Y2 space and 50 points around the steady-state equilibrium.
d.  analyze the stability of the system.
2.  A simple macro model of an economy is represented by the following equations: C = 800 + .9Y - 2r,  I = 200 - 4r, Y = C + I,  Md = .4Y - r, and Md = Ms = 3800.  The first equation expresses consumption (C) as a function of income (y) and interest rate (r); the second equation expresses investment ( I ) as a function of interest rate; the third equation is the equilibrium condition in the goods market; the fourth equation expresses demand for money (Md) as a function of income and interest rate; and the last equation is the equilibrium condition in the money market.  Suppose that the change in income in the goods market and the change in interest rate in the money market are proportional to excess demand in the corresponding markets, with a proportionality factor of .5 for the goods market and a proportionality factor of .8 for the money markets.  Given that the initial values of y and r are Y(0) = 10000 and r(0) = 8, find the time path of Y and r and analyze the stability of the model and the dynamic relationship between Y and r over time.