이 페이지의 최신 내용은 아직 번역되지 않았습니다. 최신 내용은 영문으로 볼 수 있습니다.

비선형 데이터 피팅

이 예제에서는 여러 Optimization Toolbox™ 알고리즘을 사용하여 비선형 함수를 데이터에 피팅하는 방법을 보여줍니다.

문제 설정

다음과 같은 데이터가 있다고 가정하겠습니다.

Data = ...
  [0.0000    5.8955
   0.1000    3.5639
   0.2000    2.5173
   0.3000    1.9790
   0.4000    1.8990
   0.5000    1.3938
   0.6000    1.1359
   0.7000    1.0096
   0.8000    1.0343
   0.9000    0.8435
   1.0000    0.6856
   1.1000    0.6100
   1.2000    0.5392
   1.3000    0.3946
   1.4000    0.3903
   1.5000    0.5474
   1.6000    0.3459
   1.7000    0.1370
   1.8000    0.2211
   1.9000    0.1704
   2.0000    0.2636];

이 데이터 점을 플로팅해 보겠습니다.

t = Data(:,1);
y = Data(:,2);
% axis([0 2 -0.5 6])
% hold on
plot(t,y,'ro')
title('Data points')

% hold off

다음 함수를

y = c(1)*exp(-lam(1)*t) + c(2)*exp(-lam(2)*t)

데이터에 피팅해 보겠습니다.

lsqcurvefit을 사용하는 풀이법

lsqcurvefit 함수는 이러한 유형의 문제를 손쉽게 풉니다.

먼저, 하나의 변수 x로 파라미터를 정의합니다.

x(1) = c(1)

x(2) = lam(1)

x(3) = c(2)

x(4) = lam(2)

그런 다음 곡선을 파라미터 x와 데이터 t로 구성된 함수로 정의합니다.

F = @(x,xdata)x(1)*exp(-x(2)*xdata) + x(3)*exp(-x(4)*xdata);

다음과 같이 c(1) = 1, lam(1) = 1, c(2) = 1, lam(2) = 0으로 하여 초기점 x0을 임의로 설정합니다.

x0 = [1 1 1 0];

솔버를 실행하고 결과로 생성되는 피팅을 플로팅합니다.

[x,resnorm,~,exitflag,output] = lsqcurvefit(F,x0,t,y)
Local minimum possible.

lsqcurvefit stopped because the final change in the sum of squares relative to 
its initial value is less than the value of the function tolerance.
x = 1×4

    3.0068   10.5869    2.8891    1.4003

resnorm = 0.1477
exitflag = 3
output = struct with fields:
    firstorderopt: 7.8881e-06
       iterations: 6
        funcCount: 35
     cgiterations: 0
        algorithm: 'trust-region-reflective'
         stepsize: 0.0096
          message: '...'

hold on
plot(t,F(x,t))
hold off

fminunc를 사용하는 풀이법

fminunc를 사용하여 문제를 풀려면 목적 함수를 잔차의 제곱에 대한 합으로 설정해야 합니다.

Fsumsquares = @(x)sum((F(x,t) - y).^2);
opts = optimoptions('fminunc','Algorithm','quasi-newton');
[xunc,ressquared,eflag,outputu] = ...
    fminunc(Fsumsquares,x0,opts)
Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.
xunc = 1×4

    2.8890    1.4003    3.0069   10.5862

ressquared = 0.1477
eflag = 1
outputu = struct with fields:
       iterations: 30
        funcCount: 185
         stepsize: 0.0017
     lssteplength: 1
    firstorderopt: 2.9476e-05
        algorithm: 'quasi-newton'
          message: '...'

fminunclsqcurvefit과 동일한 해를 구했지만 함수 실행 횟수가 훨씬 많았음을 알 수 있습니다. fminunc의 파라미터는 lsqcurvefit의 파라미터와 순서가 반대입니다. 즉, lam(1)이 아니라 lam(2)가 더 큰 lam입니다. 변수의 순서는 임의로 지정 가능하므로 이는 놀랍지 않은 일입니다.

fprintf(['There were %d iterations using fminunc,' ...
    ' and %d using lsqcurvefit.\n'], ...
    outputu.iterations,output.iterations)
There were 30 iterations using fminunc, and 6 using lsqcurvefit.
fprintf(['There were %d function evaluations using fminunc,' ...
    ' and %d using lsqcurvefit.'], ...
    outputu.funcCount,output.funcCount)
There were 185 function evaluations using fminunc, and 35 using lsqcurvefit.

선형 문제와 비선형 문제로 분리하기

이 피팅 문제는 파라미터 c(1)과 c(2)에서 선형 문제임을 알 수 있습니다. 이는 lam(1)과 lam(2)의 모든 값에 대해 백슬래시 연산자를 사용하여 c(1)과 c(2)의 값을 구하여 최소제곱 문제를 풀 수 있음을 의미합니다.

이제 이 문제를 2차원 문제로 다시 구성하여 lam(1)과 lam(2)에 대한 최적의 값을 찾아보겠습니다. c(1)과 c(2)의 값은 위에 설명한 대로 각 단계에서 백슬래시 연산자를 사용하여 계산됩니다.

type fitvector
function yEst = fitvector(lam,xdata,ydata)
%FITVECTOR  Used by DATDEMO to return value of fitting function.
%   yEst = FITVECTOR(lam,xdata) returns the value of the fitting function, y
%   (defined below), at the data points xdata with parameters set to lam.
%   yEst is returned as a N-by-1 column vector, where N is the number of
%   data points.
%
%   FITVECTOR assumes the fitting function, y, takes the form
%
%     y =  c(1)*exp(-lam(1)*t) + ... + c(n)*exp(-lam(n)*t)
%
%   with n linear parameters c, and n nonlinear parameters lam.
%
%   To solve for the linear parameters c, we build a matrix A
%   where the j-th column of A is exp(-lam(j)*xdata) (xdata is a vector).
%   Then we solve A*c = ydata for the linear least-squares solution c,
%   where ydata is the observed values of y.

A = zeros(length(xdata),length(lam));  % build A matrix
for j = 1:length(lam)
   A(:,j) = exp(-lam(j)*xdata);
end
c = A\ydata; % solve A*c = y for linear parameters c
yEst = A*c; % return the estimated response based on c

lsqcurvefit을 사용하여 2차원 초기점 lam(1), lam(2)부터 시작해 문제를 풉니다.

x02 = [1 0];
F2 = @(x,t) fitvector(x,t,y);

[x2,resnorm2,~,exitflag2,output2] = lsqcurvefit(F2,x02,t,y)
Local minimum possible.

lsqcurvefit stopped because the final change in the sum of squares relative to 
its initial value is less than the value of the function tolerance.
x2 = 1×2

   10.5861    1.4003

resnorm2 = 0.1477
exitflag2 = 3
output2 = struct with fields:
    firstorderopt: 4.4107e-06
       iterations: 10
        funcCount: 33
     cgiterations: 0
        algorithm: 'trust-region-reflective'
         stepsize: 0.0080
          message: '...'

2차원 풀이법과 4차원 풀이법의 효율성은 서로 비슷합니다.

fprintf(['There were %d function evaluations using the 2-d ' ...
    'formulation, and %d using the 4-d formulation.'], ...
    output2.funcCount,output.funcCount)
There were 33 function evaluations using the 2-d formulation, and 35 using the 4-d formulation.

분리된 문제가 초기 추측값에 대해 더 견고함

4개의 파라미터를 갖는 원래 문제에서 잘못된 시작점을 선택하면 전역해가 아닌 국소해를 찾습니다. 2개의 파라미터로 분리된 문제에서 똑같이 잘못된 lam(1) 값과 lam(2) 값을 갖는 시작점을 선택하는 경우에는 전역해를 찾습니다. 이를 확인하기 위해 상대적으로 나쁜 국소해를 생성하게 하는 시작점을 사용하여 원래 문제를 다시 실행한 후 결과로 생성되는 피팅을 전역해와 비교해 보겠습니다.

x0bad = [5 1 1 0];
[xbad,resnormbad,~,exitflagbad,outputbad] = ...
    lsqcurvefit(F,x0bad,t,y)
Local minimum possible.

lsqcurvefit stopped because the final change in the sum of squares relative to 
its initial value is less than the value of the function tolerance.
xbad = 1×4

  -20.2519    2.4796   25.3757    2.4795

resnormbad = 2.2173
exitflagbad = 3
outputbad = struct with fields:
    firstorderopt: 0.0036
       iterations: 31
        funcCount: 160
     cgiterations: 0
        algorithm: 'trust-region-reflective'
         stepsize: 0.0012
          message: '...'

hold on
plot(t,F(xbad,t),'g')
legend('Data','Global fit','Bad local fit','Location','NE')
hold off

fprintf(['The residual norm at the good ending point is %f,' ...
   ' and the residual norm at the bad ending point is %f.'], ...
   resnorm,resnormbad)
The residual norm at the good ending point is 0.147723, and the residual norm at the bad ending point is 2.217300.

관련 항목