Main Content

이 번역 페이지는 최신 내용을 담고 있지 않습니다. 최신 내용을 영문으로 보려면 여기를 클릭하십시오.

3차 평활화 스플라인

이 예제에서는 Curve Fitting Toolbox™의 csapsspaps 명령을 사용하여 3차 평활화 스플라인을 생성하는 방법을 보여줍니다.

CSAPS 명령

csaps 명령은 평활화 스플라인을 제공합니다. 평활화 스플라인은 잡음이 있는 데이터의 추정되는 기본 추세를 대략적으로 따르는 3차 스플라인입니다. 사용자가 선택하는 평활화 파라미터는 평활화 스플라인이 주어진 데이터를 얼마나 근접하게 따라갈지를 결정합니다. 다음 기본 정보에 문서의 내용이 요약되어 있습니다.

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

예: 3차 다항식의 잡음이 있는 데이터

다음은 몇 차례의 시험적인 실행입니다. 간단한 3차식 q(x) := x^3의 데이터로 시작해서 잡음으로 데이터를 오염시키고, 평활화 파라미터의 값을 .5로 선택합니다. 그런 다음 결과로 생성되는 평활화된 값을 기본 3차식 및 오염된 데이터와 함께 플로팅합니다.

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이 되면 (잡음이 있는) 데이터에 대한 보간을 얻게 됩니다.

실제로 csapi에서 사용하는 정식화(A Practical Guide to Splines, p.235ff)는 독립 변수의 스케일링에 매우 민감합니다. 사용된 방정식을 간단히 분석해 보면 p에 대한 민감한 범위가 약 1/(1+epsilon)임을 알 수 있습니다. 여기서 epsilon := h^3/16이고, h는 인접 지점 간 차이의 평균입니다. 구체적으로 말해 p = 1/(1+epsilon/100)이면 데이터를 밀접하게 따라가는 반면, p = 1/(1+epsilon*100)이면 평활화가 어느 정도 만족스러울 것으로 예상할 수 있습니다.

아래 플롯은 이 매직 넘버 1/(1+epsilon)에 가까운 p 값에 대한 평활화 스플라인을 보여줍니다. 여기서는, 매직 넘버 1/(1+epsilon)이 1에 매우 가까우므로 1-p를 살펴보는 것이 보다 유용합니다.

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 명령

spaps 명령이 제공하는 3차 평활화 스플라인과 csaps에서 생성되는 3차 평활화 스플라인은 스플라인이 선택되는 방식만 다릅니다. 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)).

허용오차와 평활화 파라미터 비교

csaps에 필요한 평활화 파라미터 p를 제공하는 것보다 spaps에 대해 적합한 허용오차를 제공하는 것이 더 쉬울 수 있습니다. 이전 예제에서는 구간 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' )

Figure 제목에는 이러한 데이터에 대해 이 평활화 스플라인과 정확히 똑같은 결과를 얻기 위해 csaps에서 사용할 p의 값이 표시됩니다.

다음은 평활화 파라미터가 주어지지 않은 경우 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 비교

csaps 명령과 spaps 명령은 구체적인 평활화 스플라인을 지정하는 방식(평활화 파라미터 대 허용오차)이 다릅니다. 또 다른 차이점은 spaps는 3차 평활화 스플라인 외에도 1차 또는 5차 평활화 스플라인을 제공할 수 있다는 것입니다.

5차 평활화 스플라인은 2계 도함수가 가능한 한 적게 움직이면 좋을 상황에서 3차 평활화 스플라인보다 유용합니다.