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

среда, 20 октября 2021 г.

Прогнозирование временных рядов в R с помощью lstm

Оказывается, глубокое обучение со всей его мощью также можно использовать для прогнозирования. Особенно модель LSTM (Long Short Term Memory), которая оказалась полезной при решении проблем, связанных с последовательностями с автокорреляцией.

Сети с долговременной краткосрочной памятью - это своего рода рекуррентные нейронные сети (RNN), которые способны изучать долгосрочные зависимости. LSTM позволяет сохранять долгосрочные состояния в дополнение к краткосрочным, с которыми традиционные RNN испытывают трудности. LSTM весьма полезны в задачах прогнозирования временных рядов, включающих автокорреляцию, из-за их способности поддерживать состояние и распознавать закономерности по длине ряда.

Здесь я показываю, как реализовать прогнозирующую модель LSTM с помощью языка R.

Модель LSTM доступна в пакете keras R, который работает поверх Tensorflow.

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

1
2
3
4
library(keras)
library(tensorflow)
install_keras()
install_tensorflow(version = "nightly")

Подготовка данных

Для этого примера я использовал набор данных economics, который можно найти в пакете ggplot2. Я хочу спрогнозировать безработицу на следующие 12 месяцев (столбец unemploy).

1
2
library(ggplot2)
head(economics)



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

scale_factors <- c(mean(economics$unemploy), sd(economics$unemploy))


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

1
2
3
scaled_train <- economics %>%
    dplyr::select(unemploy) %>%
    dplyr::mutate(unemploy = (unemploy - scale_factors[1]) / scale_factors[2])

Алгоритм LSTM создает прогнозы на основе запаздывающих значений. Это означает, что ему нужно оглянуться на столько предыдущих значений, сколько точек мы хотим предсказать. Поскольку мы хотим сделать прогноз на 12 месяцев, нам нужно основывать каждый прогноз на 12 точках данных.

Для демонстрации предположим, что наш ряд имеет 10 точек данных [1, 2, 3,…, 10], и мы хотим предсказать 3 значения. Тогда наши обучающие данные выглядят так:

1
2
3
4
5
[1, 2, 3] -> [4, 5, 6]
[2, 3, 4] -> [5, 6, 7]
[3, 4, 5] -> [6, 7, 8]
[4, 5, 6] -> [7, 8, 9]
[5, 6, 7] -> [8, 9, 10]

Предикторы и целевые значения имеют форму:

1
2
3
      3 4 5 6 7                6 7 8 9 10
X = [ 2 3 4 5 6 ]        Y = [ 5 6 7 8 9  ]
      1 2 3 4 5                4 5 6 7 8


Кроме того, keras LSTM ожидает определенного тензорного формата трехмерного массива [samples, timesteps, features] для предикторов (X) и для целевых (Y) значений:
  • samples указывает количество наблюдений, которые будут обрабатываться партиями.
  • timesteps сообщает нам количество временных шагов (лагов). Или, другими словами, сколько единиц назад во времени мы хотим, чтобы увидела наша сеть.
  • features указывает количество предикторов (1 для одномерного ряда и n для многомерного).
В случае предикторов, которые преобразуются в массив измерений: (nrow (data) - lag – prediction + 1, 12, 1), lag = prediction = 12.

1
2
prediction <- 12
lag <- prediction

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

3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
scaled_train <- as.matrix(scaled_train)
 
# we lag the data 11 times and arrange that into columns
x_train_data <- t(sapply(
    1:(length(scaled_train) - lag - prediction + 1),
    function(x) scaled_train[x:(x + lag - 1), 1]
  ))
 
# now we transform it into 3D form
x_train_arr <- array(
    data = as.numeric(unlist(x_train_data)),
    dim = c(
        nrow(x_train_data),
        lag,
        1
    )
)

По сути, каждый столбец является отстающей версией предыдущего - последний отстает на 11 шагов по сравнению с первым.


И финальная версия:



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

Теперь применим аналогичное преобразование для значений Y.

1
2
3
4
5
6
7
8
9
10
11
12
13
y_train_data <- t(sapply(
    (1 + lag):(length(scaled_train) - prediction + 1),
    function(x) scaled_train[x:(x + prediction - 1)]
))
 
y_train_arr <- array(
    data = as.numeric(unlist(y_train_data)),
    dim = c(
        nrow(y_train_data),
        prediction,
        1
    )
)


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

x_test<-economics$unemploy[(nrow(scaled_train)-prediction+1):nrow(scaled_train)]


Нам нужно его масштабировать и преобразовать.

4
5
6
7
8
9
10
11
12
# scale the data with same scaling factors as for training
x_test_scaled <- (x_test - scale_factors[1]) / scale_factors[2]
 
# this time our array just has one sample, as we intend to perform one 12-months prediction
x_pred_arr <- array(
    data = x_test_scaled,
    dim = c(
        1,
        lag,
        1
    )
)



Прогноз lstm

Мы можем построить модель LSTM, используя функцию keras_model_sequential и добавляя к ней слои. Первый слой LSTMer принимает требуемую форму ввода, которая имеет вид [samples, timesteps, features]. Мы устанавливаем для обоих слоев return_sequences = TRUE и stateful = TRUE. Второй уровень такой же, за исключением batch_input_shape, который нужно указать только на первом слое.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
lstm_model <- keras_model_sequential()
 
lstm_model %>%
  layer_lstm(units = 50, # size of the layer
       batch_input_shape = c(1, 12, 1), # batch size, timesteps, features
       return_sequences = TRUE,
       stateful = TRUE) %>%
  # fraction of the units to drop for the linear transformation of the inputs
  layer_dropout(rate = 0.5) %>%
  layer_lstm(units = 50,
        return_sequences = TRUE,
        stateful = TRUE) %>%
  layer_dropout(rate = 0.5) %>%
  time_distributed(keras::layer_dense(units = 1))


Вы также можете попробовать различные функции активации с помощью параметра activation (по умолчанию используется гиперболический тангенс tanh).

Также выберите функцию потерь для оптимизации, тип оптимизатора и метрику для оценки производительности модели. Информацию о различных оптимизаторах можно найти здесь.

1
2
3
4
lstm_model %>%
    compile(loss = 'mae', optimizer = 'adam', metrics = 'accuracy')
 
summary(lstm_model)


Затем мы можем подогнать нашу модель LSTM с отслеживанием состояния. Мы устанавливаем shuffle = FALSE, чтобы сохранить последовательности временных рядов.

3
4
5
6
7
8
lstm_model %>% fit(
    x = x_train_arr,
    y = y_train_arr,
    batch_size = 1,
    epochs = 20,
    verbose = 0,
    shuffle = FALSE
)

И создаем прогноз:

1
2
3
4
5
6
lstm_forecast <- lstm_model %>%
    predict(x_pred_arr, batch_size = 1) %>%
    .[, , 1]
 
# we need to rescale the data to restore the original values
lstm_forecast <- lstm_forecast * scale_factors[2] + scale_factors[1]


Модель LSTM более сложна, чем модели обычных временных рядов, потому что вы не передаете явное количество точек прогноза для прогноза. Вместо этого вам нужно спроектировать свою модель таким образом, чтобы она прогнозировала желаемое количество периодов. Если вы хотите прогнозировать больше периодов, вам необходимо предоставить дополнительные столбцы в вашем наборе prediction, содержащие значения, прогнозируемые для предыдущих периодов. Следуя этому примеру, прогнозирование на 13-й месяц вперед требует «знания» результата на 1-й месяц вперед. Другой вариант - перестроить модель для прогнозирования 13 значений вместо 12.

Объект forecast

Поскольку у нас есть прогнозируемые значения, мы можем превратить результаты в объект forecast, как если бы мы использовали пакет forecast. Это позволит использовать функцию forecast::autoplot для построения результатов прогноза. Для этого нам нужно определить несколько объектов, которые создают объект forecast.

Прогноз на обучающей выборке

1
2
fitted <- predict(lstm_model, x_train_arr, batch_size = 1) %>%
     .[, , 1]



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

1
2
3
4
5
6
7
8
9
if (dim(fitted)[2] > 1) {
    fit <- c(fitted[, 1], fitted[dim(fitted)[1], 2:dim(fitted)[2]])
} else {
    fit <- fitted[, 1]
}
 
# additionally we need to rescale the data
fitted <- fit * scale_factors[2] + scale_factors[1]
nrow(fitted) # 562

В связи с тем, что наш прогноз начинается со смещения 12 месяцев, нам необходимо предоставить искусственные (или реальные) значения для этих месяцев:

1
2
# I specify first forecast values as not available
fitted <- c(rep(NA, lag), fitted)

Прогноз в виде объекта ts

Нам нужно преобразовать предсказанные значения в объект временного ряда.

1
2
3
4
5
lstm_forecast <- timetk::tk_ts(lstm_forecast,
    start = c(2015, 5),
    end = c(2016, 4),
    frequency = 12
)



Входной ряд

Кроме того, нам необходимо преобразовать экономические данные в объект временного ряда.

1
2
3
4
input_ts <- timetk::tk_ts(economics$unemploy,
    start = c(1967, 7),
    end = c(2015, 4),
    frequency = 12)


Объект forecast

Наконец, мы можем определить объект forecast:

1
2
3
4
5
6
7
8
9
10
forecast_list <- list(
    model = NULL,
    method = "LSTM",
    mean = lstm_forecast,
    x = input_ts,
    fitted = fitted,
    residuals = as.numeric(input_ts) - as.numeric(fitted)
  )
 
class(forecast_list) <- "forecast"

Теперь мы легко можем построить график данных

forecast::autoplot(forecast_list)


lstm-прогноз с регрессорами

Работа с регрессорами в LSTM сводится к обработке ряда как многомерного, а не одномерного. Давайте создадим несколько случайных регрессоров для этого примера:

reg<-100*runif(nrow(economics))


Как и в случае с обучающим набором, нам необходимо также масштабировать регрессоры.

1
2
3
4
5
6
7
8
9
10
scale_factors_reg <- list(
     mean = mean(reg),
     sd = sd(reg)
  )
 
scaled_reg <- (reg - scale_factors_reg$mean)/scale_factors_reg$sd
 
# additionally 12 values for the forecast
scaled_reg_prediction <- (reg[(length(reg) -12): length(reg)] -
                      scale_factors_reg$mean) /scale_factors_reg$sd

Теперь мы можем добавить их к обучающим данным и преобразовать в тензоры, как раньше.

1
2
3
4
5
6
7
8
9
10
11
# combine training data with regressors
x_train <- cbind(scaled_train, scaled_reg)
x_train_data <- list()
 
# transform the data into lagged columns
for (i in 1:ncol(x_train)) {
    x_train_data[[i]] <- t(sapply(
        1:(length(x_train[, i]) - lag - prediction + 1),
        function(x) x_train[x:(x + lag - 1), i]
    ))
}

На этот раз в нашем списке осталось 2 записи - входные данные точно такие же, как и раньше, и регрессоры, следуют той же логике.



Снова преобразуем этот список в трехмерный массив.

1
2
3
4
5
6
7
8
x_train_arr <- array(
    data = as.numeric(unlist(x_train_data)),
    dim = c(
        nrow(x_train_data[[1]]),
        lag,
        2
    )
)


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

3
4
5
6
7
8
9
10
11
12
# combine the data with regressors
x_test_data <- c(x_test_scaled, scaled_reg_prediction)
 
# transform to tensor
x_pred_arr <- array(
    data = x_test_data,
    dim = c(
        1,
        lag,
        2
   )
)



Остальная часть моделирования остается прежней.

Комментариев нет:

Отправить комментарий