勾配およびヘッシアンを使った最小化

この例では、明示的な三重対角ヘッセ行列 $H (x)$ を使用して非線形最小化問題を解く方法を説明します。この問題は以下を最小化する $x$ を見つけるものです。

$f (x) = \sum_{i = 1}^{n - 1} ({(x_{i}^{2})}^{(x_{i + 1}^{2} + 1)} + {(x_{i + 1}^{2})}^{(x_{i}^{2} + 1)}),$

ここで $n$ = 1000 です。

補助関数 brownfgh (この例の終わりに掲載) が、 $f (x)$ 、その勾配 $g (x)$ 、およびそのヘッシアン $H (x)$ を計算します。fminunc ソルバーが導関数情報を使用するよう指定するには、optimoptions を使用して SpecifyObjectiveGradient オプションと HessianFcn オプションを設定します。fminunc でヘッシアンを使用するには、'trust-region' アルゴリズムを使用しなければなりません。

options = optimoptions(@fminunc,'Algorithm','trust-region',...
    'SpecifyObjectiveGradient',true,'HessianFcn','objective');

パラメーター n を 1000 に設定し、初期点 xstart として奇数の要素に対して -1 を設定し、偶数の要素に対して +1 を設定します。

n = 1000;
xstart = -ones(n,1);
xstart(2:2:n) = 1;

$f$ の最小値を求めます。

[x,fval,exitflag,output] = fminunc(@brownfgh,xstart,options);

Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.

解法および解法プロセスを検証します。

disp(fval)

   2.8709e-17

disp(exitflag)

disp(output)

         iterations: 7
          funcCount: 8
           stepsize: 0.0039
       cgiterations: 7
      firstorderopt: 4.7948e-10
          algorithm: 'trust-region'
            message: 'Local minimum found....'
    constrviolation: []

関数 $f (x)$ は、二乗の和であるため、非負となります。解 fval はほぼゼロであるため、最小値であることは明らかです。終了フラグが 1 であることも fminunc が解を見つけたことを示しています。output 構造体は、fminunc がわずか 7 回の反復で解に到達したことを示しています。

解の最大要素と最小要素を表示します。

disp(max(x))

   1.1987e-10

disp(min(x))

  -1.1987e-10

この解は、すべての要素が x = 0 となる点に極めて近いものとなります。

補助関数

次のコードは、補助関数 brownfgh を作成します。

function [f,g,H] = brownfgh(x)
%BROWNFGH  Nonlinear minimization problem (function, its gradients
% and Hessian)
% Documentation example        

%   Copyright 1990-2008 The MathWorks, Inc.

% Evaluate the function.
  n=length(x); y=zeros(n,1);
  i=1:(n-1);
  y(i)=(x(i).^2).^(x(i+1).^2+1)+(x(i+1).^2).^(x(i).^2+1);
  f=sum(y);
%
% Evaluate the gradient.
  if nargout > 1
     i=1:(n-1); g = zeros(n,1);
     g(i)= 2*(x(i+1).^2+1).*x(i).*((x(i).^2).^(x(i+1).^2))+...
             2*x(i).*((x(i+1).^2).^(x(i).^2+1)).*log(x(i+1).^2);
     g(i+1)=g(i+1)+...
              2*x(i+1).*((x(i).^2).^(x(i+1).^2+1)).*log(x(i).^2)+...
              2*(x(i).^2+1).*x(i+1).*((x(i+1).^2).^(x(i).^2));
  end
%
% Evaluate the (sparse, symmetric) Hessian matrix
  if nargout > 2
     v=zeros(n,1);
     i=1:(n-1);
     v(i)=2*(x(i+1).^2+1).*((x(i).^2).^(x(i+1).^2))+...
            4*(x(i+1).^2+1).*(x(i+1).^2).*(x(i).^2).*((x(i).^2).^((x(i+1).^2)-1))+...
            2*((x(i+1).^2).^(x(i).^2+1)).*(log(x(i+1).^2));
     v(i)=v(i)+4*(x(i).^2).*((x(i+1).^2).^(x(i).^2+1)).*((log(x(i+1).^2)).^2);
     v(i+1)=v(i+1)+...
              2*(x(i).^2).^(x(i+1).^2+1).*(log(x(i).^2))+...
              4*(x(i+1).^2).*((x(i).^2).^(x(i+1).^2+1)).*((log(x(i).^2)).^2)+...
              2*(x(i).^2+1).*((x(i+1).^2).^(x(i).^2));
     v(i+1)=v(i+1)+4*(x(i).^2+1).*(x(i+1).^2).*(x(i).^2).*((x(i+1).^2).^(x(i).^2-1));
     v0=v;
     v=zeros(n-1,1);
     v(i)=4*x(i+1).*x(i).*((x(i).^2).^(x(i+1).^2))+...
            4*x(i+1).*(x(i+1).^2+1).*x(i).*((x(i).^2).^(x(i+1).^2)).*log(x(i).^2);
     v(i)=v(i)+ 4*x(i+1).*x(i).*((x(i+1).^2).^(x(i).^2)).*log(x(i+1).^2);
     v(i)=v(i)+4*x(i).*((x(i+1).^2).^(x(i).^2)).*x(i+1);
     v1=v;
     i=[(1:n)';(1:(n-1))'];
     j=[(1:n)';(2:n)'];
     s=[v0;2*v1];
     H=sparse(i,j,s,n,n);
     H=(H+H')/2;
  end
end

勾配およびヘッシアンを使った最小化

補助関数

関連するトピック