最小二乗法による最適化をしたい場面に使えるscipyの関数を一部まとめています。原理的にはscikit-learnで回帰分析やscipyのカーブフィットを行うことと変わりないですが、扱う係数行列が巨大であったり、スパースであったり、もしくは明示的に最小二乗法をときたい場合には以下の関数が使えます。
なお、この記事では全て重みなし最小二乗法を扱っています。
最小二乗法の原理
いま、最小二乗法における目的関数を以下のように定義します。
min \Sigma_i e^2_i = min \Sigma_i (y_i - f(x_i; \beta))^2 = min |Y - X \beta|^2
ここで、$x$と$y$はそれぞれ説明変数と被説明変数の関係にあり、モデル$f(x;\beta)$はパラメータ$\beta$に関して線形です。モデルの未知パラメータ$\beta$を求めるため、残差$e$の二乗和を最小化することが目的です。web記事には基本的にこの線形回帰モデルのカーブフィット問題として説明されていることが多い印象です。
一方、上式後半は行列形式に等価変換したものです。観測数を$n$、未知パラメータの次元数を$p$とすると、$ X = ((x_{1,1},...,x_{1,p})^T,...,(x_{n,1},...,x_{n,p})^T)^T, Y=(y_1,...y_n)^T, \beta=(\beta_1,...,\beta_p)^T$となります。かなりややこしいが、各行列は以下のような形をとっています。
- $Y$：観測数$n \times 1$の行列。非説明変数の観測値を縦に並べたもの。
- $X$：観測数$n \times$未知数$p$の行列。説明変数の観測値$x_i=(x_{i,1},...,x_{i,p})^T$を並べたもの。
- $\beta$：$1 \times$未知数$p$の行列。推定すべき未知パラメータを縦に並べたもの。
ちなみに、この時の$\beta$に関する最小二乗推定量$\hat{\beta}$は以下のようになります。
\hat{\beta} = (X^TX)^{-1}X^TY
行列$X$と$Y$をnumpy.arrayなどで明示的に作成してから上記の計算をすることで$\hat{\beta}$を求めてもいいですが、計算時間やメモリ効率の観点からはあまりお勧めできません。以下の関数を使って解くようにしましょう。
#最小二乗法を解くための関数
scipy.optimizeでは以下が関数がサポートされています。
https://docs.scipy.org/doc/scipy/reference/optimize.html#equation-local-minimizers
| 関数 | 概要 |
|---|---|
| leastsq | Levenberg-Marquardt法を用いて（非線形）最小二乗解を求める |
| least_squares | 解の取りうる値に制約がある場合の（非線形）最小二乗解を求める |
| nnls | 行列形式で書き下した最小二乗問題について非負の条件付きの解を求める |
| lsq_linear | 行列形式で書き下した最小二乗問題について、解の取りうる値に制約がある場合の最小二乗解を求める |
leastsq
leastsqでは、モデル$f(x;\beta)$と解の初期値$\beta_0$を入力します。具体的には、defで定義した関数そのものを入力するため、関数の定義が必要です。
また、このメソッドでは非線形最小二乗問題によく用いられるLevenberg-Marquardt法を用います。そのため、必ずしもこれまで述べてきた線形最小二乗問題でなくてもよく、非線形最小二乗問題への適用も可能です。ただしそれゆえ、繰り返し計算のための初期値解$\beta_0$の設定が必要となります。
least_squares
下記のように解に制約がある場合の（非線形）最小二乗問題を解くためのメソッドです。
min \Sigma_i (y_i - f(x_i; \beta))^2 \text{ subject to }lb \leq\beta\leq ub
ここで、$lb$（lower bound）は解の取りうる値の最小値、$ub$（upper bound）は解の取りうる値の最大値を意味します。
least_squaresでは、leastsq同様、関数および初期値を入力する必要があります。lbとubの入力は任意です。
異なる点として、求解手法をLevenberg-MarquardtやTrust Region Reflective、dogboxから選択することが可能です。それぞれの長短があるらしいので、解こうとしている問題に応じて選択することが必要となります。
nnls
一方、nnlsでは行列形式で書き下された問題に対応しています。非負条件のもとで下記問題を解くことに相当します。
min |Y - X \beta|^2 \text{ subject to }０ \leq\beta
入力としては、$X$と$Y$に相当する行列を準備します。numpy.arrayなどで作成した行列を入力として食べさせる必要があります。行列が簡単にかける場合などは使いやすですが、巨大な係数行列を扱う際にはメモリ管理の観点から注意が必要です。
lsq_linear
下記の制約条件付き最小二乗問題を解くためのメソッドです。
min |Y - X \beta|^2 \text{ subject to }lb \leq\beta\leq ub
入力には$X$と$Y$に相当する行列を準備します。numpy.arrayのような実際の行列だけでなく、scipy.sparseで定義されるようなスパース行列を入力にすることも可能です。
lbとubの設定は任意なので、行列形式で書き下した制約なしの最小二乗問法を使いたいときはこのメソッドを用いると良いでしょう。Trust Region ReflectiveとBounded-Variable Least-Squaresから求解方法を選択することが可能です。
その他
他にも、下記のようなメソッドが用意されているようです。
- scipy.linalg.lstsq
- scipy.sparse.linalg.lsqr
- Find the least-squares solution to a large, sparse, linear system of equations
- https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.lsqr.html#scipy.sparse.linalg.lsqr
- scipy.sparse.linalg.lsmr
- Iterative solver for least-squares problems
- https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.lsmr.html#scipy.sparse.linalg.lsmr