自适应滤波算法

发布时间：March 14, 2023, 3:08 a.m.编辑：李佳生阅读（1589）

一、Wiener FIltering

二、LMS类算法

三、Kalman-Filter

四、数值优化(数学规划)思想

一、Wiener FIltering

1、原理

    维纳解是一种描述线性时不变系统(LTI System)输入和输出的理想解，目标是使系统输出的估计值与实际值具有最小的线性误差。下面简述一下推导过程:

    假设一个LTI系统存在一组输入[x(0)....x(n)]，经过系统响应后输出为[y(0)....y(n)]。这时我们希望通过一组滤波器对输出进行估计预测，则有ye(n)= Σh(m)x(n-m)，此时误差为e(n)=y(n)-ye(n)，理想滤波器解应该使整体误差最小，若写成线性最小均方误差(LMMSE)的形势则有Σe(n)*e(n)→min。

    由于我们假设系统为线性系统，可以理解为完美维纳解描述了两个信号之间线性相关的部分，等价的说法是估计的残留误差应改是线性滤波器无法描述的信号成分，也就是说误差与输入应该尽量线性不想关、互相正交，因此维纳解的求解目标也可以写为Σe(n)*x(n)→min。

    上述两种思路得到的解是一样的，最终可推导出下述由x自相关矩阵、x和y的互相关矩阵以及理想滤波器频响h组成的等式(Wiener-Hopf方程):

    其中:

    将等式中的互相关和自相关写成矩阵形式并移项即得到经典的维纳解方程:

    如果把线性系统理解为带噪信号到期望信号的线性变换，则输入可理解为带噪信号，输出可理解为期望信号，当我们知道了带噪信号和期望信号，就可以得到相应的维纳解，从带噪信号滤波得到期望信号。

    上图中d为期望信号，v为噪声，x为输入的带噪信号。求解目标为得到理想的滤波器h，能够从x中分离出d信号，若假设噪声与期望信号不想关，则根据维纳解可得到另外一个经典公式(常见的信号降噪原理):

    上述过程为维纳解的时域描述，如果转换为频域，则优化目标变为频响误差最小，推导过程与最终结果与时域类似，只需把上述公式中的互相关和自相关矩阵则理解为功率谱和互功率谱:

    对于噪声与信号不想关的带噪信号滤波，维纳解公式变为:

2、应用

    从上述原理可知，如果我们知道了线性系统的输入(带噪信号)和输出(期望信号)，就可以计算出理想的滤波器解来描述这个线性系统输入与输出的关系，从而利用输入信号预测输出信号(期望信号)。大部分同学到这里会问，这个求解有点像达文西的手电筒，如果我们事先知道了期望信号，为什么还要去求解滤波器来滤波，这里可以通过两个情况来理解:

    1、预训练: 我们可以在滤波器实际工作之前，对其系数进行训练。这时我们就可以人为的给出两种信号进行维纳解求解。因为系统我们认为是线性时不变的，因此我们将训练好的系数用到系统中，即可得到实际的期望信号。

    2、在线自适应求解: 如果我们能在线测量输入信号和期望信号，那也可以进行维纳求解。比如回声消除中的播放数字信号已知，播放后mic采集的信号也已知，就可以进行求解更新系数；另外在语音降噪中，在非语音时间点，我们可以知道噪声的特征，也可以进行维纳求解，来达到降噪的目的(从带噪语音信号中得到纯语音信号)。

    虽然维纳解是一种很好的滤波器求解方法，但是也需要注意实际情况是否符合使用条件和一些需要注意的点

    1、线性：在原理中有描述，维纳解为线性时不变系统的最优解，对于非线性部分，是当做随机误差来看待的，因此如果系统实际工作时，如果进入严重的非线性区域，那维纳解的偏差就会变大，导致预测偏离实际值。

    2、时不变系统：比较好理解，因为维纳解是根据系统某个时间段的输入、输出统计数据得到的，因此描述的也是这段时间的系统特性，只有系统特性不随着时间改变，或者说在维纳解作用期间是稳态的，这个根据过去数据得到的描述才是准确的，否则这个解就失去了意义。

    3、稳态信号过程：与时不变的理解类似，维纳解是某个时间段内输入和输出数据计算出的期望最优解，如果信号非稳态，各种统计期望值都会发生变化，最优解也就失效了。

    4、解的非因果性：上述维纳解形式只描述了输入和输出的理想线性关系，但是不能保证解的因果性，单纯看Wiener-Hopf方程，其中计算的数据是可以包含未来时刻数据的。滤波器解决的问题可分为以下几种：

    (1) filtering：根据当前和过去数据估计出当前数据的真实值，需要保证因果性。

    (2) prediciton：根据当前和过去数据估计出未来数据的真实值，需要保证因果性。

    (3) smoothing：根据当前和过去数据估计出过去数据的真实值，不需要保证因果性。这也是为什么非因果滤波器虽然物理上无法实现，但是还是有研究意义的原因，可以在软件中对已经存保存的数据进行处理，非因果数据也可以工作。

    设计因果维纳滤波器可以通过白化输入信号来实现，对于白噪声来说，因果维纳解和非因果维纳解中的因果部分完全一致。因此可以通过白化滤波器将输入转化为白噪声，然后计算白化输入与输出的非因果维纳解，去除非因果部分后得到白化信号与输出的因果维纳解，再将白化滤波器考虑在内得到最终的传递函数。参考链接

    5、滤波器种类和阶数：维纳解是一组理想的滤波器系数或频响描述，而实际使用滤波器阶数或种类(FIR、IIR)可能有限，因此只能近似的描述出理想的维纳解，其逼近的方法、描述特征的取舍和一些工程的实现则需要根据具体的使用场景来考虑。

    如果不符合上述条件，求解可能是不准确的，或者需要对计算过程进行修正、修饰。

二、LMS类算法

    LMS类算法的种类比较多，这里只讲下最常见的类型，基本所有LMS类算法的原理都与以下几种算法的原理类似，只是加上了一些变形和修饰。下述推导过程和描述只以时域为例进行介绍，所有公式也可以推广到频域，进行频域的LMS算法。

1、LMS

    LMS是比较经典的自适应滤波算法了，这里不做过多介绍，基本原理就是根据输入x(n)、期望信号d(n)、预测误差e(n)和当前的滤波器系数w(n)，通过梯度来迭代优化调整滤波器系数得到w(n+1)，最终使滤波器预测的信号输出与期望信号的误差方差J(w)最小。

    其中步长μ是一个重要的经验参数，与具体场景有关系，表示滤波器参数根据梯度向量调整幅度的大小关系，过大会导致自适应过程失稳波动或者发散，而过小又会收敛很慢无法跟踪快速变化的信号，也更容易进入局部最优解，导致最终收敛结果效果很差，因此尝试出合适的步长大小对于算法的效果与性能很重要。

2、NLMS

    如上述内容描述，LMS算法调整的时候只是确定了每个参数的调整方向(梯度)，而每个参数的调整幅度是可变的，理论上只要按照梯度方向调整就能使误差减小，参数的调整幅度存在无数种组合。实际使用中，通常使用固定步长μ来确定调整幅度，但这样存在两个问题：当x(n)很大时，梯度的噪声也会被放大很多，导致最终调整结果产生偏离，甚至失稳。

    因此引入了NLMS的概念，通过归一化、最小化干扰原则约束来解决上述问题。最小化干扰原则即假设原本的系统是稳定、次优的，这时对参数做出的调整肯定是越小越有可能保持系统参数原本的稳定性和特性。因此基于最小化干扰原则，在LMS的目标函数上可以增加一个约束，使参数调整的幅度尽可能的小：

    假设参数某次迭代的变化量为，

    若希望参数变化的模长和最小，则优化目标可以写为，

    转变为拉格朗日算子表达形式，则目标函数可以写为，

    可以得到梯度向量，

    因为J(n)可认为是二次凸函数，因此g=0时取得最小极值，可得

    联系上面的公式即可求出拉格朗日算子的解:

    将上述结果带入前面公式，并加上步长μ以及防止除数过小的修正因子δ(导致除0或极大值的数值计算问题)后，最终得到NLMS公式，

    因此可以看出，NLMS不仅避免了参数调整幅度过大，并且实现了步长自适应调整，因此也可把 NLMS 看作自适应时变步长参数的 LMS。无论对于不相关数据还是相关数据，NLMS 一般都要比标准 LMS 呈现出更快的收敛速度。上面说的LMS/NLMS都是在时域逐点进行计算的，每个点都要进行一次梯度计算和系数刷新，开销比较大，于是有了以下两种变种LMS：

    BLMS：block LMS，将数据分块，攒够一个block的点数再进行计算(可以有overlap)，这样就大大减少了计算量，对应频域即FDAF(Frequency-Domain Block Least Mean Square Adaptive Filter)。

    SBLMS：sub-band LMS，对滤波器进行子带拆分，将原滤波器拆分成若干个子带滤波器，可以降低计算复杂度，提高计算的并行性。并且每个子带内的信号比较均匀，避免了全带信号特征值较大容易发散的缺点，使得各个子带中进行自适应滤波的数据信号具有更低的相关性(接近白噪声)，从而加快了自适应滤波器的收敛速度和收敛效果。

    子带的缺点是子带滤波器设计的好坏会极大影响算法最终收敛的效果，由于子带之间存在混叠分量，会导致收敛后的稳态误差增大。虽然采用正交镜像滤波器组，可以让子带系统的混叠部分相互抵消掉，但在现实中却无法实现，目前解决子带混叠问题的方法有：

        (1) 利用采样滤波器组技术，但会增加算法复杂度。

        (2) 相邻子带间留安全频带，但会引入空白频带，降低信号质量。

        (3) 重叠子滤波器补偿法，缺点是交叉项会增加运算量、降低收敛速度。

        A_Real-Time_Wideband_Subband_LMS_Algorithm_for_Full-Duplex_Communications.pdf

3、FxLMS/FxNLMS

    Filtered-x LMS/NLMS假设系统中存在路径h，则最终的输出为滤波器w和路径H共同作用的结果：y=x*w*h。因此输入和输出误差之间多了个路径的滤波作用，自适应调整滤波器参数的时候也需要将这个作用考虑在内。可以通过将输入x先进行路径h的滤波处理得到x'，这时x'与e的关系就与传统LMS/NLMS算法一致了，可以通过相应的步骤进行计算。

    这时存在的问题是路径h的实际模型我们是不知道的。一般情况下我们可以预先离线进行次级路径的测量与建模，得到一个估计模型h'，在FxLMS/FxNLMS工作过程中，使用这个估计的模型h'来对输入进行处理得到近似的x'进行后续处理。以FxNLMS为例，参数更新的公式变为：

    虽然路径模型的估计值h'与实际值h一般是有误差的，但通过调整步长μ为合适值，与这个误差满足一定的关系，LMS就可以自动修正这个误差，使最终的输出误差e达到预期的最小值，相关证明比较多，这里不赘述了。

An analysis of the convergence condition for narrowband FxLMS.pdf

4、RLS&APA

    递归最小均方RLS(Recuisive Least Mean Square)方法利用系统过去和现在的所有数据，逐点进行递归计算，调整参数，收敛期望是针对所有点误差最小(也可以加入遗忘系数，逐步减弱很久远采样点的影响)，因此收敛比较快，但计算量很大。

    推导过程和公式很多地方已经写的很明白了，可参考文章：

        RLS.pdf

        递归最小二乘法.pdf

    对应来说，LMS只针对当前点的数据进行梯度下降搜索、调整参数，收敛速度和结果很不稳定，但计算量很小。因此在很多情况我们需要在RLS和LMS特性之间进行一个trade off，于是就出现了仿射投影算法APA，通过保留N组输入和输出数据来计算梯度，因此计算量和收敛特性介于RLS和LMS之间。A Fast APA.pdf

    用一段代码来说明：

三、Kalman-Filter

1、Bayes滤波

    贝叶斯滤波通过贝叶斯概率理论手段进行变量估计，是一种最大后验MAP(Maximum A Posteriori)的滤波方法(注意定义与最大似然估计MLE的区别)，一些经典算法比如卡尔曼滤波、粒子滤波、信息滤波等都是以贝叶斯滤波为基础实现的算法，只是对概率模型种类做了不同的假设而已。贝叶斯滤波的使用条件和基本假设为：

    (1) 系统为环境、模型时不变的稳态系统；

    (2) 观测噪声、模型噪声相互独立；

    (3) 一阶Markov性假设：系统t时刻的状态只由系统在(t-1)时刻的状态和t时刻的控制输入量决定(齐次假设)，并且t时刻的观测值仅与系统在t时刻的状态相关(观测独立假设)。

    贝叶斯滤波的基本思想和计算过程：根据之前所有时刻的数据得出状态变量的概率分布，再利用这个概率分布、前一时刻的状态、当前的控制变量(一阶Markov性假设)预测当前状态变量的分布，再根据实际测量值对分布进行递归修正，修正后的分布均值即我们想要得到的最优估计值。整体计算流程为：

    其中主要计算在3、4步，利用到了递归贝叶斯估计：

    对应贝叶斯滤波的详细过程公式为：R为观测误差、Q为状态误差

2、KF卡尔曼滤波原理

    卡尔曼滤波KF(Kalman Filter)是个很经典的适用于线性系统的自适应滤波算法，大概原理是：在一个已知模型的线性系统中，存在过程和观测噪声，因此如果我们想确认当前状态变量的真实值，无论通过模型来计算预测还是实际观察都是存在误差的。所以最好的办法就是结合预测值和估计值的特征，评估出一个最优的真实值，这就是卡尔曼滤波的基本思想。

    假设已知系统模型以及过程噪声、观测噪声特征(设定为高斯白噪声，其他类型噪声也可以通过白化滤波器转化为高斯白噪声)，于是可以根据模型和噪声特征进行系统下一时刻状态量、观测量及其方差的预测(先验值)，再结合实际测量值的特征，利用预测值和观测量对状态量的真实值进行估计(后验值)，最终得到最优值使状态量的协方差最小(与维纳滤波一致)。根据上述原理推导出的卡尔曼滤波的步骤如下图所示，推导过程和公式已经被写了无数遍，这里就不重复写了，可以参考文章，解释的比较清楚：

            卡尔曼滤波原理与公式推导

            卡尔曼滤波视频

    上述公式中有一个重要的中间变量：卡尔曼增益 K，这个变量表征了状态最优估计过程中模型预测误差(Predicted Error)与观测误差(Measurement Error)的比重关系，哪个误差小就在决定状态量最终估计值的时候比重更大。K=0表示预测误差为0，系统的状态值完全取决与预测值；而K=1表示观测误差为0，系统的状态值完全取决于量侧值。

    卡尔曼滤波性能好、收敛快，可快速准确跟踪状态变量的变化，因此广泛的用于位置轨迹跟踪、系统辨识、自适应滤波等领域，对于音频来说也可以用于降噪、回声消除等算法，将滤波器参数当做空间状态变量来跟踪，最终得到最优解，比如：基于卡尔曼滤波器的回声消除算法

    维纳滤波和卡尔曼滤波的异同：

    (1) 维纳滤波要求已知观测数据的自相关函数和观测数据与期望信号的互相关函数，通过建立滤波器进行滤波处理得到期望信号的估计值.

    (2) 卡尔曼滤波要求已知观测数据与期望数据之间的数学模型(状态方程和量测方程)，通过递推算法得到期望信号的估计值。

    (3) 维纳滤波的解是以H(z)或h(n)的形式给出，而卡尔曼滤波的解是以状态变量的估计值给出。

    (4) 维纳滤波与卡尔曼滤波都是解决最佳线性滤波问题的，并且准则都是最小均方误差，优化目标一致，因此在信号平稳的前提下，理论上二者的稳态结果是一致的。但卡尔曼滤波有一个过渡过程，在过渡期间其结果与维纳滤波计算值不完全相同。

    (5) 维纳滤波的求解和工作要求数据是平稳随机过程，而卡尔曼不需要。

3、EKF扩展卡尔曼滤波

    扩展卡尔曼滤波EKF(Extend Kalman Filter)将卡尔曼滤波器扩展到非线性系统，主要变化是状态转换模型和观测模型可以是非线性的，计算时利用一阶泰勒展开式(Jacobians矩阵)将模型线性化。因此两者的主要区别是状态估计方程和状态估计转移方程有差异，而协方差矩阵的公式推导过程与KF一致：

    预测：

    使用Jacobians矩阵更新系统线性化模型：

    利用新的线性化模型进行后续更新计算，基本与卡尔曼滤波一致：

    虽然EKF可以处理非线性系统，但存在以下缺点:

    1、Jacobians矩阵解析解计算困难

    2、利用数值计算的方法得到Jacobians开销巨大

    3、工作区间内系统函数需要严格的连续可微

    4、线性化只使用一阶泰勒展开，误差很可能比较大，导致最后的滤波结果误差增大，特别是对于高度非线性的系统，各个问题难以解决。

    5、系统在某些状态点的线性化方程可能并不能反应真实的系统模型，比如线性化点附近，梯度急剧变化。

4、UKF无迹卡尔曼滤波

    无迹卡尔曼滤波UKF(Unscented Kalman Filter)在KF的基础上加入了UT变换来估计协方差矩阵，是无损变换UT(Unscented Transform)与标准卡尔曼滤波体系的结合，通过UT变换使非线性系统方程适用于线性假设下的标准卡尔曼体系。

    因此UKF和EKF都是对KF的扩展，用于处理非线性问题。EKF是利用泰勒展开将模型线性化来近似非线性模型，而UKF则是直接寻找一个与真实分布近似的高斯分布，没有用线性表征，因此UKF避免了Jacobian matrix矩阵的计算，不要求可微，甚至不要求连续，也就不存在上述EKF介绍中存在的问题，可以用于任何非线性模型，是一种适用性更强的扩展，在实际中也经常使用。下面介绍一下UT的原理和计算过程：

    (1) 基本原理：根据当前带噪状态量的高斯分布取若干个Sigma-Point，每个点有自己的权重。这些点经过非线性模型计算后得到预测值的分布，最后根据预测值的分布和相应的权重估计出新的协方差矩阵，通过概率统计处理避免了EKF中模型线性化、Jacobbis矩阵计算等步骤。

    (2) Sigma-Point的选取：从带噪输入状态量的分布中进行采样，注意这里不是随机采样，每个采样点距离中心值的距离都是标准差σ的整数倍，因此这些采样点被称为Sigma-Point。然后经过y=h(x)非线性计算得到预测的输出分布。通常Sigma-Point选取的个数与高斯分布的维度相关，N维的高斯分布可选择2N+1个Sigma-Point(一个中心值，其它点关于中心值对称分布)。

    (3) UKF计算步骤：上面描述了UT变化的原理和Sigma-Point选取的方法，将其运用到Kalman-Filter中得到以下计算步骤，

    非线性状态方程：

    非线性观测方程：

    Sigma-Point选取：

    根据Sigma-Point计算先验预测值，并计算Sigma-Point权重、先验状态均值和方差：

    根据Sigma-Point预测先验观测值，并计算先验观测值均值和方差：

    计算协方差和卡尔曼增益：

    刷新估计值和方差矩阵：

四、数值优化(数学规划)思想

1、转化为优化问题

    在前面描述中可以看出，大部分自适应滤波过程其实也可以看做是一个数值优化问题。这样就可以把思路打开，将数值优化中的方法应用到自适应滤波中， 于是就需要将相关问题转化为优化问题中的概念，比如：

    优化目标：以LMS为例，优化目标即预测值与实际值的误差均方值MSE最小。其实在其他问题中，也可以采用对数误差等等，这和机器学习中loss的概念一样，不同loss对最后优化结果、优化侧重点和寻优方法的选择、优化效果等方面有很大影响，应该考虑实际情况(系统模型、平台、开销限制、期望效果等)选择合适的目标函数进行优化。

    单/多目标：实际问题中滤波器优化的目标可能不只是使误差最小，即单目标优化问题。还有可能是同时存在许多优化目标和约束的问题，即多目标优化问题，比如滤波器某个频段的最大增益存在约束或者对滤波器参数数值范围有约束等等，当遇到多目标问题时需要确定合适的优化顺序和处理方法。

    多目标转化为多个单目标问题：有一些方法可以将多目标问题转化为多个单目标问题，并通过约束来保持各个优化问题结果的同步性，比如ADMM(Alternating Direction Multiplier Method)等

    优化顺序：当多目标存在时，优化目标的顺序需要规定好，比如是按照权重进行trade-off还是进行字典序优化。

    Hard/Soft Constarint：表示约束是一定要强制完全符合还是尽量符合，可转化为拉格朗日算子等问题。

    我们需要根据实际情况，将问题转化为合适的优化问题，然后使用合适的方法来解决问题，才能得到比较好的效果。一个比较经典的应用场景为模型预测控制MPC(Model Predictive Control)，可以比较好的解决MIMO问题，并且可以处理非线性系统。模型预测控制简介

2、优化方法

(1) 梯度下降

    线性系统的自适应滤波器问题基本都可以转化为凸优化问题，因此可以很容易使用梯度下降寻优的一些方法，LMS其实就算是使用了牛顿梯度的一个优化问题。于是根据具体情况，可以使用一些更为复杂梯度方法来代替牛顿梯度，以达到更快的收敛效率或更好的收敛效果，比如Gauss-Newton梯度、随机梯度、dogleg等等

(2) 启发式算法

    在开销可以接受的情况下，当然也可以使用GA、PSO、鱼群、蚁群、退火等启发式算法来解决滤波器参数优化的问题。好处是定义优化目标函数时，不用过多考虑梯度计算、连续性、可导性、参数互相耦合等问题，问题也可以是非凸优化，适用性更广。但同时也存在收敛速度慢、开销大、有限时间内最终效果比较随机的问题，实际使用过程中需要很好的设计、利用很多修饰才能保证收敛速度和结果的准确性。