Материал: Matlab. Практический подход. Самоучитель

Самоучитель Matlab

и аппроксимация на основе полиномиального выражения. Коэффициенты соответствующего полинома вычисляются командой P2=polyfit(x,y,5). Аппроксимация выполняется на основе полинома пятой степени (третий аргумент функции polyfit()). Кривая на основе этого полинома отображается с помощью команды plot(z,polyval(P2,z),'k-.').

Заголовок добавляем командой title('Интерполяционный полином'), а легенду – командой legend('базовые точки','интерполяция', 'аппроксимация').Результат графических построенийпредставленна рис.8.4.

Рис. 8.4. Базовые точки, интерполяционный и аппроксимирующий полиномы

Интерполяционная кривая сплошная, а штрихпунктирная кривая соответствует полиномиальной аппроксимации. Что касается интерполяционного полинома, то по определению в узловых точках он дает табличные значения (те, на основе которых строится полином). Как отмечалось выше, в задаче аппроксимации такое условие не ставится. Достаточно, чтобы кривая проходила достаточно близко к базовым точкам. Что такое "достаточно близко" – вопрос отдельный.

На заметку

Каким бы способом ни строился интерполяционный полином, он один и тот же для фиксированного набора базовых точек. Другими словами, для данного набора базовых точек интерполяционный полином один. Как бы мы ни вычисляли этот полином, если вычисления правильные, получим один и тот же результат.

316

Глава 8. Обработка данных

При большом количестве базовых точек интерполяция полиномом может оказаться малопродуктивной. Одна из проблем связана с существенно немонотонным поведением полиномов высокой степени. Чтобы решить эту проблему, нередко используют интерполяцию сплайнами. Главная идея метода состоит в том, что диапазон интерполирования разбивается на поддиапазоны (например, в пределах между двумя соседними узловыми точками). В пределах каждого поддиапазона используются разные функции одного вида (полиномы невысокой степени). Эти функции "сшиваются" на границах поддиапазонов вместе с производными. В результате получаем кусочно-гладкую (вместе с производными нескольких первых порядков) интерполяционную зависимость.

Сплайн-интерполяцию в Matlab можно реализовать с помощью функции interp1(). Обычно аргументами функции interp1() передают набор узловых точек аргумента табулированной функции и значений функции в этих узловых точках. Это два списка, определяющие те базовые точки, на основе которых выполняется интерполяция. Третьим аргументом указывается список значений, для которых вычисляется значение интерполяционной зависимости. Этот аргумент может быть скаляром – тогда значение интерполяционной зависимости вычисляется в одной точке. Также четвертым аргументом можно указать (в одинарных кавычках) ключевое слово, определяющее тип базового сплайна (степень полинома). Для определения типа базового полинома используют следующие ключевые слова (табл. 8.1).

Табл. 8.1. Ключевые слова, определяющие тип базового полинома

Как правило, к интерполяции на основе базовых полиномов нулевой и первой степени на практике прибегают достаточно редко, поэтому рассмотрим методы интерполяции кубическими полиномами. Обратимся к рис. 8.5, на котором представлен документ с соответствующим примером.

317

Самоучитель Matlab

Рис. 8.5. Интерполяция кусочно-непрерывными функциями

Использована такая последовательность команд:

>>x=-1:0.2:1;

>>y=cos(3*pi*x);

>>z=-1:0.01:1;

>>f1=interp1(x,y,z,'spline');

>>f2=interp1(x,y,z,'pchip');

>>plot(x,y,'ks')

>>hold on

>>plot(z,f1,'r-','LineWidth',2)

>>plot(z,f2,'g--','LineWidth',2)

>>grid on

>>title('Интерполяция сплайнами')

>>legend('базовые точки','кубические сплайны','полиномы Эрмита')

В переменную x записываются узловые точки аргумента табулируемой функции. Значения табулируемой функции в узловых точках вычисляются командой y=cos(3*pi*x). В переменную z записываются значения аргумента, для которых необходимо вычислить интерполяционные зависимости. Интерполяционная зависимость на основе кубических сплайнов строится командой f1=interp1(x,y,z,'spline'). Здесь первые два аргумента (переменные x и y) задают базовые точки, на основе которых строится интерполяционная зависимость. Третий аргумент (переменная z) задает точки аргумента, для которых необходимо вычислить значения интерполяционной зависимости. Эти значения, как результат вызова функции interp1(), записываются в переменную f1. Опционный параметр 'spline' означает, что интерполяционная зависимость строится на основе базовых полиномов третьей степени. Аналогично, интерполяционная

318

Глава 8. Обработка данных

зависимость на основе полиномов Эрмита создается с помощью команды f2=interp1(x,y,z,'pchip').

На заметку

На рис. 8.6 показан результат интерполирования разными методами (обычными сплайнами и полиномами Эрмита) одной и той же табличной функции.

Рис. 8.6. Интерполяция кубическими сплайнами и полиномами Эрмита

Рассуждать вне контекста конкретной решаемой задачи о том, какой способ интерполирования лучше, особого смысла нет. Отметим лишь, каждый способ по-своему хорош.

319

Самоучитель Matlab

Аппроксимация

Это простейшая цепь рассуждений.

К/ф "Приключения Шерлока Холмса и доктора Ватсона. Знакомство"

Задача аппроксимации близка идеологически к задаче интерполяции, но все же это разные задачи. Например, предположим, что имеется набор узловых точек xk (индекс k = 1,2,...,n ) и значения некоторой функции yk в этих точках. Также есть известная функция f (x,a1,a2,...,am ), которая кроме аргумента x зависит еще и от некоторых параметров as (индекс s = 1,2,...,m ). Задача состоит в том, чтобы подобрать такие значения параметров as , что функция f (x,a1,a2,...,am ) наилучшим образом описывала бы зависимость, заданную параметрами xk и yk (индекс k = 1,2,...,n ). При этом количество m параметров as меньше (или даже намного меньше) количества узловых точек n . Поэтому добиться того, чтобы функция f (x,a1,a2,...,am ) давала "точные" значения yk в узловых точках xk , не удастся. Нужен критерий, который бы позволил определить, какая зависимость "лучше", а какая "хуже" аппроксимирует табличные значения. Обычно в качестве такого критерия используется принцип наибольшего правдоподобия, следствием которого является метод наименьших квадратов. В соответствии с этим методом параметры as (индекс s = 1,2,...,m )

n

выбираются так, чтобы сумма ∑(yk − f(xk,a1,a2,...,am ))2 принимала наи-

k =1

меньшее значение. Такая задача сводится к решению системы из m алгебраических уравнений – в общем случае нелинейных. Система решается относительно параметров as (индекс s = 1,2,...,m ). Эти уравнения име-

(для всех индексов s = 1,2,...,m ). В случае если функция f (x,a1,a2,...,am )

зависит от параметров as линейным образом, задача значительно упрощается и сводится к решению системы линейных уравнений. Действительно, в этом случае для функции f (x,a1,a2,...,am ) может быть использо-

320