如何有效地生成模拟值的数据框？答案

【问题标题】：How can I efficiently generate a dataframe of simulated values?如何有效地生成模拟值的数据框？
【发布时间】：2014-04-23 22:54:27
【问题描述】：

我正在尝试根据现有分布参数生成模拟值的数据框。我的主要数据框包含每个观察值的平均值和标准差，如下所示：

example.data <- data.frame(country=c("a", "b", "c"), 
                           score_mean=c(0.5, 0.4, 0.6), 
                           score_sd=c(0.1, 0.1, 0.2))

#   country score_mean score_sd
# 1       a        0.5      0.1
# 2       b        0.4      0.1
# 3       c        0.6      0.2

我可以使用sapply() 和自定义函数来使用 score_mean 和 score_sd 参数从正态分布中随机抽取：

score.simulate <- function(score.mean, score.sd) {
  return(mean(rnorm(100, mean=score.mean, sd=score.sd)))
}

simulated.scores <- sapply(example.data$score_mean, 
                       FUN=score.simulate, 
                       score.sd=example.data$score_sd)

# [1] 0.4936432 0.3753853 0.6267956

这将生成一轮（或一列）模拟值。但是，我想生成很多列（如 100 或 1,000）。我发现这样做的唯一方法是将我的 sapply() 函数包装在 lapply() 内的通用函数中，然后将结果列表转换为带有 ldply() 在 plyr 中的数据框：

results.list <- lapply(1:5, FUN=function(x) sapply(example.data$score_mean, FUN=score.simulate, score.sd=example.data$score_sd))

library(plyr)
simulated.scores <- as.data.frame(t(ldply(results.list)))

#           V1        V2        V3        V4        V5
# V1 0.5047807 0.4902808 0.4857900 0.5008957 0.4993375
# V2 0.3996402 0.4128029 0.3875678 0.4044486 0.3982045
# V3 0.6017469 0.6055446 0.6058766 0.5894703 0.5960403

这可行，但是 (1) 看起来确实很复杂，尤其是使用 as.data.frame(t(ldply(lapply(... FUN=function(x) sapply ...)))) 方法时，(2) 使用大量迭代或更大的数据时它真的很慢——我的实际数据集有 3,000 行，运行 1,000迭代需要 1-2 分钟。

有没有更有效的方法来创建这样的模拟值数据框？

【问题讨论】：

您可能不想从现有数据框中获取值，而是从向量中获取它们。然后创建新的向量。之后，将它们放回数据框中。这样，R 就不必花时间搜索数据框来获取值。
所以手动提取感兴趣的列，比如scores.mean <- example.data$score_mean?
你的意思是lapply( 1:nrow(df) , function(i) rnorm( 3 , df[i,'mean'] , sd=df[i,'sd'] ) )？（df 是 example.data）。要随机使用平均值和标准差，您可以将lapply 的第一个参数中的1:nrow(df) 替换为sample( 1:nrow(df) , n , TRUE ) 以使用随机行（带替换）。
附注在您的示例中，将results.list 变成data.frame 和as.data.frame( results.list ) 就足够了！！（data.frame 本质上是一个基于列表的结构本身）。
我使用Rcpp 模拟了大型n，如果你的函数足够简单，那么易于编写和加速是很重要的。

标签： r simulation

【解决方案1】：

我能想到的最快方法是利用 rnorm 内置的矢量化。 mean 和 sd 参数都是矢量化的，但是您只能为抽奖次数提供一个整数。如果您为mean 和sd 参数提供向量，R 将循环遍历它们，直到完成所需的绘制次数。因此，只需将 n 到 rnorm 的参数设为 mean 向量长度的倍数。乘数将是 data.frame 每一行的重复次数。在下面的函数中，这是n。

我想不出比单独使用base::rnorm 更重要的方法。

工作示例

#example data
df <- data.frame(country=c("a", "b", "c"), 
                           mean=c(1, 10, 100), 
                           sd=c(1, 2, 10))

#function which returns a matrix, and takes column vectors as arguments for mean and sd
normv <- function( n , mean , sd ){
    out <- rnorm( n*length(mean) , mean = mean , sd = sd )
    return( matrix( out , , ncol = n , byrow = FALSE ) )
    }

#reproducible result (note order of magnitude of rows and input sample data)
set.seed(1)
normv( 5 , df$mean , df$sd )
#           [,1]      [,2]       [,3]        [,4]        [,5]
#[1,]  0.3735462  2.595281   1.487429   0.6946116   0.3787594
#[2,] 10.3672866 10.659016  11.476649  13.0235623   5.5706002
#[3,] 91.6437139 91.795316 105.757814 103.8984324 111.2493092

【讨论】：

哦，太棒了！我不知道rnorm() 中的n 参数可以采用向量。太棒了。谢谢！

【解决方案2】：

如果您记得 rnorm(1, mean, sd) 与 rnorm(1)*sd + mean 相同，那么这可以很快完成，因此使用您的数据框 df，您可以生成 obs 观察结果的 sim 模拟，例如：

obs = nrow(df)
sim = 1000
mat = data.frame(matrix(rnorm(obs*sim), obs, sim) * df$sd + df$mean)

您可以使用 rowMeans(mat) 检查这是否具有所需的均值，并检查例如第 1 行的标准偏差为 sd(mat[1,])。

【讨论】：