Example 1

parts 1 and 2

P <- matrix(c(0, 1/3, 0, 1/3, 1/3, 0,
             1/2, 0, 0, 1/2, 0, 0,
             0, 0, 0, 1/2, 1/2, 0,
             1/4, 1/4, 1/4, 0, 1/4, 0,
             1/4, 0, 1/4, 1/4, 0, 1/4,
             0, 0, 0, 0, 1, 0), 
           nrow = 6, ncol = 6, byrow = TRUE)
matrixpower(mat = P, k = 5)
##           [,1]      [,2]       [,3]      [,4]      [,5]       [,6]
## [1,] 0.1588542 0.1425058 0.10590278 0.2484086 0.3022280 0.04210069
## [2,] 0.2137587 0.1085069 0.12999132 0.2697483 0.2035590 0.07443576
## [3,] 0.1588542 0.1299913 0.10026042 0.2615017 0.3044705 0.04492188
## [4,] 0.1863064 0.1348741 0.13075087 0.2343750 0.2540148 0.05967882
## [5,] 0.2266710 0.1017795 0.15223524 0.2540148 0.1727431 0.09255642
## [6,] 0.1263021 0.1488715 0.08984375 0.2387153 0.3702257 0.02604167
matrixpower(mat = P, k = 10)
##           [,1]      [,2]      [,3]      [,4]      [,5]       [,6]
## [1,] 0.1926234 0.1223991 0.1282373 0.2506364 0.2401568 0.06594699
## [2,] 0.1835986 0.1273149 0.1227218 0.2490602 0.2575114 0.05979302
## [3,] 0.1923560 0.1227218 0.1283517 0.2500969 0.2406491 0.06582441
## [4,] 0.1879773 0.1245301 0.1250485 0.2505551 0.2490803 0.06280872
## [5,] 0.1801176 0.1287557 0.1203246 0.2490803 0.2642060 0.05751584
## [6,] 0.1978410 0.1195860 0.1316488 0.2512349 0.2300634 0.06962590
matrixpower(mat = P, k = 25)
##           [,1]      [,2]      [,3]      [,4]      [,5]       [,6]
## [1,] 0.1874715 0.1250147 0.1249818 0.2499965 0.2500549 0.06248060
## [2,] 0.1875220 0.1249887 0.1250140 0.2500027 0.2499576 0.06251497
## [3,] 0.1874727 0.1250140 0.1249826 0.2499967 0.2500525 0.06248147
## [4,] 0.1874974 0.1250014 0.1249983 0.2499997 0.2500050 0.06249822
## [5,] 0.1875412 0.1249788 0.1250262 0.2500050 0.2499207 0.06252801
## [6,] 0.1874418 0.1250299 0.1249629 0.2499929 0.2501120 0.06246042
matrixpower(mat = P, k = 50)
##        [,1]  [,2]  [,3] [,4] [,5]       [,6]
## [1,] 0.1875 0.125 0.125 0.25 0.25 0.06250000
## [2,] 0.1875 0.125 0.125 0.25 0.25 0.06250000
## [3,] 0.1875 0.125 0.125 0.25 0.25 0.06250000
## [4,] 0.1875 0.125 0.125 0.25 0.25 0.06250000
## [5,] 0.1875 0.125 0.125 0.25 0.25 0.06250000
## [6,] 0.1875 0.125 0.125 0.25 0.25 0.06250001
matrixpower(mat = P, k = 100)
##        [,1]  [,2]  [,3] [,4] [,5]   [,6]
## [1,] 0.1875 0.125 0.125 0.25 0.25 0.0625
## [2,] 0.1875 0.125 0.125 0.25 0.25 0.0625
## [3,] 0.1875 0.125 0.125 0.25 0.25 0.0625
## [4,] 0.1875 0.125 0.125 0.25 0.25 0.0625
## [5,] 0.1875 0.125 0.125 0.25 0.25 0.0625
## [6,] 0.1875 0.125 0.125 0.25 0.25 0.0625

The rows of the matrix powers are settling down toward a single limiting row vector.

Looking at the entries of the \(P^{100}\) matrix, we find that \[P(X_{100} = 4 \mid X_0 = 1) = P^{100}_{14} = 0.25\] and \[P(X_{100} = 4 \mid X_0 = 1) = P^{100}_{34} = 0.25.\] In general, \(P(X_n = j \mid X_0 = i) = P^n_{ij}\), so \[P(X_{100} = 1 \mid X_0 = 5) = 0.1875,\quad P(X_{100} = 6 \mid X_0 = 2) = 0.0625.\]

parts 3, 4, and 5

If \(X_0 = 1\), then the initial distribution of the chain is \[\alpha = (1, 0, 0, 0, 0, 0).\] If the chain starts at a uniformly random state, then \(X_0 \sim \beta\) where \[\beta = (1/6, 1/6, 1/6, 1/6, 1/6, 1/6).\]

When \(X_0 \sim \alpha\), the distribution of \(X_n\) is \(\alpha P^n\). Similarly, when \(X_0 \sim \beta\), the distribution of \(X_n\) is \(\beta P^n\). The following commands compute the distributions of \(X_{10}\) and \(X_{100}\) for initial distributions \(\alpha\) and \(\beta\).

alpha <- c(1, 0, 0, 0, 0, 0)
beta <- c(1/6, 1/6, 1/6, 1/6, 1/6, 1/6)
alpha %*% matrixpower(P, 10)
##           [,1]      [,2]      [,3]      [,4]      [,5]       [,6]
## [1,] 0.1926234 0.1223991 0.1282373 0.2506364 0.2401568 0.06594699
beta %*% matrixpower(P, 10)
##           [,1]      [,2]      [,3]      [,4]      [,5]       [,6]
## [1,] 0.1890856 0.1242179 0.1260555 0.2501106 0.2469445 0.06358582
alpha %*% matrixpower(P, 100)
##        [,1]  [,2]  [,3] [,4] [,5]   [,6]
## [1,] 0.1875 0.125 0.125 0.25 0.25 0.0625
beta %*% matrixpower(P, 100)
##        [,1]  [,2]  [,3] [,4] [,5]   [,6]
## [1,] 0.1875 0.125 0.125 0.25 0.25 0.0625

part 6

The interpretation of what we’ve seen is as follows: for this random walk, it seems that there is a limiting distribution for \(X_n\) as \(n \to \infty\), which is given by one of the rows of \(P^n\) for large \(n\). Moreover, the chain reaches this limiting distribution regardless of the initial distribution.

Example 2

Q = matrix(c(0, 1/2, 0, 0, 0, 1/2,
             1/2, 0, 1/2, 0, 0, 0,
             0, 1/2, 0, 1/2, 0, 0,
             0, 0, 1/2, 0, 1/2, 0,
             0, 0, 0, 1/2, 0, 1/2,
             1/2, 0, 0, 0, 1/2, 0),
           nrow = 6, ncol = 6, byrow = TRUE)

matrixpower(Q, 8)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.3359375 0.0000000 0.3320312 0.0000000 0.3320312 0.0000000
## [2,] 0.0000000 0.3359375 0.0000000 0.3320312 0.0000000 0.3320312
## [3,] 0.3320312 0.0000000 0.3359375 0.0000000 0.3320312 0.0000000
## [4,] 0.0000000 0.3320312 0.0000000 0.3359375 0.0000000 0.3320312
## [5,] 0.3320312 0.0000000 0.3320312 0.0000000 0.3359375 0.0000000
## [6,] 0.0000000 0.3320312 0.0000000 0.3320312 0.0000000 0.3359375
matrixpower(Q, 9)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.0000000 0.3339844 0.0000000 0.3320312 0.0000000 0.3339844
## [2,] 0.3339844 0.0000000 0.3339844 0.0000000 0.3320312 0.0000000
## [3,] 0.0000000 0.3339844 0.0000000 0.3339844 0.0000000 0.3320312
## [4,] 0.3320312 0.0000000 0.3339844 0.0000000 0.3339844 0.0000000
## [5,] 0.0000000 0.3320312 0.0000000 0.3339844 0.0000000 0.3339844
## [6,] 0.3339844 0.0000000 0.3320312 0.0000000 0.3339844 0.0000000
matrixpower(Q, 10)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.3339844 0.0000000 0.3330078 0.0000000 0.3330078 0.0000000
## [2,] 0.0000000 0.3339844 0.0000000 0.3330078 0.0000000 0.3330078
## [3,] 0.3330078 0.0000000 0.3339844 0.0000000 0.3330078 0.0000000
## [4,] 0.0000000 0.3330078 0.0000000 0.3339844 0.0000000 0.3330078
## [5,] 0.3330078 0.0000000 0.3330078 0.0000000 0.3339844 0.0000000
## [6,] 0.0000000 0.3330078 0.0000000 0.3330078 0.0000000 0.3339844
matrixpower(Q, 11)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.0000000 0.3334961 0.0000000 0.3330078 0.0000000 0.3334961
## [2,] 0.3334961 0.0000000 0.3334961 0.0000000 0.3330078 0.0000000
## [3,] 0.0000000 0.3334961 0.0000000 0.3334961 0.0000000 0.3330078
## [4,] 0.3330078 0.0000000 0.3334961 0.0000000 0.3334961 0.0000000
## [5,] 0.0000000 0.3330078 0.0000000 0.3334961 0.0000000 0.3334961
## [6,] 0.3334961 0.0000000 0.3330078 0.0000000 0.3334961 0.0000000
matrixpower(Q, 12)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.3334961 0.0000000 0.3332520 0.0000000 0.3332520 0.0000000
## [2,] 0.0000000 0.3334961 0.0000000 0.3332520 0.0000000 0.3332520
## [3,] 0.3332520 0.0000000 0.3334961 0.0000000 0.3332520 0.0000000
## [4,] 0.0000000 0.3332520 0.0000000 0.3334961 0.0000000 0.3332520
## [5,] 0.3332520 0.0000000 0.3332520 0.0000000 0.3334961 0.0000000
## [6,] 0.0000000 0.3332520 0.0000000 0.3332520 0.0000000 0.3334961
matrixpower(Q, 50)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [2,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [3,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [4,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [5,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [6,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
matrixpower(Q, 51)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [2,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [3,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [4,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [5,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [6,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
matrixpower(Q, 52)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [2,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [3,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [4,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [5,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [6,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
matrixpower(Q, 53)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [2,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [3,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [4,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [5,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [6,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
matrixpower(Q, 54)
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [2,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [3,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [4,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333
## [5,] 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000
## [6,] 0.0000000 0.3333333 0.0000000 0.3333333 0.0000000 0.3333333

Conclusions

Looking back at the conclusions made at the end of Example 1, a natural question is whether every random walk reaches a limiting distribution and whether the limiting distribution is unique and reached regardless of initial distribution. What we see in this example is that this is not the case. There seem to be two different long term distributions, \[\mu_1 = (0, 1/3, 0, 1/3, 0, 1/3),\quad \mu_2 = (1/3, 0, 1/3, 0, 1/3, 0)\] neither of which is really reached in a limit because \(Q^n\) oscillates between two limiting matrices for large \(n\) depending on whether \(n\) is even or odd.