矩阵分析与应用-1.9-Moore-Penrose逆矩阵-Section3

前言

本文学习过程来源是《矩阵分析与应用-张贤达》一书. 可以通过 z-lib 下载.

一、非一致方程的最小范数最小二乘解

我们在之前广义逆矩阵一节中讲到了一致方程的最小范数解和非一致方程的最小二乘解. 但是我们要注意, 最小二乘解非唯一, 我们要在其中找到一个范数最小的解, 我们将其称为非一致方程的最小范数最小二乘解 ( minimun norm least squares solution ), 我们也称其为半范数 ( seminorm ) 最小二乘解

定义 1: 对于非一致方程 \(A_{m \times n} x_{n \times 1} = y_{m \times 1}\) , 如果矩阵 \(G\) 满足条件:

\[ \lVert Gy \rVert _n \leq \lVert \hat{x} \rVert _n ,其中 \forall \hat{x} \in \lbrace \hat{x}: \lVert A\hat{x}-y \rVert _m \leq \lVert Az-y \rVert _m \quad \forall y \in R^m, z \in R \rbrace \tag{1} \]

式子中, \(\lVert \cdot \rVert_n\) 和 \(\lVert \cdot \rVert_m\) 分别是在 \(R^n\) 和 \(R^m\) 空间的范数 (半范数); 花括号 \(\{ \cdot \}\) 表示 \(\hat{x}\) 是非一直方程 \(Ax=y\) 的最小二乘解, 而 \(\lVert Gy \rVert_n \leq \lVert \hat{x} \rVert_n\) 表示 \(Gy\) 是在所有的最小二乘解中具有最小范数的那个解.

定理 1: 矩阵 \(G\) 使得 \(Gy\) 是非一致方程 \(Ax=y\) 的最小范数最小二乘解, 当且仅当 \(G\) 满足条件:

\[ AGA = A, \quad (AG)^{\sharp} = AG, \quad GAG = G, \quad (GA)^{\sharp} = GA \tag{2} \]

式子中, \(A^{\sharp}\) 是 \(A\) 的伴随矩阵.

该定理也可以等价表述为: 矩阵 \(G\) 使得 \(Gy\) 是非一致方程 \(Ax = y\) 的最小范数最小二乘解, 当且仅当 \(G\) 是 \(A\) 的 \(\mathrm{Moore-Penrose}\) 逆矩阵.

二、广义逆矩阵的阶数递推计算

系统辨识? 想必又是信号处理那方面的内容. 确实看起来有点打脑壳, 也不知道它的目的是什么. 就直接提出左、右伪逆矩阵的阶数递推计算方法.

1. 左伪逆矩阵的阶数递推

考虑 \(n \times m\) 矩阵 \(F_m\), 并设 \(F^{\dagger}=(F_m^{\mathrm{H}}F_m)^{-1}F_m^{\mathrm{H}}\) 是 \(F_m\) 的左伪逆矩阵.

定理 2: 令

\[ F_m = \begin{bmatrix} F_{m-1}& f_m \end{bmatrix} \tag{3} \]

式子中, \(f_m\) 是矩阵 \(F_m\) 的第 \(m\) 列, 且 \(\mathrm{rank}(F_m)=m\), 则计算 \(F_m^{\dagger}\) 的递推公式由下可得:

\[ F_m^{\dagger} = \begin{bmatrix} F_{m-1}^{\dagger} - F_{m-1}^{\dagger}f_me_m^{\mathrm{H}}\Delta_m^{-1}\\ e_m^{\mathrm{H}}\Delta_m^{-1} \end{bmatrix} \tag{4} \]

式子中,

\[ e_m = [I_n - F_{m-1}F_{m-1}^{\dagger}] f_m \tag{5} \]

\[ \Delta_m^{-1} = [f_m^{\mathrm{H}}e_m]^{-1} \tag{6} \]

且初始值为

\[ F_1^{\dagger} = f_1^{\mathrm{H}}/(f_1^{\mathrm{H}}f_1) \tag{7} \]

2. 右伪逆矩阵的阶数递推

考虑矩阵 \(F_m \in C^{n \times m}, \ n < m\).

定理 3: 记 \(F_m = [F_{m-1}, f_m]\), 则右伪逆矩阵 \(F_m^{\dagger}=F_m^{\mathrm{H}}(F_mF_m^{\mathrm{H}})^{-1}\) 具有以下递推公式:

\[ F_m^{\dagger} = \begin{bmatrix} F_{m-1}^{\dagger} - \Delta_mF_{m-1}^{\dagger}f_mc_m\\ \Delta_m c_m^{\mathrm{H}} \end{bmatrix} \tag{8} \]

式中, \(c_m^{\mathrm{H}} = f_m^{\mathrm{H}}(I_n - F_{m-1}F_{m-1}^{\dagger}), \ \Delta_m = c_m^{\mathrm{H}}f_m\). 递推的初始值为 \(F_1^{\dagger} = f_1^{\mathrm{H}}/(f_1^{\mathrm{H}}f_1)\).

三、超定二维超越方程的求解

这部分就真的是完全看不懂. 从二维信号处理中的一个基本问题: 辨识二维实自回归 —— 移动平均 (\(\mathrm{ARMA}\)) 模型的移动平均 (\(\mathrm{MA}\)) 参数.

关于 \(\mathrm{ARMA}\) 常用作经济计量和工程预测. 本质上还是统计学, 和机器学习内核是一样的.

令 \(\mathrm{ARMA}(p_1,p_2;q_1,q_2)\) 模型为

\[ \sum_{i=0}^{p_1}\sum_{j=0}^{p_2}a_{ij}x(n_1 - i, n_2 - j) = \sum_{i=0}^{q_1}\sum_{j=0}^{q_2}b_{ij}e(n_1 - i, n_2 - j) \tag{9} \]

不失一般性, 假定 \(a_{00} = 1, q_1 \le p_1, q_2 \le p_2\), 且不可观测的实输入激励 \(e(n_1,n_2)\) 为二维正态白噪声, 其均值为零, 方差等于 \(\sigma^2\). 二维 \(\mathrm{ARMA}\) 过程 \(\{ x(n_1,x_2) \}\) 的功率谱密度由下式给出:

\[ P(\omega_1,\omega_2) = \frac{B(z_1,z_2)B(z_1^{-1},z_2^{-1})}{A(z_1,z_2)A(z_1^{-1},z_2^{-1})}\sigma^2 \bigg|_{z_1=e^{i\omega_1},z_2=e^{j\omega_2}} \tag{10} \]

式中

\[ A(z_1,z_2) = \sum_{i=0}^{p_1}\sum_{j=0}^{p_2}a_{ij}z_1^iz_2^j, \quad B=(z_1,z_2) \sum_{i=0}^{q_1}\sum_{j=0}^{q_2}b_{ij}z_1^iz_2^j \]

注意, 多项式 \(B(z_1,z_2)B(z_1^{-1},z_2^{-1})\) 中的 \(z_1^kz_2^m\) 和 \(z_1^{-k}z_2^{-m}\) 项的系数相同, \(z_1^kz_2^{-m}\) 和 \(z_1^{-k}z_2^m\) 的系数相同. 所以可以将式子 (10) 的分子写作

\[ B(z_1,z_2)B(z_1^{-1},z_2^{-1})\sigma^2=C(z_1,z_2) + D(z_1,z_2^{-1}) + D(z_1^{-1},z_2) + C(z_1^{-1},z_2^{-1}) \tag{11} \]

式中

\[ C(z_1,z_2) = \sum_{i=0}^{q_1}\sum_{j=0}^{q_2}c_{ij}z_1^iz_2^j, \quad D=(z_1,z_2^{-1}) \sum_{i=0}^{q_1}\sum_{j=0}^{q_2}d_{ij}z_1^iz_2^{-j} \]

比较式 (11) 两边相同幂次项的系数, 立即得

\[ c_{00}=\frac{\sigma^2}{2}\sum_{i=0}^{q_1}\sum_{j=0}^{q_2}b_{ij}^2 \tag{12} \]

\[ c_{km}=\sigma^2\sum_{i=0}^{q_1}\sum_{j=0}^{q_2}b_{ij}b_{i+k,j+m}, \quad k = 0,1,\dots,q_1, \ m = 0,1,\dots,q_2, (k,m) \neq (0,0) \tag{13} \]

\[ d_{km}=\sigma^2\sum_{i=0}^{q_1}\sum_{j=0}^{q_2}b_{ij}b_{i-k,j+m}, \quad k = 0,1,\dots,q_1, \ m = 1,2,\dots,q_2 \tag{14} \]

上述方程称为超定二维超越方程. 二维 \(\mathrm{ARMA}\) 模型的 \(\mathrm{MA}\) 参数辨识问题的提法是: 已知二维 \(\mathrm{MA}\) 谱参数 \(c_{km}\) 和 \(d_{km}\), 求解超定二维超越方程式 (12) ~ (13), 即得到二维 \(\mathrm{MA}\) 参数 \(b_{ij}\).

定义估计误差

\[ f_{0,0} = \sum_{i=0}^{q_1}\sum_{j=0}^{q_2}b_{ij}^2 - 2c_{0,0} \tag{15} \]

\[ f_{k,m} = \sum_{i=0}^{q_1-k}\sum_{j=0}^{q_2-m}b_{i,j}b_{i+k,j+m} - c_{k,m} \quad k = 0,1,\dots,q_1, \ m = 0,1,\dots,q_2, (k,m) \neq (0,0) \tag{16} \]

\[ f_{k,m} = \sum_{i=0}^{q_1-k}\sum_{j=0}^{q_2-m}b_{i,j}b_{i+k,j+m} - d_{k,m} \quad k = 1,2,\dots,q_1, \ m = 1,2,\dots,q_2, (k,m) \tag{17} \]

记

\[ \begin{aligned} f &= [f_{0,0} \ , \dots,f_{q_1,q_2} \ , f_{q_1+1,q_2+1} \ , \dots, f_{2q_1,2q_2}]^{\mathrm{T}} \\ b &= [b_{0,0} \ , b_{0,1} \, \dots, b_{0,q_2} \ \dots \ , b_{q_1,0} \ , b_{q_1,1} \dots, b_{q_1,q_2}]^{\mathrm{T}} \end{aligned} \]

则梯度矩阵为

\[ F = \frac{\partial f}{\partial b^{\mathrm{T}}}=\begin{bmatrix} \frac{\partial f_{0,0}}{\partial b_{0,0}}& \dots & \frac{\partial f_{0,0}}{\partial b_{0,q_2}} & \dots& \frac{\partial f_{0,0}}{\partial b_{q_1,0}}& \dots & \frac{\partial f_{0,0}}{\partial b_{q_1,q_2}}\\ \vdots& & \vdots& & \vdots& & \vdots\\ \frac{\partial f_{q_1,q_2}}{\partial b_{0,0}}& \dots & \frac{\partial f_{q_1,q_2}}{\partial b_{0,q_2}} & \dots& \frac{\partial f_{q_1,q_2}}{\partial b_{q_1,0}}& \dots & \frac{\partial f_{q_1,q_2}}{\partial b_{q_1,q_2}}\\ \frac{\partial f_{q_1+1,q_2+1}}{\partial b_{0,0}}& \dots & \frac{\partial f_{q_1+1,q_2+1}}{\partial b_{0,q_2}} & \dots& \frac{\partial f_{q_1+1,q_2+1}}{\partial b_{q_1,0}}& \dots & \frac{\partial f_{q_1+1,q_2+1}}{\partial b_{q_1,q_2}}\\ \vdots& & \vdots& & \vdots& & \vdots\\ \frac{\partial f_{2q_1,2q_2}}{\partial b_{0,0}}& \dots & \frac{\partial f_{2q_1,2q_2}}{\partial b_{0,q_2}} & \dots& \frac{\partial f_{2q_1,2q_2}}{\partial b_{q_1,0}}& \dots & \frac{\partial f_{2q_1,2q_2}}{\partial b_{q_1,q_2}}\\ \end{bmatrix} \tag{18} \]

定义矩阵

\[ B_i = \begin{bmatrix} b_{i,0}& \dots& b_{i,q_2-1}& b_{i,q_2}\\ b_{i,1}& \dots& b_{i,q_2}& 0\\ \vdots& &\vdots & \vdots\\ b_{i,q_2}& \dots& 0& 0 \end{bmatrix} \tag{19} \]

\[ \bar{B_i} = \begin{bmatrix} b_{i,0}& b_{i,1}& \dots& b_{i,q_2}\\ 0& b_{i,0}& \dots& b_{i,q_2-1}\\ \vdots& &\vdots & \vdots\\ 0& \dots& 0& b_{i,0} \end{bmatrix} \tag{20} \]

和

\[ A_i = \begin{bmatrix} 0& b_{i,0}& b_{i,1}& \dots & b_{i,q_2-1}\\ 0& 0& b_{i,0}& \dots & b_{i,q_2-2}\\ \vdots& &\vdots & \vdots & \vdots\\ 0& \dots& 0& 0 & b_{i,0} \end{bmatrix} \tag{21} \]

\[ \bar{A_i} = \begin{bmatrix} b_{i,1}& \dots& b_{i,q_2-1}& b_{i,q_2} & 0\\ b_{i,2}& \dots& b_{i,q_2}& 0 & 0\\ \vdots& &\vdots & \vdots & \vdots\\ b_{i,q_2}& 0& \dots& 0 & 0 \end{bmatrix} \tag{22} \]

矩阵 \(F\) 具有下列形式:

\[ F=\begin{bmatrix} B_0& \dots& B_{q_1-1}& B_{q_1}\\ B_1& \dots& B_{q_1}& O\\ \vdots& & \vdots& \vdots\\ B_{q_1}& \dots& O& O\\ A_1& \dots& A_{q_1}& O\\ A_2& \dots& O& O\\ \vdots& & \vdots& \vdots\\ A_{q_1}& \dots& O & O \end{bmatrix} + \begin{bmatrix} \bar{B_0}& \bar{B_1}& \dots& \bar{B_{q_1}}\\ O& \bar{B_0}& \dots& \bar{B_{q_1-1}}\\ \vdots& \vdots& & \vdots\\ O& O& \dots& \bar{B_0}\\ O& \bar{A_0}& \dots& \bar{A_{q_1-1}}\\ O& O& \dots& \bar{A_{q_1-2}}\\ \vdots& \vdots& & \vdots\\ O& O& O & \bar{A_0} \end{bmatrix} \tag{23} \]

我们就能求解, 若第 \(i\) 步迭代估计的 \(\mathrm{MA}\) 参数向量记为 \(b^{(i)}\), 则第 \((i+1)\) 步的估计 \(b^{(i+1)}\) 由下式给出:

\[ b^{(i+1)} = b^{(i)} - F^{\dagger}f \tag{24} \]

式子中, \(F^{\dagger}\) 是矩阵 \(F\) 的广义逆矩阵.

算法 1: (二维 \(\mathrm{MA}\) 参数估计)

• Step1 : 取初始值 \(b_{0,0} = \sqrt{2c_{0,0}}\) 和 \(b_{i,j}=0, \ (i,j) \neq (0,0)\)

• Step2 : 利用式子 (15)~(17) 计算拟合误差 \(f_{k,m}\), 然后利用式子 (23) 构造矩阵 \(F\)

• Step3 : 根据式子 (24) 更新二维 \(\mathrm{MA}\) 参数向量估计 \(b^{(i+1)}\)

• Step4 : 返回 Step2, 并重复以上步骤, 知道所有 \(\mathrm{MA}\) 参数估计值 \(b_{k,m}\) 都收敛为止. 例如, 当

  \[ \bigg| \frac{b_{km}^{(i+1)}-b_{km}^{(i)}}{b_{km}^{(i+1)}} \le 0.05 \bigg| \tag{25} \]

  对所有 \(k = 0,1,\dots,q_1, m=0,1,\dots,q_2\) 成立时, 即认可参数估计已经收敛.