RNN之：LSTM 长短期记忆模型-结构-理论详解-及实战（Matlab向）

0.前言

递归！循环神经网络Recurrent Neural Network

循环神经网络（又称递归神经网络，Recurrent Neural Network，RNN）。是一种用于处理序列数据的神经网络结构，具有记忆功能，能够捕捉序列中的时间依赖关系。RNN通过循环连接的方式，将当前时间步的输入和上一时间步的隐藏状态作为输入，计算当前时间步的隐藏状态和输出，从而实现对序列数据的建模和预测。

说到递归描述，各位肯定会想到马尔科夫链，循环神经网络和马尔科夫链具有一定的相似性，两者都是用来捕捉和描述序列数据的变化规律。

但两者在其原理上有根本区别：

马尔可夫链需要满足马尔可夫状态，即知道t时刻的状态时，t+1时刻的状态的概率分布明确。我们可以通过转移矩阵基于t时刻的信息计算t+1、t+2等未来时刻中不同状态的预期概率。因此未来状态仅受当前状态影响，过去发生过什么并不重要。它常常被用来描述理想状态下物理场中概率变化目标的属性。而循环神经网络则不需要满足任何前提条件，它可以基于过去任何长短时间的数据，基于统计经验给出未来的最可能状态。

但是相较于循环神经网络而言，马尔可夫链在实际应用中具有天生缺陷：

一方面是马尔可夫链的转移矩阵的概率映射是线性的，并不能很好的描述非线性的长迭代状态（比如，当转移矩阵会随时间步长发生变化时）。所以马尔可夫链在自然语言中的处理和应用非常拉跨。对于RNN而言，由于其本身结构的复合性，RNN的状态转移可以满足非线性的变化需求。
同时马尔可夫状态的设定过于理想，现实世界的应用工程中，状态空间本身，转移矩阵的精确求解困难。而对于RNN而言，无论满足满足马尔可夫状态，都可以从递归过程中挖掘统计规律。因此在一些具有统计学特性外的外部附加规则的预测任务中（比如游戏抽卡），RNN能够反应与时间步长相关的特征。

然而，对于传统的RNN，随着序列长度的增加，计算得到的梯度在反向传播过程中会逐渐消失或爆炸，导致模型难以训练。这种现象被称为“梯度消失”或“梯度爆炸”，它限制了RNN捕捉和利用长距离依赖关系的能力。这对于许多需要理解长期依赖关系的任务，如自然语言处理、语音识别等，是一个严重的挑战。

LSTM的诞生

1997年，由德国慕尼黑工业大学的计算机科学家Sepp Hochreiter与Jürgen Schmidhuber共同提出了LSTM（Long Short-Term Memory）模型。LSTM是一种特别的RNN，旨在解决传统RNN在处理长序列数据时遇到的梯度消失和梯度爆炸问题。

原文链接：Long Short-Term Memory | MIT Press Journals & Magazine | IEEE Xplore

原文PDF：nc.dvi

Sepp Hochreiter与Jürgen Schmidhuber在实验中对LSTM和之前其他RNN模型做了进一步对比：利用无噪声和有噪声的序列对不同的RNN进行训练，这些序列长度较长，对于有噪声序列而言，其中只有少数数据是重要的，而其余数据则起到干扰作用。

表：LSTM与RTRL（RealTime Recurrent Learning）、ELM（Elman nets，又称单循环网络）和RCC（Recurrent Cascade-Correlation，级联相关学习架构）在长序列无噪声数据中的表现对比

可以发现，相较于RTRL、ELM、RCC等传统模型，LSTM在长序列无噪声模型中的训练成功率更高，同时训练速度更快。

注：表中LSTM的Block为一个LSTM单元，Block为4代表4层LSTM单元进行堆叠，Size代表一个Block或LSTM单元中memory cell（记忆单元，也就是长期记忆）并列的数量，为了方便计算，这些记忆单元共享输入门和输出门和遗忘门。其实，文中所说的Size其实就是我们现在说的隐藏单元数量，由于1997年计算机技术不是很发达，GPU并行计算的运用不多，所以原文的结构写得很复杂。（具体示例结构如原文图2所示，由于是编码训练的形式，原文展示得很复杂，看不懂没关系，不影响后面理解）。

其中，原文对Size的表述是：Memory cell blocks of size S form a structure where S memory cells share the same input gate and the same output gate. These blocks facilitate information storage. However, as with conventional neural nets, it is not so easy to code a distributed input within a single cell. Since each memory cell block has as many gate units as a single memory cell, namely two, the block architecture can be even slightly more efficient (see paragraph on computational complexity in Section).

表：LSTM与RTRL、BPTT（Back-Propagation Through Time）和RCC在长序列无局部规律性数据中的表现对比

由于长期记忆（细胞记忆）的特性，在无局部规律数据中，LSTM的表现相较于RTRL、BPTT和CH模型在训练速度及跟踪能力上有质的差异，不仅训练迭代需求缩小了十倍以上，而且表现出极好的追踪能力。这也是为什么在那么多RNN模型中，LSTM可以获得成功的原因。

简单来说而言：LSTM（Long Short-Term Memory）是一种特殊的递归神经网络（RNN），旨在解决传统RNN在处理长序列数据时遇到的梯度消失和梯度爆炸问题。LSTM通过引入“记忆单元”和“门控机制”，能够高效建模长期和短期依赖关系。

1.长短期记忆网络 LSTM（Long Short-Term Memory）的基本结构

常规LSTM单元的结构如图1所示，在一些论文中常被为“记忆细胞”（Memory Cell），但在中文语境我更倾向于称之为LSTM单元，这样不容易产生混淆。

这是由于GPU的推广，目前的LSTM设计类似于多个“传统的记忆细胞”并联在一起，组成一个集合体，对应LSTM开山论文中的“Cell block”，来传输多维的记忆信息，而在开山论文中，Memory cell代表只能处理数值而非多维向量的单元，因此在当前时代，Memory cell和Cell block在很多语境下基本上是一个东西，所以这个cell究竟是那个cell？这个描述不是很严谨。

一般而言，这就是，当其在处理数据时，数据内容在内部循环流动，序列数据或任何长度的单维度数据（dim=1）在输入LSTM时，按照不同的时间步长t，依次输入模型，在这一过程中当前的细胞状态（cell sate）和隐藏状态（hidden state）在不断变化。

它包含三个关键的门控机制：输入门（Input Gate）、遗忘门（Forget Gate）和输出门（Output Gate）。其详细信息如下图所示：

图1：LSTM单元经典结构细节（一图流），以下介绍皆与该图对应

图中：

为了简化表达，不同位置weight、Bias代表不同权重乘法和偏置加法操作，尽管名字相同，但不同位置的weight、Bias的参数相互独立
“（）”括号内的数值代表在该环节传输或该参数的“张量”大小，如Weight（N*N）代表权重矩阵的张量规模是N*N。
短期记忆（隐藏状态 $h_t$ ）和长期记忆（细胞状态 $C_t$ ）的张量尺度必须一致，均为 N 维向量，N是可以自己需求设计数值，称为隐藏单元的数量

1.1.遗忘门（Forget Gate）：

遗忘门决定了前一时刻的记忆细胞状态中有多少信息需要保留或遗忘。

它接收当前的输入 $X_t$ 和上一时间步的隐藏状态 $h_{t-1}$ ，通过Sigmoid函数输出一个0到1之间的值，表示遗忘的比例。

计算公式为：

$f_t=sigmoid(W_{f1}\cdot h_{t-1}+W_{f2} \cdot X_t+b_f)$ (1)

其中 $W_{f}$ 和 $b_f$ 是遗忘门的权重和偏置的参数矩阵，为了区分对应隐藏状态和当前输入的权重矩阵相互独立，格外加了数字下表如： $W_{f1}$ ,后同。

1.2.输入门（Input Gate）：

输入门决定了当前时间步的输入信息中有多少需要被加入到细胞状态中。

如一图流所示，它包含两个部分：首先，使用Sigmoid函数决定哪些信息需要更新（Potential Memory to remember)， $i_t$ ；其次，使用Tanh函数生成一个候选长期(Potential Long-Trem menmory),也就是 $\tilde{C}_t$ 。

计算公式为：

$i_t=sigmoid(W_{i1}\cdot h_{t-1}+W_{i2} \cdot X_t+b_{i1})$ (2)

$\tilde{C}_t=tanh(W_{i3}\cdot h_{t-1}+W_{i4} \cdot X_t+b_{i2})$ (3)

1.3.细胞状态-长期记忆（Cell-State）：

Cell-State负责在序列的各个时间步长之间存储和传递信息，简单的来说，它就是长期记忆。它像是一个传输带，在整个序列处理过程中，信息在上面流动，只有少量的线性操作被应用于信息本身，信息流动相对简单。

每次迭代中，细胞状态更新公式为：

$C_t=C_{t-1}\cdot f_t +i_t\cdot \tilde{C}_t$ (4)

根据式（1）、（2）、（3），（4）可完整展开为：

$C_t=C_{t-1}\cdot sigmoid(W_{f1}\cdot h_{t-1}+W_{f2} \cdot X_t+b_f)+sigmoid(W_{i1}\cdot h_{t-1}+W_{i2} \cdot X_t+b_{i1})\cdot tanh(W_{i3}\cdot h_{t-1}+W_{i4} \cdot X_t+b_{i2})$ (5)

其中， $C_t$ 为该时间步输出的细胞状态， $C_{t-1}$ 为上一时间步输出的的细胞状态，其中 $f_t$ 由遗忘门计算而出，决定了上一刻信息保留的程度。 $i_t\cdot \tilde{C}_t$ 是选择门从隐藏状态与当前时间步或取信息产生的。总结而言：遗忘门决定了前一时刻的记忆细胞状态中有多少信息需要保留，输入门决定了当前时间步的输入信息中有多少需要被加入到细胞状态中。

1.4.隐藏状态-短期记忆（Hidden State）与输出门（Output Gate）：

输入门决定了基于当前时间步的细胞状态、输入数据和上一时间步的隐藏状态，输出当前时间步隐藏状态 $h_t$ ,对于序列预测，ht就是LSTM对于下一时间步输入 $X_{t+1}$ 的预测。如果你是以序列输出且以回归损失函数反向传播的。

它接收当前的输入xt和上一时间步的隐藏状态ht−1，通过Sigmoid函数输出一个0到1之间的值，表示输出的比例。

计算公式为：

$h_t=tanh(W_{i3}\cdot h_{t-1}+W_{i4} \cdot X_t+b_{i2})\cdot tanh(C_t)$ （6）

同样可以完整展开为：

$h_t=tanh(W_{i3}\cdot h_{t-1}+W_{i4} \cdot X_t+b_{i2})\cdot tanh(C_{t-1}\cdot sigmoid(W_{f1}\cdot h_{t-1}+W_{f2} \cdot X_t+b_f)+sigmoid(W_{i1}\cdot h_{t-1}+W_{i2} \cdot X_t+b_{i1})\cdot tanh(W_{i3}\cdot h_{t-1}+W_{i4} \cdot X_t+b_{i2}))$ （7）

2.长短期记忆网络 LSTM（Long Short-Term Memory）的具体设计细节及应用注意事项

2.1.LSTM的参数

LSTM的参数量计算

如基本结构介绍及一图流中所示，LSTM的模型参数主要包括权重和偏置，其中权重根据计算对象分为输入参数权重（Input Weight）和针对隐藏状态的权重（Recurrent Weight），偏置则对应遗忘门、输入门和输出门中的计算。具体信息可以参照一图流，里面已经标的很清楚了

我们来简单地举个例子，下图是我设计的一个简单的LSTM模型，其中输入为长度807的序列数据，LSTM隐藏单元的大小被设置成4000，因此HiddenState和Cell State的大小都为4000

图二

根据图一的展示，对应图二我们可以简单的计算得到

输入参数的权重参数数量为：

4（输入接口）*[4000*1]（单个权重矩阵大小）*[807*3]（输入数据维度=所有空间长度累乘*通道长度）=16000*2421=38736000;

隐藏状态的权重参数数量为：

4（输入接口）*[4000*4000](单个权重矩阵大小)=16000*4000=64000000;

偏置参数的数量为：

4（四个处理点位）*[4000*1]（单个偏置矩阵大小）=16000*1=16000;

总参数为：16000+64000000+38736000= 102752000

注意：在matlab中，LSTM的设计总是对象单一时间步的，比如这里，当我的输入具有空间属性时，LSTM会先将其flatten至序列再进行处理，所以输入参数的数量就是：所有空间长度累乘*通道长度。并不是你输入一个序列长度为100（S）通道为3的数据，程序就建立输入维度为3的模型，在计算中迭代100次。而是输入维度为300，不限迭代次数模型。

在设计LSTM的输入时，一般不包含空间属性，只有通道属性、批次、和时间步，也就是CBT。在计算参数的例子中，故意加入了空间属性，也是为了说明这一特性。读者可以参考后面的实战示例，加深理解。

图：当LSTM输入图像时

2.2.状态激活方法 StateActivation Fuction 与门激活方法 GateActivation Function

在LSTM（长短期记忆网络）中，状态激活方法（State Activation Function）和门激活方法（Gate Activation Function）起着至关重要的作用。它们分别用于生成候选记忆状态和控制信息的流动。虽然默认情况下，状态激活方法使用tanh函数，门激活方法使用sigmoid函数，但在实际应用中，这些激活函数可以根据具体需求进行调整和替换。

状态激活方法（State Activation Function）

默认激活函数：tanh

tanh函数（双曲正切函数）将输入映射到-1到1之间，其输出均值为0，这有助于保持梯度在传播过程中相对稳定，从而在一定程度上缓解梯度消失问题。

softsign

softsign函数是另一种非线性激活函数，其输出范围也是-1到1，但与tanh函数相比，softsign函数在输入接近0时具有更大的梯度，这有助于模型在训练初期更快地收敛。此外，softsign函数在输入绝对值较大时，其输出趋近于±1的速度更慢，这有助于避免梯度爆炸问题，以改善模型的收敛速度和稳定性。

ReLU（Rectified Linear Unit）

ReLU函数是深度学习中常用的激活函数，其输出在输入为正数时为输入本身，在输入为负数时为0。ReLU函数具有计算简单、非饱和性的特点，能够有效缓解梯度消失问题。尽管ReLU函数在LSTM中不常用作状态激活函数，但在某些特定的应用场景下，如处理稀疏特征或加速训练过程时，可以尝试使用ReLU函数。然而，需要注意的是，由于ReLU函数的输出范围不是-1到1，直接使用ReLU函数可能会破坏LSTM中记忆细胞的0中心化特性，因此在实际应用中需要进行适当的调整。

门激活方法（Gate Activation Function）

默认激活函数：sigmoid

特点：sigmoid函数将输入映射到0到1之间，其输出值可以解释为信息通过门的概率。sigmoid函数在输入较小或较大时，其输出趋近于0或1，这有助于实现门控机制，即控制信息的流动。

hard-sigmoid

hard sigmoid函数作为sigmoid函数的近似，计算速度更快。

2.3.Bi-LSTM

Bi-LSTM特别简单，其结构就是两个LSTM模块镜像并列，然后方向相反，因此其有两条流向相反的细胞状态和两条流向相反的隐藏状态，因此其参数量的计算就是简单的LSTM*2即可（隐藏单元数量相同时），读者可以对应图二和图三进行领悟，这里就没必要展开来讲了。

图三

2.4.LSTM的输入和输出（重点）

LSTM及Bi-LSTM的输入一般为序列数据，格式通常为（CBT），即一个通道度、批次度和时间度。在2.1的示例中，我们可以发现LSTM的单次输入总会铺平成向量。

对于输出长度为K的向量X时，其对应输入的权重转化矩阵W为N*K，N为隐藏单元数量，在计算权重时，W在前，输入向量X在后，进行矩阵乘法，W×X，输出的参数大小正好是N*1，对应隐藏状态的大小，使得数据对齐。

LSTM的输出比较重要，分为序列状态输出（sequence output）和最后状态输出（last output），注意！无论那种输出，它们输出的皆是隐藏状态信息，而不是细胞状态信息。

输出为隐藏状态主要是LSTM的设计考量，一是隐藏状态就是LSTM对下一个时间步信息的预测信息，二是隐藏状态在不同时间步的变化更大，含有更多信息，在序列输出时也能够保留更多的特征细节，读者可以对照式5和式7进行理解。

从应用层面而言：

序列状态输出常用于序列到序列的任务，如机器翻译、信号转化等任务，其中每个时间步的输出都对应着输出序列中的一个元素。但是记住，序列状态输出的时间序列数据，虽然它可以跟神经元连接，但输出是有时间步长的，也就是你的每个训练样本的对应学习标签要相同长度的数据。

最后状态输出用于序列到标签的任务，如文本分类、情感分析等，其中LSTM的输出用于预测整个序列的类别或标签。回归预测同样如此。需要注意的是最后状态输出没有时间步长，或者时间步长为1。

其实无论是序列状态输出还是最后状态输出，在分类和回归问题都可以使用，这主要看最后一层的LSTM是序列输出还是步长输出，比如我堆叠三层LSTM，前两层用序列输出的LSTM，最后一层用最后状态输出的LSTM。LSTM的设计是很自由的，各位只需要知道序列输出的信息更多即可。

图：三层嵌套LSTM，第三层LSTM用了最后状态输出，时间步长消去了。

图：三层嵌套LSTM，但第三层还是用sequence输出，输出始终保持是时间序列状态。

关于MinLength：

SequenceinputLayer常常作物LSTM的前置输入层，可以设置时间步长的最短长度，实际上这个操作就是检验，如果输入的数据时间步长小于MinLength的值，那么程序会将自动进行报错（不会自动填充），而且这并不会影响输出数据时间维度的大小，seqence输出的话，输入是长度K，输出就是K，无关MinLength。

2.5.LSTM特别的偏置初始化：unit-forget-gate

unit-forget-gate作为偏置初始化策略时，这通常意味着对遗忘门的偏置进行特定的初始化，以鼓励或抑制遗忘门在训练初期的行为。例如，在某些情况下，可能会将遗忘门的偏置初始化为一个较大的正值，以使得模型在训练开始时倾向于保留更多的信息（因为遗忘门的输出接近于1），这有助于模型更稳定地学习长期依赖关系。

3.LSTM实战

3.1 LSTM-涡散发动机组寿命预测（序列到序列）：

在该实战示例中，对于LSTM模型而言：单次输入是一系列多维向量（多个步长的不同属性信息）（K*P，K为时间步幅，长度任意，P为属性数量即向量维度=17），输出为多个（K个）1*1的数值，训练最小批次量为20。

也就是序列到序列，序列的长度即时间步长K，即：

输入：17（C）*20（B）*K（T）

输出：1（C）*20（B）*K（T）

因此损失函数也是基于输出（1（C）*20（B）*K（T））与标准预期（1（C）*20（B）*K（T））的MSE均值进行计算。

原代码下载：

25-LSTM-SequenceToSequence_Example.zip
链接: https://pan.baidu.com/s/1c1_GBWSi4d946T1NizQoeA?pwd=as2g提取码: as2g

代码展示及分析：

%使用深度学习进行序列到序列回归
%本示例展示了如何使用深度学习来预测发动机的剩余使用寿命（RUL）。
%为了训练一个深度神经网络，使其能够从时间序列或序列数据中预测数值，您可以使用长短期记忆（LSTM）网络。
%本示例使用了[1]中描述的涡扇发动机退化模拟数据集。示例中训练了一个LSTM网络，以根据表示发动机中各种传感器的时间序列数据，预测发动机的剩余使用寿命（预测性维护），以循环次数为单位。训练数据包含100台发动机的模拟时间序列数据。每个序列的长度各不相同，并且对应于一个从开始到故障（RTF）的完整实例。测试数据包含100个部分序列，以及每个序列结束时对应的剩余使用寿命值。
%数据集包含100个训练观测值和100个测试观测值。

%涡扇发动机退化模拟数据集中的每个时间序列代表一台不同的发动机。每台发动机在开始时具有不同程度的初始磨损和制造差异。每台发动机在时间序列开始时均正常运行，并在序列中的某个时刻出现故障。在训练集中，故障幅度逐渐增大，直至系统失效。
%在训练数据，一个cell是一个样本，cell中的矩阵行对应时间步幅，列对应不同的属性（通道）
%文件中有26列数字，用空格分隔。每一行是在单个运行周期内采集的数据快照，每一列代表一个不同的变量。各列对应属性如下：
%第1列 – 机组编号
%第2列 – 循环时间
%第3–5列 – 运行设置
%第6–26列 –21个不同的传感器测量值

%因此在这个示例中，我们只知道发动机在不同时间步的属性，以及它不同时间步对应的剩余寿命，模型所要做的，就是利用不同时间步的数据，预测不同时间步中对应的剩余寿命。
clear all
load("data.mat");
%移除具有恒定值的特征
%对于所有时间步长都保持恒定的特征可能会对训练产生负面影响。找出具有相同最小值和最大值的行（属性）数据，并将这些属性移除。
XTrainConcatenatedTimesteps = cat(1,XTrain{:});
m = min(XTrainConcatenatedTimesteps,[],1);
M = max(XTrainConcatenatedTimesteps,[],1);
idxConstant = M == m;

for i = 1:numel(XTrain)
    XTrain{i}(:,idxConstant) = [];
end
numFeatures = size(XTrain{1},2)
Normalize Training Predictors
%将训练预测变量标准化为均值为零、方差为一。按样本处理。
XTrainConcatenatedTimesteps = cat(1,XTrain{:});
mu = mean(XTrainConcatenatedTimesteps,1);
sig = std(XTrainConcatenatedTimesteps,0,1);

for i = 1:numel(XTrain)
    XTrain{i} = (XTrain{i} - mu) ./ sig;
end
%Clip Responses
%为了在发动机接近故障时从序列数据中学习更多信息，我们将数据进行裁剪，减去前面保持健康的部分，将剩余寿命（RUC）在阈值150处进行截断。也就是说，预测参数RUC的范围仅为1~150，超过了150的部分，我们令其为150，这样，数据就能在低寿命样本群体中有较好的表现。
thr = 150;
for i = 1:numel(TTrain)
    TTrain{i}(TTrain{i} > thr) = thr;
end
%此图显示了第一个观测样本及其对应的截断响应。
    
%Prepare Data for Padding
%准备数据进行填充
%为了最小化添加到小批量数据中的填充量，按序列长度对训练数据进行排序。然后，选择一个能够均匀划分训练数据的小批量大小，以减少小批量数据中的填充量。
按序列长度对训练数据进行排序。
for i=1:numel(XTrain)
    sequence = XTrain{i};
    sequenceLengths(i) = size(sequence,1);
end

[sequenceLengths,idx] = sort(sequenceLengths,"descend");
XTrain = XTrain(idx);
TTrain = TTrain(idx);
%View the sorted sequence lengths in a bar chart.
figure
bar(sequenceLengths)
xlabel("Sequence")
ylabel("Length")
title("Sorted Data")
%选择一个能均匀分割训练数据的小批量大小，以减少小批量中的填充量。该图展示了在批量大小为20时，对未排序和已排序序列所添加的填充。

%定义网络架构。创建一个LSTM网络，该网络由一个包含200个隐藏单元的LSTM层、一个大小为50的全连接层以及一个丢弃概率为0.5的丢弃层组成。
numResponses = size(TTrain{1},2);%为1
numHiddenUnits = 200;

layers = [ ...
    sequenceInputLayer(numFeatures)
    lstmLayer(numHiddenUnits,OutputMode="sequence")
    fullyConnectedLayer(50)
    dropoutLayer(0.5)
    fullyConnectedLayer(numResponses)];
%指定训练选项。使用“adam”优化器，以小批量大小为20进行60个周期的训练。设置学习率为0.01。为防止梯度爆炸，将梯度阈值设置为1。为保持序列按长度排序，将“Shuffle（洗牌）”选项设置为“never（从不）”。在图表中显示训练进度，并监控均方根误差（RMSE）指标。
maxEpochs = 60;
miniBatchSize = 20;

options = trainingOptions("adam", ...
    MaxEpochs=maxEpochs, ...
    MiniBatchSize=miniBatchSize, ...
    InitialLearnRate=0.01, ...
    GradientThreshold=1, ...
    Shuffle="never", ...
    Metrics="rmse", ...
    Plots="training-progress", ...
    Verbose=0);

Train the Network
net = trainnet(XTrain,TTrain,layers,"mse",options);

%Test the Network
%使用本示例附带的 processTurboFanDataTest 函数准备测试数据。processTurboFanDataTest 函数从 filenamePredictors 和 filenameResponses 中提取数据，并返回元胞数组 XTest 和 TTest，它们分别包含测试预测序列和响应序列。
filenamePredictors = fullfile("test_FD001.txt");
filenameResponses = fullfile("RUL_FD001.txt");
[XTest,TTest] = processTurboFanDataTest(filenamePredictors,filenameResponses);
%利用从训练数据中计算得出的 idxConstant 来移除具有常数值的特征（不随其他变量或时间变化的固定值）。接着，采用与训练数据相同的参数对测试预测数据进行标准化处理。最后，根据训练数据所使用的阈值对测试响应进行裁剪。
for i = 1:numel(XTest)
    XTest{i}(:,idxConstant) = [];
    XTest{i} = (XTest{i} - mu) ./ sig;
    TTest{i}(TTest{i} > thr) = thr;
end
%使用神经网络进行预测。若要对多个观测值进行预测，请使用 minibatchpredict 函数。minibatchpredict 函数在可用的情况下会自动使用 GPU。
%为了防止函数向数据添加填充，请将小批量大小指定为 1。若要以元胞数组的形式返回预测结果，请将 UniformOutput 设置为 false。
YTest = minibatchpredict(net,XTest,MiniBatchSize=1,UniformOutput=false);
%LSTM 网络对部分序列进行预测时，每次处理一个时间步。在每个时间步，网络仅使用当前时间步的值以及根据之前时间步计算得到的网络状态进行预测。网络在每个预测之间更新其状态。minibatchpredict 函数返回这些预测的序列。预测序列的最后一个元素对应于部分序列的预测剩余使用寿命（RUL）。
%或者，也可以使用 predict 函数并更新网络的 State 属性，来逐个时间步进行预测。这在时间步的值以流的形式到达时非常有用。
%在图表中可视化部分预测结果。
idx = randperm(numel(YTest),4);
figure
for i = 1:numel(idx)
    subplot(2,2,i)
    
    plot(TTest{idx(i)},"--")
    hold on
    plot(YTest{idx(i)},".-")
    hold off
    
    ylim([0 thr + 25])
    title("Test Observation " + idx(i))
    xlabel("Time Step")
    ylabel("RUL")
end
legend(["Test Data" "Predicted"],Location="southeast")
%对于给定的部分序列，预测的当前剩余使用寿命（RUL）是预测序列的最后一个元素。计算预测值的均方根误差（RMSE），并在直方图中可视化预测误差。
for i = 1:numel(TTest)
    TTestLast(i) = TTest{i}(end);
    YTestLast(i) = YTest{i}(end);
end
figure
rmse = sqrt(mean((YTestLast - TTestLast).^2))
histogram(YTestLast - TTestLast)
title("RMSE = " + rmse)
ylabel("Frequency")
xlabel("Error")

%References
%Saxena, Abhinav, Kai Goebel, Don Simon, and Neil Eklund. "Damage propagation modeling for aircraft engine run-to-failure simulation." In Prognostics and Health Management, 2008. PHM 2008. International Conference on, pp. 1-9. IEEE, 2008.

在示例中，为提高模型对低寿命的样本的敏感程度，对剩余寿命的大于150的预测标签修改为150。

RNN数据的预处理-填充与截断

注意！在训练模型前，我们需要将不同时间步幅的数据进行排序，以提高训练效率，这是因为在每次最小批次训练前，程序总是自动地将该批次的数据按最长或最短时间步长的数据进行对齐，即将该批次的数据步长统一提高或裁剪，使样本的K对齐。（后面的例子同理）

在训练长短期记忆网络（LSTM）或其他循环神经网络（RNN）时，处理不同长度序列的常见方法是对序列进行填充（padding）或截断（truncation），以确保每个批次中的所有序列具有相同的长度。在Matlab中，对于函数trainnet（），LSTM默认的处理方法是后填充。

填充：

确定批次中所有序列的最大长度。

将所有序列填充到这个最大长度。通常使用零填充（即在序列的末尾添加零），因为零对大多数神经网络来说是一个中性元素。

后填充（post-padding）：在序列的末尾添加零，这是最常见的方法。

前填充（pre-padding）：在序列的开头添加零，某些情况下可能更有用，但较少使用。

截断：

如果某些序列的长度超过了预设的最大长度限制，可以选择截断这些序列。

截断通常是从序列的开头或末尾去掉一些元素，具体取决于任务需求和序列的特性。

训练结果如图所示：

子函数：

function [predictors,responses] = processTurboFanDataTest(filenamePredictors,filenameResponses)
% processTurboFanDataTest 函数从 filenamePredictors 和 filenameResponses 
% 中提取数据，并返回包含测试预测序列和响应序列的元胞数组 predictors 和 responses。
% 在 filenamePredictors 中，时间序列在系统故障前的某段时间结束。
% filenameResponses 中的数据为测试数据提供了真实剩余使用寿命（RUL）值的向量。

predictors = processTurboFanDataTrain(filenamePredictors);

RULTest = dlmread(filenameResponses);

numObservations = numel(RULTest);

responses = cell(numObservations,1);
for i = 1:numObservations
    X = predictors{i};
    sequenceLength = size(X,1);
    
    rul = RULTest(i);
    responses{i} = rul+sequenceLength-1:-1:rul;
end

end


function [predictors,responses] = processTurboFanDataTrain(filenamePredictors)
% processTurboFanDataTrain 函数从 filenamePredictors 中提取数据，
% 并返回分别包含预测序列和响应序列的元胞数组 predictors 和 responses。
% 数据包含以 zip 格式压缩的文本文件，文件中有 26 列数字，各列之间用空格分隔。
% 每一行代表在一个运行周期内采集的数据快照，每一列代表一个不同的变量。
% 各列对应的内容如下：
% 1: 机组编号
% 2: 运行周期时间
% 3–5: 运行设置
% 6–26: 传感器测量值 1–17

dataTrain = dlmread(filenamePredictors);

numObservations = max(dataTrain(:,1));

predictors = cell(numObservations,1);
responses = cell(numObservations,1);
for i = 1:numObservations
    idx = dataTrain(:,1) == i;
    
    predictors{i} = dataTrain(idx,3:end);
    
    timeSteps = dataTrain(idx,2);
    responses{i} = flipud(timeSteps);
end

end

3.2 Bi-LSTM的信号分类预测（序列到标签）：

源码：
链接: https://pan.baidu.com/s/12bm6sMVHEvlD_5OOf73BCg?pwd=uk5m 提取码: uk5m

数据示例：

单批次模型训练的输入和输出：

模型输入：3（C）*64（B）*K（T）

模型输出：1（C）*64（B）*1（T）或 1（C）*64（B）

注意：由于是last输出，输出的时间步长被消掉了

%Sequence Classification Using Deep Learning
%本示例展示了如何使用长短期记忆（LSTM）网络对序列数据进行分类。
%要训练一个用于分类序列数据的深度神经网络，您可以使用 LSTM 神经网络。LSTM 神经网络允许您将序列数据输入到网络中，并基于序列数据的各个时间步进行预测。

%这个示例使用了波形数据集。本示例通过训练一个长短期记忆（LSTM）神经网络来识别给定时间序列数据的波形类型。训练数据包含四种波形的时间序列数据。每个序列具有三个通道，且长度各不相同。

%Load Sequence Data
%从 WaveformData 加载示例数据。序列数据是一个 numObservations×1 的单元数组，其中 numObservations 是序列的数量。每个序列都是一个 numTimeSteps×numChannels 的数值数组，其中 numTimeSteps 是序列的时间步数，numChannels 是序列的通道数。标签数据是一个 numObservations×1 的分类向量。
clear variables
load WaveformData 
%展示数据
numChannels = size(data{1},2);

idx = [3 4 5 12];
figure
tiledlayout(2,2)
for i = 1:4
    nexttile
    stackedplot(data{idx(i)},DisplayLabels="Channel "+string(1:numChannels))
    xlabel("Time Step")
    title("Class: " + string(labels(idx(i))))
end
classNames = categories(labels)
%留出部分数据进行测试。将数据划分为训练集和测试集，其中训练集包含 90% 的数据，测试集包含剩余的 10% 数据。要使用 trainingPartitions 函数来划分数据，该函数作为本示例的辅助文件附加提供。要访问此文件，请将示例作为实时脚本打开。
numObservations = numel(data);
[idxTrain,idxTest] = trainingPartitions(numObservations,[0.9 0.1]);
XTrain = data(idxTrain);
TTrain = labels(idxTrain);

XTest = data(idxTest);
TTest = labels(idxTest);

%Prepare Data for Padding
%数据填充工作
%在训练过程中，软件默认会将训练数据拆分为小批量，并对序列进行填充，使它们具有相同的长度。过多的填充可能会对网络性能产生负面影响。
%为了防止训练过程添加过多的填充，您可以按序列长度对训练数据进行排序，并选择一个适当的小批量大小，以便一个小批量中的序列具有相似的长度。下图显示了数据排序前后对序列进行填充的效果。

%Get the sequence lengths for each observation.
numObservations = numel(XTrain);
for i=1:numObservations
    sequence = XTrain{i};
    sequenceLengths(i) = size(sequence,1);
end
Sort the data by sequence length.
[sequenceLengths,idx] = sort(sequenceLengths);
XTrain = XTrain(idx);
TTrain = TTrain(idx);
%View the sorted sequence lengths in a bar chart.
figure
bar(sequenceLengths)
xlabel("Sequence")
ylabel("Length")
title("Sorted Data")

%Define LSTM Neural Network Architecture
%定义LSTM神经网络架构如下：
%指定输入大小：
%输入大小设定为输入数据的通道数。
%双向LSTM层：
%使用一个包含120个隐藏单元的双向LSTM层。双向LSTM层能够在每个时间步长上从序列的前向和后向两个方向进行学习。
%设置输出为序列的最后一个元素。
%全连接层：
%添加一个全连接层，其输出大小与类别数相匹配。
%softmax层：
%在全连接层之后，添加一个softmax层，以将输出转换为概率分布。
%如果你在预测时能够访问完整序列（例如，在处理静态数据集或进行批处理预测时），则可以在网络中使用双向LSTM层。双向LSTM层能够在每个时间步长上利用整个序列的信息。如果你无法在预测时访问完整序列（例如，在进行时间序列预测或逐步预测时），则应使用单向LSTM层代替。
%简而言之，这段描述定义了一个LSTM神经网络，该网络包括一个与输入数据通道数相匹配的输入层、一个双向LSTM层（或根据需要替换为单向LSTM层）、一个全连接层和一个softmax层。
numHiddenUnits = 120;
numClasses = 4;

layers = [
    sequenceInputLayer(numChannels)
    bilstmLayer(numHiddenUnits,OutputMode="last")
    fullyConnectedLayer(numClasses)
    softmaxLayer]

%Specify Training Options
%指定训练选项。选择哪个选项需要进行实证分析。为了通过实验探索不同的训练选项配置，可以使用实验管理器应用程序。
%使用Adam求解器进行训练。
%训练200个周期。
%指定学习率为0.002。
%将梯度裁剪阈值设置为1。
%为了保持序列按长度排序，禁用洗牌功能。
%在图中显示训练进度，并监控准确率。
%禁用详细输出。
options = trainingOptions("adam", ...
    MaxEpochs=200, ...
    InitialLearnRate=0.002,...
    GradientThreshold=1, ...
    Shuffle="never", ...%注意这里不要乱序训练
    Plots="training-progress", ...
    minibatchsize=64,...
    Metrics="accuracy", ...
    Verbose=false,...
    ValidationData={XTest,TTest},...
    ValidationFrequency=50);

%Train LSTM Neural Network

net = trainnet(XTrain,TTrain,layers,"crossentropy",options);

%Test LSTM Neural Network
%对测试数据进行分类，并计算预测的分类准确率。
%LSTM神经网络net是使用相似长度的序列组成的小批量数据进行训练的。请确保测试数据以相同的方式组织。按序列长度对测试数据进行排序。
numObservationsTest = numel(XTest);
for i=1:numObservationsTest
    sequence = XTest{i};
    sequenceLengthsTest(i) = size(sequence,1);
end

[sequenceLengthsTest,idx] = sort(sequenceLengthsTest);
XTest = XTest(idx);
TTest = TTest(idx);
%对测试数据进行分类，并计算预测的分类准确率。
%使用minibatchpredict函数进行预测，并使用scores2label函数将得分转换为标签。
scores = minibatchpredict(net,XTest);
YTest = scores2label(scores,classNames);
Calculate the classification accuracy. The accuracy is the percentage of correctly correctly predicted labels.
acc = mean(YTest == TTest)
%Display the classification results in a confusion chart.
figure
confusionchart(TTest,YTest)

子函数

function varargout = trainingPartitions(numObservations,splits)
%TRAININGPARTITIONS 用于拆分训练数据的随机索引
% [idx1,...,idxN] = trainingPartitions(numObservations,splits) 返回
% 随机索引向量，以帮助拆分具有指定观察数量的数据集，其中 SPLITS 是一个
% 长度为 N 的分区大小向量，且其元素之和为 1。
%
% % 示例：获取500个观察值的50%-50%训练-测试拆分的索引。
% [idxTrain,idxTest] = trainingPartitions(500,[0.5 0.5])
%
% % 示例：获取500个观察值的80%-10%-10%训练、验证、测试拆分的索引。
% [idxTrain,idxValidation,idxTest] = trainingPartitions(500,[0.8 0.1 0.1])

arguments
    numObservations (1,1) {mustBePositive}
    splits {mustBeVector,mustBeInRange(splits,0,1,"exclusive"),mustSumToOne}
end

numPartitions = numel(splits);
varargout = cell(1,numPartitions);

idx = randperm(numObservations);

idxEnd = 0;

for i = 1:numPartitions-1
    idxStart = idxEnd + 1;
    idxEnd = idxStart + floor(splits(i)*numObservations) - 1;

    varargout{i} = idx(idxStart:idxEnd);
end

% Last partition.
varargout{end} = idx(idxEnd+1:end);

end

function mustSumToOne(v)
% Validate that value sums to one.

if sum(v,"all") ~= 1
    error("Value must sum to one.")
end

end

训练结果：