0

我正在模拟从坐标（0,0）开始的随机游走。当我使用循环执行此操作时，效果很好：

require(ggplot2)
n <- 1000   #number of walks 

# first solution, w/ loop... works but is SLOOOW
coord <- data.frame (x=0, y=0, step=0) #origin
for (i in 1:n){
  dir <- sample(c("w", "e", "n", "s"), 1) #random direction
  step <- sample(1:4, 1) #how far to go in each walk
  startx <- coord[nrow(coord), 1]
  starty <- coord[nrow(coord), 2]
  endx <- ifelse (dir=="w", startx-step, ifelse(dir=="e", startx+step, startx))
  endy <- ifelse (dir=="n", starty+step, ifelse(dir=="s", starty-step, starty))
  newcoord <- data.frame (x=endx, y=endy, step=step)
  coord <- rbind(coord, newcoord)
}
rw <- ggplot(coord, aes(x=x, y=y))
rw + geom_path() + 
  ggtitle(paste(n, "walks")) + 
  geom_point(aes(x=0, y =0), color="green", size=I(5)) +
  geom_point(aes(x=endx, y =endy), color="red", size=I(5))

但是，当 n>10,000 时，它会变得非常慢，因此希望避免循环并使用某种形式的“应用”，但不知道如何添加第 n 行和第 n-1 行的坐标值。请帮忙，谢谢。

# second solution
d <- data.frame(dir=sample(c("w", "e", "n", "s"), n, replace=T), step=sample(1:4, n, replace=T))
xy <- data.frame(x=0, y=0)
x. <- data.frame(x=with(d, ifelse (dir=="w", -step, ifelse(dir=="e", step, 0))))
y. <- data.frame(y=with(d, ifelse (dir=="s", -step, ifelse(dir=="n", step, 0))))
x.y. <- cbind(x.,y.)
xy <- rbind(xy, x.y.)
head(xy)
# ... stuck here
4

4 回答 4

3

data.table对于这种问题很快......

walk.dt.f<-function(n=10000L, stepsize=1L:4L) {
  # lookup table with direction vector info
  dir.dt<-data.table(dir=c("w", "e", "n", "s"), xdir=c(-1L,1L,0L,0L), ydir=c(0L,0L,1L,-1L), key="dir")

  # initial position for random walk table
  walk.ini.dt<-data.table(rowid=0L,dir="n",step=0L)

  # generate table with random walk info
  walk.dt<-rbindlist(list(data.table(rowid=1L:n, dir=sample(dir.dt[,dir],n,replace=T), step=sample(stepsize,n,replace=T)), walk.ini.dt))

  # join the two tables, and multiply the step info by the direction vectors
  setkey(walk.dt,dir)
  walk.dt[dir.dt,c("xstep","ystep"):=list(step*xdir,step*ydir)]

  # update the key and reorder the rows
  setkey(walk.dt,rowid)

  # update the walk info table with the cumulative position
  walk.dt[,c("x","y"):=list(cumsum(xstep),cumsum(ystep))]

}

system.time(walk.dt.f(10000L))
## user  system elapsed 
## 0       0       0 

system.time(walk.dt.f(1e6L))
## user  system elapsed
## 0.25    0.00    0.25

编辑：将起始位置设置为 (0,0)

于 2013-02-14T19:51:15.030 回答
2

我想你越来越近了。如果您阅读已经发布的评论，您可以使其更快。所以我建议不要看这个：

n=10000
x.=sample(-4:4,n,rep=T)
y.=sample(-4:4,n,rep=T)
x=cumsum(x.)
y=cumsum(y.)

coord=data.frame(x,y)

然后准确地绘制你是怎么做的：

rw <- ggplot(coord, aes(x=x, y=y))
rw + geom_path() + 
  ggtitle(paste(n, "walks")) + 
  geom_point(aes(x=0, y =0), color="green", size=I(5)) +
  geom_point(aes(x=startx, y =starty), color="red", size=I(5))

更新：对于大于 10^5 的 n，绘图非常慢。也许基本图形会更快。

update2：这几乎与 joran 的响应一样慢/快。

于 2013-02-14T19:15:53.553 回答
2

呸！希望这将进一步推动我消除 R 中愚蠢的“for 循环本来就很慢”的谣言的目标，这里是您的第一个版本的重新工作，仍然使用快 40 倍以上的 for 循环。

我什至没有考虑过您实施随机游走是否有意义。我的意思只是指出如何以更快的速度实现原始代码的结果，同时仍然使用“慢”的 for 循环。

#My version
foo <- function(n){ 
    coord <- matrix(NA,nrow = n,ncol = 3) #origin
    coord[1,] <- c(0,0,0)
    dir <- sample(c("w", "e", "n", "s"), n,replace = TRUE) #random direction
    step <- sample(1:4, n,replace = TRUE) #how far to go in each walk
    for (i in 2:n){
      startx <- coord[i-1, 1]
      starty <- coord[i-1, 2]
      endx <- ifelse (dir[i]=="w", startx-step[i], ifelse(dir[i]=="e", startx+step[i], startx))
      endy <- ifelse (dir[i]=="n", starty+step[i], ifelse(dir[i]=="s", starty-step[i], starty))
      coord[i,] <- c(endx,endy,step[i])
    }
}

#Your version    
foo2 <- function(n){
    coord <- data.frame (x=0, y=0, step=0) #origin
    for (i in 1:n){
      dir <- sample(c("w", "e", "n", "s"), 1) #random direction
      step <- sample(1:4, 1) #how far to go in each walk
      startx <- coord[nrow(coord), 1]
      starty <- coord[nrow(coord), 2]
      endx <- ifelse (dir=="w", startx-step, ifelse(dir=="e", startx+step, startx))
      endy <- ifelse (dir=="n", starty+step, ifelse(dir=="s", starty-step, starty))
      newcoord <- data.frame (x=endx, y=endy, step=step)
      coord <- rbind(coord, newcoord)
    }
}


system.time(foo(10000))
   user  system elapsed 
  0.353   0.001   0.353 
> system.time(foo2(10000))
   user  system elapsed 
 11.374   2.061  13.308

我在这里所做的只是：

  1. 停止。使用。RBIND。并预先分配。
  2. 切换到矩阵。
  3. 将sample呼叫移出循环。
于 2013-02-14T19:17:30.217 回答
2

由于您正在尝试二维随机游走，因此有 4x4 可能的位移。你可以用 1 到 16 的数字对它们进行编码。但是，为了减少计算并将这些编码的数字映射到方向和位移量，我玩了一个小技巧，我没有用 1:16 编码步骤，而是用c(-7:0,4:11)

d <- sample(c(-7:0,4:11),n,replace=T)
delta <- d%%4+1
dir <- d%/%4
xd <- dir
xd[xd%%2 ==0]=0
yd <- dir
yd[xd%%2 ==1]=0
yd <- yd/2
x=c(0,xd*delta)
y=c(0,yd*delta)
x=cumsum(x)
y=cumsum(y)

coords<-data.frame(x,y)

这个版本只使用向量化操作，只有一点点开销。我认为它的性能接近于data.table之前给出的基于解决方案。

于 2013-02-14T20:11:53.097 回答