【现代信号处理】 11 -递归最小二乘

递归最小二乘 Recursive Least Square(RLS)

1. 问题引入

1.1 LMS自适应滤波回顾

  首先我们来回顾一下自适应问题。

$$\begin{gathered} Z(1),Z(2),...,Z(n) \\ d(1),d(2),...,d(n) \end{gathered}$$

  我们在不同时刻有不同的数据，同时也有不同的目标，我们希望通过某种线性变换，使得我们的数据能够去逼近我们的目标。

$$g(Z(n))\to d(n)\Rightarrow g(Z(n))={\omega }^{T}Z(n)$$

  然后我们会去寻找一个损失函数，去评估我们的估计情况，我们用均方误差作为评价标准，这样得到的结果其实就是维纳滤波。因为我们的系统是时变的，所以估计的线性系数ω也应该是个与时间有关的数字。

$$Loss:minE|d(n)-{\omega }^{T}Z(n){|}^2\Rightarrow minE|d(n)-\omega (n{)}^{T}Z(n){|}^2\text{ Adaptive}$$

  但是在实际求解的时候，我们会遇到一个问题，就是这个系统具有两种属性。首先，这个系统是时变的，每一时刻的估计都需要进行更改；同时这个系统的目标优化必须通过迭代的方法不断去寻找，这样就出现了迭代时间轴和采样时间轴的两套不同的下标体系。

  基于这个问题，我们有了自适应的思想，我们把两种时间系统进行统一，我们每次采样一次，就进行一次迭代，这样看起来似乎是种四不像的信号处理方法，但是实际上却有着不错的性能。

  自适应的继续发展就是神经网络。人们对神经网络为什么可以工作知之甚少。比人们对单层自适应算法了解的少很多，但是神经网络已经得到广泛应用了。

  在这种自适应的指导思想下，我们就有了LMS。因为每次都只采样一次，所以这个式子是没有办法求期望的，因此首先就去掉了期望符号

$$mi{n}_{\omega }E|d(n)-{\omega }^{T}Z(n){|}^2\Rightarrow mi{n}_{\omega }|d(n)-{\omega }^{T}Z(n){|}^2$$

$$\Rightarrow {\omega }_{n+1}={\omega }_n+\mu Z(n)\ast e(n)$$

  其中

$$e(n)=d(n)-{\omega }_n^{T}Z(n)$$

2. 递归最小二乘

2.1 LMS的缺陷与递归最小二乘模型建立

  今天，我们来重新审视一下LMS，我们觉得LMS是有缺陷的。

  LMS的缺陷在于，他只考虑了最新的数据点的信息，而没有考虑历史的信息。虽然马尔科夫性告诉我们，历史的可以遗忘的，但是马尔科夫性只是个假设。因此，如果我们想要利用历史的信息，我们应该怎么调整损失函数呢?

  事实上，我们可以考虑利用加法的形式，把各个时刻的损失函数加在一起

$$\sum _{k=1}^{n}E|d(k)-{\omega }_n^{T}Z(k){|}^2$$

  但是这么做，其实也有些极端了，我们应该更加重视最近的数据的影响，历史固然应该得到尊重，但也更应该与时俱进，因此我们引入遗忘因子的方式，继续修正损失函数

$$\sum _{k=1}^n\lambda (n,k)\ast E|d(k)-{\omega }_n^{T}Z(k){|}^2$$

$$\lambda (n,k)\Rightarrow \text{Forgetting Factor}$$

  最常见的遗忘因子是指数遗忘

$$\sum _{k=1}^n{\lambda }^{n-k}\ast E|d(k)-{\omega }_n^{T}Z(k){|}^2$$

  我们这里的期望仍然没法求，我们学习widrow的思想，去掉期望，得到了最终的损失函数

$$g({\omega }_n)=\sum _{k=1}^n{\lambda }^{n-k}\ast |d(k)-{\omega }_n^{T}Z(k){|}^2(1)$$

  我们对损失函数进行变形

$$\begin{gathered} g({\omega }_n)=\sum _{k=1}^n{\lambda }^{n-k}{d}^2(k)-2\sum _{k=1}^{n}d(k){\omega }_n^{T}Z(k){\lambda }^{n-k}+\sum _{k=1}^n{\lambda }^{n-k}({\omega }_n^{T}Z(k){)}^2 \\ =\sum _{k=1}^n{\lambda }^{n-k}{d}^2(k)-2{\omega }_n^T(\sum _{k=1}^n{\lambda }^{n-k}Z(k)d(k))+{\omega }_n^T(\sum _{k=1}^n{\lambda }^{n-k}Z(k)\ast Z^T(k)){\omega }_n(2) \end{gathered}$$

  这里注意一下，关于第一项，其实是个常数。

  对于第二项，因为ω^T是个行向量，Z是个列向量，二者相乘是个数字，所以可以提到前面去

  而第三项其实是个二次型

$$({\omega }_n^{T}Z(k){)}^2={\omega }_n^{T}Z(k)\ast Z^T(k){\omega }_n$$

  我们希望式子(2)能够更加简洁，因此我们令

$$\begin{gathered} Z(n)=\sum _{k=1}^n{\lambda }^{n-k}Z(k)d(k) \\ \varphi (n)=\sum _{k=1}^n{\lambda }^{n-k}Z(k)\ast Z^T(k) \end{gathered}$$

  这里我们需要记得，φ是个对称矩阵,求梯度的时候有用。

$${\varphi }^T(n)=\varphi (n)$$

  因此(2)式可以变为

$$g({\omega }_n)=C-2{\omega }_n^{T}Z(n)+{\omega }_n^T\varphi (n){\omega }_n(3)$$

  我们对我们的目标函数进行求梯度

$${\nabla }_{{\omega }_n}g({\omega }_n)=-2Z(n)+2\varphi (n){\omega }_n=0$$

$${\omega }_n={\varphi }^{-1}(n)Z(n)$$

  然后我们就得到了一个最小二乘的结果。

2.2 矩阵求逆公式

2.2.1 两个重要的中间式

  其实，我们更想得到一个自适应形式的递归最小二乘，这样的话便于求解。因为对于计算机来说，求逆是一件非常繁琐的事情。但是，首先我们要介绍一点知识，就是矩阵求逆公式。

  其实我们想推导一下这个矩阵求逆公式,我们目标是要求

$$(A+BCD{)}^{-1}=?$$

  这个矩阵求逆是个小技巧，但是这个小技巧对于递归最小二乘来说，具有本质性的作用。我们希望通过分块对角化的方法来进行求解。

  首先我们要通过对下面的矩阵进行一些操作，得到一些东西。

$$\left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)$$

  我们有这样的想法，因为对角阵求逆比较简单，所以我们希望把这个矩阵化成对角阵的形式。

  这个对角阵的转换思路有两种

• 一种是通过先消除B再消除D的方法得到对角阵
• 另一种是通过先消除D再消除B的手段得到对角阵。

  我们先来尝试第一种，我们通过左乘做行变换消除B，然后通过右乘做列变换消除D

$$\begin{gathered} \left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)\Rightarrow \left(\begin{array}{cc}I& BC\\ 0& I\end{array}\right)\ast \left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& 0\\ CD& I\end{array}\right) \\ =\left(\begin{array}{cc}A+BCD& 0\\ D& -C^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& 0\\ CD& I\end{array}\right)=\left(\begin{array}{cc}A+BCD& 0\\ 0& -C^{-1}\end{array}\right)(i) \end{gathered}$$

  我们发现，这么做完了以后，就有了A+BCD了，而且这个矩阵是个对角阵，求逆非常简单，但是前面的三个式子中，中间那一项不是对角阵，求逆还是很复杂，为了解决这个问题，我们就通过第二个处理思路，把中间矩阵通过第二种思路处理成对角阵，从而最终就能够得到A+BCD的逆了

  我们用第二种思路，先对这个矩阵进行左乘消掉D，再通过右乘消掉B，得到一个对角阵

$$\begin{gathered} \left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)\Rightarrow \left(\begin{array}{cc}I& 0\\ -D{A}^{-1}& I\end{array}\right)\ast \left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& -A^{-1}B\\ 0& I\end{array}\right) \\ =\left(\begin{array}{cc}A& B\\ 0& -D{A}^{-1}B-C^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& -A^{-1}B\\ 0& I\end{array}\right)=\left(\begin{array}{cc}A& 0\\ 0& -(D{A}^{-1}B+C^{-1})\end{array}\right)(ii) \end{gathered}$$

  于是我们得到了两个式子

$$\begin{gathered} \left(\begin{array}{cc}I& BC\\ 0& I\end{array}\right)\ast \left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& 0\\ CD& I\end{array}\right)=\left(\begin{array}{cc}A+BCD& 0\\ 0& -C^{-1}\end{array}\right)(i) \\ \left(\begin{array}{cc}I& 0\\ -D{A}^{-1}& I\end{array}\right)\ast \left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& -A^{-1}B\\ 0& I\end{array}\right)=\left(\begin{array}{cc}A& 0\\ 0& -(D{A}^{-1}B+C^{-1})\end{array}\right)(ii) \end{gathered}$$

---

2.2.2 求逆公式

  我们对得到的两个重要的中间式继续进行变换，就可以得到A+BCD的逆矩阵公式了

  我们对(1)取逆

$${\left(\begin{array}{cc}A+BCD& 0\\ 0& -C^{-1}\end{array}\right)}^{-1}={\left(\begin{array}{cc}I& 0\\ CD& I\end{array}\right)}^{-1}\ast {\left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)}^{-1}\ast {\left(\begin{array}{cc}I& BC\\ 0& I\end{array}\right)}^{-1}(iii)$$

  对上三角和下三角矩阵取逆，其实就是改变个符号

$${\left(\begin{array}{cc}I& 0\\ CD& I\end{array}\right)}^{-1}=\left(\begin{array}{cc}I& 0\\ -CD& I\end{array}\right){\left(\begin{array}{cc}I& BC\\ 0& I\end{array}\right)}^{-1}=\left(\begin{array}{cc}I& -BC\\ 0& I\end{array}\right)$$

  代入(iii)得

$${\left(\begin{array}{cc}A+BCD& 0\\ 0& -C^{-1}\end{array}\right)}^{-1}=\left(\begin{array}{cc}I& 0\\ -CD& I\end{array}\right)\ast {\left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)}^{-1}\ast \left(\begin{array}{cc}I& -BC\\ 0& I\end{array}\right)(iii)$$

  最难求的还是中间的那个逆。但是我们可以用式(ii)变形来得到

$$\left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)={\left(\begin{array}{cc}I& 0\\ -D{A}^{-1}& I\end{array}\right)}^{-1}\ast \left(\begin{array}{cc}A& 0\\ 0& -(D{A}^{-1}B+C^{-1})\end{array}\right)\ast {\left(\begin{array}{cc}I& -A^{-1}B\\ 0& I\end{array}\right)}^{-1}(iv)$$

  对(iv)取逆可得

$${\left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)}^{-1}=\left(\begin{array}{cc}I& -A^{-1}B\\ 0& I\end{array}\right)\ast \left(\begin{array}{cc}{A}^{-1}& 0\\ 0& -(D{A}^{-1}B+C^{-1}{)}^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& 0\\ -D{A}^{-1}& I\end{array}\right)(v)$$

  把(v)代入(iii)中,而且我们记住一点，我们只需要矩阵左上角的值，其他位置的值不重要

$$\begin{gathered} {\left(\begin{array}{cc}A+BCD& 0\\ 0& -C^{-1}\end{array}\right)}^{-1} \\ =\left(\begin{array}{cc}I& 0\\ -CD& I\end{array}\right)\ast {\left(\begin{array}{cc}A& B\\ D& -C^{-1}\end{array}\right)}^{-1}\ast \left(\begin{array}{cc}I& -BC\\ 0& I\end{array}\right) \\ =\left(\begin{array}{cc}I& 0\\ -CD& I\end{array}\right)\ast \left(\begin{array}{cc}I& -A^{-1}B\\ 0& I\end{array}\right)\ast \left(\begin{array}{cc}{A}^{-1}& 0\\ 0& -(D{A}^{-1}B+C^{-1}{)}^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& 0\\ -D{A}^{-1}& I\end{array}\right)\ast \left(\begin{array}{cc}I& -BC\\ 0& I\end{array}\right) \\ =\left(\begin{array}{cc}I& -A^{-1}B\\ -CD& \ast \end{array}\right)\ast \left(\begin{array}{cc}{A}^{-1}& 0\\ 0& -(D{A}^{-1}B+C^{-1}{)}^{-1}\end{array}\right)\ast \left(\begin{array}{cc}I& -BC\\ -D{A}^{-1}& \ast \end{array}\right) \\ =\left(\begin{array}{cc}{A}^{-1}& A^{-1}B(D{A}^{-1}B+C^{-1}{)}^{-1}\\ \ast & \ast \end{array}\right)\ast \left(\begin{array}{cc}I& -BC\\ -D{A}^{-1}& \ast \end{array}\right) \\ =\left(\begin{array}{cc}{A}^{-1}-A^{-1}B(D{A}^{-1}B+C^{-1}{)}^{-1}D{A}^{-1}& \ast \\ \ast & \ast \end{array}\right) \end{gathered}$$

  最终我们证明了矩阵求逆公式

$$\text{Matrix Inversion}(A+BCD{)}^{-1}=A^{-1}-A^{-1}B(D{A}^{-1}B+C^{-1}{)}^{-1}D{A}^{-1}$$

2.3 递归最小二乘的自适应形式

2.3.1 φ(n)的递推

  得到了矩阵求逆公式以后，我们就可以进行递归最小二乘自适应形式的推导了。

  我们希望能够基于n个数据得到估计的基础上，对这个估计稍作调整，就能得到n+1个数据产生的估计结果。其实这种做法，我们不论是在前后向滤波，还是卡尔曼滤波中，都是这么做的。因此，我们希望能够得到一个递推的形式，最终能简化运算。
$${\omega }_n={\varphi }^{-1}(n)Z(n)$$

  其中

$$\begin{gathered} Z(n)=\sum _{k=1}^n{\lambda }^{n-k}Z(k)d(k) \\ \varphi (n)=\sum _{k=1}^n{\lambda }^{n-k}Z(k)\ast Z^T(k) \end{gathered}$$

  为了能够对ω_n进行递推，我们首先要找到Z(n)和φ(n)的递推公式

$$\begin{gathered} \varphi (n)=\sum _{k=1}^{n-1}{\lambda }^{n-k}Z(k)\ast Z^T(k)+Z(n)Z^T(n) \\ =\lambda \ast \sum _{k=1}^{n-1}{\lambda }^{n-k-1}Z(k)\ast Z^T(k)+Z(n)Z^T(n) \\ =\lambda \varphi (n-1)+Z(n)Z^T(n) \end{gathered}$$

  

$${\varphi }^{-1}(n)=(\lambda \varphi (n-1)+Z(n)Z^T(n){)}^{-1}$$

  这个时候就用到了我们的矩阵求逆公式。

$$(A+BCD{)}^{-1}=A^{-1}-A^{-1}B(D{A}^{-1}B+C^{-1}{)}^{-1}D{A}^{-1}$$
  把φ^-1也表示成这个样子

$${\varphi }^{-1}(n)=(\lambda \varphi (n-1)+Z(n)\ast 1\ast Z^T(n){)}^{-1}$$

$$\{\begin{array}{l}A=\lambda \varphi (n-1)\\ B=Z(n)\\ C=1\\ D=Z^T(n)\end{array}$$

  因为DA^-1B是一个二次型的形式，其实最终乘积的结果是一个数字，因此C应该表示为1，而不是I，这样的话，使用矩阵求逆公式可以得到

$${\varphi }^{-1}(n)={\lambda }^{-1}{\varphi }^{-1}(n-1)-\frac{{\lambda }^{-1}{\varphi }^{-1}(n-1)Z(n)Z^T(n){\varphi }^{-1}(n-1){\lambda }^{-1}}{1+{\lambda }^{-1}{Z}^T(n){\varphi }^{-1}(n-1)Z(n)}$$

  我们发现，得到的结果中，虽然有很多逆，但是全部都是标量的逆，或者φ(n-1)^-1，而φ(n-1)^-1是递推的前一项，我们默认就可以认为这个是个已知的东西了。

  但是这个式子看起来还是非常的麻烦，我们希望能够让他更简洁，我们令

$$K(n)=\frac{{\lambda }^{-1}{\varphi }^{-1}(n-1)Z(n)}{1+{\lambda }^{-1}{Z}^T(n){\varphi }^{-1}(n-1)Z(n)}$$

$${\varphi }^{-1}(n)={\lambda }^{-1}{\varphi }^{-1}(n-1)-{\lambda }^{-1}K(n)Z^T(n){\varphi }^{-1}(n-1)[1]$$

  我们继续研究下K(n)
$$K(n)=\frac{{\lambda }^{-1}{\varphi }^{-1}(n-1)Z(n)}{1+{\lambda }^{-1}{Z}^T(n){\varphi }^{-1}(n-1)Z(n)}$$

  把分母乘到左边

$$\begin{gathered} (1+{\lambda }^{-1}{Z}^T(n){\varphi }^{-1}(n-1)Z(n))\ast K(n)={\lambda }^{-1}{\varphi }^{-1}(n-1)Z(n) \\ K(n)={\lambda }^{-1}{\varphi }^{-1}(n-1)Z(n)-{\lambda }^{-1}K(n)Z^T(n){\varphi }^{-1}(n-1)Z(n) \\ =({\lambda }^{-1}{\varphi }^{-1}(n-1)-{\lambda }^{-1}K(n)Z^T(n){\varphi }^{-1}(n-1))\ast Z(n) \end{gathered}$$

  我们发现前一项就是φ^-1(n)

$$K(n)={\varphi }^{-1}(n)\ast Z(n)[2]$$

2.3.2 Z(n)的递推

  我们再来找找Z(n)的递推公式

$$\begin{gathered} Z(n)=\sum _{k=1}^n{\lambda }^{n-k}Z(k)d(k)=\sum _{k=1}^{n-1}{\lambda }^{n-k}Z(k)d(k)+Z(n)d(n) \\ =\lambda \ast \sum _{k=1}^{n-1}{\lambda }^{n-k-1}Z(k)d(k)+Z(n)d(n)=\lambda Z(n-1)+Z(n)d(n) \end{gathered}$$

2.3.3 ω_n的递推公式与Kalman Filter

$$\begin{gathered} {\omega }_n={\varphi }^{-1}(n)Z(n) \\ =({\lambda }^{-1}{\varphi }^{-1}(n-1)-{\lambda }^{-1}K(n)Z^T(n){\varphi }^{-1}(n-1))\ast (\lambda Z(n-1)+Z(n)d(n)) \\ ={\varphi }^{-1}(n-1)\ast Z(n-1)-K(n)Z^T(n){\varphi }^{-1}(n-1)\ast Z(n-1)+{\varphi }^{-1}(n)\ast Z(n)d(n)) \\ ={\omega }_{n-1}-K(n)Z^T(n)\omega (n-1)+K(n)d(n) \\ ={\omega }_{n-1}+K(n)(d(n)-{\omega }^T(n-1)Z(n)) \end{gathered}$$

  我们发现得到的这个结果就是卡尔曼滤波。k就是卡尔曼增益

  我们现在对卡尔曼滤波有了新认识，卡尔曼滤波就是在做优化。卡尔曼滤波其实就是把跟踪的历史，全部拿来做最小二乘优化。这就是递归最小二乘

3. SVD 奇异值分解

3.1 数学表示

  下面，我们希望通过几何的方式来理解，什么是奇异值分解。

  首先来给奇异值分解做个数学上的定义。

  我们希望分解一个矩阵A，假设有n行m列，并且我们可以假设n>m(这里只是为证明而假设一个条件，其实n和m谁大谁小不影响，任何矩阵都可以做SVD)。

$$A=\left(\begin{array}{c}{a}_1^T\\ ...\\ a_n^T\end{array}\right)\in R^{n\ast m}\text{ (n>m)}$$

  我们知道A可以表示为这样的形式

$$A=U\Sigma V^T$$

3.2 证明

  我们从行空间的角度来分析，其实A就是m维空间的n个点。我们希望能够做一个优化，我们希望能够找到一个超平面，使得A中的点，到这个超平面的距离的平方和最小。

  这听起来像是在做最小二乘，但是其实并不是，因为最小二乘是把点对应到回归方程上，然后计算平方和，计算的并不是点到直线的距离。

  最小二乘的目标函数
$$min|y-cx{|}^2$$

  点到超平面的距离，我们在SVM中证明过

$$\begin{gathered} d(a,v)=\frac{a^{T}v}{||v||} \\ {d}^2(a,v)=\frac{|a^{T}v{|}^2}{||v|{|}^2} \end{gathered}$$

  我们假设超平面的法向量v是个单位向量

$$||v||=1\Rightarrow d^2(a,v)=|a^{T}v{|}^2$$

  因此，我们的目标函数是这样的

$$\sum _{k=1}^n{d}^2(a_k,v)=\sum _{k=1}^n|a_k^{T}v{|}^2=|Av{|}_2^2$$

  实际上结果就是Av的二范数的平方

$$Av=\left(\begin{array}{c}{a}_1^T\ast v\\ ...\\ a_n^T\ast v\end{array}\right)$$

$$|Av{|}_2=(|a_1^{T}v{|}^2+...+|a_n^{T}v{|}^2{)}^{\frac{1}{2}}$$

  然后要求这个超平面的法向量，就变成了一个条件极值问题

$$min|Av{|}_2^2\text{ s.t. ||v||=1}$$

  我们假设得到了法向量v₁

$$v_1=argmi{n}_{||v||=1}||Av|{|}_2^2$$

  这里，我们只往一个平面上进行投影，我们想继续找一个平面，让a到v₁和v₂构成的子空间的距离和最短

d 2 ( a , { v 1 , v 2 } ) d^2(a,\{v_1,v_2\}) d2(a,{v1​,v2​})

  因为我们已经知道a到v₁的平面的距离了，我们希望这个到子空间的距离能够有这样的表示

d ( a , { v 1 , v 2 } ) = d ( a , v 1 ) + d ( a , v 2 ) d(a,\{v_1,v_2\}) = d(a,v_1)+d(a,v_2) d(a,{v1​,v2​})=d(a,v1​)+d(a,v2​)

  只要v₁和v₂是正交的，我们就能够得到这样的式子。因此，我们实际上就是要在v₁超平面的正交补空间里面再找一个单位法向量。

$$v_2=argmi{n}_{v\perp v_1,||V||=1}||Av|{|}_2^2$$
v 3 = a r g m i n v ⊥ { v 1 , v 2 } , ∣ ∣ V ∣ ∣ = 1 ∣ ∣ A v ∣ ∣ 2 2 v_3 = argmin_{v \perp \{v_1,v_2\},||V||=1} ||Av||_2^2 v3​=argminv⊥{v1​,v2​},∣∣V∣∣=1​∣∣Av∣∣22​

  这样的向量我们一共能够找到r个，其中r是矩阵A的秩，因为A的空间维度就是r，不可能再找出更多维度的向量了

$$rank(A)=r$$
$$\Rightarrow (v_1,...v_r)$$

  实际上，因为a_k向量是m维度的，因此空间的单位正交基应该有m个，我们需要对前面的r个向量进行扩充。前面r个向量是正交向量，是对n个给定点进行最优逼近的子空间。而m个向量是A空间的正交基

$$\begin{gathered} (v_1,...,v_r)\Rightarrow (v_1,v_2,...,v_r,...v_m) \\ Orthogonal\Rightarrow \text{Orthogonal basis} \end{gathered}$$

  我们假设有向量u_k满足这样的关系,我们假设u_k是经过了归一化的

$$u_k=\frac{A{v}_k}{||A{v}_k||},k=1,2...,n$$
  令
$${\sigma }_k=||A{v}_k||$$
  则

$$A{v}_k={\sigma }_k{u}_k$$

$$A(v_1,...,v_m)=\left(\begin{array}{ccc}{u}_1& ...& u_m\end{array}\right)\ast \left(\begin{array}{ccc}{\sigma }_1& & \\ & ...& \\ & & {\sigma }_m\end{array}\right)$$

  因为v矩阵是正交矩阵，引起求逆就是求转置

$$A=\left(\begin{array}{ccc}{u}_1& ...& u_m\end{array}\right)\ast \left(\begin{array}{ccc}{\sigma }_1& & \\ & ...& \\ & & {\sigma }_m\end{array}\right)\ast \left(\begin{array}{c}{v}_1^T\\ ...\\ v_m\end{array}\right)$$

  我们得到的这个式子其实就已经是奇异值分解了。

  现在我们知道u是单位向量，但是还不确定是否是正交的。A*A^T之后，得到了对角分解的形式，能够说明u是正交矩阵

$$A{A}^T=\left(\begin{array}{ccc}{u}_1& ...& u_m\end{array}\right)\ast \left(\begin{array}{ccc}{\sigma }_1^2& & \\ & ...& \\ & & {\sigma }_m^2\end{array}\right)\ast \left(\begin{array}{c}{u}_1^T\\ ...\\ u_m\end{array}\right)$$

  按理来说，AA^T是个n*n的矩阵，但是他一定是个缺秩的，我们可以按照特征分解的方式进行，也就是对前面进行补项，通过在u的正交补空间(也就是零空间)中寻找向量，然后再做施密特正交化，就可以得到补满以后的u，相应位置的特征值是0

  而u₁到u_m是A*A^T对应的非零特征值的特征向量

$$A{A}^T=\left(\begin{array}{ccccc}{u}_1& ...& u_m& ...& u_n\end{array}\right)\ast \left(\begin{array}{ccccccc}{\sigma }_1^2& & & & & & \\ & ...& & & & & \\ & & {\sigma }_m^2& & & & \\ & & & & 0& & \\ & & & & & ...& \\ & & & & & & 0\end{array}\right)\ast \left(\begin{array}{c}{u}_1^T\\ ...\\ u_m\\ ...\\ u_n\end{array}\right)$$

  实际上，我们可以看出，奇异值分解实际上的得到扩充的矩阵。

$$\begin{gathered} A=U\Sigma V^T \\ U\in R^{n\ast n},U{U}^T=I \\ V\in R^{m\ast m},V{V}^T=I \end{gathered}$$

4. 奇异值分解与最小二乘

4.1 公式推导

  我们接下来要探究奇异值分解和最小二乘之间的关系。

  约束方程

$$mi{n}_{\omega }||\left(\begin{array}{c}d(1)\\ ...\\ d(n)\end{array}\right)\left(\begin{array}{c}{z}^T(1)\\ ...\\ z^T(n)\end{array}\right)\ast \omega |{|}_2^2\Rightarrow mi{n}_{\omega }||d-Z\omega |{|}_2^2=mi{n}_{\omega }(d-Z\omega {)}^T(d-Z\omega )$$

  如果我们要考虑遗忘因子的话，我们可以把约束方程写成这样的形式

$$\begin{gathered} mi{n}_{\omega }(d-Z\omega {)}^T\Lambda (d-Z\omega ) \\ \Lambda =diag({\lambda }^{n-1},{\lambda }^{n-2},...,1) \end{gathered}$$

  最小二乘解

$${\omega }_{LS}=(Z^{T}Z{)}^{-1}{Z}^{T}d$$

  我们知道，最小二乘是一个方程数量比未知数个数要多的超定方程，因此最小二乘是没有办法求解的，只能进行拟合。并且做最小二乘的条件是Z*Z^T必须是满秩的。如果想要做到满秩，Z本身必须是列满秩的。

  如果Z是列满秩的，我们对Z做奇异值分解有这样的式子

$$Z=U\Lambda V^T$$

  其中

$$\Lambda =\left(\begin{array}{ccc}{\sigma }_1& & \\ & ...& \\ & & {\sigma }_n\\ & 0& \end{array}\right)=\left(\begin{array}{c}\tilde{\Lambda }\\ 0\end{array}\right)$$

  因此，最小二乘解可以做这样的变形

$$\begin{gathered} {\omega }_{LS}=[(V{\Lambda }^T{U}^T)(U\Lambda V^T){]}^{-1}(V{\Lambda }^T{U}^T)d \\ =(V{\Lambda }^T\ast \Lambda V^T{)}^{-1}\ast (V{\Lambda }^T{U}^T)d(a) \end{gathered}$$

  其中

$${\Lambda }^T\Lambda =\left(\begin{array}{cc}\tilde{\Lambda }& 0\end{array}\right)\ast \left(\begin{array}{c}\tilde{\Lambda }\\ 0\end{array}\right)={\tilde{\Lambda }}^2(b)$$

  (b)代入(a)中可得

$$\begin{gathered} {\omega }_{LS}=V\ast \left(\begin{array}{ccc}{\sigma }_1^{-2}& & \\ & ...& \\ & & {\sigma }_n^{-2}\end{array}\right)V^T\ast (V\left(\begin{array}{cccc}{\sigma }_1& & & \\ & ...& & 0\\ & & {\sigma }_n& \end{array}\right)U^T)d \\ =V\left(\begin{array}{cccc}{\sigma }_1^{-1}& & & \\ & ...& & 0\\ & & {\sigma }_n^{-1}& \end{array}\right)\ast U^{T}d=Z^{+}d \end{gathered}$$

  前面一部分被称为伪逆，或者广义逆
$$\begin{gathered} \text{Pesudo Inversion} \\ \text{Gerenalized} \\ {Z}^{+}=V\left(\begin{array}{cccc}{\sigma }_1^{-1}& & & \\ & ...& & 0\\ & & {\sigma }_n^{-1}& \end{array}\right)\ast U^T \end{gathered}$$

  前面一部分如果不可逆，用符号伪逆来表示，也叫广义逆

$${\omega }_{LS}=Z^{+}d$$

4.2 分析

  现在我们对最小二乘有了新的认识，最小二乘本来是想要解方程，但是没有办法解，因为方程的数量比未知数要多。如果方程数量和未知数数量是一样的，我们就可以求解。当二者数量不一样的时候，就只能寻找二者最接近的ω

  现在，我们有了伪逆，我们可以利用伪逆，直接求得最小二乘的解。

$$Z\omega =d\Rightarrow \omega =Z^{+}d$$

$$Z=U\Lambda V^T\Rightarrow Z^{-1}=V{\lambda }^{-1}V$$

  我们看，对奇异值分解的直接求逆是不显示的，因为奇异值矩阵既不方，也可能不满秩。因此，我们处理的时候就先把他转置，然后能求逆的部分就求逆，不能求逆的部分就写0，得到的就是伪逆。伪逆的存在，让最小二乘本来是一个优化问题，现在变成了一个方程求解问题

4.3 伪逆的性质

  矩阵的逆具有这样的性质

$$Z^{-1}Z=Z{Z}^{-1}=I$$

  伪逆具有这样的性质

$$\begin{gathered} Z{Z}^{+}Z=Z \\ {Z}^{+}Z{Z}^{+}=Z^{+} \end{gathered}$$

  证明一下第一个式子

$$(U\Sigma V^T)(V{\Sigma }^{+}{U}^T)(U\Sigma V^T)=U\Sigma {\Sigma }^{+}\Sigma V^T$$

$$\Sigma =\left(\begin{array}{cc}\tilde{\Sigma }& 0\\ 0& 0\end{array}\right){\Sigma }^{+}=\left(\begin{array}{cc}{\tilde{\Sigma }}^{-1}& 0\\ 0& 0\end{array}\right)$$

$$\Rightarrow \Sigma {\Sigma }^{+}\Sigma =\left(\begin{array}{cc}\tilde{\Sigma }& 0\\ 0& 0\end{array}\right)=\Sigma$$

4.4 奇异值分解、伪逆与数值计算

  在计算机数值计算中，是不会采用求逆的方式计算最小二乘的，一般都是用奇异值分解，用伪逆，因为他的稳定性最好，性能最好。

  奇异值分解稳定性好的原因是，A^TA的条件数比A大很多，条件数多了，数值计算稳定性就比较差

$$Cond(A^{T}A)\ge Cond(A)$$

  矩阵求逆很怕矩阵的刚性大，就是大的特征值特别大，小的特征值特别小。但凡有点截断误差，逆就不容易求对，平方以后就更加脆弱了。

$$\left(\begin{array}{cc}1& 0\\ 0& 0.001\end{array}\right)\Rightarrow {\left(\begin{array}{cc}1& 0\\ 0& 0.001\end{array}\right)}^2=\left(\begin{array}{cc}1& 0\\ 0& 10^{-6}\end{array}\right)$$