**Freakonometrics - Tag - R-english**, and kindly contributed to R-bloggers)

Markov chains is a very interesting and powerful tool. Especially for parents. Because if you think about it quickly, most of the games our kids are playing at are Markovian. For instance, snakes and ladders…

It is extremely easy to write down the transition matrix, one just need to define all snakes and ladders. For the one above, we have,

n=100 M=matrix(0,n+1,n+1+6) rownames(M)=0:n colnames(M)=0:(n+6) for(i in 1:6){diag(M[,(i+1):(i+1+n)])=1/6} M[,n+1]=apply(M[,(n+1):(n+1+6)],1,sum) M=M[,1:(n+1)] starting=c(4,9,17,20,28,40,51,54,62, 64,63,71,93,95,92) ending =c(14,31,7,38,84,59,67,34,19, 60,81,91,73,75,78) for(i in 1:length(starting)){ v=M[,starting[i]+1] ind=which(v>0) M[ind,starting[i]+1]=0 M[ind,ending[i]+1]=M[ind,ending[i]+1]+v[ind]}

So, why is it important to have a Markov Chain ? Because, once you’ve noticed that you had a Markov Chain game, you can derive anything you want. For instance, you can get the distribution after some turns,

powermat=function(P,h){ Ph=P if(h>1){ for(k in 2:h){ Ph=Ph%*%P}} return(Ph)} initial=c(1,rep(0,n)) COLOR=rev(heat.colors(101)) u=1:sqrt(n) boxes=data.frame( index=1:n, ord=rep(u,each=sqrt(n)), abs=rep(c(u,rev(u)),sqrt(n)/2)) position=function(h=1){ D=initial%*%powermat(M,h) plot(0:10,0:10,col="white",axes=FALSE, xlab="",ylab="",main=paste("Position after",h,"turns")) segments(0:10,rep(0,11),0:10,rep(10,11)) segments(rep(0,11),0:10,rep(10,11),0:10) for(i in 1:n){ polygon(boxes$abs[i]-c(0,0,1,1), boxes$ord[i]-c(0,1,1,0), col=COLOR[min(1+trunc(500*D[i+1]),101)], border=NA)} text(boxes$abs-.5,boxes$ord-.5, boxes$index,cex=.7) segments(c(0,10),rep(0,2),c(0,10),rep(10,2)) segments(rep(0,2),c(0,10),rep(10,2),c(0,10))}

Here, we have the following (note that I assume that once 100 is reached, the game is over)

Assume for instance, that after 10 turns, your daughter accidentally drops her pawn out of the game. Here is the theoretical (unconditional) position of her pawn after 10 turns,

so, if she claims she was either on 58, 59 or 60, here are the theoretical probabilities to be in each cell after 10 turns,

> h=10 > (initial%*%powermat(M,h))[59:61]/ + sum((initial%*%powermat(M,h))[59:61]) [1] 0.1597003 0.5168209 0.3234788

i.e. it is more likely she was on 59 (60th cell of the vector since we start in 0). You can also look at the distribution of the number of turns (at first with only one player).

distrib=initial%*%M game=rep(NA,1000) for(h in 1:length(game)){ game[h]=distrib[n+1] distrib=distrib%*%M} plot(1-game[1:200],type="l",lwd=2,col="red", ylab="Probability to be still playing")

Once you have that survival distribution, you have the expected number of turns to finish the game,

> sum(1-game) [1] 32.16499

i.e. in 33 turns, on average, your daughter reaches the 100 cell. But in 50% of the games, it takes less than 29,

But assuming that you are playing with your daughter, and that the game is over once one player reaches the 100 cell, it is possible to get the survival distribution of the first time one of us reaches the 100 cell,

plot((1-game[1:200])^3,type="l",lwd=2,col="blue", ylab="Probability to be still playing (2 players)")

Here, the expected number of turns before ending the game is

> sum((1-game)^2) [1] 23.40439

And if you ask your son to join the game, the survival distribution function is

plot((1-game[1:200])^3,type="l",lwd=2,col="purple", ylab="Probability to be still playing (3 players)")

i.e. the expected number of turns before the end is now

> sum((1-game)^3) [1] 20.02098

**leave a comment**for the author, please follow the link and comment on his blog:

**Freakonometrics - Tag - R-english**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...