Наука о данных в R

Сколько рекламы вы смотрели сегодня?

Полный анализ временных рядов рекламных данных

Вы когда-нибудь смотрели внутриигровую рекламу, чтобы отчаянно собирать монеты и бриллианты или возродиться на уровне и попробовать еще раз? Если это так, вы не одиноки. Некоторые объявления хороши, я имею в виду, что они должны продвигать свой продукт и при этом добиваться успеха, но другие просто раздражают. Если вы достаточно смотрели YouTube, вы понимаете, о чем я говорю.

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

Библиотека

Мы будем использовать R на протяжении всей этой статьи. Во-первых, импортируйте необходимые библиотеки.

library(dplyr)          # data wrangling
library(lubridate)      # date functions
library(xts)            # time series object
library(padr)           # padding time series
library(forecast)       # forecasting time series
library(tseries)        # where is pewdiepie
library(ggplot2)        # graphing

Набор данных

Давайте прочитаем набор данных. Этот набор данных доступен на Kaggle.

ads <- read.csv('ads.csv')
head(ads)
#>                  Time    Ads
#> 1 2017-09-13T00:00:00  80115
#> 2 2017-09-13T01:00:00  79885
#> 3 2017-09-13T02:00:00  89325
#> 4 2017-09-13T03:00:00 101930
#> 5 2017-09-13T04:00:00 121630
#> 6 2017-09-13T05:00:00 116475
range(ads$Time)
#> [1] "2017-09-13T00:00:00" "2017-09-21T23:00:00"
range(ads$Ads)
#> [1]  70335 169900

У нас есть почасовые данные о том, сколько объявлений было просмотрено в мобильная игра, содержащая от 70 000 до 170 000 объявлений.

Очистка данных

Начнем с определения.

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

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

ads <- ads %>% 
  mutate(Time = ymd_hms(Time)) %>% 
  arrange(Time) %>% 
  pad()

ads %>% is.na() %>% colSums()
#> Time  Ads 
#>    0    0

Здорово! Отсутствующих значений не существует.

Теперь нам нужно создать объект временного ряда из данных. Поскольку наши данные представлены в часовых интервалах, мы можем сделать это с помощью функции xts(). Упорядочить данные по time_index, что является просто почасовым индексом, таким как Time функция в исходных данных. Кроме того, нам нужно установить частоту на 24, поскольку в каждых 24 метках времени данные имеют периодический характер.

time_index <- seq(
  from = as.POSIXct("2017-09-13 00:00:00"), 
  to = as.POSIXct("2017-09-21 23:00:00"), 
  by = "hour"
)

ads_xts <- xts(ads$Ads, order.by = time_index, frequency = 24)
ads_xts %>% autoplot() + ggtitle("Ads Watched per Hour")

Поскольку объекты временного ряда, происходящие из xts(), не могут быть легко разложены, с этого момента мы будем преобразовывать его в другой объект временного ряда, используя функцию ts(). В результате мы можем применять множество встроенных функций (таких как decompose()), хотя недостатком является то, что часовые интервалы временного индекса будут игнорироваться. Автор был бы очень признателен, если бы существовало более простое решение, напрямую использующее временные ряды, полученные из xts(). Дайте мне знать в разделе ответов ниже!

ads_ts <- ts(ads_xts, frequency = 24)

Разложение

Временной ряд представляет собой композицию из трех основных временных рядов:

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

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

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

dim(ads_ts) <- NULL
ads_ts %>% decompose() %>% autoplot()

Очевидно, что у нас довольно плоский тренд (без тренда). Тем не менее, 16–17 сентября 2017 года на кривой тренда наблюдается небольшой скачок. Это связано с тем, что эти два дня являются выходными, а люди, похоже, больше играют в мобильные игры по выходным, поэтому рекламу смотрят больше. Обратите внимание, что мы потеряли часть информации целых две половины дня в начале и в конце кривой тренда из-за реализации центральной скользящей средней при ее расчете.

Также понятно, что у нас 24-часовая сезонность. Давайте немного увеличим масштаб.

window(
  x = ads_xts, 
  start = as.POSIXct("2017-09-13 00:00:00"), 
  end = as.POSIXct("2017-09-13 23:00:00")
) %>% autoplot() + ggtitle("Ads Watched per Hour on 13 September 2017")

На приведенном выше графике представлены данные нашего временного ряда, особенно за полный день 13 сентября 2017 года (представитель для целей анализа). Количество просмотров рекламы ночью составляет примерно 80 тысяч. Это число увеличилось и достигло своего пика около 120 тысяч в 04:00. Возможно, некоторые утренние геймеры достали свой телефон и хотели поймать какие-то внутриигровые события? Затем число уменьшилось во время утренних занятий и увеличилось до своего второго пика около 160 тысяч в 12:00 во время обеда. После этого кажется, что некоторые люди закончили день и продолжают играть, в то время как другие вернулись к работе, что указывает на небольшое снижение до 16:00. Число снова увеличилось, достигнув своего третьего и последнего пика за день 160-170 тысяч в 18:00–19:00 (после рабочего дня), а затем резко упало примерно до 80 тысяч ночью.

Наконец, остальные наши данные довольно малы (что хорошо), за исключением некоторых всплесков 16–17 сентября 2017 года. Опять же, это может быть связано с нарушениями в выходные дни.

Метрики и проверка

Мы будем использовать MAPE (средняя абсолютная процентная ошибка) в качестве нашей метрики, поскольку она не зависит от величины данных, а только от процента к исходным данным, что делает удобным сравнение производительности различных моделей. Кроме того, мы разделим данные на train и test, как обычно, для проверки модели. Тестовыми данными будут последние 24 часа, которые нам нужно спрогнозировать, как того требовало условие задачи.

test_num <- 24
test <- tail(ads_ts, test_num)
train <- head(ads_ts, length(ads_ts)-length(test))

Моделирование

У нас есть несколько простых моделей на выбор, а также данные временных рядов, для которых они подходят:

  1. Простая скользящая средняя (SMA) → без тренда и сезонность
  2. Экспоненциальное сглаживание
  • Простое экспоненциальное сглаживание (SES) → отсутствие тренда и сезонность
  • Двойное экспоненциальное сглаживание (Экспоненциальное сглаживание Холта) → с трендом без сезонности
  • Тройное экспоненциальное сглаживание (Холт-Винтерс) → с трендом и сезонностью

3. Авторегрессионное интегрированное скользящее среднее (ARIMA) → с трендом без сезонности

4. Seasonal ARIMA (САРИМА) → с трендом и сезонностью

Поскольку у нас есть сезонные данные, мы будем использовать Holt-Winters и SARIMA.

Экспонента Холта-Уинтерса

Сначала постройте модель. Модель Холта-Винтерса автоматически дает нам параметры сглаживания α = 0,92, β = 0 и γ = 1 для этой задачи. Мы можем настроить эти параметры вручную, но пока будем придерживаться того, что дано.

hw_model <- HoltWinters(x = train)
hw_model
#> Holt-Winters exponential smoothing with trend and additive seasonal component.
#> 
#> Call:
#> HoltWinters(x = train)
#> 
#> Smoothing parameters:
#>  alpha: 0.9177049
#>  beta : 0
#>  gamma: 1
#> 
#> Coefficients:
#>            [,1]
#> a   116134.6370
#> b     -112.5410
#> s1  -42551.9149
#> s2  -46114.6303
#> s3  -38622.4558
#> s4  -23062.7216
#> s5     766.4552
#> s6   -2643.0073
#> s7  -10590.9356
#> s8  -16361.4525
#> s9  -18599.8027
#> s10  -9730.1904
#> s11   6773.7268
#> s12  23255.7037
#> s13  32034.0766
#> s14  29038.2208
#> s15  27572.8383
#> s16  26323.7715
#> s17  24075.6307
#> s18  26062.1503
#> s19  36996.4620
#> s20  39800.2418
#> s21  16732.8980
#> s22 -12780.5149
#> s23 -27289.9691
#> s24 -40084.6370

Теперь время прогноза…

hw_forecast <- forecast(hw_model, h = test_num)
hw_forecast
#>          Point Forecast     Lo 80     Hi 80     Lo 95     Hi 95
#> 9.000000       73470.18  64566.94  82373.42  59853.85  87086.51
#> 9.041667       69794.92  57710.82  81879.03  51313.88  88275.97
#> 9.083333       77174.56  62587.46  91761.66  54865.52  99483.60
#> 9.125000       92621.75  75902.26 109341.24  67051.50 118192.00
#> 9.166667      116338.39  97729.27 134947.50  87878.20 144798.57
#> 9.208333      112816.38  92492.58 133140.19  81733.81 143898.96
#> 9.250000      104755.91  82851.24 126660.59  71255.60 138256.23
#> 9.291667       98872.86  75493.96 122251.75  63117.92 134627.79
#> 9.333333       96521.97  71756.45 121287.48  58646.38 134397.55
#> 9.375000      105279.04  79200.53 131357.54  65395.40 145162.67
#> 9.416667      121670.41  94341.92 148998.90  79875.09 163465.73
#> 9.458333      138039.85 109516.10 166563.60  94416.54 181663.16
#> 9.500000      146705.68 117034.78 176376.58 101327.96 192083.40
#> 9.541667      143597.28 112821.97 174372.60  96530.50 190664.06
#> 9.583333      142019.36 110177.91 173860.81  93322.07 190716.65
#> 9.625000      140657.75 107784.73 173530.78  90382.80 190932.70
#> 9.666667      138297.07 104423.87 172170.27  86492.48 190101.66
#> 9.708333      140171.05 105326.37 175015.72  86880.72 193461.38
#> 9.750000      150992.82 115203.03 186782.61  96257.06 205728.58
#> 9.791667      153684.06 116973.47 190394.64  97540.06 209828.05
#> 9.833333      130504.17  92895.33 168113.01  72986.41 188021.93
#> 9.875000      100878.22  62392.08 139364.36  42018.75 159737.69
#> 9.916667       86256.22  46912.35 125600.10  26084.96 146427.49
#> 9.958333       73349.02  33165.70 113532.33  11893.94 134804.09

…и сюжет.

train %>% 
  autoplot(series = "train") +
  autolayer(test, series = "test") +
  autolayer(hw_model$fitted[,1], series = "fitted") +
  autolayer(hw_forecast$mean, series = "forecast") + 
  ggtitle("Holt-Winters Prediction")

hw_mape <- cbind(
  'train_mape' = accuracy(hw_model$fitted[,1], train)[, 'MAPE'],
  'test_mape' = accuracy(hw_forecast$mean, test)[, 'MAPE']
)

hw_mape <- as.data.frame(hw_mape)
rownames(hw_mape) <- c('Holt-Winters')
hw_mape
#>              train_mape test_mape
#> Holt-Winters   3.800135  3.691772

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

САРИМА

Использование SARIMA связано со сложностью: мы должны определить множество параметров и построить модель на основе выбранных параметров, предполагая, что данные стационарны. Чтобы подтвердить это стационарное предположение, мы можем использовать так называемый расширенный тест гипотез Дики-Фуллера (или, для краткости, тест ADF).

  • H₀: данные не являются стационарными
  • H₁: данные стационарны

Мы хотим, чтобы значение p ‹ альфа (0,05), чтобы H₀ отбрасывался, а данные были стационарными.

adf.test(train)
#> 
#>  Augmented Dickey-Fuller Test
#> 
#> data:  train
#> Dickey-Fuller = -4.9338, Lag order = 5, p-value = 0.01
#> alternative hypothesis: stationary
train %>% 
  tsdisplay()

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

train %>% 
  diff(lag = 24) %>%
  adf.test()
#> 
#>  Augmented Dickey-Fuller Test
#> 
#> data:  .
#> Dickey-Fuller = -2.4897, Lag order = 5, p-value = 0.3719
#> alternative hypothesis: stationary
train %>% 
  diff(lag = 24) %>%
  tsdisplay()

После исчезновения сезонности данные становятся нестационарными (p-значение > 0,05). Кроме того, график ACF по-прежнему показывает много значительных лагов. Чтобы преодолеть эту проблему, мы возьмем первое отличие, вычитая ряд из самого себя с запаздыванием 1.

train %>% 
  diff(lag = 24) %>%
  diff() %>%
  adf.test()
#> 
#>  Augmented Dickey-Fuller Test
#> 
#> data:  .
#> Dickey-Fuller = -6.4821, Lag order = 5, p-value = 0.01
#> alternative hypothesis: stationary
train %>% 
  diff(lag = 24) %>%
  diff() %>%
  tsdisplay()

Теперь мы в деле! Определяемся с параметрами.

  1. Поскольку мы делаем сезонную разницу один раз и первую разницу, тогда d = 1 и D = 1.
  2. Из графиков PACF и ACF мы видим, что задержки 3 и 4 значительны, поэтому p = 3 или 4 и q = 3 или 4.
  3. На графике PACF лаги 24 и 48 значительны, поэтому P = 1 или 2.
  4. Судя по графику ACF, задержка 24 значительна, а задержка 48 — нет, поэтому Q = 0 или 1.

Обобщить,

p = 3 or 4

d = 1

q = 3 or 4

P = 1 or 2

D = 1

Q = 0 or 1

Исходя из этих параметров, мы изготовим 4 модели SARIMA следующим образом.

sarima1_model <- Arima(y = train, order = c(4,1,4), seasonal = c(2,1,1))
sarima2_model <- Arima(y = train, order = c(4,1,4), seasonal = c(1,1,1))
sarima3_model <- Arima(y = train, order = c(3,1,3), seasonal = c(1,1,1))
sarima4_model <- Arima(y = train, order = c(3,1,3), seasonal = c(1,1,0))

Взгляните на AIC каждой модели. Помните, что AIC — это мера потери информации, что делает ее sarima3_model лучшей на данный момент. Однако более низкий AIC не означает, что модель работает лучше на test наборе данных.

sarima1_model$aic
#> [1] 3383.505
sarima2_model$aic
#> [1] 3382.48
sarima3_model$aic
#> [1] 3379.964
sarima4_model$aic
#> [1] 3410.026

Как и прежде, время прогноза…

sarima1_forecast <- forecast(sarima1_model, h = test_num)
sarima2_forecast <- forecast(sarima2_model, h = test_num)
sarima3_forecast <- forecast(sarima3_model, h = test_num)
sarima4_forecast <- forecast(sarima4_model, h = test_num)

…и сюжет.

train %>% 
  autoplot(series = "train") +
  autolayer(test, series = "test") +
  autolayer(sarima1_model$fitted, series = "fitted") +
  autolayer(sarima1_forecast$mean, series = "forecast") +
  ggtitle("ARIMA(4,1,4)(2,1,1)[24] Prediction")

train %>% 
  autoplot(series = "train") +
  autolayer(test, series = "test") +
  autolayer(sarima2_model$fitted, series = "fitted") +
  autolayer(sarima2_forecast$mean, series = "forecast") +
  ggtitle("ARIMA(4,1,4)(1,1,1)[24] Prediction")

train %>% 
  autoplot(series = "train") +
  autolayer(test, series = "test") +
  autolayer(sarima3_model$fitted, series = "fitted") +
  autolayer(sarima3_forecast$mean, series = "forecast") +
  ggtitle("ARIMA(3,1,3)(1,1,1)[24] Prediction")

train %>% 
  autoplot(series = "train") +
  autolayer(test, series = "test") +
  autolayer(sarima4_model$fitted, series = "fitted") +
  autolayer(sarima4_forecast$mean, series = "forecast") +
  ggtitle("ARIMA(3,1,3)(1,1,0)[24] Prediction")

Между всеми четырьмя моделями есть лишь небольшие различия. Точно так же, как Холт-Винтерс, модели с трудом предсказывают рабочие часы после обеда, даже не могут предсказать второй пик в 12:00. Мы также можем видеть, что sarima4_model лучше всего прогнозирует утренние часы.

train_mape <- rbind(
  'ARIMA(4,1,4)(2,1,1)[24]' = accuracy(sarima1_model)[, 'MAPE'],
  'ARIMA(4,1,4)(1,1,1)[24]' = accuracy(sarima2_model)[, 'MAPE'],
  'ARIMA(3,1,3)(1,1,1)[24]' = accuracy(sarima3_model)[, 'MAPE'],
  'ARIMA(3,1,3)(1,1,0)[24]' = accuracy(sarima4_model)[, 'MAPE']
)

test_mape <- rbind(
  accuracy(sarima1_forecast$mean, test)[, 'MAPE'],
  accuracy(sarima2_forecast$mean, test)[, 'MAPE'],
  accuracy(sarima3_forecast$mean, test)[, 'MAPE'],
  accuracy(sarima4_forecast$mean, test)[, 'MAPE']
)

sarima_mape <- cbind(train_mape, test_mape)
sarima_mape <- as.data.frame(sarima_mape)
colnames(sarima_mape) <- c('train_mape', 'test_mape')
result <- rbind(hw_mape, sarima_mape)
result
#>                         train_mape test_mape
#> Holt-Winters              3.800135  3.691772
#> ARIMA(4,1,4)(2,1,1)[24]   2.628641  4.185137
#> ARIMA(4,1,4)(1,1,1)[24]   2.680071  4.013859
#> ARIMA(3,1,3)(1,1,1)[24]   2.711371  3.995025
#> ARIMA(3,1,3)(1,1,0)[24]   3.341027  3.758555

Первые три модели SARIMA в некоторой степени указывают на переоснащение. С другой стороны, последняя модель определенно лучше, чем Холт-Винтерс на наборе данных train, но немного хуже, чем Холт-Винтерс на наборе данных test.

Авто АРИМА

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

auto_arima_model <- auto.arima(train)
auto_arima_model
#> Series: train 
#> ARIMA(3,0,3)(1,1,0)[24] 
#> 
#> Coefficients:
#>          ar1      ar2     ar3      ma1     ma2      ma3     sar1
#>       1.5185  -1.0912  0.4706  -0.5546  0.6366  -0.3294  -0.3786
#> s.e.  0.2028   0.2050  0.1065   0.2012  0.1330   0.1084   0.0756
#> 
#> sigma^2 estimated as 38627419:  log likelihood=-1705.08
#> AIC=3426.17   AICc=3427.07   BIC=3451.16

Мы получили ARIMA(3,0,3)(1,1,0)[24] модель, которая очень близка к нашей лучшей модели SARIMA из предыдущих, с расхождением только по параметру d. Эта новая модель не учитывает разницу в 1 лаг, как предполагалось ранее. AIC также является худшим из всех встречавшихся до сих пор моделей. Однако продолжим и посмотрим, что получится.

auto_arima_forecast <- forecast(auto_arima_model, h = test_num)
train %>% 
  autoplot(series = "train") +
  autolayer(test, series = "test") +
  autolayer(auto_arima_model$fitted, series = "fitted") +
  autolayer(auto_arima_forecast$mean, series = "forecast") + 
  ggtitle("Auto ARIMA Prediction")

auto_arima_mape <- cbind(
  'train_mape' = accuracy(auto_arima_model)[, 'MAPE'],
  'test_mape' = accuracy(auto_arima_forecast$mean, test)[, 'MAPE']
)

auto_arima_mape <- as.data.frame(auto_arima_mape)
rownames(auto_arima_mape) <- c('Auto ARIMA')
result <- rbind(result, auto_arima_mape)
result
#>                         train_mape test_mape
#> Holt-Winters              3.800135  3.691772
#> ARIMA(4,1,4)(2,1,1)[24]   2.628641  4.185137
#> ARIMA(4,1,4)(1,1,1)[24]   2.680071  4.013859
#> ARIMA(3,1,3)(1,1,1)[24]   2.711371  3.995025
#> ARIMA(3,1,3)(1,1,0)[24]   3.341027  3.758555
#> Auto ARIMA                3.332454  2.890477

Мы нашли победителя! Не только хорошо работает train с набором данных, но и Auto ARIMA действительно хорошо работает test с набором данных. Это лучшая модель, которую мы создали до сих пор.

Предположения

Мы будем проверять только предположения лучшей модели (Auto ARIMA). Предположения о данных временных рядов проверяются, чтобы определить, достаточно ли хорош остаток, полученный в результате моделирования, чтобы объяснить, что тренд и сезонные композиции уже охватывают как можно больше информации. Наилучшей моделью является модель без автокоррелированных остатков, а остатки имеют нормальное распределение.

Неавтокоррелированные остатки

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

auto_arima_model$residuals %>% 
  tsdisplay()

Как мы видим, на графиках ACF и PACF не осталось значительных лагов. Подождите, на самом деле они есть: отставание 42 на графике ACF и отставание 27 на графике PACF. Однако эти лаги слишком далеки от настоящего и, следовательно, их значение ничтожно.

Если сложно проверить предположение по сюжету, есть более объективный способ сделать это: тест Бокса-Льюнга.

  • H0: остатки не имеют автокорреляции
  • H1: остатки имеют автокорреляцию

Нам нужно значение p > альфа (0,05), чтобы H0 нельзя было отклонить, а остатки не имели автокорреляции.

Box.test(auto_arima_model$residuals, type = "Ljung-Box")
#> 
#>  Box-Ljung test
#> 
#> data:  auto_arima_model$residuals
#> X-squared = 0.044539, df = 1, p-value = 0.8329

Как мы видим, значение p > 0,05, что подтверждает предыдущее наблюдение по графикам ACF и PACF.

Нормально распределенные остатки

Чтобы проверить это, мы можем наблюдать по гистограмме остатков. Распределяются ли они нормально вокруг нуля?

hist(auto_arima_model$residuals)

Если трудно сказать, тест Шапиро-Уилка всегда будет нашим другом.

  • H0: остатки нормально распределены
  • H1: остатки не имеют нормального распределения

Нам нужно значение p > альфа (0,05), чтобы H0 нельзя было отклонить, а остатки распределялись нормально.

shapiro.test(auto_arima_model$residuals)
#> 
#>  Shapiro-Wilk normality test
#> 
#> data:  auto_arima_model$residuals
#> W = 0.94388, p-value = 0.0000008032

Поскольку p-значение ниже альфа (0,05), отклоните H0. Следовательно, остатки не имеют нормального распределения.

Вывод

result
#>                         train_mape test_mape
#> Holt-Winters              3.800135  3.691772
#> ARIMA(4,1,4)(2,1,1)[24]   2.628641  4.185137
#> ARIMA(4,1,4)(1,1,1)[24]   2.680071  4.013859
#> ARIMA(3,1,3)(1,1,1)[24]   2.711371  3.995025
#> ARIMA(3,1,3)(1,1,0)[24]   3.341027  3.758555
#> Auto ARIMA                3.332454  2.890477

Мы реализовали прогнозирование временных рядов для почасового набора данных рекламы мобильных игр. Поскольку набор данных имеет сезонность, мы используем модели Холта-Винтерса и SARIMA. Мы видим, что Auto ARIMA с параметрами ARIMA(3,0,3)(1,1,0)[24] является лучшей моделью с 3,33% MAPE в обучающем наборе данных и 2,89% MAPE в тестовом наборе данных. Однако эта модель нарушает предположение о нормально распределенных остатках данных временных рядов.

Сноска

Привет! Спасибо, что дошли до конца. Это пятая моя статья из серии Наука о данных в R. Прошу найти предыдущую здесь для любознательных, а следующую здесь:



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