Помимо автономного обучения и тестирования прогнозов

Существует множество онлайн-ресурсов об использовании машинного обучения для прогнозирования. Тем не менее, они сосредоточены на обучении в автономном режиме и тестировании прогнозов.

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

Введение

Ресурсы по прогнозированию часто упускают из виду применение моделей для оперативных прогнозов.

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

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

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

Давайте посмотрим, как мы можем построить модель и использовать ее для прогнозирования в реальном времени.

Практический пример — прогнозирование высоты волн

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

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

Этот временной ряд постоянно пополняется новыми наблюдениями. Таким образом, это идеальный пример для разработки модели прогнозирования для прогнозов в реальном времени.

Мы рассмотрим следующие шаги:

  1. Получение исторических данных;
  2. Предварительная обработка и проектирование признаков;
  3. Выбор и построение модели прогнозирования;
  4. Получение последних наблюдений и составление прогнозов.

Получение исторических данных

Начнем с получения данных.

Вы можете прочитать его непосредственно с сервера ERDAP следующим образом:

import pandas as pd

START_DATE = '2022-01-01'
URL = f'https://erddap.marine.ie/erddap/tabledap/IWaveBNetwork.csv?time%2CSignificantWaveHeight&time%3E={START_DATE}T00%3A00%3A00Z&station_id=%22AMETS%20Berth%20B%20Wave%20Buoy%22'

def reading_data(url: str) -> pd.Series:
    """
    Reading ERDAP data

    :param url: ERDAP url as string
    :return: hourly wave height time series as pd.Series
    """

    # reading data directly from erdap
    data = pd.read_csv(url, skiprows=[1], parse_dates=['time'])

    # setting time to index and getting the target series
    series = data.set_index('time')['SignificantWaveHeight']

    # transforming data to hourly and from centimeters to meters
    series_hourly = series.resample('H').mean() / 100

    return series_hourly

series = reading_data(URL)

Вот график временного ряда:

Предварительная обработка и проектирование функций

Перед моделированием может потребоваться предварительная обработка.

Для наглядности сделаем две вещи:

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

Вот код для обеих этих операций:

import numpy as np


class LogTransformation:
    """
    Log transformation and inverse transformation

    Taking the log helps stabilize the variance
    """

    @staticmethod
    def transform(x):
        xt = np.sign(x) * np.log(np.abs(x) + 1)

        return xt

    @staticmethod
    def inverse_transform(xt):
        x = np.sign(xt) * (np.exp(np.abs(xt)) - 1)

        return x

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

def feature_engineering(X: pd.DataFrame) -> pd.DataFrame:
    """
    param X: lagged observations (explanatory variables)

    :return: new features
    """

    summary_stats = {'mean': np.mean, 'sdev': np.std}

    features = {}
    for f in summary_stats:
        features[f] = X.apply(lambda x: summary_stats[f](x), axis=1)

    features_df = pd.concat(features, axis=1)
    X_feats = pd.concat([X, features_df], axis=1)

    return X_feats

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

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

Построение модели прогнозирования

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

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

from sklearn.model_selection import train_test_split
# https://github.com/vcerqueira/blog
from src.tde import time_delay_embedding

# using last 24 observations as lags, 
# and next 24 observations as the forecasting horizon
N_LAGS, HORIZON = 24, 24

train, test = train_test_split(series, test_size=0.2, shuffle=False)

X_train, Y_train = time_delay_embedding(train, n_lags=N_LAGS, horizon=HORIZON, return_Xy=True)
X_test, Y_test = time_delay_embedding(test, n_lags=N_LAGS, horizon=HORIZON, return_Xy=True)

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

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

from sklearn.model_selection import RandomizedSearchCV, TimeSeriesSplit

## apply preprocessing steps
# log transformation
X_train = LogTransformation.transform(X_train)
Y_train = LogTransformation.transform(Y_train)
# feature engineering
X_train_ext = feature_engineering(X_train)

# time series cv procedure
tscv = TimeSeriesSplit(n_splits=5, gap=50)

# defining the search space
# a simple optimization of the number of trees of a RF
model = RandomForestRegressor()
param_search = {'n_estimators': [10, 50, 100, 200],
                'criterion': ['squared_error', 'absolute_error'],
                'max_depth': [None, 2, 5, 10],
                'max_features': ['log2', 'sqrt']}

# applying CV with random search on the training data
gs = RandomizedSearchCV(estimator=model,
                        cv=tscv,
                        refit=True,
                        param_distributions=param_search,
                        n_iter=10, n_jobs=1)

gs.fit(X_train_ext, Y_train)

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

# applying preprocessing steps to test data
X_test = LogTransformation.transform(X_test)
X_test_ext = feature_engineering(X_test)

# inference on test set and evaluation
preds = gs.predict(X_test_ext)

# log forecasts
preds_log = gs.predict(X_test_ext)

# reverting the log transformation
preds = LogTransformation.inverse_transform(preds_log)

# estimating performance using r-squared
estimated_performance = r2_score(Y_test, preds)

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

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

from joblib import dump

# preparing all available data for auto-regression
X, Y = time_delay_embedding(series, n_lags=N_LAGS, horizon=HORIZON, return_Xy=True)

# applying preprocessing steps
X = LogTransformation.transform(X)
Y = LogTransformation.transform(Y)
X_ext = feature_engineering(X)

# model fitting
final_model = RandomForestRegressor(**gs.best_params_)
final_model.fit(X_ext, Y)

dump(final_model, 'random_forest_v1.joblib')

Применение модели

На этом этапе мы сделали несколько вещей:

  • Прочитать и предварительно обработать временной ряд;
  • Построена и оптимизирована модель прогнозирования с использованием перекрестной проверки;
  • Оценил его производительность с помощью тестового набора.

Теперь мы готовы применить эту модель в реальных условиях.

Начнем с получения последних наблюдений с буя.

import datetime

# setting the max history to yesterday
yesterday = datetime.date.today() - datetime.timedelta(days=1)
yesterday = yesterday.strftime('%Y-%m-%d')

LIVE_URL = f'https://erddap.marine.ie/erddap/tabledap/IWaveBNetwork.csv?time%2CSignificantWaveHeight&time%3E={yesterday}T00%3A00%3A00Z&station_id=%22AMETS%20Berth%20B%20Wave%20Buoy%22'

# reading the data from the ERDAP server
new_series = reading_data(LIVE_URL)

# getting the last 24 observations needed for the model
lags = new_series.tail(N_LAGS)

Зачем нужны эти наблюдения?

Напомним, что наша модель основана на авторегрессии. Это означает, что он использует недавние прошлые наблюдения для предсказания будущих. В нашем примере мы используем последние 24 наблюдения для прогнозирования данных на следующие 24 часа. Таким образом, входные данные для модели основаны на последних 24 недавних наблюдениях.

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

from joblib import load

# structuring the lags as a DataFrame
lags_df = pd.DataFrame(lags).T.reset_index(drop=True)
lags_df.columns = X.columns

# applying preprocessing steps
lags_df = LogTransformation.transform(lags_df)
lags_feats = feature_engineering(lags_df)

# loading the model from disk
final_model = load('random_forest_v1.joblib')

# applying the model
log_forecasts = final_model.predict(lags_feats)

# reverting the log transformation
forecasts = LogTransformation.inverse_transform(log_forecasts)

Вот как выглядят прогнозы:

Прогнозы отдельных деревьев также включены, чтобы передать неопределенность прогнозов.

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

Ключевые выводы

В этой статье мы построили модель и использовали ее для составления прогнозов в реальном сценарии.

Мы исследовали:

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

Спасибо за прочтение и до встречи в следующей истории!

Статьи по Теме

Рекомендации

[1] Irish Wave Buoys от Морского института (ID набора данных: IWaveBNetwork). Лицензия: CC BY 4.0