5.1, Part 2: Adjacency Matrices

The Game¶

That sounds abstract – what is a system, and what is long-term behavior? Let’s illustrate with an example. “The Game” refers to the annual football game between Michigan and Ohio State. Let’s suppose that:

• If Michigan wins in year $k$, there’s an 80% chance we win again in year $k+1$, and a 20% chance that Ohio State wins in year $k+1$.

• If Ohio State wins in year $k$, there’s a 70% chance they win again in year $k+1$, and a 30% chance that Michigan wins in year $k+1$.

Question: In the long-run, what percentage of the time does Michigan win? Perhaps a more approachable way of phrasing this is, what is the (unconditional) probability that Michigan wins a particular game?

This system is an example of a Markov chain, in which the probabilities of transitioning between states don’t change over time. Here, a “state” is the outcome of the game in a particular year. The probabilities described above are fixed: the chance Michigan wins in year 100 given they won in year 99 is 80%, and the chance they win in year 55555 given they won in year 55554 is also 80%.

It’s often useful to visualize Markov chains using a state diagram, which is a directed graph in which each node represents a state, and each edge represents a transition between states. Each edge is labeled with the probability of the transition.

Notice that the sum of the probabilities of the outgoing edges from each node (state) is 1. Intuitively, this says we must account for all possible transitions from a given state.

This Markov chain can be described by a $2 \times 2$ adjacency matrix, $A$:

Note that in an adjacency matrix, for our purposes:

• Each column describes the movement from a given node; the sums of the columns in an adjacency matrix are equal to 1.

• Each row describes the movement into a given node; the sums of the rows in an adjacency matrix don’t necessarily add to 1.

That is, $a_{ij}$ is the probability of transitioning from node $j$ to node $i$.

Recall, our goal is to determine the long-run fraction of games that Michigan wins. How does the adjacency matrix $A$ help us with this?

Simulating a Chain¶

Let’s use it to perform a simulation of this system. Let’s seed our simulation in the year $t = 0$, in which Michigan wins the first game, i.e.

where $\vec x_k$ is a vector with two components, the first of which is the probability that Michigan wins in year $k$, and the second of which is the probability that Ohio State wins in year $k$.

To find the distribution of wins in year 1, we can multiply $\vec x_0$ by $A$:

So, in year 1, we’d expect Michigan to win 80% of the time, and Ohio State to win 20% of the time. What about in year 2?

Or year 3?

In general, the distribution of wins in year $t$ is given by

If you look at $\vec x_1$, $\vec x_2$, and $\vec x_3$ above, they don’t seem to be exploding to some really large value. Despite repeatedly multiplying by $A$, the values seem to be stabilizing. Are these vectors converging? And, if so, to what?

Python can help verify our intuition.

It seems that the sequence of $\vec x_k$’s is converging to $\begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$, which would imply that in the long-run, Michigan wins 60% of the time.

An Eigenvalue Problem in Disguise¶

Could we have found the vector $\begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$ without simulating the system? You sure bet. This is the eigenvector of $A$ corresponding to the eigenvalue 1! Intuitively, since I’m searching for a long-run distribution, I’m really looking for a vector $\vec x$ such that

meaning that advancing $\vec x$ one step in time won’t change it. This is precisely the definition of an eigenvector with eigenvalue $\lambda = 1$. Eventually, we will prove that all adjacency matrices have an eigenvalue of 1, and that the corresponding eigenvector is the long-run distribution we’re searching for.

But for now, let’s verify that 1 is indeed an eigenvalue of $A$ for this specific $A$ by computing $A$’s eigenvalues and eigenvectors ourselves.

You might be able to spot the eigenvalues from the trace and determinant facts (they need to add to 1.5 and multiply to 0.5). If not, there’s nothing wrong with explicitly writing out the characteristic polynomial.

This tells us that $A$’s eigenvalues are $\lambda_1 = 1$ and $\lambda_2 = 0.5$.

• For $\lambda_1 = 1$, we can find an eigenvector $\vec v_1 = \begin{bmatrix} a \\ b \end{bmatrix}$ by solving $A \begin{bmatrix} a \\ b \end{bmatrix} = \begin{bmatrix} a \\ b \end{bmatrix}$.

  $$A \begin{bmatrix} a \\ b \end{bmatrix} = \begin{bmatrix} a \\ b \end{bmatrix} \implies \begin{bmatrix} 0.8a + 0.3b \\ 0.2a + 0.7b \end{bmatrix} = \begin{bmatrix} a \\ b \end{bmatrix}$$

  Both components of the equation above give the relationship

  $$0.8a + 0.3b = a \implies 0.2a = 0.3b \implies b = \frac{2}{3}a$$

  So, one eigenvector is $\vec v_1 = \begin{bmatrix} 3 \\ 2 \end{bmatrix}$. But, since we’d like to interpret the components of $\vec v_1$ as a probability distribution, we should scale it so that its components sum to 1. Since $3 + 2 = 5$, this gives $\vec v_1 = \begin{bmatrix} 3/5 \\ 2/5 \end{bmatrix} = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$, which is precisely what our simulated distributions converged to!

• For $\lambda_2 = 0.5$, notice that

  $$A \begin{bmatrix} 1 \\ -1 \end{bmatrix} = \begin{bmatrix} 0.8 & 0.3 \\ 0.2 & 0.7 \end{bmatrix} \begin{bmatrix} 1 \\ -1 \end{bmatrix} = \begin{bmatrix} 0.8 - 0.3 \\ 0.2 - 0.7 \end{bmatrix} = \begin{bmatrix} 0.5 \\ -0.5 \end{bmatrix} = 0.5 \begin{bmatrix} 1 \\ -1 \end{bmatrix}$$

  So, $\vec v_2 = \begin{bmatrix} 1 \\ -1 \end{bmatrix}$ is an eigenvector of $A$ corresponding to the eigenvalue $\lambda_2 = 0.5$. There’s no way to scale $\vec v_2$ so that its components are all positive and sum to 1, but there’s no need to, since we’re only interested in interpreting $\lambda_1 = 1$’s eigenvector as a probability distribution.

I want to revisit our Python simulation from above. Using the initial vector $\vec x_0 = \begin{bmatrix} 1 \\ 0 \end{bmatrix}$, the state vectors $\vec x_k$ converged to $\vec v_1 = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$.

But, there was nothing special about our choice of $\vec x_0$. Had we started with, say, $\vec x_0 = \begin{bmatrix} 0.5 \\ 0.5 \end{bmatrix}$ or even $\vec x_0 = \begin{bmatrix} 0 \\ 1 \end{bmatrix}$, the state vectors would have still converged to the eigenvector $\vec v_1 = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$.

What if we start with a $\vec x_0$ that isn’t a probability distribution, like $\vec x_0 = \begin{bmatrix} 50 \\ -15 \end{bmatrix}$?

The state vector is converging on $\begin{bmatrix} 21 \\ 14 \end{bmatrix}$, which is still on the line described by $\vec v_1 = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$. So, it seems that no matter where we start, the state vectors will converge to the eigenvector $\vec v_1 = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$ or some scalar multiple of it.

The Dominant Eigenvalue¶

Why is it the case that $\vec x_k = A^k \vec x_0$ converges to an eigenvector for $\lambda_1 = 1$, and not to an eigenvector for the other eigenvalue of $\lambda_2 = 0.5$?

The answer is that the dominant eigenvalue, i.e. the one with the largest magnitude, is the one that determines the long-run behavior of the system. Let’s see why this is the case.

First, recall that $A$’s eigenvalues are $\lambda_1 = 1$ and $\lambda_2 = 0.5$, and corresponding eigenvectors are $\vec v_1 = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$ and $\vec v_2 = \begin{bmatrix} 1 \\ -1 \end{bmatrix}$. Note that these two eigenvectors are linearly independent, so any other vector $\vec x \in \mathbb{R}^2$ can be written as a linear combination of $\vec v_1$ and $\vec v_2$. The fact that these eigenvectors are linearly independent is key to the argument below.

So, suppose we start with some vector $\vec x = c_1 \vec v_1 + c_2 \vec v_2$. Then, what we’re really interested in is the behavior of the vectors $A \vec x$, $A^2 \vec x$, $A^3 \vec x$, and so on. What is $A \vec x$?

Above, I used the fact that $A \vec v_1 = \lambda_1 \vec v_1$ and $A \vec v_2 = \lambda_2 \vec v_2$ by the definition of eigenvectors.

Recall from the original Chapter 5.1 that if $\mathbf{\lambda}$ is an eigenvalue of $A$ with eigenvector $\vec v$, then $A^k \vec v = \mathbf{\lambda}^k \vec v$ also. So,

and more generally, if $k$ is a positive integer,

So, $A^k \vec x$ is a linear combination of $\vec v_1$ and $\vec v_2$ also, with weights $c_1 \lambda_1^k$ and $c_2 \lambda_2^k$. Here, $\lambda_1 = 1$ and $\lambda_2 = 0.5$, so $\lambda_1^k = 1$ and $\lambda_2^k = 0.5^k$.

But, as $k$ increases, $0.5^k$ approaches 0, while $1^k = 1$. So,

This explains why the state vectors $\vec x_k$ converge to (some scalar multiple of) $\vec v_1 = \begin{bmatrix} 0.6 \\ 0.4 \end{bmatrix}$, no matter where we start!

Why do Adjacency Matrices Have an Eigenvalue of 1?¶

Let’s prove that all adjacency matrices have an eigenvalue of 1. Recall that an adjacency matrix is defined such that:

• Each column describes the movement from a given node; the sums of the columns in an adjacency matrix are equal to 1.

• Each row describes the movement into a given node; the sums of the rows in an adjacency matrix don’t necessarily add to 1.

Here’s the key fact: for any square matrix $A$, the eigenvalues of $A$ are the same as the eigenvalues of $A^T$. The eigenvectors themselves are usually different, but the eigenvalues are the same. This is because $A$ and $A^T$ both have the same characteristic polynomial. Remember from Chapter 2.9 that $\text{det}(A) = \text{det}(A^T)$. So, $A$’s characteristic polynomial is

but this is the same as $\text{det}((A - \lambda I)^T) = \text{det}(A^T - \lambda I)$, which is the characteristic polynomial of $A^T$. So, $A$ and $A^T$ have the same characteristic polynomial, and thus the same eigenvalues.

Back to the original question: why does an adjacency matrix have an eigenvalue of 1? Remember, the columns of an adjacency matrix, which describe the movement from a given node, sum to 1. If $A$ is an adjacency matrix, then

which shows that $A$ has an eigenvector (the all ones vector) with an eigenvalue of 1.

For instance, if $A = \begin{bmatrix} 0.8 & 0.3 \\ 0.2 & 0.7 \end{bmatrix}$, then $A^T = \begin{bmatrix} 0.8 & 0.2 \\ 0.3 & 0.7 \end{bmatrix}$, and $A^T \begin{bmatrix} 1 \\ 1 \end{bmatrix} = \begin{bmatrix} 1 \\ 1 \end{bmatrix}$, which is an eigenvector of $A^T$ with eigenvalue 1.

But, since $A$ and $A^T$ have the same eigenvalues, $A$ must have an eigenvalue of 1.

The followup question is why it’s guaranteed that the largest eigenvalue of an adjacency matrix is 1, because this seems to be crucial in guaranteeing that the state vectors $\vec x_k$ converge to the eigenvector corresponding to the eigenvalue 1. Proving this fact is beyond the scope of our course, but it follows from the Perron-Frobenius theorem if you’d like to read more.

Example: Non-Adjacency Matrices¶

You might wonder, if $A$’s largest eigenvalue is greater than 1, what happens to the vector $A^k \vec x$ as $k$ increases? (I won’t call $A^k \vec x$ a state vector, since it’s not a probability distribution.)

Let’s consider $A = \begin{bmatrix} 3 & 1 \\ 2 & 4 \end{bmatrix}$, which was used as an example in Chapter 5.1. $A$ has eigenvalues $\lambda_1 = 2$ and $\lambda_2 = 5$, and corresponding eigenvectors $\vec v_1 = \begin{bmatrix} 1 \\ -1 \end{bmatrix}$ and $\vec v_2 = \begin{bmatrix} 1 \\ 2 \end{bmatrix}$. These eigenvectors are also linearly independent, so any $\vec x \in \mathbb{R}^2$ can be written as a linear combination of $\vec v_1$ and $\vec v_2$. Then, following the logic from the previous example,

The issue is that both $5^k$ and $2^k$ grow without bound as $k$ increases. So, $A^k \vec x$ will not converge. But, the direction of $A^k \vec x$ is still predictable. Let me start with with $\vec x = \begin{bmatrix} -13 \\ 100 \end{bmatrix}$ (nothing special about these values) and evaluate $A^k \vec x$ for $k = 1, 2, ..., 15$.

The numbers are indeed getting bigger and bigger, but a relationship is emerging between the components of $A^k \vec x$. Let’s plot the coordinates of the vector $A^k \vec x$ for $k = 1, 2, ..., 15$.

I had to log-scale the axes, since the components of $A^k \vec x$ are growing exponentially. But, you should notice that the vectors $A^k \vec x$ all lie on the same line – the line described by the eigenvector for the dominant eigenvalue, $\vec v_2 = \begin{bmatrix} 1 \\ 2 \end{bmatrix}$! Hover over any point above to see the coordinates of the vector $A^k \vec x$ for a particular $k$, and you’ll see the $y$-coordinate tends to be roughly double the $x$-coordinate.

So, when the dominant eigenvalue (that is, the eigenvalue with the largest magnitude) is greater than 1, the vector $A^k \vec x$ doesn’t converge, but its direction converges to the direction of the eigenvector corresponding to the dominant eigenvalue. Through the lens of our most recent equation for $A^k \vec x$,

Notice that if we divide both sides by $5^k$, we see

And as $k$ increases, the term $c_1 \left(\frac{2}{5}\right)^k$ approaches 0. So, the direction of $\frac{A^k \vec x}{5^k}$ approaches the direction of $\begin{bmatrix} 1 \\ 2 \end{bmatrix}$, but $A^k \vec x$ and $\frac{A^k \vec x}{5^k}$ have the same direction (but different magnitudes).

You can view the process of computing $A^k \vec x$ for larger and larger $k$ until the direction stabilizes as being a technique for approximating the eigenvector corresponding to the dominant eigenvalue. This is called the power method.

Here’s the cliffhanger: our analysis of the behavior of $A^k \vec x$ as $k$ increases has depended on the eigenvectors of $A$ (that is, $\vec v_1$ and $\vec v_2$) being linearly independent. What if the eigenvectors of a particular matrix aren’t linearly independent? Time for Chapter 5.2 to shine.