MATLAB.Exponenta
–Û·Ë͇ Matlab&Toolboxes

Spline Toolbox

Назначение и возможности приложения Spline Tool

В состав Spline Toolbox входит приложение Spline Tool с графическим интерфейсом пользователя, предназначенное для аппроксимации данных сплайнами. Приложение Spline Tool позволяет:

  1. Указать приближаемую функцию, заданную формулой или таблично, или выбрать один из стандартных примеров для исследования (см. разд. Запуск приложения Spline Tool. Импорт данных.).
  2. Строить приближения при помощи:
    • интерполяционных кубических сплайнов с различными граничными условиями (см. разд. Приближение кубическими сплайнами);
    • интерполяционных сплайнов от 1-го до 14-го порядков;
    • сглаживающих сплайнов выбирая способ и параметр сглаживания (см. разд. Приближение сглаживающими сплайнами);
    • сплайнов от 1-го до 14-го порядков по методу наименьших квадратов.
  3. Получать быстрый доступ к информации об используемых приближениях и соответствующих терминах.
  4. Изучать распределение ошибки интерполяции и аппроксимации.
  5. Изучать поведение первой и второй производных от приближающей функции.
  6. Изменять узлы интерполирования и наблюдать за качеством получающегося приближения.
  7. Сохранять в переменных рабочей среды приближение и данные для дальнейшей работы с ними
  8. Выводить результаты в отдельное графическое окно.
  9. Генерировать m-файл c функцией, вызов которой от исходных данных приведет к повторению приближения и визуализации приближения и данных в отдельном графическом окне.

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

Для чтения материала этого раздела не требуются предварительные знания MATLAB.

Запуск приложения Spline Tool. Импорт данных.

Для запуска приложения Spline Tool достаточно в командной строке MATLAB выполнить команду

>> splinetool

При этом появляется стартовое окно, в котором предлагается

  1. импортировать собственные данные (кнопка Import your oun data);
  2. выбрать один из примеров, который загружается автоматически и появляется окно с пояснением цели исследования:
    • кнопка Titanuim heat data - приближение данных сплайнами по методу наименьших квадратов, изучение влияния узлов на качество приближения;
    • кнопка Noisy values of a smooth function - приближение данных, содержащих ошибки;
    • кнопка The function sin(x) on [0..pi/2] - изучение влияния граничных условий на качество приближения кубическими сплайнами;
    • кнопка Richard Tapia's drag race data - получение ускорения по заданной табличной зависимости перемещения от времени;
    • кнопка Improve shape of interpolant by moving knots - изучение качества интерполяции кубическими сплайнами в зависимости от расположения узлов;


Стартовое окно приложения Spline Tool

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

Нажмем кнопку Import your oun data, при этом появляется диалоговое окно Please provide data, в котором в строках ввода указываются следующие величины (по порядку сверху вниз):

  1. абсциссы табличной функции;
  2. обозначение оси ординаты табличной функции;
  3. обозначение оси.

Все эти данные можно будет задать или изменить при работе в приложении Spline Tool.

Сейчас в верхней строке ввода указано linspace(0,2*pi,31). Функция linspace(a,b,n) генерирует вектор-строку из n равноотстоящих чисел от a до b, т.е. в данном случае 31 число от 0 до 2π.

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

Time	0.1   0.4	  0.5    0.7    0.9
Temp  -1.3   1.2	  1.7    1.5   -0.6

Введем в строки ввода диалогового окна Please provide data:

  1. значения времени (в квадратных скобках через пробел, т.к. вводимое выражение должно быть массивом, а массив в MATLAB формируется при помощи квадратных скобок): [0.1 0.4 0.5 0.7 0.9]
  2. название оси абсцисс: Time
  3. значения температуры: [-1.3; 1.2; 1.7; 1.5; -0.6]
  4. название оси ординат: Temp


Ввод данных в окно Please provide data

После нажатия кнопки OK в диалоговом окне Please provide data открывается основное окно приложения Spline Tool, в котором можно приближать данные и изучать качество получающегося приближения. По умолчанию данные интерполируются кубическими сплайнами и мы сразу видим результат.

Импорт данных из массивов

Если объем данных достаточно велик, то разумно сначала создать в рабочей среде соответствующие массивы, а затем указать их имена в строках ввода диалогового окна. Предположим, что мы хотим приблизить функцию y = x²sin x на отрезке [-π;&pi/2] по 12 равноотстоящим точкам. Для этого сначала создадим соответствующие массивы, выполнив две команды из командной строки MATLAB:

>> x=linspace(-pi,pi/2,12);
>> y=x.^2.*sin(x);

Здесь для генерации 12 равноотстоящих значений абсцисс использовалась функция linspace, а для заполнения массива y поэлементные операции с точкой. Точка и крышка означает поэлементное возведение в степень. а точка и звездочка - поэлементное умножение, т.е. для i=1,2,...,12 происходит вычисление значения x(i)^2*sin(x(i)) и запись его в y(i). Поэлементные операции MATLAB очень удобны, т.к. часто позволяют избежать программирования циклов. Отметим, что команда y=x^2*sin(x) привела бы к ошибке, т.к. крышка предназначена для возведения квадратной матрицы в степень, а звездочка для матричного перемножения, что в нашем случае невозможно, т.к. x является вектором.

Далее для импорта данных в приложение Spline Tool можно действовать несколькими способами.

Если приложение Spline Tool уже запущено, то в меню File выбираем пункт Import data... . Появляется предупреждение о том, что текущие данные будут потеряны. Нажимаем кнопку OK. Появляется стартовое окно приложения Spline Tool, в котором нажимаем кнопку Import your oun data. При этом появляется диалоговое окно Please provide data, в его первую и третью строки вводим названия массивов x и y (можно также ввести названия осей, как было описано выше) и нажимаем кнопку OK. Данные импортированы.

Другой способ заключается в вызове функции splinetool со входными аргументами - массивами с данными:

>> splinetool(x,y)

При этом, если приложение Spline Tool было запущено, то появится предупреждение о том, что старые данные будут потеряны, нажимаем кнопку OK и новые данные импортированы. Если приложение Spline Tool не было запущено, то просто открывается окно приложения с импортированными данными.

Импорт данных из текстовых и двоичных файлов

Часто данные, которые требуется приблизить тем или иным способом, записаны в текстовом файле, скажем в два столбца. Тогда для записи их в массивы удобно применить команду load. Предположим, что структура файла с данными data.txt такова:

    1.3300    0.7820
    1.5600    0.9870
    1.5900    1.3400
     ...

Тогда сначала записываем данные в матрицу A (у нее будет два столбца)

>> A=load('data.txt');

и затем переписываем содержимое столбцов в две переменные

>> x=A(:,1);
>> y=A(:,2);

Дальше действуем так, как описано в разделе Импорт данных из массивов.

Если указывается неполное имя файла (т.е. без пути к нему), то важно убедиться, что MATLAB сможет найти этот файл. Для этого можно задать текущим каталогом MATLAB тот, который содержит этот файл. Для задания текущего каталога можно воспользоваться кнопкой на панели инструментов MATLAB, размещенной справа от списка Current Directory.

Если в массиве данные записаны в две строки, то после вызова функции load следует переписать содержимое строк в соответствующие массивы. Для этого выполняем команды

>> x=A(1,:);
>> y=A(2,:);

Данные могут содержаться и в двоичном файле, в MATLAB это обычно файлы с расширением mat. Предположим, что у нас есть файл data.mat. Сначала при помощи функции whos выясняем, какие данные в нем записаны

>> vars = whos('-file', 'data.mat');

Сейчас vars содержит структуру с информацией о данных, хранящихся в двоичном файле data.mat. Для доступа к их именам следует обратиться к полю name этой структуры

>> vars.name
ans =
Pres
ans =
Temp
ans =
Time

В данном примере в двоичном файле data.mat содержатся переменные Pres, Temp и Time. Если нам нужны только Temp и Time, то считываем их в рабочую среду:

>> load('data.mat','Temp','Time')

а дальше действуем так, как описано в разделе Импорт данных из массивов.

Мы разобрали импорт табличных данных в приложение Spline Tool. Если мы хотим изучить приближение сплайнами функций, заданных аналитической зависимостью, то надо либо запрограммировать собственную функцию, либо использовать подходящую из набора функций MATLAB.

Импорт функций, заданных аналитической зависимостью

Если приближаемая функция одна из стандартных математических функций MATLAB, то достаточно указать ее имя в третьей строке ввода диалогового окна Please provide data (см. Запуск приложения Spline Tool. Импорт данных.). Для того, чтобы узнать список математических функций MATLAB, достаточно выполнить команду

>> doc elfun

В открывшемся окне справочной системы содержится список функций.

Если требуется приблизить какую-либо другую функцию, например x³ + x sin x, то придется запрограммировать ее в m-файле. Для этого в меню File MATLAB следует выбрать пункт New.. и подпункт M-File. Открывается окно редактора m-файлов с новым файлом. Программируем в нем требуемую функцию с использованием поэлементных операций (про поэлементные операции написано в разд. Импорт данных из массивов):

function y=myfun(x)
y=x.^3+x.*sin(x);

Сохраняем функцию myfun в файле myfun.m в текущем каталоге для того, чтобы MATLAB смог ее найти (про текущий каталог написано в разделе Импорт данных из текстовых и двоичных файлов).

Далее при импорте данных в приложение Spline Tool в диалоговом окне Please provide data в верхней строке ввода задаем абсциссы точек, которые будут использоваться для приближения, например 10 равноотстоящих на отрезке [-2, 0.5] точек:

x=linspace(-2,0.5,10) 

а в третьей сверху строке ввода пишем имя функции myfun и нажимаем на кнопку OK. Данные импортированы.

Возникает вопрос, чем же отличается импорт данных при помощи массивов от импорта данных с заданием функции. Если вместо указания функции myfun в диалоговом окне Please provide data мы напишем имя массива, содержащего значения функции, то тогда при исследовании приближений мы не сможем получить график распределения ошибки на всем интервале приближения. График будет строиться только по значениям ошибки в узлах и при интерполяции сплайнами это будут ошибки, возникшие в результате округлений чисел с плавающей точкой. При указании функции в приложении Spline Tool будет выводится график ошибки, т.е. отклонение приближения от исходной функции.

Приближение и интерполяция данных

В этом разделе мы рассмотрим интерполяцию и приближение при помощи:

  1. интерполяционных кубических сплайнов с различными граничными условиями;
  2. сглаживающих сплайнов;
  3. интерполяционных сплайнов различных порядков;
  4. сплайнов по методу наименьших квадратов.

Приближение кубическими сплайнами

В качестве примера проинтерполируем функцию y = cos x на отрезке [0.5π/2] по 7 равноотстоящим точкам при помощи кубических сплайнов с различными граничными условиями.

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

Граничные условия при интерполяции кубическими сплайнами обязательно задаются по следующей причине. Если мы интерполируем кубическим сплайном по n точкам x1, x2, …, xn, то непрерывность кубического сплайна и его производных до второго порядка включительно во внутренних точках означает, что мы наложили 3(n - 2) условий на коэффициенты полиномов третьей степени. Кроме того, кубический сплайн интерполяционный, следовательно он должен проходить через те же точки, что и приближаемая функция, что приводит еще к условиям. Следовательно всего условий

3(n - 2) + n = 4n - 6

Всего полиномов n - 1, следовательно у нас имеется 4(n - 1) = 4n - 4 искомых коэффициентов. Получается, что условий на 2 меньше, чем неизвестных величин, Следовательно требуется наложить еще два условия для однозначного определения сплайна.

Итак, запустим приложение Spline Tool, введем собственные исходные данные (см. Запуск приложения Spline Tool. Импорт данных.). В нашем примере в диалоговом окне Please provide data в верхней строке ввода задаем абсциссы 7 равноотстоящих точек на отрезке [0.5π/2] при помощи linspace(0,5*pi/2,7), а в третьей строке ввода указываем функцию cos.

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


Приближение кубическим сплайном

Граничные условия not-a-knot (отсутствие узла)

В окне Spline Tool на панели Approximation method в списке по умолчанию выбрано Cubic Spline Interpolation (т.е. как раз кубический сплайн), а в списке End conditions по умолчанию выбрано not-a-knot, что значит граничное условие - отсутствие узла. Это граничное условие означает, что третья производная сплайна в точках x2 и xn-1 непрерывна, т.е. полиномы на отрезках [x1,x2] и [x2,x3] оказываются одинаковыми. Аналогично, полиномы на отрезках [xn-2,xn-1] и [xn-1,xn] так же одинаковы (отсюда и название not-a-knot, т.е. в переводе "отсутствие узла").

Под списком с граничными условиями расположены кнопки и области ввода для задания первой или второй производных сплайна (кнопка 1st Deriv. после щелчка по ней меняет название на 2nd Deriv.). Но области ввода сейчас недоступны и нажатие на кнопки приводит к тому, что можно только посмотреть значения производных в граничных точках. Это происходит потому, что граничные условия not-a-knot не предполагают задание значений производных кубического сплайна в граничных точках отрезка интерполирования.

Легко убедиться, что первая производная сплайна непрерывна. Для этого в меню View следует выбрать пункт Show 1st Derivative. Тогда вместо графика ошибки интерполяции будет выводиться график первой производной кубического сплайна. Удобно также нанести на график информацию о положении точек разрыва, для чего в раскрывающемся списке Data, breaks/knots, weights следует выбрать Breaks.


График первой производной кубического сплайна с граничными условиями not-a-knot.

Аналогичным образом можно получить график второй производной кубического сплайна, для чего в меню View следует выбрать пункт Show 2nd Derivative.


График второй производной кубического сплайна с граничными условиями not-a-knot.

Из графика второй производной кубического сплайна (а это кусочно-линейная функция) видно, что третья его производная будет разрывна во всех узлах интерполяции, кроме первого внутреннего и предпоследнего внутреннего узлов, т.к. в качестве граничных условий выбрано условие not-a-knot.

Внизу окна приложения Spline Tool в области вывода Bottomline выводится обращение к соответствующей функции Spline Toolbox, вызов которой приведет к аналогичному результату.


Показ соответствующей функции Spline Toolbox.

В данном случае, это функция csapi, которая как раз и предназначена для интерполирования кубическими сплайнами с граничными условиями not-a-knot. Разумеется, вызова одной этой функции недостаточно для построения графика интерполируемых данных и сплайна (по поводу функции csapi см. Функция csapi для интерполяции с условиями отсутствия узла (not-a-knot)).

Рассмотрим интерполяцию кубическими сплайнами с заданием других граничных условий в раскрывающемся списке End conditions. Для получения краткой справки по граничным условиям можно в меню Help перейти к пункту Explanation of Terms, далее Interpolation и End Conditions и в появившемся подменю нужны граничные условия (см. также Задание граничных условий в функции csape). Большинство граничных условий разумно использовать, если известна дополнительная информация об интерполируемой функции на границах отрезка интерполирования.

Граничные условия clamped и complete (условия на первую производную)

При выборе в раскрывающемся списке End conditions условий clamped (или complete) появляется возможность задания первой производной сплайна в граничных точках отрезка интерполирования. Значения производной задаются в областях ввода под соответствующей кнопкой 1st Deriv. Для изменения значений с заданным шагом и немедленным выводом графика соответствующего кубического сплайна можно также использовать кнопки со знаками "плюс" и "минус" и указать желаемый шаг изменения значений в расположенной между этими кнопками области ввода.

При нажатии на кнопку 1st Deriv она меняет название на 2nd Deriv, но область ввода для значения производной становится недоступной, можно только посмотреть значение второй производной, т.к. при выборе граничных условий clamped (или complete) мы имеем право задавать только значение первой производной.

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

Граничные условия second (условия на вторую производную)

При выборе в раскрывающемся списке End conditions условия second появляется возможность задания второй производной сплайна в граничных точках отрезка интерполирования (получающиеся значения первой производной можно лишь смотреть, но не изменять).

Граничные условия natural или 'variational' (естественные)

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

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

Граничные условия periodic (периодические)

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

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

Граничные условия Lagrange

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

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

Граничные условия custom (произвольные)

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

Приближение сглаживающими сплайнами

Интерполяция сплайнами не всегда дает хорошие результаты, особенно в случае зашумленных данных (неточно заданных данных, скажем полученных экспериментально). В этом случае часто прибегают к построению сглаживающего сплайна.

Пример зашумленных данных есть в приложении Spline Tool. Для его загрузки следует перейти в стартовое окно приложения Spline Tool (см. Запуск приложения Spline Tool. Импорт данных.) и нажать кнопку Noisy values of a smooth function. Появляется основное окно приложения Spline Tool, содержащее зашумленные данные.

В данном примере исходными данными являются возмущенные значения функции sin x в 101 точках, равномерно распределенных на отрезке [0.2π]. К значению sin x в каждой точке добавлено случайное число из интервала (-0.1, 0.1), случайные числа распределены равномерно.

Интерполяция зашумленных данных при помощи кубических сплайнов (не важно с какими условиями) не очень удовлетворительна, например зададим граничные условия not-a-knot (задание и смысл граничных условий при интерполяции кубическими сплайнами описано в разделе Приближение кубическими сплайнами).


Интерполяция зашумленных данных кубическими сплайнами.

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

  1. график сплайна должен лежать как можно ближе к заданным точкам;
  2. гладкость получающейся кривой должна быть как можно больше.

У этих требований есть количественные выражения. Обозначим через ƒ(x) сглаживающий сплайн, а через (xk, yk)k=1,2,...,n исходные данные. В качестве меры гладкости сплайна ƒ(x) выбирают интегральную характеристику его производной

где m такое число, что 2m - порядок сплайна, который для кубического сплайна равен 4 (степень минус 1), т.е. для кубического сплайна мерой гладкости является

В качестве меры близости выбирают взвешенное суммарное отклонение сплайна от исходных данных:

Здесь веса wk обычно выбираются как веса обобщенной квадратурной формулы трапеций. В случае равноотстоящих точек с шагом h:

Далее можно задаться некоторой точностью ε и искать такой сглаживающий сплайн, который минимизирует выражение

где параметр ρ выбирается из условия

Приблизим данные сглаживающим сплайном порядка 4 (т.е. кубическим), выбрав ε=0.25. Для этого в приложении Spline Tool в раскрывающемся списке на панели Approximation method выберем Smoothing Spline, убедимся, что в раскрывающемся списке Order выбрано 4 и зададим значение 0.025 в строке ввода Tolerance. В результате получаем гладкое приближение зашумленных данных.


Сглаживающий сплайн и ошибка для ε=0.025

Для ε=0.025 значение параметра ρ получилось равным 0.97064. Отметим, что в данном случае для построения сглаживающего сплайна использовалась функция spaps, входящая в Spline Toolbox, что видно из содержимого строки Bottomlone внизу окна приложения SplineTool.

При задании более малых значений ε сглаживающий сплайн будет все более походить на интерполяционный кубический сплайн и превратится в него при ε=0. Напротив, при увеличении значения ε=0.2 требование близости сплайна и данных будет играть все меньшую роль, зато будет получаться все более гладкий сплайн, например, для сплайн получится таким, как на рисунке ниже (ему соответствует ρ=0.79149).


Сглаживающий сплайн для ε=0.2

При возрастании ε до бесконечности мы получим линейную аппроксимацию данных, т.к. при минимизации выражения

по всевозможным сплайнам ƒ(x) значение ρ будет все ближе к нулю и существенная часть в этом выражении будет определятся значением показателя гладкости сплайна

Но F(ƒ) тем меньше, чем ближе функция к линейной.

Мы выбирали сглаживающий сплайн четвертого порядка. Приложение Spline Tool (равно как и использующаяся в нем функция spaps) позволяет выбрать сглаживающий сплайн шестого порядка, т.е. в качестве меры гладкости будет использоваться

Для построения сглаживающего сплайна шестого порядка в раскрывающемся списке Order следует выбрать 6.

Другой вариант построения сглаживающего сплайна четвертого порядка заключается в задании параметра ρ и поиска сплайна ƒ(x), минимизирующего выражение

Здесь параметр ρ изменяется от 0 до 1, причем ρ=1 соответствует обычному кубическому сплайну, а ρ=0 - аппроксимации данных линейной функцией. Для реализации этого способа построения кубического сплайна в Spline Toolbox следует изменять значение ρ=0 в области ввода Parameter и наблюдать за получающимся сплайном и его отклонением от данных, которое выводится в область Tolerance.

Получающиеся сплайны оказываются достаточно чувствительны к выбору сглаживающего параметра. Приблизим, например, наши зашумленные данные со значением ρ=0.999, т.е. отступив всего на 0.001 от значения, соответствующего интерполяционному кубическому сплайну.


Сглаживающий сплайн для ρ=0.999

При значении параметра ρ=0.9 мы получим уже более плохое приближение.


Сглаживающий сплайн для ρ=0.9

При таком способе построения сглаживающего сплайна использовалась функция csaps, входящая в Spline Toolbox.

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


Поиск по сайту:

Система Orphus

Яндекс.Метрика