Bagaimana cara menghitung fungsi secara numerik dari gradien berisiknya?

Saya memiliki model .$\ s(x,y)=x^2+y^2, 0 \leq x \leq 1, 0 \leq y \leq 1$

Alih-alih mengamati model secara langsung, saya mengamati turunan dari model + noise (e):

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

Dari pengukuran p (x, y dan q (x, y) saya ingin memperkirakan s (x). Katakan saya tahu itu s (0,0) = 0.

Menurut teorema gradien: $\ s(x,y) = \int_{(0,0)}^{(x,y)} [s_x,s_y]d\vec{r}$

terlepas dari jalan mana kita mengintegrasikan.

Sebagai percobaan kecil (dalam Matlab) saya menambahkan noise terdistribusi normal, N (0,1), ke p = 2x dan q = 2y. Kemudian saya terintegrasi dulu bersama x diikuti oleh sepanjang y: SXY. Selanjutnya saya terintegrasi dulu sepanjang y diikuti oleh bersama x: SYX.

Hasil menunjukkan bahwa teorema gradien tidak berlaku dalam kasus ini (karena kebisingan):

Root Mean Square Error relatif terhadap model adalah:

ErmsXY =
    0.1125
ErmsYX =
    0.0920

Bagaimana saya bisa menemukan estimasi s dari p dan q yang lebih baik (lebih sedikit kesalahan RMS dan lebih halus)?

EDIT:

Dari apa yang saya baca; menggunakan integral kurva disebut sebagai integrasi lokal. Ada juga metode integrasi global di mana seseorang malah mencoba untuk memilih S (x, y) yang meminimalkan:

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

Metode integrasi global seharusnya memberikan hasil yang lebih baik ketika gradien berisik, tetapi bagaimana saya melakukan ini dalam praktik?

EDIT 2:

Salah satu pendekatan yang saya gunakan adalah ini:

pertama-tama kami memperkenalkan operator derivasi linier: .$\ s_x = D_x*s,s_y = D_y*s$

Hasilnya adalah sistem persamaan linear berikut:

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

Selanjutnya cari solusi Least Square Error untuk persamaan ini. Solusi LSE untuk persamaan ini seharusnya setara dengan meminimalkan integral dari atas. Bagaimana ini bisa ditampilkan?

Hasilnya bagus:

Kesalahan RMS adalah sekitar 1/5 dari SXY dan SYX dan solusinya juga lebih halus.

Namun ada beberapa kelemahan dari pendekatan ini:

1. sulit untuk diterapkan; harus menggunakan perbedaan pusat dan "meratakan" matriks 2D menjadi vektor dll.

2. Matriks derivasi sangat besar dan jarang sehingga mereka dapat mengkonsumsi banyak RAM.

Pendekatan lain yang tampaknya berpotensi lebih mudah untuk dikodekan, lebih sedikit memakan RAM dan lebih cepat adalah menggunakan FFT. Di ruang Fourier pdes ini menjadi persamaan aljabar. Ini dikenal sebagai algoritme Frankot-Chellappa, tetapi sayangnya saya belum membuatnya bekerja pada data contoh saya.

Anda dapat menyaring baik gradien itu sendiri atau hasilnya, . Anda perlu mengetahui karakteristik gradien sebenarnya dengan cukup baik untuk mengetahui apa bandwidth frekuensi. Pada titik itu Anda bisa mendesain filter low-pass yang akan mempertahankan sinyal tetapi menghilangkan noise frekuensi yang lebih tinggi.$s$