基于MXNet的时序预测：概率预测与多层感知机实践-优快云博客

本文链接：https://blog.youkuaiyun.com/gitblog_00797/article/details/148862833

基于MXNet的时序预测：概率预测与多层感知机实践

mxnet-the-straight-dope An interactive book on deep learning. Much easy, so MXNet. Wow. [Straight Dope is growing up] ---> Much of this content has been incorporated into the new Dive into Deep Learning Book available at https://d2l.ai/. 项目地址: https://gitcode.com/gh_mirrors/mx/mxnet-the-straight-dope

前言

时间序列预测是机器学习领域的重要应用场景。本文将介绍如何使用MXNet框架构建概率预测模型和多层感知机(MLP)模型进行时序预测。我们将从理论基础到实践应用，完整展示一个时序预测项目的实现过程。

最大似然估计基础

从线性模型到概率模型

传统线性回归模型通常表示为： $Y_t = X_tβ + ε_t$

其中参数β通过最小化均方误差(MSE)来估计。当我们假设误差项ε_t服从高斯分布时，可以将模型重新表述为条件概率分布：

$Y_t ∼ N(X_tβ, σ^2)$

这意味着y_t服从均值为X_tβ、方差为σ²的高斯分布。

似然函数与最大似然估计

似然函数定义为： $L(β; y_t, X_t, σ) = Pr(y_t | X_t, β, σ)$

最大似然估计(MLE)就是找到使似然函数最大化的参数值(β̂,σ̂)。实践中，我们通常最小化负对数似然函数。

MXNet中的实现

在MXNet中实现MLE需要定义对数似然损失函数。例如伯努利分布的负对数似然可以这样实现：

class Bernoulli(Loss):
    def hybrid_forward(self, F, output, label, sample_weight=None):
        label = _reshape_label_as_output(F, output, label)
        loss = -1 * (1 - output) * (1-label) - output * label
        return F.mean(loss, axis=self._batch_axis, exclude=True)

对于高斯分布，MLE等价于最小化L2损失，因此可以直接使用MXNet内置的L2损失函数。

广义线性模型与计数预测

泊松分布模型

当预测计数数据时，泊松分布是合适的选择。泊松随机变量取值为0,1,2,...，其概率质量函数为：

$Pr(y_t=μ) = \frac{e^{-λ_t}λ_t^μ}{μ!}$

其中λ_t > 0既是均值也是方差。

链接函数

直接将λ_t建模为X_tβ的线性组合可能导致负值问题。因此我们使用链接函数：

$λ = f(Xβ)$

常用选择包括：

对数线性：logλ_t = X_tβ
逻辑链接：f(x) = log(1 + exp(x))

我们选择后者以避免数值溢出问题。

实现泊松对数似然

class Poisson(Loss):
    def hybrid_forward(self, F, output, label, sample_weight=None):
        label = _reshape_label_as_output(F, output, label)
        loss = logistic(output) - F.log(logistic(output)) * label
        return F.mean(loss, axis=self._batch_axis, exclude=True)

实战：海平面高度预测

数据准备

我们使用CSIRO提供的海平面高度数据，包含50个月的观测值：

df = pd.read_csv("sea-level-data.csv").set_index('Time')
ts = df.values[-50:,0]

数据呈现明显上升趋势，因此模型将包含：

线性趋势项
目标变量的滞后项

特征工程

# 标准化特征
features_mean = features.mean(axis=0)
features_std = nd.array(features.asnumpy().std(axis=0)).reshape((1,1))
features = (features - features_mean) / features_std

模型构建

使用单层全连接网络：

net = gluon.nn.Sequential()
with net.name_scope():
    net.add(gluon.nn.Dense(1))
    
# 初始化参数
net.collect_params().initialize(mx.init.Normal(sigma=0.1), ctx=ctx)

训练过程

for e in range(epochs):
    for i, (data, label) in enumerate(train_data):
        with autograd.record():
            output = net(data)
            loss = poisson_log_lik(output, label)
        loss.backward()
        trainer.step(batch_size)

预测实现

由于模型是自回归的，我们需要递归预测：

def forecast_poisson(net, last_obs, features_mean, features_std,
                     forecast_length, train_length):
    forecast = nd.empty((forecast_length,1))
    for t in range(forecast_length):
        # 构建特征
        trend = nd.array([train_length + t - 1]).reshape((1,1))
        fct_feat = nd.concat(trend,prev_obs,dim=1)
        # 预测
        fct = net(fct_feat)
        forecast[t,] = fct[0][0]
        prev_obs = inverse_logistic(fct)
    return forecast