7

假设我想使用以下命令进行模拟function

fn1 <- function(N) {
  res <- c()
  for (i in 1:N) {
    x <- rnorm(2)
    res <- c(res, x[2]-x[1])
  }
  res
}

对于非常大N的 ,计算似乎挂起。有没有更好的方法来做到这一点?

(灵感来自:https ://stat.ethz.ch/pipermail/r-help/2008-February/155591.html )

4

5 回答 5

9

在 R 中,通过使用 apply 函数可以极大地提高循环的效率,这些函数本质上是一次处理整个数据向量,而不是遍历它们。对于上面显示的循环,每次迭代期间都会发生两个基本操作:

# A vector of two random numbers is generated
x <- rnorm( 2 )

# The difference between those numbers is calculated
x[2] - x[1]

在这种情况下,适当的功能将是sapply(). sapply()对对象列表进行操作,例如循环语句生成的向量1:N并返回结果向量:

sapply( 1:N, function( i ){ x <- rnorm(2); return( x[2] - x[1] ) } )

请注意,索引值i在函数调用期间可用,并连续采用 和 之间的值1N但在这种情况下不需要它。

养成识别哪里apply可以使用的习惯for是一项非常有价值的技能——许多用于并行计算的 R 库通过函数提供即插即用的并行apply化。使用apply通常可以在重构代码的情况下在多核系统上获得显着的性能提升。

于 2009-07-23T04:27:32.077 回答
4

扩展我对 chris_dubois 回答的评论,这里有一些时间信息:

> system.time(res <- rnorm(50000) - rnorm(50000))
user  system elapsed
0.06    0.00    0.06

将此与来自同一答案的 fn3 进行对比:

> system.time(res3 <- fn3(50000))
user  system elapsed
1.33    0.01    1.36

首先要注意的是我的笔记本电脑比 chris_dubois 的机器慢。:)

第二点,也是更重要的一点是,这里非常适用的矢量方法要快一个数量级。(Richie Cotton 在对相同答案的评论中也指出)。

这让我想到了最后一点:这是一个神话apply它的朋友比 R 中的循环快得多for。在我见过的大多数测量中,它们的顺序相同。因为它们只是for幕后的循环。另见这篇文章:

http://yusung.blogspot.com/2008/04/speed-issue-in-r-computing-apply-vs.html

根据 Brian Ripley 教授的说法,“apply() 只是一个循环的包装器。” 使用 apply() 的唯一好处是它可以让你的代码更整洁!

确切地。apply如果它更具表现力,则应使用它,尤其是在您以函数式编程时。不是因为它更快。

于 2009-07-26T04:34:14.843 回答
2

R 中的 for 循环是出了名的慢,但这里还有另一个问题。预先分配结果向量 res 比在每次迭代时附加到 res 上要快得多。

下面我们可以将上述版本的速度与仅以长度为 N 的向量 res 开始并在循环期间更改第 i 个元素的版本进行比较。

fn1 <- function(N) {
  res <- c()
  for (i in 1:N) {
     x <- rnorm(2)
     res <- c(res,x[2]-x[1])
  }
  res
}
fn2 <- function(N) {
  res <- rep(0,N)
  for (i in 1:N) {
     x <- rnorm(2)
     res[i] <- x[2]-x[1]
  }
  res
}
> N <- 50000
> system.time(res1 <- fn1(N))
   user  system elapsed 
  6.568   0.256   6.826 
> system.time(res2 <- fn2(N))
   user  system elapsed 
  0.452   0.004   0.496 

此外,正如Sharpie 指出的那样,我们可以通过使用 R 函数apply(或其亲属sapplylapply)来稍微加快速度。

fn3 <- function(N) {
  sapply( 1:N, function( i ){ x <- rnorm(2); return( x[2] - x[1] ) } )
}
> system.time(res3 <- fn3(N))
   user  system elapsed 
  0.397   0.004   0.397 
于 2009-07-23T04:19:53.713 回答
2

有时不需要循环。由于 rnorm 提供 iid 样本(理论上),您将 X-Y 通过执行以下操作获得相同的结果(X 和 Y 为 N(0,1) 的采样):

res <- rnorm(N)-rnorm(N)
于 2009-07-24T07:31:41.667 回答
0

也许对您的功能最有效的替代品就是:

fn <- function(n) rnorm(N,0,sqrt(2))

这比 iid 正态变量的差异快两倍。更一般地说,如果您的目标是运行简单的模拟,则向量/数组预分配和对本机函数的调用会大大加快该过程。

如果你想为统计估计执行蒙特卡罗模拟(例如,MCMC),R 有许多本地包。对于一般随机模拟,我不知道 R 包,但您可能想尝试 Simpy ( http://simpy.sourceforge.net/ ),这非常好。

于 2009-07-27T15:00:44.440 回答