# statistical comparison of two markov chain transition matrices

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.

• "Part of the problem we're running into is that both matrices have a lot of zeroes in them" can you expand on this? Are you manually deciding where the zeros should go? Or are you doing MLE/etc, and the optimum happens to be sparse? Or are you explicitly regularizing for sparsity via, say, a lasso-like penalty? Apr 16, 2023 at 14:45
• Hi John, thanks for you're reply. We aren't doing anything as sophisticated as what you're saying, basically we're just recording two different populations of the same animal, their mating dance, and trying to understand whether it's the same dance. We've built two transition matrices basically just by the frequencies that they pass from one behaviour to another, and trying to understand whether they're doing the same dance or not. The zeroes are just there because they never passed directly between those two behaviours, also the diagonal is all zeroes. Apr 18, 2023 at 0:14
• Since you have many zero's, an interesting aspect is how many states you have and whether the states very different. Also, what do the states represent? For instance, if the states are some quantity and a transition from X to Y is very similar to a transition from X to Y+1, then this is lost in your transition matrix with the many zero's. A hypothesis test or measure of statistical significance might be having little power and possibly some parametric model might be better. Apr 19, 2023 at 7:16
• Also, are you estimating the transition matrices independently for each time, or are you averaging over the entire range of times (computing a single transition matrix)? Apr 19, 2023 at 7:20
• We are right now averaging over a period of time, though we may break up that period of time if we need to. Each state is a particular behaviour, so it wouldn't make sense to talk about Y+1. Thanks for your answer. Apr 20, 2023 at 20:48

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.

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.

• Related is the G-test, which is a comparison of the KL-divergence between the observations and and a (estimated) model. Apr 19, 2023 at 7:12