Search code examples
rgraphmarkov-models

Markov Model diagram directly from data (makovchain or deemod package?)


I want to read a bunch of factor data and create a transition matrix from it that I can visualise nicely. I found a very sweet package, called 'heemod' which, together with 'diagram' does a decent job.

For my first quick-and-dirty approach, a ran a piece of Python code to get to the matrix, then used this R sniplet to draw the graph. Note that the transition probabilities come from that undisclosed and less important Python code, but you can also just assume that I calculated it on paper.

library('heemod')
library('diagram')
mat_dim <- define_transition(
  state_names = c('State_A', 'State_B', 'State_C'),
  .18, .73, .09, 
  .22, .0, .78,
  .58, .08, .33);
plot(mat_dim)

However, I would like to integrate all in R and generate the transition matrix and the graph within R and from the sequence data directly.

This is what I have so far:

library(markovchain)
library('heemod')
library('diagram')

# the data --- this is normally read from a file
data = c(1,2,1,1,1,2,3,1,3,1,2,3,1,2,1,2,3,3,3,1,2,3,2,3,1,2,3,3,1,2,3,3,1)
fdata = factor(data)
rdata = factor(data,labels=c("State_A","State_B","State_C"))

# create transition matrix
dimMatrix = createSequenceMatrix(rdata, toRowProbs = TRUE)
dimMatrix

QUESTION: how can I transfer dimMatrix so that define_transition can process it?

mat_dim <- define_transition( ??? );
plot(mat_dim)

Any ideas? Are there better/easier solutions?


Solution

  • The input to define_transition seems to be quite awkward. Perhaps this is due to my inexperience with the heemod package but it seems the only way to input transitions is element by element.

    Here is a workaround

    library(heemod)
    library(diagram)
    

    first convert the transition matrix to a list. I used rounding on the digits which is optional. This corresponds to the ... variables in define_transition

    lis <- as.list(round(dimMatrix, 3))
    

    now add to the list all other named arguments you wish:

    lis$state_names = colnames(dimMatrix)
    

    and now pass these arguments to define_transition using do.call:

    plot(do.call(define_transition, lis))
    

    enter image description here

    Update: to the question in the comments:

    lis <- as.list(t(round(dimMatrix, 3)))
    lis$state_names = colnames(dimMatrix)
    plot(do.call(define_transition, lis))
    

    enter image description here

    The reasoning behind do.call

    The most obvious way (which does not work here) is to do:

    define_transition(dimMatrix, state_names = colnames(dimMatrix))
    

    however this throws an error since the define_transition expects each transition to be supplied as an argument and not as a matrix or a list. In order to avoid typing:

    define_transition(0.182, 0.222,..., state_names = colnames(dimMatrix))
    

    one can put all the arguments in a list and then call do.call on that list as I have done.