получить средние значения и несколько доверительных интервалов нескольких моделей в R - PullRequest
0 голосов
/ 27 февраля 2019

Я хотел бы получить доверительные интервалы 95% и 90%, а также средние значения для нескольких моделей в R.

Данные



data <- data.frame(occup = c(2,3,5,4,2,2,6,1,2,0),
                   unoccup = c(1,2,0,3,0,4,1,1,2,2),
                   month = c("feb", "feb", "feb", "feb", "feb", "mar", "mar", "mar", "mar", "mar"))

У меня есть функция

binomNLL_ratio = function(p, k, N) {
  -sum(dbinom(k, prob = p, size = N, log=TRUE))
}

Необходимые библиотеки

library(purrr)
library(bbmle)

Запустить скрипт

data %>%
   split(.$month) %>% map(~mle2(minuslog = binomNLL_ratio, start = list(p = 0.5), data = list(N = .$occup + .$unoccup, k = .$occup))) %>%
  map(confint, level = 0.95)

Это дает мне 95% доверительные интервалы каждого месяца.Я также могу заменить 0.95 на 0.9, чтобы получить 90% CI, или заменить map(confint) на map(coef), чтобы получить среднее значение для каждой месячной модели.

Тем не менее, в идеале я хотел бы получить 95% CI, 90% CI и средства каждой модели в одном и том же фрейме данных.Как я могу передать несколько функций и параметров для получения желаемых результатов?

Спасибо за вашу помощь.

1 Ответ

0 голосов
/ 27 февраля 2019

Это был хитрый маленький жучок, и я потратил слишком много времени, пытаясь заставить invoke_map() работать.В конце концов, я согласился на немного хакерский способ ведения дел, но он возвращает тибл со значениями confint() и coef():

library(tidyverse)
library(bbmle)
#> Loading required package: stats4
#> 
#> Attaching package: 'bbmle'
#> The following object is masked from 'package:dplyr':
#> 
#>     slice

data <- data.frame(occup = c(2,3,5,4,2,2,6,1,2,0),
                   unoccup = c(1,2,0,3,0,4,1,1,2,2),
                   month = c("feb", "feb", "feb", "feb", "feb", "mar", "mar", "mar", "mar", "mar"))

binomNLL_ratio = function(p, k, N) {
  -sum(dbinom(k, prob = p, size = N, log=TRUE))
}

data %>%
  split(.$month) %>% 
  map(~mle2(minuslog = binomNLL_ratio, start = list(p = 0.5), 
            data = list(N = .$occup + .$unoccup, k = .$occup))) %>% 
  imap_dfr(~ {tibble(confints = confint(.x, level = .95), 
                coefs = coef(.x),
                month = .y)})
#> Warning in dbinom(k, prob = p, size = N, log = TRUE): NaNs produced

#> Warning in dbinom(k, prob = p, size = N, log = TRUE): NaNs produced

#> Warning in dbinom(k, prob = p, size = N, log = TRUE): NaNs produced

#> Warning in dbinom(k, prob = p, size = N, log = TRUE): NaNs produced
#> # A tibble: 4 x 3
#>   confints coefs month
#>      <dbl> <dbl> <chr>
#> 1    0.523 0.727 feb  
#> 2    0.881 0.727 feb  
#> 3    0.317 0.524 mar  
#> 4    0.725 0.524 mar

Создано в 2019-02-26 представительный пакет (v0.2.1)

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