Javaでのニュートン法による最小二乗法の実装

Javaでのニュートン法による最小二乗法の実装

本記事では、Javaを用いてニュートン法を利用した最小二乗法の実装について解説します。数式はMathJaxを用いて記述し、具体例を多く含めています。以下のリンクから各セクションにジャンプできます。

最小二乗法とは

最小二乗法は、観測データに対してモデルをフィッティングする方法の一つです。具体的には、次のような式で定義されます。

観測データ \( (x_i, y_i) \) に対し、モデル \( f(x; \boldsymbol{\theta}) \) が与えられるとします。このとき、最小二乗法は次の目的関数を最小化します:

\[ S(\boldsymbol{\theta}) = \sum_{i=1}^n \left( y_i – f(x_i; \boldsymbol{\theta}) \right)^2 \]

\( \boldsymbol{\theta} \) はモデルパラメータであり、これを最適化することで誤差を最小化します。

ニュートン法とは

ニュートン法は、非線形方程式の解を反復的に求める方法です。関数 \( g(\boldsymbol{\theta}) \) の根を見つける場合、次の更新式を用います:

\[ \boldsymbol{\theta}_{k+1} = \boldsymbol{\theta}_k – \mathbf{H}^{-1} \nabla g(\boldsymbol{\theta}_k) \]

ここで、\( \nabla g(\boldsymbol{\theta}) \) は勾配ベクトル、\( \mathbf{H} \) はヘッセ行列(2次微分行列)です。

最小二乗法においては、目的関数の勾配とヘッセ行列を次のように計算します:

  • 勾配: \[ \nabla S(\boldsymbol{\theta}) = -2 \sum_{i=1}^n \left( y_i – f(x_i; \boldsymbol{\theta}) \right) \nabla f(x_i; \boldsymbol{\theta}) \]
  • ヘッセ行列: \[ \mathbf{H} = 2 \sum_{i=1}^n \nabla f(x_i; \boldsymbol{\theta}) \nabla f(x_i; \boldsymbol{\theta})^\top – 2 \sum_{i=1}^n \left( y_i – f(x_i; \boldsymbol{\theta}) \right) \nabla^2 f(x_i; \boldsymbol{\theta}) \]

Javaによる実装

ここでは、Javaでニュートン法を実装し、最小二乗法を用いてモデルのパラメータを最適化する方法を示します。


// 必要なライブラリのインポート
import java.util.Arrays;

public class NewtonLeastSquares {
    // モデル関数 f(x; theta)
    public static double model(double x, double[] theta) {
        return theta[0] * x + theta[1]; // 例: 線形モデル y = theta0 * x + theta1
    }

    // 目的関数 S(theta) の勾配を計算
    public static double[] gradient(double[][] data, double[] theta) {
        double[] grad = new double[theta.length];
        Arrays.fill(grad, 0.0);
        for (double[] point : data) {
            double x = point[0];
            double y = point[1];
            double error = y - model(x, theta);
            grad[0] -= 2 * error * x;
            grad[1] -= 2 * error;
        }
        return grad;
    }

    // ヘッセ行列を計算
    public static double[][] hessian(double[][] data, double[] theta) {
        double[][] hess = new double[theta.length][theta.length];
        for (double[] point : data) {
            double x = point[0];
            double y = point[1];
            hess[0][0] += 2 * x * x;
            hess[0][1] += 2 * x;
            hess[1][0] += 2 * x;
            hess[1][1] += 2;
        }
        return hess;
    }

    // ニュートン法による更新
    public static double[] newtonStep(double[][] data, double[] theta) {
        double[] grad = gradient(data, theta);
        double[][] hess = hessian(data, theta);
        double[] step = new double[theta.length];
        // ヘッセ行列の逆行列を計算(簡略化のため固定サイズ)
        double det = hess[0][0] * hess[1][1] - hess[0][1] * hess[1][0];
        double[][] invHess = {
            { hess[1][1] / det, -hess[0][1] / det },
            { -hess[1][0] / det, hess[0][0] / det }
        };
        for (int i = 0; i < theta.length; i++) {
            for (int j = 0; j < theta.length; j++) {
                step[i] -= invHess[i][j] * grad[j];
            }
        }
        for (int i = 0; i < theta.length; i++) {
            theta[i] += step[i];
        }
        return theta;
    }

    public static void main(String[] args) {
        // 観測データ (x, y)
        double[][] data = {
            {1, 2}, {2, 4}, {3, 6}, {4, 8}, {5, 10}
        };
        // 初期パラメータ
        double[] theta = {0.0, 0.0}; // [theta0, theta1]
        // 反復計算
        for (int i = 0; i < 10; i++) {
            theta = newtonStep(data, theta);
            System.out.println("Iteration " + (i + 1) + ": " + Arrays.toString(theta));
        }
    }
}
    

具体例

上記のコードを実行すると、次のような結果が得られます。これは、線形モデル \( y = \theta_0 x + \theta_1 \) に対してニュートン法を適用した例です。

観測データが \( \{(1, 2), (2, 4), (3, 6), (4, 8), (5, 10)\} \) の場合、最終的なパラメータは次のようになります:

  • \( \theta_0 \approx 2.0 \)
  • \( \theta_1 \approx 0.0 \)

この結果は、観測データが完全に直線 \( y = 2x \) に一致することを意味します。

まとめ

本記事では、Javaを用いてニュートン法による最小二乗法を実装する方法について解説しました。ニュートン法を利用することで、高速かつ精度の高いパラメータ推定が可能です。ただし、ヘッセ行列の逆行列計算や収束性に注意が必要です。これを応用すれば、より複雑なモデルへのフィッティングも実現できます。

コメントを残す

メールアドレスが公開されることはありません。 が付いている欄は必須項目です