Ответ
Да, работал с ARIMA (AutoRegressive Integrated Moving Average) для прогнозирования временных рядов. Вот мой практический опыт:
1. Проект: Прогнозирование ежедневных продаж ритейлера
Данные: Ежедневные продажи за 2 года (730 точек) Задача: Прогноз на 30 дней вперед для управления запасами
2. Полный пайплайн анализа и прогнозирования
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.stattools import adfuller, acf, pacf
from statsmodels.tsa.arima.model import ARIMA
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
from sklearn.metrics import mean_absolute_error, mean_squared_error
import warnings
warnings.filterwarnings('ignore')
# Загрузка и подготовка данных
df = pd.read_csv('daily_sales.csv', parse_dates=['date'], index_col='date')
df = df.asfreq('D') # Установка ежедневной частоты
sales = df['sales']
# Визуализация ряда
plt.figure(figsize=(12, 6))
plt.plot(sales)
plt.title('Daily Sales Time Series')
plt.xlabel('Date')
plt.ylabel('Sales')
plt.grid(True)
plt.show()
3. Проверка стационарности
# Тест Дики-Фуллера
def test_stationarity(timeseries):
print('Results of Dickey-Fuller Test:')
dftest = adfuller(timeseries, autolag='AIC')
dfoutput = pd.Series(dftest[0:4],
index=['Test Statistic', 'p-value',
'#Lags Used', 'Number of Observations Used'])
for key, value in dftest[4].items():
dfoutput[f'Critical Value ({key})'] = value
print(dfoutput)
return dftest[1] # p-value
p_value = test_stationarity(sales)
print(f"Series is {'stationary' if p_value < 0.05 else 'non-stationary'}")
# Если ряд нестационарный - применяем дифференцирование
if p_value >= 0.05:
sales_diff = sales.diff().dropna()
p_value_diff = test_stationarity(sales_diff)
print(f"After 1st difference, p-value: {p_value_diff:.4f}")
4. Определение параметров ARIMA (p, d, q)
# ACF и PACF для определения параметров
fig, axes = plt.subplots(1, 2, figsize=(12, 4))
plot_acf(sales_diff if 'sales_diff' in locals() else sales,
lags=40, ax=axes[0])
plot_pacf(sales_diff if 'sales_diff' in locals() else sales,
lags=40, ax=axes[1], method='ywm')
plt.show()
# На основе графиков:
# p (AR) - где PACF обрезается (lag ~ 7)
# q (MA) - где ACF обрезается (lag ~ 3)
# d - порядок дифференцирования (1 в нашем случае)
5. Подбор модели ARIMA
# Разделение на train/test
split_idx = int(len(sales) * 0.8)
train = sales[:split_idx]
test = sales[split_idx:]
# Модель ARIMA(7,1,3)
model = ARIMA(train, order=(7, 1, 3))
model_fit = model.fit()
print(model_fit.summary())
# Диагностика остатков
residuals = pd.DataFrame(model_fit.resid)
fig, axes = plt.subplots(2, 2, figsize=(12, 8))
residuals.plot(ax=axes[0, 0])
axes[0, 0].set_title('Residuals over Time')
residuals.plot(kind='kde', ax=axes[0, 1])
axes[0, 1].set_title('Residuals Distribution')
plot_acf(residuals, lags=40, ax=axes[1, 0])
plot_pacf(residuals, lags=40, ax=axes[1, 1], method='ywm')
plt.tight_layout()
plt.show()
6. Прогнозирование и оценка
# Прогноз на тестовой выборке
forecast_steps = len(test)
forecast = model_fit.forecast(steps=forecast_steps)
forecast_index = test.index
# Прогноз на 30 дней вперед (для production)
future_forecast = model_fit.forecast(steps=30)
# Визуализация
plt.figure(figsize=(12, 6))
plt.plot(train.index, train, label='Train', color='blue')
plt.plot(test.index, test, label='Actual Test', color='green')
plt.plot(forecast_index, forecast, label='Forecast', color='red', linestyle='--')
plt.fill_between(forecast_index,
forecast - 1.96*model_fit.resid.std(),
forecast + 1.96*model_fit.resid.std(),
alpha=0.2, color='red')
plt.title('ARIMA Forecast vs Actual')
plt.xlabel('Date')
plt.ylabel('Sales')
plt.legend()
plt.grid(True)
plt.show()
# Метрики качества
mae = mean_absolute_error(test, forecast)
rmse = np.sqrt(mean_squared_error(test, forecast))
mape = np.mean(np.abs((test - forecast) / test)) * 100
print(f"MAE: {mae:.2f}")
print(f"RMSE: {rmse:.2f}")
print(f"MAPE: {mape:.2f}%")
7. Автоматический подбор параметров
# Использование auto_arima для автоматического подбора
from pmdarima import auto_arima
# Автоматический подбор лучших параметров
auto_model = auto_arima(train,
start_p=0, start_q=0,
max_p=10, max_q=10,
d=None, # Автоопределение
seasonal=False, # Несезонная модель
trace=True,
error_action='ignore',
suppress_warnings=True,
stepwise=True)
print(f"Best ARIMA parameters: {auto_model.order}")
print(f"Best AIC: {auto_model.aic:.2f}")
# Прогноз по лучшей модели
auto_forecast = auto_model.predict(n_periods=len(test))
8. Практические сложности и решения:
Проблема 1: Пропущенные значения
# Интерполяция пропущенных значений
sales_interpolated = sales.interpolate(method='time')
# Или forward fill для дневных данных
sales_ffill = sales.fillna(method='ffill')
Проблема 2: Выбросы
# Обнаружение и обработка выбросов
from scipy import stats
z_scores = np.abs(stats.zscore(sales))
outliers = sales[z_scores > 3]
# Замена выбросов медианным значением
sales_clean = sales.mask(z_scores > 3, sales.median())
Проблема 3: Сезонность (требуется SARIMA)
from statsmodels.tsa.statespace.sarimax import SARIMAX
# SARIMA для недельной сезонности (s=7)
seasonal_model = SARIMAX(train,
order=(7, 1, 3),
seasonal_order=(1, 1, 1, 7),
enforce_stationarity=False,
enforce_invertibility=False)
seasonal_fit = seasonal_model.fit(disp=False)
9. Альтернативы ARIMA в моей практике:
- Prophet: Для рядов с сильной сезонностью и праздничными эффектами
- LSTM/нейросети: Для сложных нелинейных зависимостей
- Ensemble методы: Комбинация нескольких моделей
- Скользящее окно: Для online-прогнозирования
10. Ключевые выводы:
- ARIMA хорошо работает для стационарных рядов без сложной сезонности
- Важно правильно определить параметры (p, d, q) через ACF/PACF
- Необходима тщательная диагностика остатков
- Для production систем часто комбинирую ARIMA с другими методами
- Автоматический подбор параметров экономит время, но ручной анализ дает лучшее понимание данных