我怎样才能将这个模拟循环 1000 次并将其中一列组合成一个矩阵?

how can i loop this simulation1000 times and make one of the column together to be a matrix?

library(tidyverse)
library(broom)

# create a tibble with an id column for each simulation and x wrapped in list()
sim <- tibble(id = 1:1000,
               x = list(rbinom(1000,1,0.5))) %>% 
# to generate z, pr, y, k use map and map2 from the purrr package to loop over the list column x
# `~ ... ` is similar to `function(.x) {...}`
# `.x` represents the variable you are using map on
          mutate(z  = map(x, ~ log(1.3) * .x), 
                 pr = map(z, ~ 1 / (1 + exp(-.x))),
                 y  = map(pr, ~ rbinom(1000, 1, .x)),
                 k  = map2(x, y, ~ glm(.y ~ .x, family="binomial")),
# use broom::tidy to get the model summary in form of a tibble
                 sum = map(k, broom::tidy)) %>% 
# select id and sum and unnest the tibbles
  select(id, sum) %>% 
  unnest(cols = c(sum)) %>% 
simAll <- sim %>% 
  filter(term !="(Intercept)")

simAll 就像:

   id  term  estimate       std.error   statistic   p.value
1   1   .x  0.4058039189    0.1275272   3.182096892 1.462129e-03
2   2   .x  0.2515178701    0.1276719   1.970033693 4.883451e-02
3   3   .x  0.2464097082    0.1274321   1.933654251 5.315565e-02
4   4   .x  0.2308803864    0.1273598   1.812819663 6.985964e-02
5   5   .x  0.3029238760    0.1271623   2.382182816 1.721035e-02
6   6   .x  0.2452264719    0.1270829   1.929657417 5.364930e-02
7   7   .x  0.2390919312    0.1270123   1.882430831 5.977754e-02
8   8   .x  0.2437134055    0.1271373   1.916930426 5.524677e-02
9   9   .x  0.4372744410    0.1274612   3.430646232 6.021453e-04
10  10  .x  0.2915176118    0.1272609   2.290708545 2.198028e-02
11  11  .x  0.3373491310    0.1271283   2.653612132 7.963531e-03
12  12  .x  0.1991820874    0.1269380   1.569128570 1.166180e-01
13  13  .x  0.3437529981    0.1272502   2.701394595 6.904936e-03
14  14  .x  0.2229632179    0.1269851   1.755822253 7.911876e-02
15  15  .x  0.2606269011    0.1271385   2.049944533 4.036984e-02
Showing 1 to 15 of 1,000 entries, 6 total columns

这里的问题是,这里的估计列是x的值,我想像simALL一样有1000个类似的表(好像把整个模拟重复1000次),然后就有1000 * 1000 x,我想把它们做成矩阵(1000 * 1000),怎么办?

您可以为要重复的代码编写函数,并且 return 仅 estimate 列,因为这是您唯一需要的信息。

library(tidyverse)
run_fun <- function() {

sim <- tibble(id = 1:1000,
              x = list(rbinom(1000,1,0.5))) %>% 
          mutate(z  = map(x, ~ log(1.3) * .x), 
                 pr = map(z, ~ 1 / (1 + exp(-.x))),
                 y  = map(pr, ~ rbinom(1000, 1, .x)),
                 k  = map2(x, y, ~ glm(.y ~ .x, family="binomial")),
                 sum = map(k, broom::tidy)) %>%
           select(id, sum) %>% 
           unnest(cols = c(sum)) %>%
           filter(term !="(Intercept)") %>% 
           pull(estimate)
  return(sim)
}

您可以调用此函数 n 次:

data <- replicate(1000, run_fun())