Чистый синтаксис оценки в пользовательской функции R

Я хочу определить общую функцию

func_boxplot2 <- function(tmp, xvar, yvar, groupvar)
{
  xvar <- enquo(xvar)
  yvar <- enquo(yvar)
  groupvar <- enquo(groupvar)

  # If variable yield exists, put concentrations to NA for all yields < annual_yield_thres
  if( "yield" %in% names(tmp) )
  {
    tmp <- tmp %>%
      mutate_at(vars(!!yvar), ~ifelse(round(yield, 0) < 85, NA, .))
  }
  
  # Compute IQR for each year
  tmp_iqr <- tmp %>%
    group_by(!!groupvar) %>%
    summarise(iqr=IQR(!!yvar, na.rm = TRUE))
  
  p <- ggplot(data = tmp %>% mutate_at(vars(!!yvar), ~ifelse(tmp_iqr[which(tmp_iqr[[!!groupvar]] %in% (!!xvar)),]$iqr == 0, . + runif(1, -0.01, 0.01), . )), aes(x = !!xvar, y = !!yvar))
  p <- p + stat_boxplot(aes(group = !!groupvar), na.rm = TRUE, coef = 10000)   # Trick (large unrealistic coef value) so whiskers end at min(y) & max(y)
  p <- p + geom_boxplot(na.rm = TRUE, outlier.shape = NA)

  return(p)
}

который способен отображать усы коробчатой диаграммы, простирающиеся до минимума / максимума, даже когда IQR равен 0. Я пытаюсь добиться этого, добавляя крошечные случайные числа (ниже уровня значимости) к инкриминируемым данным, чтобы избежать IQR = 0.

Однако я, должно быть, пропустил что-то в синтаксисе о quosure, потому что запуск этой функции

func_boxplot2(data, date, days, date)

с набором данных

structure(list(date = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 2L, 
2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 
4L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 
7L, 7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 9L, 9L, 10L, 
10L, 10L, 10L, 10L, 10L), .Label = c("2010", "2011", "2012", 
"2013", "2014", "2015", "2016", "2017", "2018", "2019"), class = c("ordered", 
"factor")), station = c("41B011", "41MEU1", "41N043", "41R001", 
"41R012", "41WOL1", "41B011", "41MEU1", "41N043", "41R001", "41R012", 
"41WOL1", "41B011", "41MEU1", "41N043", "41R001", "41R012", "41WOL1", 
"41B011", "41MEU1", "41N043", "41R001", "41R012", "41WOL1", "41B011", 
"41MEU1", "41N043", "41R001", "41R012", "41WOL1", "41B011", "41MEU1", 
"41N043", "41R001", "41R012", "41WOL1", "41B011", "41MEU1", "41N043", 
"41R001", "41R012", "41WOL1", "41B011", "41MEU1", "41N043", "41R001", 
"41R012", "41WOL1", "41B011", "41MEU1", "41N043", "41R001", "41R012", 
"41WOL1", "41B011", "41MEU1", "41N043", "41R001", "41R012", "41WOL1"
), days = c(16, 15, 45, 26, 14, 14, 32, 7, 87, 42, 24, 23, 25, 
25, 55, 29, 29, 16, 11, 14, 58, 21, 19, 10, 10, 14, 33, 18, 10, 
7, 9, 10, 19, 7, 8, 7, 1, 5, 15, 8, 1, 4, 5, 6, 14, 6, 5, 5, 
3, 5, 19, 8, 4, 5, 3, 4, 16, 3, 1, 3), yield = c(98.4817351598173, 
49.4520547945205, 95.8561643835616, 97.6712328767123, 98.2648401826484, 
95.1598173515982, 97.8767123287671, 27.9109589041096, 98.310502283105, 
98.972602739726, 97.203196347032, 96.2100456621005, 98.7818761384335, 
96.7554644808743, 97.4954462659381, 98.8046448087432, 98.747723132969, 
98.3037340619308, 99.0525114155251, 96.1986301369863, 97.1004566210046, 
96.4954337899543, 96.3698630136986, 98.2077625570776, 96.62100456621, 
98.3675799086758, 95.6963470319635, 96.8835616438356, 93.5844748858447, 
87.8196347031963, 91.2328767123288, 92.5570776255708, 81.5182648401827, 
82.7739726027397, 90.1826484018265, 87.1461187214612, 87.2153916211293, 
92.9986338797814, 94.6948998178506, 85.5760473588342, 92.3611111111111, 
96.2204007285975, 86.3698630136986, 86.4269406392694, 87.796803652968, 
93.2762557077626, 96.6438356164384, 95.6164383561644, 71.3812785388128, 
93.7442922374429, 96.3698630136986, 97.2602739726027, 95.7876712328767, 
94.7146118721461, 87.6141552511416, 43.0821917808219, 88.6872146118722, 
92.6826484018265, 90.365296803653, 86.541095890411), environ = structure(c(5L, 
4L, 6L, 3L, 5L, 3L, 5L, 4L, 6L, 3L, 5L, 3L, 5L, 4L, 6L, 3L, 5L, 
3L, 5L, 4L, 6L, 3L, 5L, 3L, 5L, 4L, 6L, 3L, 5L, 3L, 5L, 4L, 6L, 
3L, 5L, 3L, 5L, 4L, 6L, 3L, 5L, 3L, 5L, 4L, 6L, 3L, 5L, 3L, 5L, 
4L, 6L, 3L, 5L, 3L, 5L, 4L, 6L, 3L, 5L, 3L), .Label = c("Urbain avec très forte influence du trafic", 
"Urbain avec forte influence du trafic", "Urbain avec influence modérée du trafic", 
"Urbain avec faible influence du trafic", "Urbain avec très faible influence du trafic", 
"Industriel avec influence modérée du trafic"), class = "factor")), row.names = c(NA, 
-60L), class = c("tbl_df", "tbl", "data.frame"))

дает мне следующие ошибки

 Error: Problem with `mutate()` input `days`.
x Must extract column with a single valid subscript.
x Subscript `date` has size 60 but must be size 1.
ℹ Input `days` is `(structure(function (..., .x = ..1, .y = ..2, . = ..1) ...`.

Что не так в моем синтаксисе, пожалуйста?

Большое спасибо,

======== ОБНОВЛЕНИЕ ==========

Использование предлагаемой обновленной функции

func_boxplot2 <- function(tmp, xvar, yvar, groupvar)
{
  # If variable yield exists, put concentrations to NA for all yields < annual_yield_thres
  if("yield" %in% names(tmp)) {
    tmp <-
      tmp %>%
      mutate(across({{yvar}}, ~ifelse(round(yield, 0) < 85, NA, .)))
  }
  
  tmp <-
    tmp %>%
    group_by({{groupvar}}) %>%
    mutate(
      across({{yvar}}, function (x) {
        ifelse(
          IQR({{yvar}}, na.rm = TRUE) == 0,
          x + runif(1, -0.01,0.01),
          x
        )
      })
    )

  ggplot(tmp, aes(x = {{xvar}}, y = {{yvar}})) +
    stat_boxplot(aes(group = {{groupvar}}), na.rm = TRUE, coef = 10000) +
    geom_boxplot(na.rm = TRUE, outlier.shape = NA)
}

приводит к следующему сюжету

Как указано в моем комментарии, кажется, что обработка вывода tmp для всех строк того же года является первым значением yvar этого года, которое объясняет сюжет. Действительно, комментируя этот блок, мы получаем следующий рисунок.

Alessandro 13.11.2020 источник

Ответы (2)

arrow_upward
0
arrow_downward

Условие в ifelse не соответствует длине. Вы можете изменить на это аргумент data ggplot.

data =
  tmp %>%
    group_by(!!groupvar) %>%
    mutate_at(
        vars(!!yvar),
        if (IQR(., na.rm = TRUE) == 0) {
          . + runif(1, -0.01,0.01)
        } else {
          .
        }
    )

Вы используете quosure и !! правильно, однако вы должны использовать более свежий _ 7_ оператор.

Это обновленная функция

func_boxplot2 <- function(tmp, xvar, yvar, groupvar)
{
  # If variable yield exists, put concentrations to NA for all yields < annual_yield_thres
  if("yield" %in% names(tmp)) {
    tmp <-
      tmp %>%
      mutate(across({{yvar}}, ~ifelse(round(yield, 0) < 85, NA, .)))
  }
  
  tmp <-
    tmp %>%
    group_by({{groupvar}}) %>%
    mutate(
      across({{yvar}}, function (x) {
        if (IQR(x, na.rm = TRUE) == 0) {
          x + runif(length(x), -0.01, 0.01)
        } else {
          x
        }
      })
    )

  ggplot(tmp, aes(x = {{xvar}}, y = {{yvar}})) +
    stat_boxplot(aes(group = {{groupvar}}), na.rm = TRUE, coef = 10000) +
    geom_boxplot(na.rm = TRUE, outlier.shape = NA)
}

func_boxplot2(data, date, days, date)

Paul 13.11.2020

comment

Спасибо ! Я переработал старый код, но конвертирую его в синтаксис {{. Однако есть еще кое-что странное, потому что коробочные сюжеты не работают. group_by и ifelse заменяют каждый год все значения первым. - Alessandro; 13.11.2020

comment

Я обновил свой вопрос, так как сюжет, полученный в результате этой функции, все еще не работает (см. Рисунки). - Alessandro; 15.11.2020

comment

@Alessandro, см. Мой обновленный ответ. Я заменил ifelse на _2 _ / _ 3_, поскольку IQR(x, na.rm = TRUE) == 0 возвращает только одно значение. - Paul; 15.11.2020

arrow_upward
0
arrow_downward

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

Обычно статистика равна StatBoxplot, и важной частью этого ggproto объекта является метод compute_group. При использовании в графике он возвращает фрейм данных с одной строкой, содержащий столбцы.

[1] "ymin"        "lower"       "middle"      "upper"      
[5] "ymax"        "outliers"    "notchupper"  "notchlower" 
[9] "x"           "relvarwidth" "flipped_aes"

В основном они имеют очевидное значение; единственный неочевидный - это outliers, столбец в режиме списка, содержащий единственный числовой вектор, содержащий выбросы.

Итак, чтобы полностью избавиться от построения графика выбросов, вы можете создать унаследованную статистику, которая похожа на StatBoxplot, но изменяет результат compute_group:

NoOutlierStatBoxplot <- 
  ggproto("NoOutlierStatBoxplot", ggplot2::StatBoxplot,
          compute_group = function(..., self) {
            res <- ggproto_parent(StatBoxplot, self)$compute_group(...)
            res$ymin <- min(c(res$ymin, res$outliers[[1]]))
            res$ymax <- max(c(res$ymax, res$outliers[[1]]))
            res$outliers <- list(numeric())
            res
          })

(На самом деле это не идентично тому, что вы делали: он по-прежнему вычисляет верхний и нижний квартили после удаления выбросов. Если это важно для вас, вам может потребоваться более обширная модификация.)

С помощью этой модификации вы можете отбросить большую часть кода из func_boxplot2, включая удаление stat_boxplot():

func_boxplot3 <- function(tmp, xvar, yvar, groupvar)
{
  # If variable yield exists, put concentrations to NA for all yields < annual_yield_thres
  if("yield" %in% names(tmp)) {
    tmp <-
      tmp %>%
      mutate(across({{yvar}}, ~ifelse(round(yield, 0) < 85, NA, .)))
  }
  
  ggplot(tmp, aes(x = {{xvar}}, y = {{yvar}})) +
    geom_boxplot(na.rm = TRUE, outlier.shape = NA, 
                 aes(group = {{groupvar}}), 
                 stat = NoOutlierStatBoxplot)
}

func_boxplot3(mydf, date, days, date)

user2554330 15.11.2020

comment

Действительно, это наиболее элегантный подход вместо изменения данных. Однако ваше решение выходило за рамки моих нынешних познаний в программировании на R, но оно определенно улучшит мою точку зрения. Спасибо. - Alessandro; 16.11.2020

comment

Виньетка ggplot2 Расширение ggplot2 дает довольно читаемый обзор такого рода вещей. - user2554330; 16.11.2020

Чистый синтаксис оценки в пользовательской функции R

Ответы (2)

Вопросы по теме