Этот пример показывает, как использовать csaps
и команды spaps
от Curve Fitting Toolbox™, чтобы создать кубические сплайны сглаживания.
Команда csaps
обеспечивает сплайн сглаживания. Это - кубический сплайн, который более или менее следует за предполагаемым базовым трендом в шумных данных. Параметр сглаживания, чтобы быть выбранным вами, определяет, как тесно сплайн сглаживания следует за определенными данными. Вот основная информация, сокращенная версия документации:
CSAPS Cubic smoothing spline.
VALUES = CSAPS(X, Y, P, XX)
Returns the values at XX of the cubic smoothing spline for the
given data (X,Y) and depending on the smoothing parameter P, chosen
from the interval [0 .. 1]. This smoothing spline f minimizes
P * sum_i W(i)(Y(i) - f(X(i)))^2 + (1-P) * integral (D^2 f)^2
Вот некоторые пробные прогоны. Мы запускаем с данных из простого кубического, q(x) := x^3
, загрязняем значения некоторым шумом и выбираем значение параметра сглаживания, чтобы быть.5. Затем постройте получившиеся сглаживавшие значения, наряду с базовым кубическим, и загрязненные данные.
xi = (0:.05:1); q = @(x) x.^3; yi = q(xi); randomStream = RandStream.create( 'mcg16807', 'Seed', 23 ); ybad = yi+.3*(rand(randomStream, size(xi))-.5); p = .5; xxi = (0:100)/100; ys = csaps(xi,ybad,p,xxi); plot(xi,yi,':',xi,ybad,'x',xxi,ys,'r-') title('Clean Data, Noisy Data, Smoothed Values') legend( 'Exact', 'Noisy', 'Smoothed', 'Location', 'NorthWest' )
Сглаживание является путем, преувеличенным здесь. Путем выбора параметра сглаживания p
ближе к 1, мы получаем сплайн сглаживания ближе к определенным данным. Мы пробуем p = .6, .7, .8, .9, 1
и строим получившиеся сплайны сглаживания.
yy = zeros(5,length(xxi)); p = [.6 .7 .8 .9 1]; for j=1:5 yy(j,:) = csaps(xi,ybad,p(j),xxi); end hold on plot(xxi,yy); hold off title('Smoothing Splines for Various Values of the Smoothing Parameter') legend({'Exact','Noisy','p = 0.5','p = 0.6','p = 0.7','p = 0.8', ... 'p = 0.9', 'p = 1.0'}, 'Location', 'NorthWest' )
Мы видим, что сплайн сглаживания может быть очень чувствителен к выбору параметра сглаживания. Даже для p
= 0.9, сплайн сглаживания все еще далек от базового тренда, в то время как для p
= 1, мы получаем interpolant к (шумным) данным.
На самом деле формулировка, используемая csapi
(p.235ff Практического Руководства по Сплайнам), очень чувствительна к масштабированию независимой переменной. Простой анализ используемых уравнений показывает, что чувствительная область значений для p
вокруг 1/(1+epsilon)
с epsilon := h^3/16
и h
средняя разница между соседними сайтами. А именно, вы ожидали бы завершение после данных когда p = 1/(1+epsilon/100)
и некоторое удовлетворительное сглаживание когда p = 1/(1+epsilon*100)
.
График ниже показов сплайн сглаживания для значений p
около этого магического числа 1/(1+epsilon)
. Для этого случая это более информативно, чтобы посмотреть на 1-p
, поскольку магическое число, 1/(1+epsilon)
, очень близко к 1.
epsilon = ((xi(end)-xi(1))/(numel(xi)-1))^3/16; 1 - 1/(1+epsilon)
ans = 7.8124e-06
plot(xi,yi,':',xi,ybad,'x') hold on labels = cell(1,5); for j=1:5 p = 1/(1+epsilon*10^(j-3)); yy(j,:) = csaps(xi,ybad,p,xxi); labels{j} = ['1-p= ',num2str(1-p)]; end plot(xxi,yy) title('Smoothing Splines for Smoothing Parameter Near Its ''Magic'' Value') legend( [{'Exact', 'Noisy'}, labels], 'Location', 'NorthWest' ) hold off
В этом примере сплайн сглаживания очень чувствителен к изменению параметра сглаживания около магического числа. Одно самое дальнее от 1 кажется лучшим выбором от них, но можно предпочесть тот кроме того.
p = 1/(1+epsilon*10^3); yy = csaps(xi,ybad,p,xxi); hold on plot( xxi, yy, 'y', 'LineWidth', 2 ) title( sprintf( 'The Smoothing Spline For 1-p = %s is Added, in Yellow', num2str(1-p) ) ) hold off
Можно также предоставить csaps
ошибочные веса, чтобы уделить больше внимания некоторым точкам данных, чем другие. Кроме того, если вы не предоставляете сайты оценки xx
, затем csaps
возвращает ppform сплайна сглаживания.
Наконец, csaps
может также обработать данные с векторным знаком и даже многомерные, данные с координатной сеткой.
Кубический сплайн сглаживания, обеспеченный командой, spaps
отличается от того, созданного в csaps
только в способе, которым это выбрано. Вот сокращенная версия документации для spaps
:
SPAPS Smoothing spline.
[SP,VALUES] = SPAPS(X,Y,TOL) returns the B-form and, if asked,
the values at X, of a cubic smoothing spline f for the given
data (X(i),Y(:,i)), i=1,2, ..., n.
The smoothing spline f minimizes the roughness measure
F(D^2 f) := integral ( D^2 f(t) )^2 dt on X(1) < t < X(n)
over all functions f for which the error measure
E(f) := sum_j { W(j)*( Y(:,j) - f(X(j)) )^2 : j=1,...,n }
is no bigger than the given TOL. Here, D^M f denotes the M-th
derivative of f. The weights W are chosen so that E(f) is the
Composite Trapezoid Rule approximation for F(y-f).
f is constructed as the unique minimizer of
rho*E(f) + F(D^2 f),
with the smoothing parameter RHO so chosen so that E(f) equals
TOL. Hence, FN2FM(SP,'pp') should be (up to roundoff) the same
as the output from CPAPS(X,Y,RHO/(1+RHO)).
Может быть легче предоставить подходящий допуск к spaps
, чем параметр сглаживания p
, требуемый csaps
. В нашем более раннем примере мы добавили равномерно распределенный случайный шум от интервала 0.3* [-0.5.. 0.5]. Следовательно, мы можем оценить рыночную стоимость для tol
, когда значение ошибки измеряет e
в таком шуме.
tol = sum((.3*(rand(randomStream, size(yi))-.5)).^2);
Этот график показывает получившийся сплайн сглаживания, созданный spaps
. Обратите внимание на то, что ошибочные веса заданы, чтобы быть универсальными, который является их значением по умолчанию в csaps
.
[sp,ys,rho] = spaps(xi,ybad,tol,ones(size(xi))); plot(xi,yi,':',xi,ybad,'x',xi,ys,'r-') title( sprintf( 'Clean Data, Noisy Data, Smoothed Values (1-p = %s )', num2str(1/(1+rho)) ) ); legend( {'Exact','Noisy','Smoothed'}, 'location', 'NorthWest' )
Заголовок фигуры показывает значение p
, который вы использовали бы в csaps
, чтобы получить точно этот сплайн сглаживания для этих данных.
Здесь, кроме того, сплайн сглаживания, обеспеченный csaps
если не, учитывая параметр сглаживания. В этом случае csaps
выбирает параметр определенной оперативной процедурой, которая пытается определить местоположение области, где сплайн сглаживания является самым чувствительным к параметру сглаживания (подобный более раннему обсуждению).
hold on plot(xxi,fnval(csaps(xi,ybad),xxi),'-') title('Clean Data, Noisy Data, Smoothed Values') legend({'Exact' 'Noisy' 'spaps, specified tolerance' ... 'csaps, default smoothing parameter'}, 'Location', 'NorthWest' ) hold off
csaps
и команды spaps
отличаются по пути, которым вы задаете конкретный сплайн сглаживания через параметр сглаживания по сравнению с допуском. Другое различие - то, что spaps
может обеспечить линейное или quintic, сглаживающий сплайн, в дополнение к кубическому сплайну сглаживания.
quintic, сглаживающий сплайн, лучше, чем кубический сплайн сглаживания в ситуации, когда вы хотели бы, чтобы вторая производная переместилась как можно меньше.