Skip to main content

MATLABによる最適化手法と応用

概要​

本文書では、MATLABの最適化機能を活用した数値計算手法について説明します。特に境界制約付き最適化問題とShooting法による微分方程式境界値問題の解法に焦点を当てます。

境界制約付き最適化​

基本概念​

境界制約付き最適化は、変数が特定の範囲内に制限された最適化問題です。実世界の多くの問題では、解は物理的な制約内に収まる必要があります。

fminbnd関数の使用​

MATLABのfminbnd関数は、単変数関数の境界制約付き最適化を効率的に実行します。

例：二次関数の最小化​

% 目的関数を定義
f = @(x) -x.^2 + 5*x - 6;

% 境界制約 [0, 4] 内で最小値を求める
[x_opt, f_opt] = fminbnd(f, 0, 4);

fprintf('最適解: x = %.4f\n', x_opt);
fprintf('最適値: f(x) = %.4f\n', f_opt);
数値解析のポイント

fminbnd関数は黄金分割法（Golden Section Search）に基づいており、関数の単峰性を仮定しています。関数が複数の局所最小値を持つ場合は、初期範囲の設定が重要になります。

アルゴリズムの理解​

黄金分割法は以下の手順で動作します：

  1. 初期区間 [a, b] を設定
  2. 黄金比 φ = (1+√5)/2 を使用して内分点を計算
  3. 内分点での関数値を比較
  4. 収束条件まで区間を縮小
% 黄金分割法の実装例
function [x_min, f_min] = golden_section(f, a, b, tol)
    phi = (1 + sqrt(5)) / 2;
    resphi = 2 - phi;
    
    x1 = a + resphi * (b - a);
    x2 = a + (1 - resphi) * (b - a);
    
    f1 = f(x1);
    f2 = f(x2);
    
    while abs(b - a) > tol
        if f1 > f2
            a = x1;
            x1 = x2;
            f1 = f2;
            x2 = a + (1 - resphi) * (b - a);
            f2 = f(x2);
        else
            b = x2;
            x2 = x1;
            f2 = f1;
            x1 = a + resphi * (b - a);
            f1 = f(x1);
        end
    end
    
    x_min = (a + b) / 2;
    f_min = f(x_min);
end

境界値問題への応用​

以下のような境界値問題を考えます：

% 最適化結果の可視化
x_vals = linspace(0, 4, 100);
y_vals = -x_vals.^2 + 5*x_vals - 6;

figure;
plot(x_vals, y_vals, 'b-', 'LineWidth', 2);
hold on;
plot(x_opt, f_opt, 'ro', 'MarkerSize', 10, 'MarkerFaceColor', 'r');
xlabel('x');
ylabel('f(x)');
title('境界制約付き最適化結果');
legend('f(x) = -x² + 5x - 6', '最適解', 'Location', 'best');
grid on;

微分方程式最適化への応用​

線形回帰問題​

最小二乗法による線形回帰は最適化問題の典型例です。

データフィッティング​

% サンプルデータの生成
x_data = [1, 2, 3, 4, 5];
y_data = [2.1, 3.9, 6.1, 8.2, 9.8];

% 線形モデル y = ax + b の係数を求める
design_matrix = [x_data(:), ones(length(x_data), 1)];
coefficients = design_matrix \ y_data(:);

a = coefficients(1);
b = coefficients(2);

fprintf('回帰直線: y = %.4fx + %.4f\n', a, b);

最小二乗法の数学的背景​

最小二乗法は以下の目的関数を最小化します：

J(a,b)=∑i=1n(yi−axi−b)2J(a,b) = \sum_{i=1}^{n} (y_i - ax_i - b)^2

正規方程式は：

[∑xi2∑xi∑xin][ab]=[∑xiyi∑yi]\begin{bmatrix} \sum x_i^2 & \sum x_i \\ \sum x_i & n \end{bmatrix} \begin{bmatrix} a \\ b \end{bmatrix} = \begin{bmatrix} \sum x_i y_i \\ \sum y_i \end{bmatrix}

Shooting法による境界値問題​

問題設定​

以下の連立微分方程式系を考えます：

{x1′=x2x2′=−x1−3x2+e−5tx3′=x4x4′=2x1−3x2−3x3−4x4−sin⁡(t)\begin{cases} x_1' = x_2 \\ x_2' = -x_1 - 3x_2 + e^{-5t} \\ x_3' = x_4 \\ x_4' = 2x_1 - 3x_2 - 3x_3 - 4x_4 - \sin(t) \end{cases}

境界条件：

  • x1(0)=1,x2(0)=2x_1(0) = 1, x_2(0) = 2
  • x3(10)=−0.021677,x4(10)=0.15797x_3(10) = -0.021677, x_4(10) = 0.15797

Shooting法の実装​

Shooting法は境界値問題を初期値問題に変換する手法です：

function z = c3mode(y, f, x0, xn)
    % y: 未知の初期条件 [x3(0), x4(0)]
    % f: 微分方程式の右辺
    % x0: 既知の初期条件 [x1(0), x2(0)]
    % xn: 終端境界条件 [x3(10), x4(10)]

    x0_full = [x0(1:2); y];
    [t, x] = ode45(f, [0, 10], x0_full);
    z = norm([x(end, 3:4)]' - xn);
end

最適化による解法​

% 微分方程式の定義
f = @(t,x) [x(2); 
           -x(1) - 3*x(2) + exp(-5*t); 
            x(4); 
            2*x(1) - 3*x(2) - 3*x(3) - 4*x(4) - sin(t)];

% 既知の境界条件
x0 = [1; 2]; 
xn = [-0.021677; 0.15797];

% 目的関数の設定
g = @(x) c3mode(x, f, x0, xn);

% 初期推定値
x2 = rand(2, 1);

% 最適化の実行
x3 = fminunc(g, x2);

fprintf('最適解: x3(0) = %.4f, x4(0) = %.4f\n', x3(1), x3(2));
計算上の注意

Shooting法では初期推定値の選択が重要です。不適切な初期値は収束しない場合があります。また、数値積分の精度も解の品質に大きく影響します。

解の可視化と検証​

% 最適解による微分方程式の解
[t, x] = ode45(f, [0, 10], [x0; x3]);

% 結果のプロット
figure;
plot(t, x, 'LineWidth', 2);
legend('x_1', 'x_2', 'x_3', 'x_4', 'Location', 'best');
xlabel('時間 t');
ylabel('状態変数');
title('Shooting法による境界値問題の解');
grid on;

% 境界条件の確認
fprintf('終端値: x3(10) = %.6f, x4(10) = %.6f\n', x(end,3), x(end,4));
fprintf('目標値: x3(10) = %.6f, x4(10) = %.6f\n', xn(1), xn(2));

アルゴリズムの改良​

より高精度な解を得るために以下の改良が可能です：

  1. 適応的ステップサイズ制御
  2. 高次のRunge-Kutta法の使用
  3. 複数の初期推定値からの実行
% 改良版Shooting法
function [x_opt, residual] = improved_shooting(f, x0, xn, t_span, n_trials)
    best_residual = inf;
    x_opt = [];
    
    for i = 1:n_trials
        % ランダム初期値
        x_init = randn(2, 1);
        
        % 最適化オプションの設定
        options = optimoptions('fminunc', ...
            'Display', 'off', ...
            'TolFun', 1e-12, ...
            'TolX', 1e-12);
        
        % 目的関数
        objective = @(x) shooting_objective(x, f, x0, xn, t_span);
        
        try
            [x_trial, residual_trial] = fminunc(objective, x_init, options);
            
            if residual_trial < best_residual
                best_residual = residual_trial;
                x_opt = x_trial;
                residual = residual_trial;
            end
        catch
            continue;
        end
    end
end

function residual = shooting_objective(y, f, x0, xn, t_span)
    x0_full = [x0; y];
    [~, x] = ode45(f, t_span, x0_full);
    residual = norm(x(end, 3:4)' - xn);
end

実践的な応用例​

制御系設計への応用​

Shooting法は制御系の最適軌道生成にも応用できます：

% 最適制御問題の例
% 状態方程式: dx/dt = Ax + Bu
% 性能指標: J = ∫(x'Qx + u'Ru)dt

A = [0 1; -2 -3];
B = [0; 1];
Q = eye(2);
R = 1;

% 終端条件を満たす制御入力の設計
% （詳細な実装は制御理論の知識が必要）

数値解析の精度向上​

解の精度を向上させるための技法：

  1. Richardson外挿法
  2. 適応的メッシュ細分化
  3. 誤差推定と制御
% Richardson外挿法の例
function improved_solution = richardson_extrapolation(f, x0, t_span, h_coarse)
    % 粗いステップサイズでの解
    [t1, x1] = ode45(f, t_span, x0, ...
        odeset('MaxStep', h_coarse));
    
    % 細かいステップサイズでの解
    [t2, x2] = ode45(f, t_span, x0, ...
        odeset('MaxStep', h_coarse/2));
    
    % Richardson外挿（高次精度解の推定）
    % 実装の詳細は数値解析理論に基づく
end

まとめ​

本文書では以下の内容を学習しました：

  1. 境界制約付き最適化の基本概念とfminbnd関数の使用法
  2. 黄金分割法のアルゴリズムと実装
  3. Shooting法による境界値問題の解法
  4. 最小二乗法による線形回帰の最適化
  5. 数値解析における精度向上技法
学習のポイント
  • 最適化問題の定式化能力を身につける
  • 数値計算アルゴリズムの背景理論を理解する
  • 実問題への応用力を養う
  • 計算精度と効率のバランスを考慮する

参考文献​

  1. Nocedal, J., & Wright, S. J. (2006). Numerical Optimization. Springer.
  2. Ascher, U. M., & Petzold, L. R. (1998). Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations. SIAM.
  3. MathWorks Documentation: Optimization Toolbox
  4. Press, W. H., et al. (2007). Numerical Recipes: The Art of Scientific Computing. Cambridge University Press.

演習問題​

  1. 異なる境界条件でShooting法を適用し、解の変化を観察せよ
  2. 多変数最適化問題に拡張した場合の計算複雑度を分析せよ
  3. ノイズを含むデータに対する最小二乗法の頑健性を評価せよ
  4. 非線形境界値問題への応用例を考案し、実装せよ
発展学習

より高度な最適化手法（遺伝的アルゴリズム、粒子群最適化、シミュレーテッドアニーリングなど）への展開も考慮してください。