Наука о данных в 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)
Разложение
Временной ряд представляет собой композицию из трех основных временных рядов:
- тенденция: общее направление данных, независимо от того, является ли оно плоским, возрастающим или убывающим. Если в тренде все еще есть колебания, это означает, что в данных есть какая-то закономерность, не учтенная в декомпозициях.
- сезонный: периодический шаблон, который происходит в пределах фиксированной метки времени.
- остаток: любой паттерн, не учитываемый в тренде и сезонности.
Под композицией мы подразумеваем здесь просто сложение или умножение этих трех разложений. Временной ряд аддитивного типа имеет относительно постоянную дисперсию все время, тогда как мультипликативный аналог имеет тенденцию к увеличению или уменьшению дисперсии. В качестве иллюстрации, пожалуйста, взгляните на картинку ниже.
Мультипликативный временной ряд может быть логарифмически преобразован для создания аддитивного временного ряда. Очевидно, что наши данные уже являются аддитивным временным рядом. Понаблюдаем за его разложениями.
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))
Моделирование
У нас есть несколько простых моделей на выбор, а также данные временных рядов, для которых они подходят:
- Простая скользящая средняя (SMA) → без тренда и сезонность
- Экспоненциальное сглаживание
- Простое экспоненциальное сглаживание (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()
Теперь мы в деле! Определяемся с параметрами.
- Поскольку мы делаем сезонную разницу один раз и первую разницу, тогда d = 1 и D = 1.
- Из графиков PACF и ACF мы видим, что задержки 3 и 4 значительны, поэтому p = 3 или 4 и q = 3 или 4.
- На графике PACF лаги 24 и 48 значительны, поэтому P = 1 или 2.
- Судя по графику 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. Ваш членский взнос напрямую поддерживает меня и других писателей, которых вы читаете.