逻辑回归分析二元目标变量与其预测变量之间的关系,以估计因变量取值1的概率。在存在时间数据的情况下,随着时间的推移,观测值不是独立的,模型的误差会随着时间产生相关性,而引入自回归特征或滞后项可以捕捉时间依赖性并增强逻辑回归的预测能力。


NHITS 的输入包括静态外部变量 x(s)\mathbf{x}^{(s)}、历史外部变量 x[:t](h)\mathbf{x}^{(h)}_{[:t]}、预测时可用的外部变量 x[:t+H](f)\mathbf{x}^{(f)}_{[:t+H]} 和自回归特征 y[:t]\mathbf{y}_{[:t]},其中每个输入都进一步分解为分类变量和连续变量。该网络使用多分位数回归来模拟以下条件概率:P(y[t+1:t+H]  y[:t],  x[:t](h),  x[:t+H](f),  x(s))\mathbb{P}(\mathbf{y}_{[t+1:t+H]}|\;\mathbf{y}_{[:t]},\; \mathbf{x}^{(h)}_{[:t]},\; \mathbf{x}^{(f)}_{[:t+H]},\; \mathbf{x}^{(s)})

在本 Notebook 中,我们将展示如何使用 NeuralForecast 方法进行二元序列回归。我们将: - 安装 NeuralForecast。- 加载二元序列数据。- 拟合和预测时间分类器。- 绘制和评估预测结果。

您可以使用 Google Colab 的 GPU 运行这些实验。

1. 安装 NeuralForecast

!pip install neuralforecast
import numpy as np
import pandas as pd
from sklearn import datasets

import matplotlib.pyplot as plt
from neuralforecast import NeuralForecast
from neuralforecast.models import MLP, NHITS, LSTM
from neuralforecast.losses.pytorch import DistributionLoss, Accuracy

2. 加载二元序列数据

core.NeuralForecast 类包含共享的 fitpredict 和其他方法,它们接受列为 ['unique_id', 'ds', 'y'] 的 pandas DataFrames 作为输入,其中 unique_id 标识数据集中的各个时间序列,ds 是日期,而 y 是目标二元变量。

在此示例中,我们将 8x8 数字图像转换为 64 长度的序列,并定义一个分类问题,以识别像素何时超过特定阈值。我们声明一个长格式的 pandas dataframe,以匹配 NeuralForecast 的输入。

digits = datasets.load_digits()
images = digits.images[:100]

plt.imshow(images[0,:,:], cmap=plt.cm.gray, 
           vmax=16, interpolation="nearest")

pixels = np.reshape(images, (len(images), 64))
ytarget = (pixels > 10) * 1

fig, ax1 = plt.subplots()
ax2 = ax1.twinx()
ax1.plot(pixels[10])
ax2.plot(ytarget[10], color='purple')
ax1.set_xlabel('Pixel index')
ax1.set_ylabel('Pixel value')
ax2.set_ylabel('Pixel threshold', color='purple')
plt.grid()
plt.show()

# We flat the images and create an input dataframe
# with 'unique_id' series identifier and 'ds' time stamp identifier.
Y_df = pd.DataFrame.from_dict({
            'unique_id': np.repeat(np.arange(100), 64),
            'ds': np.tile(np.arange(64)+1910, 100),
            'y': ytarget.flatten(), 'pixels': pixels.flatten()})
Y_df
unique_iddsypixels
00191000.0
10191100.0
20191205.0
301913113.0
40191409.0
6395991969114.0
6396991970116.0
639799197103.0
639899197200.0
639999197300.0

3. 拟合和预测时间分类器

拟合模型

使用 NeuralForecast.fit 方法,您可以根据数据集训练一组模型。您可以定义预测 horizon(本例中为 12),并修改模型的超参数。例如,对于 NHITS,我们更改了编码器和解码器的默认隐藏大小。

请参阅 NHITSMLP 模型文档

警告

目前,基于循环的模型族尚不支持使用 Bernoulli 分布输出。这影响以下方法 LSTMGRUDilatedRNNTCN。此功能正在开发中。

# %%capture
horizon = 12

# Try different hyperparmeters to improve accuracy.
models = [MLP(h=horizon,                           # Forecast horizon
              input_size=2 * horizon,              # Length of input sequence
              loss=DistributionLoss('Bernoulli'),  # Binary classification loss
              valid_loss=Accuracy(),               # Accuracy validation signal
              max_steps=500,                       # Number of steps to train
              scaler_type='standard',              # Type of scaler to normalize data
              hidden_size=64,                      # Defines the size of the hidden state of the LSTM
              #early_stop_patience_steps=2,         # Early stopping regularization patience
              val_check_steps=10,                  # Frequency of validation signal (affects early stopping)
              ),
          NHITS(h=horizon,                          # Forecast horizon
                input_size=2 * horizon,             # Length of input sequence
                loss=DistributionLoss('Bernoulli'), # Binary classification loss
                valid_loss=Accuracy(),              # Accuracy validation signal                
                max_steps=500,                      # Number of steps to train
                n_freq_downsample=[2, 1, 1],        # Downsampling factors for each stack output
                #early_stop_patience_steps=2,        # Early stopping regularization patience
                val_check_steps=10,                 # Frequency of validation signal (affects early stopping)
                interpolation_mode="nearest",
                )             
          ]
nf = NeuralForecast(models=models, freq=1)
Y_hat_df = nf.cross_validation(df=Y_df, n_windows=1)
# By default NeuralForecast produces forecast intervals
# In this case the lo-x and high-x levels represent the 
# low and high bounds of the prediction accumulating x% probability
Y_hat_df
unique_iddscutoffMLPMLP-medianMLP-lo-90MLP-lo-80MLP-hi-80MLP-hi-90NHITSNHITS-medianNHITS-lo-90NHITS-lo-80NHITS-hi-80NHITS-hi-90y
00196219610.1730.00.00.01.01.00.7611.00.00.01.01.00
10196319610.7841.00.00.01.01.00.5711.00.00.01.01.01
20196419610.0420.00.00.00.00.00.0090.00.00.00.00.00
30196519610.0720.00.00.00.01.00.0540.00.00.00.01.00
40196619610.0590.00.00.00.01.00.0000.00.00.00.00.00
119599196919610.5511.00.00.01.01.00.6971.00.00.01.01.01
119699197019610.6621.00.00.01.01.00.4650.00.00.01.01.01
119799197119610.3690.00.00.01.01.00.3820.00.00.01.01.00
119899197219610.0560.00.00.00.01.00.0000.00.00.00.00.00
119999197319610.0000.00.00.00.00.00.0000.00.00.00.00.00
# Define classification threshold for final predictions
# If (prob > threshold) -> 1
Y_hat_df['NHITS'] = (Y_hat_df['NHITS'] > 0.5) * 1
Y_hat_df['MLP'] = (Y_hat_df['MLP'] > 0.5) * 1
Y_hat_df
unique_iddscutoffMLPMLP-medianMLP-lo-90MLP-lo-80MLP-hi-80MLP-hi-90NHITSNHITS-medianNHITS-lo-90NHITS-lo-80NHITS-hi-80NHITS-hi-90y
001962196100.00.00.01.01.011.00.00.01.01.00
101963196111.00.00.01.01.011.00.00.01.01.01
201964196100.00.00.00.00.000.00.00.00.00.00
301965196100.00.00.00.01.000.00.00.00.01.00
401966196100.00.00.00.01.000.00.00.00.00.00
1195991969196111.00.00.01.01.011.00.00.01.01.01
1196991970196111.00.00.01.01.000.00.00.01.01.01
1197991971196100.00.00.01.01.000.00.00.01.01.00
1198991972196100.00.00.00.01.000.00.00.00.00.00
1199991973196100.00.00.00.00.000.00.00.00.00.00

4. 绘制和评估预测结果

最后,我们绘制了两个模型的预测结果与实际值的对比。并评估 MLPNHITS 时间分类器的准确性。

plot_df = Y_hat_df[Y_hat_df.unique_id==10]

fig, ax = plt.subplots(1, 1, figsize = (20, 7))
plt.plot(plot_df.ds, plot_df.y, label='target signal')
plt.plot(plot_df.ds, plot_df['MLP'] * 1.1, label='MLP prediction')
plt.plot(plot_df.ds, plot_df['NHITS'] * .9, label='NHITS prediction')
ax.set_title('Binary Sequence Forecast', fontsize=22)
ax.set_ylabel('Pixel Threshold and Prediction', fontsize=20)
ax.set_xlabel('Timestamp [t]', fontsize=20)
ax.legend(prop={'size': 15})
ax.grid()

def accuracy(y, y_hat):
    return np.mean(y==y_hat)

mlp_acc = accuracy(y=Y_hat_df['y'], y_hat=Y_hat_df['MLP'])
nhits_acc = accuracy(y=Y_hat_df['y'], y_hat=Y_hat_df['NHITS'])

print(f'MLP Accuracy: {mlp_acc:.1%}')
print(f'NHITS Accuracy: {nhits_acc:.1%}')
MLP Accuracy: 77.8%
NHITS Accuracy: 74.5%

参考文献