Сглаживающий сплайн

Сглаживающий сплайн (англ. smoothing spline) — оценка функции [math]\displaystyle{ \hat f(x) }[/math], полученная из набора зашумлённых наблюдений [math]\displaystyle{ y_i }[/math] за исходными данными [math]\displaystyle{ f(x_i) }[/math] и используемая в дальнейших вычислениях для балансировки адекватности модели функции [math]\displaystyle{ \hat f(x_i) }[/math] к [math]\displaystyle{ y_i }[/math] с основанной на производной мере кривизной функции [math]\displaystyle{ \hat f(x) }[/math]. Иными словами, сглаживающий сплайн является важным средством при работе с зашумленными данными типа [math]\displaystyle{ x_i }[/math], [math]\displaystyle{ y_i }[/math]. Наиболее известным видом сглаживающего сплайна является кубический сплайн.

Определение кубического сплайна

Пусть [math]\displaystyle{ (x_i,Y_i);x_1\lt x_2\lt \dots\lt x_n, i \in \mathbb{Z} }[/math] — последовательность наблюдений, порождённых выражением [math]\displaystyle{ Y_i = \mu(x_i) }[/math]. Приближение сглаживающими сплайнами [math]\displaystyle{ \hat\mu }[/math] функции [math]\displaystyle{ \mu }[/math] определяется как функция (в классе дважды дифференцируемых функций), минимизирующая

[math]\displaystyle{ \sum_{i=1}^n (Y_i - \hat\mu(x_i))^2 + \lambda \int_{x_1}^{x_n} \hat\mu''(x)^2 \,dx. }[/math]

Замечания:

1. [math]\displaystyle{ \lambda \ge 0 }[/math] параметр сглаживания, контролирующий соотношение между точностью воспроизведения данных и «неровностью» аппроксимирующей функции.
2. интеграл вычисляется по всему диапазону [math]\displaystyle{ x_i }[/math].
3. при [math]\displaystyle{ \lambda\to 0 }[/math] (нет сглаживания), сглаживающий сплайн превращается в интерполяционный сплайн.
4. при [math]\displaystyle{ \lambda\to\infty }[/math] (бесконечное сглаживание), штраф за неровность становится преобладающим и аппроксимация превращается в линейную МНК аппроксимацию.
5. наиболее часто в современной статистической литературе используется штраф за неровность на основе второй производной, однако метод может быть легко адаптирован к использованию штрафов на основе других производных.
6. в ранней литературе, с равноудалёнными [math]\displaystyle{ x_i }[/math], для вычисления штрафа вместо производной использовались конечные разности второго и третьего порядка.
7. если сумму квадратов отклонений сплайна от исходных данных (первый член функционала) заменить на логарифм функции правдоподобия, получим оценку максимального правдоподобия со штрафной функцией. В такой постановке обычный сглаживающий сплайн представляет собой специальный случай, когда правдоподобие рассчитывается исходя из нормального распределения погрешности.

Вывод кубического сглаживающего сплайна

Разделим нахождение выражений, описывающих сглаживающий сплайн, на два этапа:

1. Сначала найдём значения [math]\displaystyle{ \hat\mu(x_i);i=1,\ldots,n }[/math].
2. Из этих значений найдём [math]\displaystyle{ \hat\mu(x) }[/math] для всех x.

Начнём со второго этапа:

Дан вектор [math]\displaystyle{ \hat{m} = (\hat\mu(x_1),\ldots,\hat\mu(x_n))^T }[/math] «подогнанных» значений; сумма квадратов в критерии сплайна — константа. Требуется только минимизировать [math]\displaystyle{ \int \hat\mu''(x)^2 \, dx }[/math], и минимизация — натуральный кубический сплайн, интерполирующий точки [math]\displaystyle{ (x_i,\hat\mu(x_i)) }[/math]. Данный интерполяционный сплайн — линейный оператор — может быть представлен в виде:

[math]\displaystyle{ \hat\mu(x) = \sum_{i=1}^n \hat\mu(x_i) f_i(x) }[/math],

где [math]\displaystyle{ f_i(x) }[/math] — набор базисных сплайн-функций. В результате штраф за отсутствие у функции признака гладкости имеет форму

[math]\displaystyle{ \int \hat\mu''(x)^2 dx = \hat{m}^T A \hat{m}. }[/math]

где элементы A — [math]\displaystyle{ \int f_i''(x) f_j''(x)dx }[/math]. Базисные функции и матрица A зависят от конфигурации независимых переменных [math]\displaystyle{ x_i }[/math], но не от [math]\displaystyle{ Y_i }[/math] или [math]\displaystyle{ \hat m }[/math].

Возвращаясь к первому этапу, взвешенная сумма квадратов может быть записана так:

[math]\displaystyle{ \|Y - \hat m\|^2 + \lambda \hat{m}^T A \hat m, }[/math]

где [math]\displaystyle{ Y=(Y_1,\ldots,Y_n)^T }[/math]. минимизация по [math]\displaystyle{ \hat m }[/math] даёт

[math]\displaystyle{ \hat m = (I + \lambda A)^{-1} Y. }[/math]

Создание многомерных сплайнов

Из приведённого ограничения на формулу из определения [math]\displaystyle{ x_1\lt x_2\lt \dots \lt x_n }[/math] следует, что алгоритм не работает для произвольного набора данных. Если планируется использование алгоритма для произвольного набора точек в многомерном пространстве необходим алгоритм, в котором нет таких ограничений. Возможное решение заключается во введении параметра таким образом, что входные данные могут быть представлены как одномерные функции, зависящие от данного параметра; после можно применить сглаживание для каждой функции. В двумерном пространстве решение состоит в параметризации [math]\displaystyle{ x }[/math] и [math]\displaystyle{ y }[/math] как [math]\displaystyle{ x(t) }[/math] and [math]\displaystyle{ y(t) }[/math] где [math]\displaystyle{ t_1\lt t_2\lt \dots \lt t_n }[/math]. Подходящее решение для [math]\displaystyle{ t }[/math] это накопленное расстояние [math]\displaystyle{ t_{i+1}=t_{i}+\sqrt{(x_{i+1}-x_{i})^2+(y_{i+1}-y_{i})^2} }[/math] где [math]\displaystyle{ t_1=0 }[/math].

Более детальный анализ параметризации выполнен E.T.Y Lee.

Связанные методы

Сглаживающие сплайны имеют отношение, но отличаются от:

• Регрессионных сплайнов (англ. Regression Splines). Метод, при использовании которого данные аппроксимируются с помощью набора базисных сплайн-функций с уменьшенным количеством узлов, в большинстве случаев при помощи метода наименьших квадратов. При этом в случае отсутствия у функции признака гладкости штрафы не используются.
• Штрафных сплайнов, сплайнов со штрафами (англ. Penalized Splines). Сочетают уменьшенное количество узлов регрессионных сплайнов со штрафом за отсутствие у функций сглаживающих сплайнов признака гладкости.
• Метод упругой карты. Метод, сочетающий штрафы по методу наименьших квадратов для ошибки аппроксимации со штрафами за кривизну и растяжение аппроксимирующего множества и использующий крупный шаг дискретизации для оптимизации проблемы.

Исходный код

Исходный код для сглаживающих сплайнов может быть взят из примеров к книге Carl de Boor’s A Practical Guide to Splines. Примеры написаны на Фортране. Обновлённые исходные коды также доступны на официальном сайте Carl de Boor’s [1].

Литература

• Wahba, G. (1990). Spline Models for Observational Data. SIAM, Philadelphia.
• Green, P. J. and Silverman, B. W. (1994). Nonparametric Regression and Generalized Linear Models. CRC Press.
• De Boor, C. (2001). A Practical Guide to Splines (Revised Edition). Springer.
• Березовский, М. В. Сглаживающие изогеометрические и робастные сплайны: методы и алгоритмы. Диссертация.