Time is continuous time. I have a 3D state space, and transition rates across all of these.

Using the transition rates, I can compute a generator matrix A associated with the markov process. My A satisfies all requirements on a generator matrix:

  1. A_ii < 0 and A_ij >=0 for all i and j
  2. Also, rows sum up to zero.

I want to compute the implied stationary distribution g given by

dot(A.T,  g) = z
sum(g) = 1

where z is a zero vector. To ensure the second condition, I set z_i = 1 for some i, A_ii = 1, A_ij = 0 for all j != i.

Now, my problem is that I don't always get to invert A, despite it satisfying all the conditions for a generator matrix. Python just gives me

scipy.sparse.linalg.dsolve.linsolve.MatrixRankWarning: Matrix is exactly singular

I have a suspect: My grid size. Let's call one of the three dimensions of the state space u. For a fixed u.min(), u.max(), if I increase the grid size enough, I can usually invert A. For the order of magnitudes: If the transition rate across grid points is around 100, inverting A fails. When it's only around 10, it works.

Now, I can't always throw more grid points at the problem. Also, I'm not sure of the fundamental problem behind why the grid size matters where. Could anyone shed some light on what's happening here, and whether there's an alternative way forward?


1 Answer 1


It turns out that there was some floating point error. The key insight to get around this was that the stationary distribution is invariant to the speed of convergence.

So then I just rescaled A2 = A/x, for x > 0. Depending which way the floating point errors went, x was either above or below 1. Inverting A2 was then possible without problems.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.