线性逆问题中的非正态统计范数求解方法

### 线性逆问题中的非正态统计范数求解方法 #### 1. 混合确定问题的求解方法 混合确定问题可以通过多种方法解决。类比L2方法，有两种主要途径： - **方法一**：选择一些先验模型参数并最小化E + L。这会导致一个线性规划问题，具有更多的变量（5M + 3N）、等式约束（2M + 2N）和不等式约束（5M + 2N）。目标是最小化： \[ z = \sum_{i = 1}^{M} \frac{\alpha_i}{\sigma_{m_i}} + \sum_{i = 1}^{N} \frac{\alpha'_i}{\sigma_{m_i}} \] 约束条件如下： \[ \begin{cases} m' - m'' + x - \alpha = \langle m \rangle \\ m' - m'' - x' + \alpha = \langle m \rangle \\ G(m' - m'') + x'' - \alpha' = d_{obs} \\ G(m' - m'') - x''' + \alpha' = d_{obs} \\ m' \geq 0, m'' \geq 0, \alpha \geq 0, \alpha' \geq 0 \\ x \geq 0, x' \geq 0, x'' \geq 0, x''' \geq 0 \end{cases} \] - **方法二**：首先使用奇异值分解来确定G的零空间。解的形式为： \[ m_{est} = \sum_{i = 1}^{p} a_i v^{(i)}_p + \sum_{i = p + 1}^{M} b_i v^{(i)}_0 = V_p a + V_0 b \] 其中，v是特征向量，a和b是未知系数向量。只有向量a会影响预测误差，因此使用超定算法确定a： \[ \text{找到最小化 } E = \left\lVert d_{obs} - Gm \right\rVert_1 = \sum_{i = 1}^{N} \frac{\left| d_{obs} - U_p \Lambda_p a \right|_i}{\sigma_{d_i}} \text{ 的 } a \] 其中奇异值分解为\(G = U_p \Lambda_p V_p^T\)。然后使用欠定算法确定b： \[ \text{找到最小化 } L = \left\lVert \langle m \rangle - m \right\rVert_1 = \sum_{i = 1}^{M} \frac{\left| [V_0 b]_i - (\langle m_i \rangle - \Lambda_p a) \right|}{\sigma_{m_i}} \text{ 的 } b \] 此外，还可以实现基本的欠定和超定L1算法，避免显式计算许多额外变量，从而大大减少存储和计算时间，使这些算法适用于解决中等规模（如M = 1000）的逆问题。 #### 2. 通过重加权L2最小化求解L1范数问题 之前通过将L1最小化问题转化为等价的线性规划问题来求解，但这种转化引入了许多新变量，增加了计算机内存需求。实际上，L1问题也可以转化为等价的L2问题，并使用标准最小二乘法求解。 首先考虑Ln范数，目标是使\(\left\lVert v \right\rVert_n^n\)看起来像加权L2范数\(\left\lVert W[v] \right\rVert_2^2\)的平方。其中： \[ \left\lVert v \right\rVert_n^n \equiv \sum_{k} |v_k|^n \] \[ \left\lVert W[v] \right\rVert_2^2 = v^T W v = \sum_{k} w_k |v_k|^2 \] 选择权重\(w_k = (|v_k|^{\gamma} + \delta^{\gamma})^{(n - 2) / \gamma}\)，其中\(0 < \delta \ll 1\)且\(1 \leq \gamma \leq 2\)。这样的选择使两个范数相等。为避免舍入误差，推荐使用更稳定的公式： \[ w_k = \begin{cases} \delta^{n - 2} \exp \left( \frac{n - 2}{\gamma} \log_{1p} \left( \frac{|v_k|}{\delta} \right)^{\gamma} \right) & \text{如果 } |v_k| \leq \delta \\ |v_k|^{n - 2} \exp \left( \frac{n - 2}{\gamma} \log_{1p} \left( \frac{\delta}{|v_k|} \right)^{\gamma} \right) & \text{如果 } |v_k| > \delta \end{cases} \] 对于L1范数（n = 1），权重矩阵是v的函数，这给最小化问题带来了困难。可以通过逐次逼近的方法解决： 1. 初始权重矩阵近似为\(W^{(0)} = I\)，使用标准最小二乘法确定模型参数的初始近似\(m^{(0)}\)。 2. 迭代更新权重矩阵\(W^{(j)} = W(m^{(j - 1)})\)，然后使用最小二乘法更新解\(m^{(j)}\)。 3. 当连续迭代之间的解不再变化时，终止重加权过程，此时\(m_{est} = m^{(j)}\)。 例如，对于找到单个模型参数\(m_1\)以最小化误差\(\left\lVert e \right\rVert_1\)的问题，其中\(e = d_{obs} - Gm\)，\(d_{obs} = Gm_{true} + n\)，\(n\)是指数分布噪声，\(G = [1, \cdots, 1]\)。L1问题的解是\(m_{1_{est}} = m_{1_{median}} = \text{median}(d_{obs})\)。通过迭代公式： \[ m^{(j)} = (G^T W^{(j - 1)} G)^{-1} G^T W^{(j - 1)} d_{obs} \] 迭代通常在约25次后收敛到预期解。 迭代算法还可用于解决具有先验信息\(Hm = h_{pri}\)的问题。目标是最小化： \[ \left\lVert e \right\rVert_2^2 + \mu \left\lVert h_{pri} - Hm \right\rVert_1 \] 解通过迭代公式计算： \[ m^{(j)} = (G^T G + \mu H^T W^{(j - 1)} H)^{-1} (G^T d_{obs} + \mu H^T W^{(j - 1)} h_{pri}) \] 在大多数情况下，重加权算法比线性规划算法更可取，主要例外是当存在不等式约束形式的先验信息时，线性规划算法能更自然地处理这些约束。 #### 3. 通过转化为线性规划问题求解L∞范数问题 L1范数对“坏”数据的权重小于L2范数，而L∞范数对其权重更大。目标是最小化： \[ \minimize E + L \equiv \left\lVert e \right\rVert_{\infty} + \left\lVert l \right\rVert_{\infty} = \max_i \frac{|e_i|}{\sigma_{d_i}} + \max_i \frac{|l_i|}{\sigma_{m_i}} \] L∞形式主要用于提供模型参数的“最坏情况”估计，用于与基于其他范数的估计进行比较。 一般线性方程\(Gm = d\)可以通过转化为线性规划问题在L∞意义下求解。 - **欠定问题**：引入新变量\(m_i'\)、\(m_i''\)、\(x_i\)和\(x_i'\)，每个长度为M，以及单个参数\(\alpha\)（共4M + 1个变量）。线性规划问题为： \[ \minimize \alpha \text{ ，约束条件为：} \begin{cases} G[m' - m''] = d_{obs} \\ m'_i - m''_i + x_i - \alpha \sigma_{m_i} = \langle m_i \rangle \\ m'_i - m''_i - x'_i + \alpha \sigma_{m_i} = \langle m_i \rangle \\ m' \geq 0, m'' \geq 0, \alpha \geq 0, x \geq 0, x' \geq 0 \end{cases} \] - **超定问题**： \[ \minimize \alpha \text{ ，约束条件为：} \begin{cases} \sum_{j = 1}^{M} G_{ij} [m'_j - m''_j] + x_i - \alpha \sigma_{d_i} = d_{obs_i} \\ \sum_{j = 1}^{M} G_{ij} [m'_j - m''_j] - x'_i + \alpha \sigma_{d_i} = d_{obs_i} \\ m'