#
# Written by:
# -- 
# John L. Weatherwax                2009-04-21
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
#-----

set.seed(0)

Markov = function( N=100, initial.value=1, P ){
  X = numeric(N)
  X[1] = initial.value + 1 # convert from dollars to states i.e. States 0:5 <=> subscripts 1:6 (0 = no money; 5 = all money)
  n = nrow(P)
  for( i in 2:N ){
    X[i] = sample( 1:n, size=1, prob=P[X[i-1],] )
  }
    X - 1 # return the dollar amount 
}

P = matrix( 0, nrow=6, ncol=6 )
P[1,1] = 1.;
P[2,1] = 0.5; P[2,3] = 0.5;
P[3,2] = 0.5; P[3,4] = 0.5;
P[4,3] = 0.5; P[4,5] = 0.5;
P[5,4] = 0.5; P[5,6] = 0.5;
P[6,6] = 1.;

# Simulate our Markov chain several times: 
for( ii in 1:10 ){ 
  print( Markov( N=15, initial.value=2, P ) )
}