Прогнозирование временных рядов с помощью ARIMA в Python 3

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

В этом руководстве мы сосредоточимся на использовании модели ARIMA, одного из наиболее часто применяемых подходов в области анализа временных рядов. Мы подробно рассмотрим процесс использования модели в Python 3 — от начальных этапов по загрузке и обработке данных до конечного этапа — прогнозирования. Мы также изучим, как определить и интерпретировать параметры модели ARIMA и как оценить ее качество.

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

cloud

Подготовка рабочей среды для анализа данных на Python 

Установка Python

В первую очередь, вам необходимо установить сам Python — язык программирования, который мы будем использовать для анализа данных. Можно скачать его с официального сайта python.org, следуя представленным там инструкциям по установке.

После завершения установки откройте командную строку (на Windows) или терминал (на MacOS/Linux) и введите:

python --version

Если все было сделано правильно, вы увидите номер установленной версии Python.

Установка программной среды

Для работы с Python вы можете выбрать удобную для вас среду разработки (IDE). В данном руководстве мы будем работать с Jupyter Notebook, также большой популярностью среди аналитиков данных пользуются PyCharm, Visual Studio Code, Spyder. 

Для установки Jupyter Notebook введите в командной строке: 

pip install jupyter

Установка необходимых пакетов Python

Вместе с Python всегда ставятся некоторые базовые библиотеки. Это крайне полезные инструменты, но для более глубокого анализа данных вам могут понадобиться дополнительные инструменты. В этом руководстве мы будем использовать: pandas (для работы с данными в табличной форме), numpy (для научных вычислений), matplotlib (для визуализации данных) и statsmodels (библиотека статистических моделей). 

Вы можете установить эти библиотеки при помощи команды pip3 install в терминале или командой строке:

pip3 install pandas numpy matplotlib statsmodels

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

Для проверки установленных пакетов используйте команду: 

pip list

В результате вы получите список всех установленных модулей и их версии: 

Создание рабочего каталога

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

cd path_to_your_directory
mkdir Your_Project_Name
cd Your_Project_Name

Здесь path_to_your_directory — путь к месту, где будет создана папка проекта, а Your_Project_Name — это название вашего проекта. 

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

Загрузка и обработка данных

Запуск Jupyter Notebook

Начнем с запуска Jupyter Notebook, нашего основного инструмента для написания и тестирования кода Python. В командной строке (или терминале) перейдите в ваш рабочий каталог и введите следующую команду: 

jupyter notebook

В вашем браузере откроется новая вкладка с интерфейсом Jupyter Notebook. Чтобы создать новый документ, выберете вкладку «New» в правом верхнем углу окна и в выпадающем меню выберите «Python 3». Вы будете автоматически перенаправлены в новую вкладку, где и будет создан ваш блокнот.

Импорт библиотек

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

import warnings
import itertools
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm

Для запуска кода используйте Shift+Enter. Теперь эти библиотеки доступны в вашем проекте, и вы можете использовать их функционал для выполнения различных задач анализа данных. 

При нажатии Shift+Enter текущая ячейка кода будет выполнена, а фокус переместится на следующую ячейку, где вы можете продолжать писать свой код. 

Загрузка данных

Для прогнозирования временных рядов в Python мы будем использовать набор данных Airline Passengers. Этот датасет представляет собой месячное отслеживание количества пассажиров в международных авиалиниях, выраженное в тысячах, за период с 1949 по 1960 годы. Эти данные можно найти здесь.

Для загрузки данных из csv-файла по URL используем библиотеку pandas: 

url = "https://raw.githubusercontent.com/jbrownlee/Datasets/master/airline-passengers.csv"
time_series = pd.read_csv(url)

Если же ваш csv-файл хранится локально на вашем компьютере, то для загрузки используйте следующую строку кода:

time_series = pd.read_csv('airline_passengers.csv')

Теперь данные сохранены в DataFrame, обозначенном как time_series. DataFrame — это основная структура данных в библиотеке pandas. Она представляет собой двумерную таблицу, где каждая строка — отдельное наблюдение, а столбцы — различные признаки или переменные этого наблюдения.

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

print(time_series.head())

Этот код выдаст первые пять строк загруженного набора данных, позволяя быстро проверить, правильно ли они были загружены и выглядят ли они так, как ожидалось. По умолчанию метод head() выводит именно пять строк, но для просмотра другого количества, можно указать нужно число в скобках метода. 

Также можно просмотреть последние строки датафрейма:

print(time_series.tail())

Обработка данных

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

1. Проверка пропущенных значений

Обработка отсутствующих значений — важный этап предварительной обработки в анализе динамических рядов. Отсутствующие значения могут вызывать проблемы при анализе и искажать результаты прогнозирования. Для проверки пропущенных значений можно использовать метод isnull() из библиотеки pandas: 

print(time_series.isnull().sum())

Если для всех колонок указан 0, это означает отсутствие пропусков в данных. Однако, если в процессе выполнения будут обнаружены пропущенные значения, их следует обработать.

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

time_series = time_series.fillna(time_series.mean())

Для замены пропущенных значений только в определенных столбцах используйте эту команду:

time_series['Столбец 1'] = time_series['Столбец 1'].fillna(time_series['Столбец 1'].mean())

2. Преобразование типа данных

Каждый столбец в датафрейме имеет определенный тип данных. В динамических рядах особенно важен тип данных DataTime, который специально предназначен для хранения дат и времени. Наш датафрейм в pandas по умолчанию считывает информацию как текст. Даже если в столбце содержатся даты, pandas будет воспринимать их как обычные строки. В нашем случае необходимо преобразовать столбец 'Month' в DataTime, чтобы мы могли работать с временными данными: 

time_series['Month'] = pd.to_datetime(time_series['Month'])

3. Указание DataTime в качестве индекса

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

Это позволит легко выбирать и анализировать данные за определенные периоды времени. 

В нашем случае используем столбец 'Month' в качестве индекса:

time_series.set_index('Month', inplace=True)

4. Изменение масштаба данных

Еще один важный шаг в предварительной обработке данных — это проверка необходимости изменения масштаба данных. Если размах ваших данных слишком велик (например, значение 'Passengers' колеблется от тысяч до миллионов), вам может потребоваться преобразование данных. 

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

Пример стандартизации данных при большом размахе: 

from sklearn.preprocessing import StandardScaler

scaler = StandardScaler()
time_series[['Passengers']] = scaler.fit_transform(time_series[['Passengers']])

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

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

Визуализация данных

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

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

Чтобы нарисовать график, используйте следующие строки кода: 

plt.figure(figsize=(15,8))
plt.plot(time_series['Passengers'])
plt.title('График временного ряда')
plt.xlabel('Дата')
plt.ylabel('Пассажиры')
plt.show()

Для наших данных мы получаем такой график: 

Стационарность данных

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

Модель ARIMA способна адаптировать ряд к стационарности «самостоятельно», за это свойство ответственен специальный параметр модели (d). Но все же понимание того, стационарен ли ваш начальный временной ряд, помогает лучше разобраться в работе ARIMA. 

Существует несколько методов для проверки ряда на стационарность. 

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

• Среднее значение: колебания этого показателя во времени могут сигнализировать о том, что временной ряд не стационарен.

• Дисперсия: если дисперсия меняется со временем, это также свидетельствует о нестационарности ряда.

• Тренд: если на графике есть видимый тренд, это еще одно уведомление о нестационарности.

• Сезонность: если на графике присутствуют сезонные колебания, это может указывать на нестационарность. 

Во-вторых, можно провести статистический анализ, например, воспользоваться тестом Дики-Фуллера. Этот метод представляет количественную оценку стационарности временного ряда. Основная гипотеза теста предполагает, что ряд динамики не стационарен. Если значение p меньше уровня значимости 0.05 в результате теста, то основная гипотеза отвергается, и ряд может считаться стационарным.

Проведение теста Дики-Фуллера для наших данных может выглядеть так: 

from statsmodels.tsa.stattools import adfuller

print('Результат теста:')
df_result = adfuller(time_series['Passengers'])
df_labels = ['ADF Test Statistic', 'p-value', '#Lags Used', 'Number of Observations Used']
for result_value, label in zip(df_result, df_labels):
    print(label + ' : ' + str(result_value))

if df_result[1] <= 0.05:
    print("Сильные доказательства против нулевой гипотезы, ряд является стационарным.")
else:
    print("Слабые доказательства против нулевой гипотезы, ряд не является стационарным.")

Наш временной ряд не является стационарным, но в следующих разделах мы займемся автоматизированным поиском параметров для модели ARIMA и найдем необходимые параметры, которые приведут ряд к стационарности.

Несмотря на то, что нам не надо самостоятельно приводить ряд к стационарности, полезно знать, какими методами можно это сделать. Существует множество способов, некоторые из них:

1. Дифференцирование: один из самых распространенных методов, именно с помощью дифференцирования ARIMA приводит ряд к стационарности. Этот процесс основывается на вычислении разницы между последовательными наблюдениями во временном ряду. 

2. Сезонное дифференцирование: вариация обычного дифференцирования, которая применяется для данных с сезонной компонентой

3. Логарифмирование: создание логарифма данных может помочь снизить вариативность ряда и сделать его более стационарным. 

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

Если вы решили преобразовать ряд, то после всех процессов нужно снова проверить ряд на стационарность, чтобы убедиться, что преобразование было успешным. Это можно сделать, повторно выполнив тест Дики-Фуллера. 

Модель ARIMA

ARIMA (AutoRegressive Integrated Moving Average), Авторегрессионное интегрированное скользящее среднее — это статистическая модель, которая используется для анализа и прогнозирования данных динамических рядов. 

1. Авторегрессионная (AR) часть относится к использованию зависимости между текущим наблюдением и определенным количеством предыдущих наблюдений. Например, чтобы предсказать погоду на завтра, мы смотрим погоду в предыдущие дни. 

2. Интегрированная (I) часть обозначает, что данные временного ряда преобразуются с целью сделать ряд стационарным. Под стационарностью подразумевается такое свойство временного ряда, при котором его статистические характеристики не меняются во времени. Например, среднее и дисперсия остаются постоянными со временем. Обычно это достигается путем применения операции разности к ряду один или несколько раз. 

3. Скользящее среднее (MA) относится к моделированию ошибки между реальным наблюдением и прогнозируемым значением, используя комбинацию значений ошибок предыдущих шагов. 

Модель ARIMA обычно обозначается как ARIMA (p, d, q), где p, d и q — параметры модели: 

• p — порядок авторегрессии, который описывает количество предыдущих наблюдений, принимаемых во внимание в модели. 

• d — порядок интегрирования, который указывает, сколько раз нужно взять разность временного ряда для достижения стационарности. 

• q — порядок скользящего среднего, который связан с количеством предыдущих ошибок, используемых в модели.

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

Сезонная модель ARIMA

Сезонная модель ARIMA или SARIMA (Seasonal Autoregressive Integrated Moving Average) — это расширение модели ARIMA, добавляющее в нее понятие сезонности. В ряде ситуаций динамические ряды отражают явные сезонные колебания. Например, продажи мороженого, вероятно, увеличатся в летние месяцы и уменьшатся в зимние. Это тенденция, которая повторяется каждый год. Здесь и приходит на помощь модель SARIMA — она позволяет улавливать и использовать эти сезонные колебания. 

Сезонная модель ARIMA обычно обозначается как SARIMA (p, d, q) (P, D, Q) m, где p, d, q — это несезонные параметры, а P, D, Q — сезонные:

• p, d, q работают также, как и в ARIMA.

• P — порядок сезонной авторегрессии, и он основан на количестве предыдущих сезонов, которые влияют на текущий сезон.

• D — порядок сезонной интеграции, то есть сколько раз нужно продифференцировать данные, чтобы убрать сезонные тренды.

• Q — порядок сезонного скользящего среднего, то есть как много предыдущих сезонных ошибок прогноза используется в модели.

• m — это длина сезонного периода. Например, если данные имеют годовую сезонность и данные собираются ежемесячно, то m будет равно 12. 

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

Модель ARIMA, в частности сезонная ARIMA, может показаться сложной, поскольку необходимо тщательно подбирать многочисленные параметры. Чтобы упростить эту задачу, необходимо автоматизировать этот процесс. 

Определение параметров модели

Первый шаг в настройке модели ARIMA — это определение оптимальных значений параметров для нашего конкретного набора данных. 

Для настройки параметров ARIMA мы будем использовать «поиск по сетке» (grid search). Суть метода заключается в том, что он проходит через все возможные комбинации параметров из заранее определенной сетки значений и обучает модель на каждой из этих комбинаций.

После обучения модели на каждой из комбинаций выбирается модель с лучшими показателями. 

Чем больше различных значений параметров, тем больше комбинаций придется проверить и тем дольше будет идти процесс. Для нашего случая мы будем использовать только два возможных значения (0 и 1) для каждого из параметров, то есть всего 8 комбинаций для параметров ARIMA и 8 для сезонной части (при длине сезонного периода = 12). Таким образом, общее количество комбинаций, которые подлежат проверке, равно 64, что приведет к относительно быстрому выполнению.

Важно помнить, что цель — найти баланс между временем, затраченным на «поиск по сетке», и качеством исходной модели, то есть найти такие значения параметров, при которых качество будет наибольшим, а затраты времени минимальными. 

Импорт необходимых пакетов

from statsmodels.tsa.arima.model import ARIMA
from statsmodels.tsa.stattools import acf, pacf
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf

Statsmodels предоставляет нам методы для построения моделей ARIMA, а itertools (мы импортировали ранее) используется для создания комбинаций из возможных значений параметров. 

Игнорирование предупреждений

При работе с большими наборами данных и сложными вычислениями, такими как статистический анализ или машинное обучение, библиотеки и функции могут генерировать предупреждения о возможных проблемах или неоптимальности. Однако, эти предупреждения могут быть в большинстве случаев незначимыми или нерелевантными для вашего конкретного кейса. Поэтому мы устанавливаем фильтр предупреждений на «ignore»: 

warnings.filterwarnings("ignore")

Создание диапазона параметров для настройки модели

Для определения параметров модели определим функцию search_optimal_sarima. 

def search_optimal_sarima(time_series, seasonal_cycle):
    order_vals = diff_vals = ma_vals = range(0, 2)
    pdq_combinations = list(itertools.product(order_vals, diff_vals, ma_vals))
    seasonal_combinations = [(combo[0], combo[1], combo[2], seasonal_cycle) for combo in pdq_combinations]
       
    smallest_aic = float("inf")
    optimal_order_param = optimal_seasonal_param = None

    for order_param in pdq_combinations:
        for seasonal_param in seasonal_combinations:
            try:
                sarima_model = sm.tsa.statespace.SARIMAX(time_series,
                                                         order=order_param,
                                                         seasonal_order=seasonal_param,
                                                         enforce_stationarity=False,
                                                         enforce_invertibility=False)

                model_results = sarima_model.fit()
                if model_results.aic < smallest_aic:
                    smallest_aic = model_results.aic
                    optimal_order_param = order_param
                    optimal_seasonal_param = seasonal_param
            except:
                continue

    print('ARIMA{}x{} - AIC:{}'.format(optimal_order_param, optimal_seasonal_param, smallest_aic))

seasonal_cycle_length = 12
search_optimal_sarima(time_series, seasonal_cycle_length)

За создание диапазона параметров в нашей функции отвечают первые три строчки кода. 

Как мы уже знаем, в модели ARIMA есть три основных параметра, p, d, q. В коде выше p, d и q являются диапазонами от 0 до 2, то есть они могут принимать значения 0 или 1. Функция itertools.product() генерирует все возможные комбинации этих трех параметров. 

Примеры комбинаций включают (0, 0, 0), (0, 0, 1), (0, 1, 1) и так далее. 

Затем мы создаем дополнительные комбинации, добавляя сезонный период к каждой из комбинаций pdq. Это позволяет модели учитывать сезонные влияния на ряд динамики. 

Поиск лучших параметров для модели

Теперь нам необходимо применить те параметры, которые уже были определены ранее, для автоматической настройки моделей ARIMA. 

Когда мы работаем с моделями прогнозирования, нашей задачей является выбор той модели, которая лучше всего подходит для объяснения и предсказания данных. Однако, выбор лучшей модели не всегда прост. На помощь приходит критерий информативности AIC (Akaike Information Criterion), который позволяет нам сравнить разные модели и определить, какая из них лучше. AIC помогает оценить, насколько модель соответствует данным, учитывая ее сложность. Итак, цель — найти модель с наименьшим значением AIC.

Описанный выше код проводит итерацию по всевозможным комбинациям параметров и применяет функцию SARIMAX для построения сезонной модели ARIMA. Параметр 'order' устанавливает основные параметры (p, d, q), а 'seasonal_order' задает сезонные параметры модели (P, D, Q, S).

Для наших данных мы получим такой результат:

ARIMA(0, 1, 1)x(1, 1, 1, 12) - AIC:920.3192974989254

Построение и оценка модели SARIMAX

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

Для начала определяем модель SARIMAX с найденными ранее параметрами: 

from statsmodels.tsa.statespace.sarimax import SARIMAX
model = SARIMAX(time_series, order=(0, 1, 1), seasonal_order=(1, 1, 1, 12))

Затем обучаем модель: 

results = model.fit()

После чего можно выводить сводку по модели, которая предоставляет информацию о качестве подбора параметров:

print(results.summary())

Вывод модели широко используется для оценки ее качества. Основные пункты, на которые следуют обратить внимание:

1. Коэффициенты: они должны быть статистически значимыми. Проверьте p-значения коэффициентов (P>|z|); они должны быть меньше 0.05. 

2. Информационный критерий (AIC): меньшее значение AIC указывает на лучшую подгонку модели.

3. Ljung-Box (L1) (Q): это p-значение для Ljung-Box Q-статистики. Если значение больше 0.05, остатки случайны, что хорошо.

4. Jarque-Bera (JB): это тест на нормальность остатков. Если Prob(JB) больше 0.05, остатки нормальны, что хорошо. 

5. Heteroskedasticity (H): это тест на гетероскедастичность остатков. Если Prob(H) (two-sided) больше 0.05, остатки гетероскедастичны, что хорошо. Гетероскедастичность — это ситуация, когда разброс ошибок ваших прогнозов меняется в зависимости от того, в какой точке вы находитесь. Более простыми словами, это когда у вас имеется неоднородность в ваших данных.

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

Для нашей модели мы получили следующий вывод:

Следующий код служит для визуализации диаграммы диагностики модели из statsmodels:

results.plot_diagnostics(figsize=(12, 8))
plt.show()

Берется заполненная ARIMA-модель и создается четыре графика:

1. Кривая остатков — это график остатков модели по времени. Если модель хороша, остатки будут случайными, и график будет выглядеть как белый шум. 

2. График нормального квантиля (Q-Q plot) — это график, который сравнивает распределение остатков с идеальным нормальным распределением. Если точки следуют диагональной линии, это означает, что остатки нормально распределены.

3. График коррелограммы — это график автокорреляции остатков. Если модель хороша, остатки не будут коррелировать друг с другом. Отсутствие голубых полосок за пределами голубого диапазона шума указывает на это. 

4. Гистограмма остатков — это гистограмма распределения остатков. Если модель хороша, остатки должны быть нормально распределены, и эта гистограмма будет напоминать форму колокола. 

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

Наши графики выглядят следующим образом: 

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

Статическое и динамическое прогнозирование

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

Статический прогноз

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

st_pred = results.get_prediction(start=pd.to_datetime('1955-12-01'), dynamic=False) 
forecast_values = st_pred.predicted_mean

Здесь используются данные с декабря 1955 года в качестве примера, вы можете поменять это значение в соответствии со своими потребностями. 

Теперь у нас есть прогнозные значения, которые мы можем сравнить с реальными данными временного ряда. Мы будем использовать среднеквадратичную ошибку (MSE) как нашу метрику оценки точности прогноза: 

actual_values = time_series['1955-12-01':]['Passengers']
forecast_mse = ((forecast_values - actual_values) ** 2).mean()
print('Среднеквадратичная ошибка прогноза составляет {}'.format(round(forecast_mse, 2)))

MSE является общепринятой метрикой для оценки качества моделей прогнозирования. Чем ниже значение MSE, тем точнее модель. Разумеется, нет безупречной модели, и всегда будет некоторое отклонение между прогнозами и реальными данными. В нашем случае среднеквадратичная ошибка прогноза составляет 170.37.

Наконец, визуализируем результаты. Это позволяет нам визуально проверить точность наших прогнозов относительно реальных данных:

plt.figure(figsize=(15,8))

plt.plot(actual_values.index, actual_values, label='Реальные значения', color='blue')

plt.plot(forecast_values.index, forecast_values, label='Спрогнозированные значения', color='red')

plt.title('Реальные и cпрогнозированные значения')
plt.xlabel('Дата')
plt.ylabel('Пассажиры')
plt.legend()

plt.show()

Этот код строит график реальных и прогнозируемых значений количества пассажиров в течение времени. Красная линия обозначает прогнозируемые значения, а синяя линия — реальные данные. 

Эта визуализация может помочь вам понять, насколько хорошо модель предсказывает данные. 

Динамический прогноз

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

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

Для реализации динамического прогноза необходимо изменить параметр dynamic на True:

dyn_pred = results.get_prediction(start=pd.to_datetime('1955-12-01'), dynamic=True) 
dynamic_forecast_values = dyn_pred.predicted_mean

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

mse_dynamic_forecast = ((dynamic_forecast_values - actual_values) ** 2).mean()
print('Среднеквадратичная ошибка динамического прогноза составляет {}'.format(round(mse_dynamic_forecast, 2)))

И построить график реальных и прогнозируемых значений: 

plt.figure(figsize=(15,8))

plt.plot(actual_values.index, actual_values, label='Реальные значения', color='blue')

plt.plot(dynamic_forecast_values.index, dynamic_forecast_values, label='Динамический прогноз', color='green')

plt.title('Реальные и динамически спрогнозированные значения')
plt.xlabel('Дата')
plt.ylabel('Пассажиры')
plt.legend()

plt.show()

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

Создание и визуализация прогноза

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

Чтобы выполнить прогнозирование на определенное количество шагов вперед, можно использовать метод get_forecast модели results:

pred_future = results.get_forecast(steps=12)

Мы используем обученную модель (results) для получения прогноза на следующие 12 периодов времени. Поскольку наши данные содержали информацию до декабря 1960 года, здесь мы получим предсказание о количестве пассажиров каждый месяц на 1961 год. 

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

print(f'Средние прогнозируемые значения:\n\n{pred_future.predicted_mean}')
print(f'\nДоверительные интервалы:\n\n{pred_future.conf_int()}')

Также мы можем визуализировать наш прогноз: 

fig = plt.figure()
plt.plot(pred_future.predicted_mean, label='Средние прогнозируемые значения')
plt.fill_between(pred_future.conf_int().index,
                 pred_future.conf_int().iloc[:, 0],
                 pred_future.conf_int().iloc[:, 1], color='k', alpha=.2)
plt.legend()  
plt.show()

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

Подготовили для вас выгодные тарифы на облачные серверы

Заключение 

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

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

Теперь вы можете автоматизировано прогнозировать динамические ряды при помощи модели ARIMA и языка программирования Python. Советуем потренироваться и проделать это руководство еще раз, но уже с другими данными. 

# Проверка на отсутствующие значения
print(time_series.isnull().sum())

# Визуализация данных
time_series.plot(figsize=(10, 6))
plt.show()

def search_optimal_sarima(time_series, seasonal_cycle):
    order_vals = range(0, 3)  # Увеличен диапазон параметров
    diff_vals = range(0, 3)
    ma_vals = range(0, 3)
    pdq_combinations = list(itertools.product(order_vals, diff_vals, ma_vals))
    seasonal_combinations = [(combo[0], combo[1], combo[2], seasonal_cycle) for combo in pdq_combinations]
       
    smallest_aic = float("inf")
    optimal_order_param = optimal_seasonal_param = None

    for order_param in pdq_combinations:
        for seasonal_param in seasonal_combinations:
            try:
                sarima_model = sm.tsa.statespace.SARIMAX(time_series,
                                                         order=order_param,
                                                         seasonal_order=seasonal_param,
                                                         enforce_stationarity=False,
                                                         enforce_invertibility=False)

                model_results = sarima_model.fit()
                if model_results.aic < smallest_aic:
                    smallest_aic = model_results.aic
                    optimal_order_param = order_param
                    optimal_seasonal_param = seasonal_param
            except Exception as e:
                print(f"Exception for parameters {order_param}x{seasonal_param}: {e}")
                continue

    if optimal_order_param is not None and optimal_seasonal_param is not None:
        print('ARIMA{}x{} - AIC:{}'.format(optimal_order_param, optimal_seasonal_param, smallest_aic))
    else:
        print('No suitable model found.')

# Указание длины сезонного цикла
seasonal_cycle_length = 12

# Поиск оптимальных параметров
search_optimal_sarima(time_series, seasonal_cycle_length)

try:
    sarima_model = sm.tsa.statespace.SARIMAX(time_series,
                                             order=order_param,
                                             seasonal_order=seasonal_param,
                                             enforce_stationarity=False,
                                             enforce_invertibility=False)

    model_results = sarima_model.fit()
    if model_results.aic < smallest_aic:
        smallest_aic = model_results.aic
        optimal_order_param = order_param
        optimal_seasonal_param = seasonal_param
except Exception as e:
    print(f"Exception for parameters {order_param}x{seasonal_param}: {e}")
    continue

time_series_diff = time_series.diff().dropna()