时频分析的时频谱分解

截至目前为止发展的几种时频分析方法中,科恩系列分布是最为有力的转换之一，但是其计算复杂度高于短时距傅立叶变换和小波分析，因此在应用上有所受限。 透过将时频分布的核分解为不同时频谱的线性组合，每一份时频谱都是由简单的短时距傅立叶变换所计算而得，可以协助我们有效降低分析信号时频特性的计算量。一般来说当我们使用哈尔小波转换来分解时频分部的核的时候，我们可以将STFT的计算量降至只需少量的短时距傅立叶变换即可。

基本概念

根据Cunningham和Williams最早提出的想法，首先我们将时频分布的核分解为不同短时距傅立叶变换的结果

$C(n\,,k)=\sum _{i=1}^{M}\lambda _{i}|STFT_{i}(n,k)|^{2}$

然后将短时距傅立叶变换的结果以时频谱的窗函数$w_{i}(n)$和信号$x(n)$加以分解而得

$STFT_{i}(n,k)=\sum _{m=0}^{N-1}x(m-n)w_{i}(m)e^{-j2\pi \,mk/N}$

则我们能顺利将核分解为不同时频谱的窗函数$w_{i}(t)$，但是尽管不同的时频谱之间彼此正交，不同的时频谱之间并没有重要性的顺序关联性让我们能依据逼近精细度的不同而选择要保留那些项，因此还是必须计算完整的分解以决定所有的时频谱。

时频谱分解

一般来说，一个离散双线性的时频分布计算由此时频分布的核$\psi$所决定，并可表达为内积的形式

$X(n,\omega \,,\psi )=\sum _{n_{1}}\sum _{n_{2}}x(n+n_{1})e^{-j\omega \,(n+n_{1})}\psi (-{\frac {n_{1}+n_{2}}{2}},n_{1}-n_{2})[x(n+n_{2})e^{-j\omega \,(n+n_{2})}]^{*}=<,{\tilde {\psi }}{\bar {S}}_{-n}{\bar {M}}_{-\omega }x,{\bar {S}}_{-n}{\bar {M}}_{-\omega }x>$

其中${\bar {S}}_{-n}$和${\bar {M}}_{-\omega }$分别代表时间位移和频率位移的算子，而${\tilde {\psi }}$则是代表核$\psi$的作用算子,透过将算子${\tilde {\psi }}$以特征值分解，我们可以将核分解为线性成分，因而时频分布本身被拆解为时频谱的和。具体的来说和可以被分解为

$\sum _{k}\sum _{l}\lambda _{k,l}w_{k}(n_{1})w_{l}(n_{2})$

其中$w_{k}$和$w_{l}$是窗函数，而$\lambda _{k,l}$是特征值，因此我们重新改写时频分析的计算为

$X(n,\omega \,,\psi )=\sum _{n_{1}}x(n+n_{1})e^{-j\omega \,(n+n_{1})}\sum _{n_{2}}[x(n+n_{2})e^{-j\omega \,(n+n_{2})}]^{*}\sum _{k}\sum _{l}\lambda _{k,l}w_{k}(n_{1})w_{l}(n_{2})$

则$\sum _{n_{1}}w_{k}(n_{1})x(n+n1)e^{-j\omega \,(n+n_{1})}$可以被看作信号$x(n)$的短时距傅立叶变换$STFT_{k}(n,\omega )$,因而若令$SP_{k,l}(n,\omega )=STFT_{k}(n,\omega )STFT_{l}^{*}(n,\omega )$

$X(n,\omega \,,\psi )=\sum _{k}\sum _{l}\lambda _{k,l}SP_{k,l}(n,\omega )=\sum _{k}\sum _{l}\lambda _{k,l}STFT_{k}(n,\omega )STFT_{l}^{*}(n,\omega )$

则时频分布拆解为不同的时频谱的和，系数$\lambda _{k,l}$代表权重，并可被简单的由特征值分解计算而得

降低STFT的计算量

$8\times 8$的离散维格纳分布的核矩阵为

${\begin{bmatrix}0&0&0&0&0&0&0&1\\0&0&0&0&0&0&1&0\\0&0&0&0&0&1&0&0\\0&0&0&0&1&0&0&0\\0&0&0&1&0&0&0&0\\0&0&1&0&0&0&0&0\\0&1&0&0&0&0&0&0\\1&0&0&0&0&0&0&0\end{bmatrix}}$

选取8点的哈尔小波作为基底，可计算出特征值的矩阵

$(\lambda _{k,l})={\begin{bmatrix}1&0&0&0&0&0&0&0\\0&-1&0&0&0&0&0&0\\0&0&0&-1&0&0&0&0\\0&0&-1&0&0&0&0&0\\0&0&0&0&0&0&0&-1\\0&0&0&0&0&0&-1&0\\0&0&0&0&0&-1&0&0\\0&0&0&0&-1&0&0&0\end{bmatrix}}$

可见一开始需要计算$8\times 8=64$次的时频谱,透过特征值分解只有少数的非零系数，因此可大大减低计算量。

参考资料

• G.S. Cunningham, W.J. Williams, “Kernel decomposition of time-frequency distributions”, IEEE Trans. Signal Process. 42 (6) (1994) 1425-1442.
• W.J. Williams, T.-H. Sang, J.C. O’Neill, E.J. Zalubas, “Wavelet windowed time-frequency distribution decompositions”, in: Proc. SPIE: Advanced Signal Processing Algorithms, Architectures, and Implementations VII, vol. 3162, Soc. of Photo-Optical Instrumentation Engineers, 1997, pp. 149-160.
• Boualem Boashash, “Time frequency signal analysis and processing”, 2016