Markov Chains

by Teresa Anderson

Final Project
Using Mathematica in the Mathematics Classroom 
June, 1998

This unit is designed for use in any Algebra class where students are familiar with matrices, matrix multiplication, and properties of matrices.

Initialization

Before proceeding, activate the following code:

<< "Graphics`Colors`"

<< "Graphics`Arrow`"

$TextStyle = {FontFamily -> "Geneva", FontSize -> 12}

PD  =  
      Graphics[{RGBColor[1, 0, 0], Text["D", {0, 0}, {0, 0}]}, 
        AspectRatio -> Automatic];

\(\(PDE\  = \ \ 
      Graphics[{RGBColor[1, 0, 0], Text[".3", {1.3, 1}, {\(-1\), 0}]

          	}, AspectRatio -> Automatic]; \)\)

PDW  =   
      Graphics[{RGBColor[1, 0, 0], Text[".2", {-1.3, 1}, {1, 0}]}, 
        AspectRatio -> Automatic];

PDD  =   
      Graphics[{RGBColor[1, 0, 0], 
          Text[".5", {0, -0.8}, {0, -1}]}, 
        AspectRatio -> Automatic];

DD  =  
      Graphics[{RGBColor[1, 0, 0], 
          Circle[{0, -0.25}, 0.25, {3  π/4, 9  π/4}]}, 
        AspectRatio -> Automatic];

DE  =  
      Graphics[{RGBColor[1, 0, 0], Arrow[{0.1, -0.1}, {2.1, 1.85}]}, 
        AspectRatio -> Automatic];

DW  =  
      Graphics[{RGBColor[1, 0, 0], 
          Arrow[{-0.1, -0.1}, {-2.1, 1.85}]}, 
        AspectRatio -> Automatic];

PE  =   
      Graphics[{RGBColor[0, 0, 1], Text["E", {2, 2.1}, {0, 0}]}, 
        AspectRatio -> Automatic];

PW  =  
      Graphics[{RGBColor[1, 0, 1], Text["W", {-2, 2.1}, {0, 0}]}, 
        AspectRatio -> Automatic];

Downtown  = {PD, PDE, PDW, PDD, DD, DE, DW, PE, PW};

\(\(PED\  = \ \ 
      Graphics[{RGBColor[0, 0, 1], Text[".6", { .8, 1}, {0, 0}]
	}, 
        AspectRatio -> Automatic]; \)\)

PEW  =   
      Graphics[{RGBColor[0, 0, 1], Text[".2", {0, 2.15}, {0, 0}]}, 
        AspectRatio -> Automatic];

PEE  =   
      Graphics[{RGBColor[0, 0, 1], Text[".2", {2.5, 2.1}, {0, 0}]}, 
        AspectRatio -> Automatic];

EE  =  
      Graphics[{RGBColor[0, 0, 1], 
          Circle[{2.1, 2.1}, 0.25, {5  π/4, 11  π/4}]}, 
        AspectRatio -> Automatic];

ED  =  
      Graphics[{RGBColor[0, 0, 1], Arrow[{1.9, 1.9}, {0.1, 0.1}]}, 
        AspectRatio -> Automatic];

EW  =  
      Graphics[{RGBColor[0, 0, 1], Arrow[{1.8, 2}, {-1.8, 2}]}, 
        AspectRatio -> Automatic];

Eastside  = {PD, PED, PEW, PEE, ED, EE, EW, PE, PW};

\(\(PWD\  = \ \ 
      Graphics[{RGBColor[1, 0, 1], Text[".3", {\(- .8\), 1}, {0, 0}]

          	}, AspectRatio -> Automatic]; \)\)

PWW  =   
      Graphics[{RGBColor[1, 0, 1], Text[".3", {-2.5, 2.1}, {0, 0}]}, 
        AspectRatio -> Automatic];

PWE  =   
      Graphics[{RGBColor[1, 0, 1], Text[".4", {0, 1.75}, {0, 0}]}, 
        AspectRatio -> Automatic];

WW  =  
      Graphics[{RGBColor[1, 0, 1], 
          Circle[{-2.1, 2.1}, 0.25, {π/4, 7  π/4}]}, 
        AspectRatio -> Automatic];

WD  =  
      Graphics[{RGBColor[1, 0, 1], Arrow[{-1.9, 1.9}, {-0.1, 0.1}]}, 
        AspectRatio -> Automatic];

WE  =  
      Graphics[{RGBColor[1, 0, 1], Arrow[{-1.9, 1.9}, {1.9, 1.9}]}, 
        AspectRatio -> Automatic];

Westside  = {PD, PWD, PWW, PWE, WD, WE, WW, PE, PW};

T  =  {{ .5,  .3,  .2}, { .6,  .2,  .2}, { .3,  .4,  .3}};

Taxis  =  {{100, 100, 100}};

WA = 
      ListPlot[{{1, 3}, {0, 3}, {0, 7}, {3, 7}, {3, 6}}, 
        PlotStyle -> {PointSize[0.01]}, PlotJoined -> True, 
        Axes -> False, DisplayFunction -> Identity];

WB = 
      ListPlot[{{2, 3}, {3, 3}, {3, 5}}, PlotStyle -> {PointSize[0.01]}, 
        PlotJoined -> True, Axes -> False, 
        DisplayFunction -> Identity];

WC = 
      ListPlot[{{3, 4}, {4, 4}}, PlotStyle -> {PointSize[0.01]}, 
        PlotJoined -> True, Axes -> False, 
        DisplayFunction -> Identity];

WF = 
      ListPlot[{{3, 7}, {6, 7}, {6, 4}, {5, 4}}, 
        PlotStyle -> {PointSize[0.01]}, PlotJoined -> True, 
        Axes -> False, DisplayFunction -> Identity];

WG = 
      ListPlot[{{6, 7}, {9, 7}, {9, 3}, {8, 3}}, 
        PlotStyle -> {PointSize[0.01]}, PlotJoined -> True, 
        Axes -> False, DisplayFunction -> Identity];

WH = 
      ListPlot[{{7, 3}, {6, 3}}, PlotStyle -> {PointSize[0.01]}, 
        PlotJoined -> True, Axes -> False, 
        DisplayFunction -> Identity];

WI = 
      ListPlot[{{6, 4}, {6, 2}}, PlotStyle -> {PointSize[0.01]}, 
        PlotJoined -> True, Axes -> False, 
        DisplayFunction -> Identity];

WJ = 
      ListPlot[{{6, 1}, {6, 0}}, PlotStyle -> {PointSize[0.01]}, 
        PlotJoined -> True, Axes -> False, 
        DisplayFunction -> Identity];

WK = 
      ListPlot[{{9, 3}, {9, 0}, {0, 0}, {0, 3}}, 
        PlotStyle -> {PointSize[0.01]}, PlotJoined -> True, 
        Axes -> False, DisplayFunction -> Identity];

RA =  Graphics[Text["1", {1.5, 5}]];

RB =  Graphics[Text["2", {4.5, 5.5}]];

RC =  Graphics[Text["3", {7.5, 5}]];

RD =  Graphics[Text["4", {7.5, 1.5}]];

RE =  Graphics[Text["5", {3.5, 1.5}]];

Maze = {WA, WB, WC, WF, WG, WH, WI, WJ, WK};

Titles = {RA, RB, RC, RD, RE};

Introduction

Markov chains use probability theory and matrix multiplication  to model a variety of situations. Let's examine the following situation to see how Markov chains work.

A taxi company with 300 cabs serves a city from three dispatch centers:  Downtown, Eastside, and Westside.  The taxi drivers' average pickup interval is 20 minutes (the time between picking up one fare and the next) and the demand for taxis is the same at all times of day and from day to day.  Based on records from previous years, the head dispatcher has found that, on average, 50% of the fares picked up in the Downtown region go to another Downtown location, 30% go to an Eastside address, and the remaining 20% go to the Westside.  Enter the command which follows for a diagram of this  information.

Show[Downtown];

Eastside fares follow a different pattern:  20% stay on the Eastside, 60% go Downtown, and 20% go to the Westside.  For a diagram, enter the command below.

Show[Eastside];

30% of Westside fares stay on the Westside, 40% go to the Eastside, and 30% go Downtown as shown in the diagram produced by entering the command below.

Show[Westside];

By combining the data from all three areas, the dispatcher produces diagram of taxi movement through the city; enter the following command to see it.

Show[{Downtown, Eastside, Westside}];

The data can also be displayed in the table.  Entering the following command will produce the table.

TableForm[T, 
      TableHeadings -> {{"Downtown", "Eastside", "Westside"}, {
            "Downtown", "Eastside", "Westside"}}, 
      TableAlignments -> Center]

The table can condensed into the matrix T, which will be shown below after entering the command.

MatrixForm[T]

Matrix T is called a transition matrix because each entry represents the probability of transition from one region to another.  For example, (Subscript[T, 21]) is the probability of moving from region 2 to region 1; in the context of our taxi model, this is the probability that a taxi picking up a fare on the Eastside will drop off that fare Downtown.  Notice that T does not include the identification of the regions; when working with transition matrices it is useful have a labeled table (like the one above) at hand for reference.

Our taxi model describes a system with a finite number of states - the three regions of the city - and probabilities for transitions from each state to every other state.  These probabilities are constant and independent of previous behavior.  A transition - picking up a new fare - occurs each time the system is observed,  and observations occur at regular intervals, every 20 minutes in this model.  Systems with these characteristics are called Markov chains or Markov processes.

Our transition matrix demonstrates some properties which are common to all transition matrices for Markov processes:

1.  A transition matrix must be square.  (The number of rows and number of columns are both the number of states in the system.)

2.  All the entries are between 0 and 1 inclusive.  (The entries are probabilities which must lie in that interval.)

3.  The sum of the entries in any row must be 1.  (Each row represents probabilities of moving from a particular state to any of the states in the system.  Since the system is closed and the movement must take place, the sum of the probabilities must be 1.)

4.  The entries in the matrix are constant.  (The Markov chain model depends upon the assumption that the transition matrix does not change throughout the process.)

Multiple Transitions

When the morning drivers go off duty at 2:00 PM, 300 afternoon drivers are waiting at the three regional dispatch centers to pick up a cab.  Because each taxi has picked up 24 fares during the eight hour morning shift, there may no longer be 100 taxis in each region; therefore the dispatcher must determine how the cabs are distributed at the shift change to optimize the assignment of afternoon shift drivers between the three dispatch centers.

The information presented in the transition matrix, T, allows the dispatcher to approximate the distribution of taxis throughout the city.  Assume the company starts each morning at 6:00 AM with 100 cabs in each of the three regions of the city.  This can be represented by the matrix Taxis which will be shown below after entering the command.

MatrixForm[Taxis]

From the dispatcher's records, it is known that on average 50% of the 100 fares which begin Downtown will remain there, 60% of the 100 Eastside taxis will pickup fares going Downtown, and 30% of the 100 Westside fares are going Downtown.  The number of taxis Downtown after the first twenty minutes (or after one fare) can be projected by the computation below; click on the computation and Enter to complete the calculation.

.5*100 +  .6*100 +  .3*100

This computation is the result of multiplying the first (and only) row in the matrix Taxis by the first column in the transition matrix (which contains the probabilities of a fare going Downtown from each of the three regions).  The estimated number of taxis on the Eastside after one fare can be computed in a similar manner using the second column of the transition matrix.  Thus, the approximate distribution of cabs across the city after twenty minutes is the product of the matrix Taxis and the transition matrix.  The command below will show the results of that computation.

MatrixForm[Taxis . T]

The exact distribution of taxis will vary from day to day; this matrix product shows the long term average distribution to be 140 taxis downtown, 90 on the Eastside, and 70 on the Westside at 6:20 AM, after one fare.  

Multiplying the 6:20 AM distribution matrix by the transition matrix, T, will give the approximate distribution of taxis after 40 minutes.  Enter the computation below to see the results.

MatrixForm[{{140, 90, 70}} . T]

At 6:40 AM, after two fares, there are an estimated 145 taxis Downtown, 88 on the Eastside, and 67 cabs on the Westside.  Continuing this multiplication process through all 24 morning fares will produce the projected distribution of cabs at 2:00 PM.  If we examine this process more closely, a less cumbersome method of computation appears.

Using only 6 fares to streamline the notation, the necessary matrix multiplication is:

         {{{{{{Taxis . T} . T} . T} . T} . T} . T}
         
Application of the associative law of matrix multiplication allows this to be rewritten as:         

	{Taxis . {T . {T . {T . {T . {T . T}}}}}} = {Taxis . (T^6)}

Extending this concept to all 24 fares throughout the morning, the projected distribution of taxis throughout the city can be found by computing {Taxis . (T^24)} which will be found by entering the command below.

MatrixForm[Taxis . MatrixPower[T, 24]]

Although the dispatcher knows that there will not be exactly 145.455 taxis near the Downtown center, the resulting distribution matrix does provide a basis upon which to make driver assignments.  One scheme the dispatcher might use is to pre-assign 145 drivers Downtown, 87 to the Eastside, 66 to the Westside, and leave 3 drivers as unassigned "floaters" who would not be given an assignment until the actual distribution of taxis at 2:00 PM is known.

Once the morning shift is over, the taxi company is responsible for delivering the drivers to the dispatch center where they began the shift.  To arrange adequate transportation for the drivers, the dispatcher needs to know where the drivers at each center began the shift.   Looking only at the transition matrix portion of the computation used to determine cab distribution at 2:00 PM as written below, can provide the necessary information.  Enter the command to see the results.

MatrixForm[MatrixPower[T, 24]]

Approximately 48% of the drivers who started the shift in any region will end up Downtown, 29% will end the shift on the Eastside, and 22% will be on the Westside at 2:00 PM.  Therefore, the dispatcher can estimate that transportation will be needed Downtown to take approximately 48 drivers (48% of 100) to the Eastside and 48 to the Westside, approximately 29 drivers will need to go from the Eastside to Downtown and 29 to the Westside, and the Westside dispatch center will deliver approximately 22 drivers Downtown and 22 to the Eastside.

Long Term Behavior

In computing the distribution of taxis at various times during the morning shift, we saw a significant change of distribution from the 100 taxis per region at 6:00 AM to that after one fare.  Enter the command below to review those results.

MatrixForm[Taxis . T]

The next command shows what happens between 6:20 AM and 6:40 AM; notice the distribution doesn't change as dramatically as in the first twenty minutes.

MatrixForm[Taxis . MatrixPower[T, 2]]

The command below shows the difference in distribution between 6:40 AM and 2:00 PM to be very small.

MatrixForm[Taxis . MatrixPower[T, 24]]

Let's examine what happens to taxi distribution across the city between 6:00 AM and 2:00 PM in more detail.  The command below will compute the approximate distribution of taxis every twenty minutes from 6 AM to 2 PM (0600 to 1400 hours in military time).

\(\(For[i = 0, i < 25, 
      Print[600 + 100*IntegerPart[i/3] + 60*FractionalPart[i/3], 
        MatrixForm[Taxis . MatrixPower[T, i]]]; \(i++\)]
\)\)

As you can see, after 8:00 AM the distribution of taxis becomes constant and doesn't change with subsequent transitions.  When this happens in a system represented by a Markov process, the system is said to stabilize and the resulting distribution matrix (or vector, since it has only one row) is called the stable state vector of the Markov process.

Exercises

A rat is placed in the maze shown below.  During  a fixed time interval, the rat randomly chooses one of the doors available to it (depending upon which room it is in) and moves through that door to the next room - it always moves to a new room and never stays where it is.  Each movement of the rat is taken as a transition in a Markov process in which a state is identified with the room the rat is in.

Show[Maze, Titles, DisplayFunction -> $DisplayFunction];

The first row of the transition matrix is

TableForm[{{0, 1/2, 0, 0, 1/2}}, TableHeadings -> Automatic, 
      TableAlignments -> Center]

1.	Construct the entire transition matrix for this process.

2.	If the rat starts in room 1, what is the probability that it is in room 3 after two transitions?  After three transitions?

3.	Determine the stable state vector.

4.	After a large number of transitions, what is the probability that the rat is in room 4?

5.	In the long run, what percentage of the time will the rat spend in rooms 2 or 3?

Answers

1.	The transition matrix is:

MatrixForm[{{0, 1/2, 0, 0, 1/2}, {1/2, 0, 0, 0, 1/2}, {0, 0, 0, 1, 0}, {
          0, 0, 1/2, 0, 1/2}, {1/3, 1/3, 0, 1/3, 0}}]

2.  	0; 0.083

MatrixForm[
      N[{{1, 0, 0, 0, 0}} . 
          MatrixPower[{{0, 1/2, 0, 0, 1/2}, {1/2, 0, 0, 0, 1/2}, {0, 0, 0, 1, 
                0}, {0, 0, 1/2, 0, 1/2}, {1/3, 1/3, 0, 1/3, 0}}, 2], 3]]

MatrixForm[
      N[{{1, 0, 0, 0, 0}} . 
          MatrixPower[{{0, 1/2, 0, 0, 1/2}, {1/2, 0, 0, 0, 1/2}, {0, 0, 0, 1, 
                0}, {0, 0, 1/2, 0, 1/2}, {1/3, 1/3, 0, 1/3, 0}}, 3], 3]]

3.	The stable state vector is:

MatrixForm[
      N[{{1, 0, 0, 0, 0}} . 
          MatrixPower[{{0, 1/2, 0, 0, 1/2}, {1/2, 0, 0, 0, 1/2}, {0, 0, 0, 1, 
                0}, {0, 0, 1/2, 0, 1/2}, {1/3, 1/3, 0, 1/3, 0}}, 50], 3]]

4.	0.2

5.	The stable state vector shows the long run probability of being in room 2 to be 20% and in room 3 to be 10%.  So in the long run, the rat will spend 20% + 10%, or 30%, of the time in rooms 2 or 3.