Time Series: Прогнозирование продаж

Time Series: Прогнозирование продаж

1. Анализ временного ряда (EDA)

Декомпозиция временного ряда:

import pandas as pd
from statsmodels.tsa.seasonal import seasonal_decompose

# Загрузка данных
df = pd.read_csv('sales.csv', parse_dates=['date'], index_col='date')

# Декомпозиция
decomposition = seasonal_decompose(df['sales'], model='additive', period=12)
trend = decomposition.trend
seasonal = decomposition.seasonal
residual = decomposition.resid

decomposition.plot()

Компоненты:

• Тренд (Trend): долгосрочное движение вверх/вниз
• Сезонность (Seasonality): повторяющиеся паттерны (месячные, годовые)
• Остатки (Residuals): шум, нерегулярные колебания

Статистический анализ:

• Стационарность: проверка ADF-тестом (Augmented Dickey-Fuller)
• Если не стационарный: применить differencing

from statsmodels.tsa.stattools import adfuller

result = adfuller(df['sales'])
print(f'ADF статистика: {result[0]:.4f}')
print(f'p-value: {result[1]:.4f}')

if result[1] > 0.05:
    print('Ряд не стационарный - нужна дифференциация')
    df['sales_diff'] = df['sales'].diff()

Выявление сезонности:

import matplotlib.pyplot as plt
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf

# ACF (автокорреляция) и PACF (частичная автокорреляция)
fig, (ax1, ax2) = plt.subplots(2, 1)
plot_acf(df['sales'], lags=40, ax=ax1)
plot_pacf(df['sales'], lags=40, ax=ax2)
plt.show()

Пики в ACF/PACF указывают на сезонность (e.g., на отставании 12 месяцев).

2. Модели прогнозирования

Вариант A: ARIMA (AutoRegressive Integrated Moving Average)

Параметры (p, d, q):

• p — порядок AR (авторегрессия)
• d — порядок differencing (интеграция)
• q — порядок MA (скользящее среднее)

from statsmodels.tsa.arima.model import ARIMA

# Автоматический подбор параметров
from pmdarima import auto_arima

auto_model = auto_arima(df['sales'], seasonal=True, m=12,
                         stepwise=True, trace=True)
print(auto_model.summary())

# Обучение
model = ARIMA(df['sales'], order=(1, 1, 1), seasonal_order=(1, 1, 1, 12))
fitted = model.fit()

# Прогноз
forecast = fitted.get_forecast(steps=12)
forecast_df = forecast.conf_int()
print(forecast_df)

Плюсы: интерпретируемость, быстрота, мало данных нужно Минусы: предполагает линейность, требует стационарности

Вариант B: Prophet (Facebook)

from prophet import Prophet

# Подготовка данных
df_prophet = df.reset_index()
df_prophet.columns = ['ds', 'y']

# Создание модели
model = Prophet(
    growth='linear',  # или 'logistic'
    yearly_seasonality=True,
    weekly_seasonality=False,
    daily_seasonality=False,
    interval_width=0.95
)

# Добавление праздников
holidays = pd.DataFrame({
    'holiday': 'New Year',
    'ds': pd.to_datetime(['2023-01-01', '2024-01-01']),
    'lower_window': 0,
    'upper_window': 1,
})
model.add_seasonality(name='monthly', period=30, fourier_order=5)
model.add_country_holidays(country_name='RU')

# Обучение
model.fit(df_prophet)

# Прогноз
future = model.make_future_dataframe(periods=12, freq='MS')
forecast = model.predict(future)
forecast[['ds', 'yhat', 'yhat_lower', 'yhat_upper']].tail(12)

# Визуализация
model.plot(forecast)

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

Вариант C: LSTM (Long Short-Term Memory)

import numpy as np
from tensorflow.keras.models import Sequential
from tensorflow.keras.layers import LSTM, Dense, Dropout
from sklearn.preprocessing import MinMaxScaler

# Подготовка данных
scaler = MinMaxScaler(feature_range=(0, 1))
scaled_data = scaler.fit_transform(df[['sales']])

# Создание sequence данных
lookback = 12  # используем 12 месяцев для прогноза следующего
X, y = [], []
for i in range(lookback, len(scaled_data)):
    X.append(scaled_data[i-lookback:i, 0])
    y.append(scaled_data[i, 0])

X, y = np.array(X), np.array(y)
X = X.reshape(X.shape[0], X.shape[1], 1)

# Обучение
model = Sequential([
    LSTM(50, activation='relu', input_shape=(lookback, 1)),
    Dropout(0.2),
    LSTM(50, activation='relu'),
    Dropout(0.2),
    Dense(25, activation='relu'),
    Dense(1)
])

model.compile(optimizer='adam', loss='mse')
model.fit(X, y, epochs=100, batch_size=32, validation_split=0.2)

# Прогноз
last_sequence = scaled_data[-lookback:].reshape(1, lookback, 1)
forecast_scaled = model.predict(last_sequence)
forecast = scaler.inverse_transform(forecast_scaled)

Плюсы: учитывает сложные нелинейные паттерны, гибкая архитектура Минусы: требует много данных (>1000 точек), нужен GPU, сложнее интерпретировать

Вариант D: Ensemble методы

Комбинирование предсказаний нескольких моделей:

# Простое усреднение
ensemble_forecast = (arima_pred + prophet_pred + lstm_pred) / 3

# Взвешенное усреднение (по исторической точности)
ensemble_forecast = 0.4*arima_pred + 0.3*prophet_pred + 0.3*lstm_pred

Рекомендация по выбору:

1. Быстрый старт (< 1 месяца): ARIMA или Prophet
2. Много данных (> 2 лет) и нужна максимальная точность: LSTM
3. Production система: Ensemble (Prophet + ARIMA)

3. Учет внешних факторов

Праздники и событ:

# В Prophet
holidays = pd.DataFrame({
    'holiday': ['New Year', 'May Day', 'Black Friday'],
    'ds': pd.to_datetime(['2023-01-01', '2023-05-01', '2023-11-24']),
    'lower_window': -1,  # эффект за день до
    'upper_window': 1,   # эффект за день после
})
model.add_country_holidays(country_name='RU')

# В LSTM — добавить как экзогенная переменная
X_exog = holidays_binary + promo_flags  # one-hot encoding

Акции и промоции:

# Бинарная переменная для периодов с акциями
df['is_promo'] = 0
df.loc[df['date'].isin(promo_dates), 'is_promo'] = 1

# Для ARIMA (ARIMAX) — добавить в модель
from statsmodels.tsa.arima.model import ARIMA

model = ARIMA(df['sales'], 
              exog=df[['is_promo', 'price', 'competitor_price']],
              order=(1, 1, 1))

Экзогенные переменные:

# Погода, конкурентные цены, макроэкономические индексы
exog_vars = df[['temperature', 'competitor_sales', 'oil_price']]

# Нормализация
from sklearn.preprocessing import StandardScaler
scaler_exog = StandardScaler()
exog_scaled = scaler_exog.fit_transform(exog_vars)

# LSTM с экзогенными переменными
X_combined = np.concatenate([X, exog_scaled[:-1]], axis=1)

4. Оценка качества прогноза

Метрики:

• MAE (Mean Absolute Error) = среднее абсолютное отклонение
  • Интерпретируемо (те же единицы, что и исходные данные)

from sklearn.metrics import mean_absolute_error
mae = mean_absolute_error(y_test, forecast)

• RMSE (Root Mean Squared Error) = корень из среднеквадратической ошибки
  • Штрафует большие ошибки сильнее

from sklearn.metrics import mean_squared_error
rmse = np.sqrt(mean_squared_error(y_test, forecast))

• MAPE (Mean Absolute Percentage Error) = средняя абсолютная процентная ошибка
  • Нормализована, легко сравнивать разные ряды

mape = np.mean(np.abs((y_test - forecast) / y_test)) * 100

• SMAPE (Symmetric MAPE) — более справедливая версия MAPE

smape = 2.0 * np.mean(np.abs(forecast - y_test) / (np.abs(y_test) + np.abs(forecast))) * 100

Backtesting (проверка на исторических данных):

# Разбиение на train/test
train_size = int(len(df) * 0.8)
train, test = df[:train_size], df[train_size:]

# Обучение на train
model.fit(train)

# Предсказание на test периоде
forecast = model.predict(len(test))

# Расчет ошибки
mae = mean_absolute_error(test['sales'], forecast)
rmse = np.sqrt(mean_squared_error(test['sales'], forecast))
mape = np.mean(np.abs((test['sales'] - forecast) / test['sales'])) * 100

print(f'MAE: {mae:.2f}')
print(f'RMSE: {rmse:.2f}')
print(f'MAPE: {mape:.2f}%')

Cross-validation для временных рядов:

# Time series split (не стандартный cross-val!)
from sklearn.model_selection import TimeSeriesSplit

tscv = TimeSeriesSplit(n_splits=5)
scores = []

for train_idx, test_idx in tscv.split(df):
    train, test = df.iloc[train_idx], df.iloc[test_idx]
    model.fit(train)
    pred = model.predict(len(test))
    mape = np.mean(np.abs((test['sales'] - pred) / test['sales'])) * 100
    scores.append(mape)

print(f'Средний MAPE: {np.mean(scores):.2f}% +/- {np.std(scores):.2f}%')

Проверка стабильности:

• Сравнить точность на разных периодах
• Если MAPE сильно различается → модель нестабильна
• Добавить регуляризацию или пересмотреть features

Доверительные интервалы:

# Prophet автоматически возвращает доверительные интервалы
forecast[['ds', 'yhat', 'yhat_lower', 'yhat_upper']]

# ARIMA
results = model.get_forecast(steps=12)
conf_int = results.conf_int()  # 95% доверительный интервал

Итоговая рекомендация

Roadmap:

1. Месяц 1: ARIMA/Prophet — baseline (MAPE ~15-20%)
2. Месяц 2: добавить экзогенные переменные (MAPE ~10-15%)
3. Месяц 3: LSTM ансамбль (MAPE ~8-12%)
4. Production: Prophet + ARIMA ensemble с автоматическим переобучением каждый день

Мониторинг:

• Отслеживать MAPE в production
• Алерт, если MAPE > 25%
• Переобучение модели 1 раз в неделю
• A/B тест: compare прогнозов vs actual спросом для оптимизации inventory