基于深度学习的振动信号自适应模态分解与瞬时频率估计
郭一冰
a
,
b
,
c
a
,
b
,
c
^(a,b,c) { }^{a, b, c}
a
,
b
,
c
,
∗
a
,
b
,
c
,
∗
^(a,b,c,**) { }^{a, b, c, *}
a,b,c
a,b,c
^("a,b,c ") { }^{\text {a,b,c }}
d,e
d,e
^("d,e ") { }^{\text {d,e }}
a
a
^("a ") { }^{\text {a }}
b
b
^(b) { }^{\mathrm{b}}
c
c
^("c ") { }^{\text {c }}
d
d
^("d ") { }^{\text {d }}
e
e
^("e ") { }^{\text {e }}
文章信息 关键词： 时频分析 非线性和非平稳 自适应模态分解 单分量 瞬时频率 深度学习
摘要 在时变环境和载荷下，结构的振动监测数据呈现多频率成分的非平稳特性，极端载荷下还可能表现出显著非线性。时频分析方法非常适合分析这类非线性和非平稳信号。然而，大型结构的振动信号通常包含多种频率成分，使得现有自适应模态分解方法难以有效实现模态分离，且易受虚假单分量干扰。虚假单分量会严重降低瞬时频率的估计精度。针对模态分解中的模态混叠和虚假单分量问题，本文提出了一种基于深度学习的信号自适应模态分解与瞬时频率估计数据驱动方法。通过将时频分析问题转化为深度神经网络的优化任务，设计了内置冗余自适应滤波器与稀疏正则化机制，从而能准确估计本征模态函数和时变频率。 所提方法通过数值非平稳信号及长缆索模型实验振动信号验证。结果表明，该方法能自适应地将信号分解为单分量且无虚假分量，同时实现对应瞬时频率的高精度估计。在频谱重叠的无噪声信号案例中，瞬时频率估计精度达
98.52
%
98.52
%
98.52% 98.52 \%
4.90
%
4.90
%
4.90% 4.90 \%
1. 引言 结构健康监测（SHM）是理解工程结构演化规律、保障结构安全最直接先进的方法，已在众多大型结构中得到广泛应用与研究[1-3]。监测数据的有效分析是结构参数识别、损伤检测与状态评估的基础[4]。在时...
术语表
A
自适应滤波器
B
傅里叶变换矩阵
C
中心矩阵
D
带宽矩阵
e
e
e e 相对估计误差
ε
ε
epsi \varepsilon 少量
f
f
f f 合成信号
f
(
t
)
f
(
t
)
f(t) \mathbf{f}(t) 结构振动信号
f
^
r
e
c
(
ω
)
f
^
r
e
c
(
ω
)
widehat(f)_(rec)(omega) \widehat{\mathbf{f}}_{r e c}(\omega) 重构信号频谱
F
F
F \mathscr{F} 傅里叶变换
F
非线性变换
F
实测索力
F
^
F
^
widehat(F) \widehat{\mathrm{F}} 信号的频谱矩阵
F
^
IMF
F
^
IMF
widehat(F)_("IMF ") \widehat{\mathbf{F}}_{\text {IMF }} 本征模态函数的频谱矩阵
F
^
red
F
^
red
widehat(F)_("red ") \widehat{\mathbf{F}}_{\text {red }} 冗余单分量矩阵
H
滤波器组
IF
估计缆索力
I
r
,
k
I
r
,
k
I_(r,k) I_{r, k} 正交值
K
滤波器数量
L
L
L L 损失函数
L
c
L
c
L_(c) L_{c} 电缆长度
m
m
m m 质量密度
M
本征模态函数数量
N
N
N N 数据点数量
r
(
t
)
r
(
t
)
r(t) \mathbf{r}(t) 残差
u
k
(
t
)
u
k
(
t
)
u_(k)(t) \mathbf{u}_{k}(t) 本征模态函数
ω
~
(
t
)
ω
~
(
t
)
widetilde(omega)(t) \widetilde{\boldsymbol{\omega}}(t) 瞬时频率
W
稀疏权重矩阵
z
k
(
t
)
z
k
(
t
)
z_(k)(t) \mathrm{z}_{\mathrm{k}}(t) 解析信号
Z
解析信号矩阵
Z
^
Z
^
widehat(Z) \widehat{\mathbf{Z}} 解析信号的频谱矩阵
Nomenclature
A adaptive filters
B Fourier transform matrix
C center matrix
D bandwidth matrix
e relative estimation error
epsi a small number
f synthetic signal
f(t) structural vibration signal
widehat(f)_(rec)(omega) spectrum of reconstructed signal
F Fourier transform
F nonlinear transformation
F measured cable force
widehat(F) spectrum matrix of signal
widehat(F)_("IMF ") spectrum matrix of intrinsic mode functions
widehat(F)_("red ") redundant mono-component matrix
H filter set
IF estimated cable force
I_(r,k) orthogonal value
K number of filters
L loss function
L_(c) cable length
m mass density
M number of intrinsic mode functions
N number of data points
r(t) residual
u_(k)(t) intrinsic mode function
widetilde(omega)(t) instantaneous frequency
W sparse weights matrix
z_(k)(t) analytic signal
Z analytic signal matrix
widehat(Z) spectrum matrix of analytic signal | Nomenclature | |
| :---: | :---: |
| A | adaptive filters |
| B | Fourier transform matrix |
| C | center matrix |
| D | bandwidth matrix |
| $e$ | relative estimation error |
| $\varepsilon$ | a small number |
| $f$ | synthetic signal |
| $\mathbf{f}(t)$ | structural vibration signal |
| $\widehat{\mathbf{f}}_{r e c}(\omega)$ | spectrum of reconstructed signal |
| $\mathscr{F}$ | Fourier transform |
| F | nonlinear transformation |
| F | measured cable force |
| $\widehat{\mathrm{F}}$ | spectrum matrix of signal |
| $\widehat{\mathbf{F}}_{\text {IMF }}$ | spectrum matrix of intrinsic mode functions |
| $\widehat{\mathbf{F}}_{\text {red }}$ | redundant mono-component matrix |
| H | filter set |
| IF | estimated cable force |
| $I_{r, k}$ | orthogonal value |
| K | number of filters |
| $L$ | loss function |
| $L_{c}$ | cable length |
| $m$ | mass density |
| M | number of intrinsic mode functions |
| $N$ | number of data points |
| $\mathbf{r}(t)$ | residual |
| $\mathbf{u}_{k}(t)$ | intrinsic mode function |
| $\widetilde{\boldsymbol{\omega}}(t)$ | instantaneous frequency |
| W | sparse weights matrix |
| $\mathrm{z}_{\mathrm{k}}(t)$ | analytic signal |
| Z | analytic signal matrix |
| $\widehat{\mathbf{Z}}$ | spectrum matrix of analytic signal |
在多变环境与荷载作用下，结构振动监测数据呈现多频率成分的非平稳特性。某些复杂环境效应和荷载条件下，振动数据往往具有非线性特征，例如重载车辆作用下的桥梁拉索振动数据及强震作用下受损结构的振动数据[5-7]。这些信号的瞬时频率为结构健康状态评估与损伤检测提供了丰富信息，因此精确的时频分析估计对结构健康监测(SHM)至关重要[8-10]。此类非线性非平稳 SHM 信号的时频分析关键在于精准的自适应模态分解过程，该过程可将振动信号分解为一组称为本征模态函数(IMFs)的单分量项[11]。每个 IMF 对应一个独特的时变频率成分，进而可通过基于希尔伯特变换的解析信号
[
12
,
13
]
[
12
,
13
]
[12,13] [12,13] 传统的基于基展开的模式分解方法，如傅里叶变换，采用形式已知的正交基对信号进行分解，存在灵活性差和适应性不足的缺点[14,15]，无法为非线性非平稳信号提供足够的精度。为此，Huang 等人[11]提出了经验模态分解(EMD)，其具有自适应性和操作简便的优势，已被广泛应用于结构健康监测(SHM)中的模态参数识别
[
16
,
17
]
[
16
,
17
]
[16,17] [16,17] 这些虚假分解的单分量存在于 EEMD 和 LMD[35]中。预先选定的分解模式会导致 VMD[27,36]产生过分割或欠分割的 IMF，而不恰当的初始值会影响 ASTFA 的收敛性并导致虚假单分量。然而，对于大型结构而言，振动信号通常包含丰富的模态信息。在周围环境存在显著背景噪声的情况下，很难从振动信号中精确分解这些模态，现有方法容易获得虚假单分量，从而为有效的瞬时频率估计带来挑战。
自适应模态分解和瞬时频率估计的时频分析理论上是受非线性目标函数控制的优化问题。因此，深度学习网络（DNNs）等机器学习方法因其更强大的全局最优性和逼近任意非线性模型的特性而得以应用[37-39]。在结构健康监测（SHM）领域，已有研究将机器学习与结构模态识别、系统辨识、损伤检测等任务相结合[40-42]。Liu 等[43]提出了一种结合独立特性的机器学习方法用于结构模态参数识别；Guo 等[44]在损失函数中嵌入物理信息的长短期记忆神经网络，用于振动台系统辨识；Hesser 等[45]利用二维卷积神经网络进行特征提取，并通过迁移学习训练数据，实现了声发射源的精准识别。前期自适应稀疏时频分析的研究工作已表明机器学习在精确瞬时频率估计方面的良好潜力[46]。 然而，该方法中目标 IMF 的手动筛选及基于傅里叶的字典使用，对于非三角信号的自适应性存在明显不足。
为提升自适应性并克服模态混叠与虚假单分量问题，本文提出一种基于深度学习的信号自适应稀疏模态分解与瞬时频率估计新方法。通过将时频分析转化为深度学习中的优化问题，构建了内置自适应冗余滤波器及稀疏正则化的深度神经网络。借助 DNN 强大的优化能力，实现了 IMFs 的精确分解与高分辨率瞬时频率获取。
2. 面向时频分析的深度学习方法 2.1 理论基础 考虑一个非线性非平稳的结构振动信号
f
(
t
)
f
(
t
)
f(t) \mathbf{f}(t)
N
N
N N
M
M
M M
{
minimize
M
subject to :
f
(
t
)
=
∑
k
=
1
M
u
k
(
t
)
+
r
(
t
)
,
t
∈
R
,
M
∈
R
,
f
(
t
)
∈
R
N
minimize
M
subject to :
f
(
t
)
=
∑
k
=
1
M
u
k
(
t
)
+
r
(
t
)
,
t
∈
R
,
M
∈
R
,
f
(
t
)
∈
R
N
{[minimize M],[" subject to : "f(t)=sum_(k=1)^(M)u_(k)(t)+r(t)","t inR","M inR","f(t)inR^(N)]:} \left\{\begin{array}{l}
\operatorname{minimize} M \\
\text { subject to : } \mathbf{f}(t)=\sum_{k=1}^{M} \mathbf{u}_{k}(t)+\mathbf{r}(t), t \in \mathbb{R}, M \in \mathbb{R}, \mathbf{f}(t) \in \mathbb{R}^{N}
\end{array}\right. 其中
f
(
t
)
f
(
t
)
f(t) \mathbf{f}(t)
M
M
M M
u
k
(
t
)
u
k
(
t
)
u_(k)(t) \mathbf{u}_{k}(t)
r
(
t
)
r
(
t
)
r(t) \mathbf{r}(t)
k
k
k k 模态确定后，通过频谱方法获取每个 IMF 的解析形式，仅选取其傅里叶变换的正频率成分[8,9]。因此，每个 IMF 解析信号的频谱为：
z
^
k
(
ω
)
=
{
u
^
k
(
ω
)
,
for
ω
>
0
0
,
otherwise
z
^
k
(
ω
)
=
u
^
k
(
ω
)
,
for
ω
>
0
0
,
otherwise
widehat(z)_(k)(omega)={[ widehat(u)_(k)(omega)","," for "omega > 0],[0","," otherwise "]:} \widehat{\mathbf{z}}_{k}(\omega)= \begin{cases}\widehat{\mathbf{u}}_{k}(\omega), & \text { for } \omega>0 \\ 0, & \text { otherwise }\end{cases} 其中
u
^
k
(
ω
)
=
F
[
u
k
(
t
)
]
u
^
k
(
ω
)
=
F
u
k
(
t
)
widehat(u)_(k)(omega)=F[u_(k)(t)] \widehat{\mathbf{u}}_{k}(\omega)=\mathscr{F}\left[\mathbf{u}_{k}(t)\right]
F
F
F \mathscr{F}
z
k
(
t
)
=
F
−
1
[
z
^
k
(
ω
)
]
=
z
k
,
r
e
(
t
)
+
j
⋅
z
k
,
i
m
(
t
)
ω
~
k
(
t
)
=
d
d
t
(
arctan
z
k
,
i
m
(
t
)
z
k
,
r
e
(
t
)
)
z
k
(
t
)
=
F
−
1
z
^
k
(
ω
)
=
z
k
,
r
e
(
t
)
+
j
⋅
z
k
,
i
m
(
t
)
ω
~
k
(
t
)
=
d
d
t
arctan
z
k
,
i
m
(
t
)
z
k
,
r
e
(
t
)
{:[z_(k)(t)=F^(-1)[ widehat(z)_(k)(omega)]=z_(k,re)(t)+j*z_(k,im)(t)],[ widetilde(omega)_(k)(t)=(d)/(dt)(arctan((z_(k,im)(t))/(z_(k,re)(t))))]:} \begin{aligned}
& \mathbf{z}_{k}(t)=\mathscr{F}^{-1}\left[\widehat{\mathbf{z}}_{k}(\omega)\right]=\mathbf{z}_{k, r e}(t)+j \cdot \mathbf{z}_{k, i m}(t) \\
& \widetilde{\boldsymbol{\omega}}_{k}(t)=\frac{d}{d t}\left(\arctan \frac{\mathbf{z}_{k, i m}(t)}{\mathbf{z}_{k, r e}(t)}\right)
\end{aligned} 其中
z
k
,
r
e
(
t
)
z
k
,
r
e
(
t
)
z_(k,re)(t) \mathbf{z}_{k, r e}(t)
z
k
,
i
m
(
t
)
z
k
,
i
m
(
t
)
z_(k,im)(t) \mathbf{z}_{k, i m}(t)
k
k
k k
F
−
1
F
−
1
F^(-1) \mathscr{F}^{-1}
ω
~
k
(
t
)
ω
~
k
(
t
)
widetilde(omega)_(k)(t) \widetilde{\boldsymbol{\omega}}_{k}(t)
k
k
k k 上述问题的核心在于如何有效分解模态而不产生虚假单分量信号，因为分解结果直接影响瞬时频率的准确性。每个模态中应仅呈现单一特定频率成分。式(1)中的模态分解问题可转化为如式(5)所示的非凸最小二乘优化问题：
minimize
M
,
subject to :
‖
f
(
t
)
−
∑
k
=
1
M
u
k
(
t
)
‖
2
2
⩽
ε
minimize
M
,
subject to :
f
(
t
)
−
∑
k
=
1
M
u
k
(
t
)
2
2
⩽
ε
minimize M," subject to : "||f(t)-sum_(k=1)^(M)u_(k)(t)||_(2)^(2) <= epsi \operatorname{minimize} M, \text { subject to : }\left\|\mathbf{f}(t)-\sum_{k=1}^{M} \mathbf{u}_{k}(t)\right\|_{2}^{2} \leqslant \varepsilon 其中
ε
ε
epsi \varepsilon
F
F
F F
ω
~
(
t
)
=
F
(
f
(
t
)
)
,
f
(
t
)
∈
R
N
,
t
∈
[
0
,
T
]
ω
~
(
t
)
=
F
(
f
(
t
)
)
,
f
(
t
)
∈
R
N
,
t
∈
[
0
,
T
]
widetilde(omega)(t)=F(f(t)),quadf(t)inR^(N),quad t in[0,T] \widetilde{\boldsymbol{\omega}}(t)=F(\mathbf{f}(t)), \quad \mathbf{f}(t) \in \mathbb{R}^{N}, \quad t \in[0, T] 上述内容构成了基于深度学习的数据驱动型信号自适应稀疏模态分解与瞬时频率估计技术的基石，其程序流程图如图 1 所示。首先，将结构健康监测系统（SHM）采集的振动信号
f
(
t
)
f
(
t
)
f(t) \mathbf{f}(t)
F
^
F
^
widehat(F) \widehat{\mathbf{F}}
A
A
A \mathbf{A}
W
W
W \mathbf{W}
F
^
I
M
F
F
^
I
M
F
widehat(F)_(IMF) \widehat{\mathbf{F}}_{I M F}
Z
^
Z
^
widehat(Z) \widehat{\mathbf{Z}}
ω
~
(
t
)
ω
~
(
t
)
widetilde(omega)(t) \widetilde{\boldsymbol{\omega}}(t) 图1. 程序流程图
2.3 深度神经网络构建 式(5)中的优化问题类似于深度神经网络的目标函数。所设计的八层全连接 DNN 优化器细节如图 2 所示，其实现了时频分析中的两个关键步骤：从原始信号到本征模态函数(IMFs)的自适应稀疏模态分解，以及 IMFs 的瞬时频率估计。
自适应稀疏模态分解在设计的 DNN 前四层完成。时间序列信号首先输入至输入层，随后通过傅里叶基矩阵转换至频域。其次，在第二层将信号分离为实部与虚部后，引入具有可变中心频率和带宽的自适应滤波器组，并筛选单分量以实现自适应滤波。为避免虚假单分量，在第三层输出端施加稀疏正则化（亦称
ℓ
1
ℓ
1
ℓ_(1) \ell_{1} 固有模态函数（IMFs）的瞬时频率估计通过第 5、6、7 和 8 层的希尔伯特变换实现。首先，对解析信号的频谱进行正频率滤波。接着，乘以逆傅里叶矩阵并合并结果，从解析信号的相位中计算出瞬时频率。
在该设计的深度神经网络中，三类权重至关重要。前两类是不同自适应滤波器的峰值中心
C
C
C \mathbf{C}
D
D
D \mathbf{D}
W
W
W \mathbf{W} 具体细节如下。对于设计的深度神经网络的第二层，执行如式(7)所示的傅里叶变换，其中
F
^
F
^
widehat(F) \widehat{\mathbf{F}}
f
(
t
)
f
(
t
)
f(t) \mathbf{f}(t)
B
B
B \mathbf{B}
F
^
=
f
(
t
)
⋅
B
=
[
f
^
1
(
ω
)
,
f
^
2
(
ω
)
,
…
,
f
^
N
(
ω
)
]
∈
R
N
×
N
B
=
(
b
r
,
k
)
r
=
1
,
…
,
N
,
k
=
1
,
…
,
N
F
^
=
f
(
t
)
⋅
B
=
f
^
1
(
ω
)
,
f
^
2
(
ω
)
,
…
,
f
^
N
(
ω
)
∈
R
N
×
N
B
=
b
r
,
k
r
=
1
,
…
,
N
,
k
=
1
,
…
,
N
{:[ widehat(F)=f(t)*B=[ widehat(f)_(1)(omega), widehat(f)_(2)(omega),dots, widehat(f)_(N)(omega)]inR^(N xx N)],[B=(b_(r,k))_(r=1,dots,N,k=1,dots,N)]:} \begin{aligned}
& \widehat{\mathbf{F}}=\mathbf{f}(t) \cdot \mathbf{B}=\left[\widehat{\mathbf{f}}_{1}(\omega), \widehat{\mathbf{f}}_{2}(\omega), \ldots, \widehat{\mathbf{f}}_{N}(\omega)\right] \in \mathbb{R}^{N \times N} \\
& \mathbf{B}=\left(b_{r, k}\right)_{r=1, \ldots, N, k=1, \ldots, N}
\end{aligned} 在设计的 DNN 第三层中，引入了冗余自适应滤波器以提取频谱单分量。此步骤中，冗余滤波器提取的单分量数量多于实际固有模态函数（IMFs）。滤波是一种提取有意义信号成分并去除不理想部分的方法[47]。本方法的一个主要差异在于采用了一组具有不同中心频率和带宽的自适应滤波器，其设计基于波数域中的 Meyer 尺度函数，如式(8)所示：
Φ
=
{
1
,
|
l
|
⩽
1
1
2
(
1
−
cos
π
l
)
,
1
<
|
l
|
<
2
0
,
otherwise
Φ
=
1
,
|
l
|
⩽
1
1
2
(
1
−
cos
π
l
)
,
1
<
|
l
|
<
2
0
,
otherwise
Phi={[1",",|l| <= 1],[(1)/(2)(1-cos pi l)",",1 < |l| < 2],[0","," otherwise "]:} \Phi=\left\{\begin{array}{lr}
1, & |l| \leqslant 1 \\
\frac{1}{2}(1-\cos \pi l), & 1<|l|<2 \\
0, & \text { otherwise }
\end{array}\right. 其中
l
l
l l
K
(
K
>
M
)
K
(
K
>
M
)
K(K > M) K(K>M)
A
A
A \mathbf{A} 图 2. 所提方法的 DNN 结构示意图。
A
=
m
1
m
2
2
(
1
−
cos
|
l
|
−
C
∣
−
D
D
π
)
+
(
1
−
m
1
)
m
2
A
=
[
a
^
1
(
ω
)
,
a
^
2
(
ω
)
,
…
,
a
^
K
(
ω
)
]
∈
R
N
×
K
,
C
=
[
c
1
,
c
2
,
…
,
c
K
]
∈
R
K
,
D
=
[
d
1
,
d
2
,
…
,
d
K
]
∈
R
K
A
=
m
1
m
2
2
1
−
cos
|
l
|
−
C
∣
−
D
D
π
+
1
−
m
1
m
2
A
=
a
^
1
(
ω
)
,
a
^
2
(
ω
)
,
…
,
a
^
K
(
ω
)
∈
R
N
×
K
,
C
=
c
1
,
c
2
,
…
,
c
K
∈
R
K
,
D
=
d
1
,
d
2
,
…
,
d
K
∈
R
K
{:[A=(m_(1)m_(2))/(2)(1-cos((|l|-C∣-D)/(D)pi))+(1-m_(1))m_(2)],[A=[ widehat(a)_(1)(omega), widehat(a)_(2)(omega),dots, widehat(a)_(K)(omega)]inR^(N xx K)","C=[c_(1),c_(2),dots,c_(K)]inR^(K)","D=[d_(1),d_(2),dots,d_(K)]inR^(K)]:} \begin{aligned}
& \mathbf{A}=\frac{\mathbf{m}_{1} \mathbf{m}_{2}}{2}\left(1-\cos \frac{|l|-\mathbf{C} \mid-\mathbf{D}}{\mathbf{D}} \pi\right)+\left(1-\mathbf{m}_{1}\right) \mathbf{m}_{2} \\
& \mathbf{A}=\left[\widehat{\mathbf{a}}_{1}(\omega), \widehat{\mathbf{a}}_{2}(\omega), \ldots, \widehat{\mathbf{a}}_{K}(\omega)\right] \in \mathbb{R}^{N \times K}, \mathbf{C}=\left[c_{1}, c_{2}, \ldots, c_{K}\right] \in \mathbb{R}^{K}, \mathbf{D}=\left[d_{1}, d_{2}, \ldots, d_{K}\right] \in \mathbb{R}^{K}
\end{aligned} 式中
a
^
i
(
ω
)
,
i
=
1
,
…
,
K
a
^
i
(
ω
)
,
i
=
1
,
…
,
K
widehat(a)_(i)(omega),i=1,dots,K \widehat{\mathbf{a}}_{i}(\omega), i=1, \ldots, K
i
i
i i
l
l
l l
m
1
=
sgn
(
l
l
−
C
∣
)
+
1
2
m
1
=
sgn
(
l
l
−
C
∣
)
+
1
2
m_(1)=(sgn(ll-C∣)+1)/(2) \mathbf{m}_{1}=\frac{\operatorname{sgn}(l l-\mathbf{C} \mid)+1}{2}
m
2
=
sgn
(
D
−
|
|
l
|
−
C
|
)
+
1
2
m
2
=
sgn
(
D
−
|
|
l
|
−
C
|
)
+
1
2
m_(2)=(sgn(D-||l|-C|)+1)/(2) \mathbf{m}_{2}=\frac{\operatorname{sgn}(\mathbf{D}-||l|-\mathbf{C}|)+1}{2}
C
C
C \mathbf{C}
D
D
D \mathbf{D}
K
K
K K 冗余单分量矩阵
F
^
red
F
^
red
widehat(F)_("red ") \widehat{\mathbf{F}}_{\text {red }}
f
^
i
(
ω
)
f
^
i
(
ω
)
widehat(f)_(i)(omega) \widehat{\mathbf{f}}_{i}(\omega)
i
i
i i
i
i
i i
F
^
red
=
A
⋅
F
^
=
[
f
^
1
(
ω
)
,
f
^
2
(
ω
)
,
…
,
f
^
K
(
ω
)
]
∈
C
N
×
K
F
^
red
=
A
⋅
F
^
=
f
^
1
(
ω
)
,
f
^
2
(
ω
)
,
…
,
f
^
K
(
ω
)
∈
C
N
×
K
widehat(F)_("red ")=A* widehat(F)=[ widehat(f)_(1)(omega), widehat(f)_(2)(omega),dots, widehat(f)_(K)(omega)]inC^(N xx K) \widehat{\mathbf{F}}_{\text {red }}=\mathbf{A} \cdot \widehat{\mathbf{F}}=\left[\widehat{\mathbf{f}}_{1}(\omega), \widehat{\mathbf{f}}_{2}(\omega), \ldots, \widehat{\mathbf{f}}_{K}(\omega)\right] \in \mathbb{C}^{N \times K} 在第四层应用稀疏矩阵以确保信号重建所需 IMF 数量最小化，此步骤将剔除冗余不必要的单分量。重建信号
f
^
rec
(
ω
)
f
^
rec
(
ω
)
widehat(f)_("rec ")(omega) \widehat{\mathbf{f}}_{\text {rec }}(\omega)
F
^
I
M
F
=
F
^
red
⋅
W
f
^
r
e
c
(
ω
)
=
∑
k
=
1
M
F
^
I
M
F
(
r
,
k
)
,
f
^
rec
(
ω
)
∈
R
N
F
^
I
M
F
=
F
^
red
⋅
W
f
^
r
e
c
(
ω
)
=
∑
k
=
1
M
F
^
I
M
F
(
r
,
k
)
,
f
^
rec
(
ω
)
∈
R
N
{:[ widehat(F)_(IMF)= widehat(F)_("red ")*W],[ widehat(f)_(rec)(omega)=sum_(k=1)^(M) widehat(F)_(IMF)(r","k)"," widehat(f)_("rec ")(omega)inR^(N)]:} \begin{aligned}
& \widehat{\mathbf{F}}_{I M F}=\widehat{\mathbf{F}}_{\text {red }} \cdot \mathbf{W} \\
& \widehat{\mathbf{f}}_{r e c}(\omega)=\sum_{k=1}^{M} \widehat{\mathbf{F}}_{I M F}(r, k), \widehat{\mathbf{f}}_{\text {rec }}(\omega) \in \mathbb{R}^{N}
\end{aligned} 其中
W
W
W \mathbf{W}
F
^
IMF
F
^
IMF
widehat(F)_("IMF ") \widehat{\mathbf{F}}_{\text {IMF }}
r
r
r r
k
k
k k
F
^
IMF
F
^
IMF
widehat(F)_("IMF ") \widehat{\mathbf{F}}_{\text {IMF }} 本方法的核心在于从原始信号中获取自适应稀疏 IMF 分量，进而精确识别瞬时频率。该优化过程通过最小化重建时序信号（即 IMF 分量之和）与输入原始信号之间的差异来实现，其机器学习任务的最优化问题可表述为式(13)，精确表达式见式(14)：
min
(
‖
f
(
t
)
−
f
rec
(
t
)
‖
2
2
+
‖
W
‖
1
)
min
A
∈
C
N
×
K
,
W
∈
R
N
‖
f
−
A
⋅
F
^
⋅
W
⋅
B
−
1
‖
2
2
+
‖
W
‖
1
min
f
(
t
)
−
f
rec
(
t
)
2
2
+
‖
W
‖
1
min
A
∈
C
N
×
K
,
W
∈
R
N
f
−
A
⋅
F
^
⋅
W
⋅
B
−
1
2
2
+
‖
W
‖
1
{:[ min(||f(t)-f_("rec ")(t)||_(2)^(2)+||W||_(1))],[min_(AinC^(N xx K),WinR^(N))||f-A*( widehat(F))*W*B^(-1)||_(2)^(2)+||W||_(1)]:} \begin{aligned}
& \min \left(\left\|\mathbf{f}(t)-\mathbf{f}_{\text {rec }}(t)\right\|_{2}^{2}+\|\mathbf{W}\|_{1}\right) \\
& \min _{\mathbf{A} \in \mathbb{C}^{N \times K}, \mathbf{W} \in \mathbb{R}^{N}}\left\|\mathbf{f}-\mathbf{A} \cdot \widehat{\mathbf{F}} \cdot \mathbf{W} \cdot \mathbf{B}^{-1}\right\|_{2}^{2}+\|\mathbf{W}\|_{1}
\end{aligned} 其中
‖
f
(
t
)
−
f
r
e
c
(
t
)
‖
2
2
=
∑
r
=
1
N
(
f
(
t
r
)
−
f
r
e
c
(
t
r
)
)
2
f
(
t
)
−
f
r
e
c
(
t
)
2
2
=
∑
r
=
1
N
f
t
r
−
f
r
e
c
t
r
2
||f(t)-f_(rec)(t)||_(2)^(2)=sum_(r=1)^(N)(f(t_(r))-f_(rec)(t_(r)))^(2) \left\|\mathbf{f}(t)-\mathbf{f}_{r e c}(t)\right\|_{2}^{2}=\sum_{r=1}^{N}\left(f\left(t_{r}\right)-f_{r e c}\left(t_{r}\right)\right)^{2}
‖
W
‖
1
=
∑
i
=
1
K
|
w
i
|
‖
W
‖
1
=
∑
i
=
1
K
w
i
||W||_(1)=sum_(i=1)^(K)|w_(i)| \|\mathbf{W}\|_{1}=\sum_{i=1}^{K}\left|w_{i}\right| 为实现设计 DNN 对信号的精确模态分解，式(15)提出了一种新颖的双重损失函数：
L
=
1
N
∑
r
=
1
N
(
f
(
t
r
)
−
F
−
1
[
∑
k
=
1
M
A
⋅
F
^
⋅
W
(
r
,
k
)
]
)
2
+
∑
i
=
1
K
|
w
i
|
L
=
1
N
∑
r
=
1
N
f
t
r
−
F
−
1
∑
k
=
1
M
A
⋅
F
^
⋅
W
(
r
,
k
)
2
+
∑
i
=
1
K
w
i
L=(1)/(N)sum_(r=1)^(N)(f(t_(r))-F^(-1)[sum_(k=1)^(M)A*( widehat(F))*W(r,k)])^(2)+sum_(i=1)^(K)|w_(i)| L=\frac{1}{N} \sum_{r=1}^{N}\left(f\left(t_{r}\right)-\mathscr{F}^{-1}\left[\sum_{k=1}^{M} \mathbf{A} \cdot \widehat{\mathbf{F}} \cdot \mathbf{W}(r, k)\right]\right)^{2}+\sum_{i=1}^{K}\left|w_{i}\right| 其中第一项用于最小化输入信号与重构信号之间的重建误差，第二项则用于在冗余分量上施加最稀疏[48]的权重，从而在倒数第二层获得最稀疏的 IMF；
N
N
N N
k
k
k k
Im
[
⋅
]
Im
[
⋅
]
Im[*] \operatorname{Im}[\cdot] 所提出的方法为数据驱动型，依赖于数据训练，但其模型架构具有通用性，可适用于不同场景。对于该显式神经网络模型，训练过程同时也是获取解的优化过程。在优化过程中，我们采用随机梯度下降法（SGD）更新权重。由于损失函数是非凸函数，不仅存在一个全局最小值，还可能存在多个局部极小值，因此每一步参数更新仅使用单个样本。这种方式下，梯度与传统梯度下降法相比可能存在轻微偏差。
在完成本征模态函数（IMFs）的自适应稀疏分解后，瞬时频率估计流程被嵌入后续无需优化权重的深度神经网络层中。通过如式(16)所示在频域构建解析信号并提取正频率分量，在获得最优自适应滤波器组
A
A
A \mathbf{A}
W
W
W \mathbf{W}
Z
^
=
A
⋅
F
^
⋅
W
⋅
H
=
[
z
^
1
(
ω
)
,
z
^
2
(
ω
)
,
…
,
z
^
M
(
ω
)
]
∈
C
N
×
M
Z
^
=
A
⋅
F
^
⋅
W
⋅
H
=
z
^
1
(
ω
)
,
z
^
2
(
ω
)
,
…
,
z
^
M
(
ω
)
∈
C
N
×
M
widehat(Z)=A* widehat(F)*W*H=[ widehat(z)_(1)(omega), widehat(z)_(2)(omega),dots, widehat(z)_(M)(omega)]inC^(N xx M) \widehat{\mathbf{Z}}=\mathbf{A} \cdot \widehat{\mathbf{F}} \cdot \mathbf{W} \cdot \mathbf{H}=\left[\widehat{\mathbf{z}}_{1}(\omega), \widehat{\mathbf{z}}_{2}(\omega), \ldots, \widehat{\mathbf{z}}_{M}(\omega)\right] \in \mathbb{C}^{N \times M} 其中
Z
^
Z
^
widehat(Z) \widehat{\mathbf{Z}}
H
H
H \mathbf{H}
H
=
{
f
^
r
e
c
(
ω
)
,
for
ω
=
0
2
f
^
r
e
c
(
ω
)
,
for
1
⩽
ω
⩽
N
2
−
1
f
^
r
e
c
(
N
2
)
,
for
ω
=
N
2
0
,
for
N
2
+
1
⩽
ω
⩽
N
−
1
H
=
f
^
r
e
c
(
ω
)
,
for
ω
=
0
2
f
^
r
e
c
(
ω
)
,
for
1
⩽
ω
⩽
N
2
−
1
f
^
r
e
c
N
2
,
for
ω
=
N
2
0
,
for
N
2
+
1
⩽
ω
⩽
N
−
1
H={[ widehat(f)_(rec)(omega)","," for "omega=0],[2 widehat(f)_(rec)(omega)","," for "1 <= omega <= (N)/(2)-1],[ widehat(f)_(rec)((N)/(2))","," for "omega=(N)/(2)],[0","," for "(N)/(2)+1 <= omega <= N-1]:} \mathbf{H}= \begin{cases}\widehat{\mathbf{f}}_{r e c}(\omega), & \text { for } \omega=0 \\ 2 \widehat{\mathbf{f}}_{r e c}(\omega), & \text { for } 1 \leqslant \omega \leqslant \frac{N}{2}-1 \\ \widehat{\mathbf{f}}_{r e c}\left(\frac{N}{2}\right), & \text { for } \omega=\frac{N}{2} \\ 0, & \text { for } \frac{N}{2}+1 \leqslant \omega \leqslant N-1\end{cases} 频谱
Z
^
Z
^
widehat(Z) \widehat{\mathbf{Z}}
B
−
1
B
−
1
B^(-1) \mathbf{B}^{-1}
Z
Z
Z \mathbf{Z}
Z
=
Z
^
⋅
B
−
1
ω
~
(
t
)
=
1
2
π
d
d
t
(
arctan
Z
i
m
Z
r
e
)
Z
=
Z
^
⋅
B
−
1
ω
~
(
t
)
=
1
2
π
d
d
t
arctan
Z
i
m
Z
r
e
{:[Z= widehat(Z)*B^(-1)],[ widetilde(omega)(t)=(1)/(2pi)(d)/(dt)(arctan((Z_(im))/(Z_(re))))]:} \begin{aligned}
& \mathbf{Z}=\widehat{\mathbf{Z}} \cdot \mathbf{B}^{-1} \\
& \widetilde{\boldsymbol{\omega}}(t)=\frac{1}{2 \pi} \frac{d}{d t}\left(\arctan \frac{\mathbf{Z}_{i m}}{\mathbf{Z}_{r e}}\right)
\end{aligned} 其中
Z
i
m
Z
i
m
Z_(im) \mathbf{Z}_{i m}
Z
r
e
Z
r
e
Z_(re) \mathbf{Z}_{r e}
Z
Z
Z \mathbf{Z} 3. 示例 3.1. 示例1 具有幅频调制的合成信号 以一个包含五个调幅或调频本征模态函数(IMF)的非平稳信号为例。每个 IMF 的解析表达式均已给出，原始输入信号为
f
f
f f 图3. 示例1中用于合成信号的深度神经网络结构。
f
1
=
sin
(
500
π
t
+
cos
(
4
π
t
)
)
,
f
2
=
(
0.5
+
cos
(
2
π
t
)
)
⋅
cos
(
380
π
t
)
,
f
3
=
1
1.5
+
sin
(
2
π
t
)
⋅
cos
(
240
π
t
+
sin
(
16
π
t
)
)
,
f
4
=
sin
(
150
π
t
+
2
sin
(
6
π
t
)
)
,
f
5
=
cos
(
30
π
t
+
10
sin
(
2
π
t
)
)
,
f
1
=
sin
(
500
π
t
+
cos
(
4
π
t
)
)
,
f
2
=
(
0.5
+
cos
(
2
π
t
)
)
⋅
cos
(
380
π
t
)
,
f
3
=
1
1.5
+
sin
(
2
π
t
)
⋅
cos
(
240
π
t
+
sin
(
16
π
t
)
)
,
f
4
=
sin
(
150
π
t
+
2
sin
(
6
π
t
)
)
,
f
5
=
cos
(
30
π
t
+
10
sin
(
2
π
t
)
)
,
{:[f_(1)=sin(500 pi t+cos(4pi t))","],[f_(2)=(0.5+cos(2pi t))*cos(380 pi t)","],[f_(3)=(1)/(1.5+sin(2pi t))*cos(240 pi t+sin(16 pi t))","],[f_(4)=sin(150 pi t+2sin(6pi t))","],[f_(5)=cos(30 pi t+10 sin(2pi t))","]:} \begin{gathered}
f_{1}=\sin (500 \pi t+\cos (4 \pi t)), \\
f_{2}=(0.5+\cos (2 \pi t)) \cdot \cos (380 \pi t), \\
f_{3}=\frac{1}{1.5+\sin (2 \pi t)} \cdot \cos (240 \pi t+\sin (16 \pi t)), \\
f_{4}=\sin (150 \pi t+2 \sin (6 \pi t)), \\
f_{5}=\cos (30 \pi t+10 \sin (2 \pi t)),
\end{gathered}
f
=
f
1
+
f
2
+
f
3
+
f
4
+
f
5
f
=
f
1
+
f
2
+
f
3
+
f
4
+
f
5
f=f_(1)+f_(2)+f_(3)+f_(4)+f_(5) f=f_{1}+f_{2}+f_{3}+f_{4}+f_{5} 根据图2所示的网络结构，本例采用的深度神经网络如图3所示，其预训练网络参数列于表1。神经网络超参数基于经验选取。信号总长度为
N
N
N N
=
1024
=
1024
=1024 =1024 网络重构的时频信号如图 4 所示。根据设计的 DNN 框架，输出的重构信号是各 IMF 分量的总和。IMF 分量的提取结果如图 5 和图 6 所示，表明所提方法的分解结果更接近 IMF 的解析表达式。图 5(b)和图 6(b)展示了 EMD 的相应波形，其中提取的 IMF 分量含有额外的固有模态。对于 EMD 而言，虚假的 IMF 分量几乎不可避免，因为它们是通过三次样条插值从包络中提取的[11]。此外，EMD 的另一个局限是当某些 IMF 分量的频率过于接近时（即两个 IMF 分量的频率比高于 0.75[34]），它们难以被分离。如图 6(b)所示，
f
1
f
1
f_(1) f_{1}
f
2
(
f
1
/
f
2
=
0.76
>
0.75
)
f
2
f
1
/
f
2
=
0.76
>
0.75
f_(2)(f_(1)//f_(2)=0.76 > 0.75) f_{2}\left(f_{1} / f_{2}=0.76>0.75\right) 在详细说明 DNN 中自适应稀疏滤波器如何提取 IMF 时，图 7 展示了频域中每个 IMF 对应的滤波器波形。当训练轮次为 0 时，中心频率和带宽的初始值设置如表 2 和表 3 所示。每个滤波器的初始带宽均为 20，其中五个滤波器的初始中心频率经过粗略调整接近解析值，其余则随机选择。图 7 表明自适应滤波器能同步进行滤波操作，每个滤波器独立提取一个 IMF，且滤波器的最终中心频率接近解析解，与初始值相差较大。此外，尽管初始赋予相同的滤波器带宽，但根据各滤波器提取的 IMF 带宽特性，其最终带宽值会动态调整。这些特性验证了所提方法具有良好的自适应能力。
合成信号
f
f
f f
IMF
3
IMF
3
IMF_(3) \mathrm{IMF}_{3}
IMF
1
IMF
1
IMF_(1) \mathrm{IMF}_{1}
IMF
3
IMF
3
IMF_(3) \mathrm{IMF}_{3}
IMF
1
IMF
1
IMF_(1) \mathrm{IMF}_{1}
99.95
%
99.95
%
99.95% 99.95 \%
98.55
%
98.55
%
98.55% 98.55 \% 冗余自适应滤波器的稀疏系数如图 9(b)所示。在 30 个训练周期时，几乎所有权重已完成优化且方差极小，表明优化接近完成。当周期达到 100 时，权重呈现最稀疏分布。仅五个滤波器的最终权重约为 1，说明 IMFs 的最稀疏数量为 5 个，同时...
表1 合成信号
f
f
f f
参数 数量
初始滤波器数量 20
初始滤波器系数数量 20
学习率 1
动量 0.9
最大训练轮数 100
批次大小 64
Parameter Number
Initial filter numbers 20
Initial filter coefficient numbers 20
Learning rate 1
Momentum 0.9
Max epoch 100
Batch size 64 | Parameter | Number |
| :--- | :--- |
| Initial filter numbers | 20 |
| Initial filter coefficient numbers | 20 |
| Learning rate | 1 |
| Momentum | 0.9 |
| Max epoch | 100 |
| Batch size | 64 |
图 4. 输入合成信号
f
f
f f 表2 合成信号
f
f
f f
迭代次数 固有模态函数 Epochs
IMF | Epochs |
| :--- |
| IMF | 25
50
75
1
270.00
259.74
256.60
254.82
2
175.00
182.86
185.19
186.61
3
135.00
121.10
120.15
119.98
4
90.00
80.60
78.16
78.01
5
5.00
14.99
14.99
14.99
"Epochs
IMF" 25 50 75
1 270.00 259.74 256.60 254.82
2 175.00 182.86 185.19 186.61
3 135.00 121.10 120.15 119.98
4 90.00 80.60 78.16 78.01
5 5.00 14.99 14.99 14.99 | Epochs <br> IMF | | 25 | 50 | 75 |
| :--- | ---: | ---: | ---: | ---: |
| 1 | 270.00 | 259.74 | 256.60 | 254.82 |
| 2 | 175.00 | 182.86 | 185.19 | 186.61 |
| 3 | 135.00 | 121.10 | 120.15 | 119.98 |
| 4 | 90.00 | 80.60 | 78.16 | 78.01 |
| 5 | 5.00 | 14.99 | 14.99 | 14.99 |
表3 合成信号
f
f
f f
训练周期 0
25
50
75
IMF
1
20.00
27.46
28.52
28.94
2
20.00
24.25
24.88
29.14
3
20.00
29.35
31.43
33.03
4
20.00
27.31
28.58
29.23
5
20.00
29.59
32.32
34.27
Epochs 0 25 50 75
IMF
1 20.00 27.46 28.52 28.94
2 20.00 24.25 24.88 29.14
3 20.00 29.35 31.43 33.03
4 20.00 27.31 28.58 29.23
5 20.00 29.59 32.32 34.27 | Epochs | 0 | 25 | 50 | 75 |
| :--- | :--- | :--- | :--- | :--- |
| IMF | | | | |
| 1 | 20.00 | 27.46 | 28.52 | 28.94 |
| 2 | 20.00 | 24.25 | 24.88 | 29.14 |
| 3 | 20.00 | 29.35 | 31.43 | 33.03 |
| 4 | 20.00 | 27.31 | 28.58 | 29.23 |
| 5 | 20.00 | 29.59 | 32.32 | 34.27 |
图8. 合成信号
f
f
f f 其他值接近于0。 图 10 对比了不同方法下分解得到的 IMF 分量之间的正交性。每个方块代表通过如下公式[11]计算的正交性结果：
I
O
r
,
k
=
|
⟨
IMF
r
,
IMF
k
⟩
|
‖
IMF
r
‖
2
⋅
‖
IMF
k
‖
2
I
O
r
,
k
=
IMF
r
,
IMF
k
IMF
r
2
⋅
IMF
k
2
IO_(r,k)=(|(:IMF_(r),IMF_(k):)|)/(||IMF_(r)||_(2)*||IMF_(k)||_(2)) I O_{r, k}=\frac{\left|\left\langle\mathrm{IMF}_{r}, \mathrm{IMF}_{k}\right\rangle\right|}{\left\|\mathrm{IMF}_{r}\right\|_{2} \cdot\left\|\mathrm{IMF}_{k}\right\|_{2}} 其中下标
r
r
r r
k
k
k k
I
O
r
,
k
I
O
r
,
k
IO_(r,k) I O_{r, k} 在获取本征模态函数(IMFs)后，相关瞬时频率通过深度神经网络(DNN)识别并输出(图 11(a))。作为对比方法，图 11(b)展示了结合经验模态分解(EMD)与希尔伯特变换的希尔伯特-黄变换(HHT)处理结果。图 11(c)则显示对变分模态分解(VMD)获得的 IMFs 进行希尔伯特变换实现的瞬时频率识别。通过计算各瞬时频率与其解析解之间的相对识别误差(
e
e
e e 图 11. 合成信号
f
f
f f
e
e
e e
a
k
(
k
=
1
,
…
,
5
a
k
(
k
=
1
,
…
,
5
a_(k)(k=1,dots,5 \mathbf{a}_{k}(k=1, \ldots, 5
ω
~
k
(
k
=
1
,
…
,
5
)
ω
~
k
(
k
=
1
,
…
,
5
)
widetilde(omega)_(k)(k=1,dots,5) \widetilde{\boldsymbol{\omega}}_{k}(k=1, \ldots, 5) 3.2 示例2 不同信噪比下频率重叠的合成信号分析 为验证所提方法的抗噪性能，分析了不同信噪比(SNR)下的非线性非平稳信号。相较于前一示例，该信号的单分量具有频谱重叠特性且频带更宽。原始输入信号
s
s
s s
s
1
=
{
4
t
⋅
sin
(
(
50
+
50
t
)
⋅
2
π
t
)
,
0
⩽
t
⩽
0.25
−
4
3
(
t
−
1
)
⋅
sin
(
(
50
+
50
t
)
⋅
2
π
t
)
,
0.25
⩽
t
⩽
1
s
2
=
{
2
t
⋅
sin
(
(
120
+
100
t
)
⋅
2
π
t
)
,
0
⩽
t
⩽
0.5
−
2
(
t
−
1
)
⋅
sin
(
(
120
+
100
t
)
⋅
2
π
t
)
,
0.5
⩽
t
⩽
1
s
3
=
(
2
+
cos
(
8
π
t
)
)
⋅
cos
(
300
π
(
t
+
1
)
2
)
,
s
1
=
4
t
⋅
sin
(
(
50
+
50
t
)
⋅
2
π
t
)
,
0
⩽
t
⩽
0.25
−
4
3
(
t
−
1
)
⋅
sin
(
(
50
+
50
t
)
⋅
2
π
t
)
,
0.25
⩽
t
⩽
1
s
2
=
2
t
⋅
sin
(
(
120
+
100
t
)
⋅
2
π
t
)
,
0
⩽
t
⩽
0.5
−
2
(
t
−
1
)
⋅
sin
(
(
120
+
100
t
)
⋅
2
π
t
)
,
0.5
⩽
t
⩽
1
s
3
=
(
2
+
cos
(
8
π
t
)
)
⋅
cos
300
π
(
t
+
1
)
2
,
{:[s_(1)={[4t*sin((50+50 t)*2pi t)",",0 <= t <= 0.25],[-(4)/(3)(t-1)*sin((50+50 t)*2pi t)",",0.25 <= t <= 1]:}],[s_(2)={[2t*sin((120+100 t)*2pi t)",",0 <= t <= 0.5],[-2(t-1)*sin((120+100 t)*2pi t)",",0.5 <= t <= 1]:}],[s_(3)=(2+cos(8pi t))*cos(300 pi(t+1)^(2))","]:} \begin{aligned}
& s_{1}= \begin{cases}4 t \cdot \sin ((50+50 t) \cdot 2 \pi t), & 0 \leqslant t \leqslant 0.25 \\
-\frac{4}{3}(t-1) \cdot \sin ((50+50 t) \cdot 2 \pi t), & 0.25 \leqslant t \leqslant 1\end{cases} \\
& s_{2}= \begin{cases}2 t \cdot \sin ((120+100 t) \cdot 2 \pi t), & 0 \leqslant t \leqslant 0.5 \\
-2(t-1) \cdot \sin ((120+100 t) \cdot 2 \pi t), & 0.5 \leqslant t \leqslant 1\end{cases} \\
& s_{3}=(2+\cos (8 \pi t)) \cdot \cos \left(300 \pi(t+1)^{2}\right),
\end{aligned} 图 12. 合成信号
s
s
s s
20
dB
,
10
dB
,
5
dB
20
dB
,
10
dB
,
5
dB
20dB,10dB,5dB 20 \mathrm{~dB}, 10 \mathrm{~dB}, 5 \mathrm{~dB}
s
=
s
1
+
s
2
+
s
3
s
=
s
1
+
s
2
+
s
3
s=s_(1)+s_(2)+s_(3) s=s_{1}+s_{2}+s_{3} 三个单分量的瞬时频率均为线性单调递增函数，频率区间分别为[50 150]、[120 320]和[300 600]。原始时间序列信号及其在添加白噪声（
SNR
=
20
dB
,
10
dB
SNR
=
20
dB
,
10
dB
SNR=20dB,10dB \mathrm{SNR}=20 \mathrm{~dB}, 10 \mathrm{~dB} 预训练的 DNN 参数列于表 4。为验证噪声干扰下的识别精度，考虑了
SNR
20
dB
,
10
dB
SNR
20
dB
,
10
dB
SNR20dB,10dB \mathrm{SNR} 20 \mathrm{~dB}, 10 \mathrm{~dB} 如图 13 和图 14 所示，与 VMD（图 17 和图 18）相比，所提方法的单分量在时域和频域上分离效果更佳，而 VMD 在处理高度非平稳信号时失效[35]。然而，EEMD 产生了超出预期的单分量数量，且部分分量并非 IMF（图 16）。
初始中心频率和带宽值在 epoch 为 0 时按表 5 和表 6 设定，各滤波器初始带宽均为 50。最终 IMFs 的中心频率分别为
91.84
Hz
,
217.29
Hz
91.84
Hz
,
217.29
Hz
91.84Hz,217.29Hz 91.84 \mathrm{~Hz}, 217.29 \mathrm{~Hz}
80.22
Hz
,
104.70
Hz
80.22
Hz
,
104.70
Hz
80.22Hz,104.70Hz 80.22 \mathrm{~Hz}, 104.70 \mathrm{~Hz} 图 22 对比了所提方法、EEMD 结合 HT 及 VMD 结合 HT 对三个单分量瞬时频率的识别结果。时变频率识别考虑了两种典型场景：无噪声原始信号及 5dB 信噪比环境。各方法与解析解之间的识别误差
(
e
)
(
e
)
(e) (e) 表4 合成信号的 DNN 相关参数
s
s
s s
参数 数量
初始滤波器数量 10
初始滤波器系数数量 10
学习率 1
动量 0.9
最大训练轮数 1000
批次大小 256
Parameter Number
Initial filter numbers 10
Initial filter coefficient numbers 10
Learning rate 1
Momentum 0.9
Max epoch 1000
Batch size 256 | Parameter | Number |
| :--- | :--- |
| Initial filter numbers | 10 |
| Initial filter coefficient numbers | 10 |
| Learning rate | 1 |
| Momentum | 0.9 |
| Max epoch | 1000 |
| Batch size | 256 |
图 13. 所提方法在不同噪声水平下对合成信号
s
s
s s
SNR
=
20
dB
SNR
=
20
dB
SNR=20dB \mathrm{SNR}=20 \mathrm{~dB}
SNR
=
10
SNR
=
10
SNR=10 \mathrm{SNR}=10
(
d
)
SNR
=
5
dB
(
d
)
SNR
=
5
dB
(d)SNR=5dB (\mathrm{d}) \mathrm{SNR}=5 \mathrm{~dB} 图 14. 所提方法在不同噪声水平下对合成信号
s
s
s s
=
20
dB
=
20
dB
=20dB =20 \mathrm{~dB}
=
10
dB
=
10
dB
=10dB =10 \mathrm{~dB}
=
5
dB
=
5
dB
=5dB =5 \mathrm{~dB} 图 15. EEMD 方法在不同噪声水平下对合成信号
s
s
s s
SNR
=
20
dB
SNR
=
20
dB
SNR=20dB \mathrm{SNR}=20 \mathrm{~dB}
SNR
=
10
dB
SNR
=
10
dB
SNR=10dB \mathrm{SNR}=10 \mathrm{~dB}
SNR
=
5
dB
SNR
=
5
dB
SNR=5dB \mathrm{SNR}=5 \mathrm{~dB} 图 16. EEMD 方法在不同噪声水平下对合成信号
s
s
s s
SNR
=
20
dB
SNR
=
20
dB
SNR=20dB \operatorname{SNR}=20 \mathrm{~dB}
SNR
=
10
dB
SNR
=
10
dB
SNR=10dB \mathrm{SNR}=10 \mathrm{~dB}
SNR
=
5
dB
SNR
=
5
dB
SNR=5dB \mathrm{SNR}=5 \mathrm{~dB} 图 17. 合成信号
s
s
s s
=
20
dB
=
20
dB
=20dB =20 \mathrm{~dB}
SNR
=
10
dB
SNR
=
10
dB
SNR=10dB \mathrm{SNR}=10 \mathrm{~dB}
=
5
dB
=
5
dB
=5dB =5 \mathrm{~dB} 实例中，结果显示无噪声和 5dB 信噪比时，所提方法的数值分别为
1.48
%
1.48
%
1.48% 1.48 \%
4.90
%
4.90
%
4.90% 4.90 \%
40.43
%
40.43
%
40.43% 40.43 \%
53.12
%
53.12
%
53.12% 53.12 \%
30.37
%
30.37
%
30.37% 30.37 \%
36.56
%
36.56
%
36.56% 36.56 \% 本例中，考虑
30
%
30
%
30% 30 \% 图 18. 合成信号
s
s
s s
SNR
=
20
dB
SNR
=
20
dB
SNR=20dB \mathrm{SNR}=20 \mathrm{~dB}
SNR
=
10
dB
SNR
=
10
dB
SNR=10dB \mathrm{SNR}=10 \mathrm{~dB}
SNR
=
5
dB
SNR
=
5
dB
SNR=5dB \mathrm{SNR}=5 \mathrm{~dB} 表5 合成信号
s
s
s s
训练周期 0 250 IMF Epochs 0 250
IMF | Epochs | 0 | 250 |
| :--- | :--- | :--- |
| IMF | | | 500
750
1
300
369.65
439.90
448.69
2
200
213.51
216.64
216.01
3
100
90.94
90.91
91.41
"Epochs 0 250
IMF " 500 750
1 300 369.65 439.90 448.69
2 200 213.51 216.64 216.01
3 100 90.94 90.91 91.41 | Epochs 0 250 <br> IMF | | 500 | 750 | |
| :--- | :--- | :--- | ---: | ---: |
| 1 | 300 | 369.65 | 439.90 | 448.69 |
| 2 | 200 | 213.51 | 216.64 | 216.01 |
| 3 | 100 | 90.94 | 90.91 | 91.41 |
表6 合成信号
s
s
s s
(
Hz
)
(
Hz
)
(Hz) (\mathrm{Hz})
本征模态函数(IMF)的迭代周期 Epochs
IMF | Epochs |
| :--- |
| IMF | 0
250
500
750
1
50
130.61
228.24
297.38
2
50
71.68
85.14
104.53
3
50
62.84
69.93
79.50
"Epochs
IMF" 0 250 500 750
1 50 130.61 228.24 297.38
2 50 71.68 85.14 104.53
3 50 62.84 69.93 79.50 | Epochs <br> IMF | 0 | 250 | 500 | 750 |
| :--- | :--- | ---: | ---: | ---: |
| 1 | 50 | 130.61 | 228.24 | 297.38 |
| 2 | 50 | 71.68 | 85.14 | 104.53 |
| 3 | 50 | 62.84 | 69.93 | 79.50 |
图 19. 合成信号
s
s
s s 3.3 实例3 电缆模型实验示例 如图 23 所示进行缆索实验以验证所提方法的瞬时频率估计精度，并分别采用所提方法、EEMD 和 VMD 对缆索振动信号的 IMF 分量进行分析。测试缆索置于缆索制造厂区坑道内，一侧安装液压千斤顶加载，另一侧设置紧固机构。缆索质量密度为
22
kg
⋅
m
−
1
22
kg
⋅
m
−
1
22kg*m^(-1) 22 \mathrm{~kg} \cdot \mathrm{~m}^{-1} 图20. 合成信号
s
s
s s 图21. 合成信号
s
s
s s 试验中，缆索分三步张拉。首先进行预张拉阶段，施加 1250 kN 的力以测量缆索实际长度，随后释放。第二步，将缆索复位至初始状态后重新加载至 1250 kN，并逐步增加荷载至
1500
kN
,
2000
kN
1500
kN
,
2000
kN
1500kN,2000kN 1500 \mathrm{kN}, 2000 \mathrm{kN} 由 1 号加速度计监测的缆索振动数据波形与傅里叶谱如图 25 所示。所提方法重构信号与原输入信号在时域和频域上差异均较小。采用与图 3 相同的 DNN 结构，相关参数列于表 7。
在此实例中，15 个滤波器的初始带宽为 0.4，中心频率大致按 0.88 的倍数设置。冗余滤波器的系数初始值设为 1。信号长度为 4096 个监测离散点。图 26 和图 27 展示了所提方法、EEMD 和 VMD 分解得到的 IMF 分量，分别对应 9 个、4 个和 9 个 IMF。此处 VMD 的 IMF 数量人为设置为与所提方法相同，因为若预设模态数不准确，该方法会失效。值得注意的是，所提方法能良好分离 IMF 分量，频谱中未出现虚假单分量，而 EEMD 无法实现有效分离（图 27(b)）。对于 VMD 而言，该方法未能捕获频率在
2.5
Hz
,
5.5
Hz
2.5
Hz
,
5.5
Hz
2.5Hz,5.5Hz 2.5 \mathrm{~Hz}, 5.5 \mathrm{~Hz} 提取出的 IMF 优化自适应滤波器如图 28 所示，神经网络模型的损失曲线展示于图 29。经过 1000 轮神经网络训练后，自适应滤波器从 15 个冗余滤波器缩减至 9 个，完全覆盖了 IMFs 的频谱范围，同时相关滤波器系数被稀疏优化，其中 9 个滤波器权重值接近 1（图 30）。 图 31. 缆索振动信号 IMFs 的正交矩阵：(a)本文方法，(b)EEMD 及(c)VMD。 图 31 检验了 IMFs 间的正交性，显示本文方法的
I
O
r
,
k
I
O
r
,
k
IO_(r,k) I O_{r, k} 接下来，在获得电缆信号的自适应稀疏 IMF 后，通过深度神经网络识别相关的瞬时频率。首先对电缆振动信号进行短时傅里叶变换的频谱分析，图 32 揭示了 9 个频率成分，其中能量最大的频率集中在 6 至 8 赫兹之间。选取了五个离散时间点（标记为
t
1
∼
t
5
t
1
∼
t
5
t_(1)∼t_(5) t_{1} \sim t_{5}
6
∼
8
Hz
6
∼
8
Hz
6∼8Hz 6 \sim 8 \mathrm{~Hz}
I
F
=
4
m
L
c
2
(
ω
~
k
(
t
)
2
π
k
)
2
I
F
=
4
m
L
c
2
ω
~
k
(
t
)
2
π
k
2
IF=4mL_(c)^(2)(( widetilde(omega)_(k)(t))/(2pi k))^(2) \mathbf{I F}=4 m L_{c}^{2}\left(\frac{\widetilde{\boldsymbol{\omega}}_{k}(t)}{2 \pi k}\right)^{2} 式中
L
c
L
c
L_(c) L_{c}
m
m
m m
ω
~
k
(
t
)
ω
~
k
(
t
)
widetilde(omega)_(k)(t) \widetilde{\boldsymbol{\omega}}_{k}(t)
k
k
k k 图32. 电缆振动信号的频谱图。右图显示了五个时间点的绝对振幅。 图 33. 由 1 号、3 号和 5 号加速度计监测信号识别的时变电缆力及其识别误差(e)。 图 33 对比了加速度计 1 号、3 号和 5 号的数据。识别出的索力按 IMF 分量顺序排列，其中最低固有频率对应最大阶数。各索力估计值与实测真实索力之间的相对误差通过以下公式计算：
e
=
‖
I
F
k
−
F
‖
2
‖
F
‖
2
×
100
%
e
=
I
F
k
−
F
2
‖
F
‖
2
×
100
%
e=(||IF_(k)-F||_(2))/(||F||_(2))xx100% e=\frac{\left\|\mathbf{I} \mathbf{F}_{k}-\mathbf{F}\right\|_{2}}{\|\mathbf{F}\|_{2}} \times 100 \% 其中
I
F
k
I
F
k
IF_(k) \mathbf{I F}_{k}
k
k
k k
F
F
F \mathbf{F}
2.80
%
,
1.18
%
,
1.16
%
,
0.98
%
,
1.93
%
,
1.52
%
,
1.70
%
,
2.43
%
2.80
%
,
1.18
%
,
1.16
%
,
0.98
%
,
1.93
%
,
1.52
%
,
1.70
%
,
2.43
%
2.80%,1.18%,1.16%,0.98%,1.93%,1.52%,1.70%,2.43% 2.80 \%, 1.18 \%, 1.16 \%, 0.98 \%, 1.93 \%, 1.52 \%, 1.70 \%, 2.43 \%
2.18
%
2.18
%
2.18% 2.18 \% 如图所示，
I
F
3
I
F
3
IF_(3) \mathbf{I F}_{3}
I
F
4
I
F
4
IF_(4) \mathbf{I F}_{4}
6
∼
8
Hz
6
∼
8
Hz
6∼8Hz 6 \sim 8 \mathrm{~Hz} 为了说明 IMF 能量对瞬时频率估计精度的影响，图 34 展示了不同时间点 IMF 的振幅及相对误差。识别出的索力相对误差与频谱图中五个不同时间点电缆振动信号的绝对振幅相关联。结果显示，对于
I
F
3
I
F
3
IF_(3) \mathbf{I F}_{3}
I
F
4
I
F
4
IF_(4) \mathbf{I F}_{4} 4. 结论 本文提出了一种基于深度学习的自适应模态分解与瞬时频率估计方法，用于结构健康监测中的非线性和非平稳信号分析。为实现无模态混叠和虚假分量的精确本征模态函数提取以及高精度的瞬时频率估计，构建了一种新型深度神经网络框架。在该网络中，采用自适应滤波器进行单分量频谱提取，并将希尔伯特变换嵌入网络层中以实现 图34. 索缆振动信号频谱图幅值与相对误差的关系。 瞬时频率估计。通过考虑频谱重叠和不同噪声干扰水平的数值信号，以及具有变化张力工况的长索模型实验振动数据，验证了所提方法的有效性。 该方法能够自适应地分解具有紧密间隔单分量的合成信号。此外，在示例中实现了
1.48
%
1.48
%
1.48% 1.48 \%
4.90
%
4.90
%
4.90% 4.90 \%
53.12
%
53.12
%
53.12% 53.12 \%
36.56
%
36.56
%
36.56% 36.56 \%
10
−
5
10
−
5
10^(-5) 10^{-5}
10
−
6
10
−
6
10^(-6) 10^{-6} 从缆索振动数据的结果来看，时变缆力估计值与实测值之间的识别误差小于
2.80
%
2.80
%
2.80% 2.80 \%
0.98
%
0.98
%
0.98% 0.98 \% 基金资助 本研究得到国家自然科学基金资助[项目编号：51978216、52192664]。
CRediT 作者贡献声明 郭一兵：概念化、方法论、软件、验证、形式分析、初稿撰写。鲍跃全：概念化、方法论、文稿审阅与编辑、监督指导。李惠：概念化、方法论、文稿审阅与编辑、监督指导。张玉峰：资源提供。
利益冲突声明 作者声明，对于可能影响本研究报告的已知竞争性经济利益或个人关系不存在。
数据可用性 作者无权限共享数据。
参考文献 [1] 鲍勇，陈志，魏松，徐毅，唐震，李华，结构健康监测中的数据科学与工程研究现状，《工程》5 (2019)
234
−
242
234
−
242
234-242 234-242 [2] 李华楠，任力，贾志国，易天昊，李大胜，大型复杂土木基础设施结构健康监测研究进展，《土木结构健康监测杂志》6 (2016) 3-16。 [3] 欧进萍, 李惠, 中国大陆结构健康监测研究综述与未来趋势, 结构健康监测 9 (2010) 219-231. [4] 侯锐, 夏勇, 土木工程结构基于振动的损伤识别新进展综述: 2010-2019, 声与振动杂志 491 (2021) 115741. [5] 杨勇, 李素娟, S. Nagarajaiah, 李惠, 周培, 基于复杂度追踪的加速度实时输出唯一识别时变索力, 结构工程杂志 142 (2016) 04015083. [6] B. Bhowmik, T. Tripura, B. Hazra, V. Pakrashi, 基于递归典型相关分析的鲁棒线性和非线性结构损伤检测, 机械系统与信号处理 136 (2020) 106499. [7] S. Panda, T. Tripura, B. Hazra, 基于近实时包络提取的地震激励结构在线损伤检测, 结构健康监测. 21 (2022) 298-319. [8] L. Li, H. Zhou, H. Liu, C. Zhang, J. Liu, 一种耦合经验模态分解与长短期记忆网络的混合方法用于预测 SHM 系统缺失测量信号数据, 结构健康监测. 20 (2021) 1778-1793. [9] L. Lian, W. Ren, S. Wang, 分数阶傅里叶变换：时频表征与结构瞬时频率识别, 机械系统与信号处理. 178 (2022) 109305. [10] D. Goyal, B. Pabla, 结构健康监测中的振动监测方法与信号处理技术综述, 计算机方法工程档案. 23 (2016) 585-594. [11] 黄锷、沈正、龙慎仪等，《基于经验模态分解与希尔伯特谱的非线性非平稳时间序列分析方法》，《皇家学会会刊 A 辑-数学物理》454 卷（1998 年）903-995 页。 [12] 丹尼斯·盖博，《通信理论》，《电气工程师学会杂志》93卷（1946年）42-57页。 [13] 让·阿莱克西·维尔，《解析信号理论及其应用》，《电缆与传输》2卷（1948年）61-74页。 [14] 斯特凡·马莱，《小波基》，载于《信号处理的小波之旅：稀疏方法》（第3版），爱思唯尔出版社，伯灵顿，2009年，263-376页。 [15] W.J. Fitzgerald，《傅里叶分析》，载于：《电信工程师参考手册》，Butterworth-Heinemann 出版社，牛津，1993 年，第 2-17 页。 [16] H. Aied, A. González, D. Cantero，《基于集合经验模态分解的桥梁加速度响应中移动荷载引起的刚度突变识别》，《机械系统与信号处理》66 卷（2016 年）314-338 页。 [17] M. Lazhari, AJJoS Sadhu，《采用自适应经验模态分解法的结构分散模态识别》，《声音与振动》447 卷（2019 年）20-41 页。 [18] P. Singh, M. Keyvanlou, A. Sadhu，《一种基于有限传感器的结构状态评估改进时变经验模态分解方法》，《工程结构》232 卷（2021 年）111882 页。 [19] C. Delprete, E. Brusa, C. Rosso, F. Bruzzone, 基于正交经验模态分解的轴承健康监测, 《振动与冲击》2020 年 (2020). 
