How to better my code?

I have codes below:

transitionProbability(mc^2, "1", "0")
transitionProbability(mc^2, "1", "1")
transitionProbability(mc^2, "1", "2")
transitionProbability(mc^2, "1", "3")
transitionProbability(mc^2, "1", "4")

transitionProbability(mc^2, "2", "0")
transitionProbability(mc^2, "2", "1")
transitionProbability(mc^2, "2", "2")
transitionProbability(mc^2, "2", "3")
transitionProbability(mc^2, "2", "4")

transitionProbability(mc^2, "3", "0")
transitionProbability(mc^2, "3", "1")
transitionProbability(mc^2, "3", "2")
transitionProbability(mc^2, "3", "3")
transitionProbability(mc^2, "3", "4")

transitionProbability(mc^2, "4", "0")
transitionProbability(mc^2, "4", "1")
transitionProbability(mc^2, "4", "2")
transitionProbability(mc^2, "4", "3")
transitionProbability(mc^2, "4", "4")

Then I need to do mc^10, 50, 100, 200 again......

I recommend this material.
"21 Iteration | R for Data Science" https://r4ds.had.co.nz/iteration.html

This material doesn't help me. I want to use apply family function to transitionProbability.

Iteration is the general concept, of which apply functions are one of several examples.

I will quote from the linked material

The apply family of functions in base R ( apply() , lapply() , tapply() , etc) solve a similar problem, but purrr is more consistent and thus is easier to learn.

If anyone can help me write down the codes, I will appreciate it. I am learning purr and map function now, but I still get stuck. I still prefer apply family function, since I used before and it worked out. But in this new problem I face, the function I want to deal with is not what I create. 'transitionProbability' is from the package. Furthermore, this function has 3 arguments that I need to fill in.
First, mc^n, here n = 2, 10, 50, 100, 200
Second is 0,1,2,3,4
Third is also 0,1,2,3,4.

Previously, when I use apply, the function I used has only 1 argument.

But if someone use other iteration method, I can accept!

purrr::map2 allows you to map a function over two data structures. I am not aware of any purrr functions that go beyond two, so you probably have to add another iteration/for-loop etc. to iterate across your different powers. Here's some examples to guide you along:

mc <- 0.5

# Using an anonymous function 
# here I am multiplying the 3 things together
map2_dbl(c(5, 4, 3, 2, 1), c(6, 5, 4, 3, 2), function(x , y)(mc^2 * x * y))

# Using a pre-defined function "mean" 
# in your case, it's transitionProbability
map2_dbl(c(5, 4, 3, 2, 1), c(6, 5, 4, 3, 2), function(x, y)(mean(mc^2 + x + y)))

A cumbersome for-loop:

n.vec <- c(2, 10, 50, 100, 200)

# to store the values
mat <- matrix(ncol = length(n.vec), nrow = 5)

for(i in n.vec){
mat[, which(n.vec == i)] <- map2_dbl(c(5, 4, 3, 2, 1), c(6, 5, 4, 3, 2), 
                                     function(x, y)(mean(mc^i * x * y)))
}

mat

     [,1]        [,2]         [,3]         [,4]         [,5]
[1,]  7.5 0.029296875 2.664535e-14 2.366583e-29 1.866905e-59
[2,]  5.0 0.019531250 1.776357e-14 1.577722e-29 1.244603e-59
[3,]  3.0 0.011718750 1.065814e-14 9.466331e-30 7.467618e-60
[4,]  1.5 0.005859375 5.329071e-15 4.733165e-30 3.733809e-60
[5,]  0.5 0.001953125 1.776357e-15 1.577722e-30 1.244603e-60
1 Like

Thank you so much! Your code works very close to what my desired. Here I got 5 numbers, because 5x6, 4x5,... What I want is a 5 by 5 table. The first row is 5x6, 5x5, 5x4...

I don't know what markov chains are but I tried to recreate something as close as possible using the example given by the markovchain package.

library(markovchain)

# The markovchain object from their CRAN
statesNames <- c("1", "2", "3")
markovB <- new("markovchain", states = statesNames, transitionMatrix =
                 matrix(c(0.2, 0.5, 0.3, 0, 1, 0, 0.1, 0.8, 0.1), nrow = 3,
                        byrow = TRUE, dimnames=list(statesNames,statesNames)),
               name = "A markovchain Object" 
)    

# I was unable to create a 5x5 markovchain object so I use a shorter version of what you want to do

n.vec <- c(2, 10, 50)

# to store the values
mat <- matrix(ncol = length(n.vec), nrow = 3)

for(i in n.vec){
  mat[, which(n.vec == i)] <- map2_dbl(c(1, 1, 1), c(0, 1, 2), 
                                       function(x, y)(transitionProbability(markovB^i, x, y)))
}

mat # if you want to populate the matrix row-wise, use [which(n.vec == i), ] 

     [,1]         [,2]      [,3]
[1,]   NA 7.000000e-02 0.8400000
[2,]   NA 9.864400e-06 0.9999773
[3,]   NA 5.613498e-25 1.0000000

@Tung, I have forgotten most of my Markov chain knowledge, but I think if you have the transition matrix (mc in your code) , raising it up to the number of steps would give you the transition probabilities already. I am not following why you need to call that function multiple times and create a matrix out of the results. Isn't mc^2 enough to give you all the probabilities?

Anyway, this is probably one more possibility of using outer.

outer(c("0", "1", "2", "3", "4"), c("0", "1", "2", "3", "4"), Vectorize(function(x, y) transitionProbability(mc^2, x, y)))

But as far as I can see, it'll be same as (mc^2)@transitionMatrix.


@bayesian, briefly answering two of your points, in case it helps.

  1. purrr::pmap is supposed to go beyond two. It has support for ...1, ...2 etc. I almost never used it, so can't give an example. Sorry.
  2. Markov chain is a type of stochastic process. In the simplest case, given present future is independent of past (conditional probability of next step given the state history up to current state is same as conditional probability of next step given only the current step).

Edit (in reply to post #13)

I mentioned this previously:

You should have used this only. There is no need to create a matrix using outer, map2, pmap, for or any other method, as what you have is enough to check convergence. Probably some wrapper like this will be enough, but to be fair, probably unnecessary:

f <- function(markov_object, number_of_steps)
{
    probabilities_after_steps <- markov_object ^ number_of_steps
    return(probabilities_after_steps@transitionMatrix)
}

Then call it as f(mc, 10), f(mc, 50), etc. to get matrices for that many steps.


Please note that this is not a site to solve your problem for you. I gave you a solution and its explanation (here and in previous thread), and from that you should have been able to make a function to repeat for multiple steps. At least, you are supposed to show us some of your effort where you may or may not have failed.

And, for your future threads, please provide a minimal reproducible example, and ask only one concise question. Ideally, each thread should be about only one question, without multiple follow up questions. Asking explanations about a solution is fine, but not a new question. Please take a look at #meta for FAQ's.

1 Like

If you don't want to use matrix, probably simply double loop will do the job:

for (i in 1:4)
   for (j in 0:4)
     {
      x <- as.character(i)
      y <- as.character(j)
      transitionProbability(mc^2, x,  y)
     }

should work properly. Now you can make function from these 2 loops and just run it 4 times:

probabilityCalc <- function(z = 2)
{
for (i in 1:4)
   for (j in 0:4)
     {
      x <- as.character(i)
      y <- as.character(j)
      transitionProbability(mc^z, x,y)
     }
}

probabilityCalc(2)
probabilityCalc(10)
probabilityCalc(50)
probabilityCalc(100)
probabilityCalc(200)

You are right. After forming a Markov chain (say, the name is mc.) Then mc^k means the transition probability after k steps.

The reason I want to multiple many times is to show simulation and convergence. After c(2, 10, 50, 100, 200) steps, they should be convergent. (It's kind of like law of large numbers, but they are the same thing. But the idea is similar.)

> map2_dbl("1", c("0", "1", "2", "3", "4"), function(x , y)(transitionProbability(mc^2, x, y)))
[1] 0.19 0.39 0.14 0.14 0.14

Thank you for @bayesian. This code means given my initial state is 1, after 2 steps (mc^2), what my transition probability matrix would be. Thus 0.19 means from state 1 to state 0 after 2 steps; 0.39 means from state 1 to state 1 after 2 steps, etc.

Now the above code just gives me a vector, since the initial state is fixed at 1. I have 5 states (0,1,2,3,4). Thus, what I want is a 5 by 5 matrix. For example, the (3,3) entry should be transitionProbability(mc^2, "2", "2"). If it could be done, then this 5 by 5 matrix is for markov chain after 2 steps.

Then, I need to do the same 5 by 5 matrix after these steps: c(2, 10, 50, 100, 200) to show the "convergence" (or steady state)

Thank you. Your code makes sense to me. But after running your code, I don't get any output.

@Yarnabrina, Your code works! One minor questions is how to change the column name and row name to 0-4, currently it shows 1-5? The second question is what I asked before, I still need to do mc^2, 10, 50, 100, 200. In this case, what iteration method should I use?

> outer(c("0", "1", "2", "3", "4"), c("0", "1", "2", "3", "4"), Vectorize(function(x, y) transitionProbability(mc^2, x, y)))
     [,1] [,2] [,3] [,4] [,5]
[1,] 1.00 0.00 0.00 0.00 0.00
[2,] 0.19 0.39 0.14 0.14 0.14
[3,] 0.19 0.14 0.39 0.14 0.14
[4,] 0.19 0.14 0.14 0.39 0.14
[5,] 0.19 0.14 0.14 0.14 0.39

I don't know what output you are looking at. If transitionProbability function gives values which u need, you can add matrix to keep track of outputs i.e.


probabilityCalc <- function(z  =  2)
{
 output <- matrix(1:20, nrow = 4)

for (i in 1:4)
   for (j in 0:4)
     {
      x <- as.character(i)
      y <- as.character(j)
      output[i,j+1] <- transitionProbability(mc^z, x,y) #matrix is numereted from 1
     }

    return(output)
}

Output1 <- probabilityCalc(2)
View(Output1)
...

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.

If you have a query related to it or one of the replies, start a new topic and refer back with a link.