Полигармонический сплайн - Polyharmonic spline

Полигармонические сплайны используются для аппроксимация функции и данные интерполяция. Они очень полезны для интерполяции и подгонки разрозненных данных во многих измерениях. Особые случаи включают шлицы тонкой пластины^[1]^[2] и естественные кубические шлицы в одном измерении.^[3]

Определение

Полигармонический сплайн - это линейная комбинация полигармонических радиальные базисные функции (RBFs) обозначается ${displaystyle phi}$ плюс полиномиальный член:

{displaystyle f (mathbf {x}) = sum _ {i = 1} ^ {N} w_ {i} phi (| mathbf {x} -mathbf {c} _ {i} |) + mathbf {v} ^ { extrm {T}} {egin {bmatrix} 1 mathbf {x} end {bmatrix}}}

(1)

куда

Полигармонические базисные функции

${displaystyle mathbf {x} = [x_ {1} x_ {2} cdots x_ {d}] ^ {extrm {T}}}$ ( ${displaystyle {extrm {T}}}$ обозначает транспонирование матрицы, что означает ${displaystyle mathbf {x}}$ вектор-столбец) является вещественным вектором ${displaystyle d}$ независимые переменные,
${displaystyle mathbf {c} _ {i} = [c_ {i, 1} c_ {i, 2} cdots c_ {i, d}] ^ {extrm {T}}}$ находятся ${displaystyle N}$ векторы того же размера, что и ${displaystyle mathbf {x}}$ (часто называемые центрами), которые кривая или поверхность должна интерполировать,
${displaystyle mathbf {w} = [w_ {1} w_ {2} cdots w_ {N}] ^ {extrm {T}}}$ являются ${displaystyle N}$ веса RBFs,
${displaystyle mathbf {v} = [v_ {1} v_ {2} cdots v_ {d + 1}] ^ {extrm {T}}}$ являются ${displaystyle d + 1}$ веса полинома.

Многочлен с коэффициентами ${displaystyle mathbf {v}}$ повышает точность подгонки полигармонических сглаживающих сплайнов, а также улучшает экстраполяцию вдали от центров ${displaystyle mathbf {c} _ {i}.}$ См. Рисунок ниже для сравнения шлицев с полиномиальным членом и без полиномиального члена.

Полигармонические RBF имеют вид:

{displaystyle {egin {matrix} phi (r) = {egin {case} r ^ {k} & {mbox {with}} k = 1,3,5, точки, r ^ {k} ln (r) & {mbox {with}} k = 2,4,6, точки на конце {case}} [5mm] r = | mathbf {x} -mathbf {c} _ {i} | = {sqrt {(mathbf {x} -mathbf {c} _ {i}) ^ {T}, (mathbf {x} -mathbf {c} _ {i})}}. конец {матрица}}}

Другие значения показателя степени ${displaystyle k}$ бесполезны (например, ${displaystyle phi (r) = r ^ {2}}$ ), потому что решение проблемы интерполяции может не существовать. Чтобы избежать проблем при ${displaystyle r = 0}$ (поскольку ${displaystyle log (0) = - infty}$ ) полигармонические RBF с натуральным логарифмом могут быть реализованы как:

{displaystyle phi (r) = {egin {case} r ^ {k-1} ln (r ^ {r}) & {mbox {for}} r <1 r ^ {k} ln (r) & {mbox {for}} rgeq 1end {case}}.}

Веса ${displaystyle w_ {i}}$ и ${displaystyle v_ {j}}$ определены так, что функция интерполирует ${displaystyle N}$ данные баллы ${displaystyle (mathbf {c} _ {i}, f_ {i})}$ (за ${displaystyle i = 1,2, точки, N}$ ) и выполняет ${displaystyle d + 1}$ условия ортогональности

{displaystyle sum _ {i = 1} ^ {N} w_ {i} = 0, ;; sum _ {i = 1} ^ {N} w_ {i} mathbf {c} _ {i} = mathbf {0} .}

В совокупности эти ограничения эквивалентны симметричной линейной системе уравнений

{displaystyle {egin {bmatrix} A&B B ^ {extrm {T}} & 0end {bmatrix}}; {egin {bmatrix} mathbf {w} mathbf {v} end {bmatrix}}; =; {egin {bmatrix} mathbf {f} mathbf {0} конец {bmatrix}} ;;;;}

(2)

куда

{displaystyle A_ {i, j} = phi (| mathbf {c} _ {i} -mathbf {c} _ {j} |), quad B = {egin {bmatrix} 1 & 1 & cdots & 1 mathbf {c} _ {1 } & mathbf {c} _ {2} & cdots & mathbf {c} _ {N} end {bmatrix}} ^ {extrm {T}}, quad mathbf {f} = [f_ {1}, f_ {2}, cdots, f_ {N}] ^ {extrm {T}}.}

Чтобы эта система уравнений имела единственное решение, ${displaystyle B}$ должен быть в полном звании. ${displaystyle B}$ является полным рангом для очень мягких условий во входных данных. Например, в двух измерениях три центра, образующие невырожденный треугольник, гарантируют, что ${displaystyle B}$ имеет полный ранг, а в трех измерениях четыре центра, образующие невырожденный тетраэдр, гарантируют, что B имеет полный ранг. Как поясняется позже, линейное преобразование, являющееся результатом ограничения области линейного преобразования ${displaystyle A}$ к пустое пространство из ${displaystyle extstyle {B ^ {extrm {T}}}}$ положительно определен. Это означает, что если ${displaystyle B}$ полного ранга, система уравнений (2) всегда имеет уникальное решение, и его можно решить с помощью Разложение Холецкого после подходящего преобразования. Вычисленные веса позволяют оценить сплайн для любых ${displaystyle mathbf {x} в mathbb {R} ^ {d}}$ используя уравнение (1). Многие практические детали реализации и использования полигармонических сплайнов объясняются в Фасхауэре.^[4] В Иске^[5] полигармонические сплайны рассматриваются как частные случаи других методов множественного разрешения при моделировании рассеянных данных.

Причина названия "полигармоническая"

Полигармоническое уравнение - это уравнение в частных производных формы ${displaystyle Delta ^ {m} f = 0}$ для любого натурального числа ${displaystyle m}$ , куда ${displaystyle Delta}$ это Оператор Лапласа. Например, бигармоническое уравнение является ${displaystyle Delta ^ {2} f = 0}$ и уравнение тригармонии ${displaystyle Delta ^ {3} f = 0}$ . Все полигармонические радиальные базисные функции являются решениями полигармонического уравнения (или, точнее, модифицированного полигармонического уравнения с Дельта-функция Дирака справа, а не 0). Например, радиальная базисная функция тонкой пластины является решением модифицированного 2-мерного бигармонического уравнения.^[6] Применяя двумерный оператор Лапласа ( ${displaystyle Delta = частично _ {xx} + частично _ {yy}}$ ) к радиальной базисной функции тонкой пластины ${displaystyle f_ {ext {tp}} (x, y) = (x ^ {2} + y ^ {2}) log {sqrt {x ^ {2} + y ^ {2}}}}$ либо вручную, либо с помощью система компьютерной алгебры показывает, что ${displaystyle Delta f_ {ext {tp}} = 4 + 4log r}$ . Применение оператора Лапласа к ${displaystyle Delta f_ {ext {tp}}}$ (это ${displaystyle Delta ^ {2} f_ {ext {tp}}}$ ) дает 0. Но 0 не совсем правильно. Чтобы увидеть это, замените ${displaystyle r ^ {2} = x ^ {2} + y ^ {2}}$ с ${displaystyle ho ^ {2} = x ^ {2} + y ^ {2} + h ^ {2}}$ (куда ${displaystyle h}$ некоторое небольшое число, стремящееся к 0). Оператор Лапласа применяется к ${displaystyle 4log ho}$ дает ${displaystyle Delta ^ {2} f_ {ext {tp}} = 8h ^ {2} / ho ^ {4}}$ . За ${displaystyle (x, y) = (0,0),}$ правая часть этого уравнения стремится к бесконечности при ${displaystyle h}$ приближается к 0. Для любого другого ${displaystyle (x, y)}$ , правая часть стремится к 0 при ${displaystyle h}$ приближается к 0. Это означает, что правая часть является дельта-функцией Дирака. Система компьютерной алгебры покажет, что

{displaystyle lim _ {h o 0} int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} 8h ^ {2} / (x ^ {2} + y ^ {2} + h ^ {2}) ^ {2}, dx, dy = 8pi.}

Таким образом, радиальная базисная функция тонкой пластины является решением уравнения ${displaystyle Delta ^ {2} f_ {ext {tp}} = 8pi delta (x, y)}$ .

Применяя трехмерный лапласиан ( ${displaystyle Delta = частично _ {xx} + частично _ {yy} + частично _ {zz}}$ ) к бигармоническому RBF ${displaystyle f_ {ext {bi}} (x, y, z) = {sqrt {x ^ {2} + y ^ {2} + z ^ {2}}}}$ дает ${displaystyle Delta f_ {ext {bi}} = 2 / r}$ и применяя 3D ${displaystyle Delta ^ {2}}$ оператор к тригармоническому RBF ${displaystyle f_ {ext {tri}} (x, y, z) = (x ^ {2} + y ^ {2} + z ^ {2}) ^ {3/2}}$ дает ${displaystyle Delta ^ {2} f_ {ext {tri}} = 24 / r}$ . Сдача ${displaystyle ho ^ {2} = x ^ {2} + y ^ {2} + z ^ {2} + h ^ {2}}$ и вычисления ${displaystyle Delta (1 / ho) = - 3h ^ {2} / ho ^ {5}}$ снова указывает на то, что правая часть PDE для бигармонических и тригармонических RBF являются дельта-функциями Дирака. С

{displaystyle lim _ {h o 0} int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} -3h ^ {2} / (x ^ { 2} + y ^ {2} + z ^ {2} + h ^ {2}) ^ {5/2}, dx, dy, dz = -4pi,}

точные PDE, удовлетворяющие бигармоническим и тригармоническим RBF, равны ${displaystyle Delta ^ {2} f_ {ext {bi}} = - дельта 8pi (x, y, z)}$ и ${displaystyle Delta ^ {3} f_ {ext {tri}} = - дельта 96 пикселей (x, y, z)}$ .

Полигармонические сглаживающие сплайны

Полигармонические сплайны минимизировать

{displaystyle sum _ {i = 1} ^ {N} (f (mathbf {c} _ {i}) - f_ {i}) ^ {2} + lambda int limits _ {{mathcal {B}} subset mathbb { R} ^ {d}} | abla ^ {m} f | ^ {2}, dmathbf {x}}

(3)

куда ${displaystyle {mathcal {B}}}$ какая-то коробка в ${displaystyle mathbb {R} ^ {d}}$ содержащая окрестность всех центров, ${displaystyle lambda}$ - некоторая положительная константа, а ${displaystyle abla ^ {m} f}$ вектор всех ${displaystyle m}$ частные производные порядка ${displaystyle f.}$ Например, в 2D ${displaystyle abla ^ {1} f = (f_ {x} f_ {y})}$ и ${displaystyle abla ^ {2} f = (f_ {xx} f_ {xy} f_ {yx} f_ {yy})}$ и в 3D ${displaystyle abla ^ {2} f = (f_ {xx} f_ {xy} f_ {xz} f_ {yx} f_ {yy} f_ {yz} f_ {zx} f_ {zy} f_ {zz})}$ . В 2D ${displaystyle | abla ^ {2} f | ^ {2} = f_ {xx} ^ {2} + 2f_ {xy} ^ {2} + f_ {yy} ^ {2},}$ превращая интеграл в упрощенный тонкая пластина энергетический функционал.

Чтобы показать, что полигармонические сплайны минимизируют уравнение (3), аппроксимирующий член должен быть преобразован в интеграл, используя определение дельта-функции Дирака:

{displaystyle sum _ {i = 1} ^ {N} (f (mathbf {c} _ {i}) - f_ {i}) ^ {2} = int limits _ {mathcal {B}} sum _ {i = 1} ^ {N} (f (mathbf {x}) -f_ {i}) ^ {2} delta (mathbf {x} -mathbf {c} _ {i}), dmathbf {x}.}

Итак, уравнение (3) можно записать как функционал

{displaystyle J [f] = int limits _ {mathcal {B}} F ​​(mathbf {x}, f, частично ^ {alpha _ {1}} f, частично ^ {alpha _ {2}} f, точки частично ^ {alpha _ {n}} f), dmathbf {x} = int limits _ {mathcal {B}} left [sum _ {i = 1} ^ {N} (f (mathbf {x}) -f_ {i} ) ^ {2} delta (mathbf {x} -mathbf {c} _ {i}) + lambda | abla ^ {m} f | ^ {2} ight], dmathbf {x}.}

куда ${displaystyle alpha _ {i}}$ это мультииндекс пробегает все частные производные порядка ${displaystyle m}$ за ${displaystyle mathbb {R} ^ {d}.}$ Чтобы применить Уравнение Эйлера-Лагранжа для одной функции нескольких переменных и производных более высокого порядка величины

{displaystyle {частичное F вместо частичного f} = 2sum _ {i = 1} ^ {N} (f (mathbf {x}) -f_ {i}) delta (mathbf {x} -x_ {i})}

и

{displaystyle sum _ {i = 1} ^ {m ^ {2}} partial ^ {alpha _ {i}} {частичное F над частичным (частичное ^ {alpha _ {i}} f)} = 2 Дельта лямбда ^ {m } f}

необходимы. Подставление этих величин в уравнение E-L показывает, что

{displaystyle left (sum _ {i = 1} ^ {N} (f (mathbf {x}) -f_ {i}) delta (mathbf {x} -mathbf {c} _ {i}) ight) + (- 1) ^ {m} лямбда-дельта ^ {m} f = 0.}

(4)

Слабое решение ${displaystyle f (mathbf {x})}$ из (4) удовлетворяет

${displaystyle int ограничивает _ {mathcal {B}} left (sum _ {i = 1} ^ {N} (f (mathbf {x}) -f_ {i}) delta (mathbf {x} -mathbf {c} _ {i}) ight) g (mathbf {x}) + (- 1) ^ {m} lambda (Delta ^ {m} f) g (mathbf {x}), dmathbf {x} = 0}$

(5)

для всех гладких тестовых функций ${displaystyle g}$ которые исчезают за пределами ${displaystyle {mathcal {B}}.}$ Слабое решение уравнения (4) по-прежнему будет минимизировать (3), избавившись от дельта-функции путем интеграции.^[7]

Позволять ${displaystyle f}$ - полигармонический сплайн, как определено уравнением (1). Следующие расчеты покажут, что ${displaystyle f}$ удовлетворяет (5). Применяя ${displaystyle Delta ^ {m}}$ оператор к уравнению (1) дает

{displaystyle Delta ^ {m} f = sum _ {i = 1} ^ {M} w_ {i} C_ {m, d} delta (mathbf {x} -mathbf {c} _ {i})}

куда ${displaystyle C_ {2,2} = 8pi,}$ ${displaystyle C_ {2,3} = - 8pi,}$ и ${displaystyle C_ {3,3} = - 96pi.}$ Так (5) эквивалентно

{displaystyle {egin {align} int limits _ {mathcal {B}} & sum _ {i = 1} ^ {N} delta (mathbf {x} -mathbf {c} _ {i}) (f (mathbf {x} ) -f_ {i} + (- 1) ^ {m} лямбда C_ {m, d} w_ {i}) g (mathbf {x}), dmathbf {x} & = sum _ {i = 1} ^ {N} (f (mathbf {c} _ {i}) - f_ {i} + (- 1) ^ {m} лямбда C_ {m, d} w_ {i}) g (mathbf {c} _ {i }) & = 0.end {выровнено}}}

(6)

Единственно возможное решение (6) для всех тестовых функций ${displaystyle g}$ является

{displaystyle f (mathbf {c} _ {j}) - f_ {j} + (- 1) ^ {m} lambda C_ {m, d} w_ {j} = 0quad {ext {for}} j = 1, 2, точки, N}

(7)

(что подразумевает интерполяцию, если ${displaystyle lambda = 0}$ ). Объединяя определение ${displaystyle f}$ в уравнении (1) с уравнением (7) приводит к почти той же линейной системе, что и уравнение (2) за исключением того, что матрица ${displaystyle A}$ заменяется на ${displaystyle A + (- 1) ^ {m} C_ {m, d} лямбда I}$ куда ${displaystyle I}$ это ${displaystyle N imes N}$ единичная матрица. Например, для трехмерных тригармонических RBFs, ${displaystyle A}$ заменяется на ${displaystyle A + 96pi lambda I.}$

Объяснение дополнительных ограничений

В (2) нижняя половина системы уравнений ( ${displaystyle B ^ {extrm {T}} mathbf {w} = 0}$ ) дается без объяснения причин. Для объяснения сначала требуется получить упрощенную форму ${displaystyle extstyle {int _ {mathcal {B}} | abla ^ {m} f | ^ {2}, dmathbf {x}}}$ когда ${displaystyle {mathcal {B}}}$ все из ${displaystyle mathbb {R} ^ {d}.}$

Во-первых, потребуйте, чтобы ${displaystyle extstyle {sum _ {i = 1} ^ {N} w_ {i} = 0.}}$ Это гарантирует, что все производные порядка ${displaystyle m}$ и выше ${displaystyle extstyle {f (mathbf {x}) = sum _ {i = 1} ^ {N} w_ {i} phi (| mathbf {x} -mathbf {c} _ {i} |)}}$ исчезают в бесконечности. Например, пусть ${displaystyle m = 3}$ и ${displaystyle d = 3}$ и ${displaystyle phi}$ быть тригармоническим RBF. потом ${displaystyle phi _ {zzy} = 3y (x ^ {2} + y ^ {2}) / (x ^ {2} + y ^ {2} + z ^ {2}) ^ {3/2}}$ (учитывая ${displaystyle phi}$ как отображение из ${displaystyle mathbb {R} ^ {3}}$ к ${displaystyle mathbb {R}}$ ). Для данного центра ${displaystyle mathbf {P} = (P_ {1}, P_ {2}, P_ {3}),}$

{displaystyle phi _ {zzy} (mathbf {x} -mathbf {P}) = {frac {3 (y-P_ {2}) ((y-P_ {2}) ^ {2} + (x-P_ { 1}) ^ {2})} {((x-P_ {1}) ^ {2} + (y-P_ {2}) ^ {2} + (z-P_ {3}) ^ {2}) ^ {3/2}}}.}

На линии ${displaystyle mathbf {x} = mathbf {a} + tmathbf {b}}$ для произвольной точки ${displaystyle mathbf {a}}$ и единичный вектор ${displaystyle mathbf {b},}$

{displaystyle phi _ {zzy} (mathbf {x} -mathbf {P}) = {frac {3 (a_ {2} + b_ {2} t-P_ {2}) ((a_ {2} + b_ {2 } t-P_ {2}) ^ {2} + (a_ {1} + b_ {1} t-P_ {1}) ^ {2})} {((a_ {1} + b_ {1} t- P_ {1}) ^ {2} + (a_ {2} + b_ {2} t-P_ {2}) ^ {2} + (a_ {3} + b_ {3} t-P_ {3}) ^ {2}) ^ {3/2}}}.}

Разделив числитель и знаменатель этого на ${displaystyle t ^ {3}}$ показывает, что ${displaystyle extstyle {lim _ {t o infty} phi _ {zyy} (mathbf {x} -mathbf {P}) = 3b_ {2} (b_ {2} ^ {2} + b_ {1} ^ {2} ) / (b_ {1} ^ {2} + b_ {2} ^ {2} + b_ {3} ^ {2}) ^ {3/2}},}$ количество, не зависящее от центра ${displaystyle mathbf {P}.}$ Итак, в данной строке

{displaystyle lim _ {tightarrow infty} f_ {zyy} (mathbf {x}) = lim _ {tightarrow infty} sum _ {i = 1} ^ {N} w_ {i} phi _ {zyy} (mathbf {x} -mathbf {c} _ {i}) = left (сумма _ {i = 1} ^ {N} w_ {i} ight) 3b_ {2} (b_ {2} ^ {2} + b_ {1} ^ { 2}) / (b_ {1} ^ {2} + b_ {2} ^ {2} + b_ {3} ^ {2}) ^ {3/2} = 0.}

Недостаточно требовать, чтобы ${displaystyle extstyle {sum _ {i = 1} ^ {N} w_ {i} = 0,}}$ потому что в дальнейшем это необходимо для ${displaystyle f_ {alpha} g_ {eta}}$ исчезнуть на бесконечности, где ${displaystyle alpha}$ и ${displaystyle eta}$ мультииндексы такие, что ${displaystyle | alpha | + | eta | = 2m-1.}$ Для тригармонии ${displaystyle phi,}$ ${displaystyle w_ {i} u_ {j} phi _ {alpha} (mathbf {x} -mathbf {c} _ {i}) phi _ {eta} (mathbf {x} -mathbf {d} _ {j}) }$ (куда ${displaystyle u_ {j}}$ и ${displaystyle mathbf {d} _ {j}}$ веса и центры ${displaystyle g}$ ) всегда является суммой полиномов 5-й степени от ${displaystyle x,}$ ${displaystyle y,}$ и ${displaystyle z}$ делится на квадратный корень из общего многочлена степени 8. Рассмотрим поведение этих терминов в строке ${displaystyle mathbf {x} = mathbf {a} + tmathbf {b}}$ в качестве ${displaystyle t}$ приближается к бесконечности. Числитель представляет собой многочлен 5-й степени от ${displaystyle t.}$ Разделив числитель и знаменатель на ${displaystyle t ^ {4}}$ оставляет в числителе члены степени 4 и 5 и функцию ${displaystyle mathbf {b}}$ только в знаменателе. Срок 5 степени, разделенный на ${displaystyle t ^ {4}}$ это продукт пяти ${displaystyle b}$ координаты и ${displaystyle t.}$ В ${displaystyle extstyle {сумма w = 0}}$ (и ${displaystyle extstyle {sum u = 0}}$ ) ограничение заставляет это исчезать везде на линии. Срок 4 степени, разделенный на ${displaystyle t ^ {4}}$ является либо произведением четырех ${displaystyle b}$ координаты и ${displaystyle a}$ координата или произведение четырех ${displaystyle b}$ координаты и единый ${displaystyle c_ {i}}$ или же ${displaystyle d_ {j}}$ координировать. В ${displaystyle extstyle {сумма w = 0}}$ ограничение приводит к тому, что член первого типа исчезает везде на линии. Дополнительные ограничения ${displaystyle extstyle {sum _ {i = 1} ^ {N} w_ {i} mathbf {c} _ {i} = 0}}$ заставит исчезнуть второй тип члена.

Теперь определите внутренний продукт двух функций. ${displaystyle f, g: mathbb {R} ^ {d} ightarrow mathbb {R}}$ определяется как линейная комбинация полигармонических RBFs ${displaystyle phi _ {m, d}}$ с ${displaystyle extstyle {сумма w = 0}}$ и ${displaystyle extstyle {sum wmathbf {c}} = 0}$ в качестве

{displaystyle langle f, gangle = int limits _ {mathbb {R} ^ {d}} (abla ^ {m} f) cdot (abla ^ {m} g), dmathbf {x}.}

Интеграция по частям показывает, что

${displaystyle langle f, gangle = (- 1) ^ {m} int limits _ {mathbb {R} ^ {d}} f (Delta ^ {m} g), dmathbf {x}.}$

(8)

Например, пусть ${displaystyle m = 2}$ и ${displaystyle d = 2.}$ потом

${displaystyle langle f, gangle = int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} (f_ {xx} g_ {xx} + 2f_ {xy} g_ {xy} + f_ {yy} g_ {yy}), dx, dy.}$

(9)

Интегрирование первого члена этого по частям один раз дает

{displaystyle int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} f_ {xx} g_ {xx}, dx, dy = int _ {- infty} ^ {infty} f_ {x} g_ {xx} {ig |} _ {- infty} ^ {infty}, dy-int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} f_ {x} g_ {xxx}, dx, dy = -int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} f_ {x} g_ {xxx}, dx, dy}

поскольку ${displaystyle f_ {x} g_ {xx}}$ исчезает на бесконечности. Снова интеграция по частям приводит к ${displaystyle extstyle {int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} fg_ {xxxx}, dx, dy.}}$

Таким образом, интегрируя по частям дважды для каждого члена (9) дает

{displaystyle langle f, gangle = int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} f (g_ {xxxx} + 2g_ {xxyy} + g_ {yyyy}), dx, dy = int _ {- infty} ^ {infty} int _ {- infty} ^ {infty} f (Delta ^ {2} g), dx, dy.}

С ${displaystyle extstyle {(Delta ^ {m} f) (mathbf {x}) = sum _ {i = 1} ^ {N} w_ {i} C_ {m, d} delta (mathbf {x-mathbf {c} _ {i}}),}}$ (8) показывает, что

{displaystyle {egin {align} langle f, fangle & = (- 1) ^ {m} int limits _ {mathbb {R} ^ {d}} f (mathbf {x}) sum _ {i = 1} ^ { N} w_ {i} (- 1) ^ {m} C_ {m, d} delta (mathbf {x-mathbf {c} _ {i}}), dmathbf {x} = (- 1) ^ {m} C_ {m, d} sum _ {i = 1} ^ {N} w_ {i} f (mathbf {c} _ {i}) & = (- 1) ^ {m} C_ {m, d} sum _ {i = 1} ^ {N} сумма _ {j = 1} ^ {N} w_ {i} w_ {j} phi (mathbf {c} _ {i} -mathbf {c} _ {j}) = (-1) ^ {m} C_ {m, d} mathbf {w} ^ {extrm {T}} Amathbf {w} .end {выровнено}}}

Так что если ${displaystyle extstyle {сумма w = 0}}$ и ${displaystyle extstyle {sum wmathbf {c} = 0},}$

${displaystyle int ограничивает _ {mathbb {R} ^ {d}} | abla ^ {m} f | ^ {2}, dmathbf {x} = (- 1) ^ {m} C_ {m, d} mathbf {w } ^ {extrm {T}} Amathbf {w}.}$

(10)

Теперь происхождение ограничений ${displaystyle B ^ {extrm {T}} mathbf {w} = 0}$ можно объяснить. Здесь ${displaystyle B}$ является обобщением ${displaystyle B}$ определенное выше, чтобы по возможности включать мономы до степени ${displaystyle m-1.}$ Другими словами,

${displaystyle B = {egin {bmatrix} 1 & 1 & dots & 1 mathbf {c} _ {1} & mathbf {c} _ {2} & dots & mathbf {c} _ {N} vdots & vdots & dots & vdots mathbf {c} _ {1 } ^ {m-1} & mathbf {c} _ {2} ^ {m-1} & dots & mathbf {c} _ {N} ^ {m-1} end {bmatrix}} ^ {extrm {T}}}$

куда ${displaystyle mathbf {c} _ {i} ^ {j}}$ вектор-столбец любой степени ${displaystyle j}$ мономы координат ${displaystyle mathbf {c} _ {i}.}$ Верхняя половина (2) эквивалентно ${displaystyle Amathbf {w} + Bmathbf {v} -mathbf {f} = 0.}$ Итак, чтобы получить сглаживающий сплайн, нужно минимизировать скалярное поле ${displaystyle F: mathbb {R} ^ {N + d + 1} ightarrow mathbb {R}}$ определяется

{displaystyle F (mathbf {w}, mathbf {v}) = | Amathbf {w} + Bmathbf {v} -mathbf {f} | ^ {2} + lambda Cmathbf {w} ^ {extrm {T}} Amathbf { w}.}

Уравнения

{displaystyle {frac {partial F} {partial w_ {i}}} = 2A_ {i *} (Amathbf {w} + Bmathbf {v} -mathbf {f}) + 2lambda CA_ {i *} mathbf {w} = 0quad {extrm {for}} i = 1,2, точки, N}

и

{displaystyle {frac {partial F} {partial v_ {i}}} = 2B_ {i *} ^ {extrm {T}} (Amathbf {w} + Bmathbf {v} -mathbf {f}) = 0quad {extrm { for}} i = 1,2, точки, d + 1}

(куда ${displaystyle A_ {i *}}$ обозначает строку ${displaystyle i}$ из ${displaystyle A}$ ) эквивалентны двум системам линейных уравнений ${displaystyle A (Amathbf {w} + Bmathbf {v} -mathbf {f} + lambda Cmathbf {w}) = 0}$ и ${displaystyle B ^ {extrm {T}} (Amathbf {w} + Bmathbf {v} -mathbf {f}) = 0.}$ С ${displaystyle A}$ обратима, первая система эквивалентна ${displaystyle Amathbf {w} + Bmathbf {v} -mathbf {f} + lambda Cmathbf {w} = 0.}$ Итак, первая система подразумевает, что вторая система эквивалентна ${displaystyle B ^ {extrm {T}} mathbf {w} = 0.}$ Как и в предыдущем вычислении коэффициента сглаживающего сплайна, верхняя половина (2) становится ${displaystyle (A + lambda CI) mathbf {w} + Bmathbf {v} = mathbf {f}.}$

Этот вывод системы уравнений полигармонического сглаживающего сплайна не предполагал ограничений, необходимых для гарантии того, что ${displaystyle extstyle {int _ {{mathcal {mathbb {R}}} ^ {d}} | abla ^ {m} f | ^ {2}, dmathbf {x}} = Cw ^ {extrm {T}} AW. }$ Но ограничения, необходимые для этого, ${displaystyle extstyle {сумма w = 0}}$ и ${displaystyle extstyle {sum wmathbf {c} = 0},}$ являются подмножеством ${displaystyle B ^ {extrm {T}} w = 0}$ что верно для критической точки ${displaystyle w}$ из ${displaystyle F.}$ Так ${displaystyle extstyle {int _ {{mathcal {mathbb {R}}} ^ {d}} | abla ^ {m} f | ^ {2}, dmathbf {x}} = Cw ^ {extrm {T}} Aw}$ верно для ${displaystyle f}$ формируется из решения системы уравнений полигармонического сглаживающего сплайна. Поскольку интеграл положителен для всех ${displaystyle weq 0,}$ линейное преобразование, возникающее в результате ограничения области линейного преобразования ${displaystyle A}$ к ${displaystyle w}$ такой, что ${displaystyle B ^ {T} w = 0}$ должно быть положительно определенным. Этот факт позволяет преобразовать систему полигармонических сглаживающих сплайновых уравнений в симметричную положительно определенную систему уравнений, которую можно решить вдвое быстрее с помощью разложения Холецкого.^[6]

Примеры

На следующем рисунке показана интерполяция через четыре точки (отмеченные «кружками») с использованием различных типов полигармонических сплайнов. «Кривизна» интерполированных кривых растет с порядком сплайна, и экстраполяция на левой границе (x <0) разумна. На рисунке также представлены радиальные базисные функции phi = exp (-r²), что также дает хорошую интерполяцию. Наконец, на рисунке есть также неполигармонический сплайн phi = r² чтобы продемонстрировать, что эта радиальная базисная функция не может проходить через заранее определенные точки (линейное уравнение не имеет решения и решается методом наименьших квадратов).

Интерполяция с помощью различных полигармонических сплайнов, которые должны проходить через 4 предопределенные точки, отмеченные кружком (интерполяция с phi = r² бесполезен, так как система линейных уравнений интерполяционной задачи не имеет решения; решается методом наименьших квадратов, но не проходит центров)

На следующем рисунке показана та же интерполяция, что и на первом рисунке, за исключением того, что точки, которые должны быть интерполированы, масштабируются с коэффициентом 100 (и случай phi = r² больше не входит). Поскольку phi = (scale * r)^k = (масштаб^k)*р^k, коэффициент (масштаб^k) можно извлечь из матрицы А системы линейных уравнений и, следовательно, масштабирование не влияет на решение. Это отличается для логарифмической формы сплайна, хотя масштабирование не оказывает большого влияния. Этот анализ отражен на рисунке, где интерполяция не показывает больших различий. Обратите внимание, что для других радиальных базисных функций, таких как phi = exp (-k * r²) с k = 1 интерполяция больше нецелесообразна и необходимо адаптировать k.

Та же интерполяция, что и на первом рисунке, но точки для интерполяции масштабируются на 100

На следующем рисунке показана та же интерполяция, что и на первом рисунке, за исключением того, что полиномиальный член функции не принимается во внимание (и случай phi = r² больше не входит). Как видно из рисунка, экстраполяция для x <0 уже не такая «естественная», как на первом рисунке для некоторых базисных функций. Это означает, что полиномиальный член полезен при экстраполяции.

Та же интерполяция, что и на первом рисунке, но без полиномиального члена

Обсуждение

Основное преимущество полигармонической сплайн-интерполяции состоит в том, что обычно очень хорошие результаты интерполяции получаются для разрозненных данных без выполнения какой-либо «настройки», поэтому возможна автоматическая интерполяция. Это не относится к другим радиальным базисным функциям. Например, функция Гаусса ${displaystyle e ^ {- kcdot r ^ {2}}}$ нужно настроить, чтобы ${displaystyle k}$ выбирается в соответствии с базовой сеткой независимых переменных. Если эта сетка неоднородна, правильный выбор ${displaystyle k}$ добиться хорошего результата интерполяции сложно или невозможно.

Основные недостатки:

Чтобы определить веса, необходимо решить плотную линейную систему уравнений. Решение плотной линейной системы становится непрактичным, если размер ${displaystyle N}$ большой, так как требуется память ${displaystyle O (N ^ {2})}$ а количество необходимых операций равно ${displaystyle O (N ^ {3}).}$
Оценка вычисленной полигармонической сплайн-функции на ${displaystyle M}$ точки данных требует ${displaystyle O (MN)}$ операции. Во многих приложениях (например, обработка изображений), ${displaystyle M}$ намного больше, чем ${displaystyle N,}$ и если оба числа большие, это непрактично.

Недавно были разработаны способы преодоления вышеупомянутых трудностей. Например, Beatson et al.^[8] представить метод интерполяции полигармонических сплайнов в одной точке в трех измерениях в ${displaystyle O (log N)}$ операции вместо ${displaystyle O (N)}$ операции.

Смотрите также

Взвешивание обратных расстояний
Радиальная базисная функция
Подразделение поверхности (новая альтернатива сплайн-поверхностям)
Сплайн