随着现代工业的发展，工业过程日益繁复，设备集成度越来越高，随之而来的故障种类也愈发多样^[1-2]. 得益于数据采集、传输、存储技术的发展，可以收集大量工业过程的测量数据，因此基于数据驱动的故障检测和诊断成为研究热点^[3-4]. 多元统计分析方法作为典型的数据驱动方法，在化工、石油、炼钢等行业的故障检测和诊断中得到了广泛的应用.

基于多元统计分析的故障诊断方法包括2类，即贡献图法和故障重构法，其中贡献图法通过计算过程变量对统计量的贡献程度，进而定位故障变量^[5-6]. Westerhuis等^[7]研究贡献图的涂抹效应，并引入贡献的控制限；Li等^[8]提出基于全潜结构投影的贡献图方法，用统一的方式描述所有变量对监测统计量的贡献；Alcala等^[9-10]相继提出基于重构贡献和相对贡献的故障诊断方法. 相较于贡献图法，故障重构是一种监督的故障诊断方法，通过对历史故障数据的分析，提取故障子空间，其包含了故障的特征^[5]. Dunia等^[11]利用故障方向信息通过重构进行故障识别；Li等^[12]提出与输出相关的故障子空间提取方法和约简故障子空间提取方法，有效降低了故障子空间的维数；冯晓伟等^[13]采用广义主成分分析方法剔除正常建模数据和故障建模数据之间的公共信息，构造更精确的故障子空间.

由于物料、能量及信息间的传递和反馈控制作用的存在，过程变量高度耦合^[14]. 故障子空间之间存在相似性，一个故障子空间可以重构多个类型的故障，造成故障子空间与故障类型之间存在一对多的问题. 在现有研究中，故障子空间提取高度依赖故障样本，往往假设各故障模式下都具有足够的理想标签样本；然而在实际生产过程中，故障样本往往难以获取，由于不同类型的故障发生概率不同，不同模式的故障样本之间存在严重不平衡^[15-17]. 此外，大多现有方法致力于提高故障子空间的高效性，而忽略了故障子空间之间的可分辨性. 这使得故障子空间的故障表征能力和分辨能力严重下降. 因此开展基于数据增强的故障分类研究成为工程实践中的迫切需要.

目前，基于数据增强的虚拟样本生成方法已经被广泛应用到图像分类领域，如生成对抗网络(generative adversarial nets, GAN)^[18-19]、自编码器(variational autoencoder, VAE)^[20-21]，但这些方法难以有效捕获数据在时间上的相关关系及其动态特性. Arjovsky等^[22]使用Wasserstein距离代替JS散度来衡量真实样本和生成样本之间的分布距离，提出Wasserstein距离生成对抗网络(Wasserstein GAN, WGAN). 随后，Gulrajani等^[23]增加了梯度惩罚项，提出基于梯度惩罚的WGAN(gradient penalty WGAN, WGAN-GP). 此外，合成少数过采样技术(synthetic minority oversampling technique, SMOTE)作为一种重采样方法，利用随机线性插值在一个少数类样本和它最近的邻居样本之间产生新的样本^[24]. 这2种数据增强方法被广泛应用到故障诊断领域以解决故障样本不平衡问题. 时间序列生成对抗网络(time-series generative adversarial networks, TimeGAN)是专门用于时间序列数据增强的生成模型，能够有效捕获时间序列数据中的时间相关性^[25-27]. 因此，本研究提出基于TimeGAN数据增强的非线性复杂过程故障分类方法，主要创新包括以下4点：1)采用TimeGAN网络对历史故障数据进行数据增强，克服故障样本少、不平衡的问题；2)采用马氏距离对生成虚拟样本进行质量评估，并剔除不可信样本，避免其对模型训练产生负面影响；3)在核空间而非原始数据空间提取故障子空间，提高子空间故障表征能力和分辨能力；4)设计基于故障重构的故障分类策略以实现在线故障样本分类，形成完整的基于数据增强的故障分类方案.

时间序列数据的时间特征对传统的数据增强方法提出了严峻的挑战. 时间序列数据的数据增强不仅要处理一个时间点内的特征，而且要捕捉这些特征在时间范围内的复杂属性.

TimeGAN是由Yoon等^[25]提出的时间序列生成框架. 它有效结合了监督训练的控制性和无监督训练的灵活性. TimeGAN模型结构如图1所示.

自编码器网络降低了输入数据的维度，利用嵌入模块将原始特征空间映射到包含静态特征和时态特征的潜空间，映射过程如下：

(1)$ {{\boldsymbol{h}}_{{\mathrm{s}}}} = {e_{\text{S}}}({\boldsymbol{s}}),\; {{\boldsymbol{h}}_{t}} = {e_{\text{X}}}({{\boldsymbol{h}}_{{\mathrm{s}}}},{{\boldsymbol{h}}_{{t - }{\text{1}}}},{{\boldsymbol{x}}_{t}}). $

式中：${\boldsymbol{s}}$和${\boldsymbol{x}}$分别表示静态和时态特征数据，${{e}_{\text{S}}}$和${{e}_{\text{X}}}$分别表示静态和时态特征数据到潜空间的映射函数，${{\boldsymbol{h}}_{{\mathrm{s}}}}$和$ {{\boldsymbol{h}}_{t}} $则对应${\boldsymbol{s}}$和${\boldsymbol{x}}$在潜空间的低维映射. 嵌入模块主要是对输入数据进行降维，而复现模块则是通过训练逆映射函数完成数据从低维到高维的恢复，逆映射过程如下：

(2)$ {\boldsymbol{\tilde s}} = {{r}_{\text{S}}}({{\boldsymbol{h}}_{{\mathrm{s}}}}),\; {{{\boldsymbol{\tilde x}}}_{t}} = {r_{\text{X}}}({{\boldsymbol{h}}_{t}}{\text{)}}. $

式中：${\boldsymbol{\tilde s}}$和${{\boldsymbol{\tilde x}}_{t}}$表示经过逆映射后恢复的静态和时态特征数据. 自编码器实现了特征空间与潜空间的可逆映射.

为了提高自编码器中嵌入模块和复现模块对输入数据分布的学习能力，使用数据重构损失函数${L_{\text{R}}}$对自编码器进行优化.

(3)$ {L_{\text{R}}} = {E_{{\boldsymbol{s}},\;{{\boldsymbol{x}}_{1{\text{ }}:\;{{ T}}}} \sim {{p}}}}\left[ {{{\left\| {{\boldsymbol{s}} - {\boldsymbol{\tilde s}}} \right\|}_2}+\sum\limits_t {{{\left\| {{{\boldsymbol{x}}_{t}} - {{{\boldsymbol{\tilde x}}}_{t}}} \right\|}_2}} } \right]. $

表中：T为时间序列的长度，p为数据的分布.

生成对抗网络是为了实现生成器和判别器之间的平衡. 生成器捕获实时序列的潜在分布并生成新的合成序列. 判别器在真实序列和合成序列之间提供正确的分类. 其生成函数可以定义为

(4)$ {{{\boldsymbol{\hat h}}}_{{\mathrm{s}}}} = {g_{\text{S}}}({{\boldsymbol{z}}_{{\mathrm{s}}}}),\; {{{\boldsymbol{\hat h}}}_{t}} = {g_{\text{X}}}({{{\boldsymbol{\hat h}}}_{{\mathrm{s}}}},{{{\boldsymbol{\hat h}}}_{{t - }{\text{1}}}},{{\boldsymbol{z}}_{t}}). $

式中：$ {{\boldsymbol{z}}_{{\mathrm{s}}}} $和$ {{\boldsymbol{z}}_{t}} $表示随机静态特征和时态特征，$ {{\boldsymbol{\hat h}}_{{\mathrm{s}}}} $和$ {{\boldsymbol{\hat h}}_{t}} $则表示潜空间中合成静态特征和时态特征. 为了提高生成器对真实数据时序相关性学习能力和潜空间表征能力，引入有监督误差损失函数：

(5)$ {L_{\text{S}}} = {E_{{\boldsymbol{s}}{\text{,}}{{\boldsymbol{x}}_{1{\text{ }}:\;{{ T}}}} \sim {{p}}}}\left[ {\sum\limits_t {{{\left\| {{{\boldsymbol{h}}_{t}} - {g_{\text{X}}}({{\boldsymbol{h}}_{{\mathrm{s}}}},{{\boldsymbol{h}}_{{t - }{\text{1}}}},{{\boldsymbol{z}}_{t}})} \right\|}_2}} } \right]. $

判别器模块采用双向递归神经网络(recurrent neural network, RNN)区分真实样本和合成样本，判别模型如下：

(6)$\tilde{y}_{\mathrm{s}}=d_{\mathrm{S}}\left(\tilde{\boldsymbol{h}}_{\mathrm{s}}\right) ,\;\tilde{y}_t=d_{\mathrm{X}}\left({\boldsymbol{u}}_t^{\mathrm{b}}, {\boldsymbol{u}}_t^{\mathrm{f}}\right) . $

式中：${{d}_{\text{S}}}$和${{d}_{\text{X}}}$分别表示静态和时态特征判别函数，$\tilde{\boldsymbol{h}}_{\mathrm{s}} $表示${\boldsymbol{h}}_{\mathrm{s}} $或$\hat{\boldsymbol{h}}_{\mathrm{s}} $，$ {\tilde y_{\mathrm{s}}} $和$ {\tilde y_t} $表示合成序列静态和时态特征的判别器标签，$ {\boldsymbol{u}}_t^{\mathrm{b}} $和$ {\boldsymbol{u}}_t^{\mathrm{f}} $分别表示反向和前向RNN的隐藏状态. 定义无监督的对抗损失函数${L_{\text{U}}}$反映生成器和判别器之间的竞争关系：

(7)$ \begin{split} {L_{\text{U}}} =& {E_{{\boldsymbol{s}},{{\boldsymbol{x}}_{1:\;{{ T}}}} \sim {\tilde{{p}}}}}\left[ {\ln\; {y_{\mathrm{s}}}+\sum\limits_t {\ln\; {y_t}} } \right]+ \\& {E_{{\boldsymbol{s}},{{\boldsymbol{x}}_{1:\;{{ T}}}} \sim {\tilde{{p}}}}}\left[ {\ln\; (1 - {{\tilde y}_{\mathrm{s}}})+\sum\limits_t {\ln\; (1 - {{\tilde y}_t})} } \right]. \end{split} $

式中：$ {y_{\mathrm{s}}} $和$ {y_t} $表示真实序列静态和时态特征的标签.

尽管TimeGAN的目标是生成与真实数据分布相似的样本，但生成的虚拟样本可能仍存在与真实数据分布不完全一致的问题，或者生成的虚拟样本存在噪声或不可信，这些样本会对模型训练产生负面影响进而影响故障子空间的故障表征能力. 因此，采用合理的方法对生成的虚拟样本进行过滤筛选是必要的，可以保证虚拟样本质量. 在评估虚拟样本质量时，常见的距离准则包括马氏距离、欧式距离及三角隶属度函数法^[28]等. 马氏距离考虑了特征之间的相关性和协方差矩阵，因此适用于高维数据和具有相关性的数据^[16]. 欧氏距离简单、直观、易于理解和实现，对于线性关系较强的数据有较好效果，但不考虑特征之间的相关性和尺度差异，可能受到特征尺度影响. 三角隶属度函数法具有较高的灵活性，可以根据实际情况调整隶属度函数的形状和参数，适应不同的数据分布和异常情况，但隶属度函数的选择和参数的确定需要一定的经验，可能存在主观性和不确定性，且计算复杂度较高. 在具体应用中，须根据实际数据的特点选择合适的距离准则. 考虑到工业过程数据往往呈现非线性、多维度、高度相关等特点，选择马氏距离来评估虚拟样本质量.

假设系统共有$N$类故障，其中$M$类为小子样故障，第$\lambda $类小子样故障矩阵$ {{\boldsymbol{X}}^\lambda } \in {{{\bf{R}}}^{{\text{ }}{\alpha _\lambda } \times m}} $，其中${\alpha _\lambda }$和$m$分别表示真实故障样本数和变量数，则样本${\boldsymbol{x}}_i^\lambda $与${{\boldsymbol{X}}^\lambda }$之间的距离可以表示为

(8)$ {\text{Dis}}({\boldsymbol{x}}_i^\lambda ,{{\boldsymbol{X}}^\lambda }) = \left[{({\boldsymbol{x}}_i^\lambda - {{{\boldsymbol{\bar x}}}^\lambda }){{({{\boldsymbol{S}}^\lambda })}^{ - 1}}{{({\boldsymbol{x}}_i^\lambda - {{{\boldsymbol{\bar x}}}^\lambda })}^{\text{T}}}}\right]^{1/2} . $

式中：${{\boldsymbol{\bar x}}^\lambda }$为样本均值；${{\boldsymbol{S}}^\lambda }$为协方差矩阵. ${\boldsymbol{S}}^\lambda $必须可逆，但在真实数据中，各变量之间高度耦合，会导致协方差矩阵奇异. 为了克服这一问题，采用主成分分析(principal component analysis, PCA)进行去相关，并在主元空间中计算样本之间的距离. 首先对${{\boldsymbol{S}}^\lambda }$进行特征值分解：

(9)$ {\sigma _j}{{\boldsymbol{\mu }}_j} = {{\boldsymbol{S}}^\lambda }{{\boldsymbol{\mu }}_j}. $

式中：${{\boldsymbol{\sigma}} _j} $、${{\boldsymbol{\mu }}_j}{\text{ }}(j = 1,2, \cdots ,m)$分别为第$j$个特征值和特征向量. 将所有特征值降序排列，选取前${A_\lambda }$个特征值对应的特征向量组成投影矩阵${{\boldsymbol{U}}^\lambda }, {\boldsymbol{U}}^\lambda = \left[ {{\boldsymbol{\mu }}_1^{},{\boldsymbol{\mu }}_2^{}, \cdots ,{\boldsymbol{\mu }}_{{A_\lambda }}^{}} \right] \in {{{\bf{R}}}^{{\text{ }}m \times {A_\lambda }}}$，其中主元个数${A_\lambda }$由累计方差贡献率确定. 将${{\boldsymbol{X}}^\lambda }$投影到主元空间，其得分可以表示为

(10)$ {{\boldsymbol{T}}^\lambda } = [{\boldsymbol{t}}_1^\lambda ,{\boldsymbol{t}}_2^\lambda , \cdots ,{\boldsymbol{t}}_{{\alpha _\lambda }}^\lambda ] . $

式中：$ {\boldsymbol{t}}_i^\lambda = {\boldsymbol{x}}_i^\lambda {{\boldsymbol{U}}^\lambda } , {{\boldsymbol{T}}^\lambda } \in {{\bf{R}}^{{\text{ }}{\alpha _\lambda } \times {A_\lambda }}}$. 计算样本得分均值$ {{\boldsymbol{\bar T}}^\lambda } = \sum\nolimits_{i = 1}^{{\alpha _\lambda }} {{\boldsymbol{t}}_i^\lambda } $，则样本${\boldsymbol{x}}_i^\lambda $与${{\boldsymbol{X}}^\lambda }$之间的距离可由式(8)转换为

(11)$ {\text{Dis}}({\boldsymbol{x}}_i^\lambda ,{{\boldsymbol{X}}^\lambda }) = \left[{({\boldsymbol{t}}_i^\lambda - {{{\boldsymbol{\bar T}}}^\lambda }){{({\boldsymbol{S}}_{\text{T}}^\lambda )}^{ - 1}}{{({\boldsymbol{t}}_i^\lambda - {{{\boldsymbol{\bar T}}}^\lambda })}^{\text{T}}}}\right]^{1/2} . $

假设第$\lambda $类小子样故障经过1.1节TimeGAN网络生成${\beta _\lambda }$个虚拟故障样本$ {\boldsymbol{x}}_i^{*\lambda }{\text{ (}}{i}{\text{ = 1,2,}} \cdots {\text{,}}\;{\beta _\lambda }{\text{)}} $，对其进行基于马氏距离的虚拟数据筛选的具体步骤总结如算法1所示.

算法1　基于马氏距离的虚拟数据筛选

输入：${M}$类小子样故障真实数据${{\boldsymbol{X}}^\lambda }$ $(\lambda = 1,2, \cdots ,{M})$; $({N} - {M})$类大子样故障真实数据；生成的虚拟故障样本集${{\boldsymbol{X}}^{*\lambda }} = ({\boldsymbol{x}}_1^{*\lambda },{\boldsymbol{x}}_2^{*\lambda }, \cdots ,{\boldsymbol{x}}_{{\beta _\lambda }}^{*\lambda })$$(\lambda = 1,2, \cdots ,{M})$.

1 for $\lambda = 1,2, \cdots ,{M}$ do

2 　定义第$\lambda $类经过滤后的故障数据集为空集，即${{\boldsymbol{X}}^{**\lambda }} = \emptyset $;

3 　对${{\boldsymbol{X}}^\lambda }$进行中心化，并计算协方差矩阵${{\boldsymbol{S}}^\lambda }$；

4 　按照式(9)对${{\boldsymbol{S}}^\lambda }$进行特征值分解，求得投影矩阵${{\boldsymbol{U}}^\lambda }$；

5 　按照式（11）分别计算每个样本${\boldsymbol{x}}_{i}^\lambda $与${{\boldsymbol{X}}^\lambda }$之间的距离$ {\text{Dis}}({\boldsymbol{x}}_{i}^\lambda ,{{\boldsymbol{X}}^\lambda }){\text{ }}({i} = 1,2, \cdots ,{\alpha _\lambda }) $；

6 　求各样本与${{\boldsymbol{X}}^\lambda }$之间距离的最大值，$\max\; ({\text{Dis}} ({\boldsymbol{x}}_{i}^\lambda ,{{\boldsymbol{X}}^\lambda }))$；

7　 对虚拟故障数据集${{\boldsymbol{X}}^{*\lambda }}$进行中心化，即${\boldsymbol{x}}_{i}^{*\lambda } = {\boldsymbol{x}}_{i}^{*\lambda } - {{\boldsymbol{\bar x}}^\lambda }$；

8 　for ${i} = 1,2, \cdots ,{\beta _\lambda }$ do

9 　　按照式(11)求得样本${\boldsymbol{x}}_{i}^{*\lambda }$与${{\boldsymbol{X}}^\lambda }$之间的距离$ {\text{Dis}}({\boldsymbol{x}}_{i}^{*\lambda },{{\boldsymbol{X}}^\lambda }) $；

10 　　if $ {\text{Dis}}({\boldsymbol{x}}_{i}^{*\lambda },{{\boldsymbol{X}}^\lambda }) < \max\;({\text{Dis}}({\boldsymbol{x}}_{i}^\lambda ,{{\boldsymbol{X}}^\lambda })) $

11　　　将样本${\boldsymbol{x}}_{i}^{*\lambda }$添加到${{\boldsymbol{X}}^{**\lambda }}$中，即${{\boldsymbol{X}}^{**\lambda }} = [{{\boldsymbol{X}}^{**\lambda }};{\boldsymbol{x}}_{i}^{*\lambda }]$；

12 　　else

13 　　　舍弃样本${\boldsymbol{x}}_{i}^{*\lambda }$；

14 　　endif

15 　endfor

16 　求得经数据增强后的${N}$类故障的最小故障样本数$\gamma $；

17 　从${{\boldsymbol{X}}^{**\lambda }}$中随机选择$\gamma $个样本作为最终故障数据矩阵，记作${{\boldsymbol{X}}^{**\lambda }}$.

18 endfor

输出：${N}$类经过滤后的故障数据集${{\boldsymbol{X}}^{**\lambda }}(\lambda = 1,2, \cdots ,N)$.

核独立成分分析(kernel independent component analysis, KICA)^[29]作为独立成分分析的非线性扩展，已经成功地应用于非线性非高斯过程的故障检测中，其故障检测原理如下. 设有非线性系统包含$L$个样本，$m$个变量，数据标准化后记作$ {\boldsymbol{X}} = {[{\boldsymbol{x}}_1^{},{\boldsymbol{x}}_2^{}, \cdots ,{\boldsymbol{x}}_{L}^{}]^{\text{T}}} \in {{\bf{R}}^{{ L} \times {m}}} $. KICA的核心思想是将非线性数据${\boldsymbol{X}}$经过转换函数$\phi$映射到高维线性空间，即${\boldsymbol{\varPhi }} = {[\phi ({\boldsymbol{x}}_1^{}),\phi ({\boldsymbol{x}}_2^{}), \cdots ,\phi ({\boldsymbol{x}}_{L}^{})]^{\text{T}}}$. 由于$\phi$难以显式表示，作为一种替代方法，核函数法常常被用来计算Gram矩阵：

(12)$ {\boldsymbol{K}} = {\boldsymbol{\varPhi }}{{\boldsymbol{\varPhi }}^{\text{T}}} \in {{\bf{R}}^{{\text{ }}{L} \times {L}}}. $

通常采用径向基作为核函数，即$ {{{K}}_{{i},{j}}} = \exp \;\left( - {\left\| {{\boldsymbol{x}}_{i}^{} - {\boldsymbol{x}}_{j}^{}} \right\|^2}/{c}\right) $，$c$为核参数. 对核矩阵${\boldsymbol{K}}$进行中心化和缩放：

(13)$\left. \begin{array}{l} {{\boldsymbol{K}}^*} = ({{\boldsymbol{I}}_{L}} - \dfrac{{\text{1}}}{{L}}{{\boldsymbol{1}}_{L}}{\boldsymbol{1}}_{L}^{\text{T}}){\boldsymbol{K}}{({{\boldsymbol{I}}_{L}} - \dfrac{{\text{1}}}{{L}}{{\boldsymbol{1}}_{L}}{\boldsymbol{1}}_{L}^{\text{T}})^{\text{T}}}, \\ {\boldsymbol{\bar K}} = \dfrac{{{{\boldsymbol{K}}^*}}}{{{{\mathrm{trace}}}\;({{\boldsymbol{K}}^*})/{L}}}. \end{array} \right\}$

式中：${{\boldsymbol{I}}_{L}}$为$L$维单位矩阵，${{\boldsymbol{1}}_{L}}$表示元素全为1的$L$维列向量. 对${\boldsymbol{\bar K}}$进行特征值分解$ {{\eta }_{i}}{{\boldsymbol{v}}_i} = {\boldsymbol{\bar K}}{{\boldsymbol{v}}_i} $，取前$p$个最大特征值对应的特征向量构建白化矩阵：

(14)$ {\boldsymbol{P}} = \sqrt {L} {{\boldsymbol{\bar \varPhi }}^{\text{T}}}{\boldsymbol{V}}{{\boldsymbol{\varLambda }}^{ - 1}}. $

式中：${\boldsymbol{\varLambda }} = {\text{diag}}\;({\eta _{\text{1}}},{\eta _{\text{2}}}, \cdots ,{\eta _{p}})$,$ {\boldsymbol{V}} = [{{\boldsymbol{v}}_{\text{1}}},{{\boldsymbol{v}}_{\text{2}}}, \cdots ,{{\boldsymbol{v}}_{p}}] $；${\boldsymbol{\bar \varPhi }} $为${\boldsymbol{\varPhi }}$标准化后的矩阵，$ {\boldsymbol{\bar \varPhi }} = {[\bar \phi ({\boldsymbol{x}}_1^{}),\bar \phi ({\boldsymbol{x}}_2^{}), \cdots ,\bar \phi ({\boldsymbol{x}}_{L}^{})]^{\text{T}}} $. 则${\boldsymbol{\bar \varPhi }}$经过白化后可得正交矩阵：

(15)$ {\boldsymbol{Z}} = {\boldsymbol{\bar \varPhi P}} = \sqrt {L} {\boldsymbol{\bar \varPhi }}{{\boldsymbol{\bar \varPhi }}^{\text{T}}}{\boldsymbol{V}}{{\boldsymbol{\varLambda }}^{ - 1}} = \sqrt {L} {\boldsymbol{\bar KV}}{{\boldsymbol{\varLambda }}^{ - 1}}. $

假设存在单位正交矩阵${\boldsymbol{C}} \in {{\bf{R}}^{{\text{ }}{p} \times {d}}}$，使得

(16)$ {\boldsymbol{S}} = {\boldsymbol{Z}}{\boldsymbol{C}}{{\boldsymbol{D}}^{1/2}} = \sqrt {L} \bar {\boldsymbol{K}}{\boldsymbol{V}}{{\boldsymbol{\varLambda }}^{ - 1}}{\boldsymbol{C}}{{\boldsymbol{D}}^{1/2}} = \bar {\boldsymbol{K}}{\boldsymbol{W}}. $

式中：$ {\boldsymbol{W}} = \sqrt {L} {\boldsymbol{V}}{{\boldsymbol{\varLambda }}^{ - 1}}{\boldsymbol{C}}{D^{1/2}} $，$ {\boldsymbol{D}} = {\text{diag}}\;\left( {{\eta _1},{\eta _2}, \cdots ,{\eta _{d}}}\right) $，$d$为独立成分的个数. 通过最大化${\boldsymbol{S}}$的非高斯性可求得矩阵${\boldsymbol{C}}$，具体的推导和求解过程可参考文献[29].

一般情况下，KICA模型采用${{{I}}^2}$和${\text{SPE}}$统计量作为监测指标. 对于给定样本${{\boldsymbol{x}}_{{\text{new}}}}$及其核向量$ {{\boldsymbol{\bar k}}_{{\text{new}}}} $，则对应的统计量可以计算为

(17)$ {{{I}}^2}({{\boldsymbol{x}}_{{\text{new}}}}) = {{\boldsymbol{\bar k}}_{{\text{new}}}}{\boldsymbol{W}}{{\boldsymbol{W}}^{\text{T}}}{\boldsymbol{\bar k}}_{{\text{new}}}^{\text{T}}, $

(18)$ {\text{SPE}}({{\boldsymbol{x}}_{{\text{new}}}}){\text{ = }}{{\boldsymbol{\bar k}}_{{\text{new}}}}L{\boldsymbol{V}}{{\boldsymbol{\varLambda }}^{ - 1}}({\boldsymbol{I}} - {\boldsymbol{C}}{{\boldsymbol{C}}^{\text{T}}}){{\boldsymbol{\varLambda }}^{ - 1}}{{\boldsymbol{V}}^{\text{T}}}{\boldsymbol{\bar k}}_{{\text{new}}}^{\text{T}}. $

为了方便后续表述，将式(17)、(18)统一表示为

(19)$ {\text{Index}}({{\boldsymbol{x}}_{{\text{new}}}}) = {{\boldsymbol{\bar k}}_{{\text{new}}}}{\boldsymbol{G\bar k}}_{{\text{new}}}^{\text{T}}, $

(20)$ \boldsymbol{G}= \left\{\begin{array}{ll}\boldsymbol{W} \boldsymbol{W}^{\mathrm{T}}, & {I}^2 \text { 统计量; } \\ L \boldsymbol{V} \boldsymbol{{\boldsymbol{\varLambda}}}^{-1}\left(\boldsymbol{I}-\boldsymbol{C}^{\mathrm{T}}\right) \boldsymbol{{\boldsymbol{\varLambda}}}^{-1} \boldsymbol{V}^{\mathrm{T}}, & \text { SPE统计量. }\end{array}\right. $

由于复杂系统的变量分布形式未知，难以确保过程数据遵循特定的分布假设. 为了解决这个问题，核密度估计(kernel density estimation, KDE)通常被用来确定控制限，KDE是一种非参数经验密度估计技术，其不依赖于任何特定的分布假设，适合于估计过程数据的概率密度函数^[30-31]. 因此，本研究采用KDE获取控制限${ I}_{{\text{limit}}}^2$和$ {\text{SP}}{{\text{E}}_{{\text{limit}}}} $. 对于在线监测样本${{\boldsymbol{x}}_{{\text{new}}}}$，按照如下诊断逻辑判断系统运行状态：如果$ {{ I}^2}({{\boldsymbol{x}}_{{\text{new}}}}) > { I}_{{\text{limit}}}^2 $或$ {\text{SPE}}({{\boldsymbol{x}}_{{\text{new}}}}) > $$ {\text{SP}}{{\text{E}}_{{\text{limit}}}} $，则系统发生了故障；如果$ {{{I}}^2}({{\boldsymbol{x}}_{{\text{new}}}}) \leqslant {{I}}_{{\text{limit}}}^2 $且$ {\text{SPE}}({{\boldsymbol{x}}_{{\text{new}}}}) \leqslant $$ {\text{SP}}{{\text{E}}_{{\text{limit}}}} $，则系统正常.

假设主子空间中故障方向为${{\boldsymbol{F}}_{\text{I}}}$，故障样本${{\boldsymbol{x}}_t}$的核向量在主子空间的分量可以表示为

(21)$ {\boldsymbol{k}}_t^{{\mathrm{f}}}{\boldsymbol{WA}} = {\boldsymbol{k}}_t^{{\mathrm{n}}}{\boldsymbol{WA}}+{\boldsymbol{f}}_t^{\text{I}}{{\boldsymbol{F}}_{\text{I}}}. $

式中：$ {\boldsymbol{k}}_t^{\text{f}} $为故障样本的核向量，$ {\boldsymbol{k}}_t^{\text{n}} $为修正后的核向量，A为混合矩阵，${\boldsymbol{A}} = {L^{ - 1/2}}{{\boldsymbol{D}}^{ - 1/2}}{{\boldsymbol{C}}^{\text{T}}}{\boldsymbol{\varLambda }}{{\boldsymbol{V}}^{\text{T}}}$；由于$ {\boldsymbol{k}}_t^{\mathrm{n}} $是零均值的，$ {\boldsymbol{k}}_t^{\text{f}}{\boldsymbol{WA}} $可近似表示为$ {\boldsymbol{k}}_t^{\text{f}}{\boldsymbol{WA}} \approx {\boldsymbol{f}}_t^{\text{I}}{{\boldsymbol{F}}_{\text{I}}} $；${\boldsymbol{f}}_t^{\text{I}} $为故障方向，$ {\boldsymbol{f}}_{t}^{\text{I}} \in {{\bf{R}}^{{\text{ }}1 \times {L}}} $. 假设有$q$个故障样本，则其故障数据的核矩阵在主子空间的分量可以表示为

(22)$ {{\boldsymbol{K}}^{\text{f}}}{\boldsymbol{WA}} = {[{\boldsymbol{f}}_1^{\text{I}},{\boldsymbol{f}}_2^{\text{I}}, \cdots ,{\boldsymbol{f}}_q^{\text{I}}]^{\text{T}}}{{\boldsymbol{F}}_{\text{I}}}. $

对其进行奇异值分解：

(23)$ {({{\boldsymbol{K}}^{\text{f}}}{\boldsymbol{WA}})^{\text{T}}} = {{\boldsymbol{U}}_{\text{I}}}{{\boldsymbol{D}}_{\text{I}}}{\boldsymbol{V}}_{\text{I}}^{\text{T}}. $

式中：$ {{\boldsymbol{U}}_{\text{I}}}$、${{\boldsymbol{V}}_{\text{I}}} $分别为左、右奇异矩阵，${{\boldsymbol{D}}_{\text{I}}} $为奇异值矩阵.

将${{\boldsymbol{D}}_{\text{I}}}$中的奇异值按降序排列，则${{\boldsymbol{F}}_{\text{I}}} = {{\boldsymbol{U}}_{\text{I}}}$.

同理，假设残差子空间中故障方向为${{\boldsymbol{F}}_{{\text{SPE}}}}$,故障数据的核矩阵在残差子空间的分量可以表示为

(24)$ {{\boldsymbol{K}}^{\text{f}}}({\boldsymbol{I}} - {\boldsymbol{WA}}) = {[{\boldsymbol{f}}_1^{{\text{SPE}}},{\boldsymbol{f}}_2^{{\text{SPE}}}, \cdots ,{\boldsymbol{f}}_q^{{\text{SPE}}}]^{\text{T}}}{{\boldsymbol{F}}_{{\text{SPE}}}}. $

式中：${\boldsymbol{f}}_i^{{\text{SPE}}} $为第i个样本在残差子空间的故障幅度.

对式（24）进行奇异值分解：

(25)$ {({{\boldsymbol{K}}^{\text{f}}}({\boldsymbol{I}} - {\boldsymbol{WA}}))^{\text{T}}} = {{\boldsymbol{U}}_{{\text{SPE}}}}{{\boldsymbol{D}}_{{\text{SPE}}}}{\boldsymbol{V}}_{{\text{SPE}}}^{\text{T}}. $

将${{\boldsymbol{D}}_{{\text{SPE}}}}$中的奇异值按降序排列，则${{\boldsymbol{F}}_{{\text{SPE}}}} = {{\boldsymbol{U}}_{{\text{SPE}}}}$.

至此可得故障在主子空间和残差子空间的故障方向（即故障子空间）${{\boldsymbol{F}}_{\text{I}}}$和${{\boldsymbol{F}}_{{\text{SPE}}}}$，在实际应用中，故障子空间的维度为使得重构后的故障样本的统计量低于控制限的最小值.

对于已知故障方向${\boldsymbol{F}}$的故障样本，假设沿着故障方向${\boldsymbol{F}}$对其核向量进行修正：

(26)$ {\boldsymbol{k}}_t^{\text{n}} = {\boldsymbol{k}}_t^{\text{f}} - {{\boldsymbol{f}}_t}{\boldsymbol{F}}. $

式中：$ {{\boldsymbol{f}}_t} $为待估计的故障方向. 由式(19)，修正后的统计量可以表示为

(27)$ {\text{Index}}({\boldsymbol{k}}_t^{\text{n}}) = ({\boldsymbol{k}}_t^{\text{f}} - {{\boldsymbol{f}}_t}{\boldsymbol{F}}){\boldsymbol{G}}{({\boldsymbol{k}}_t^{\text{f}} - {{\boldsymbol{f}}_t}{\boldsymbol{F}})^{\text{T}}}. $

为了使修正后的样本回到正常状态，则将估计故障方向转化为如下优化问题：

(28)$ \mathop {{\min}}\limits_{{{\boldsymbol{f}}_t}}\; \{ {\text{Index}}({\boldsymbol{k}}_t^{\text{n}})\} . $

这是一个无约束最小二乘问题，其解析解为

(29)$ {{\boldsymbol{f}}_t} = {\boldsymbol{k}}_t^{\text{f}}{{\boldsymbol{G}}^{\text{T}}}{\boldsymbol{F}}{({{\boldsymbol{F}}^{\text{T}}}{{\boldsymbol{G}}^{\text{T}}}{\boldsymbol{F}})^\dagger }. $

式中：$\dagger $表示伪逆运算符.

上述为针对某一特定故障方向对故障样本进行重构，进而求得对应故障幅度的过程. 在实际应用中，分别对每类故障提取故障子空间${\boldsymbol{F}}_{\text{I}}^i$和${\boldsymbol{F}}_{{\text{SPE}}}^i$，并将这些故障子空间存到故障库中. 当在线监测到某一样本${{\boldsymbol{x}}_{{\text{new}}}}$的${{{I}}^2}$统计量超出控制限${{I}}_{{\text{limit}}}^2$时，分别用${\boldsymbol{F}}_{\text{I}}^i$$(i = 1,2, \cdots ,N)$故障子空间对其进行重构，按照式(29)，其对应的故障幅度和修正后的核向量为

(30)$ \boldsymbol{f}_{\text {new }}^i=\boldsymbol{k}_{\text {new }} \boldsymbol{G}^{\mathrm{T}} \boldsymbol{F}_{\mathrm{I}}^i\left(\boldsymbol{F}_{\mathrm{I}}^{i \mathrm{~T}} \boldsymbol{G}^{\mathrm{T}} \boldsymbol{F}_{\mathrm{I}}^i\right)^\dagger , $

(31)$ {\boldsymbol{k}}_{{\text{new}}}^i = {{\boldsymbol{k}}_{{\text{new}}}} - {\boldsymbol{f}}_{{\text{new}}}^i{\boldsymbol{F}}_{\text{I}}^i. $

将其代入式(17)得到修正后的${{{I}}^2}$统计量为

(32)$ {{I}}_{{\text{new}}}^{2\_i} = {\boldsymbol{k}}_{{\text{new}}}^i{\boldsymbol{G}}{({\boldsymbol{k}}_{{\text{new}}}^i)^{\text{T}}}. $

如此可得故障样本${{\boldsymbol{x}}_{{\text{new}}}}$在$N$个${\boldsymbol{F}}_{\rm{I}}^i$故障子空间下重构后的统计量$ {{I}}_{{\text{new}}}^{2\_1} $，$ {{I}}_{{\text{new}}}^{2\_2}, \cdots ,{{I}}_{{\text{new}}}^{2\_N} $. 类似地，如果样本${{\boldsymbol{x}}_{{\text{new}}}}$的${\text{SPE}}$统计量超出控制限${\text{SP}}{{\text{E}}_{{\text{limit}}}}$时，可以求得其在$N$个${\boldsymbol{F}}_{{\text{SPE}}}^i$故障子空间下重构后的统计量$ {\text{SPE}}_{{\text{new}}}^1,{\text{SPE}}_{{\text{new}}}^2, \cdots ,{\text{SPE}}_{{\text{new}}}^N $. 通过对比这些统计量，可以将其最小值所对应的故障子空间编号判断为样本${{\boldsymbol{x}}_{{\text{new}}}}$的故障类型. 故障子空间提取及在线故障分类具体步骤总结如算法2所示.

算法2　故障子空间提取及在线故障分类

输入：正常工况数据${\boldsymbol{X}};\;{N}$类经过滤后的故障数据集${{\boldsymbol{X}}^{**\lambda }}(\lambda = 1,2, \cdots ,{N})$；在线样本${{\boldsymbol{x}}_{{\text{new}}}}$；重构参数$\eta $.

1) 故障子空间提取

(a) 对${\boldsymbol{X}}$进行标准化，并建立KICA模型，求得混合矩阵${\boldsymbol{A}}$,系数矩阵${\boldsymbol{\varSigma }}$及混解矩阵${\boldsymbol{W}}$.

(b) for $\lambda = 1,2, \cdots ,{N}$ do

(c) 　对${{\boldsymbol{X}}^{**\lambda }}$进行非线性映射，求得核矩阵${\boldsymbol{K}}_\lambda ^{\text{f}}$；

(d) 　按照式(23)对核矩阵的主子空间分量进行奇异值分解，求得${{\boldsymbol{U}}_{\text{I}}}$；

(e) 　for ${i} = 1,2, \cdots ,{L}$ do

(f) 　　　令${\boldsymbol{F }}_{\text{I}}^\lambda = {{\boldsymbol{U}}_{\text{I}}}(:,1:{i})$；

(g) 　　　按照式(30)~(32)对${{\boldsymbol{X}}^{**\lambda }}$中所有样本进行故障重构，并计算重构恢复率$ \varsigma ={{{\boldsymbol{X}}}^{**\lambda }中重构后{{I}}^{2}统计量恢复到{{I}}_{\text{limit}}^{2}以下的样本数}/ \quad\;\; {{{\boldsymbol{X}}}^{**\lambda }中所有故障样本数} $

(h) 　　　if $\varsigma \geqslant \eta $

(i) 　　　　break；

(j) 　　　 endif

(k) 　endfor

(l) 　同理，类似于步骤(d)~(k)，求得第$\lambda $类故障残差分量的故障子空间${\boldsymbol{F}}_{{\text{SPE}}}^{\lambda }$；

(m) endfor

2) 在线故障分类

(a) 在线样本${{\boldsymbol{x}}_{{\text{new}}}}$，对其进行标准化并求得核向量$ {{\boldsymbol{\bar k}}_{{\text{new}}}} $，按照式(17)、(18)，求得其$ {{{I}}^2}({{\boldsymbol{x}}_{{\text{new}}}}) $及$ {\text{SPE}}({{\boldsymbol{x}}_{{\text{new}}}}) $统计量；

(b) if $ {{{I}}^2}({{\boldsymbol{x}}_{{\text{new}}}}) > {{I}}_{{\text{limit}}}^2 $

(c) 　分别以$F_{{I}}^1,F_{{I}}^2, \cdots ,F_{{I}}^{N}$为故障方向进行故障重构并计算重构后的$ {{I}}_{{\text{new}}}^{2\_\lambda } $；

(d) 　对比统计量$ {{I}}_{\text{new}}^{2\_\lambda }(\lambda =\text{1,2},\mathrm{ \cdots },{N}) $，将其最小值所对应的故障子空间编号判断为样本${{\boldsymbol{x}}_{{\text{new}}}}$的故障类型；

(e) else if $ {{{I}}^2}({{\boldsymbol{x}}_{{\text{new}}}}) \leqslant {{I}}_{{\text{limit}}}^2 $且$ {\text{SPE}}({{\boldsymbol{x}}_{{\text{new}}}}) > {\text{SP}}{{\text{E}}_{{\text{limit}}}} $

(f) 　　类似于步骤(c)~(g)，在残差空间判断故障类型；

(g) else

(h) 　　${{\boldsymbol{x}}_{{\text{new}}}}$为正常样本；

(i) endif

输出：故障库${\boldsymbol{F}}_{\text{I}}^{} = [{\boldsymbol{F}}_{\text{I}}^1,{\boldsymbol{F}}_{\text{I}}^2, \cdots ,{\boldsymbol{F}}_{\text{I}}^{N}]$和$\;{\boldsymbol{F}}_{{\text{SPE}}}^{} = [{\boldsymbol{F}}_{{\text{SPE}}}^1, {\boldsymbol{F}}_{{\text{SPE}}}^2, \cdots ,{\boldsymbol{F}}_{{\text{SPE}}}^{N}];$${{\boldsymbol{x}}_{{\text{new}}}}$的运行状态及其故障类型.

综上所述，本研究所提基于数据增强的复杂过程故障分类方法主要包括KICA模型建立、故障样本数据增强、故障子空间提取及在线故障分类几个部分，算法整体流程如图2所示.

为了验证所提故障分类方法的性能，采用田纳西伊士曼过程 (Tennessee Eastman process, TEP)^[32]进行实验验证. TEP仿真平台包括正常工况和21种故障工况，正常工况共包括1460个观测样本，随机使用500个样本进行建模，剩余样本用于验证模型并计算控制限，每种故障工况下有2个故障数据文件，例如对于故障1，文件d01中包含480个故障样本，而d01_te中包含960个样本，其中前160个为正常样本，后800个为故障样本.

选取故障1、 2、 6、 7、 12、 14作为已知故障，52个变量作为观测建模变量，具体实验数据设置如表1所示. 其中，故障1作为大子样故障，不对其进行数据增强，其他故障均视作小子样故障进行数据增强. 表中，N₁、N₂分别为历史故障数据中故障样本数量和增强后故障样本数量. 以此验证所提算法的数据增强性能，并讨论多模式故障分类问题.

在通常情况下，一个性能良好的故障分类算法应该具备以下几点：1)故障子空间能够精确表征故障属性；2)模型对正样本具有较高的识别能力，对负样本具有较强的区分能力. 具体表现为：如果在线样本属于某一特定故障类型时，只有用该故障类型对应的故障子空间进行故障重构，才可以将统计量恢复到控制限以下；相反，如果使用其他故障类型对应的故障子空间进行故障重构，则无法将统计量恢复到控制限以下. 因此，分别采用精确率P、召回率R和F1得分3个指标来定量评价故障分类性能，其中精确率用于衡量模型对负样本的区分能力，精确率越高，则区分能力越高；召回率用于衡量模型对正样本的识别能力，召回率越高，识别能力越强；F1得分则为精确率和召回率的调和平均值，体现模型的稳健性. 此外，本研究定义重构恢复率(reconstruction recovery rate, R3)指标以衡量故障子空间表征故障属性的能力. 各性能指标具体定义如下：

(33)$ {{P}} = \frac{{{\text{TP}}}}{{{\text{TP}}+{\text{FP}}}}, $

(34)$ {{R}} = \frac{{{\text{TP}}}}{{{\text{TP + FN}}}}, $

(35)$ {\text{F}}1 = \frac{{2 {{PR}} }}{{{{P + R}}}}, $

(36)$ {\text{R3 = }}\frac{{{\text{NU}}{{\text{M}}_{{\text{FN}}}}}}{{{\text{NU}}{{\text{M}}_{\text{F}}}}}. $

式中：TP、FN分别表示正样本被诊断为正样本、负样本的样本数，FP、TN分别表示负样本被诊断为正样本、负样本的样本数，$ {\text{NU}}{{\text{M}}_{{\text{FN}}}} $表示经重构后统计量恢复到控制限以下的故障样本数，$ {\text{NU}}{{\text{M}}_{\text{F}}} $表示故障样本总数.

分别将小子样故障2、 6、 7、 12、 14的历史故障数据作为TimeGAN网络的输入进行模型训练，迭代10 000次后，生成虚拟故障样本，TimeGAN网络的部分参数设置如表2所示.

为了验证1.2节所述基于马氏距离的虚拟样本质量评估效果，分别采用欧氏距离、三角隶属度函数法及马氏距离作为距离准则，对生成的虚拟故障数据进行筛选，将不可信样本剔除，最终在每种故障模式下形成容量为400的故障样本集. 为了直观展示不同故障数据的分布情况，采用PCA将数据降至三维空间，绘制各故障模式的故障样本的三维得分分量，如图3所示. 图3(b)~3(d)分别表示采用欧式距离、三角隶属度函数法及马氏距离作为距离准则的筛选结果，图3(e)为不同故障模式下的全部真实故障样本，即文件d01、d02、d06、d07、d12和d14的全部故障样本.

由图3可知，使用TimeGAN网络生成的虚拟数据与真实故障样本的分布基本吻合. 在对虚拟样本进行基于欧氏距离和三角隶属函数的筛选后，仍然存在一些离群点，如图3(b)和3(c)中黑色椭圆所示. 然而，使用马氏距离对虚拟样本进行筛选后，相同故障模式的样本分布较为集中.

按照3.1节所示故障数据设置情况，大子样故障1的历史故障样本为400个，而小子样故障的历史故障样本甚至低至20个. 在不进行数据增强的情况下直接提取故障子空间$ {F}_{\text{I}}^{1}、{F}_{\text{I}}^{2}、{F}_{\text{I}}^{6}、{F}_{\text{I}}^{7}、{F}_{\text{I}}^{12}和{F}_{\text{I}}^{14} $，并对故障样本的${{{I}}^2}$统计量进行重构，可得到如图4所示的重构恢复率混淆矩阵，混淆矩阵中的元素表示纵坐标对应的子空间重构横坐标对应的故障样本时所得的重构恢复率，比如，故障7与子空间12所对应的值为0.3450，表示用故障12的子空间重构故障7，可以将34.50%故障样本的统计量恢复到控制限以下.

由图4可知，在故障样本不平衡条件下进行故障分类主要存在以下2个问题.

1)在部分故障模式下，用自身的子空间重构故障样本，重构恢复率较低. 例如故障12和故障14的重构恢复率仅为30.13%和46.50%，具体重构情况如图5所示. 图中，蓝色线条代表样本原始统计量，红色粗线条代表相应的控制限，绿色点代表重构后的统计量. 为了清晰展示故障样本（即样本161~960）经重构后的统计量与控制限之间的关系，对图中局部细节进行放大，用紫色虚线框表示. 可以看出，重构后依旧有大量样本的统计量高于控制限. 这说明该故障子空间对正样本的识别能力较弱，即故障子空间表征系统属性的能力有限，本质原因为真实历史故障样本较少，难以表示故障的全部特性.

2)故障子空间与故障样本之间存在一对多的关系，即故障样本可同时被多个故障子空间恢复. 例如图4中，故障7中31%、34.5%的样本可被子空间$ F_{\text{I}}^1、F_{\text{I}}^{12} $恢复到控制限以下，具体重构情况如图6所示，说明故障子空间对负样本的区分能力不足.

为了克服不平衡故障数据对故障分类性能的影响，采用TimeGAN对小子样故障数据进行数据增强，将所有故障样本增至400个，分别在400个故障样本中提取故障子空间$ {F}_{\text{I}}^{1}、{F}_{\text{I}}^{2}、{F}_{\text{I}}^{6}、{F}_{\text{I}}^{7}、{F}_{\text{I}}^{12}和{F}_{\text{I}}^{14} $，并对${{{I}}^2}$统计量进行重构，可得如图7(a)所示的重构恢复率混淆矩阵. 为了对比数据增强故障分类性能，从各故障模式下全部真实故障样本中提取故障子空间，即从d01中全部480个样本提取故障子空间$F_{\text{I}}^1$，依此类推，分别提取$ {F}_{\text{I}}^{2}、{F}_{\text{I}}^{6}、{F}_{\text{I}}^{7}、 {F}_{\text{I}}^{12}和{F}_{\text{I}}^{14} $，同样对${{{I}}^2}$统计量进行重构，重构恢复率如图7(b)混淆矩阵所示.

对比图4和图7(a)，由于故障1为大子样故障，并没有对故障样本进行数据增强，因此图7(a)中子空间1对应的重构恢复率与图4中保持一致；其他故障模式经过数据增强后，自重构恢复率均有不同幅度的提高，尤其是故障12和故障14，分别增高了38.74个百分点和45.63个百分点，具体故障重构情况如图8所示. 对比图7(a)与图7(b)，除故障12外，数据增强子空间与真实故障数据子空间的自重构恢复率相差较小，说明经过TimeGAN网络生成的虚拟故障数据与真实数据分布相似，可以准确表征故障特性.

为了验证本研究所提TimeGAN数据增强方法的有效性，采用WGAN-GP、SMOTE作为对比算法，对表1的历史故障样本进行数据增强，然后用故障重构方法进行故障分类，并将分类结果与采用不平衡历史故障样本和全部真实故障样本提取子空间的分类结果进行对比. 分别计算6种故障模式下的精确率、召回率和F1，结果如表3所示. 可以看出，相比直接采用不平衡历史故障样本提取子空间方法，使用数据增强方法可以有效提升故障分类性能. 在3种数据增强方法中，TimeGAN方法在大多数故障模式下可以获得优于其他2种数据增强方法的故障分类性能，尤其在故障7、12和14中，其精确率、召回率和F1得分均高于其他数据增强方法的. 在使用TimeGAN进行数据增强后，故障样本的故障分类性能指标与使用全部真实故障样本的基本相似，只有部分故障的诊断性能略逊于使用全部真实故障样本的，如故障7的召回率和F1得分. 但总体而言，各性能指标之间的差异较小，均低于4个百分点. 为了展示各算法的分类性能，鉴于F1得分能够综合体现模型对正负样本的识别能力和区分能力，将各算法在不同故障模式下的F1得分绘制如图9所示.

为了更加直观地展示故障分类性能，将各模式下的验证故障样本串联组成在线故障样本数据集，共由4 800个样本组成，采用基于TimeGAN数据增强的故障分类方法进行诊断，分别采用6种模式的故障子空间进行重构，可以得到如图10所示的故障重构图. 按照2.3节所述故障分类策略可以得到如图11所示诊断结果图.

由图10可知，对于所有故障类型，使用该故障对应子空间进行重构后的统计量达到最小. 为了清晰表示，对图中局部进行放大显示. 例如对故障1（即样本1~800），最上方的蓝色点状线表示原始统计量，而红色实线表示用不同故障子空间重构后的统计量，可以看出只有用子空间1进行故障重构可以将大多数样本的统计量恢复至控制限以下，结果与图7(a)所示一致. 故障样本的在线分类结果如图11所示. 图中，${t}$表示采样时刻，${k}$表示分类标签. 可以看出，2.3节所述故障分类策略对于大多数样本都可以正确判断. 对于故障7（即样本2401~3200），存在部分误判，将样本诊断为故障12，表明子空间7对正样本的识别能力和子空间12对负样本的区分能力较弱，该结果也间接解释了表3中故障7的召回率和故障12的精确率较低的原因. 整体而言，本研究所提基于数据增强的复杂过程故障分类方法具有较高的诊断性能.

针对实际工业过程中普遍存在的故障样本稀缺、样本不平衡问题，提出基于TimeGAN数据增强的非线性复杂过程故障分类方法. 研究结果显示，相较于直接使用不平衡数据进行故障重构，采用TimeGAN网络数据增强方法能显著提高诊断准确性. 与WGAN-GP和SMOTE数据增强方法相比，TimeGAN方法在提升诊断精度方面表现出色. 此外，将原始数据映射到核空间后，在核空间提取的故障子空间具有更强的故障特征表征和分辨能力.

由于故障重构是基于多元统计分析的故障诊断方法，诊断精度依赖于故障监测模型的泛化能力. 本研究重点关注非线性非高斯过程的故障诊断，然而实际工业过程通常具有混合特性. 因此，建立更具泛化能力的故障诊断模型是下一步的研究计划.