Временные ряды в R : измерение ошибки прогноза

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

Рассмотрим часто используемые показатели производительности.
Введем следующие обозначения:
$y_{t}$ - представляет фактическое значение в момент времени t.
$f_{t}$ - представляет собой прогнозируемое значение (расчетное значение) в момент времени t
$e_{t}=y_{t}-f_{t}$ - представляет собой ошибку прогноза (остаточную) в момент времени t
n - количество наблюдений

Средняя ошибка прогноза (mean forecast error - MFE)

$$MFE=\frac{1}{n} \sum_{1}^{n}e_{t}$$
- среднее значение всех ошибок прогноза, величина и знак показывает смещение прогноза относительно фактических значений, если равно нулю, смещения нет
Зависит от масштаба измерений и преобразования данных. 
Не наказывает экстремальные значения ошибок.

Средняя абсолютная ошибка прогноза ( mean absolute error - MAE)

$$MAE=\frac{1}{n}\sum_{1}^{n}\left | e_{t} \right |$$
среднее значение всех ошибок прогноза по абсолютной величине, показывает на величину общей погрешности
Желательны небольшие значения MAE.
Зависит от масштаба измерений и преобразования данных.

Не наказывает экстремальные значения ошибок.

Средняя процентная ошибка (mean percentage error  - MPE)

$$MPE=\frac{1}{n} \sum_{1}^{n}\left ( \frac{e_{t}}{y_{t}} \right )*100$$
средний процент ошибок, величина и знак показывает смещение прогноза относительно фактических значений в процентах

Желательно, чтобы MPE был близок к нулю.

Средняя абсолютная процентная ошибка (mean absolute percentage error - MAPE)

$$MAPE=\frac{1}{n} \sum_{1}^{n}\left | \frac{e_{t}}{y_{t}} \right |*100$$
средняя абсолютная ошибка в процентах без знака ошибки.
Зависит от преобразование данных, но не от масштаба измерений

Не наказывает экстремальные значения ошибок.

Среднеквадратичная ошибка (mean squared error - MSE)

$$MSE=\frac{1}{n}\sum_{1}^{n}e_{t}^{2}$$
средняя квадратичная ошибка, подчеркивает большие ошибки за счет возведения каждой ошибки в квадрат.

Зависит от масштаба измерений и преобразования данных.

Сумма квадратов ошибок (sum of squared estimate of errors - SSE)

$$SSE=\sum_{1}^{n}e_{t}^{2}$$

сумма квадратов ошибок прогноза.

Свойства такие же, как и у MSE.


Среднеквадратичная ошибка со знаком (signed mean squared error - SMSE)

$$SMSE=\frac{1}{n} \sum_{1}^{n}\left ( \frac{e_{t}}{\left | e_{t} \right |} \right )*e_{t}^{2}$$

тоже, что и MSE только с учетом знака.


U-статистика Тейла

$$U=\frac{\sqrt{\frac{1}{n}\sum_{1}^{n}e_{t}^{2}}}{\sqrt{\frac{1}{n}\sum_{1}^{n}f_{t}^{2}}+\sqrt{\frac{1}{n}\sum_{1}^{n}y_{t}^{2}}}$$

Диапазоны 0 ≤ U ≤ 1
U = 0 означает идеальный прогноз.
Зависит от масштаба измерений и преобразования данных.
Индекс Тейла показывает степень схожести временных рядов  и чем ближе он к нулю, тем ближе сравниваемые ряды.

Корень среднеквадратичной ошибки (root mean squared error - RMSE)

$$RMSE=\sqrt{MSE}=\sqrt{\frac{1}{n}\sum_{1}^{n}e_{t}^{2}}$$

Такие же свойства как MSE

Нормированная среднеквадратичная ошибка (normalized mean square error-NMSE)

$$NMSE=\frac{MSE}{\sigma ^{2}}=\frac{1}{\sigma ^{2}n}\sum_{1}^{n}e_{t}^{2}$$

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


Проверим эти метрики на примере смоделированных прогнозов. Для этого берем датафрейм с месячными показателями посещаемости магазина одежды за пять лет с января 2015 по декабрь 2019 годы.

head(df)

# A tibble: 6 x 3
   year month visit
  <dbl> <dbl> <dbl>
1  2015     1 14328
2  2015     2  6974
3  2015     3 10174
4  2015     4  8994
5  2015     5 10715

6  2015     6  9381

str(df)

Classes ‘tbl_df’, ‘tbl’ and 'data.frame': 60 obs. of  3 variables:
 $ year : num  2015 2015 2015 2015 2015 ...
 $ month: num  1 2 3 4 5 6 7 8 9 10 ...

 $ visit: num  14328 6974 10174 8994 10715 

С имитируем три три прогноза с помощью трех последовательностей мультипликативной случайной составляющей, представляющей собой случайную величину с нормальным распределением с показателями 𝛔=0.95 и sd=0.05 (прогноз с занижением),𝛔=1.05 и sd=0.05 (прогноз с завышением),𝛔=1.0 и sd=0.05 (идеальный прогноз)

df$visit_n <- df$visit*rnorm(60,0.95,0.05)
df$visit_p <- df$visit*rnorm(60,1.05,0.05)
df$visit_z <- df$visit*rnorm(60,1.00,0.05)


round(colSums(df[,c(3:6)]),0)

  visit      visit_n    visit_p visit_z 

 592714  563092  626292  594431 

Как видим, суммарное количество посетителей по полученным временным рядам получилось, как и было задумано. Сделаем из полученного датафрейма объект класса ts и выведем его график.

ts_month <- ts(data = df_month[,c(3:6)],
                start = c(2015, 1),
                frequency = 12)

autoplot(ts_month)

 Можно переходить к расчету ошибок прогноза.

mfe <- df %>%
  summarise(mfe_n = mean(visit_n-visit),
            mfe_p = mean(visit_p-visit),
            mfe_z = mean(visit_z-visit))

mae <- df %>%
  summarise(mae_n = mean(abs(visit_n-visit)),
            mae_p = mean(abs(visit_p-visit)),
            mae_z = mean(abs(visit_z-visit)))

mpe <- df %>%
  summarise(mpe_n = mean((visit_n-visit)/visit*100),
            mpe_p = mean((visit_p-visit)/visit*100),
            mpe_z = mean((visit_z-visit)/visit*100))

mape <- df %>%
  summarise(mape_n = mean(abs(visit_n-visit)/visit*100),
            mape_p = mean(abs(visit_p-visit)/visit*100),
            mape_z = mean(abs(visit_z-visit)/visit*100))
mse <- df %>%
  summarise(mse_n = mean((visit_n-visit)^2),
            mse_p = mean((visit_p-visit)^2),
            mse_z = mean((visit_z-visit)^2))

sse <- df %>%
  summarise(sse_n = sum((visit_n-visit)^2),
            sse_p = sum((visit_p-visit)^2),
            sse_z = sum((visit_z-visit)^2))

smse <- df %>%
  summarise(smse_n = mean(sign(visit_n-visit)*(visit_n-visit)^2),
            smse_p = mean(sign(visit_p-visit)*(visit_p-visit)^2),
            smse_z = mean(sign(visit_z-visit)*(visit_z-visit)^2))

ut <- df %>%
  summarise(ut_n = sqrt(mean((visit_n-visit)^2))/(sqrt(mean(visit_n^2))+sqrt(mean(visit^2))),
            ut_p = sqrt(mean((visit_p-visit)^2))/(sqrt(mean(visit_p^2))+sqrt(mean(visit^2))),
            ut_z = sqrt(mean((visit_z-visit)^2))/(sqrt(mean(visit_z^2))+sqrt(mean(visit^2))))

rmse <- df %>%
  summarise(rmse_n = sqrt(mean((visit_n-visit)^2)),
            rmse_p = sqrt(mean((visit_p-visit)^2)),
            rmse_z = sqrt(mean((visit_z-visit)^2)))

nmse <- df %>%
  summarise(nmse_n = mean((visit_n-visit)^2)/var(visit_n),
            nmse_p = mean((visit_p-visit)^2)/var(visit_p),
            nmse_z = mean((visit_z-visit)^2)/var(visit_z))

df_err <- data.frame(mfe=c(mfe$mfe_n,mfe$mfe_p,mfe$mfe_z),

                  mae=c(mae$mae_n,mae$mae_p,mae$mae_z),

                  mpe=c(mpe$mpe_n,mpe$mpe_p,mpe$mpe_z),

                  mape=c(mape$mape_n,mape$mape_p,mape$mape_z),

                  mse=c(mse$mse_n,mse$mse_p,mse$mse_z),

                  sse=c(sse$sse_n,sse$sse_p,sse$sse_z),

                  smse=c(smse$smse_n,smse$smse_p,smse$smse_z),

                  ut=c(ut$ut_n,ut$ut_p,ut$ut_z),

                  rmse=c(rmse$rmse_n,rmse$rmse_p,rmse$rmse_z),

                  nmse=c(nmse$nmse_n,nmse$nmse_p,nmse$nmse_z),

                  row.names = c('n','p','z'))

round(df_err,2)

      mfe       mae         mpe    mape      mse            sse            smse          ut      rmse       nmse
n -493.876  538.623   -5.19   5.627   462593.4   27755605   -454568.44   0.034   680.142   0.102
p  660.635  700.996    6.38   6.840   760386.2   45623174    746613.43    0.042   872.001   0.128
z  -24.242   394.060    0.02   3.943   274906.5   16494390    -68096.34     0.026   524.315   0.068


Оценим величину разброса значений оценок через коэффициент вариации

round(sapply(df_err, function(x){sd(x)/mean(x)}),2)

  mfe   mae   mpe   mape   mse   sse   smse    ut     rmse  nmse 

12.22  0.28   14.35   0.27     0.49   0.49   8.21    0.23   0.25    0.30 

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