我有模型$\ s(x,y)=x^2+y^2, 0 \leq x \leq 1, 0 \leq y \leq 1$.

我不是直接观察模型，而是观察模型的导数+一些噪声（e）：

$\ p(x,y)=s_x+e, q(x,y)=s_y+e$

根据 p(x,y 和 q(x,y) 的测量值，我想估计 s(x)。假设我知道 s(0,0)=0。

根据梯度定理： $\ s(x,y) = \int_{(0,0)}^{(x,y)} [s_x,s_y]d\vec{r}$

无论我们沿着哪条路径整合。

作为一个小实验（在 Matlab 中），我将正态分布噪声 N(0,1) 添加到 p=2x 和 q=2y。然后我首先沿 x 集成，然后沿 y 集成：SXY。接下来，我首先沿 y 集成，然后沿 x 集成：SYX。

结果表明梯度定理在这种情况下不成立（因为噪声）：

与模型相关的均方根误差为：

ErmsXY =
    0.1125
ErmsYX =
    0.0920

如何从 p 和 q 中找到更好的（更少的 RMS 误差和更平滑的）估计 s？

编辑：

从我读到的；使用曲线积分称为局部积分。还有一些全局积分方法，其中尝试选择最小化的 S(x,y)：

$\ \int_0^1\int_0^1 [|S_x - P|^2 + |S_y - Q|^2] dxdy$

当梯度嘈杂时，全局积分方法应该会给出更好的结果，但在实践中我该如何做到这一点？

编辑2：

我使用的一种方法是：

首先我们引入线性求导算子：$\ s_x = D_x*s,s_y = D_y*s$.

结果是以下线性方程组：

$\ D_x*s=p, D_y*s=q$

接下来找到这些方程的最小二乘误差解。这些方程的 LSE 解应该等效于从上面最小化积分。 这怎么能显示出来？

结果很好：

RMS 误差约为 SXY 和 SYX 的 1/5，求解也更平滑。

但是，这种方法有一些缺点：

1. 难以实施；必须使用中心差异并将二维矩阵“展平”为向量等。

2. 推导矩阵非常大且稀疏，因此它们可能会消耗大量 RAM。

另一种看起来可能更容易编码、更少 RAM 消耗和更快的方法是使用 FFT。在傅立叶空间中，这些 pdes 成为代数方程。这被称为 Frankot-Chellappa 算法，但不幸的是我还没有让它在我的示例数据上工作。