Вычисление верхнего и нижнего доверительных интервалов по группам в dplyr - PullRequest
Новая
       73

Вычисление верхнего и нижнего доверительных интервалов по группам в dplyr

0 голосов
/ 22 января 2020

Я пытаюсь создать таблицу, которая показывает N (количество наблюдений), процентную частоту (ответов> 0), а также нижний и верхний доверительные интервалы для процентной частоты, и я хочу сгруппировать это по типу.

Пример данных 

dat <- data.frame(
  "type" = c("B","B","A","B","A","A","B","A","A","B","A","A","A","B","B","B"),
  "num" = c(3,0,0,9,6,0,4,1,1,5,6,1,3,0,0,0)
)

Ожидаемый результат (с заполненными значениями): 

Type   N   Percent   Lower 95% CI   Upper 95% CI
A
B

Попытка 

library(dplyr)
library(qwraps2)

table<-dat %>%
  group_by(type) %>%
  summarise(N=n(),
            mean.ci = mean_ci(dat$num),
            "Percent"=n_perc(num > 0))

Это сработало, чтобы получить N и процентную частоту, но вернуло ошибку: «Столбец должен иметь длину 1 (итоговое значение), а не 3», когда я добавил в mean_ci

Второй код, который я пробовал, нашел здесь : 

table2<-dat %>%
  group_by(type) %>%
  summarise(N.num=n(),
            mean.num = mean(dat$num),
            sd.num = sd(dat$num),
            "Percent"=n_perc(num > 0)) %>%
  mutate(se.num = sd.num / sqrt(N.num),
         lower.ci = 100*(mean.num - qt(1 - (0.05 / 2), N.num - 1) * se.num),
         upper.ci = 100*(mean.num + qt(1 - (0.05 / 2), N.num - 1) * se.num))

# A tibble: 2 x 8
#  type  N.num mean.num sd.num Percent        se.num lower.ci upper.ci
# <fct> <int>    <dbl>  <dbl> <chr>           <dbl>    <dbl>    <dbl>
#1 A         8     2.44   2.83 "6 (75.00\\%)"   1.00     7.35     480.
#2 B         8     2.44   2.83 "4 (50.00\\%)"   1.00     7.35     480.

Это дало мне вывод, но доверительные интервалы не логичны.

Ответы [ 2 ]

5 голосов
/ 22 января 2020

Вывод mean_ci представляет собой вектор длины 3. Это может быть неожиданно, потому что пакет добавил метод печати, так что, когда вы видите это в консоли, он выглядит как однозначное значение, а не цифра c длина> 1 вектор Но вы можете увидеть базовую структуру данных, посмотрев на str. 

mean_ci(dat$num) %>% str
 # 'qwraps2_mean_ci' Named num [1:3] 2.44 1.05 3.82
 # - attr(*, "names")= chr [1:3] "mean" "lcl" "ucl"
 # - attr(*, "alpha")= num 0.05

В итоге, каждый элемент каждого столбца выходных данных должен иметь длину 1, поэтому для суммирования необходим объект длины 3 положить в одну "ячейку" (элемент столбца) приводит к ошибке. Обходной путь - поместить вектор длины 3 в список, чтобы теперь он был списком длины 1. Затем вы можете использовать unnest_wider, чтобы разделить его на 3 столбца (и, следовательно, сделать таблицу «шире») 

library(tidyverse)

dat %>%
  group_by(type) %>%
  summarise( N=n(),
            mean.ci = list(mean_ci(num)),
            "Percent"= n_perc(num > 0)) %>% 
  unnest_wider(mean.ci)
# # A tibble: 2 x 6
#   type      N  mean   lcl   ucl Percent       
#   <fct> <int> <dbl> <dbl> <dbl> <chr>         
# 1 A         8  2.25 0.523  3.98 "6 (75.00\\%)"
# 2 B         8  2.62 0.344  4.91 "4 (50.00\\%)"
1 голос
/ 02 марта 2020

IceCreamToucan очень хороший ответ. Я публикую этот ответ, чтобы предложить другой способ представления информации. 

library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(qwraps2)

dat <- data.frame("type" = c("B","B","A","B","A","A","B","A","A","B","A","A","A","B","B","B"),
                  "num"  = c(3,0,0,9,6,0,4,1,1,5,6,1,3,0,0,0))

При построении вызова dplyr::summarize вы можете использовать вызов qwraps2::frmtct для форматирования вывода qwraps2::mean_ci в символьная строка длиной один.

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

dat %>%
  dplyr::group_by(type) %>%
  dplyr::summarize(N = n(),
                   mean.ci = qwraps2::frmtci(qwraps2::mean_ci(.data$num)),
                   Percent = qwraps2::n_perc(.data$num > 0))
#> # A tibble: 2 x 4
#>   type      N mean.ci           Percent       
#>   <fct> <int> <chr>             <chr>         
#> 1 A         8 2.25 (0.52, 3.98) "6 (75.00\\%)"
#> 2 B         8 2.62 (0.34, 4.91) "4 (50.00\\%)"

Создано в 2020- 03-02 представьте пакет (v0.3.0)

...