No way to extract the transition matrix from a Markov chain?

Issue #44 resolved
Paul Moore
created an issue

I'm experimenting with using Lea to analyze Markov chains. As an example, rules for rolling dice in certain games can be modelled as a Markov chain (a simplified example - roll 2 dice, rerolling doubles, is an absorbing chain where (n, m) is an absorbing state if n != m, and transitions to a new rolled pair if n == m).

One thing I'd like to do is extract the transition matrix from a Chain object. I can do this using the private _next_state_lea_per_state attribute, and deconstructing the Alea object for each starting state - but that's a little messy, and uses private information.

Is there a better way that I've missed? And if not, could one be added?

(Actually, Lea in general seems to make it hard to extract the underlying representation of an object - the best way I've found of getting the probability values for an Alea object is list(zip(obj.support, obj.ps)) which is a little non-obvious (neither "support" nor "ps" are attribute names I would have immediately thought of).

Comments (20)

  1. Pierre Denis repo owner

    Hi Paul,

    Happy to see your interest in Markov chain, a topic, which is for sure not much elaborated in Lea!

    Let me first answer on the last point. Your need is addressed by methods described in the third part of the tutorial. You could for example replace your non-obvious expression by obj.pmf_tuple (well, you get a tuple, not a list!). I think that I should move this section in the first part of the wiki because it's seemingly a basic need of many users.

    Now, about the main question, indeed, you didn’t miss anything. There is nothing out-of-the-box to get the transition matrix easily. Python is hopefully tolerant on accessing private attribute : as you have discovered already, you can get the transition NxN matrix as N tuple of N probabilities by

    tuple(a.ps for (_,a) in your_markov._next_state_lea_per_state)
    

    and get the associated state names with

    tuple(s for (s,_) in your_markov._next_state_lea_per_state)
    

    Definitely, I’m willing to add such access methods in Lea package. Now, is it exactly what you want? What method names would be meaningful? BTW, would it be sensible to rename _next_state_lea_per_state and promote it as public (and duly document it)?

    Yes,… It’s now time to give your wish list. :)

  2. Paul Moore reporter

    I'm only really playing with Markov chains. There are some specific things I'm interested in, but they may not be particularly appropriate for Lea (and once I have the transition matrix, they aren't that hard to produce manually).

    Specifically:

    1. Is a Markov chain absorbing? (The ones I'm working with are, by construction, but checking that fact before I continue is not a bad idea).
    2. Extract the number of transient states and the number of absorbing states.
    3. Put the transition matrix into canonical form (trying to lay this out messes up Bitbucket's formatting, but basically 4 quadrants, Q in the upper left, R in the upper right, a zero matrix in the lower left and an identity matrix in the lower right - you probably already know this) and extract Q and R.
    4. Calculate the fundamental matrix (I-Q)**(-1)

    The matrix arithmetic really needs numpy (at least, it seems silly to avoid numpy) so would probably have to be optional in Lea anyway. So (4) probably isn't suitable for Lea.

    Suggestions:

    # The name here is problematic as currently the states method does
    # something different. But I can't think of a better name for a list of
    # the states (as Python objects, not probability distributions)
    all_states = chain.states()
    # Boolean
    chain.is_absorbing()
    # A subset of all_states
    absorbing = chain.absorbing_states()
    transition_matrix = chain.transition_matrix()
    # This is a bit long, not sure if it can be shortened?
    q, r = chain.transition_matrix_canonical_decomposition()
    
  3. Pierre Denis repo owner

    OK. That sounds interesting and general enough to be integrated in Lea. Let me work on this. Hopefully, I'll commit soon something that you could review.

  4. Pierre Denis repo owner

    Hi Paul,

    Sorry for the long delay! I’ve eventually done something on the topic of absorbing MC (I know that you were able to manage, by hard ways, to get the information you need).

    Everything is currently on the branch dev_upgrade_markov_chains of lea (to be merged on the main branch). There are 3 new methods on markov.Chain : absorbing_mc_info, matrix and reachable_states. The first is the main one, the last two are auxiliary but these could be helpful as such. These methods are documented in the code. Here is an example of use, for an absorbing MC:

    mc = markov.chain_from_matrix(
                   ( 'A'  , 'B', 'C' ),
            ( 'A', ( 0.5 , 0.3 , 0.2 )),
            ( 'B', ( 0.0 , 1.0 , 0.0 )),
            ( 'C', ( 0.3 , 0.5 , 0.2 )))
    print(mc.states)
    # -> ('A', 'B', 'C')
    print(mc.matrix())
    # -> ((0.5, 0.3, 0.2), (0.0, 1.0, 0.0), (0.3, 0.5, 0.2))
    print(mc.absorbing_mc_info())
    # -> (True, ('A', 'C'), ('B',), ((0.5, 0.2), (0.3, 0.2)), ((0.3,), (0.5,)))
    

    Notes:

    1. Unlike your suggestions, I’ve preferred to return all data in one single method (absorbing_mc_info) because these are closely interdependent.
    2. So far, numpy is not used in the methods; also, the fundamental matrix is not calculated. However, as I’ve done for other dependencies, I may envision an optional support to numpy, for example by adding an optional Boolean arg on absorbing_mc_info for asking to get Q/R matrices as numpy arrays and provide (I-Q)**(-1) at the end of the output tuple.
    3. Of course, there remain work to do, at least: adding tests and documenting the new features on the wiki.

    I would appreciate if you can give some feedback on these points before I go on.

  5. Paul Moore reporter

    Hi, The above looks good to me. I haven't done much more since my original experiments - even though I had to play with internal attributes, that was enough to solve the problem I was working on at the time. But my solution would have been much cleaner with the new methods you're adding here. Optional numpy-based calculation of the fundamental matrix sounds like a reasonable way to go (although of course, it's possible to do that myself using the data you currently return in absorbing_mc_info).

  6. Log in to comment