【现代信号处理】09 - 自适应与LMS

自适应 Adaptive

1. 最小二乘的两个难题

1.1 概述

  我们从最小二乘开始。最小二乘有这样的基本模型。

  假设我们有这样的数据

Z = ( z 1 , z 2 , . . . , z n ) T Z = (z_1,z_2,...,z_n)^T Z=(z1,z2,...,zn)T

  我们希望用这组数据去估计我们的目标。我们希望能够找到一组数据,来极小化我们的目标函数。最终我们得到的解就是最小二乘解。

ω = ( ω 1 , ω 2 , . . . , ω n ) T m i n E ( d − Z T ω ) 2 ⇒ ω L S \omega = (\omega_1,\omega_2,...,\omega_n)^T \\ min E(d-Z^T \omega)^2 \Rightarrow \omega_{LS} ω=(ω1,ω2,...,ωn)TminE(dZTω)2ωLS
  其实我们得到最小二乘解的时候,没有遇到过多的问题。但是实际上,我们可能会面临两方面的困难。

  首先我们要面对的问题是,目标函数可能具有自身的复杂性。因为我们是以均方误差作为目标函数的,并且我们采用的是线性估计。均方和线性使得我们在求解目标函数的时候没有遇到很多问题。我们也可以通过引入对角加载项的方法,让最小二乘进行正则化。但是,我们之前做的种种努力,都没有脱离均方和线性的范围。

  但是实际情况是,我们的目标函数可能并不是均方误差。这就使得,我们的目标函数求梯度,可能没有解析解。

J ( ω ) ⇒ ∇ J ( ω ) = 0 ⇒ ω J(\omega) \Rightarrow \nabla J(\omega) = 0 \cancel \Rightarrow \omega J(ω)J(ω)=0 ω

  其次,我们可能会面临非稳态的环境。这种非稳态来源于两部分。首先,可能我们的目标函数本身是在变的。其次,我们的数据也是时变的,其统计特性会随着时间的改变而发生改变。这样的话,我们就不能一劳永逸的求出估计了。我们就需要在每一个时刻,都进行重新求解。

  一旦引入了这两个问题,我们之前做的努力,就有待商榷了。

1.2 Complexity of Objective Function

  首先我们来看待第一个问题。如果目标函数得不得解析解怎么办。

  这个时候我们可以引入数值优化的方法。

1.2.1 最速下降法

  最经典的就是最速下降法。Steepest Descent

  最速下降法是一种迭代的方法,我们从某一个初值开始,逐次进行迭代,每一次迭代都进行某种计算。我们希望我们的数据能够收敛到最优值。

ω 0 → ω 1 → ω 2 → . . . → ω n \omega_0 \rightarrow \omega_1 \rightarrow \omega_2 \rightarrow... \rightarrow \omega_n ω0ω1ω2...ωn

ω n → ω o p t ⇒ ∇ J ( ω o p t ) = 0 \omega_n \rightarrow \omega_{opt} \Rightarrow \nabla J(\omega_{opt}) = 0 ωnωoptJ(ωopt)=0

  迭代的话,就会有数值方程,我们假设已经有了ωk了,我们想知道ωk和ωk+1之间有什么联系。

ω k + 1 = g ( ω k ) \omega_{k+1} = g(\omega_k) ωk+1=g(ωk)

  我们可以通过泰勒展开来寻找这种关系

J ( ω k + 1 ) = J ( ω k ) + ( ∇ J ( ω k ) ) T ( ω k + 1 − ω k ) + o ( ∣ ∣ ω k + 1 − ω k ∣ ∣ 2 ) J(\omega_{k+1}) = J(\omega_k) + (\nabla J(\omega_k))^T(\omega_{k+1} - \omega_k) +o(||\omega_{k+1}- \omega_k||^2) J(ωk+1)=J(ωk)+(J(ωk))T(ωk+1ωk)+o(ωk+1ωk2)

  我们假定高阶无穷小可以忽略,只考虑线性的部分。

J ( ω k + 1 ) = J ( ω k ) + ( ∇ J ( ω k ) ) T ( ω k + 1 − ω k ) J(\omega_{k+1}) = J(\omega_k) + (\nabla J(\omega_k))^T(\omega_{k+1} - \omega_k) J(ωk+1)=J(ωk)+(J(ωk))T(ωk+1ωk)

  我们想让这个式子尽可能的小,也就是

m i n ( ∇ J ( ω k ) ) T ( ω k + 1 − ω k ) min(\nabla J(\omega_k))^T(\omega_{k+1} - \omega_k) min(J(ωk))T(ωk+1ωk)

  根据柯西不等式等号成立条件
ω k + 1 − ω k = − μ ∇ ( ω k ) \omega_{k+1} - \omega_k = -\mu \nabla(\omega_k) ωk+1ωk=μ(ωk)

  我们就到了迭代关系
ω k + 1 = ω k − μ ∇ ( ω k ) ( 1 ) \omega_{k+1} = \omega_k-\mu \nabla(\omega_k) \quad\quad(1) ωk+1=ωkμ(ωk)(1)

1.2.2 牛顿法

  最速下降法是一个一阶的优化方法,如果我们不满足于一阶的优化方法,我们想继续做二阶的。我们可以进行二阶展开

J ( ω k + 1 ) = J ( ω k ) + ( ∇ J ( ω k ) ) T ( ω k + 1 − ω k ) + 1 2 ( ω k + 1 − ω k ) T ∇ 2 J ( ω k ) ( ω k + 1 − ω k ) J(\omega_{k+1}) = J(\omega_k) + (\nabla J(\omega_k))^T(\omega_{k+1} - \omega_k) +\frac{1}{2} (\omega_{k+1} - \omega_k)^T \nabla^2 J(\omega_k) (\omega_{k+1} - \omega_k) J(ωk+1)=J(ωk)+(J(ωk))T(ωk+1ωk)+21(ωk+1ωk)T2J(ωk)(ωk+1ωk)

  其中二阶项是海森矩阵

∇ 2 J ( ω k ) ⇒  Hessian Matrix \nabla^2 J(\omega_k) \Rightarrow \text{ Hessian Matrix} 2J(ωk) Hessian Matrix

  一阶展开只有线性的部分。二阶展开以后,我们得到了一个二次型。现在,如果我们想要让目标函数尽可能小,我们只要对后面的部分进行求导即可

m i n [ ( ∇ J ( ω k ) ) T ( ω k + 1 − ω k ) + 1 2 ( ω k + 1 − ω k ) T ∇ 2 J ( ω k ) ( ω k + 1 − ω k ) ] min[(\nabla J(\omega_k))^T(\omega_{k+1} - \omega_k) +\frac{1}{2} (\omega_{k+1} - \omega_k)^T \nabla^2 J(\omega_k) (\omega_{k+1} - \omega_k)] min[(J(ωk))T(ωk+1ωk)+21(ωk+1ωk)T2J(ωk)(ωk+1ωk)]

∇ ω k + 1 ( ( ∇ J ( ω k ) ) T ( ω k + 1 − ω k ) + 1 2 ( ω k + 1 − ω k ) T ∇ 2 J ( ω k ) ( ω k + 1 − ω k ) ) = 0 ∇ J ( ω k ) + 2 ∗ 1 2 ∇ 2 J ( ω k ) ∗ ( ω k + 1 − ω k ) = 0 ω k + 1 = ω k − μ ∇ ( ∇ 2 J ( ω k ) ) − 1 J ( ω k ) ( 2 ) \nabla \omega_{k+1}((\nabla J(\omega_k))^T(\omega_{k+1} - \omega_k) +\frac{1}{2} (\omega_{k+1} - \omega_k)^T \nabla^2 J(\omega_k) (\omega_{k+1} - \omega_k)) =0\\ \nabla J(\omega_k) + 2*\frac{1}{2}\nabla^2 J(\omega_k)*(\omega_{k+1} - \omega_k) =0 \\ \omega_{k+1} = \omega_k-\mu \nabla (\nabla^2 J(\omega_k))^{-1}J(\omega_k) \quad\quad(2) ωk+1((J(ωk))T(ωk+1ωk)+21(ωk+1ωk)T2J(ωk)(ωk+1ωk))=0J(ωk)+2212J(ωk)(ωk+1ωk)=0ωk+1=ωkμ(2J(ωk))1J(ωk)(2)

  这里需要注意一下

∇ x ( x T A x ) = ( A + A T ) x \nabla_x(x^T A x) = (A +A^T)x x(xTAx)=(A+AT)x

  这里求导出来了两倍是因为hessian矩阵是对称的,转置等于其本身。

  这个二阶优化的方法叫做牛顿法。因此,如果我们的目标函数解不出解析解,我们就只能通过数值计算的方法进行逼近了。

1.3 Non-Stationary Environment

  我们再来说说第二个问题,我们让下标k代表时间。现在我们的数据和目标都是随时间变化的了。并且求解的系数也是个时变的参数了。

k → T i m e Z ( k ) = ( z 1 ( k ) , . . . , z n ( k ) ) T d → d ( k ) ω → ω ( k ) k \rightarrow Time \\ Z(k) = (z_1(k),...,z_n(k))^T \\ d \rightarrow d(k) \\ \omega \rightarrow \omega(k) kTimeZ(k)=(z1(k),...,zn(k))Tdd(k)ωω(k)

  现在研究的对象全都是与时间相关的了,我们重写一下目标函数,然后进行求解就行

m i n E ( d ( k ) − Z T ( k ) ω ( k ) ) 2 ⇒ ω o p t ( k ) min E(d(k) - Z^T(k) \omega(k))^2 \Rightarrow \omega_{opt}(k) minE(d(k)ZT(k)ω(k))2ωopt(k)

2. 最小均方自适应滤波器LMS

2.1 问题引入

  我们重新审视一下我们提出的有关最小二乘的两个问题。

  假设我们的研究对象是非平稳的,如果优化函数还是均方的,我们要求目标函数可以简单的用最小二乘解。因为虽然ω是随着时间变化的,但是我们计算最小二乘解只需要一步就能求解完成。

  但是,如果我们的目标函数不是均方的,很复杂,而且还是非平稳的,目标函数会与时间相关,这个时候怎么办呢?因为,首先,我们的最优解是一步解不出来的,需要通过不断的迭代。而且,最优解又是与时间相关的,会随时变化。即使,我们可以假设,最优解的时间变化是个慢过程,计算是个快过程,计算机计算性能足够高,能够在最优解状态改变之前得到结果。如果这个假设成立,会对计算机产生极大的压力,而且这样的算法也不具备普适性,没有办法得到广泛的应用。

  因此,我们需要用其他办法来解决这个既不平稳,目标函数又复杂的问题。事实上,如果我们能够让迭代的尺度和时间的尺度联系起来,这个问题就能得到解决了。

  这就是我们今天要思考的核心问题。

  我们一开始做的事情是这样的。

在这里插入图片描述

  在一个状态下,对采样得到的数据,通过不断的迭代,获得最优的参数估计。只要没找到最优解,我们的迭代就不会结束。然后进入下一个状态,继续采样,对采样数据不断迭代,得到新状态的最优估计参数,然后依次类推。

  现在,我们不这么做了。我们每次获得一个采样点就迭代一次。采样一次以后,目标函数就变了,我们在新的目标函数下继续迭代一次,是一种阶梯状的迭代方式。这么做其实更有利于对变化的环境进行追踪,因为虽然我们没有在原系统中得到最优解,但是系统马上就发生改变了,得到最优解也没有什么意义。这就是今天要探讨的自适应问题。

在这里插入图片描述

2.2 LMS

2.2.1 模型建立

  现在,我们要让我们的算法去适应环境的变化。我们写一下新的目标函数

J ( ω ( k ) ) = E ( d ( k ) − Z T ( k ) ω ( k ) ) 2 J(\omega(k)) = E(d(k) - Z^T(k)\omega (k))^2 J(ω(k))=E(d(k)ZT(k)ω(k))2

  我们来求一下最速下降
∇ ω ( J ( ω ( k ) ) ) = − 2 Z ( k ) ∗ ( d ( k ) − Z T ( k ) ω ( k ) ) = − 2 Z ( k ) e ( k ) \nabla_{\omega}(J(\omega(k))) = -2Z(k)*(d(k) - Z^T(k)\omega (k)) = -2Z(k)e(k) ω(J(ω(k)))=2Z(k)(d(k)ZT(k)ω(k))=2Z(k)e(k)

  其中e(k)表示估计误差

e ( k ) = d ( k ) − Z T ( k ) ω ( k ) ⇒ error e(k) = d(k) - Z^T(k)\omega (k) \Rightarrow \text{error} e(k)=d(k)ZT(k)ω(k)error

  目前我们使用的符号k还都是表示时间,现在我们混用符号了。
ω k + 1 = ω k − μ ∇ J ( ω k ) \omega_{k+1} = \omega_{k} - \mu \nabla J(\omega_k) ωk+1=ωkμJ(ωk)

  现在我们要把计算得到的梯度进行代入

ω k + 1 = ω k − μ Z ( k ) e ( k ) \omega_{k+1} = \omega_{k} - \mu Z(k)e(k) \\ ωk+1=ωkμZ(k)e(k)

  注意,原本上,ω的角标k是迭代次数,其他的k表示时间,在经过我们的符号统一之后,我们给k赋予了时间的含义,也就是每次采样一次,就发生一次迭代,然后新的迭代是用来估计新状态的最优参数的。其实,我们每次采样就迭代一次的想法是正确的,因为每次采样之后,系统就发生改变了,就算我们把最优参数算出来,也没有多大的意义。因为系统马上就变了,约等于四九年入国军。

ω k + 1 = ω k − μ Z ( k ) ( d ( k ) − Z T ( k ) ω ( k ) ) ( i ) \omega_{k+1}= \omega_{k} - \mu Z(k)(d(k) - Z^T(k)\omega (k)) \quad\quad(i) ωk+1=ωkμZ(k)(d(k)ZT(k)ω(k))(i)

  (i)式就是非常有名的LMS自适应滤波器。(LMS = Least Mean Square)

2.2.2 LMS细节剖析

  其实自适应的目的就是为了把迭代和自身的发展实现统一。

  这里其实需要注意一下,因为看起来有些地方有些问题。就是我们在求梯度的时候,并没有写成期望的形式

J ( ω ( k ) ) = E ( d ( k ) − Z T ( k ) ω ( k ) ) 2 J(\omega(k)) = E(d(k) - Z^T(k)\omega (k))^2 J(ω(k))=E(d(k)ZT(k)ω(k))2

  我们求梯度是这么写的,也就是没有加期望。
∇ ω ( J ( ω ( k ) ) ) = − 2 Z ( k ) ∗ ( d ( k ) − Z T ( k ) ω ( k ) ) = − 2 Z ( k ) e ( k ) \nabla_{\omega}(J(\omega(k))) = -2Z(k)*(d(k) - Z^T(k)\omega (k)) = -2Z(k)e(k) ω(J(ω(k)))=2Z(k)(d(k)ZT(k)ω(k))=2Z(k)e(k)

  这样是不是错了呢?实际上,并没有,因为我们在做LMS的时候,每次只采用1次,也就是说Z只代表了一个数据。一个数据是没有办法求期望的,因此就索性不写期望符号了。

  这里我们使用了一个小技巧,就是令瞬时值代替了期望值。 instant value -> expectation value

  LMS是个革命性的算法,里面有两个非常重要的思想

  • LMS首先将迭代和时间的角标完成了统一,把采样和优化完成了统一
  • 其次LMS使用瞬时值取代了平均值

  后面人们发现,这两个做法是非常聪明的,这个算法也得到了非常广泛的应用。因为实现极其简单,并且效果出人意料的好。

  这个算法来自与widrow,这个人做了第一个感知器,而且用模拟电路实现了。对神经网络的发展功不可没。而且这个人谦虚的说,自己不懂数学。

2.3 LMS与学习曲线

2.3.1 学习曲线的建立

  所谓学习曲线,就是误差随着时间的变化规律,变化规律可能会有很多形式,比如线性的,指数的等等

在这里插入图片描述

  我们希望知道LMS的误差曲线是什么情况。我们来计算两件事情。

( 1 ) J ( k ) = E ( d ( k ) − Z T ( k ) ω k ) 2 ( 2 ) ξ ( k ) = E ( ω k − ω 0 ) 2 (1) \quad\quad J(k) = E(d(k) - Z^T(k) \omega_k)^2 \\ (2) \quad\quad\quad\quad\quad\quad \xi(k) = E(\omega_k - \omega_0)^2 (1)J(k)=E(d(k)ZT(k)ωk)2(2)ξ(k)=E(ωkω0)2

  其中,我们假设
d ( k ) = Z T ( k ) ω 0 + e 0 ( k ) ( i i ) d(k) = Z^T(k) \omega_0 +e_0(k) \quad\quad(ii) d(k)=ZT(k)ω0+e0(k)(ii)

e 0 ( k ) : White Noise E ( e 0 ( k ) e 0 ( n ) ) = δ k n e_0(k):\text{White Noise} \\ E(e_0(k)e_0(n)) = \delta_{kn} e0(k):White NoiseE(e0(k)e0(n))=δkn
  这是一个简化的模型,我们假设目标值d(k)就是最优估计和白噪声。

  通常来说,我们的目标值d(k)应该具有这样的形式

d ( k ) = Z T ( k ) ω ( k ) d(k) = Z^T(k) \omega(k) d(k)=ZT(k)ω(k)

  这个模型其实是比较难处理的,我们就认为线性系数有一个最优的值可以去追寻。但是,如果直接把线性系数认为是不变的显然不太合理,于是我们又通过增加白噪声的方式,增加了一个加性扰动。

  因此我们的两个计算式子的计算意义如下

  • 第一个式子是在计算估计值和目标值相差多少
  • 第二个式子考察的就是,随着时间k的变化,我们的估计系数与最优的系数之间的距离是怎么减少的。这是一种跟踪的性能,我们跟踪系数,克服白噪声的影响,不断的逼近我们的最优系数。
2.3.2 估计系数误差的递推

  我们假设对估计的系数和最优的系数之间的距离是e,并且k时刻之间的距离是e(k)。这个距离就是误差
e ( k ) = ω k − ω 0 ( i i i ) e(k) = \omega_k - \omega_0 \quad\quad(iii) e(k)=ωkω0(iii)

  现在我们想寻求k+1时刻的距离,并建立递推关系

e ( k + 1 ) = ω k + 1 − ω 0 ( i v ) e(k+1) = \omega_{k+1} - \omega_0 \quad\quad(iv) e(k+1)=ωk+1ω0(iv)

  我们引入LMS,也就是把(i)式代入(iv)式

e ( k + 1 ) = ω k − μ Z ( k ) ( d ( k ) − Z T ( k ) ω k ) − ω 0 ( v ) e(k+1)= \omega_{k} - \mu Z(k)(d(k) - Z^T(k)\omega_k) - \omega_0 \quad\quad(v) e(k+1)=ωkμZ(k)(d(k)ZT(k)ωk)ω0(v)

  我们再引入d(k)的近似模型,就是把(ii)代入(v)

e ( k + 1 ) = ω k − ω 0 − μ Z ( k ) ( Z T ( k ) ω 0 + e 0 ( k ) − Z T ( k ) ω k ) = ω k − ω 0 + μ Z ( k ) Z ( k ) T ( ω k − ω 0 ) − μ Z ( k ) e 0 ( k ) = ( I + μ Z ( k ) Z ( k ) T ) ( ω k − ω 0 ) − μ Z ( k ) e 0 ( k ) ( v i ) e(k+1)= \omega_{k}- \omega_0 - \mu Z(k)( Z^T(k) \omega_0 +e_0(k) - Z^T(k)\omega_k) \\ =\omega_{k}- \omega_0 + \mu Z(k)Z(k)^T(\omega_k - \omega_0) - \mu Z(k)e_0(k) \\ = (I+\mu Z(k)Z(k)^T)(\omega_k - \omega_0)- \mu Z(k)e_0(k) \quad\quad(vi) e(k+1)=ωkω0μZ(k)(ZT(k)ω0+e0(k)ZT(k)ωk)=ωkω0+μZ(k)Z(k)T(ωkω0)μZ(k)e0(k)=(I+μZ(k)Z(k)T)(ωkω0)μZ(k)e0(k)(vi)

  再把(iii)代入(vi)可得

e ( k + 1 ) = ( I + μ Z ( k ) Z ( k ) T ) e ( k ) − μ Z ( k ) e 0 ( k ) ( v i i ) e(k+1) = (I+\mu Z(k)Z(k)^T)e(k) - \mu Z(k)e_0(k) \quad\quad(vii) e(k+1)=(I+μZ(k)Z(k)T)e(k)μZ(k)e0(k)(vii)

  我们就得到了误差的递推式子(vii)

  但是通过观测,我们发现这个式子中,每一个变量都是随机变量。如果递推系数本身是个确定值的话,递推比较容易进行,但是现在,递推系数也是个随机变量。我们希望能够把某些变量确定下来

  当μ足够小的时候,有这样的式子成立

R ≈ Z ( k ) Z T ( k ) μ < < 1 R \approx Z(k)Z^T(k) \\ \mu << 1 RZ(k)ZT(k)μ<<1

  这样做是有意义的。因为μ是步长,如果步长足够小,就可以对Z进行近似。

ω k + 1 = ω k − μ ∇ J \omega_{k+1} = \omega_k - \mu \nabla J ωk+1=ωkμJ

  加入近似条件以后(vii)就变成了

e ( k + 1 ) = ( I + μ R ) e ( k ) − μ Z ( k ) e 0 ( k ) e(k+1) = (I+ \mu R)e(k) - \mu Z(k) e_0(k) e(k+1)=(I+μR)e(k)μZ(k)e0(k)

  但是这个矩阵结构仍然非常复杂,我们可以对R做特征分解

R = Q T Λ Q R = Q^T \Lambda Q R=QTΛQ

e ( k + 1 ) = ( I + μ Q T Λ Q ) e ( k ) − μ Z ( k ) e 0 ( k ) e(k+1) = (I +\mu Q^T \Lambda Q)e(k) - \mu Z(k) e_0(k) \\ e(k+1)=(I+μQTΛQ)e(k)μZ(k)e0(k)

  因为特征矩阵是正交矩阵,满足

Q T ∗ Q = Q ∗ Q T = I Q^T*Q = Q*Q^T = I QTQ=QQT=I

  所以,我们又可以把式子写成

e ( k + 1 ) = Q T ( I + μ Λ ) Q e ( k ) − μ Z ( k ) e 0 ( k ) Q e ( k + 1 ) = ( I + μ Λ ) Q e ( k ) − μ Q Z ( k ) e 0 ( k ) [ 1 ] e(k+1) = Q^T (I +\mu\Lambda) Qe(k) - \mu Z(k) e_0(k) \\ Qe(k+1) = (I +\mu\Lambda)Qe(k) - \mu QZ(k) e_0(k) \quad\quad[1] e(k+1)=QT(I+μΛ)Qe(k)μZ(k)e0(k)Qe(k+1)=(I+μΛ)Qe(k)μQZ(k)e0(k)[1]

  这里做一次变量代换,令

v ( k ) = Q e ( k ) ϕ ( k ) = Q Z ( k ) e 0 ( k ) [ 2 ] v(k) = Q e(k) \\ \phi(k) = QZ(k)e_0(k) \quad\quad[2] v(k)=Qe(k)ϕ(k)=QZ(k)e0(k)[2]

  则[1]可以写成

v ( k + 1 ) = ( I + μ Λ ) v ( k ) − μ ϕ ( k ) [ 3 ] v(k+1) = (I +\mu \Lambda)v(k) - \mu \phi(k) \quad\quad[3] v(k+1)=(I+μΛ)v(k)μϕ(k)[3]

  由于v(k)和φ(k)可以写成列向量的形式

v ( k ) = ( v 1 ( k ) , . . . , v n ( k ) ) T ϕ ( k ) = ( ϕ 1 ( k ) , . . . , ϕ n ( k ) ) T v(k) = (v_1(k),...,v_n(k))^T \\ \phi(k) = (\phi_1(k),...,\phi_n(k))^T v(k)=(v1(k),...,vn(k))Tϕ(k)=(ϕ1(k),...,ϕn(k))T

  所以有
v i ( k + 1 ) = ( 1 + μ λ i ) v i ( k ) − μ ϕ i ( k ) [ 4 ] v_i(k+1) = (1+ \mu \lambda_i)v_i(k) - \mu \phi_i(k) \quad\quad[4] vi(k+1)=(1+μλi)vi(k)μϕi(k)[4]
  我们可以进行递推了

v i ( k ) = ( 1 + μ λ i ) k v i ( 0 ) − ∑ j = 0 k − 1 ( 1 + μ λ i ) k − 1 − j ( μ ϕ i ( j ) ) [ 5 ] v_i(k) = (1+ \mu \lambda_i)^k v_i(0) - \sum_{j=0}^{k-1}(1+ \mu \lambda_i)^{k-1-j}(\mu \phi_i(j)) \quad\quad[5] vi(k)=(1+μλi)kvi(0)j=0k1(1+μλi)k1j(μϕi(j))[5]
  这里的这个v是对我们的误差进行某种正交变换的结果,并不影响误差的模

  我们发现误差随时间的规律是个指数变化,如果我们要让误差减小,必须让底数小于1。我们来计算一下这个步长的条件

∣ 1 + μ λ i ∣ < 1 |1+ \mu \lambda_i| <1 1+μλi<1

⇒ − 1 < 1 + μ λ i < 1 ⇒ − 2 < μ λ i < 0 \Rightarrow -1<1+\mu\lambda_i <1 \\ \Rightarrow -2 < \mu \lambda_i <0 1<1+μλi<12<μλi<0

  因为正交阵的特征值一定大于0,所以我们得到了

⇒ − 2 λ i < μ < 0 ∣ μ ∣ < 2 λ i , ∀ 1 ≤ i ≤ n ⇒ ∣ μ ∣ < 2 λ m a x \Rightarrow -\frac{2}{\lambda_i} <\mu <0 \\ |\mu| < \frac{2}{\lambda_i},\forall 1\leq i \leq n \\ \Rightarrow |\mu| < \frac{2}{\lambda_{max}} λi2<μ<0μ<λi2,1inμ<λmax2

  因此,μ一定小于最大的特征值分之2

2.3.3 估计目标误差的递推

  我们研究的误差是两个。第一个是估计距离估计目标还有多远,第二个是估计系数距离真正的系数还有多远。我们刚才已经完成了第二个误差的计算,下面我们要开始计算第一个误差了。我们希望把这个两个误差能够联系在一起

d ( k ) − Z T ( k ) ω k = d ( k ) − Z T ( k ) ( ω k − ω 0 + ω 0 ) = d ( k ) − Z T ( k ) ω 0 − Z T e ( k ) = d 0 − Z T e ( k ) d(k) - Z^T(k) \omega_k = d(k) - Z^T(k)(\omega_k - \omega_0 + \omega_0) \\ = d(k)- Z^T(k) \omega_0 - Z^T e(k) \\ =d_0 - Z^T e(k) d(k)ZT(k)ωk=d(k)ZT(k)(ωkω0+ω0)=d(k)ZT(k)ω0ZTe(k)=d0ZTe(k)

  其中,d0表示最优的系数产生的估计误差,也就是白噪声e0(k)

d 0 = d ( k ) − Z T ( k ) ω 0 = e 0 ( k ) d_0 = d(k)- Z^T(k) \omega_0 = e_0(k) d0=d(k)ZT(k)ω0=e0(k)

  因此,目标误差可以写做

E ( d ( k ) − Z T ( k ) ω k ) 2 = E ( e 0 ( k ) − Z T ( k ) e ( k ) ) 2 = E [ e 0 ( k ) 2 ] − 2 ∗ E [ e 0 ( k ) ∗ Z T ( k ) ∗ e ( k ) ] + E [ Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ] < 1 > E(d(k) - Z^T(k) \omega_k)^2 = E(e_0(k) - Z^T(k) e(k))^2 \\ = E[e_0(k)^2] - 2 *E[e_0(k) * Z^T(k) * e(k)] + E[Z^T(k) e(k) * e(k)^T Z(k)] \quad\quad<1> E(d(k)ZT(k)ωk)2=E(e0(k)ZT(k)e(k))2=E[e0(k)2]2E[e0(k)ZT(k)e(k)]+E[ZT(k)e(k)e(k)TZ(k)]<1>

  我们发现这个式子可以写成三部分。

  • 第一部分:是最优系数产生的误差,也就是最小误差,是个常数
  • 第二部分:是交叉项,因为白噪声和后面的乘积都不相关,这个交叉项实际是0
  • 第三部分:是我们需要认真研究的部分

  重写一下式<1>

E ( d ( k ) − Z T ( k ) ω k ) 2 = J m i n + E [ Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ] < 2 > E(d(k) - Z^T(k) \omega_k)^2 = J_{min} + E[Z^T(k) e(k) * e(k)^T Z(k)] \quad\quad<2> E(d(k)ZT(k)ωk)2=Jmin+E[ZT(k)e(k)e(k)TZ(k)]<2>

  我们把第三项单独拿出来分析

E [ Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ] < 3 > E[Z^T(k) e(k) * e(k)^T Z(k)] \quad\quad<3> E[ZT(k)e(k)e(k)TZ(k)]<3>

  我们发现,ZTe和eTZ都是个标量。标量乘以标量还是个标量。由于我们知道,标量的迹就是它本身,所以我们可以把式子变成迹的样子。

E [ Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ] = E [ T r ( Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ) ] < 4 > E[Z^T(k) e(k) * e(k)^T Z(k)] =E[Tr(Z^T(k) e(k) * e(k)^T Z(k))] \quad\quad<4> E[ZT(k)e(k)e(k)TZ(k)]=E[Tr(ZT(k)e(k)e(k)TZ(k))]<4>

  矩阵的迹就是矩阵对角线主元之和,也等于特征值之和

T r ( A ) = ∑ k = 1 n A k k = ∑ k = 1 n λ k Tr(A)= \sum_{k=1}^n A_{kk} = \sum_{k=1}^n \lambda_{k} Tr(A)=k=1nAkk=k=1nλk

  迹具有这样的性质,求迹的矩阵可以交换乘积顺序

T r ( A B ) = T r ( B A ) Tr(AB) = Tr(BA) Tr(AB)=Tr(BA)

  因此<4>式可以进行这样的变形

E [ Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ] = E [ T r ( Z T ( k ) e ( k ) e ( k ) T ∗ Z ( k ) ) ] = E [ T r ( Z ( k ) ∗ Z T ( k ) e ( k ) e ( k ) T ) ] = T r [ E ( Z ( k ) Z T ( k ) ∗ e ( k ) e ( k ) T ) ] E[Z^T(k) e(k) * e(k)^T Z(k)] =E[Tr(Z^T(k) e(k) e(k)^T *Z(k))] \\ = E[Tr(Z(k)*Z^T(k) e(k) e(k)^T)] \\ = Tr[E(Z(k)Z^T(k) *e(k) e(k)^T)] E[ZT(k)e(k)e(k)TZ(k)]=E[Tr(ZT(k)e(k)e(k)TZ(k))]=E[Tr(Z(k)ZT(k)e(k)e(k)T)]=Tr[E(Z(k)ZT(k)e(k)e(k)T)]

  相当于前三项作为一个整体,与第四项交换了乘积顺序。

  下面我们要引入玄妙的一步。我们假设Z和e是独立的。也就是前后两部分能分开算。这个假设其实是比较离谱的,但是最终得到的结果却是比较贴近的。

E [ Z T ( k ) e ( k ) ∗ e ( k ) T Z ( k ) ] = T r [ E ( Z ( k ) Z T ( k ) ∗ e ( k ) e ( k ) T ) ] = T r [ E ( ( Z ( k ) Z T ( k ) ) ∗ E ( e ( k ) e ( k ) T ) ] = T r ( R E ( k ) ) < 5 > E[Z^T(k) e(k) * e(k)^T Z(k)] = Tr[E(Z(k)Z^T(k) *e(k) e(k)^T)] \\ = Tr[E((Z(k)Z^T(k))*E(e(k) e(k)^T)] \\ = Tr(RE(k)) \quad\quad <5> E[ZT(k)e(k)e(k)TZ(k)]=Tr[E(Z(k)ZT(k)e(k)e(k)T)]=Tr[E((Z(k)ZT(k))E(e(k)e(k)T)]=Tr(RE(k))<5>
  其中
R = E ( ( Z ( k ) Z T ( k ) ) E ( k ) = E ( e ( k ) e ( k ) T ) R=E((Z(k)Z^T(k)) \\ E(k) = E(e(k) e(k)^T) R=E((Z(k)ZT(k))E(k)=E(e(k)e(k)T)

  把<5>代入<2>中得

J ( k ) = J m i n + T r ( R E ( k ) ) < 6 > J(k) = J_{min} + Tr(RE(k)) \quad\quad<6> J(k)=Jmin+Tr(RE(k))<6>

  我们希望把J和v联系在一起

T r ( R E ( k ) ) = T r [ E ( Z ( k ) Z T ( k ) e ( k ) ∗ e ( k ) T ) ] = T r [ E ( e ( k ) T ∗ Z ( k ) Z T ( k ) e ( k ) ) ] = T r [ E ( e ( k ) T ∗ E ( Z ( k ) Z T ( k ) ) ∗ e ( k ) ) ] = T r [ E ( e ( k ) T ∗ R ∗ e ( k ) ) ] = T r [ E ( e ( k ) T ∗ Q T Λ Q ∗ e ( k ) ) ] = T r [ E ( v T ( k ) ∗ Λ ∗ v ( k ) ) ] = T r [ E ( ∑ i = 1 k ( λ i ∗ ∣ v i ( k ) ∣ 2 ) ) ] Tr(RE(k)) = Tr[E(Z(k)Z^T(k)e(k) *e(k)^T)] \\ = Tr[E(e(k)^T*Z(k)Z^T(k)e(k))] \\ =Tr[E(e(k)^T*E(Z(k)Z^T(k))*e(k))] \\ = Tr[E(e(k)^T *R *e(k))] \\ =Tr[E(e(k)^T *Q^T \Lambda Q *e(k))] \\ = Tr [E(v^T(k)*\Lambda* v(k))] \\ = Tr[E(\sum_{i=1}^k(\lambda_i* |v_i(k)|^2))] Tr(RE(k))=Tr[E(Z(k)ZT(k)e(k)e(k)T)]=Tr[E(e(k)TZ(k)ZT(k)e(k))]=Tr[E(e(k)TE(Z(k)ZT(k))e(k))]=Tr[E(e(k)TRe(k))]=Tr[E(e(k)TQTΛQe(k))]=Tr[E(vT(k)Λv(k))]=Tr[E(i=1k(λivi(k)2))]

  因为已经是标量相加的结果了,迹可以去掉了

T r [ E ( ∑ i = 1 k ( λ i ∗ ∣ v i ( k ) ∣ 2 ) ) ] = ∑ i = 1 k ( λ i ∗ E ( ∣ v i ( k ) ∣ 2 ) ) < 7 > Tr[E(\sum_{i=1}^k(\lambda_i* |v_i(k)|^2))] = \sum_{i=1}^k(\lambda_i* E(|v_i(k)|^2)) \quad\quad<7> Tr[E(i=1k(λivi(k)2))]=i=1k(λiE(vi(k)2))<7>

  <7>代入<6>中

J ( k ) = J m i n + T r ( R E ( k ) ) = J m i n + ∑ i = 1 k ( λ i ∗ E ( ∣ v i ( k ) ∣ 2 ) ) < 8 > J(k) = J_{min}+ Tr(RE(k)) = J_{min} + \sum_{i=1}^k(\lambda_i* E(|v_i(k)|^2)) \quad\quad<8> J(k)=Jmin+Tr(RE(k))=Jmin+i=1k(λiE(vi(k)2))<8>

  我们对v继续进行处理,代入[5]式的结果,因为我们在计算系数误差的时候已经得到了v的具体表达式了

E ( ∣ v i ( k ) ∣ 2 ) = E [ ( 1 + μ λ i ) k v i ( 0 ) − ∑ j = 0 k − 1 ( 1 + μ λ i ) k − 1 − j ( μ ϕ i ( j ) ) ] 2 = ( 1 + μ λ i ) 2 k ∣ v i ( 0 ) ∣ 2 − 2 ∗ ( 1 + μ λ i ) k v i ( 0 ) ∗ E [ ∑ j = 0 k − 1 ( 1 + μ λ i ) k − 1 − j ( μ ϕ i ( j ) ) ] + E [ ∑ j = 0 k − 1 ∑ s = 0 k − 1 ( 1 + μ λ i ) k − 1 − j ( 1 + μ λ i ) k − 1 − s ( μ ϕ i ( j ) ) ( μ ϕ i ( s ) ) ] < 9 > E(|v_i(k)|^2) = E[(1+ \mu \lambda_i)^k v_i(0) - \sum_{j=0}^{k-1}(1+ \mu \lambda_i)^{k-1-j}(\mu \phi_i(j)) ]^2 \\ = (1+ \mu \lambda_i)^{2k} |v_i(0)|^2 - 2*(1+ \mu \lambda_i)^k v_i(0)*E[\sum_{j=0}^{k-1}(1+ \mu \lambda_i)^{k-1-j}(\mu \phi_i(j))] \\+ E[\sum_{j=0}^{k-1}\sum_{s=0}^{k-1}(1+ \mu \lambda_i)^{k-1-j}(1+ \mu \lambda_i)^{k-1-s}(\mu \phi_i(j))(\mu \phi_i(s))] \quad\quad<9> E(vi(k)2)=E[(1+μλi)kvi(0)j=0k1(1+μλi)k1j(μϕi(j))]2=(1+μλi)2kvi(0)22(1+μλi)kvi(0)E[j=0k1(1+μλi)k1j(μϕi(j))]+E[j=0k1s=0k1(1+μλi)k1j(1+μλi)k1s(μϕi(j))(μϕi(s))]<9>

  首先,我们来看一下交叉项,由于e0与其他的东西不相关,求期望是0,所以交叉项是0

E ( μ ϕ i ( j ) ) = E ( Q μ Z ( j ) e 0 ( j ) ) = 0 < 10 > E(\mu \phi_i(j)) = E(Q \mu Z(j)e_0(j)) = 0 \quad\quad<10> E(μϕi(j))=E(QμZ(j)e0(j))=0<10>

  后面部分中

E ( ϕ i ( j ) ϕ i ( s ) ) = 0 , j = s E(\phi_i(j) \phi_i(s)) = 0, j \cancel= s E(ϕi(j)ϕi(s))=0,j= s
  因此

E ( μ 2 ϕ i ( j ) ϕ i ( s ) ) = μ 2 Q T ∗ Q E ( e 0 ( j ) ) 2 E ( Z i ( j ) ) 2 = μ 2 E ( e 0 ( j ) ) 2 E ( Z i ( j ) ) 2 < 11 > E(\mu^2\phi_i(j) \phi_i(s)) = \mu^2 Q^T*Q E(e_0(j))^2 E(Z_i(j))^2 \\ = \mu^2 E(e_0(j))^2 E(Z_i(j))^2 \quad\quad<11> E(μ2ϕi(j)ϕi(s))=μ2QTQE(e0(j))2E(Zi(j))2=μ2E(e0(j))2E(Zi(j))2<11>

  这里面由于Q是正交阵,乘积为1

  把<10>和<11>代入 <9>中

E ( ∣ v i ( k ) ∣ 2 ) = ( 1 + μ λ i ) 2 k ∣ v i ( 0 ) ∣ 2 + ∑ j = 0 k − 1 ( 1 + μ λ i ) 2 ( k − 1 − j ) ∗ μ 2 E ( e 0 ( j ) ) 2 E ( Z i ( j ) ) 2 < 12 > E(|v_i(k)|^2) = (1+ \mu \lambda_i)^{2k} |v_i(0)|^2 + \sum_{j=0}^{k-1}(1+ \mu \lambda_i)^{2(k-1-j)}*\mu^2 E(e_0(j))^2 E(Z_i(j))^2 \quad\quad<12> E(vi(k)2)=(1+μλi)2kvi(0)2+j=0k1(1+μλi)2(k1j)μ2E(e0(j))2E(Zi(j))2<12>

  其中

J m i n = E [ ( e 0 ( j ) ) 2 ] < 13 > J_{min} = E[(e_0(j))^2] \quad\quad<13> Jmin=E[(e0(j))2]<13>
  我们继续做个假设,这个假设影响不是非常大,就是结果看起来比较简洁

E [ ( Z i ( j ) ) 2 ] = λ i < 14 > E[(Z_i(j))^2] = \lambda_i \quad\quad<14> E[(Zi(j))2]=λi<14>

  把<14>和<13>代入<12>中可得

E ( ∣ v i ( k ) ∣ 2 ) = ( 1 + μ λ i ) 2 k ∣ v i ( 0 ) ∣ 2 + μ 2 J m i n λ i ∑ j = 0 k − 1 ( 1 + μ λ i ) 2 ( k − 1 − j ) < 15 > E(|v_i(k)|^2) = (1+ \mu \lambda_i)^{2k} |v_i(0)|^2 +\mu^2 J_{min} \lambda_i \sum_{j=0}^{k-1}(1+ \mu \lambda_i)^{2(k-1-j)} \quad\quad<15> E(vi(k)2)=(1+μλi)2kvi(0)2+μ2Jminλij=0k1(1+μλi)2(k1j)<15>

  我们发现<15>后面部分是个等比数列求和问题,因此,我们可以把<15>继续写成

E ( ∣ v i ( k ) ∣ 2 ) = ( 1 + μ λ i ) 2 k ∣ v i ( 0 ) ∣ 2 + μ 2 ∗ λ i ∗ J m i n ∗ 1 − ( 1 + μ λ i ) 2 k 1 − ( 1 + μ λ i ) 2 = ( 1 + μ λ i ) 2 k [ ∣ v i ( 0 ) ∣ 2 − μ 2 ∗ λ i ∗ J m i n 1 − ( 1 + μ λ i ) 2 ] + μ 2 ∗ λ i ∗ J m i n 1 − ( 1 + μ λ i ) 2 < 16 > E(|v_i(k)|^2) = (1+ \mu \lambda_i)^{2k} |v_i(0)|^2 +\mu^2 *\lambda_i*J_{min} * \frac{1-(1+\mu \lambda_i)^{2k}}{1-(1+\mu \lambda_i)^{2}} \\ = (1+ \mu \lambda_i)^{2k}[|v_i(0)|^2 - \frac{\mu^2 *\lambda_i*J_{min}}{1-(1+\mu \lambda_i)^{2}}] + \frac{\mu^2 *\lambda_i*J_{min}}{1-(1+\mu \lambda_i)^{2}} \quad\quad<16> E(vi(k)2)=(1+μλi)2kvi(0)2+μ2λiJmin1(1+μλi)21(1+μλi)2k=(1+μλi)2k[vi(0)21(1+μλi)2μ2λiJmin]+1(1+μλi)2μ2λiJmin<16>

  其中

μ 2 ∗ λ i ∗ J m i n 1 − ( 1 + μ λ i ) 2 = μ 2 ∗ λ i ∗ J m i n − 2 μ λ i − μ 2 λ i 2 = μ J m i n − 2 − μ λ i < 17 > \frac{\mu^2 *\lambda_i*J_{min}}{1-(1+\mu \lambda_i)^{2}} = \frac{\mu^2 *\lambda_i*J_{min}}{-2\mu\lambda_i - \mu^2 \lambda_i^2} = \frac{\mu J_{min}}{-2-\mu \lambda_i } \quad\quad<17> 1(1+μλi)2μ2λiJmin=2μλiμ2λi2μ2λiJmin=2μλiμJmin<17>

   <17>代入<16>得

E ( ∣ v i ( k ) ∣ 2 ) = [ ∣ v i ( 0 ) ∣ 2 + μ J m i n 2 + μ λ i ] ( 1 + μ λ i ) 2 k − μ J m i n 2 + μ λ i E(|v_i(k)|^2) = [|v_i(0)|^2 +\frac{\mu J_{min}}{2+\mu \lambda_i } ] (1+ \mu \lambda_i)^{2k} - \frac{\mu J_{min}}{2+\mu \lambda_i} E(vi(k)2)=[vi(0)2+2+μλiμJmin](1+μλi)2k2+μλiμJmin

  最终我们得到了对目标估计的学习曲线方程

{ E ( ∣ v i ( k ) ∣ 2 ) = [ ∣ v i ( 0 ) ∣ 2 + μ J m i n 2 + μ λ i ] ( 1 + μ λ i ) 2 k − μ J m i n 2 + μ λ i J ( k ) = J m i n + ∑ i = 1 k ( λ i ∗ E ( ∣ v i ( k ) ∣ 2 ) ) \begin{cases} E(|v_i(k)|^2) = [|v_i(0)|^2 +\frac{\mu J_{min}}{2+\mu \lambda_i } ] (1+ \mu \lambda_i)^{2k} - \frac{\mu J_{min}}{2+\mu \lambda_i} \\ J(k) = J_{min} + \sum_{i=1}^k(\lambda_i* E(|v_i(k)|^2)) \end{cases} {E(vi(k)2)=[vi(0)2+2+μλiμJmin](1+μλi)2k2+μλiμJminJ(k)=Jmin+i=1k(λiE(vi(k)2))

3. 小结

  我们为了求得LMS自适应滤波器的学习曲线,做了很多很多的假设,看似LMS非常不靠谱,因为做的假设往往很离谱。但是实际上,LMS效果非常好,并且学习曲线与我们做了很多假设的理论计算基本吻合。

  自适应这个思想,成功解决了迭代优化和非稳态两件事情,实现了采用时间和迭代时间两个时间尺度的统一。我们要构造一个自适应算法可能并不是很难。但是要注意,我们做自适应算法的时候,要做两件事情,首先要做仿真,然后要做理论上的分析,计算这个自适应算法的学习曲线。在计算学习曲线的时候,要敢于假设,把最后学习曲线的计算公式变得非常简洁漂亮才行。 算法+仿真+数学分析才是一个完整的研究。

  • 0
    点赞
  • 3
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值