Погрешность численного дифференцирования аппроксимацией МНК

[worm2]

Задана сеточная функция $y_i=f(x_i)$ на равномерной сетке (для простоты пусть сетка бесконечная в обе стороны, чтобы не возиться с краевыми эффектами).
Задана она с некоторой постоянной погрешностью $\delta$. Требуется её численно продифференцировать в узлах сетки через аппроксимацию прямой $y=Ax+B$, построенной стандартным методом наименьших квадратов по фиксированному (симметричному) числу точек слева и справа.

Ну стандартнейшая же задача, согласитесь? Все программы цифровой обработки данных так делают, и не стесняются этого.
А вы попробуйте найти оценки погрешности для этого метода! Я не нашёл. Я, правда, смотрел только популярные вузовские учебники по численным методам, но должно же было хоть что-то найтись!

Только здесь, на странице 83, есть такой пассаж:

Цитата:

Физики издавна употребляют (без строгого обоснования)...

Ах они, физики, такие-сякие! Употребляют, понимаешь, а строго не обосновывают!

Цитата:

...ещё один способ регуляризации — дифференцирование предварительно сглаженной кривой, причём сглаживание обычно выполняют методом наименьших квадратов.

Да, это оно, то, что мне надо. Ну, думаю, сейчас математики ка-а-ак покажут этим физикам мастер-класс, ка-а-ак строго обоснуют, заодно и я формулку себе урву. А вот фигушки! Дальше идёт вывод формулы приближённой производной и эмпирическое правило, как выбирать число точек в сглаживающем окне. Так я тоже могу! А формула-то оценки погрешности где? Нетути
Ну вот зачем они так, а? Это и есть вопрос для обсуждения.

Перехожу к разделу "самостоятельные попытки решения". Пришлось скрипеть застарелыми мозгами, гуглить-вспоминать, какие там в этом МНК формулы, выводить. С горем пополам что-то получилось, поэтому в раздел "Помогите решить/разобраться" вопрос не попал. Я подумал: не пропадать же добру, вдруг кому-нибудь ещё приспичит, а поди-ка найди справочник с этой формулой!

Итак, ограничиваемся формулировкой выше, простейший случай, сетка равномерная. Понятно, что не нарушая общности можно считать производную в нуле. Обозначим через $h$ шаг сетки, через $n$ число точек слева и справа, которые берём для МНК (всего с центральной точкой, в которой ищется производная, будет $2n+1$ точек). Абсциссы точек $x_i=ih$, $i=-n, \dots, -1, 0, 1, \dots, n$. Ординаты $f(x_i)=y_i$. В этом случае МНК для искомой производной даёт такую простую формулу:

$$f'(0) \approx A=\frac{\sum\limits_{i=-n}^n x_iy_i}{\sum\limits_{i=-n}^n x_i^2}=\frac{3\sum\limits_{i=-n}^n iy_i}{n(n+1)(2n+1)h}.$$

Сумма в числителе преобразуется в сумму разностей от равноудалённых точек, которые, в свою очередь, являются стандартными центральными аппроксимациями производной, взятыми с некими коэффициентами:

$$\sum\limits_{i=-n}^n iy_i=2h\sum\limits_{i=1}^n i^2\left(\frac{y_i-y_{-i}}{x_i-x_{-i}}\right).$$

Нас интересует погрешность суммы в скобках, которая состоит из двух компонентов:
1) Погрешность замены производной разностным аналогом, которая, к счастью, есть буквально везде: $M_3(hi)^2/6$, где $M_3=\max |f'''|$ (позвольте мне не уточнять, по какому интервалу берётся максимум).
2) Погрешность, связанная с погрешностью $\delta$ исходной функции. Она легко считается, это $\delta/(hi)$
Итого, общая погрешность вычисления производной:

$$\Delta A=|A-f'(0)|\leqslant \frac{\sum\limits_{i=1}^n \left(M_3h^3i^4+6\delta i\right)}{n(n+1)(2n+1)h}$$

Формулу для суммы четвёртых степеней тоже не пришлось самостоятельно выводить, гуглится легко:
$\sum\limits_{i=1}^n i^4=\frac{n(n+1)(2n+1)(3n^2+3n-1)}{30}$. Итого имеем:

$$\Delta A\leqslant \frac{M_3h^2(3n^2+3n-1)}{30}+\frac{3\delta}{(2n+1)h},$$

где первое слагаемое — погрешность, обусловленная разностной аппроксимацией, второе обусловлено погрешностью самой функции. Спасибо за внимание!

-- Ср апр 19, 2017 22:04:41 --

Возможно, более наглядно было бы привести выражение в терминах ширины окна $w=2nh$, там ещё можно было сделать грубую оценку и избавиться от $n$. Но это совсем простое упражнение.

-- Ср апр 19, 2017 22:20:49 --

Ну и стоит, наверное, добавить, что оценка погрешности — это одно, а настоящая погрешность может быть значительно меньше из-за взаимной компенсации ошибок и т.п.

[Евгений Машеров]

А тупо взять матрицу ковариаций ошибок коэффициентов
$C=\sigma^2(X^TX)^{-1}$
и вектор коэффициентов формулы дифференцирования b, выраженной через коэффициенты подгоночного полинома
и найти
$\sigma^2_D=\sigma^2b^T(X^TX)^{-1}b$?

[worm2]

Т.е. рассмотреть погрешности, как случайные величины, в простейшем случае независимые, нормально распределённые с $\sigma=\delta$ ?
Интересный подход. И, наверное, оптимальный с практической точки зрения. В учебниках по численным методам такому подходу, если и уделяют место, то мало: больше интересуются погрешностью в наихудшем случае.

Кажется, тогда второе слагаемое превращается в:

$$\frac{\delta}{h}\sqrt{\frac{3}{n(n+1)(2n+1)}},$$

а всё остальное остаётся как есть.

[Евгений Машеров]

Просто если мы используем МНК, то уже принимаем, что ошибки независимые нормально распределённые с одинаковой дисперсией.

[INGELRII]

Первое слагаемое можно грубо представить в виде $\alpha M_3 w^2$. Я правильно понимаю: выходит, что при фиксированном шаге сетки $h$ чем больше мы берем эту ширину $w$, тем потенциально выше будет погрешность?

[worm2]

Конечно. Тут та же ситуация, что и с простейшей формулой $f'(x) \approx (f(x+w)-f(x-w))/(2w)$: чем больше $w$, тем больше погрешность.