Мой канал - https://t.me/tidy_mng
Данная заметка - это любительский перевод статьи Julia Silge.
Эта статья взята из блога Julia Silge, которая демонстрирует использование пакетов tidymodels. В сегодняшней заметке будет показана относительно новая функция из пакета rsample - reg_intervals. Данная функция разработана для быстрого поиска доверительных интервалов bootstrap.
Данные: набор #TidyTuesday о рекламных роликах суперкубка.
Знакомимся с данными
Наша цель выяснить, как с течением времени менялись характеристики рекламных роликов суперкубка. В наших данных не так много наблюдений, что позволяет нам использовать метод bootstrap для демонстрации.
Прочитаем наши данные:
library(tidyverse) youtube <- read_csv("https://raw.githubusercontent.com/rfordatascience/tidytuesday/master/data/2021/2021-03-02/youtube.csv")
Давайте сделаем некоторые преобразование и посмотрим, как ведут себя данные с течением времени.
youtube %>% select(year, funny:use_sex) %>% pivot_longer(funny:use_sex) %>% group_by(year, name) %>% summarise(prop = mean(value)) %>% ungroup() %>% ggplot(aes(year, prop, color = name)) + geom_line(size = 1.2, show.legend = FALSE) + facet_wrap(~name) + scale_y_continuous(labels = scales::percent) + labs(x = NULL, y = "% of commercials")
Сделаем простую модель
Не смотря на то, что данные не выглядят идеально линейными, мы можем использовать линейную модель для оценки наших характеристик.
simple_mod <- lm(year ~ funny + show_product_quickly + patriotic + celebrity + danger + animals + use_sex, data = youtube) summary(simple_mod)
Call: lm(formula = year ~ funny + show_product_quickly + patriotic + celebrity + danger + animals + use_sex, data = youtube) Residuals: Min 1Q Median 3Q Max -12.5254 -4.1023 0.1456 3.9662 10.1727 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 2011.0838 0.9312 2159.748 < 2e-16 *** funnyTRUE -2.8979 0.8593 -3.372 0.00087 *** show_product_quicklyTRUE 0.7706 0.7443 1.035 0.30160 patrioticTRUE 2.0455 1.0140 2.017 0.04480 * celebrityTRUE 2.4416 0.7767 3.144 0.00188 ** dangerTRUE 0.4814 0.7846 0.614 0.54007 animalsTRUE 0.1082 0.7330 0.148 0.88274 use_sexTRUE -2.4041 0.8175 -2.941 0.00359 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 5.391 on 239 degrees of freedom Multiple R-squared: 0.178, Adjusted R-squared: 0.1539 F-statistic: 7.393 on 7 and 239 DF, p-value: 4.824e-08
Мы получили статистические свойства линейной модели, но мы можем использовать bootstrap для того, чтобы оценить дисперсию наших величин. В этом нам поможет пакет rsample.
library(rsample) bootstraps(youtube, times = 1000)
# Bootstrap sampling # A tibble: 1,000 × 2 splits id <list> <chr> 1 <split [247/87]> Bootstrap0001 2 <split [247/85]> Bootstrap0002 3 <split [247/88]> Bootstrap0003 4 <split [247/83]> Bootstrap0004 5 <split [247/87]> Bootstrap0005 6 <split [247/84]> Bootstrap0006 7 <split [247/92]> Bootstrap0007 8 <split [247/84]> Bootstrap0008 9 <split [247/91]> Bootstrap0009 10 <split [247/86]> Bootstrap0010 # … with 990 more rows
Нам бы пришлось выполнить достаточно много предварительных шагов, чтобы найти доверительные интервалы. Поэтому в пакете rsample предусмотрена функция reg_intervals(), которая находит доверительные интервалы для моделей типа: lm(), glm().
set.seed(123) youtube_intervals <- reg_intervals(year ~ funny + show_product_quickly + patriotic + celebrity + danger + animals + use_sex, data = youtube, type = "percentile", keep_reps = TRUE) youtube_intervals
# A tibble: 7 × 7 term .lower .estimate .upper .alpha .method .replicates <chr> <dbl> <dbl> <dbl> <dbl> <chr> <list<tibble> 1 animalsTRUE -1.22 0.144 1.51 0.05 percentile [2,001 × 2] 2 celebrityTRUE 0.828 2.46 4.06 0.05 percentile [2,001 × 2] 3 dangerTRUE -1.01 0.515 2.09 0.05 percentile [2,001 × 2] 4 funnyTRUE -4.58 -2.91 -1.26 0.05 percentile [2,001 × 2] 5 patrioticTRUE 0.112 2.05 3.88 0.05 percentile [2,001 × 2] 6 show_product_quicklyTRUE -0.839 0.740 2.23 0.05 percentile [2,001 × 2] 7 use_sexTRUE -4.04 -2.43 -0.952 0.05 percentile [2,001 × 2]
Посмотри на результаты bootstrap
Мы можем визуализировать наши результаты. Обратим внимание, что если бы мы не установил параметр keep_reps = TRUE, то у нас были бы только интервалы и мы бы не смогли построить график №3.
youtube_intervals %>% mutate(term = str_remove(term, "TRUE"), term = fct_reorder(term, .estimate)) %>% ggplot(aes(.estimate, term)) + geom_vline(xintercept = 0, size = 1.5, lty = 2, color = "gray80") + geom_errorbarh(aes(xmin = .lower, xmax = .upper), size = 1.5, alpha = 0.5, color = "midnightblue") + geom_point(size = 3, color = "midnightblue") + labs(x = "Increase in year for each commercial characteristic", y = NULL)
outube_intervals %>% mutate(term = str_remove(term, "TRUE"), term = fct_reorder(term, .estimate)) %>% unnest(.replicates) %>% ggplot(aes(estimate, fill = term)) + geom_vline(xintercept = 0, size = 1.5, lty = 2, color = "gray50") + geom_histogram(alpha = 0.8, show.legend = FALSE) + facet_wrap(vars(term))
Таким образом, мы видим, что рекламные ролики содержат меньше юмора и сексуального контента, и больше рубрик селебрити и патриотических тем.
