使用 for 循环保存模拟结果答案

【问题标题】：Saving results of a simulation using a for loop使用 for 循环保存模拟结果
【发布时间】：2021-12-24 09:19:19
【问题描述】：

我正在尝试运行更复杂的模拟并将到达复制的结果保存到数据框中。我正在尝试做的一个非常简单的抽象是：

sim.res <- data.frame(mn1 = mn1, mn2 = mn2, rep = rep)

for (k in 1:10){
  
  set.seed(k*8)
  mn1 <- mean(rnorm(25, 1, 1))
  mn2 <- mean(rnorm(25, 2, 1))

  sim.res$mn1 <- mn1
  sim.res$mn2 <- mn2
  sim.res$rep <- k
}

输出应该是一个包含三列的数据框，在本例中为 10 行，每个复制的结果存储在每一行中。我得到的是第 10 次复制。我哪里错了？

谢谢。

【问题讨论】：

标签： r loops simulation

【解决方案1】：

感谢您的回复。我还找到了一个解决方案，它给出了我所追求的结果：

sim.res <- list() # List container

for (k in 1:10){
  mn1 <- mean(rnorm(25, 1, 1))
  mn2 <- mean(rnorm(25, 2, 1))
  
  tmp <- cbind(mn1, mn2, k)
  
  sim.res[[k]] <- tmp
  
}

（我知道这不是最有效的）

最好，莱斯利

【讨论】：

【解决方案2】：

在 R 中，最基本的对象仍然是向量，所以这根本不需要循环。

当您设置sim.res$mn1 <- mn1 时，您将 mn1 放入 data.frame，但由于 mn1 的长度仅为 1 且 sim.res 有 10 行，R 重复该值 10 次并将其放入您的 data.frame。

同样 set.seed 只需使用一次即可重现性。

set.seed(42)
sim.res <- data.frame(mn1 = replicate(10, mean(rnorm(25, 1, 1))),
                      mn2 = replicate(10, mean(rnorm(25, 2, 1))), 
                      rep = 1:10)
sim.res
         mn1      mn2 rep
1  1.0818928 1.934759   1
2  1.0153544 2.185476   2
3  0.9898056 1.523066   3
4  0.9132660 2.229027   4
5  1.1003520 1.754576   5
6  0.7468452 2.211500   6
7  0.8480360 2.303661   7
8  1.2965247 2.014611   8
9  0.7511059 1.801408   9
10 1.0231177 1.924155  10

replicate 只是 sapply 的一个包装器，它使代码更易于理解。
也许从技术上讲，它使用了一个循环，但它肯定看起来不像一个。

【讨论】：

【解决方案3】：

每次遍历循环时，您都会覆盖每列的整个内容，因为默认情况下 R 是矢量化的（在列上运行）。此外，每次迭代循环时都无需重置种子。试试这样的。

library(tidyverse)

set.seed(123)

lapply(
  1:10,
  function(k) {
    data.frame(
      mn1 = mean(rnorm(25, 1, 1)), 
      mn2 = mean(rnorm(25, 2, 1)), 
      rep = k
    )
  }
) %>%
bind_rows()

         mn1      mn2 rep
1  0.9666697 2.102137   1
2  1.0102410 2.282576   2
3  0.7231910 1.769008   3
4  1.2152427 1.862371   4
5  1.1074981 1.875461   5
6  1.1727220 2.326180   6
7  0.9789770 2.025131   7
8  1.0981718 1.752828   8
9  0.9402713 1.928519   9
10 1.2714378 2.283175  10

或者，利用 n 的平均值从 N(mu, sd) 得出的事实是 N(mu, sd/sqrt(n))，

data.frame(
  rep=rep(1:10),
  mn1=rnorm(10, 1, sqrt(1/25)),
  mn2=rnorm(10, 2, sqrt(1/25))
)

   rep       mn1      mn2
1    1 1.2053570 2.176493
2    2 1.1502123 2.041120
3    3 0.6981667 1.876713
4    4 0.9809705 1.853040
5    5 0.8208104 1.973639
6    6 0.5858498 2.062003
7    7 1.0300240 1.792064
8    8 0.9841577 1.963138
9    9 0.9805261 2.193453
10  10 1.0432305 1.978344

【讨论】：