2

我对 R 及其所有智慧相对较新,我正在努力提高我的脚本效率。我正在使用循环来模拟动物如何在不同地点之间移动。我遇到的问题是,当我增加站点数量或更改初始参数(基于移动或停留在同一站点的固定概率)时,我会以一个非常复杂的循环结束。如果我必须使用不同的参数运行多个不同的模拟,我更喜欢可以适应不同情况的更有效的循环或函数。第一个循环将根据初始概率填充矩阵,第二个循环将累积概率矩阵与值列表中的随机数(本例中为 10)进行比较,并决定该个体的命运(留下或离开到新网站)

这是我的代码的简化:

N<-4 # number of sites
sites<-LETTERS[seq(from=1,to=N)]

p.stay<-0.45
p.move<-0.4

move<-matrix(c(0),nrow=N,ncol=N,dimnames=list(c(sites),c(sites)))
from<-array(0,c(N,N),dimnames=list(c(sites),c(sites)))
to<-array(0,c(N,N),dimnames=list(c(sites),c(sites)))

# 以固定概率填充矩阵#

for(from in 1:N){
  for(to in 1:N){
     if(from==to){move[from,to]<-p.stay} else {move[from,to]<-p.move/(N-1)}
     }
 }

move
cumsum.move<-cumsum(data.frame(move))

steps<-100
result<-as.character("") # for storing results
rand<-sample(random,steps,replace=TRUE) 
time.step<-data.frame(rand)
colnames(time.step)<-c("time.step")
time.step$event<-""
to.r<-(rbind(sites))
j<-sample(1:N,1,replace=T) # first column to select (random number)
k<-sample(1:N,1,replace=T) # site selected after leaving and coming back

# 较长循环的开始#

for(i in 1:steps){
   if (time.step$time.step[i]<cumsum.move[1,j]){time.step$event[i]<-to.r[1]} else 
     if (time.step$time.step[i]<cumsum.move[2,j]){time.step$event[i]<-to.r[2]} else
        if (time.step$time.step[i]<cumsum.move[3,j]){time.step$event[i]<-to.r[3]} else
           if (time.step$time.step[i]<cumsum.move[4,j]){time.step$event[i]<-to.r[4]} else
              if (time.step$time.step[i]<(0.95)){time.step$event[i]<-NA} else
                 if (time.step$time.step[i]<1.0) break # break the loop             

 result[i]<-time.step$event[i]
 j<-which(to.r==result[i])
 if(length(j)==0){j<-k} # for individuals the leave and come back later

 }

 time.step

 result

这个循环是一个更大的循环的一部分,它将在一系列模拟之后模拟和存储结果。任何关于如何提高此循环效率的想法或评论,以便我可以轻松地修改站点数量或更改初始概率参数,而无需重复或不必对循环进行重大编辑,我们将不胜感激。

4

1 回答 1

0

我不确定我是否抓住了代码的精髓,但这比 for 循环快。一旦我们开始超过几千步,这就会开始具有优势。我将“随机”替换为均匀分布的样本 ( runif())

system.time(
time.step$event <- sapply(
  time.step$time.step, 
  function(x) rownames(
    cumsum.move[which(cumsum.move[,j] > x),])[[1]]
  )
)

这是我在 10,000 步时的结果。我正在使用笔记本电脑,因此使用 for 循环的 100,000 没有在 1 分钟内计算,但 sapply 在 14 秒内完成了。

> system.time(
+ time.step$event <- sapply(
+ time.step$time.step,
+ function(x) rownames(
+ cumsum.move[which(cumsum.move[,j] > x),])[[1]]
+ )
+ )
   user  system elapsed 
  1.384   0.000   1.387 
> head(time.step)
  time.step event
1 0.2787642     C
2 0.3098240     C
3 0.9079045     D
4 0.9904031     D
5 0.3754330     C
6 0.6984415     C
> system.time(
+ for(i in 1:steps){
+ if (time.step$time.step[i]<cumsum.move[1,j]){time.step$event[i]<-to.r[1]} else
+ if (time.step$time.step[i]<cumsum.move[2,j]){time.step$event[i]<-to.r[2]} else
+ if (time.step$time.step[i]<cumsum.move[3,j]){time.step$event[i]<-to.r[3]} else
+ if (time.step$time.step[i]<cumsum.move[4,j]){time.step$event[i]<-to.r[4]}
+ result[i]<-time.step$event[i]
+ }
+ )
   user  system elapsed 
  3.137   0.000   3.143 
> head(time.step)
  time.step event
1 0.2787642     C
2 0.3098240     C
3 0.9079045     D
4 0.9904031     D
5 0.3754330     C
6 0.6984415     C
于 2012-10-24T04:30:19.140 回答