Библиотека Statsmodels предоставляет широкий спектр инструментов для работы с временными рядами, от базовых методов анализа до продвинутых эконометрических моделей. Она полезна не только для академических исследований, но и для практического применения в трейдинге, риск-менеджменте и финансовом планировании.
Что такое Statsmodels?
Statsmodels — это Python-библиотека, предоставляющая классы и функции для оценки различных статистических моделей, проведения статистических тестов и исследования данных. В контексте финансовых временных рядов она особенно ценна благодаря модулю tsa (Time Series Analysis).
Перед погружением в сложные модели важно понять базовую архитектуру библиотеки и ее возможности. Statsmodels следует последовательному подходу к моделированию: определение модели, подгонка под данные, исследование результатов и проведение диагностики. Этот процесс обеспечивает гибкость и контроль на каждом этапе анализа.
Ключевое преимущество Statsmodels перед другими инструментами — это его сосредоточенность на статистической теории. В отличие от библиотек машинного обучения, которые часто работают как «черные ящики», Statsmodels позволяет получить детальное представление о статистических свойствах модели. Это критически важно для работы с финансовыми данными, где интерпретируемость моделей часто имеет приоритет перед чистой производительностью.
Установка и настройка окружения для работы с Statsmodels
Прежде чем приступить к анализу, необходимо настроить рабочее окружение. Установка Statsmodels проста и может быть выполнена с помощью менеджера пакетов pip или conda:
pip install statsmodels
Или через conda:
conda install -c conda-forge statsmodels
Кроме самого Statsmodels, для эффективной работы с финансовыми временными рядами понадобятся дополнительные библиотеки: numpy, pandas, matplotlib, yfinance, alpha_vantage.
После установки необходимых пакетов можно создать базовый шаблон для импорта и настройки окружения:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.api as sm
import statsmodels.tsa.api as tsa
import yfinance as yf
# Настройка для графиков
plt.style.use('seaborn-v0_8-darkgrid')
plt.rcParams['figure.figsize'] = (12, 6)
plt.rcParams['font.size'] = 12
# Отключение предупреждений для более чистого вывода
import warnings
warnings.filterwarnings("ignore")
Этот базовый шаблон я использую часто проектах по анализу финансовых временных рядов. Он обеспечивает доступ к необходимым инструментам и настраивает визуализацию для лучшего восприятия результатов.
Загрузка и предварительная обработка финансовых данных
Анализ временных рядов начинается с получения качественных данных. Для финансовых временных рядов я обычно использую либо yfinance для данных с Yahoo Finance, либо Alpha Vantage API для более детальной информации.
Вот пример загрузки и предварительной обработки данных акций с использованием yfinance:
# Загрузка исторических данных Nikkei-225 за последние 3 года
ticker = "^N225"
start_date = "2022-05-14"
end_date = "2025-05-14"
# Получение данных
data = yf.download(ticker, start=start_date, end=end_date)
# Фокус на ценах закрытия
prices = data['Close']
# Базовая обработка
prices = prices.asfreq('B') # Установка бизнес-дневной частоты
prices = prices.fillna(method='bfill') # Заполнение пропущенных значений
# Расчет логарифмических доходностей
returns = np.log(prices / prices.shift(1)).dropna()
print(f"Загружено {len(prices)} точек данных для {ticker}")
print(f"Период: с {prices.index[0].date()} по {prices.index[-1].date()}")
# Визуализация цен и доходностей
fig, ax = plt.subplots(2, 1, figsize=(14, 10))
ax[0].plot(prices)
ax[0].set_title(f'{ticker} Цены закрытия')
ax[0].set_ylabel('Цена (USD)')
ax[1].plot(returns)
ax[1].set_title(f'{ticker} Логарифмические доходности')
ax[1].set_ylabel('Доходность')
plt.tight_layout()
plt.show()
# Базовая статистика
print("\nСтатистика цен:")
print(prices.describe())
print("\nСтатистика доходностей:")
print(returns.describe())
Рис. 1: График динамики индекса Nikkei-225 и его ежедневной логарифмической доходности
Загружено 782 точек данных для ^N225
Период: с 2022-05-16 по 2025-05-13
Статистика цен:
Ticker ^N225
count 782.000000
mean 33354.286862
std 4882.863918
min 25716.859375
25% 28006.930664
50% 33206.689453
75% 38382.238281
max 42224.019531
Статистика доходностей:
Ticker ^N225
count 781.000000
mean 0.000465
std 0.013649
min -0.132341
25% -0.005356
50% 0.000430
75% 0.007010
max 0.097366
После выполнения этого кода мы получаем не только данные, но и их базовую визуализацию, а также статистические характеристики. Это позволяет сразу оценить общие свойства ряда: тренды, волатильность, наличие выбросов.
Важный аспект предварительной обработки финансовых временных рядов — работа с пропущенными данными. Биржевые данные часто содержат пропуски из-за выходных дней и праздников. В примере выше я использовал метод прямого заполнения (backward fill), что является распространенной практикой для финансовых данных. Однако в зависимости от конкретной задачи могут быть более подходящие стратегии заполнения пропусков.
Исследовательский анализ финансовых временных рядов
Перед применением сложных моделей необходимо провести тщательное исследование данных. Это помогает выявить ключевые характеристики временного ряда и выбрать подходящие методы моделирования.
Анализ стационарности временных рядов
Стационарность — фундаментальное понятие при анализе временных рядов. Стационарный ряд имеет постоянное среднее значение, дисперсию и автокорреляционную структуру во времени. Большинство моделей временных рядов требуют стационарности данных, поэтому проверка и обеспечение стационарности являются критическими шагами.
Для проверки стационарности я обычно использую тест Дики-Фуллера (Augmented Dickey-Fuller test) и визуальный анализ:
from statsmodels.tsa.stattools import adfuller
def test_stationarity(timeseries, title=''):
# Определение скользящей средней и стандартного отклонения
rolling_mean = timeseries.rolling(window=252).mean() # 252 торговых дня в году
rolling_std = timeseries.rolling(window=252).std()
# Построение графиков
plt.figure(figsize=(14, 7))
plt.title(f'Анализ стационарности для {title}')
plt.plot(timeseries, label='Исходные данные')
plt.plot(rolling_mean, label='Скользящее среднее (252 дня)')
plt.legend()
plt.tight_layout()
plt.show()
# Проведение теста Дики-Фуллера
print("Результаты теста Дики-Фуллера:")
result = adfuller(timeseries.dropna())
labels = ['ADF Статистика', 'p-значение', '# Использованных лагов', '# Наблюдений']
for value, label in zip(result[:4], labels):
print(f'{label}: {value}')
critical_values = result[4]
for key, value in critical_values.items():
print(f'Критическое значение ({key}): {value}')
if result[1] <= 0.05:
print("\nРяд СТАЦИОНАРЕН при уровне значимости 5%")
else:
print("\nРяд НЕСТАЦИОНАРЕН при уровне значимости 5%")
# Проверка стационарности цен
test_stationarity(prices, title=f"{ticker} Цены закрытия")
# Проверка стационарности доходностей
test_stationarity(returns, title=f"{ticker} Логарифмические доходности")
Рис. 2: Анализ стационарности цен закрытия индекса Nikkei-225
Результаты теста Дики-Фуллера:
ADF Статистика: -1.391018604764019
p-значение: 0.5865742578333997
# Использованных лагов: 3
# Наблюдений: 778
Критическое значение (1%): -3.438783171038672
Критическое значение (5%): -2.865262118650577
Критическое значение (10%): -2.568752018688748
Ряд НЕСТАЦИОНАРЕН при уровне значимости 5%
Рис. 3: Анализ стационарности временного ряда логарифмических доходностей Nikkei-225
Результаты теста Дики-Фуллера:
ADF Статистика: -17.502441296076615
p-значение: 4.378188200626437e-30
# Использованных лагов: 2
# Наблюдений: 778
Критическое значение (1%): -3.438783171038672
Критическое значение (5%): -2.865262118650577
Критическое значение (10%): -2.568752018688748
Ряд СТАЦИОНАРЕН при уровне значимости 5%
В финансовом анализе редко можно встретить стационарные ценовые ряды. Обычно цены активов демонстрируют тренды и изменяющуюся во времени волатильность. Однако ряды доходностей (особенно логарифмических) часто проявляют свойства стационарности, что делает их более подходящими для моделирования.
Если ряд нестационарен, его можно преобразовать к стационарному виду с помощью дифференцирования, удаления тренда или других методов трансформации. Statsmodels предоставляет инструменты для таких преобразований.
# Дифференцирование ряда
diff_prices = prices.diff().dropna()
# Удаление тренда
from statsmodels.tsa.filters.hp_filter import hpfilter
cycle, trend = hpfilter(prices, lamb=1600 * 365**2) # Модифицированный параметр для дневных данных
# Проверка стационарности дифференцированного ряда
test_stationarity(diff_prices, title=f"{ticker} Дифференцированные цены")
Анализ автокорреляции и частичной автокорреляции
Автокорреляционные функции (ACF) и частичные автокорреляционные функции (PACF) — мощные инструменты для изучения структуры зависимостей во временном ряде. Они помогают определить подходящие параметры для моделей ARIMA и выявить сезонные компоненты.
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
def plot_acf_pacf(series, lags=40, title=''):
"""
Построение графиков ACF и PACF
"""
fig, axes = plt.subplots(2, 1, figsize=(14, 10))
# ACF
plot_acf(series, lags=lags, ax=axes[0])
axes[0].set_title(f'Автокорреляционная функция для {title}')
# PACF
plot_pacf(series, lags=lags, ax=axes[1])
axes[1].set_title(f'Частичная автокорреляционная функция для {title}')
plt.tight_layout()
plt.show()
# Анализ автокорреляции для доходностей
plot_acf_pacf(returns, title=f"{ticker} Логарифмические доходности")
# Анализ автокорреляции для квадратов доходностей (волатильность)
plot_acf_pacf(returns**2, title=f"{ticker} Квадраты доходностей (волатильность)")
Рис. 4: Графики автокорреляционной и частичной автокорреляционной функции логарифмических доходностей индекса
Анализ ACF и PACF для логарифмических доходностей и их квадратов дает ценную информацию о структуре временного ряда. В частности:
- ACF и PACF доходностей показывают, есть ли линейная зависимость между текущими и прошлыми значениями;
- ACF и PACF квадратов доходностей выявляют кластеризацию волатильности — типичное явление для финансовых рядов, когда периоды высокой волатильности группируются вместе.
Интерпретация графиков ACF и PACF требует определенного опыта. В общем случае, быстро затухающая ACF и PACF с выбросами на определенных лагах указывают на то, что ряд может быть хорошо описан моделью ARMA с соответствующими параметрами.
Декомпозиция временных рядов
Декомпозиция временного ряда — это процесс разделения ряда на составляющие компоненты: тренд, сезонность и остаток. Это особенно полезно для анализа долгосрочных финансовых данных, которые могут содержать сезонные паттерны.
from statsmodels.tsa.seasonal import seasonal_decompose
# Для декомпозиции лучше использовать данные с регулярной частотой
# Поэтому убедимся, что у нас нет пропусков
regular_prices = prices.asfreq('B').fillna(method='bfill')
# Декомпозиция ряда
decomposition = seasonal_decompose(regular_prices, model='multiplicative', period=252) # 252 торговых дня в году
# Визуализация результатов
fig, axes = plt.subplots(4, 1, figsize=(14, 16))
decomposition.observed.plot(ax=axes[0], title='Наблюдаемый ряд')
decomposition.trend.plot(ax=axes[1], title='Тренд')
decomposition.seasonal.plot(ax=axes[2], title='Сезонность')
decomposition.resid.plot(ax=axes[3], title='Остатки')
plt.tight_layout()
plt.show()
# Проверка стационарности остатков
residuals = decomposition.resid.dropna()
test_stationarity(residuals, title="Остатки после декомпозиции")
Рис. 5: Декомпозиция временного ряда на тренд, сезонность и остатки
Результаты теста Дики-Фуллера:
ADF Статистика: -2.8655080272603715
p-значение: 0.04950847898161726
# Использованных лагов: 1
# Наблюдений: 528
Критическое значение (1%): -3.4427957890025533
Критическое значение (5%): -2.867029512430173
Критическое значение (10%): -2.5696937122646926
Ряд СТАЦИОНАРЕН при уровне значимости 5%
Декомпозиция помогает увидеть «скрытые» компоненты финансового временного ряда. Например, на графике сезонности можно обнаружить годовые циклы, которые могут быть связаны с экономическими или отраслевыми факторами.
Остатки после декомпозиции должны быть стационарны и близки к белому шуму, если модель декомпозиции хорошо подходит для данных. Если это не так, возможно, необходимо использовать другую модель декомпозиции или учесть дополнительные факторы.
Моделирование финансовых временных рядов с помощью Statsmodels
После тщательного исследования данных можно переходить к моделированию. Statsmodels предлагает богатый набор моделей для временных рядов, от классических авторегрессионных до более сложных моделей с учетом гетероскедастичности.
ARIMA и SARIMAX модели для прогнозирования цен
Модели ARIMA (AutoRegressive Integrated Moving Average) и их сезонные расширения SARIMAX (Seasonal ARIMA with eXogenous variables) являются классическими инструментами для анализа временных рядов. Хотя они не всегда оптимальны для финансовых рынков, понимание этих моделей создает основу для более продвинутых подходов.
!pip install "numpy==1.23.5"
!pip install pmdarima==2.0.4 # Библиотека для автоматического подбора параметров ARIMA
import pmdarima as pm
from statsmodels.tsa.statespace.sarimax import SARIMAX
# Для демонстрации будем использовать лог-доходности
returns_data = returns.copy()
# Разделение на обучающую и тестовую выборки
train_size = int(len(returns_data) * 0.8)
train, test = returns_data[:train_size], returns_data[train_size:]
# Автоматический подбор параметров ARIMA
auto_model = pm.auto_arima(
train,
start_p=0, start_q=0,
max_p=5, max_q=5,
m=5, # Частота для сезонной компоненты
seasonal=True,
d=None, # Позволить модели определить порядок интеграции
trace=True,
error_action='ignore',
suppress_warnings=True,
stepwise=True
)
print(auto_model.summary())
# Построение SARIMAX модели с подобранными параметрами
best_order = auto_model.order
best_seasonal_order = auto_model.seasonal_order
model = SARIMAX(
train,
order=best_order,
seasonal_order=best_seasonal_order,
enforce_stationarity=False,
enforce_invertibility=False
)
model_fit = model.fit(disp=False)
print(model_fit.summary())
# Прогнозирование на тестовый период
forecast = model_fit.get_forecast(steps=len(test))
forecast_mean = forecast.predicted_mean
forecast_ci = forecast.conf_int()
# Визуализация результатов
plt.figure(figsize=(14, 7))
plt.title(f'SARIMAX Прогноз для {ticker} доходностей')
plt.plot(train.index, train, label='Обучающие данные')
plt.plot(test.index, test, label='Тестовые данные')
plt.plot(test.index, forecast_mean, label='Прогноз')
plt.fill_between(test.index,
forecast_ci.iloc[:, 0],
forecast_ci.iloc[:, 1],
color='gray', alpha=0.2, label='95% доверительный интервал')
plt.legend()
plt.tight_layout()
plt.show()
# Оценка качества прогноза
from sklearn.metrics import mean_squared_error, mean_absolute_error
import math
rmse = math.sqrt(mean_squared_error(test, forecast_mean))
mae = mean_absolute_error(test, forecast_mean)
print(f"Средняя квадратическая ошибка (RMSE): {rmse:.6f}")
print(f"Средняя абсолютная ошибка (MAE): {mae:.6f}")
Performing stepwise search to minimize aic
ARIMA(0,0,0)(1,0,1)[5] intercept : AIC=-3629.641, Time=0.56 sec
ARIMA(0,0,0)(0,0,0)[5] intercept : AIC=-3632.558, Time=0.14 sec
ARIMA(1,0,0)(1,0,0)[5] intercept : AIC=-3630.707, Time=0.26 sec
ARIMA(0,0,1)(0,0,1)[5] intercept : AIC=-3630.693, Time=0.54 sec
ARIMA(0,0,0)(0,0,0)[5] : AIC=-3633.250, Time=0.09 sec
ARIMA(0,0,0)(1,0,0)[5] intercept : AIC=-3631.266, Time=0.31 sec
ARIMA(0,0,0)(0,0,1)[5] intercept : AIC=-3631.179, Time=0.54 sec
ARIMA(1,0,0)(0,0,0)[5] intercept : AIC=-3632.041, Time=0.17 sec
ARIMA(0,0,1)(0,0,0)[5] intercept : AIC=-3632.102, Time=0.17 sec
ARIMA(1,0,1)(0,0,0)[5] intercept : AIC=-3631.411, Time=0.48 sec
Best model: ARIMA(0,0,0)(0,0,0)[5]
Total fit time: 3.323 seconds
SARIMAX Results
==============================================================================
Dep. Variable: y No. Observations: 624
Model: SARIMAX Log Likelihood 1817.625
Date: Thu, 15 May 2025 AIC -3633.250
Time: 22:24:14 BIC -3628.814
Sample: 05-17-2022 HQIC -3631.526
- 10-04-2024
Covariance Type: opg
==============================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------
sigma2 0.0002 2.88e-06 60.040 0.000 0.000 0.000
===================================================================================
Ljung-Box (L1) (Q): 1.49 Jarque-Bera (JB): 12115.32
Prob(Q): 0.22 Prob(JB): 0.00
Heteroskedasticity (H): 2.62 Skew: -1.29
Prob(H) (two-sided): 0.00 Kurtosis: 24.43
===================================================================================
Warnings:
[1] Covariance matrix calculated using the outer product of gradients (complex-step).
SARIMAX Results
==============================================================================
Dep. Variable: ^N225 No. Observations: 624
Model: SARIMAX Log Likelihood 1814.264
Date: Thu, 15 May 2025 AIC -3626.529
Time: 22:24:15 BIC -3622.094
Sample: 05-17-2022 HQIC -3624.805
- 10-04-2024
Covariance Type: opg
==============================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------
sigma2 0.0002 2.88e-06 59.995 0.000 0.000 0.000
===================================================================================
Ljung-Box (L1) (Q): 1.51 Jarque-Bera (JB): 12055.94
Prob(Q): 0.22 Prob(JB): 0.00
Heteroskedasticity (H): 2.60 Skew: -1.29
Prob(H) (two-sided): 0.00 Kurtosis: 24.40
===================================================================================
Рис. 6: SARIMAX Прогноз доходностей Nikkei-225
Средняя квадратическая ошибка (RMSE): 0.015490
Средняя абсолютная ошибка (MAE): 0.009464
В этом примере я использую библиотеку pmdarima для автоматического подбора параметров ARIMA, что упрощает процесс моделирования. Затем на основе найденных параметров строится модель SARIMAX, которая учитывает возможные сезонные компоненты.
Важно отметить, что для финансовых временных рядов модели ARIMA часто имеют ограниченную предсказательную способность. Это связано с тем, что финансовые рынки сложны, нелинейны и подвержены влиянию множества факторов, не учитываемых в простых авторегрессионных моделях.
Тем не менее, ARIMA-модели могут служить полезным бенчмарком и базой для более сложных моделей.
VAR модели для многомерных временных рядов
Векторные авторегрессионные (VAR) модели позволяют анализировать взаимосвязи между несколькими финансовыми временными рядами. Это особенно полезно для изучения взаимного влияния различных активов или макроэкономических показателей.
Код ниже строит модель VAR между индексом Nikkei-225 и 4 акциями топовых японских корпораций: Toyota Motor, Sony, Mitsubishi, Canon.
from statsmodels.tsa.vector_ar.var_model import VAR
# Загрузим данные для нескольких тикеров
tickers = ["^N225", "TM", "SONY", "MUFG", "CAJFF"]
start_date = "2022-05-14"
end_date = "2025-05-14"
# Получение данных
multi_data = yf.download(tickers, start=start_date, end=end_date)['Close']
# Рассчитаем логарифмические доходности
multi_returns = np.log(multi_data / multi_data.shift(1)).dropna()
# Разделение на обучающую и тестовую выборки
train_size = int(len(multi_returns) * 0.8)
train, test = multi_returns[:train_size], multi_returns[train_size:]
# Подбор оптимального порядка VAR модели
max_lags = 15
model = VAR(train)
results = {}
for i in range(1, max_lags + 1):
result = model.fit(i)
results[i] = result.aic
best_lag = min(results, key=results.get)
print(f"Оптимальный лаг для VAR модели: {best_lag}")
# Построение VAR модели с оптимальным порядком
var_model = VAR(train)
var_results = var_model.fit(best_lag)
print(var_results.summary())
# Прогнозирование на тестовый период
forecast = var_results.forecast(train.values[-best_lag:], steps=len(test))
forecast_df = pd.DataFrame(forecast, index=test.index, columns=test.columns)
# Визуализация результатов для одного из активов (например, AAPL)
plt.figure(figsize=(14, 7))
plt.title(f'VAR Прогноз доходностей для {tickers[0]}')
plt.plot(train.index, train[tickers[0]], label='Обучающие данные')
plt.plot(test.index, test[tickers[0]], label='Тестовые данные')
plt.plot(test.index, forecast_df[tickers[0]], label='Прогноз')
plt.legend()
plt.tight_layout()
plt.show()
# Оценка качества прогноза для всех активов
for ticker in tickers:
rmse = math.sqrt(mean_squared_error(test[ticker], forecast_df[ticker]))
mae = mean_absolute_error(test[ticker], forecast_df[ticker])
print(f"{ticker} - RMSE: {rmse:.6f}, MAE: {mae:.6f}")
Оптимальный лаг для VAR модели: 1
Summary of Regression Results
==================================
Model: VAR
Method: OLS
Date: Thu, 15, May, 2025
Time: 22:32:27
--------------------------------------------------------------------
No. of Equations: 5.00000 BIC: -41.6684
Nobs: 513.000 HQIC: -41.8191
Log likelihood: 7141.96 FPE: 6.25138e-19
AIC: -41.9163 Det(Omega_mle): 5.89829e-19
--------------------------------------------------------------------
Results for equation CAJFF
===========================================================================
coefficient std. error t-stat prob
---------------------------------------------------------------------------
const 0.000814 0.000716 1.137 0.256
L1.CAJFF -0.012926 0.044305 -0.292 0.770
L1.MUFG -0.032653 0.041560 -0.786 0.432
L1.SONY 0.031572 0.047649 0.663 0.508
L1.TM 0.097837 0.057081 1.714 0.087
L1.^N225 -0.011038 0.056007 -0.197 0.844
===========================================================================
Results for equation MUFG
===========================================================================
coefficient std. error t-stat prob
---------------------------------------------------------------------------
const 0.000727 0.000847 0.857 0.391
L1.CAJFF 0.071298 0.052418 1.360 0.174
L1.MUFG 0.015485 0.049170 0.315 0.753
L1.SONY -0.102523 0.056375 -1.819 0.069
L1.TM 0.105939 0.067533 1.569 0.117
L1.^N225 0.101775 0.066263 1.536 0.125
===========================================================================
Results for equation SONY
===========================================================================
coefficient std. error t-stat prob
---------------------------------------------------------------------------
const 0.000237 0.000766 0.310 0.757
L1.CAJFF 0.028404 0.047393 0.599 0.549
L1.MUFG -0.012228 0.044457 -0.275 0.783
L1.SONY -0.041218 0.050970 -0.809 0.419
L1.TM 0.061582 0.061059 1.009 0.313
L1.^N225 0.054088 0.059910 0.903 0.367
===========================================================================
Results for equation TM
===========================================================================
coefficient std. error t-stat prob
---------------------------------------------------------------------------
const 0.000174 0.000688 0.253 0.800
L1.CAJFF 0.007966 0.042533 0.187 0.851
L1.MUFG -0.069511 0.039898 -1.742 0.081
L1.SONY 0.000722 0.045744 0.016 0.987
L1.TM 0.105215 0.054798 1.920 0.055
L1.^N225 0.020301 0.053767 0.378 0.706
===========================================================================
Results for equation ^N225
===========================================================================
coefficient std. error t-stat prob
---------------------------------------------------------------------------
const 0.000276 0.000551 0.500 0.617
L1.CAJFF -0.015819 0.034099 -0.464 0.643
L1.MUFG 0.082118 0.031986 2.567 0.010
L1.SONY 0.154335 0.036673 4.208 0.000
L1.TM 0.182452 0.043932 4.153 0.000
L1.^N225 -0.162754 0.043106 -3.776 0.000
===========================================================================
Correlation matrix of residuals
CAJFF MUFG SONY TM ^N225
CAJFF 1.000000 0.047916 0.009617 0.033420 -0.011984
MUFG 0.047916 1.000000 0.301288 0.437885 0.139969
SONY 0.009617 0.301288 1.000000 0.472900 0.254494
TM 0.033420 0.437885 0.472900 1.000000 0.319147
^N225 -0.011984 0.139969 0.254494 0.319147 1.000000
Рис. 7: Прогноз доходностей для Nikkei-225 через построение VAR модели
^N225 - RMSE: 0.017561, MAE: 0.011213
TM - RMSE: 0.019036, MAE: 0.013708
SONY - RMSE: 0.021317, MAE: 0.015089
MUFG - RMSE: 0.023186, MAE: 0.016014
CAJFF - RMSE: 0.029898, MAE: 0.014763
VAR модели позволяют учитывать не только автокорреляцию каждого ряда, но и кросс-корреляцию между различными рядами. Это может быть особенно полезно при анализе тесно связанных финансовых инструментов, например, акций компаний из одного сектора.
Определение оптимального порядка лага — важная часть построения VAR модели. В примере использован информационный критерий Акаике (AIC), но могут применяться и другие критерии, такие как BIC или FPE.
Моделирование волатильности с использованием Statsmodels
Финансовые временные ряды характеризуются изменчивой волатильностью — периоды высокой волатильности чередуются с периодами относительного спокойствия. В Statsmodels есть инструменты для анализа и моделирования этой важной характеристики финансовых данных.
В модуле statsmodels.tsa есть несколько методов для оценки и анализа волатильности:
from statsmodels.tsa.stattools import acf, pacf
from statsmodels.stats.diagnostic import het_arch
import numpy as np
import matplotlib.pyplot as plt
# Будем работать с квадратами доходностей как мерой волатильности
volatility = returns**2
# Визуализация волатильности
plt.figure(figsize=(14, 7))
plt.title(f'Волатильность {ticker} (квадраты логарифмических доходностей)')
plt.plot(volatility)
plt.ylabel('Квадраты доходностей')
plt.tight_layout()
plt.show()
# Тест на наличие ARCH эффектов
lm_stat, lm_pval, f_stat, f_pval = het_arch(returns)
print("Тест на наличие ARCH эффектов (гетероскедастичности):")
print(f"LM статистика: {lm_stat:.4f}, p-значение: {lm_pval:.4f}")
print(f"F статистика: {f_stat:.4f}, p-значение: {f_pval:.4f}")
if lm_pval < 0.05:
print("Обнаружен эффект ARCH (гетероскедастичность) при уровне значимости 5%")
else:
print("Эффект ARCH не обнаружен при уровне значимости 5%")
# Анализ автокорреляции волатильности
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(14, 10))
acf_vals = acf(volatility.dropna(), nlags=40)
pacf_vals = pacf(volatility.dropna(), nlags=40, method='ols')
# Построение ACF
ax1.stem(range(len(acf_vals)), acf_vals)
ax1.set_title('ACF квадратов доходностей (волатильность)')
ax1.axhline(y=0, linestyle='-', color='black')
ax1.axhline(y=1.96/np.sqrt(len(volatility)), linestyle='--', color='gray')
ax1.axhline(y=-1.96/np.sqrt(len(volatility)), linestyle='--', color='gray')
# Построение PACF
ax2.stem(range(len(pacf_vals)), pacf_vals)
ax2.set_title('PACF квадратов доходностей (волатильность)')
ax2.axhline(y=0, linestyle='-', color='black')
ax2.axhline(y=1.96/np.sqrt(len(volatility)), linestyle='--', color='gray')
ax2.axhline(y=-1.96/np.sqrt(len(volatility)), linestyle='--', color='gray')
plt.tight_layout()
plt.show()
# Расчет скользящей волатильности
window_sizes = [21, 63, 252] # Окна для дневной, квартальной и годовой волатильности
rolling_vol = pd.DataFrame(index=returns.index)
for window in window_sizes:
rolling_vol[f'{window}-дневная'] = returns.rolling(window).std() * np.sqrt(252) * 100 # Аннуализированная волатильность в %
# Визуализация скользящей волатильности
plt.figure(figsize=(14, 7))
plt.title(f'Скользящая волатильность {ticker} (аннуализированная, %)')
for column in rolling_vol.columns:
plt.plot(rolling_vol[column], label=column)
plt.legend()
plt.tight_layout()
plt.show()
Рис. 8: Моделирование волатильности методом CAJFF (квадраты логарифмических доходностей)
Тест на наличие ARCH эффектов (гетероскедастичности):
LM статистика: 222.0012, p-значение: 0.0000
F статистика: 30.7325, p-значение: 0.0000
Обнаружен эффект ARCH (гетероскедастичность) при уровне значимости 5%
Рис. 9: График модели скользящей аннуализированной волатильности CAJFF
Анализ автокорреляционной функции для квадратов доходностей позволяет выявить паттерны волатильности во временном ряде. Высокие значения ACF на многих лагах свидетельствуют о стойкой кластеризации волатильности — типичном явлении для финансовых рынков.
Для сложных задач моделирования волатильности Statsmodels также предлагает реализацию условно гетероскедастичных моделей в модуле statsmodels.tsa.regime_switching.markov_regression. Хотя эти модели более сложны в использовании, они позволяют учитывать изменения режимов волатильности на рынке.
Расчет скользящей волатильности с разными окнами показывает, как изменяется волатильность актива с течением времени, что чрезвычайно важно для оценки рисков и принятия инвестиционных решений.
Диагностика и оценка качества моделей временных рядов
Построение модели — только часть работы с финансовыми временными рядами. Не менее важна тщательная диагностика модели и оценка качества ее прогнозов. Statsmodels предоставляет широкий набор инструментов для проверки адекватности моделей и выбора наиболее подходящей для конкретной задачи.
Первым шагом в диагностике модели временного ряда является анализ остатков. Остатки хорошей модели должны представлять собой белый шум — случайный процесс без автокорреляции и с постоянной дисперсией.
import matplotlib.pyplot as plt
import seaborn as sns
from statsmodels.stats.diagnostic import acorr_ljungbox
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
from scipy import stats
import numpy as np
from statsmodels.tsa.arima.model import ARIMA
# Для демонстрации будем использовать логарифмические доходности
returns_data = returns.copy()
# Разделение на обучающую и тестовую выборки
train_size = int(len(returns_data) * 0.8)
train, test = returns_data[:train_size], returns_data[train_size:]
# Построение модели ARIMA(1,0,1)
arima_model = ARIMA(train, order=(1, 0, 1))
arima_result = arima_model.fit()
print(arima_result.summary())
# Получаем остатки
residuals = arima_result.resid
# Статистический анализ остатков
print("\n=== Анализ остатков ===")
print("Среднее остатков:", residuals.mean())
print("Стандартное отклонение остатков:", residuals.std())
print("Асимметрия (Skewness):", stats.skew(residuals))
print("Эксцесс (Kurtosis):", stats.kurtosis(residuals))
# Тест Льюнга-Бокса на автокорреляцию остатков
lb_test = acorr_ljungbox(residuals, lags=10, return_df=True)
print("\nТест Льюнга-Бокса (p-values):")
print(lb_test)
# Визуализация остатков
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# График остатков
axes[0, 0].plot(residuals)
axes[0, 0].set_title('Остатки модели')
# Гистограмма остатков
sns.histplot(residuals, kde=True, ax=axes[0, 1])
axes[0, 1].set_title('Гистограмма остатков')
# Q-Q Plot
stats.probplot(residuals, dist="norm", plot=axes[1, 0])
axes[1, 0].set_title('Q-Q Plot остатков')
# ACF и PACF остатков
plot_acf(residuals, ax=axes[1, 1], lags=20)
axes[1, 1].set_title('ACF остатков')
plt.tight_layout()
plt.show()
SARIMAX Results
==============================================================================
Dep. Variable: ^N225 No. Observations: 624
Model: ARIMA(1, 0, 1) Log Likelihood 1819.680
Date: Thu, 15 May 2025 AIC -3631.361
Time: 22:53:49 BIC -3613.616
Sample: 05-17-2022 HQIC -3624.465
- 10-04-2024
Covariance Type: opg
==============================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------
const 0.0006 0.001 1.124 0.261 -0.000 0.002
ar.L1 0.6318 0.358 1.766 0.077 -0.069 1.333
ma.L1 -0.6812 0.348 -1.958 0.050 -1.363 0.001
sigma2 0.0002 2.97e-06 57.707 0.000 0.000 0.000
===================================================================================
Ljung-Box (L1) (Q): 0.00 Jarque-Bera (JB): 12876.81
Prob(Q): 0.95 Prob(JB): 0.00
Heteroskedasticity (H): 2.59 Skew: -1.60
Prob(H) (two-sided): 0.00 Kurtosis: 25.02
===================================================================================
Warnings:
[1] Covariance matrix calculated using the outer product of gradients (complex-step).
=== Анализ остатков ===
Среднее остатков: 3.803258336713272e-06
Стандартное отклонение остатков: 0.01311112427849995
Асимметрия (Skewness): -1.60479231044423
Эксцесс (Kurtosis): 22.021559904926043
Тест Льюнга-Бокса (p-values):
lb_stat lb_pvalue
1 0.003850 0.950522
2 0.087379 0.957251
3 0.675196 0.879022
4 1.290128 0.863048
5 1.700045 0.888894
6 1.702979 0.944891
7 1.797455 0.970196
8 1.871963 0.984689
9 3.766020 0.926127
10 5.517184 0.854067
Рис. 10: Визуализация остатков модели, гистограммы, Q-Q plot и ACF
Помимо классических моделей ARIMA и VAR, библиотека Statsmodels предлагает более сложные и гибкие инструменты для работы с финансовыми временными рядами. Эти модели позволяют учесть различные особенности финансовых данных: структурные сдвиги, нелинейность, асимметрию реакции на положительные и отрицательные шоки и другие.
Модели с марковским переключением режимов
Финансовые рынки часто демонстрируют разные режимы поведения — например, режим высокой волатильности во время кризиса и режим низкой волатильности в спокойные периоды. Модели с марковским переключением режимов позволяют учесть эти различия и адаптировать прогнозы к текущему состоянию рынка.
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.tsa.regime_switching.markov_autoregression import MarkovAutoregression
# Используем ряд логарифмических доходностей
returns_data = returns.copy()
# Разделение на обучающую и тестовую выборки
train_size = int(len(returns_data) * 0.8)
train, test = returns_data[:train_size], returns_data[train_size:]
# Построение модели с двумя режимами
markov_model = MarkovAutoregression(train, k_regimes=2, order=1, switching_variance=True)
markov_result = markov_model.fit()
print(markov_result.summary())
# Получение вероятностей режимов
filtered_probs = markov_result.filtered_marginal_probabilities[0]
smoothed_probs = markov_result.smoothed_marginal_probabilities[0]
# Визуализация вероятностей режимов
fig, axes = plt.subplots(3, 1, figsize=(14, 12))
# Оригинальный ряд
axes[0].plot(train.values, label='Доходности')
axes[0].set_title('Логарифмические доходности ^N225')
axes[0].grid(True)
# Отфильтрованные вероятности режимов
axes[1].plot(train.index[1:], filtered_probs, label='Режим 1', color='blue')
axes[1].set_title('Отфильтрованные вероятности режимов')
axes[1].set_ylim(0, 1)
axes[1].legend()
axes[1].grid(True)
# Сглаженные вероятности режимов
axes[2].plot(train.index[1:], smoothed_probs, label='Режим 1', color='green')
axes[2].set_title('Сглаженные вероятности режимов')
axes[2].set_ylim(0, 1)
axes[2].legend()
axes[2].grid(True)
plt.tight_layout()
plt.xticks(rotation=45)
plt.show()
# Вывод характеристик режимов
print("\nДоступные параметры модели:")
print(markov_result.params.index.tolist())
print("\nХарактеристики режимов:")
for i in range(2):
mean = markov_result.params[f'const[{i}]']
variance = markov_result.params[f'sigma2[{i}]']
ar_coeff = markov_result.params[f'ar.L1[{i}]']
print(f"\nРежим {i+1}:")
print(f" Среднее: {mean:.6f}")
print(f" Дисперсия: {variance:.6f}")
print(f" AR(1) коэффициент: {ar_coeff:.6f}")
print(f" Волатильность (годовая, %): {np.sqrt(variance) * np.sqrt(252) * 100:.2f}%")
# Расчет матрицы переходов
p_00 = markov_result.params['p[0->0]'] # Вероятность остаться в режиме 0
p_10 = markov_result.params['p[1->0]'] # Вероятность перейти из 1 в 0
p_01 = 1 - p_00 # Вероятность перейти из 0 в 1
p_11 = 1 - p_10 # Вероятность остаться в режиме 1
transition_matrix = np.array([[p_00, p_01], [p_10, p_11]])
print("\nМатрица переходов между режимами:")
print(transition_matrix)
# Ожидаемая продолжительность режимов
duration_0 = 1 / (1 - p_00)
duration_1 = 1 / (1 - p_11)
print(f"\nОжидаемая продолжительность режима 0: {duration_0:.2f} дней")
print(f"Ожидаемая продолжительность режима 1: {duration_1:.2f} дней")
Markov Switching Model Results
================================================================================
Dep. Variable: ^N225 No. Observations: 623
Model: MarkovAutoregression Log Likelihood 1908.125
Date: Thu, 15 May 2025 AIC -3800.250
Time: 23:06:32 BIC -3764.774
Sample: 05-17-2022 HQIC -3786.463
- 10-04-2024
Covariance Type: approx
Regime 0 parameters
==============================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------
const 0.0009 0.000 2.002 0.045 1.84e-05 0.002
sigma2 0.0001 7.13e-06 16.158 0.000 0.000 0.000
ar.L1 -0.0049 0.039 -0.127 0.899 -0.081 0.071
Regime 1 parameters
==============================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------
const -0.0221 0.021 -1.046 0.296 -0.063 0.019
sigma2 0.0031 0.002 1.926 0.054 -5.54e-05 0.006
ar.L1 -0.4607 0.561 -0.821 0.411 -1.560 0.639
Regime transition parameters
==============================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------
p[0->0] 0.9936 0.004 228.764 0.000 0.985 1.002
p[1->0] 0.4117 0.238 1.730 0.084 -0.055 0.878
==============================================================================
Рис. 11: Моделирование доходности индекса марковским переключением режимов
Доступные параметры модели:
['p[0->0]', 'p[1->0]', 'const[0]', 'const[1]', 'sigma2[0]', 'sigma2[1]', 'ar.L1[0]', 'ar.L1[1]']
Характеристики режимов:
Режим 1:
Среднее: 0.000881
Дисперсия: 0.000115
AR(1) коэффициент: -0.004948
Волатильность (годовая, %): 17.04%
Режим 2:
Среднее: -0.022089
Дисперсия: 0.003143
AR(1) коэффициент: -0.460698
Волатильность (годовая, %): 89.00%
Матрица переходов между режимами:
[[0.99364626 0.00635374]
[0.41173915 0.58826085]]
Ожидаемая продолжительность режима 0: 157.39 дней
Ожидаемая продолжительность режима 1: 2.43 дней
Модели с марковским переключением режимов обеспечивают более глубокое понимание динамики финансовых рынков. Они позволяют идентифицировать различные состояния рынка и оценить вероятность перехода между ними. Это особенно ценно для управления рисками и оптимизации портфеля, так как позволяет адаптировать стратегии к различным рыночным условиям.
В приведенном примере модель определяет два режима рынка, которые можно интерпретировать как «спокойный» и «волатильный». Анализ матрицы переходов показывает вероятность перехода из одного состояния в другое, а расчет ожидаемой продолжительности режимов дает представление о типичных периодах стабильности и турбулентности на рынке.
Модели с коинтеграцией и механизмом коррекции ошибок
Коинтеграция — важное понятие при анализе взаимосвязанных финансовых временных рядов. Она описывает ситуацию, когда два или более нестационарных ряда имеют долгосрочную равновесную связь. Это часто наблюдается в парах финансовых инструментов, например, в акциях компаний одной отрасли или валютных парах.
Поскольку в данной статье мы исследуем индекс Nikkei и входящие в него акции, то будет интересно посмотреть есть ли коинтеграция между акциями Sony и Panasonic.
from statsmodels.tsa.vector_ar.vecm import VECM, coint_johansen
# Загрузим данные для пары связанных акций (например, Sony и Panasonic)
tickers = ["SONY", "PCRFF"]
start_date = "2022-05-14"
end_date = "2025-05-14"
# Получение данных
pair_data = yf.download(tickers, start=start_date, end=end_date)['Close']
pair_data = pair_data.asfreq('B').fillna(method='ffill')
# Визуализация цен
plt.figure(figsize=(14, 7))
plt.title('Цены акций Sony и Panasonic')
plt.plot(pair_data.index, pair_data['SONY'], label='Sony (SONY)')
plt.plot(pair_data.index, pair_data['PCRFF'], label='Panasonic (PCRFF)')
plt.legend()
plt.tight_layout()
plt.show()
# Нормализация для лучшего сравнения
normalized_data = pair_data / pair_data.iloc[0]
plt.figure(figsize=(14, 7))
plt.title('Нормализованные цены акций Sony и Panasonic')
plt.plot(normalized_data.index, normalized_data['SONY'], label='Sony (SONY)')
plt.plot(normalized_data.index, normalized_data['PCRFF'], label='Panasonic (PCRFF)')
plt.legend()
plt.tight_layout()
plt.show()
# Тест на коинтеграцию Йохансена
johansen_results = coint_johansen(pair_data, det_order=0, k_ar_diff=1)
print("Результаты теста на коинтеграцию Йохансена:")
print(f"Собственные значения: {johansen_results.eig}")
print(f"Критические значения (90%, 95%, 99%):")
print(johansen_results.cvt)
print(f"Тестовая статистика:")
print(johansen_results.lr1)
# Определение ранга коинтеграции
crit_vals = johansen_results.cvt[:, 1] # 95% критические значения
test_stats = johansen_results.lr1
rank = sum(test_stats > crit_vals)
print(f"Ранг коинтеграции: {rank}")
if rank > 0:
print("Ряды коинтегрированы! Можно применить модель VECM.")
# Построение модели VECM
model = VECM(pair_data, k_ar_diff=2, coint_rank=rank, deterministic='ci')
vecm_results = model.fit()
print(vecm_results.summary())
# Получение коинтеграционного вектора
if rank == 1:
coint_vector = vecm_results.beta
print(f"Коинтеграционный вектор: {coint_vector}")
# Построение коинтеграционного спреда
spread = np.dot(pair_data, coint_vector)
# Визуализация спреда
plt.figure(figsize=(14, 7))
plt.title('Коинтеграционный спред')
plt.plot(pair_data.index, spread)
plt.axhline(y=spread.mean(), color='r', linestyle='-', label='Среднее значение')
plt.axhline(y=spread.mean() + 2*spread.std(), color='g', linestyle='--', label='+2 станд. отклонения')
plt.axhline(y=spread.mean() - 2*spread.std(), color='g', linestyle='--', label='-2 станд. отклонения')
plt.legend()
plt.tight_layout()
plt.show()
# Стратегия парного трейдинга на основе коинтеграции
z_score = (spread - spread.mean()) / spread.std()
# Сигналы для входа в позицию
threshold = 2.0 # Стандартное значение для z-score
signals = pd.DataFrame(index=z_score.index)
signals['z_score'] = z_score
signals['positions'] = 0
# Когда z-score выше threshold, шортим спред (шорт SONY, лонг PCRFF)
signals.loc[signals['z_score'] > threshold, 'positions'] = -1
# Когда z-score ниже -threshold, лонгуем спред (лонг SONY, шорт PCRFF)
signals.loc[signals['z_score'] < -threshold, 'positions'] = 1
# Закрываем позицию, когда z-score возвращается к среднему
signals['positions'] = signals['positions'].shift(1).fillna(0)
signals.loc[np.abs(signals['z_score']) < 0.5, 'positions'] = 0
# Визуализация сигналов
plt.figure(figsize=(14, 10))
plt.subplot(2, 1, 1)
plt.title('Z-score коинтеграционного спреда')
plt.plot(signals.index, signals['z_score'])
plt.axhline(y=threshold, color='r', linestyle='--', label=f'+{threshold} станд. отклонения')
plt.axhline(y=-threshold, color='g', linestyle='--', label=f'-{threshold} станд. отклонения')
plt.axhline(y=0, color='k', linestyle='-')
plt.legend()
plt.subplot(2, 1, 2)
plt.title('Позиции на основе стратегии парного трейдинга')
plt.plot(signals.index, signals['positions'])
plt.axhline(y=0, color='k', linestyle='-')
plt.tight_layout()
plt.show()
else:
print("Ряды не коинтегрированы при уровне значимости 5%.")
Рис. 12: Сравнительный график динамики цен акций Sony и Panasonic за последние 3 года
На первый взгляд выглядит, что акции Sony и стоят дороже, и растут быстрее. Однако это не так. Если мы нормализуем данные, то получим следующую картину.
Рис. 13: Динамика нормализованных цен акций Sony и Panasonic
На графике выше видно, что несмотря на то, что обе компании представляют по сути одну отрасль и производят одну и ту же продукцию, динамика их акций не всегда совпадает: есть периоды когда котировки движутся синхронно, а есть когда Sony существенно отстает от своего конкурента.
После визуального анализа давайте посмотрим что скажут статистические показатели:
Результаты теста на коинтеграцию Йохансена:
Собственные значения: [0.01065511 0.00158335]
Критические значения (90%, 95%, 99%):
[[13.4294 15.4943 19.9349]
[ 2.7055 3.8415 6.6349]]
Тестовая статистика:
[9.59156969 1.23598855]
Ранг коинтеграции: 0
Ряды не коинтегрированы при уровне значимости 5%.
Мы видим, что конкретно эти активы и в этом временном диапазоне не являются коинтегрированными. На практике же, есть немало примеров, когда коинтеграция есть и это является прекрасной возможностью для исследования и торговли долгосрочных отношений между финансовыми активами. Этот подход особенно полезен для:
- Разработки статистических арбитражных стратегий (парный трейдинг);
- Понимания долгосрочного равновесия между связанными активами;
- Создания портфелей с учетом структурной взаимосвязи активов;
- Анализа эффективности ценообразования на различных рынках.
Модель с механизмом коррекции ошибок (VECM) не только учитывает краткосрочную динамику рядов, но и корректирует прогнозы с учетом долгосрочного равновесия, что делает ее особенно ценной для анализа финансовых данных.
Заключение
Анализ финансовых временных рядов с использованием библиотеки Statsmodels открывает широкие возможности для исследователей и практиков.
Эта библиотека для Python сочетает в себе академическую строгость с практической применимостью, предлагая инструменты для всех этапов работы — от первичной обработки данных до сложного моделирования. В отличие от «черных ящиков» машинного обучения, Statsmodels фокусируется на прозрачности моделей. Это критически важно для финансовой аналитики, где понимание факторов риска и волатильности часто важнее точечных прогнозов.
Также стоит отметить универсальной данной библиотеки: она поддерживает как классические подходы (ARIMA, VAR), так и продвинутые методы (модели с марковскими переключениями, коинтеграция).
Statsmodels не просто предоставляет алгоритмы — она формирует методологию работы с данными. Для риск-менеджеров: это в первую очередь инструменты анализа волатильности и стресс-тестирования моделей в различных рыночных режимах. Для аналитиков: готовые решения для декомпозиции рядов, визуализации автокорреляции и диагностики моделей, что ускоряет исследовательский процесс.