Введение

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

Наиболее часто для оценки эффективности прогнозирования используется относительный (процентный) показатель ошибки - MAPE, однако у данного показателя есть ряд недостатков, которые проявляются в случае, когда временной ряд содержит нулевые или очень малые наблюдения. Даже когда нет проблемы нулевых значений, качество прогнозирования сложно оценить лишь с использованием MAPE, поскольку для каждого конкретного объекта прогнозирования действует свой набор факторов (новизна продукта, маркетинговая активность, количество и тип клиентов, степень агрегирования), который определяет пределы точности прогнозирования. Поэтому нельзя ориентироваться на усредненные, “отраслевые” значения ошибки прогнозирования и необходимо сравнивать подходы к прогнозированию в одинаковых условиях. С этой целью будут рассмотрены показатели, основанные на сопоставлении ошибки прогноза применяемого для прогнозирования метода с показателями ошибки простых, “наивных” методов прогнозирования на тех же данных.

Также будут рассмотрены методы оценки качества моделей прогнозирования, основанные на разделении данных на “обучающий” и “тестовый” периоды, а также скользящем контроле.

Материал составлен на основе раздела 3.1 и 3.4 книги: Hyndman R, Athanasopoulos G. Forecasting Principles and Practice.

library(tidyverse) # визуализация и трансформация данных
library(forecast) # анализ временных рядов и прогнозирование
library(fpp) # Примеры временных рядов
library(sophisthse) # Загрузка временных рядов из базы Sophist

Простые методы прогнозирования

Описание

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

Метод усреднения (Average Method)

В этом методе в качестве прогноза на будущие периоды берется среднее значение ряда, вычисленное на всем историческом периоде:

\[ \hat{y}_{T+h} = \frac{1}{T} \sum_{t=1}^{T} y_t \]

Метод также может применяться и для неупорядоченных данных (берется среднее по выборке).

В R для реализации этого метода используется функция: forecast::meanf(y, h), где y - временной ряд, а h - горизонт прогнозирования.

Наивный метод (Naïve method)

Наивный прогноз - это простое повторение последнего наблюдения в качестве прогноза на будущие периоды. В R метод реализован двумя функциями: forecast::naive(y, h) и forecast::rwf(y, h).

Сезонный наивный метод (Seasonal naïve method)

Этот метод использует тот же принцип, однако повторяется не последнее наблюдение, а последнее наблюдение для одноименного периода. Например, прогноз на декабрь 2018 года, сделанный в 2017 году, будет равен значению за декабрь 2017 года.

В R этот метод реализован функцией forecast::snaive(y, h)

Метод дрейфа (Drift Method)

Дрейф - это медленное изменение уровня ряда. В методе дрейфа вычисляется средний прирост уровня ряда по всем смежным периодам. Затем это значение используется для прогнозирования будущих значений.

\[ \hat{y}_{T+h} = y_T + \frac{h}{T-1} \sum_{t=2}^{T} (y_t - y_{t-1})\]

Можно показать, что данное выражение эквивалентно:

\[ \hat{y}_{T+h} = y_T + h \frac{y_T - y_1}{T-1} \]

Т.е. в методе дрейфа прогнозом является продолжение прямой линии, проведенной через первое и последнее наблюдения.

В R метод дрейфа реализован функцией: forecast::rwf(y, h, drift=TRUE).

Примеры

# Отбор данных
beer2 <- ausbeer %>% 
  window(start = 1992, end = c(2005, 4))

beer_label <- 'Поквартальный прогноз производства пива в Австралии'
beer_unit <- 'млн. л'

# Прогнозирование
beer_m <- meanf(beer2, h = 11) # по среднему значению ряда
beer_n <- naive(beer2, h = 11) # наивный прогноз
beer_sn <- snaive(beer2, h = 11) # сезонный наивный прогноз
beer_d <- rwf(beer2, h = 11, drift = T) # метод дрейфа

Для визуализации прогнозов можно воспользоваться функциями autoplot() и autolayer():

# Сравнение нескольких прогнозов
autoplot(beer2, series = 'История') +
  autolayer(beer_m, PI = FALSE, series = 'Среднее ряда') +
  autolayer(beer_n, PI = FALSE, series = 'Наивный прогноз') +
  autolayer(beer_sn, PI = FALSE, series = 'С. наивный прогноз') +
  autolayer(beer_d, PI = FALSE, series = 'Метод дрейфа') +
  labs(title = beer_label, x = NULL, color = NULL, y = beer_unit)

Пакет forecast рассчитывает и доверительные интервалы для прогноза, которые легко получить с в виде графика или таблицы. В предыдущем примере мы отключили вывод доверительных интервалов для прогноза с помощью параметра PI = FALSE, чтобы не загромождать график.

# График
autoplot(beer_sn) + 
  labs(title = beer_label, y = beer_unit, 
       x = NULL, color = NULL)

# Таблица
beer_sn
##         Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 2006 Q1            416 393.9027 438.0973 382.2051 449.7949
## 2006 Q2            403 380.9027 425.0973 369.2051 436.7949
## 2006 Q3            408 385.9027 430.0973 374.2051 441.7949
## 2006 Q4            482 459.9027 504.0973 448.2051 515.7949
## 2007 Q1            416 384.7497 447.2503 368.2068 463.7932
## 2007 Q2            403 371.7497 434.2503 355.2068 450.7932
## 2007 Q3            408 376.7497 439.2503 360.2068 455.7932
## 2007 Q4            482 450.7497 513.2503 434.2068 529.7932
## 2008 Q1            416 377.7264 454.2736 357.4655 474.5345
## 2008 Q2            403 364.7264 441.2736 344.4655 461.5345
## 2008 Q3            408 369.7264 446.2736 349.4655 466.5345

Доверительные интервалы позволяют оценить степень надежности прогноза. Чем шире интервал, тем выше неопределенность.

Работа с компонентами прогноза

Объект, который возвращают функции для прогнозирования, имеет следующую структуру:

names(beer_d)
##  [1] "method"    "model"     "lambda"    "x"         "fitted"   
##  [6] "residuals" "series"    "mean"      "level"     "lower"    
## [11] "upper"

Наиболее полезными компонентами являются mean (точечный прогноз), fitted (прогноз в историческом периоде) и residuals - остатки прогноза. С назначением остальных компонентов можно познакомиться в справке: ?forecast.

Для выделения прогноза в историческом периоде и остатков удобно использовать функции fitted() и residuals() соответственно.

# Подогнанные значения
autoplot(beer2, series = 'История') +
  autolayer(fitted(beer_sn), series = 'Прогноз') +
  labs(title = beer_label, y = beer_unit, 
       x = NULL, color = NULL)
## Warning: Removed 4 rows containing missing values (geom_path).

# Остатки
autoplot(residuals(beer_sn)) +
labs(title = 'Ошибки прогноза объема производства пива', 
     y = beer_unit, x = NULL, color = NULL)

Пример прогнозирования показателя из ЕАЭСД

Рассчитаем прогноз показателя из ЕАЭСД (sophist): Ввод в действие жилых домов, млн.кв.метров (CONSTR_Y_NAT)

constr <- sophisthse('CONSTR_Y_NAT')

#sophisthse_metadata(constr)

constr_label <- 'Ежегодный ввод в действие жилых домов (CONSTR_Y_NAT)'
constr_unit <- 'млн.кв.м'

# Оставляем данные с 2005 года
constr <- constr %>% 
  window(start = 2005)

# Визуализация
autoplot(constr) +
  labs(title = constr_label, y = constr_unit, x = NULL)

Прогнозирование ввода новых жилых домов с помощью рассмотренных методов.

constr_m <- meanf(constr, h = 10)
constr_n <- naive(constr, h = 10)
constr_d <- rwf(constr, h = 10, drift = T)


autoplot(constr, series = 'История') +
  autolayer(constr_m, PI = FALSE, series = 'Среднее ряда') +
  autolayer(constr_n, PI = FALSE, series = 'Наивный прогноз') +
  autolayer(constr_d, PI = FALSE, series = 'Метод дрейфа') +
  labs(title = constr_label, x = NULL, color = NULL, y = constr_unit)

Оценка точности прогнозирования

Показатели ошибки прогнозирования

Ошибка прогноза в периоде \(t\) - это разность фактических и прогнозных значений:

\[ e_t = y_t -\hat{y}_t \]

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

Все показатели ошибки прогнозирования можно разделить на две группы:

  • абсолютные - зависят от масштаба пронозируемой величины (scale-dependent) и, как правило, измеряются в тех же единицах, что и эта величина;

  • относительные - представляют собой безразмерные отношения и не зависят от масштаба (scale-independent). Показатели этой группы удобно использовать для сравнения точности прогнозирования для разных объектов.

Еще одним критерием для классификации показателей может служить их назначение. Одни показатели характеризуют разброс прогнозных и фактических значений, другие - систематическую ошибку, т.е. смещение прогноза.

Средняя абсолютная ошибка и стандартная ошибка (MAE и RMSE)

Наиболее популярными абсолютными показателями являются средняя абсолютная ошибка (Mean Absolute Error, MAE) и стандартная ошибка (Root Mean Squared Error, RMSE):

\[ MAE = mean(|e_t|) \]

\[ RMSE = \sqrt{mean(e_t^2)} \]

Оба показателя позволяют сравнивать точность прогнозирования при работе с одним и тем же рядом. MAE легче вычисляется и интерпретируется, а RMSE применяется для построение доверительных интервалов и расчета страхового запаса.

Средняя ошибка (ME)

Кроме того, к абсолютным относится еще один показатель - средняя ошибка (Mean Error, ME), который позволяет оценить систематическую ошибку (смещение) прогноза:

\[ ME = mean(e_t) \]

Средняя абсолютная ошибка и средняя ошибка в процентах (MAPE и MPE)

Наиболее популярным относительным показателем является средняя абсолютная ошибка в процентах (Mean Absolute Percentage Error, MAPE). Этот показатель вычисляется путем усреднения относительных ошибок прогноза:

\[ p_t = 100 \frac{e_t} {y_t} \]

\[ MAPE = mean(|p_t|) \]

Для оценки систематической ошибки прогноза применяется показатель средняя ошибка в процентах (Mean Percentage Error, MPE):

\[ MPE = mean(p_t) \]

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

Недостатком этих показателей является невозможность их вычисления, если встречаются нулевые значения ряда (например, при прогнозировании спроса на малооборачиваемые или сезонные товары). Величина MAPE получается очень большой, когда фактические значения спроса \(y_t < 1\).

MAPE и MPE нельзя применять, когда у показателя нет абсолютного нуля (например, при прогнозировании температуры в шкале Цельсия), поскольку для таких данных не имеет смысла само понятие отношения (деления).

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

Нормированная ошибка (MASE)

Нормированная ошибка (scaled error) для несезонного ряда:

\[ q_j = \frac{e_j} {\frac{1}{T - 1} \sum_{t=2}^T |y_t - y_{t-1}|} \]

Нормированная ошибка для сезонного ряда:

\[ q_j = \frac{e_j} {\frac{1}{T - m} \sum_{t=m+1}^T |y_t - y_{t-m}|} \]

Для неупорядочынных (перекрестных) данных нормированная ошибка вычисляется относительно стандартного отклонения прогнозируемой величины:

\[ q_j = \frac{e_j} {\frac{1}{N} \sum_{i=1}^N |y_i - \bar{y}|} \]

Средняя абсолютная нормированная ошибка (Mean Absolute Scaled Error) вычисляется путем усреднения нормированных ошибок для всех периодов:

\[ MASE = mean(|q_j|) \]

Для вычисления показателей ошибки в R используется функция forecast::accuracy().

Пример - ошибка прогноза объема производства пива

Ряд сезонный, поэтому наилучший результат среди простых методов должен давать сезонный наивный прогноз.

# Метод усреднения
autoplot(beer2, series = 'История') +
  autolayer(fitted(beer_m), series = 'Среднее ряда') +
  labs(title = beer_label, x = NULL, color = NULL, y = beer_unit)

# Наивный прогноз
autoplot(beer2, series = 'История') +
  autolayer(fitted(beer_n), series = 'Наивный прогноз') +
  labs(title = beer_label, x = NULL, color = NULL, y = beer_unit)
## Warning: Removed 1 rows containing missing values (geom_path).

# Сезонный наивный прогноз
autoplot(beer2, series = 'История') +
  autolayer(fitted(beer_sn), series = 'С. наивный прогноз') +
  labs(title = beer_label, x = NULL, color = NULL, y = beer_unit)
## Warning: Removed 4 rows containing missing values (geom_path).

# Метод дрейфа
autoplot(beer2, series = 'История') +
  autolayer(fitted(beer_d), series = 'Метод дрейфа') +
  labs(title = beer_label, x = NULL, color = NULL, y = beer_unit)
## Warning: Removed 1 rows containing missing values (geom_path).

Показатели ошибки прогноза автоматически рассчитываются функцией forecast::accuracy():

accuracy(beer_m)
##                        ME    RMSE      MAE        MPE     MAPE     MASE
## Training set 8.121418e-15 44.1763 35.91135 -0.9510944 7.995509 2.444228
##                    ACF1
## Training set -0.1256697

Для удобства сравнения прогнозов, объединим результаты по всем прогнозам в одну таблицу:

rbind(accuracy(beer_m),
      accuracy(beer_n),
      accuracy(beer_sn),
      accuracy(beer_d)) %>%
  as_tibble() %>%
  round(2) %>%
  mutate(`Метод` = c('Среднее', 'Наивный', 'С. наивный', 'Дрейф')) %>%
  select(`Метод`, MASE, MAPE, everything()) %>%
  arrange(MASE)
## # A tibble: 4 x 8
##   Метод       MASE  MAPE    ME  RMSE   MAE   MPE  ACF1
##   <chr>      <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 С. наивный  1      3.4 -1.85  17.2  14.7 -0.48 -0.34
## 2 Среднее     2.44   8    0     44.2  35.9 -0.95 -0.13
## 3 Наивный     3.77  12.3  0.71  66.6  55.4 -0.9  -0.25
## 4 Дрейф       3.78  12.3  0     66.6  55.5 -1.06 -0.25

Показатель ошибки MASE для сезонного наивного прогноза равен 1, поскольку для сезонного ряда (когда задана частота) именно этот прогноз используется как основа для сравнения. Усреднение в ~2.5 раза хуже по этому показателю, чем сезонный наивный прогноз.

Примечание: по умолчанию выводится слишком большое количество значащих цифр. Для упрощения интерпретации показателей можно округлить все значения в таблице с помощью функции round(). Обратите внимание, что это надо делать до того, как в таблицу будет добавлен текстовый столбец с названиями методов.

Прогнозирование Ex-Post

При выборе метода прогнозирования аналитика интересует не то, как эта модель объясняла прошлое, а то, как она предсказывает будущее. Поэтому показатели ошибки, вычисленные на историческом периоде могут вводить в заблуждение. Можно построить такую модель, которая “предсказывает” уже известное прошлое идеально точно.

Поэтому распространенным приемом для выбора модели является разделение ряда на два периода - обучающий (training set) и тестовый (test set). При этом модель строится на обучающем периоде, а показатели ошибки вычисляются на тестовом. Этот способ называется прогнозирование Ex-Post.

Прогнозирование Ex-Post

Обычно в качестве тестового множества выбирается порядка 20% наблюдений, однако это ориентировочная величина, которую можно скорректировать исходя из объема имеющихся данных и горизонта прогнозирования. Рекомендуется делать тестовый период как минимум такой же длины, как и горизонт прогнозирования.

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

autoplot(ausbeer, series = 'Полный ряд') +
  autolayer(beer2, series = 'Сокращенный ряд')

Прогнозы в тестовом периоде можно сравнить с фактическими данными визуально.

# Сравнение нескольких прогнозов
autoplot(window(ausbeer, start = 2004), series = 'Факт') +
  autolayer(beer_m, PI = FALSE, series = 'Среднее ряда') +
  autolayer(beer_n, PI = FALSE, series = 'Наивный прогноз') +
  autolayer(beer_sn, PI = FALSE, series = 'С. наивный прогноз') +
  autolayer(beer_d, PI = FALSE, series = 'Метод дрейфа') +
  labs(title = beer_label, x = NULL, color = NULL, y = beer_unit)

Функция accuracy() автоматически вычисляет показатели ошибки на тестовом периоде, если предоставить данные тестового периода, или просто весь ряд целиком.

# Метод усреднения
accuracy(beer_m, ausbeer)
##                         ME     RMSE      MAE        MPE     MAPE     MASE
## Training set  8.121418e-15 44.17630 35.91135 -0.9510944 7.995509 2.444228
## Test set     -1.718344e+01 38.01454 33.77760 -4.7345524 8.169955 2.298999
##                     ACF1 Theil's U
## Training set -0.12566970        NA
## Test set     -0.08286364 0.7901651
# Наивный прогноз
accuracy(beer_n, ausbeer)
##                       ME     RMSE      MAE         MPE     MAPE     MASE
## Training set   0.7090909 66.60207 55.43636  -0.8987351 12.26632 3.773156
## Test set     -62.2727273 70.90647 63.90909 -15.5431822 15.87645 4.349833
##                     ACF1 Theil's U
## Training set -0.25475212        NA
## Test set     -0.08286364  1.428524
# Сезонный наивный прогноз
accuracy(beer_sn, ausbeer)
##                     ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -1.846154 17.24261 14.69231 -0.4803931 3.401224 1.0000000
## Test set     -2.545455 12.96849 11.27273 -0.7530978 2.729847 0.7672537
##                    ACF1 Theil's U
## Training set -0.3408329        NA
## Test set     -0.1786912   0.22573
# Метод дрейфа
accuracy(beer_d, ausbeer)
##                         ME     RMSE      MAE        MPE     MAPE     MASE
## Training set  1.240308e-14 66.59830 55.50083  -1.062644 12.29043 3.777543
## Test set     -6.652727e+01 74.83196 67.64793 -16.567964 16.79620 4.604310
##                     ACF1 Theil's U
## Training set -0.25475212        NA
## Test set     -0.07101826  1.509823

При использовании тестового периода функция вычисляет еще один показатель - U-статистику Тейла (Theil’s U).

Этот показатель похож на MASE, и в нем используется для нормализации показатель ошибки наивного прогноза. Однако, в отличие от MASE, используются не абсолютные а процентные ошибки. Поэтому данный показатель страдает от тех же недостатков при малых или нулевых значениях ряда, что и показатели MAPE и MPE.

Суммы квадратов процентных ошибок (Sum of Squared Percent Errors, SSPE) для модели и наивного прогноза вычисляются по формулам:

\[ SSPE(модель) = \sum_{t=2}^N \left(\frac{e_t} {y_{t-1}} \right)^2 \]

\[ SSPE(\text{наивный прогноз}) = \sum_{t=1}^N \left( \frac{y_t - y_{t-1}} {y_{t-1}} \right)^2 \]

\[ U= \sqrt{ \frac{SSPE(модель)}{SSPE(\text{наивный прогноз}) }} \]

Интерпретация U-статистики Тейла аналогична интерпретации показателя MASE:

  • при \(U < 1\) используемый метод работает лучше, чем наивный прогноз;
  • при \(U >= 1\) не имеет смысла использовать этот метод, т.к. он не лучше, чем наивный прогноз.

Как сравнить модели по ошибке в тестовом периоде

Чтобы сравнить модели на основе ошибки в тестовом периоде, можно собрать все показатели ошибки в одну таблицу, с которой можно работать, используя функции пакета dplyr:

rbind(accuracy(beer_m, ausbeer)['Test set',],
      accuracy(beer_n, ausbeer)['Test set',],
      accuracy(beer_sn, ausbeer)['Test set',],
      accuracy(beer_d, ausbeer)['Test set',]) %>%
  as_tibble() %>%
  round(2) %>%
  mutate(`Метод` = c('Среднее', 'Наивный', 'С. наивный', 'Дрейф')) %>%
  select(`Метод`, MASE, MAPE, everything()) %>%
  arrange(MASE)
## # A tibble: 4 x 9
##   Метод       MASE  MAPE     ME  RMSE   MAE    MPE  ACF1 `Theil's U`
##   <chr>      <dbl> <dbl>  <dbl> <dbl> <dbl>  <dbl> <dbl>       <dbl>
## 1 С. наивный  0.77  2.73  -2.55  13.0  11.3  -0.75 -0.18        0.23
## 2 Среднее     2.3   8.17 -17.2   38.0  33.8  -4.73 -0.08        0.79
## 3 Наивный     4.35 15.9  -62.3   70.9  63.9 -15.5  -0.08        1.43
## 4 Дрейф       4.6  16.8  -66.5   74.8  67.6 -16.6  -0.07        1.51

Метод скользящего контроля (cross-validation)

Более сложный и полезный метод оценки качества прогноза основан на скользящем контроле (cross-validation, или rolling forecasting origin). При этом в качестве обучающего периода используются \(n\) начальных наблюдений, а в качестве тестового - одно наблюдение, отстоящее от обучающего периода на величину горизонта прогнозирования (\(n + h\)). Значение \(n\) постепенно увеличивается, насколько позволяют данные. Ошибки прогноза на \(h\) шагов вперед усредняются. Принцип работы метода показан на рисунке:

Скользящий контроль для оценки ошибки прогноза

В это примере горизонт прогноза \(h = 4\). Синие точки - обучающий период, красные - тестовый.

В R скользящий контроль реализован в функции forecast::tsCV(). Эта функция возвращает ошибки прогноза на горизонте \(h\) для каждого периода.

beer3 <- ausbeer %>% window(start = 2003)
tsCV(beer3, meanf, h = 4)
##               h=1        h=2        h=3        h=4
## 2003 Q1 -55.00000 -14.000000  55.000000   0.000000
## 2003 Q2  13.50000  82.500000  27.500000 -17.500000
## 2003 Q3  78.00000  23.000000 -22.000000   0.000000
## 2003 Q4   3.50000 -41.500000 -19.500000  22.500000
## 2004 Q1 -42.20000 -20.200000  21.800000 -16.200000
## 2004 Q2 -13.16667  28.833333  -9.166667 -22.166667
## 2004 Q3  30.71429  -7.285714 -20.285714 -15.285714
## 2004 Q4 -11.12500 -24.125000 -19.125000  54.875000
## 2005 Q1 -22.88889 -17.888889  56.111111  12.111111
## 2005 Q2 -15.60000  58.400000  14.400000 -37.600000
## 2005 Q3  59.81818  15.818182 -36.181818 -17.181818
## 2005 Q4  10.83333 -41.166667 -22.166667  63.833333
## 2006 Q1 -42.00000 -23.000000  63.000000  -1.000000
## 2006 Q2 -20.00000  66.000000   2.000000 -42.000000
## 2006 Q3  67.33333   3.333333 -40.666667 -29.666667
## 2006 Q4  -0.87500 -44.875000 -33.875000  45.125000
## 2007 Q1 -44.82353 -33.823529  45.176471  -7.823529
## 2007 Q2 -31.33333  47.666667  -5.333333 -35.333333
## 2007 Q3  49.31579  -3.684211 -33.684211 -13.684211
## 2007 Q4  -6.15000 -36.150000 -16.150000         NA
## 2008 Q1 -35.85714 -15.857143         NA         NA
## 2008 Q2 -14.22727         NA         NA         NA
## 2008 Q3        NA         NA         NA         NA

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

RMSE <- function(errors) {
  errors^2 %>% mean(na.rm = T) %>% sqrt()
}

# Усреднение
beer3 %>% tsCV(meanf, h = 4) %>% RMSE() %>% round(1)
## [1] 34.5
# Наивный прогноз
beer3 %>% tsCV(naive, h = 4) %>% RMSE() %>% round(1)
## [1] 48.7
# Сезонный наивный прогноз
beer3 %>% tsCV(snaive, h = 4) %>% RMSE() %>% round(1)
## [1] 23.1
# Метод дрейфа
beer3 %>% tsCV(rwf, drift = TRUE, h = 4) %>% RMSE() %>% round(1)
## [1] 71.2

Исследуем, как зависит стандартная ошибка прогноза от горизонта прогнозирования.

beer_err <- tibble(horizon = 1:12) %>%
  mutate(error = map_dbl(horizon, 
                         function(x) tsCV(beer3, naive, h = x) %>% RMSE()))

ggplot(beer_err, aes(horizon, error)) +
  geom_line() +
  labs(title = 'RMSE прогноза производства пива по наивной модели',
       x = 'h', y = beer_unit) +
  scale_x_continuous(breaks = beer_err$horizon)