Как сохранить множество brms в виде отдельных объектов R или файла Rds/Rda?

Я пытаюсь подогнать байесовскую модель, используя brms::brm(), заменяя априор одного параметра в модели один за другим на purrr::map2() (у меня есть 63 априора для этого параметра). Я могу «теоретически» сохранить каждую подогнанную модель как отдельный объект в глобальной среде (т. моего предыдущего вопроса. К list2env() объекты (brmsfits) будут сохранены после выполнения всех 63 итераций. Однако, когда я запускаю весь код, я всегда получаю сообщение об ошибке, в котором говорится, что Error in scan(con, nlines = 1, sep = ",", quiet = TRUE) : could not allocate memory (0 Mb) in C function 'R_AllocStringBuffer' и объект brmsfit не хранится в глобальной среде, хотя подгонка модели, кажется, выполняется 63 раза, столько раз, сколько существует мой предыдущий.

Поэтому я хотел бы сохранить каждый объект brmsfit как объект R и как файл .Rds или .Rda сразу после его итерации, чтобы избежать такой проблемы с памятью. Но что мне сделать, чтобы это осознать?

MWE

Обратите внимание, что следующая команда — это всего лишь схематический пример с общедоступными данными того, что я пытаюсь сделать. Это будет работать без проблемы, о которой я упоминал выше, поскольку данные и модель в brm() в этом примере намного проще, а количество априорных значений намного меньше, чем то, что я решаю.

library(lme4)
library(tidyverse)
library(magrittr)
library(brms)
library(cmdstanr)
library(rstan)

## Parallelize the chains using all the cores:
options(mc.cores = parallel::detectCores())

prior_test <- data.frame(
  sd = c(0.01, 0.01, 0.01),
  mean = c(50, -50, 0)
) %>% 
  mutate(
    id = row_number()
  )

list2env(
  purrr::map2(
    prior_test$sd,
    prior_test$mean,
    function(psd, pm){
          gc()
          gc()
          cbfm1 <- brm(
            Reaction ~ 0 + Intercept + Days + (0 + Intercept + Days|Subject), 
            data = sleepstudy,
            family = "normal",
            prior =
              c(
                prior(normal(0, 1), class = b, coef = Intercept),
                set_prior(
                  paste0(
                    "normal(",
                    pm,
                    ", ",
                    psd,
                    ")"
                  ),
                  class = "b",
                  coef = "Days"
                ),
                prior(normal(0, 1), class = sd),
                prior(lkj(2),       class = cor)    
              ),
            save_pars = save_pars(all = TRUE),
            backend = "cmdstanr"
          )
        }
  ) %>%
    setNames(paste0('model', prior_test$id)),
  .GlobalEnv
)

Сообщение об ошибке с обратной трассировкой

Error in scan(con, nlines = 1, sep = ",", quiet = TRUE) : 
  could not allocate memory (0 Mb) in C function 'R_AllocStringBuffer' 
15.
scan(con, nlines = 1, sep = ",", quiet = TRUE) 
14.
rstan::read_stan_csv(out$output_files()) 
13.
.fit_model(model, ...) 
12.
.fun(model = .x1, sdata = .x2, algorithm = .x3, backend = .x4, 
    iter = .x5, warmup = .x6, thin = .x7, chains = .x8, cores = .x9, 
    threads = .x10, inits = .x11, exclude = .x12, control = .x13, 
    future = .x14, seed = .x15, silent = .x16) 
11.
eval(expr, envir, ...) 
10.
eval(expr, envir, ...) 
9.
eval2(call, envir = args, enclos = parent.frame()) 
8.
do_call(fit_model, fit_args) 
7.
brm(voice ~ 0 + Intercept + agent + patient + REGION1 + REGION2 + 
    agent:patient + agent:REGION1 + agent:REGION2 + patient:REGION1 + 
    patient:REGION2 + agent:patient:REGION1 + agent:patient:REGION2 + 
    (0 + Intercept + agent + patient + agent:patient | subj) +  ... 
6.
.f(.x[[i]], .y[[i]], ...) 
5.
map2(R1data$prior_sd, R1data$prior_mean, function(psd, pm) {
    gc()
    gc()
    brm(voice ~ 0 + Intercept + agent + patient + REGION1 + REGION2 +  ... 
4.
eval(lhs, parent, parent) 
3.
eval(lhs, parent, parent) 
2.
map2(R1data$prior_sd, R1data$prior_mean, function(psd, pm) {
    gc()
    gc()
    brm(voice ~ 0 + Intercept + agent + patient + REGION1 + REGION2 +  ... 
1.
list2env(map2(R1data$prior_sd, R1data$prior_mean, function(psd, 
    pm) {
    gc()
    gc() ... 

person Carlos Luis Rivera    schedule 29.10.2020    source источник


Ответы (2)


Я думаю, вы могли бы использовать что-то вроде следующего (это более простая версия, которая избегает того, откуда поступают ваши данные, имена переменных, спецификация модели и все такое)


dir.create("models")

# This function generates a file path, runs the `brm` function and saves the object
# returned by `brm()` as a .rds file
my_fun = function(prior_sd, prior_mean) {
  file_name = paste0("models/model_sd_", prior_sd, "_mean_", prior_mean, ".rds")
  brms_object = brm(...) # pass formula, prior, data, options, etc.
  saveRDS(brms_object, file_name) # you can control compression too if you want
}

sds = c(0.01, 0.01, 0.01)
means = c(50, -50, 0)
purrr::walk2(sds, means, my_fun)

Вы должны убедиться, что фрейм данных, который вы собираетесь использовать, существует в глобальной среде.

person Tomas Capretto    schedule 29.10.2020
comment
Спасибо за отличный ответ! Я принял и проголосовал за ваш, но также проверьте мой ответ, где я сохраняю результаты вычислений как объекты R в глобальной среде и избегайте использования - в именах объектов R. - person Carlos Luis Rivera; 30.10.2020

Я немного изменил функцию в принятом ответе, чтобы (1) мы могли сохранить объекты brmsfit в глобальной среде, используя assign(..., envir = .GlobalEnv) и (2) мы можем избежать имени объекта с помощью -, которое R интерпретирует как оператор вычитания, когда мы используем отрицательные значения в качестве априорного среднего.

data <- data.frame(
  prior_sd = prior_sd <- c(0.01, 0.01, 0.01),
  prior_mean = prior_mean <- c(50, -50, 0),
  id = case_when(
    prior_mean == 0 ~ paste0("model_mean_", prior_mean, "_sd_", prior_sd),
    prior_mean < 0 ~ paste0("model_mean_m", abs(prior_mean), "_sd_", prior_sd),
    prior_mean > 0 ~ paste0("model_mean_p", prior_mean, "_sd_", prior_sd)
    )
)

computation = function(prior_mean, prior_sd) {
  file_name = case_when(
    prior_mean == 0 ~ paste0("path-to-your-directory/model_mean_", prior_mean, "_sd_", prior_sd, ".rds"),
    prior_mean < 0 ~ paste0("path-to-your-directory/model_mean_m", abs(prior_mean), "_sd_", prior_sd, ".rds"),
    prior_mean > 0 ~ paste0("path-to-your-directory/model_mean_p", prior_mean, "_sd_", prior_sd, ".rds")
  )
  # pass formula, prior, data, options, etc.
  brms_object = brm(
    Reaction ~ 0 + Intercept + Days + (0 + Intercept + Days|Subject), 
    data = sleepstudy,
    family = "normal",
    prior =
      c(
        prior(normal(0, 1), class = b, coef = Intercept),
        set_prior(
          paste0(
            "normal(",
            prior_mean,
            ", ",
            prior_sd,
            ")"
          ),
          class = "b",
          coef = "Days"
        ),
        prior(normal(0, 1), class = sd),
        prior(lkj(2),       class = cor)    
      ),
    save_pars = save_pars(all = TRUE),
    backend = "cmdstanr"
  )
  # save the computation results with object name
  assign(
    case_when(
      prior_mean == 0 ~ paste0("model_mean_",  prior_mean, "_sd_", prior_sd),
      prior_mean <  0 ~ paste0("model_mean_m", abs(prior_mean), "_sd_", prior_sd),
      prior_mean >  0 ~ paste0("model_mean_p", prior_mean, "_sd_", prior_sd)
    ), 
    brms_object,
    envir = .GlobalEnv)
  # you can control compression too if you want
  saveRDS(brms_object, file_name) 
}

purrr::walk2(data$prior_mean, data$prior_sd, computation)
person Carlos Luis Rivera    schedule 30.10.2020