Одномерная сплайновая интерполяция и аппроксимация

При небольшом числе узловых точек (менее 10) линейная интерполя­ция оказывается довольно грубой. При ней даже первая производная функции аппроксимации испытывает резкие скачки в узловых точках. Для целей экстраполяции функция linterp нe предназначена и за пре­делами области определения может вести себя непредсказуемо. Гораздо лучшие результаты дает сплайновая интерполяция. При ней исходная функция заменяется отрезками кубических полиномов, про­ходящих через три смежные узловые точки. Триады точек перемеща­ются по оси абсцисс, что создает набор полиномов. Коэффициенты по­линомов рассчитываются так, чтобы непрерывными были первая и вторая производные. Линия, которую описывает сплайн-функция, на­поминает по форме гибкую линейку, закрепленную в узловых точках. Для осуществления сплайновой интерполяции система Mathcad пред­лагает четыре встроенные функции. Три из них служат для получения векторов вторых производных сплайн - функций при различном виде интерполяции:

• cspline(VX, VY) — возвращает вектор VS вторых производных при приближении в опорных точках к кубическому полиному;

• pspline(VX, VY) — возвращает вектор VS вторых производных при приближении в опорных точках к параболической кривой;

• lspline(VX, VY) — возвращает вектор VS вторых производных при приближении в опорных точках к прямой.

Рис. 1

Четвертая функция — interp(VS, VX, VY, x) — возвращает значение у(х) для заданных векторов VS, VX, VY и заданного значения х. Таким образом, сплайн - интерполяция проводится в два этапа. На пер­вом с помощью функции сsplinе, рsplinе или lspline отыскивается вектор вторых производных функции у(х), заданной векторами VX и VY ее значений (абсцисс и ординат). Затем на втором этапе для каж­дой искомой точки вычисляется значение у(х) с помощью функции interp. Дополнительно следует использовать функции:

• csort(A, n) - перестановка строк матрицы А таким образом, чтобы отсортированным в порядке возрастания значений элементов оказался п- й столбец;

• length(v) – число элементов в векторе v.

На рис. 1 показан фрагмент документа Mathcad, иллюстрирующий применение описанных функций для линейной и сплайновой интер­поляций. Нетрудно заметить, что при линейной интерполяции график оказы­вается слишком грубым — отчетливо видны точки перегибов. В то же время сплайновая интерполяция, несмотря на малое число точек в этом примере (их всего 10), дает прекрасные результаты: график функции оказывается плавным и точки его перегиба вообще незаметны. Впрочем, стоит отметить, что при неудачном расположении узлов вы­бег сплайновой интерполяции может оказаться весьма значительным.

Регрессия

Еще одной широко распространенной задачей обработки данных явля­ется представление их совокупности некоторой функцией у(х). Зада­ча регрессии заключается в получении параметров этой функции таки­ми, чтобы функция приближала «облако» исходных точек (заданных векторами VX и VY) с наименьшей среднеквадратичной погрешностью. Mathcad позволяет выполнять также многомерную регрессию. Самый типичный случай ее использования — приближение поверхностей в трехмерном пространстве. Их можно описать, задав массив значений высот z, соответствующих двухмерному массиву Мху координат точек (х, у) на горизонтальной плоскости. Новых функций для этого не задано. Используются уже описанные ранее функции, но в несколько иной форме:

• regress (Мху ,Vz, n) — возвращает вектор, запрашиваемый функ­цией interp(VS, Мху ,Vz , V) для вычисления многочлена n-й степе­ни, который наилучшим образом приближает точки множества Мху и Vz (Мху — матрица размера m ´ 2, содержащая координаты х и у, Vz — m-мерный вектор, содержащий z-координаты, соответствую­щие m точкам, указанным в Мху);

• loess(Mxy, Vz, span) —аналогична loess(VX, VY, span), но в мно­гомерном случае;

• interp(VS, Mxy ,Vz, V) — возвращает значение z по заданным век­торам VS (создается функциями regress или loess) и Мху, Vz и V (вектор координат х и у заданной точки, для которой находится z).

Линейная регрессия

Чаще всего используется линейная регрессия, при которой функция у(х) описывает отрезок прямой и имеет вид: у(х) = а + bх. К линейной регрессии можно свести многие виды нелинейной регрес­сии при зависимостях вида у(х). Для проведения линейной регрессии в систему встроен ряд приведен­ных ниже функций:

• corr(VX, VY) — возвращает скаляр — коэффициент корреляции Пирсона;

• intercrpt(VX, VY) — возвращает значение параметра а (смещение линии регрессии по вертикали);

• slope(VX, VY) — возвращает значение параметра b (угловой коэф­фициент линии регрессии).

Рис. 2

На рис. показан фрагмент документа Mathcad с примером проведе­ния линейной регрессии для данных, представленных значениями эле­ментов в векторах VX и VY. Как видно из рисунка, линия регрессии проходит в «облаке» исходных точек с максимальным среднеквадратичным приближением к ним. Чем ближе коэффициент корреляции к 1, тем точнее представленная ис­ходными точками зависимость приближается к линейной.