rw.sim <- function( k, p, theta.start, n.sim, seed ) { set.seed( seed ) theta <- rep( 0, n.sim + 1 ) theta[ 1 ] <- theta.start for ( i in 1:n.sim ) { theta[ i + 1 ] <- move( k, p, theta[ i ] ) } return( table( theta ) ) } move <- function( k, p, theta ) { repeat { increment <- sample( x = c( -1, 0, 1 ), size = 1, prob = p ) theta.next <- theta + increment if ( abs( theta.next ) <= k ) { return( theta.next ) break } } } p <- c( 1, 1, 1 ) / 3 k <- 5 theta.start <- 0 seed <- 9626954 rw.sim( k, p, theta.start, 10, seed ) # -2 -1 0 # 3 4 4 rw.sim( k, p, theta.start, 100, seed ) # -5 -4 -3 -2 -1 0 1 2 3 # 15 22 18 11 11 8 2 9 5 rw.sim( k, p, theta.start, 1000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 100 113 87 65 66 68 83 95 130 124 70 rw.sim( k, p, theta.start, 10000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 768 1052 1054 1013 980 909 877 916 939 906 587 rw.sim( k, p, theta.start, 100000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 6663 10111 10201 9864 9841 9587 9380 9305 9535 9337 6177 rw.sim( k, p, theta.start, 1000000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 65232 97285 97862 97707 96845 96007 96624 96173 96415 96317 63534 rw.sim( k, p, 5, 100000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 6661 10109 10197 9859 9839 9586 9382 9307 9539 9341 6181 p <- c( 0.2, 0.3, 0.5 ) rw.sim( k, p, 0, 10, seed ) # -2 -1 0 # 5 5 1 rw.sim( k, p, 0, 100, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 1 7 13 9 13 4 4 6 9 18 17 rw.sim( k, p, 0, 1000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 1 7 13 9 16 17 55 88 150 296 349 rw.sim( k, p, 0, 10000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 1 7 17 25 61 117 257 564 1395 3370 4187 rw.sim( k, p, 0, 100000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 13 28 64 140 382 946 2236 5509 13945 34454 42284 rw.sim( k, p, 0, 1000000, seed ) # -5 -4 -3 -2 -1 0 1 2 3 4 5 # 77 247 584 1413 3579 8991 22214 54773 137010 342069 429044