R help - Dec 2009 - cycling k times a realization of a random walk.....problems..

hello R-masters.

i have an R-issue here that i don't know if you'd wish to help me? about
it:

briefly i'd like to generate many (say hundred) realizations of a random
walk, execute a few operations on each of them (mean time of return), and graph
each realization on the same plot.
IN OTHER WORDS I'D LIKE TO IMPOSE A LOOPING CYCLE TO THE COMMAND NOT THE
ARGUMENT OF THE COMMAND.

for some of these questions i have already a partial answer: my main problem
here is automatizing the process for 100 times.

my random walk generating function is:

rwalk <- function(n,p, x0=0)
? ? ? ???{
? ? ? ? ? x <- rbinom(n,1,p)
? ? ? ? ? x[which(x==0)] <- -1
? ? ? ? ? y<-c(x0,x)
? ? ? ? ? y <- cumsum(y)
? ? ? ? ???}

THIS IS THE CODE THAT I'D LIKE TO REITERATE

a<- rwalk(n,p,x0=0) #whish to generate 100 of those 
#and for each calculate the following
o<-which(a==0)
N<-length(o) # number of times the walk returns to 0
Nn<-(o[-1]-N) # number of steps to get to 0
V<-mean(Nn)? # mean time of return to 0


par(mfrow=c(1,1))
plot(0:n, a, type= "l", main="...", xlab="tempo",
ylab="stato",
ylim=c(...), lwd=3, col="red", lty=1)

par(new=T) 

plot(0:n, b, type= "l", main="", xlab="",
ylab="",
ylim=c(), lwd=3, col="green", lty=2)

#technically i've tryed to fuse all those parts in a function reiterating
k=100 times with a for cycle, but with no succes, and i dare there must be
something really disturbing in the code below ....

function(k){
for (i in 1:k){
a[i]<-rwalk(1e+04,0.5,0)
o[i]<-which(a[i]==0)
N[i]<-length(o[i]) 
Nn[i]<-(o[i][-1]-N[i]) 
V[i]<-mean(Nn[i]) 
w<-c(V[1]:V[k])
M<-mean(w) #mean over all the realizations
par(mfrow=c(1,1))
plot(0:n, a[1], type= "l", main="...",
xlab="tempo", ylab="stato",
ylim=c(...), lwd=3, col="red", lty=1)

par(new=T) 

plot(0:n, a[i], type= "l", main="", xlab="",
ylab="",
ylim=c(...), lwd=3, col=c(1:i), lty=2)

}
}???#NOTE: I'VE ALSO TRIED TO REITERATE THE WALK WITH THE rep()
? ? # t<-rep(rwalk(1e+04,0.5,0),100)
? ???# then of course i have the problem of splicing this ridicolously???#one
milion long vector in 100 bits of tenthousand, and still goig through all the
calculations and plotting i meant before....

i wonder what......

THANK YOU FOR YOUR PRECIUOS ATTENTION!!! 

?

On Fri, 2009-12-04 at 05:54 -0800, Federico Bonofiglio
wrote:
> hello R-masters.
> 
> i have an R-issue here that i don't know if you'd wish to help me
> about it:
> 
> briefly i'd like to generate many (say hundred) realizations of a
> random walk, execute a few operations on each of them (mean time of
> return), and graph each realization on the same plot.
> IN OTHER WORDS I'D LIKE TO IMPOSE A LOOPING CYCLE TO THE COMMAND NOT
> THE ARGUMENT OF THE COMMAND.
> 
> for some of these questions i have already a partial answer: my main
> problem here is automatizing the process for 100 times.
> 
> my random walk generating function is:
> 
> rwalk <- function(n,p, x0=0)
>          {
>           x <- rbinom(n,1,p)
>           x[which(x==0)] <- -1
>           y<-c(x0,x)
>           y <- cumsum(y)
>            }
<snip />
 #NOTE: I'VE ALSO TRIED TO REITERATE THE WALK WITH THE
rep()
>     # t<-rep(rwalk(1e+04,0.5,0),100)
>      # then of course i have the problem of splicing this
> ridicolously   #one milion long vector in 100 bits of tenthousand, and
> still goig through all the calculations and plotting i meant
> before....
One way of doing away with the loop, is to just get all the rbinoms you
want and format it appropriately, before cumsum: Try this:

rwalk <- function(n, p, x0 = 0, sim = 1)
         {
x <- rbinom(n * sim, 1, p)
x[!x] <- -1
x <- matrix(x, nrow = sim, byrow = TRUE)
x <- cbind(rep(x0, sim), x)
y <- t(apply(x, 1, cumsum))
return(y)
}

sim is the extra argument, indicating the number of walks you want. Set
it to 100 to get what you need. Each row in the returned matrix is a
random walk of the form you describe.

You can then use the apply() function to compute your statistics on each
row (walk) and aggregate the results, e.g.
> set.seed(123)
> system.time(walks <- rwalk(10000, p = 0.5, x0 = 0, sim = 100))
   user  system elapsed 
  0.747   0.010   0.759 
> tmp[1:10, 1:10] # first 10 walks, first 10 data in each
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
 [1,]    0    1    0    1    0   -1    0   -1   -2    -3
 [2,]    0   -1   -2   -3   -4   -3   -4   -5   -6    -7
 [3,]    0    1    2    1    0   -1    0   -1   -2    -3
 [4,]    0    1    0   -1    0   -1    0   -1    0    -1
 [5,]    0    1    0   -1    0   -1   -2   -3   -4    -3
 [6,]    0   -1   -2   -3   -4   -5   -4   -5   -4    -3
 [7,]    0   -1    0    1    0   -1   -2   -1    0    -1
 [8,]    0    1    0   -1   -2   -3   -2   -3   -2    -1
 [9,]    0   -1   -2   -3   -4   -3   -4   -5   -4    -3
[10,]    0   -1   -2   -3   -2   -3   -4   -5   -4    -3
> apply(tmp[1:10, ], 1, function(x) which(x == 0))
[[1]]
[1] 1 3 5 7

[[2]]
[1] 1

[[3]]
[1] 1 5 7

[[4]]
[1] 1 3 5 7 9

[[5]]
[1] 1 3 5

[[6]]
[1] 1

[[7]]
[1]  1  3  5  9 11 17 23 27 29

[[8]]
[1] 1 3

[[9]]
[1] 1

[[10]]
[1] 1

That might make it easier to work with as you only need to generate the
walks once and then develop each summary function independently. If you
start with a small example, say n = 10 and sim = 10, you should be able
to follow the computations and check they are correct.

HTH

G
-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson             [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,          [f] +44 (0)20 7679 0565
 Pearson Building,             [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London          [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT.                 [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%