随机梯度下降算法简介

之前的章节中介绍了利用最速下降算法可以实现维纳滤波器的最优解（LMMSE），其最优解的形式为：
$$w_0=R^{-1}P$$

它基于两个假设：

1. 环境的联合平稳，即输入$u(n)$以及期望输出$d(n)$是联合平稳的。

2. 环境的有关信息，即$u(n)$构成的自相关矩阵$\mathbf{R}$，以及由$u(n)$和$d(n)$构成的互相关矩阵$\mathbf{P}$需要已知。

遗憾的是$\mathbf{R}$和$\mathbf{P}$所包含的环境信息不可获得或不可用，因此我们应该寻找一类具有适应未知环境统计变化能力的新算法，导出这些算法由以下两种不同的方法：

1. 一种是基于随机梯度下降的方法，将在此章节进行讨论。

2. 另一种是基于最小二乘的方法。

随机梯度下降算法表示从某一次自适应循环到下一次自适应循环的自适应滤波算法中梯度的"随机选择"。随机梯度下降算法和最速下降算法都是局部滤波方法。

严格意义上不把随机梯度看成一种优化算法，理由如下：

因为梯度上由随机噪声的存在，它永远也无法达到凸优化问题所期望的最优解。相反，一旦它到达最优解的附近，它将用一种随机游走的方式持续在该解周围徘徊，从而无法停留在一个平衡点上。

故随机梯度下降算法是次优的。

LMS算法结构推导

最小均方LMS算法作为随机梯度算法的第一个应用，其独特特征可以概括如下：

1. LMS算法简易，复杂度和FIR滤波器的维数呈线性相关。

2. 不用要知道P和R

3. 在位置环境干扰下，具有一定的鲁棒性。

4. 不需要相关矩阵求逆，所以它比RLS更加简单。

基于以上特征，LMS是目前最流行的自适应滤波算法之一，下面是LMS算法结构图：

它包含了3个部分：

1. FIR滤波器，实现估计输出的计算：
   $$\hat{d}(n)=\sum _{k=0}^{M-1}u(n-k)\ast {w_k}^{\ast }(n)$$
2. 比较器，实现估计误差的计算

$$e(n)=d(n)-\hat{d}(n)$$

3. 自适应权值控制机制，实现滤波器系数的更新

$$w_k(n+1)=w_k(n)+\mu \ast u(n-k)\ast e^{\ast }(n),k=0,1,...,M-1$$
下面给出详细的推导过程。

回顾之前定义的代价函数（均方误差），期望算子对n时刻大量-统计-独立实现的瞬时平方估计，误差$|e(n){|}^2$进行集平均。
$$J(n)=E\left(|e(n){|}^2\right)$$

遗憾的是，在自适应滤波器的实际应用中，上述方式使用集平均是不可行的。

我们如此说是因为使用自适应滤波的动机是基于估计误差$e(n)$（随时间n而变）的单一实现，以在线方式来自适应未知环境的统计变化。这些做恰好和随机梯度下降算法是一回事，所以我们移除上面的期望算子，得到了新的代价函数：
$$J_s(n)=|e(n){|}^2=e(n)\ast e(n{)}^{\ast }$$

$J_s(n)$的下标s用来表示有别于集平均表达式中的$J(n)$。估计误差$e(n)$是随机过程的样本函数，其代价函数$J_s(n)$本身就是一个随机过程的样值。因此$J_s(n)$相对于FIR滤波器第k个抽头权值$w_k(n)$的导数也是随机的，下面会给出证明。代价函数$J_s(n)$对${w_k}^{\ast }(n)$求偏导可表示为：
$$\nabla J_{s,k}(n)=\frac{\partial J_s(n)}{\partial {w_k}^{\ast }(n)}=\frac{\partial (e(n)\ast e(n{)}^{\ast })}{\partial {w_k}^{\ast }(n)}=\frac{\partial ((d(n)-\hat{d}(n))\ast (d(n)-\hat{d}(n){)}^{\ast })}{\partial {w_k}^{\ast }(n)}$$
$$=\frac{\partial \left[\left(d(n)-{\sum }_{k=0}^{M-1}u(n-k)\ast {w_k}^{\ast }(n)\right)\ast {\left(d(n)-{\sum }_{k=0}^{M-1}u(n-k)\ast {w_k}^{\ast }(n)\right)}^{\ast }\right]}{\partial {w_k}^{\ast }(n)}$$
$$=\frac{\partial \left[\left(d(n)-{\sum }_{k=0}^{M-1}u(n-k)\ast {w_k}^{\ast }(n)\right)\ast \left(d^{\ast }(n)-{\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast w_k(n)\right)\right]}{\partial {w_k}^{\ast }(n)}$$
$$=\frac{\partial \left(\begin{array}{c}-d^{\ast }(n){\sum }_{k=0}^{M-1}u(n-k)\ast {w_k}^{\ast }(n)-d(n){\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast w_k(n)\\ +{\sum }_{k=0}^{M-1}u(n-k)\ast {w_k}^{\ast }(n)\ast {\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast w_k(n)\end{array}\right)}{\partial {w_k}^{\ast }(n)}$$

$$=\frac{\partial \left(\begin{array}{c}-d^{\ast }(n){\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)-d(n){\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\\ +{\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)\ast {\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\end{array}\right)}{\partial \left(a_k(n)-j{b}_k(n)\right)}$$

$$=\begin{array}{c}\frac{\partial \left(\begin{array}{c}-d^{\ast }(n){\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)-d(n){\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\\ +{\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)\ast {\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\end{array}\right)}{\partial a_k(n)}+\\ j\ast \frac{\partial \left(\begin{array}{c}-d^{\ast }(n){\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)-d(n){\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\\ +{\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)\ast {\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\end{array}\right)}{\partial b_k(n)}\end{array}$$
$$=\begin{array}{c}\begin{array}{c}-d^{\ast }(n)u(n-k)-d(n)u^{\ast }(n-k)+u(n-k)\ast {\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\\ +u^{\ast }(n-k){\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)\end{array}\\ +j\ast \left(\begin{array}{c}j\ast d^{\ast }(n)u(n-k)-j\ast d(n)u^{\ast }(n-k)-j\ast u(n-k)\ast {\sum }_{k=0}^{M-1}{u}^{\ast }(n-k)\ast \left(a_k(n)+j{b}_k(n)\right)\\ +j\ast u^{\ast }(n-k)\ast {\sum }_{k=0}^{M-1}u(n-k)\ast \left(a_k(n)-j{b}_k(n)\right)\end{array}\right)\end{array}$$
$$=\begin{array}{c}-d^{\ast }(n)u(n-k)-d(n)u^{\ast }(n-k)+u(n-k){\hat{d}}^{\ast }(n)+u^{\ast }(n-k)\hat{d}(n)\\ -d^{\ast }(n)u(n-k)+d(n)u^{\ast }(n-k)+u(n-k){\hat{d}}^{\ast }(n)-u^{\ast }(n-k)\hat{d}(n)\end{array}$$
$$=-2d^{\ast }(n)u(n-k)+2u(n-k){\hat{d}}^{\ast }(n)$$
$$=-2u(n-k){\left(d(n)-\hat{d}(n)\right)}^{\ast }$$
$$=-2u(n-k)\ast e^{\ast }(n)$$

重写上面的结论，LMS算法的代价函数的梯度算子是：
$$\nabla J_{s,k}(n)=-2u(n-k)\ast e^{\ast }(n),k=0,1,2,...,M-1$$

MSE（均方误差）算法的代价函数的梯度算子是：
$${\nabla }_{k}J=-2\ast E\left[u(n-k)\ast e^{\ast }(n)\right],k=0,1,2,...,M-1$$

发现两者的区别仅在于是否有期望算子，那是因为新的代价函数$J_s(n)$取消了集平均，没有了期望算子，所以对系数的偏微分也就是没有期望算子。也正是没有了这个期望算子，再次重复上面的一个结论，${\mathbf{J}}_{\mathbf{s}}(\mathbf{n})$相对于FIR滤波器第k个抽头权值${\mathbf{w}}_{\mathbf{k}}(\mathbf{n})$的导数也是随机的（随机噪声没有被期望所平均）。

基于上面的随机梯度$\nabla J_{s,k}(n)$，可得LMS算法的更新规则如下：
$$w_k(n+1)=w_k(n)-\frac{1}{2}\mu \nabla J_{s,k}(n),k=0,1,...,M-1$$

将随机梯度代入可得：
$$w_k(n+1)=w_k(n)+\mu \ast u(n-k)\ast e^{\ast }(n),k=0,1,...,M-1$$

简化成矩阵形式的

1. 系数更新规则：
   $$\mathbf{w}(n+1)=\mathbf{w}(n)+\mu \mathbf{u}(n)e^{\ast }(n),k=0,1,...,M-1$$

2. 估计误差：
   $$e(n)=d(n)-\hat{d}(n)$$

3. 估计输出
   $$\hat{d}(n)={\mathbf{w}}^H(n)\mathbf{u}(n)$$

其中
$$\mathbf{w}(n)=\left[\begin{array}{c}{w}_0(n)\\ w_1(n)\\ ...\\ w_{M-1}(n)\end{array}\right],\mathbf{u}(n)=\left[\begin{array}{c}u(n)\\ u(n-1)\\ ...\\ u(n-M+1)\end{array}\right]$$

同样的，这里的$\mu$的取值将会影响LMS算法中系数的稳定性和瞬态响应，从稳定性的角度出发，$\mu$应当取一个较小的正值。

小结

不管LMS算法和维纳滤波之间可能存在何种联系，一旦从维纳滤波器导出LMS算法时删除了期望算子E，他们之间的联系就完全破坏了。

最优性考虑

考虑LMS算法的最优性有两种办法：

1. 局部最优性

2. 全局次优性

局部最优性

之前提及随机梯度下降法是一种局部方法，它表明搜索的一个梯度不仅是随机的，而且是以一种局部的方式进行的。前者说梯度是随机的，可有梯度的表达式看出来，重写如下：

$$\nabla J_{s,k}(n)=-2u(n-k)\ast e^{\ast }(n),k=0,1,2,...,M-1$$

上式的估计误差就是随机量，所以梯度也是随机的。

只要从一次自适应循环到下一次自适应循环由数据引起的扰动充分小，那么且不管其随机性，LMS算法在欧氏意义下呈现局部最优性，为了解释，我们考虑以下两种场景

1. 当前估计误差表达式如下$$e(n)=d(n)-{\mathbf{w}}^H(n)\mathbf{u}(n)$$
2. 后验估计误差表达式如下$$r(n)=d(n)-{\mathbf{w}}^H(n+1)\mathbf{u}(n)$$
   注意后验估计误差有别于一步预测${\mathbf{w}}^H(n)\mathbf{u}(n+1)$

将后验估计误差减去当前估计误差可得（中间利用了$\mathbf{w}(n+1)=\mathbf{w}(n)+\mu \mathbf{u}(n)e^{\ast }(n)$ ）：

$$r(n)-e(n)=\left({\mathbf{w}}^H(n)-{\mathbf{w}}^H(n+1)\right)\mathbf{u}(n)={\left(\mathbf{w}(n)-\mathbf{w}(n+1)\right)}^H\mathbf{u}(n)={\left(-\mu \mathbf{u}(n)e^{\ast }(n)\right)}^H\mathbf{u}(n)=-\mu e(n){\mathbf{u}}^H(n)\mathbf{u}(n)=-\mu e(n){\Vert \mathbf{u}(n)\Vert }^2$$

移项重新整理可得：
$$r(n)=(1-\mu {\Vert \mathbf{u}(n)\Vert }^2)e(n)$$

为了让后验估计误差小于当前估计误差，那么就必须满足下列约束条件：
$$-1<1-\mu {\Vert \mathbf{u}(n)\Vert }^2<1$$

或者是
$$0<\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2<1$$

在此约束条件下，进一步寻找$\mathbf{w}(n+1)$ 与过去值$\mathbf{w}(n)$
的平方欧氏范数达到最小的$\mathbf{w}(n+1)$的最优值。因此局部最优化实际上就是一个约束最优化问题。

之前已经推得$r(n)-e(n)+\mu e(n){\Vert \mathbf{u}(n)\Vert }^2=0$，将后面的拆分成两个式子分别放在式子的左右，可得
$$r(n)-e(n)+\frac{1}{2}\mu e(n){\Vert \mathbf{u}(n)\Vert }^2=-\frac{1}{2}\mu e(n){\Vert \mathbf{u}(n)\Vert }^2$$

设$g(n)$ 为
$$g(n)=r(n)-e(n)+\frac{1}{2}\mu e(n){\Vert \mathbf{u}(n)\Vert }^2=-\frac{1}{2}\mu e(n){\Vert \mathbf{u}(n)\Vert }^2$$
$g(n)$ 进一步表达为：
$$g(n)=r(n)-\left(1-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)=-\frac{1}{2}\mu e(n){\Vert \mathbf{u}(n)\Vert }^2$$

重写之前的约束条件$0<\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2<1$，因此可得
$$-1<g(n)=r(n)-\left(1-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)=-\frac{1}{2}\mu e(n){\Vert \mathbf{u}(n)\Vert }^2<0,或-1<g(n)<0$$

接下来应用著名的拉格朗日乘子法，首先引入拉格朗日函数，其标准的表达式为
$$L(n)=f(n)+\lambda (n)g(n)$$

其中$f(n)$是需要求极值（一般是极小值）的代价函数，
$g(n)$是约束条件，$\lambda (n)$是拉格朗日乘子。在LMS算法里，我们定义的拉格朗日函数为：

$$L(n)=\frac{1}{2}{\Vert \mathbf{w}(n+1)-\mathbf{w}(n)\Vert }^2+\lambda (n)\left[r(n)-\left(1-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)\right]$$

其中$\lambda (n)$ 是时变的拉格朗日乘子。

此时将$L(n)$ 对${\mathbf{w}}^H(n+1)$求导数，得到：
$$\frac{\partial L(n)}{\partial {\mathbf{w}}^H(n+1)}=\frac{\partial \left\{\frac{1}{2}{\Vert \mathbf{w}(n+1)-\mathbf{w}(n)\Vert }^2+\lambda (n)\left[r(n)-\left(1-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)\right]\right\}}{\partial {\mathbf{w}}^H(n+1)}=\frac{\partial \left[\frac{1}{2}{\left(\mathbf{w}(n+1)-\mathbf{w}(n)\right)}^H\left(\mathbf{w}(n+1)-\mathbf{w}(n)\right)+\lambda (n)\left(d(n)-{\mathbf{w}}^H(n+1)\mathbf{u}(n)\right)\right]}{\partial {\mathbf{w}}^H(n+1)}=\frac{\partial \left[\frac{1}{2}\left(\mathbf{w}(n+1{)}^H-\mathbf{w}(n{)}^H\right)\left(\mathbf{w}(n+1)-\mathbf{w}(n)\right)+\lambda (n)\left(d(n)-{\mathbf{w}}^H(n+1)\mathbf{u}(n)\right)\right]}{\partial {\mathbf{w}}^H(n+1)}=\frac{\partial \left[\frac{1}{2}\left({\mathbf{w}}^H(n+1)-{\mathbf{w}}^H(n)\right)\left(\mathbf{w}(n+1)-\mathbf{w}(n)\right)\right]}{\partial {\mathbf{w}}^H(n+1)}+\frac{\partial \left[\lambda (n)\left(d(n)-{\mathbf{w}}^H(n+1)\mathbf{u}(n)\right)\right]}{\partial {\mathbf{w}}^H(n+1)}=\frac{\partial \left[\frac{1}{2}\left({\mathbf{w}}^H(n+1)\mathbf{w}(n+1)-{\mathbf{w}}^H(n+1)\mathbf{w}(n)-{\mathbf{w}}^H(n)\mathbf{w}(n+1)+{\mathbf{w}}^H(n)\mathbf{w}(n)\right)\right]}{\partial {\mathbf{w}}^H(n+1)}-\frac{\partial \left[\lambda (n){\mathbf{w}}^H(n+1)\mathbf{u}(n)\right]}{\partial {\mathbf{w}}^H(n+1)}=\frac{\frac{1}{2}\partial {\mathbf{w}}^H(n+1)}{\partial {\mathbf{w}}^H(n+1)}\mathbf{w}(n+1)+\frac{\frac{1}{2}\partial \mathbf{w}(n+1)}{\partial {\mathbf{w}}^H(n+1)}{\mathbf{w}}^H(n+1)-\frac{\frac{1}{2}\partial {\mathbf{w}}^H(n+1)}{\partial {\mathbf{w}}^H(n+1)}\mathbf{w}(n)-\frac{\frac{1}{2}\partial \mathbf{w}(n+1)}{\partial {\mathbf{w}}^H(n+1)}{\mathbf{w}}^H(n)-\lambda (n)\mathbf{u}(n)$$

应用了$\frac{\partial \mathbf{f}}{\partial {\mathbf{w}}^H}={\left(\frac{\partial \mathbf{f}}{\partial \mathbf{w}}\right)}^H$后，可以进一步化简为下式
$$\frac{\partial L(n)}{\partial {\mathbf{w}}^H(n+1)}=\mathbf{w}(n+1)-\mathbf{w}(n)-\lambda (n)\mathbf{u}(n)$$

为使上式为0，可得
$$\mathbf{w}(n+1)=\mathbf{w}(n)+\lambda (n)\mathbf{u}(n)$$

将此代入拉格朗日方程使其为0，求解$\lambda (n)$ 的表达式。
$$L(n)=\frac{1}{2}{\Vert \mathbf{w}(n+1)-\mathbf{w}(n)\Vert }^2+\lambda (n)\left[r(n)-\left(1-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)\right]=\frac{1}{2}{\Vert \lambda (n)\mathbf{u}(n)\Vert }^2+\lambda (n)\left[\left(1-\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)-\left(1-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^2\right)e(n)\right]=\frac{1}{2}{\Vert \lambda (n)\mathbf{u}(n)\Vert }^2+\lambda (n)\left[-\frac{1}{2}\mu {\Vert \mathbf{u}(n)\Vert }^{2}e(n)\right]=\frac{1}{2}{\lambda }^2(n){\Vert \mathbf{u}(n)\Vert }^2-\frac{1}{2}\lambda (n)\mu {\Vert \mathbf{u}(n)\Vert }^{2}e(n)=0$$

故可得：
$$\lambda (n)=\mathbf{u}(n)e^{\ast }(n)$$

所以我们再次得到了系数更新规则：

$$\mathbf{w}(n+1)=\mathbf{w}(n)+\mu \mathbf{u}(n)e^{\ast }(n)$$

该式和上一节的结论是一致的。因此我们用拉格朗日乘子法得到了同样的最优系数更新规则，这也证明了LMS法在局部是最优的，或者说是约束最优化。

全局次优性

为维纳滤波器相比，LMS算法的相应学习曲线接近一个最终值$J(\infty )$，$J(\infty )$和维纳解下的最小代价$J_{\min }$存在一个额外的均方误差，因此和维纳解相比，LMS是次优的：

1. 在最速下降法中，在计算学习曲线之前要进行集平均，所以最速下降法的学习曲线在性质上是确定性的。

2. LMS法没有上述的集平均，所以梯度会存在"噪声"，因此LMS的学习曲线在性质上是统计的。

梯度自适应格型算法

其实现结构如下图所示：

因为目前应用中该算法结构没有被使用，这里不再展开讨论。

其实随机梯度下降算法还有很多其他的应用，这里不再展开。