Корректировки стационарности временного ряда, улучшение точности прогнозирования

Цель - повышение точности прогнозирования. Корректировка стационраности - такое изменение временого хода исследуемого набора данных, при котором его среднее значение не изменяется во времени (Параметры математического ожидания и дисперсии постоянны).

Блок корректировки стационарности временного ряда

#Копирование исходного дф data в дф data_correction для дальнейшего изменения

data_correction = data[data_columns]

#Формирование новой колонки 'trend_residuals' дф data_correction в который заносится

#разность максимального значения годового тренда и каждого отсчета данного тренда

#с целью определения отклонения значений данного временного ряда

(data_365d_rol['Consumption'])

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

#для выравнивания уровней амплитуды исследуемого временного хода

(data['Consumption']).

data_correction['trend_residuals'] = data_365d_rol['Consumption'].max() - data_365d_rol['Consumption']

#Формирование скорректированного временного хода (data['Consumption']). В

колонке 'Consumption_eq'

#(eq - equalized (выровненый)) заносится скорректированный ряд

data_correction['Consumption'],

#коррекция которого выполнена путем прибавления вычесленных отклонения на основе линии тренда.

data_correction['Consumption_eq'] = data_correction['Consumption'] + data_correction['trend_residuals']

#Формирование тренда, в соответствие с которым выполнены поправки амплитуды

#временного хода (data['Consumption']) с целью графической визуализации поправок.

data_correction['trend_eq'] = data_365d_rol['Consumption'] + data_correction['trend_residuals'] * 2

#Фомирование тренда по выровненному временному ходу

data_correction['Consumption_eq']. data_correction['Consumption_eq_roll_7_day'] = data_correction['Consumption_eq'].rolling(window = 7, center = True).mean()

#Формирование дф с постоянным знаяением для гафичекой визуализации. data_correction['constant'] = 1600

#Построение графиков

18

fig, ax = plt.subplots(figsize = (28,10))

# Исходный временной ряд ax.plot(data_correction['Consumption'], marker='.', markersize=3,linestyle='None', label='Consumption')

# Выровненный временной ряд ax.plot(data_correction['Consumption_eq'], marker='.', markersize=2, color='0',linestyle='None', label='Consumption equalized')

# Выровненный тренд с усреднением по 7 дней ax.plot(data_correction['Consumption_eq_roll_7_day'], linewidth=2, label='Trend (7-d Rolling Mean) equalized')

# Исходный тренд с усреднением по 7 дней ax.plot(data_7d_rol['Consumption'], linewidth=2, label='Trend (7-d Rolling Mean)')

# Выровненный тренд с усреднением по 365 дней ax.plot(data_correction['trend_eq'], color='0.2', linewidth=3, label='Trend (365-d Rolling Mean) equalized')

#ax.plot(data_correction['constant'], linewidth=2, label='7-d Rolling Mean') ax.plot(data_correction['constant'], linewidth=2, label='constant')

#Добавление легенды на график

ax.legend()

#Название оси абсцисс ax.set_xlabel('Year')

#Название оси ординат

ax.set_ylabel('Consumption (GWh)')

Text(0, 0.5, 'Consumption (GWh)')

Функция генерирования дополнительных временных признаком

Функция генерирования дополнительных признаком для временного ряда средствами функции .dt, которая позволяет выделять временные диапазоны с различным шагом (час, день, неделя, месяц, квартал, год и т.д.)

def Create_time_predictors(df):

df['datetime'] = df.index # Копирование индексов дф в колонку df['hour'] = df.datetime.dt.hour # Час

df['dayofweek'] = df.datetime.dt.dayofweek # День недели

19

df['quarter'] = df.datetime.dt.quarter # Квартал df['month'] = df.datetime.dt.month # Месяц df['dayofyear'] = df.datetime.dt.dayofyear # День года df['day'] = df.datetime.dt.day # День

df['weekofyear'] = df.datetime.dt.isocalendar().week # Неделя df = df.drop('datetime', axis = 1) # Удаление из дф колонки

# Операция возвращения функцией дф return df

import warnings warnings.filterwarnings('ignore')

#Реализация прогнозирования на основе скорректированного датасета data_correction

#Формирмирование дополнительных временных признаков средствами функции

Create_time_predictors()

data_correction = Create_time_predictors(data_correction)

#Удаление из дф неиспользуемых колонок

data_correction = data_correction.drop(['Consumption', 'Wind', 'Solar', 'Wind+Solar',

'trend_residuals', 'trend_eq', 'Consumption_eq_roll_7_day', 'constant'], axis = 1)

data_correction = data_correction.dropna()

X_train = data_correction[:'2016'].drop(['Consumption_eq'], axis = 1) y_train = data_correction.loc[:'2016', 'Consumption_eq']

X_test = data_correction['2017'].drop(['Consumption_eq'], axis = 1) y_test = data_correction.loc['2017', 'Consumption_eq']

#Объявление алгоритмов МО и аккумуляция их в массив

#для удобной работы в цикле

models = [] # Лист для хранения моделей models.append(('LR', LinearRegression()))

models.append(('NN', MLPRegressor(solver = 'lbfgs'))) # solver - определение ядра принятия решений

models.append(('KNN', KNeighborsRegressor())) models.append(('RF', RandomForestRegressor(n_estimators = 50))) # n_estimators - число деревьев в ансамбле

models.append(('SVR', SVR(gamma='auto', kernel='rbf'))) # kernel - ядро по которому формируется опорный вектор

# Рассчетный цикл

results = [] # Лист для хранения результата оценки модели names = [] # Лист для хранения имени модели

for name, model in models:

# TimeSeriesSplit определяет количество частей на которые будет разбит датасет,

20

#переменная tscv используется как атрибут функции cross_val_score(), т.е.

#определяет количество частей датасета по которым будет проведен

перекрестный

# анализ (валидация результатов по метрике R^2) tscv = TimeSeriesSplit(n_splits=10)

cv_results = cross_val_score(model, X_train, y_train, cv=tscv, scoring='r2')

# Добавление результатов в лист results и names results.append(cv_results)

names.append(name)

# Вывод результатов

print('%s: %f (%f)' % (name, cv_results.mean(), cv_results.std()))

#График сравнения используемых алгоритмов plt.boxplot(results, labels=names) plt.title('Algorithm Comparison') plt.show()

#Цикл выполнения прогноза и сравнения реального временого хода

#со спрогнозированный результатом

for name, model in models:

#Функция .fit() выполняет обучение модели (подгонку)

#на основе тренировочного набора данных

model = model.fit(X_train, y_train)

#Функция .predict() выполняет прогнозирование

#на основе тестового промежутка времени.

y_pred = model.predict(X_test)

#Далее выполняется передача результатов

#прогнозирования y_pred и реального временного

#хода y_test в функцию regression_results(), с целью

#оценки точности итоговой модели.

regression_results(y_test, y_pred, model)

# Построение графика спрогнозированного и реальных временных ходов. fig, ax = plt.subplots(figsize = (18,8))

ax.plot(X_test.index, y_pred, linewidth=2, label='Prediction') ax.plot(X_test.index, y_test, linewidth=2, label='Real') ax.set_title('Data Comparison')

ax.legend() ax.set_xlabel('Year')

ax.set_ylabel('Consumption (GWh)') plt.show()

LR: 0.421818 (0.058571)

NN: 0.291276 (0.132530)

KNN: 0.673180 (0.159450)

RF: 0.851859 (0.020055)

SVR: 0.016285 (0.026550)

21

LinearRegression() r2: 0.4433

MAE: 98.5827 MSE: 15043.5814 RMSE: 122.6523

MLPRegressor(solver='lbfgs') r2: 0.2768

MAE: 112.2256

MSE: 19542.9866

RMSE: 139.7962

22

KNeighborsRegressor() r2: 0.8756

MAE: 39.9038

MSE: 3362.0734

RMSE: 57.9834

RandomForestRegressor(n_estimators=50) r2: 0.9361

MAE: 30.2583

MSE: 1725.9785

RMSE: 41.5449

23