statistical comparison of two markov chain transition matrices

Question

I'm working on a time-dependent dataset, where basically I have two different populations and we're building Markov chains to describe their behaviors. What I'm trying to do is compare the transition matrices associated to the two populations, and decide whether the transition matrices are really different, or whether they could be due to chance. So my question is, is there a good way to compare two transition matrices in order to determine whether they are statistically different? Part of the problem we're running into is that both matrices have a lot of zeroes in them, but not always in the same places. Any suggestions are welcome.

Answer 1

The answer by user @krkeane is interesting.

Following is a different kind of answer, which considers the "effective behavior" of the Markov chain transition matrices.

For a given $n$-by-$n$ Markov chain transition matrix, the most important things are:

1. Does the Stochastic matrix converge? In other words, are the conditions of the Perron–Frobenius theorem satisfied? The original question mentioned the existence of many zeros, which is natural in practice; so to "encourage", but not guarantee, that a typical stochastic matrix with zeros complies with the Perron–Frobenius theorem, it is typical to replace zeros with very small positive values, and re-normalize row sums to 1, prior to the subsequent numerical processing (see next 2 items).
2. What is the steady-state final probability of residing in each state, i.e., what is the distribution values of the leading eigenvector? (The leading eigenvector $\mathbf{v}_1$ should contain only real non-negative numbers which sum to 1, and whose associated eigenvalue should have a real value $\lambda_1 = 1$.) Thus, the steady-state distribution of 2 stochastic matrices $M_1$ and $M_2$ can be compared simply by considering the leading eigenvector $\mathbf{v}_1(M_1)$ and $\mathbf{v}_1(M_2)$ from each of the pair of stochastic matrices, and the 2 vectors compared as any 2 discrete distributions of length $n$, with for example the Kullback-Leibler (KL) divergence, etc. (Note that this is a comparison of just a pair of length $n$ vectors, in contrast to @krkeane answer which requires $n$ times comparisons of length $n$ pairs of vectors, corresponding to the complete $n \times n$ matrices.)
3. How fast does the convergence to the steady-state final probabilities occur? I.e., the approximate dynamic behavior. This is typically the ratio of the magnitude of the $2^{nd}$-largest eigenvalue (which is typically complex valued) to the largest eigenvalue (the leading eigenvalue), which should have real value 1: $|\lambda_2|/|\lambda_1| = |\lambda_2|/1 = |\lambda_2|$. The smaller this ratio (actually simply the value of $|\lambda_2|$), the faster the convergence to steady-state. Thus, the 2 stochastic matrices $M_1$ and $M_2$ can be compared by comparing $|\lambda_2(M_1)|$ to $|\lambda_2(M_2)|$

Items 2 and 3 above give 2 comparison conditions, both of which should be passed in order for a pair of stochastic matrices to be considered "effectively equivalent" in terms of their steady-state and approximate dynamic behavior. Item 1 (replacing zeros with very small positive values) is a practical pre-processing for being able to apply items 2 and 3.

Note: Typical convention for a Stochastic matrix is that each of its rows sum to 1. Therefore, the leading eigenvector $\mathbf{v}_1$ mentioned in item 2 is the row eigenvector more commonly known as left eigenvector. However, typical numerical analysis libraries customarily compute the right eigenvector(s), i.e., column eigenvectors. Therefore, in order to use a library which computes the right eigenvector(s), it's simply required to input to it the transpose $M^T$ of the stochastic matrix $M$. Of course, $M^T$ has columns which sum to 1, instead of rows.

Answer 2

The empirical state transition matrices may be thought of as a collection of (row) histograms of a counting process given an initial state. You could compare the histograms using Kullback-Leibler divergence for the Dirichlet distribution. The row divergences could be combined using estimated state probabilities. This aggregate metric would allow you to rank the similarity of mating dances.

---

Edit: after further thought, just vectorize the transition count matrices, and do a single KL-divergence computation. Sample MATLAB code is available here for function KL_dirichlet().

Zeros will remain "problematic" without a prior for each matrix. To resolve this, add pseudocounts to your observed counts, for instance add $\frac{1}{m \times n}$ to each element of a state transition matrix of size $m \times n$. Hopefully you are okay with a Bayesian framework.

---

the diagonal is all zeroes

Given this constraint (apparently $m=n$), omit the diagonal elements from the vectorized representation of the observed state transition (count) matrices. No sense modeling and comparing the probability of something that can't happen.