【翻译】具有时间结构的特定信号的鲁棒提取（上）

具有时间结构的特定信号的鲁棒提取

《Robust extraction of specific signals with temporal structure 》

1、摘要

基于二阶统计量的盲或半盲源提取算法具有计算量小、处理速度快的特点，有发展的趋势。Barros和Cichocki做了重要和初步的工作，他们提出了一种基于时间延迟的提取算法。该算法简单、快速，但性能不理想。本文扩展了他们的工作，提出了一种基于多个延迟协方差矩阵特征值分解的鲁棒算法。理论分析和计算机仿真验证了该算法速度快、性能好。

盲源提取(BSE)[6]是一种与盲源分离(BSS)密切相关的强大技术[6,3,5]。BSE的基本任务是对观测数据中线性组合的部分源信号进行估计。与BSS相比，BSE具有许多优势，在生物医学信号处理[12,13,1]和语音处理[6]等各个领域得到了广泛关注

在BSE中，我们观察到一个$n$维随机信号向量$x$，它被认为是$n$维零均值单位方差源向量$s$的线性变换，即$x=As$，其中$A$是一个未知的混合矩阵。BSE的目标是找到一个矢量$w$，使$y=w^{T}x=w^{T}As$是一个估计的源信号，最大为标量。为了处理病态情况，并使算法更简单、更快，通常使用一种称为预白化的线性变换将观测信号$x$变换为$\tilde{x}=Vx$，令 E { x ~ x ~ T } = I E\{\tilde{x}\tilde{x}^{T}\}=I E{x~x~T}=I，其中$V$为预白化矩阵。为方便起见，下面我们假设$x$是预白化后的观测信号

许多源提取算法[6]利用稀疏性[14]和高阶统计量[12]等先验信息来提取特定的信号作为第一输出。但利用稀疏性或高阶统计量的算法往往具有较高的计算量。因此，利用信号时间结构的先验知识开发基于二阶统计量的提取算法是一种趋势。Barros和Cichocki[1]做了一项重要的工作。他们提出了一种快速而简单的算法，该算法要求精确估计与期望源信号的某些时间结构相对应的最佳时间延迟。遗憾的是，该算法的性能严重依赖于时间延迟;较小的时延估计误差往往导致性能较差。

在Barros和Cichocki的基础上，本文提出了一种基于特征值分解的鲁棒提取算法。理论分析和计算机仿真验证了该算法速度快，且只要误差不大，其性能不受时延估计误差的影响

2、提出算法

假设期望的源信号$s_i$是时间相关的，对于特定的时间延迟${\tau }^{\ast }$满足以下关系
{ E { s i ( k ) s i ( k − τ ∗ ) } > 0 E { s i ( k ) s j ( k − τ ∗ ) } = 0 E { s j ( k ) s l ( k − τ ∗ ) } = 0 , ∀ j ≠ i , l ≠ i ( 1 ) \left\{\begin{matrix} E\{s_{i}(k)s_{i}(k-\tau^{*})\}>0 \\ E\{s_{i}(k)s_{j}(k-\tau^{*})\}=0 \\ E\{s_{j}(k)s_{l}(k-\tau^{*})\}=0, \forall j\neq i,l\neq i \end{matrix}\right. \\(1) ⎩⎨⎧​E{si​(k)si​(k−τ∗)}>0E{si​(k)sj​(k−τ∗)}=0E{sj​(k)sl​(k−τ∗)}=0,∀j​=i,l​=i​(1)
$k$是时间下标，${\tau }^{\ast }$是一个整数延迟，可以是正的或负的。在不失一般性的情况下，下面假定时滞为正。注意(1)与[1]中的假设非常相似。在约束$||w||=1$下，最大化
J ( w ) = E { y ( k ) y ( k − τ ∗ ) } = w T E { x ( k ) x ( k − τ ∗ ) } w ( 2 ) J(w)=E\{y(k)y(k-\tau^{*})\}=w^{T}E\{x(k)x(k-\tau^{*})\}w\\(2) J(w)=E{y(k)y(k−τ∗)}=wTE{x(k)x(k−τ∗)}w(2)
导出所需的源信号。这里，$y(k)=w^{T}x(k)$是输出信号。当$J(w)$达到一个最大值时，$y(k)$估计所需的信号$s_i$达到一个标量。这个建议的原因是，对于期望的源信号，这种自相关应该有一个很高的值，而对于其他源信号，这个值应该非常小。为了方便起见，我们假定所需的源信号是周期性的。但这并不意味着本文提出的算法仅限于周期信号的提取。事实上，其他非周期信号也可以提取出来，只要它们满足假设(1)，并且可以估计出相应的时延。

根据Barros和Cichocki[1]的思想，根据假设(1)，我们可以很容易地得到Barros算法:
{ w + = E { x ( k ) x ( k − τ ∗ ) T } w w = w + / ∣ ∣ w + ∣ ∣ ( 3 ) \left\{\begin{matrix} w^{+}=E\{x(k)x(k-\tau^{*})^{T}\}w \\ w=w^{+}/||w^{+}|| \end{matrix}\right. \\(3) {w+=E{x(k)x(k−τ∗)T}ww=w+/∣∣w+∣∣​(3)
如果想要的信号是周期的，基本周期为${\tau }_0$，那么${\tau }^{\ast }$可以被设为${\tau }^{\ast }=r{\tau }_0$，其中r为非零整数。

现在再次考虑目标函数(2)。我们有
J ( w ) = 1 2 J ( w ) + 1 2 J ( w ) T = 1 2 w T E { x ( k ) x ( k − τ ∗ ) T } w + 1 2 w T E { x ( k − τ ∗ ) x ( k ) T } w = 1 2 w T ( R x ( τ ∗ ) + R x ( τ ∗ ) T ) w ( 4 ) J(w)=\frac{1}{2}J(w)+\frac{1}{2}J(w)^{T}\\ =\frac{1}{2}w^{T}E\{x(k)x(k-\tau^{*})^{T}\}w+\frac{1}{2}w^{T}E\{x(k-\tau^{*})x(k)^{T}\}w\\ =\frac{1}{2}w^{T}(R_{x}(\tau^{*})+R_{x}(\tau^{*})^{T})w \\(4) J(w)=21​J(w)+21​J(w)T=21​wTE{x(k)x(k−τ∗)T}w+21​wTE{x(k−τ∗)x(k)T}w=21​wT(Rx​(τ∗)+Rx​(τ∗)T)w(4)
推在约束$||w||=1$下使(2)最大化，相当于求实对称矩阵$R=R_x({\tau }^{\ast })+R_x({\tau }^{\ast }{)}^T$最大特征值对应的特征向量。其中 R x ( τ ∗ ) = E { x ( k ) x ( k − τ ∗ ) T } R_{x}(\tau^{*})=E\{x(k)x(k-\tau^{*})^{T}\} Rx​(τ∗)=E{x(k)x(k−τ∗)T}。因此，我们有以下算法:
{ R x ( τ ∗ ) = E { x ( k ) x ( k − τ ∗ ) T } w = E I G ( R x ( τ ∗ ) + R x ( τ ∗ ) T ) ( 5 ) \left\{\begin{matrix} R_{x}(\tau^{*})=E\{x(k)x(k-\tau^{*})^{T}\} \\ w=EIG(R_{x}(\tau^{*})+R_{x}(\tau^{*})^{T}) \end{matrix}\right. \\(5) {Rx​(τ∗)=E{x(k)x(k−τ∗)T}w=EIG(Rx​(τ∗)+Rx​(τ∗)T)​(5)
其中$EIG(T)$是计算实对称矩阵$T$的最大特征值对应的归一化特征向量的算子。与Barros算法(3)相比，算法(5)由于有效的特征值分解技术[8]，速度更快。

但是，应该考虑一些实际问题。一个重要的问题是有限样本的影响。虽然信号是相互不相关的，实际上有限样本上计算的源信号的互相关联值通常是非零的，这是因为期望算子 E { z ( k ) } E\{z(k)\} E{z(k)}被数学平均${\sum }_{k=1}^{N}z(k)/N$ ($N$是数据长度)所取代。即使可用的样本数量非常大，在实践中它们往往被分成连续的块，因为批处理算法的处理类型不同。所以在这种情况下互相关值仍然是非零的。另一个关键问题是时延估计误差。在许多情况下，这些错误是无法避免的。我们将在后面看到，可用样本的限制和估计误差极大地影响了算法(3)和(5)的性能。

因此，我们将算法(5)修改为(原因在第3节)
{ R x ( τ ∗ ) = E { x ( k ) x ( k − τ ∗ ) T } w = E I G ( ∑ i = 1 P ( R x ( i τ ∗ ) ) + ( R x ( i τ ∗ ) ) T ) ( 6 ) \left\{\begin{matrix} R_{x}(\tau^{*})=E\{x(k)x(k-\tau^{*})^{T}\} \\ w=EIG(\sum_{i=1}^{P}(R_{x}(i\tau^{*}))+(R_{x}(i\tau^{*}))^{T}) \end{matrix}\right. \\(6) {Rx​(τ∗)=E{x(k)x(k−τ∗)T}w=EIG(∑i=1P​(Rx​(iτ∗))+(Rx​(iτ∗))T)​(6)
最大化新的目标函数
J ( w ) = ∑ i = 1 P E { y ( k ) y ( k − i τ ∗ ) } = w T ( ∑ i = 1 E { x ( k ) x ( k − i τ ∗ ) } w ( 7 ) J(w)=\sum_{i=1}^{P}E\{y(k)y(k-i\tau^{*})\}=w^{T}(\sum_{i=1}E\{x(k)x(k-i\tau^{*})\}w \\(7) J(w)=i=1∑P​E{y(k)y(k−iτ∗)}=wT(i=1∑​E{x(k)x(k−iτ∗)}w(7)
约束$||w||=1$，其中$P$为正整数，${\tau }^{\ast }$是期望源信号的基本周期。