Прогнозирование временных рядов в R; построение графиков событий и создание новых графиков прогнозов с указанным диапазоном дат после первоначального прогноза

Я создал функцию, которая позволяет мне выполнять прогнозирование временных рядов с помощью пакета fable. Идея функции заключалась в том, чтобы анализировать наблюдаемые и прогнозируемые значения после определенной даты / события. Вот фиктивный фрейм данных, который генерирует столбец дат: -

set.seed(1)
df <- data.frame(Date = sort(sample(seq(as.Date('2018/01/01'), as.Date('2020/09/17'), by="day"),1368883, replace = T)))

А вот и созданная мной функция. Вы указываете данные, затем дату события, затем период прогноза в днях и, наконец, заголовок для вашего графика.

event_analysis<-function(data,eventdate,period,title){
  require(dplyr)
  require(tsibble)
  require(fable)
  require(fabletools)
  require(imputeTS)
  require(ggplot2)
  data_count<-data%>%
    group_by(Date)%>%
    summarise(Count=n())
  
  data_count<-as_tsibble(data_count)
  data_count<-na_mean(data_count)
  
  
  train <- data_count %>%
    #sample_frac(0.8)
    filter(Date<=as.Date(eventdate))
  
  fit <- train %>%
    model(
      ets = ETS(Count),
      arima = ARIMA(Count),
      snaive = SNAIVE(Count)
    ) %>%
    mutate(mixed = (ets + arima + snaive) / 3)
  
  
  fc <- fit %>% forecast(h = period)
  
  forecastplot<-fc %>%
    autoplot(data_count, level = NULL)+ggtitle(title)+
    geom_vline(xintercept = as.Date(eventdate),linetype="dashed",color="red")+
    labs(caption = "Red dashed line = Event occurrence")
                                                                 
  
  fc_accuracy<-accuracy(fc,data_count)
  
  #obs<-data_count
  #colnames(obs)[2]<-"Observed"
  #obs_pred<-merge(data_count,fc_accuracy, by="Date")
  return(list(forecastplot,fc_accuracy,fc))
}

И за один прогон я указываю df, дату события, количество дней, которые я хочу спрогнозировать (3 недели), а затем заголовок: -

event_analysis(df, "2020-01-01",21,"Event forecast")

Что напечатает этот результат и построит график: -

введите описание изображения здесь

введите описание изображения здесь

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

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

  1. для простоты интерпретации
  2. Я хочу видеть количество дней N до и количество дней N после даты события (N представляет период прогноза, т.е. 21).

Итак, дополнительный график (вместе с исходным полным прогнозом), который будет выглядеть так, возможно, в стиле одного вывода и нескольких графиков: -

введите описание изображения здесь

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

По сути, это две дополнительные вещи, которые я хочу добавить к своей функции, но я не уверен, как это сделать. Любая помощь очень ценится :).


person Robin Turkington    schedule 16.10.2020    source источник
comment
с coord_cartesian вы сможете увеличивать масштаб вашего графика. В качестве альтернативы вы можете преобразовать свою диаграмму в диаграмму plotly, и вы можете увеличивать масштаб так, как вам нужно.   -  person Edo    schedule 16.10.2020
comment
для записи, прикрепление библиотек внутри функции - это вторичный эффект, который обычно не принимается во внимание. Использование должно предоставить сообщение о том, что вы подключаете новые библиотеки, вы должны сделать это извне, или, если вы добавляете свою функцию в пакет, вы можете импортировать каждую отдельную функцию, которая вам нужна.   -  person Edo    schedule 16.10.2020
comment
данные вашего воспроизводимого примера можно переписать так: df <- data.frame(Date = sort(sample(seq(as.Date('2018/01/01'), as.Date('2020/09/17'), by="day"),1368883, replace = T))). На данный момент выдает ошибку   -  person Edo    schedule 16.10.2020
comment
Привет @Edu, спасибо, что указали, что воспроизводимый пример не работает. С тех пор я заменил данные вашим решением. Что касается прикрепления библиотек к функции, я не знал, что это создает какой-либо вторичный эффект. Хотя позже в него будут внесены соответствующие поправки, меня беспокоит только то, что он работает, и в некоторой степени работает.   -  person Robin Turkington    schedule 16.10.2020


Ответы (1)


Я полагаю, вы могли бы переписать это так. Я сделал несколько поправок, чтобы помочь вам.

set.seed(1)
df <- data.frame(Date = sort(sample(seq(as.Date('2018/01/01'), as.Date('2020/09/17'), by="day"),1368883, replace = T)))

event_analysis <- function(data, eventdate, period, title){
 
 # in the future, you may think to move them out
 library(dplyr)
 library(tsibble)
 library(fable)
 library(fabletools)
 library(imputeTS)
 library(ggplot2)
 
 # convert at the beginning
 eventdate <- as.Date(eventdate)
 
 # more compact sintax
 data_count <- count(data, Date, name = "Count")
 
 # better specify the date variable to avoid the message
 data_count <- as_tsibble(data_count, index = Date)
 
 # you need to complete missing dates, just in case
 data_count <- tsibble::fill_gaps(data_count)
 
 
 data_count <- na_mean(data_count)
 
 
 train <- data_count %>%
  filter(Date <= eventdate)
 
 test <- data_count %>%
  filter(Date > eventdate, Date <= (eventdate+period))
 
  fit <- train %>%
  model(
   ets    = ETS(Count),
   arima  = ARIMA(Count),
   snaive = SNAIVE(Count)
  ) %>%
  mutate(mixed = (ets + arima + snaive) / 3)
 
 
 fc <- fit %>% forecast(h = period)
 

 # your plot
 forecastplot <- fc %>%
  autoplot(data_count, level = NULL) + 
  ggtitle(title) +
  geom_vline(xintercept = as.Date(eventdate), linetype = "dashed", color = "red") +
  labs(caption = "Red dashed line = Event occurrence")
 
 
 # plot just forecast and test
 zoomfcstplot <- autoplot(fc) + autolayer(test, .vars = Count)
 
 fc_accuracy <- accuracy(fc,data_count)
 

 ### EDIT: ###

 # results vs test
 res <- fc %>% 
  as_tibble() %>% 
  select(-Count) %>% 
  tidyr::pivot_wider(names_from = .model, values_from = .mean) %>% 
  inner_join(test, by = "Date")

 ##############
 

 return(list(forecastplot = forecastplot,
             zoomplot     = zoomfcstplot,
             accuracy     = fc_accuracy,
             forecast     = fc,
             results      = res))
}


event_analysis(df, 
               eventdate = "2020-01-01",
               period    = 21,
               title     = "Event forecast")


Выход:

#> 
#> 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
#> Carico il pacchetto richiesto: fabletools
#> Registered S3 method overwritten by 'quantmod':
#>   method            from
#>   as.zoo.data.frame zoo
#> $forecastplot

введите описание изображения здесь

#> 
#> $zoomplot

введите описание изображения здесь

#> 
#> $accuracy
#> # A tibble: 4 x 9
#>   .model .type    ME  RMSE   MAE   MPE  MAPE  MASE    ACF1
#>   <chr>  <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>   <dbl>
#> 1 arima  Test  -16.8  41.8  35.2 -1.31  2.61 0.791  0.164 
#> 2 ets    Test  -16.8  41.8  35.2 -1.31  2.61 0.791  0.164 
#> 3 mixed  Test  -21.9  44.7  36.8 -1.68  2.73 0.825 -0.0682
#> 4 snaive Test  -32.1  57.3  46.6 -2.43  3.45 1.05  -0.377 
#> 
#> $forecast
#> # A fable: 84 x 4 [1D]
#> # Key:     .model [4]
#>    .model Date               Count .mean
#>    <chr>  <date>            <dist> <dbl>
#>  1 ets    2020-01-02 N(1383, 1505) 1383.
#>  2 ets    2020-01-03 N(1383, 1505) 1383.
#>  3 ets    2020-01-04 N(1383, 1505) 1383.
#>  4 ets    2020-01-05 N(1383, 1505) 1383.
#>  5 ets    2020-01-06 N(1383, 1505) 1383.
#>  6 ets    2020-01-07 N(1383, 1505) 1383.
#>  7 ets    2020-01-08 N(1383, 1505) 1383.
#>  8 ets    2020-01-09 N(1383, 1505) 1383.
#>  9 ets    2020-01-10 N(1383, 1505) 1383.
#> 10 ets    2020-01-11 N(1383, 1505) 1383.
#> # ... with 74 more rows
#>
#> $results
#> # A tibble: 21 x 6
#>    Date         ets arima snaive mixed Count
#>    <date>     <dbl> <dbl>  <dbl> <dbl> <int>
#>  1 2020-01-02 1383. 1383.   1386 1384.  1350
#>  2 2020-01-03 1383. 1383.   1366 1377.  1398
#>  3 2020-01-04 1383. 1383.   1426 1397.  1357
#>  4 2020-01-05 1383. 1383.   1398 1388.  1415
#>  5 2020-01-06 1383. 1383.   1431 1399.  1399
#>  6 2020-01-07 1383. 1383.   1431 1399.  1346
#>  7 2020-01-08 1383. 1383.   1350 1372.  1299
#>  8 2020-01-09 1383. 1383.   1386 1384.  1303
#>  9 2020-01-10 1383. 1383.   1366 1377.  1365
#> 10 2020-01-11 1383. 1383.   1426 1397.  1328
#> # ... with 11 more rows 
person Edo    schedule 16.10.2020
comment
это здорово, спасибо, что нашли время предоставить решение (а также привести в порядок мой код!). То, что я имею в виду под вышеизложенным, было бы похоже на идею, когда вы печатаете forecast, но я бы хотел, чтобы значения Count в тестовом периоде (или периоде прогноза, т.е. h = 21) были столбцом. Итак, один столбец с датой прогноза, один с количеством, а затем еще четыре столбца (каждый представляет модели) и их прогнозируемые значения. Я надеюсь это имеет смысл? - person Robin Turkington; 16.10.2020
comment
Я добавил новую часть. Это то, что вы искали? - person Edo; 16.10.2020
comment
Это именно то, чего я хотел, еще раз большое спасибо! - person Robin Turkington; 16.10.2020