如何优化 R 中重复测量模拟的内存和速度？

问题描述 投票：0回答：1

我正在尝试进行一个模拟，其中有 257 个分区以及每个分区内的村庄列表，每个分区都有已知的患有某种疾病的人口百分比。我想迭代五个村庄（{5, 10, 20 ,30 ,50} 中的

n_village
）和每个村庄五个人（{10, 15, 20, 25, 30} 中的
n_person
）对于每个我想重复调查 1000 次。这意味着 5 * 5 * 257 * 1000 = 6,425,000 个结果。此外，我的网站列表有 4,999 个，其中一个国家/地区有 3836 个，另一个国家/地区有 1163 个。

所有分区：

  1. 随机选择
    n_village
    村庄。
  2. 在每个村庄中，随机选择 
    n_person
    个人并计算有多少人患有这种情况。
  3. 计算估计的地区级患有该疾病的人的百分比，以及该百分比的置信区间。

我的第一次尝试涉及使用 mapply 迭代所有内容并将结果随机采样的村庄保存到列表中，但这花费了非常非常长的时间。

我使用嵌套 

tibbles
和 
pmap
取得了很多进展，但我现在遇到了内存问题和性能下降。

下面是我的方法的 10 个分区的近似值。我还将示例限制为 10 次模拟。

library(tidyverse)
library(tictoc)

admin_unit_df <- 
  tibble(
    shapeName_0 = "Country 1", 
    shapeName_1 = "Region 1",
    shapeName_2 = "District 1",
    shapeName_3 = paste0("Sub-district ", 1:10)
  )

sampling_frame_df <- 
  tibble(
    site_id = 1:500,
    shapeName_3 = sample(admin_unit_df$shapeName_3, 500, replace = T),
    site.prev = runif(500)
  ) %>%
  group_by(
    shapeName_3
  ) %>%
  nest()

sim_parms_df <-
  expand.grid(
    "n_village" = c(2, 5, 10, 15, 20),
    "n_person"  = c(10, 20, 30, 40, 50),
    "sim_id"    = 1:10
  ) %>%
  arrange(n_village, n_person, sim_id) %>%
  as_tibble() %>%
  rowwise() 

simulations_df <-
  admin_unit_df %>%
  mutate(
    data = sampling_frame_df$data,
    
    sim_parms = list(sim_parms_df)
  )

simulations_df

simulations_df$sim_parms[[1]]

然后我从数据中随机抽取村庄

tic()
simulations_df %<>%
  mutate(
    sim_parms = 
      pmap(
        .l = 
          list(
            "sim_parms" = sim_parms,
            "data" = data
          ),
        .f = 
          function(sim_parms, data){
            n_site_avail = nrow(data)
            
            sim_parms %>%
              mutate(
                samples =
                  list(
                    data[sample(x = 1:n_site_avail, size = min(n_village, n_site_avail), replace = F), ] %>%
                      as_tibble() %>%
                      rowwise()
                  )
              )
          }
    )
  )
toc()

此后，我使用 pmap 迭代并随机抽样每个村庄的病情呈阳性的人数。

tic()
simulations_df %<>%
  mutate(
    sim_parms = 
      pmap(
        .l = 
          list(
            "sim_parms" = sim_parms
          ),
        .f = 
          function(sim_parms){
            # cat(names(sim_parms), "\n")
            
            sim_parms %>%
              ungroup() %>%
              mutate(
                samples = 
                  pmap(
                    .l = list("samples" = samples, "n_person" = n_person),
                    .f = function(samples, n_person){
                      # cat(names(samples), "\n")
                      samples %>%
                        mutate(
                          n_pos = rbinom(n = 1, size = n_person, p = site.prev),
                          p_hat = n_pos / n_person
                        )
                    }
                  )
              )
            
          }
    )
  )
toc()

这里的事情开始变慢并爆炸，但我已经成功完成了。

我想做的下一件事是使用 

binom.test()
计算估计患病率的置信区间下限和上限。

这是一个非常基础的分析，但需要相当长的时间。在某个时候，我的目标是使用 lme4 对所有 1000 次模拟进行多级模型分析。

关于如何使这项任务易于处理，有什么建议吗？有没有办法循环并保存结果而不保存所有中间值？

顺便说一句——保存到 RDS 并读取 RDS 需要非常非常长的时间。有什么建议可以加快速度吗？

r performance loops memory multidimensional-array
1个回答
0
投票

通过操作大 

data.frame
而不是长 
list
/具有小 
data.frame
条目的嵌套列，您可以获得大量性能。分组的 
data.frame
通常操作起来更快。

另外，对于不需要分组的操作，尽量避免携带分组数据（

rowwise()
也是分组）。 
.by/by
动词的新 
dplyr
参数非常棒，因为它允许您仅针对当前操作应用分组。

这是对您的代码的快速重写。这个想法是限制 

map()
函数的使用，并在使用时将结果绑定到 
data.frame
，然后再应用下一个操作。

（我只是简单地检查了我的代码的结果，它似乎产生了与您相同的结果。请仔细检查以确保它符合您的规范。）

library(tidyverse)
library(tictoc)

tic()

set.seed(1)
admin_unit_df <-
  tibble(
    shapeName_0 = "Country 1",
    shapeName_1 = "Region 1",
    shapeName_2 = "District 1",
    shapeName_3 = paste0("Sub-district ", 1:10)
  )

sampling_frame_df <-
  tibble(
    site_id = 1:500,
    shapeName_3 = sample(admin_unit_df$shapeName_3, 500, replace = T),
    site.prev = runif(500)
  )

sim_parms_df <-
  expand.grid(
    "n_village" = c(2, 5, 10, 15, 20),
    "n_person"  = c(10, 20, 30, 40, 50),
    "sim_id"    = 1:10
  ) |>
  arrange(n_village, n_person, sim_id) |>
  as_tibble()

samples_df <-
  sim_parms_df |>
  mutate(samples =
           pmap(list(n_village,
                     n_person),
                \(villages, persons) {
                  slice_sample(sampling_frame_df,
                               n = min(villages, persons),
                               by = shapeName_3)
                })) |>
  unnest(samples)

res <-
  samples_df |>
  group_by(shapeName_3, site_id) |>
  mutate(n_pos = rbinom(n = 1, size = n_person, p = site.prev),
         p_hat = n_pos / n_person)

toc()
#> 0.364 sec elapsed

res
#> # A tibble: 24,500 × 8
#> # Groups:   shapeName_3, site_id [500]
#>    n_village n_person sim_id site_id shapeName_3    site.prev n_pos p_hat
#>        <dbl>    <dbl>  <int>   <int> <chr>              <dbl> <int> <dbl>
#>  1         2       10      1     254 Sub-district 9   0.00573     0   0  
#>  2         2       10      1     155 Sub-district 9   0.971      10   1  
#>  3         2       10      1      70 Sub-district 4   0.0764      0   0  
#>  4         2       10      1     451 Sub-district 4   0.0389      0   0  
#>  5         2       10      1     278 Sub-district 7   0.188       3   0.3
#>  6         2       10      1     157 Sub-district 7   0.0392      0   0  
#>  7         2       10      1     258 Sub-district 1   0.167       4   0.4
#>  8         2       10      1     259 Sub-district 1   0.764       9   0.9
#>  9         2       10      1       5 Sub-district 2   0.820       9   0.9
#> 10         2       10      1     496 Sub-district 2   0.305       4   0.4
#> # ℹ 24,490 more rows
最新问题
© www.soinside.com 2019 - 2024. All rights reserved.