Am Mittwoch, den 27.07.2011, 19:59 -0400 schrieb R. Michael Weylandt : > Some more skilled folks can help with the curve fitting, but the general > answer is yes -- R will handle this quite ably.
Great to read that. > Consider the following code: > > <<-------------------------------------->> > n = 1e5 > length = 1e5 > > R = matrix(sample(c(-1,1),length*n,replace=T),nrow=n) > R = apply(R,1,cumsum) ## this applies cumsum `row-wise' to R and will make > your life INFINITELY better > R = cbind(rep(0,n),R) ## Now each row is a random walk as you desired. > > <<---------------------------------------->> > > There are actually even faster ways to do what you are asking for, but this > introduces you to some useful R architecture, above all the apply function. Thank you very much. I realized the the 0 column is not need when summing this up. Additionally I posted the wrong example code and I actually am only interested how long it stays negative from the beginning. > To see how long the longest stretch of negatives in each row is, the > following is a little sneaky but works pretty well: > > countNegative = apply(R,1,function(x){which.max(table(cumsum(x>=0))}) > > then you can study these random numbers to do whatever with them. > > The gist of this code is that it counts how many positive number have been > seen up to each point: for any stretch this doesn't increase, you must be > negative, so this identifies the longest such stretch on each row and > records the length. (It may be off by one so check it on a smaller R matrix. That is a great example. It took me a while what `table()` does here but thanks to your explanation I finally understood it. […] > So all together: > > <<-------------------------------------->> > n = 1e3 > length = 1e3 > > R = matrix(sample(c(-1,1),length*n,replace=T),nrow=n) > R = apply(R,1,cumsum) ## this applies cumsum `row-wise' to R and will make > your life INFINITELY better > R = cbind(rep(0,n),R) ## Now each row is a random walk as you desired. > fTemp <- function(x) { > return(max(table(cumsum(x>=0)))) > } > countNegative = apply(R,1,fTemp) > mu = mean(countNegative) > sig = sd(countNegative)/sqrt(length(countNegative)) > > <<---------------------------------------->> > > This runs pretty fast on my laptop, but you'll need to look into the > memory.limit() function if you want to increase the simulation parameters. > There are much faster ways to handle the simulation as well, but this should > get you off to a nice start with R. > > Hope this helps, It did. Thank you again for the kind and elaborate introduction. Trying to run your example right away froze my system using `n = 1000` and `length = 1e5` [1]. So I really need to be careful how big such a matrix can get. One thing is to use integers as suggested in [2]. My current code looks like the following. -------- 8< -------- code -------- >8 -------- f4 <- function(n = 100000, # number of simulations length = 100000) # length of iterated sum { R = matrix(sample(c(-1L,1L),length*n,replace=T),nrow=n) R = apply(R,1,cumsum) ## this applies cumsum `row-wise' to R and will make your life INFINITELY better fTemp <- function(x) { if (x[1] >= 0 ) { return(1) } for (i in 1:length-1) { if (x[i] < 0 && x[i + 1] >= 0) { return(as.integer(i/2 + 2)) # simple random walks only hit 0 on even “times” } } } countNegative = apply(R,2,fTemp) tabulate(as.vector(countNegative), length) } -------- 8< -------- code -------- >8 -------- 1.I could actually avoid `cumsum()` half the time, when the first entry is already positive. So I am still looking for a way to speed that up in comparison to a simple two loops scenario. 2. The counting of the length how long the walk stayed negative is probably also inefficient and I should find a better way on how to return the values. I am still thinking about both cases, but to come up with vectoriazations of the problem is quite hard. So I welcome any suggestions. ;-) Thanks, Paul [1] http://bugs.debian.org/cgi-bin/bugreport.cgi?bug=635832 [2] http://bugs.debian.org/cgi-bin/bugreport.cgi?bug=635832#10
signature.asc
Description: This is a digitally signed message part
______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.