时序模型：线性条件随机场模型 ( Linear - CRF )

1. 线性条件随机场（Linear-CRF）

线性条件随机场（linear chain conditional random field，Linear-CRF）是对隐马尔科夫模型（hidden Markov model，HMM）在更普遍的实际问题上的推广¹。

隐马尔科夫（HMM）模型假设序列数据具有齐次马尔可夫性和观测独立性，这是对实际问题的两种近似简化。线性条件随机场（Linear-CRF）模型取消了观测独立性假设，并削弱了齐次马尔可夫性假设；Linear-CRF考虑一个模型每时刻隐藏状态，都是受其相邻时刻隐藏状态和全部时刻观测值共同影响的情况假设。即Linear-CRF模型中：
$$P(y_i|X,y_1,\cdots ,y_{i-1},y_{i+1},\cdots ,y_n)=P(y_i|X,y_{i-1},y_{i+1})\text{\hspace{1em}(1.1)}$$

线性条件随机场（Linear-CRF）是一种概率无向图模型，在实际中一般假设输入序列X和标签序列Y具有相同的序列长度，其模型结构如下所示：

图1 X和Y具有相同图结构的Linear-CRF模型
图中每个节点都表示一个随机变量，每条边表示所连接的随机变量间的概率依赖关系。且无向边表示相连变量间有概率依赖关系，但并不一定是因果关系。

此无向图中，最大团（maximal clique）由每时刻隐藏状态、其前一时刻隐藏状态和同时刻观测值组成。所以根据 Hammersley-Clifford 定理，Linear-CRF模型的条件概率分布可表示为：

$$P(y|X)=\frac{exp[{\sum }_{i\in N}{\sum }_{l\in K}{\mu }_l^{(i)}{s}_l(y_i,x_i)+{\sum }_{i\in N}{\sum }_{k\in K^2}{\lambda }_k{t}_k(y_{i-1},y_i,x_i)]}{{\sum }_{y\in Y}exp[{\sum }_{i\in N}{\sum }_{l\in K}{\mu }_l^{(i)}{s}_l(y_i,x_i)+{\sum }_{i\in N}{\sum }_{k\in K^2}{\lambda }_k{t}_k(y_{i-1},y_i,x_i)]}\text{\hspace{1em}(1.2)}$$
式中，分子是状态序列 $y$ 出现的可能性分数（下文简称：状态序列分数），分母是状态序列出现可能性分数的归一化项（下文简称：归一化项），二者相除等于状态序列出现的条件概率。$N$ 是序列长度，$K$ 是所有可能的状态数量。

$s_l(y_i,x_i)$ 是状态特征，是定义在节点上的特征函数，它依赖于当前时刻隐藏状态和观测值，表示是否由当前时刻观察值发射（emit）到第 $l$ 种隐藏状态（是：状态特征为1，否：状态特征为0）；$t_k(y_{i-1},y_i,x_i)$ 是转移特征，是定义在边上的特征函数，它依赖于前一时刻及当前时刻的隐藏状态，表示当前时刻是否由前一时刻隐藏状态转到第 $k$ 种隐藏状态（是：转移特征为1，否：转移特征为0）。

${\mu }_l^{(i)}$ 是发射分数（emission score），它表示由当前时刻观测值发射（emit）到第 $l$ 种隐藏状态的可能性权重；${\lambda }_k$ 是转移分数（transfer score），它表示由上一时刻隐藏状态转移到第 $k$ 种隐藏状态的可能性权重；发射分数与转移分数之和，决定隐藏状态 $y_i$ 出现的条件概率。

2. 线性条件随机场（Linear-CRF）的矩阵形式

2.1 模型定义

为便于并行计算，根据上述定义可将 Linear-CRF 改写为矩阵形式：

$$P(y|X)=\frac{exp[E(X)\cdot M_E(y)+T\cdot M_T(y)]}{{\sum }_{y\in Y}exp[E(X)\cdot M_E(y)+T\cdot M_T(y)]}\text{\hspace{1em}(2.1.1)}$$

式中，$X$ 是输入特征序列，$y$ 是标签序列，$Y$ 是所有可能的标签序列的集合；$E(X)$ 是发射分数矩阵，$M_E(y)$ 是发射分数掩码矩阵；$T$ 是转移分数矩阵，$M_T(y)$ 是转移分数掩码矩阵。

发射分数矩阵（emission matrix）$E(X)$ 由输入特征的线性变换得到，其形状为$[N,K]$，它表示输入序列每时刻观测值发射到每种标签状态的可能性权重，表达式如下所示：

$$E(X)=X{W}_E+b_E\text{\hspace{1em}(2.1.2)}$$

式中 $W_E$ 和 $b_E$ 是可学习的模型参数。转移分数矩阵（transfer matrix）$T$ 由可学习的模型参数得到，其形状为$[K,K]$，它表示每种标签类别转移到任意另一种类别的可能性权重。

模型训练时：发射分数掩码矩阵（emission mask matrix）$M_E(y)$ 由标签序列的 one-hot 编码得到，其形状为 $[N,K]$，它表示每时刻生效的发射分数类别。转移分数掩码矩阵（transfer mask matrix）$M_T(y)$ 由标签序列相邻时刻状态的二维 one-hot 编码得到，其形状为 $[N,K,K]$，它表示每时刻生效的转移分数类别。

模型预测时：根据发射分数矩阵和转移分数矩阵，直接由维特比算法解码（decoder）出预测标签序列，不需要计算发射分数掩码矩阵和转移分数掩码矩阵。

2.2 模型参数优化

Linear-CRF是一种判别式概率模型，其参数优化方法是针对目标函数的极大似然估计。判别模型的目的是最大化对所有样本标签预测正确的概率，所以Linear-CRF模型的目标函数是 $max{\prod }_{i=1}^{m}P(y_i|X_i)$，其中 $m$ 是小批量样本的批量大小。由此我们可简单地设定：Linear-CRF模型的损失函数等于模型预测样本标签出现的条件概率的负数。

那么根据极大似然估计法，模型参数的估计为：

$$\begin{array}{c}\hat{\theta }=argmin-{\prod }_{i=1}^{m}P(y_i|X_i)\\ \\ \hat{\theta }=argmin-log({\prod }_{i=1}^{m}P(y_i|X_i))\\ \\ \hat{\theta }=argmin-{\sum }_{i=1}^{m}log(P(y_i|X_i))\\ \\ \hat{\theta }=argmin-{\sum }_{i=1}^{m}log(\frac{exp[E(X_i)\cdot M_E(y_i)+T\cdot M_T(y_i)]}{{\sum }_{y_i\in Y}exp[E(X_i)\cdot M_E(y_i)+T\cdot M_T(y_i)]})\end{array}\text{\hspace{1em}(2.2.1)}$$

设：一种隐藏状态序列出现的可能性权重分数（Score）为 $S_i^{(y_i)}=E(X_i)\cdot M_E(y_i)+T\cdot M_T(y_i)$。将其带入上式后，Linear-CRF模型参数估计的表达式可进一步简化为：
$$\hat{\theta }=argmin\sum _{i=1}^m[log(\sum _{y_i\in Y}{e}^{S_i^{(y_i)}})-S_i^{(y_i)}]\text{\hspace{1em}(2.2.2)}$$

进而模型的损失函数为：
$$Loss=\frac{1}{m}\sum _{i=1}^m[log(\sum _{y_i\in Y}{e}^{S_i^{(y_i)}})-S_i^{(y_i)}]\text{\hspace{1em}(2.2.3)}$$

式中，归一化项 $log({\sum }_{y_i\in Y}{e}^{S_i^{(y_i)}})$ 不能直接穷举计算，因为 $S_i^{(y_i)}$ 共有 $K^N$ 种可能的情况，N为序列长度。本文设计了一种基于动态规划的高效计算方法，借助矩阵运算，可将归一化项计算的时间复杂度从 $O(K^n)$ 降低到 $O(n)$。

本方法利用指数函数性质，将累加所有可能的状态序列分数（即：${\sum }_{y_i\in Y}{e}^{S_i^{(y_i)}}$）转化为一个动态规划问题。

首先，我们设：

• 第 $i$ 个样本序列 $t$ 时刻，发射到状态 $k$ 的发射分数记作 $emi{t}_{(i,t,k)}$
• 第 $i$ 个样本序列 $t$ 时刻，从状态 $j$ 转移到状态 $k$ 的转移路径分数记作 $R_{(i,t,j,k)}$
• 第 $i$ 个样本序列 $t$ 时刻，通向状态 $k$ 的所有状态路径的分数的指数之和记作 $S_{(i,t,k)}$
• 任意样本序列的任意时刻，从状态 $l$ 转移到状态 $k$ 的转移分数记作 $tran{s}_{(l,k)}$

下图中的每一个指向箭头（黑色或灰色），都是一条可能的状态路径（route）。

图2 Linear-CRF的状态路径
所以动态规划算法的初始状态和状态转移方程为：
$$\{\begin{array}{llc}初始状态：& S_1=E_{t=1}=\left[\begin{array}{c}\left[\begin{array}{ccc}emi{t}_{(1,1,1)}& \cdots & emi{t}_{(1,1,K)}\end{array}\right]\\ ⋮\\ \left[\begin{array}{ccc}emi{t}_{(M,1,1)}& \cdots & emi{t}_{(M,1,K)}\end{array}\right]\end{array}\right]=\left[\begin{array}{c}\left[\begin{array}{ccc}{S}_{(1,1,1)}& \cdots & S_{(1,1,K)}\end{array}\right]\\ ⋮\\ \left[\begin{array}{ccc}{S}_{(M,1,1)}& \cdots & S_{(M,1,K)}\end{array}\right]\end{array}\right],& S_1\in R^{M\times 1\times K}\\ & & \\ & & \\ 状态路径集合：& R_t=T+E_t=\left[\begin{array}{c}\left[\begin{array}{ccc}tran{s}_{(1,1)}+emi{t}_{(1,t,1)}& \cdots & tran{s}_{(1,K)}+emi{t}_{(1,t,K)}\\ ⋮& \ddots & ⋮\\ tran{s}_{(K,1)}+emi{t}_{(1,t,1)}& \cdots & tran{s}_{(K,K)}+emi{t}_{(1,t,K)}\end{array}\right]\\ ⋮\\ \left[\begin{array}{ccc}tran{s}_{(1,1)}+emi{t}_{(M,t,1)}& \cdots & tran{s}_{(1,K)}+emi{t}_{(M,t,K)}\\ ⋮& \ddots & ⋮\\ tran{s}_{(K,1)}+emi{t}_{(M,t,1)}& \cdots & tran{s}_{(K,K)}+emi{t}_{(M,t,K)}\end{array}\right]\end{array}\right]=\left[\begin{array}{c}\left[\begin{array}{ccc}{R}_{(1,t,1,1)}& \cdots & R_{(1,t,1,K)}\\ ⋮& \ddots & ⋮\\ R_{(1,t,K,1)}& \cdots & R_{(1,t,K,K)}\end{array}\right]\\ ⋮\\ \left[\begin{array}{ccc}{R}_{(M,t,1,1)}& \cdots & R_{(M,t,1,K)}\\ ⋮& \ddots & ⋮\\ R_{(M,t,K,1)}& \cdots & R_{(M,t,K,K)}\end{array}\right]\end{array}\right],& R_t\in R^{M\times K\times K},t=2,\cdots ,N\\ & & \\ & & \\ 状态转移方程：& S_t=log[{\sum }_{k\in K}{e}^{(S_{t-1}+R_t^T)}]=log(\sum \left[\begin{array}{c}\left[\begin{array}{ccc}{e}^{(S_{(1,t-1,1)}+R_{(1,t,1,1)})}& \cdots & e^{(S_{(1,t-1,K)}+R_{(1,t,K,1)})}\\ ⋮& \ddots & ⋮\\ e^{(S_{(1,t-1,1)}+R_{(1,t,1,K)})}& \cdots & e^{(S_{(1,t-1,K)}+R_{(1,t,K,K)})}\end{array}\right]\\ ⋮\\ \left[\begin{array}{ccc}{e}^{(S_{(M,t-1,1)}+R_{(M,t,1,1)})}& \cdots & e^{(S_{(M,t-1,K)}+R_{(M,t,K,1)})}\\ ⋮& \ddots & ⋮\\ e^{(S_{(M,t-1,1)}+R_{(M,t,1,K)})}& \cdots & e^{(S_{(M,t-1,K)}+R_{(M,t,K,K)})}\end{array}\right]\end{array}\right]),& S_t\in R^{M\times 1\times K},t=2,\cdots ,N\end{array}$$
归一化项等于：
$$log(\sum _{y_i\in Y}{e}^{S_i^{(y_i)}})=log(\sum _{k\in K}{S}_{(i,N,k)})\text{\hspace{1em}(2.2.5)}$$

综上所述，损失函数表达式为：
$$Loss=\frac{1}{m}\sum _{i=1}^m[log(\sum _{k\in K}{S}_{(i,N,k)})-S_i^{(y_i)}]\text{\hspace{1em}(2.2.3)}$$

2.3 模型预测

每时刻发射分数为节点权重，以转移分数为节点连接边权重，使用维特比算法找出分数最高的状态路径（route）。如下图所示：

图3 Linear-CRF解码过程
设：第 $i$ 个样本序列截至 $t$ 时刻的最佳状态路径（出现可能性最大的状态路径）分数记作 $s_{(i,t)}$。则，Linear-CRF解码过程的初始状态和状态转移方程为：
$$\{\begin{array}{llc}初始状态：& {\hat{y}}_1^{(i)}=\underset{k\in K}{\mathrm{argmax}}emi{t}_{(i,1,k)},s_{(i,1)}=\underset{k\in K}{\max }emi{t}_{(i,1,k)}& \\ & & \\ 状态转移方程：& {\hat{y}}_t^{(i)}=\underset{k\in K}{\mathrm{argmax}}{s}_{(i,t-1)}+tran{s}_{({\hat{y}}_{t-1}^{(i)},k)}+emi{t}_{(i,t-1,k)},& t=2,\dots ,N\end{array}\text{\hspace{1em}(2.2.6)}$$

3. 线性条件随机场（Linear-CRF）的工程实现

Linear-CRF模型的结构如下图所示：

图4 Linear-CRF模型结构图

3.1 TensorFlow框架的代码实现

"""
V2.0 此版本在CRF的发射分数计算和损失函数计算中加入了mask机制，可有效去除padding值引起的噪声干扰。 2022.4.10
V1.0 Baseline版本 2022.3.28
"""
# %%
import numpy as np
from tqdm import tqdm
import tensorflow as tf
from tensorflow.keras import layers, losses
from tensorflow.keras.optimizers import SGD, Adam
import tensorflow_addons as tfa

gpus = tf.config.experimental.list_physical_devices(device_type="GPU")
tf.config.experimental.set_visible_devices(devices=gpus[0], device_type="GPU")


class LinerCRFlayer(layers.Layer):
    def __init__(self, num_labels, **kwargs):
        super().__init__(**kwargs)
        self.num_labels = num_labels
        self.Emit = layers.Dense(num_labels)

    def build(self, batch_input_shape):
        self.transform = self.add_weight(
            name='transform',
            shape=[self.num_labels, self.num_labels],
            initializer='glorot_normal'
        )
        super().build(batch_input_shape)

    def call(self, inputs, mask, training):
        """
        :param:
        inputs is the feature sequence, its' shape is [batch_size, time_steps, hidden_num].
        mask is the mask matrix of inputs, its' shape is [batch_size, time_steps].
        emit is the emit score matrix, its' shape is [batch_size, time_steps, num_labels].
        transform is the transform score matrix, its' shape is [num_labels, num_labels].
        transform_i is the transform score vector of label state at time i-1, its' shape is [batch_size, 1, num_labels].

        :return:
        if training is True, return the emit score matrix of inputs and transfer score matrix;
        if training is False, return the label sequence decoded by Viterbi algorithm.
        """
        # 发射矩阵的padding噪声不用在输入端滤除，
        # 而是在训练阶段的损失函数计算和预测阶段的维特比结果选择中滤除。
        emit = self.Emit(inputs)
        if training:
            return [emit, self.transform]
        else:
            """维特比寻找最佳标注路径"""
            y = []
            # 初始状态赋值
            emit_i = emit[:, :1, :]
            y.append(tf.argmax(emit_i, axis=-1))
            Score = tf.reduce_max(emit_i, axis=-1)
            # 状态转移方程
            for i in range(1, inputs.shape[1], 1):
                y_onehot = tf.one_hot(y[i - 1], self.num_labels)
                transform_i = tf.matmul(y_onehot, self.transform)
                Score = tf.expand_dims(Score, axis=-1) + transform_i + emit[:, i:i + 1, :]
                # 记录本轮状态
                y.append(tf.argmax(Score, axis=-1))
                Score = tf.reduce_max(Score, axis=-1)
            return tf.concat(y, axis=-1)


#
class LinerCRF(tf.keras.Model):
    def __init__(self, input_block, hidden_block, num_labels):
        super().__init__()
        self.num_labels = num_labels
        # model define
        self.Input = input_block
        self.Hiddens = hidden_block
        self.Output = LinerCRFlayer(self.num_labels)
        # compile
        self.optimizer = None
        self.loss_fun = None
        self.metric = None

    def call(self, X, mask, num_mask, training):
        """
        mask 是布尔型掩码矩阵，元素值为 True 或 False
        num_mask 是浮点型掩码矩阵，元素值为 1. 或 0.
        """
        hiddens = [self.Input(X, num_mask)['last_hidden_state']]
        for Hidden in self.Hiddens:
            hiddens = Hidden(hiddens[0], mask=mask, training=training)
        outputs = self.Output(hiddens[0], num_mask, training=training)
        return outputs

    def compile(self, optimizer=Adam(), loss=None, metrics=None, **kwargs):
        if loss is None:
            raise ValueError("Please give the loss function!")
        self.optimizer = optimizer
        self.loss_fun = loss

    def train(self, train_generator, epochs=10, verbose=100,  gard_theta=1, pad_id=0):
        """
        tensorflow原生LSTM/GRU层的mask接口要求是布尔矩阵（True/False），而其他层的mask接口要求是01矩阵（0./1.）;
        所以有 num_mask 转换。
        """
        loss_pre_epoch = []
        for epoch in range(epochs):
            loss_pre_batch = []
            for i, train_data in enumerate(tqdm(list(train_generator), desc='step')):
                X, y, mask, num_mask = train_data
                with tf.GradientTape(persistent=True) as tape:
                    emit, transform = self(X, mask, num_mask, training=True)
                    loss = self.loss_fun(emit, transform, y, num_mask)
                # 参数更新
                train_varb = self.trainable_variables
                del train_varb[197:199]    # 删除不使用的BERT模型pool层参数
                grads = tape.gradient(loss, train_varb)
                grads = self.grad_clipping(grads, gard_theta)  # 梯度裁剪
                self.optimizer.apply_gradients(zip(grads, train_varb))
                # record
                loss_pre_batch.append(loss.numpy())
                if i % verbose == 0:
                    print('step ', i, ' ', 'loss:', np.mean(loss_pre_batch))
            loss_pre_epoch.append(np.mean(loss_pre_batch))
            print('epoch ', epoch, ' ', 'loss:', np.mean(loss_pre_batch), '\n')
        return loss_pre_epoch

    #
    def predict(self, X):
        y_hat = self(X, training=False)
        return y_hat

    # 梯度裁剪
    def grad_clipping(self, grads, theta):
        """
        梯度裁剪
        """
        theta = tf.constant(theta, dtype=tf.float32)
        new_grads = [tf.convert_to_tensor(grad) if isinstance(grad, tf.IndexedSlices) else grad for grad in grads]
        norm = tf.math.sqrt(sum(tf.reduce_sum(grad ** 2).numpy() for grad in new_grads))
        norm = tf.cast(norm, tf.float32)
        if tf.greater(norm, theta):
            new_grads = [(grad * theta) / norm for grad in new_grads]
        return new_grads

    def CRFLossFun(self, emit, transform, y, num_mask):
        """
        **Param
        :emit:          emit is input sequences, and its' shape is [batch_size, time_steps, num_labels]
        :transform:     transform is transform score, and its' shape is [num_labels, num_labels]
        :y:             y is label sequences, and its' shape is [batch_size, time_steps]
        :y_len:         y_len is one-dim tensor and its' shape is [batch_size]
        **备注
        输入序列padding噪声包含在emit矩阵的对映时间步位置中；
        标签序列padding噪声包含在num_mask矩阵( MEy & MTy )的对映时间步位置中；
        经过矩阵数乘，两类padding噪声积累在 Ey 和 Ty 矩阵的对映时间步位置中，通过 ME 和 MT 掩码矩阵的过滤，可全完去除对映时间步位置中的padding噪声（两类噪声都被除去）。
        y 矩阵 time_steps 的尾端，包含标签序列的 padding 噪声。
        **Output
        :return:
        """
        ## 计算标签序列的分数
        # 计算标签序列的发射分数
        MEy = tf.one_hot(y, self.num_labels)   # 引入标签序列padding噪声，包含在num_mask矩阵(MEy)对映时间步位置的元素中
        Ey = tf.reduce_sum(emit * MEy, axis=[-1])
        # 计算标签序列的转移分数
        yts1 = tf.expand_dims(tf.one_hot(y[:, :y.shape[-1] - 1], self.num_labels), axis=-2)
        yt = tf.expand_dims(tf.one_hot(y[:, 1:], self.num_labels), axis=-2)
        MTy = tf.matmul(tf.transpose(yts1, [0, 1, 3, 2]), yt)
        Ty = tf.reduce_sum(MTy * transform, axis=[-1, -2])  # 引入标签序列padding噪声，包含在num_mask矩阵(MTy)对映时间步位置的元素中
        # 汇总标签序列分数，并去除标签序列中的padding噪声;
        # ME = self.generate_num_mask(Ey, y_len);
        # MT = self.generate_num_mask(Ty, y_len - 1).
        yScore = tf.reduce_sum(Ey * num_mask, axis=-1) + tf.reduce_sum(Ty * num_mask[:, 1:], axis=-1)

        ## 输入序列所有可能的状态序列的分数之和
        allScore = tfa.text.crf_log_norm(emit, tf.reduce_sum(num_mask, axis=-1), transform)

        ## 计算损失函数值
        return tf.reduce_mean(allScore - yScore)

---

1. 详细演进推导过程可参考文章Linear-Chain-CRF 的前世今生 ↩︎