「スモールデータ解析と機械学習」を読む　～第3章回帰分析と最小二乗法 ③擬似逆行列～

はじめに

藤原幸一著「スモールデータ解析と機械学習」は，サンプル数が少ない状況でのモデル構築や予測精度の確保，統計的な不確実性の扱い方など，現場で直面する悩みに直結する理論と実践がコンパクトにまとまっている。製造業のデータ解析において「少ないデータだから仕方ない」と諦めるのではなく，「少ないからこそできる工夫」を身につけるためにこの本を読み，その学びをブログで共有しようと思う。

本記事は，「第3章回帰分析と最小二乗法」における，擬似逆行列に関する読書メモである。

3.9 擬似逆行列を用いる方法

3.7節「多重共線性の問題」や3.8節「サンプル数が入力変数の数よりも少ない場合」において，最小二乗法の問題点について説明していた。いずれの場合も原因は，最小二乗推定量 $\hat{\mathbfit{\beta}} = (\mathbfit{X}^T \mathbfit{X})^{-1} \mathbfit{X}^T \mathbfit{y}$ に含まれる $\mathbfit{X}^T \mathbfit{X}$ の逆行列計算に起因するものであった。

本節ではこの問題への対策として，擬似逆行列を用いる方法を説明している。

擬似逆行列(ムーア-ペンローズの擬似逆行列)

最小二乗法において，逆行列が問題になるのであれば，逆行列と似た働きをする都合の良い行列があれば便利である。

行列 $\mathbfit{A}$ に対して，

$\mathbfit{A} \mathbfit{A}^{+} \mathbfit{A} = \mathbfit{A}$
$\mathbfit{A}^{+} \mathbfit{A} \mathbfit{A}^{+} = \mathbfit{A}^{+}$
$(\mathbfit{A} \mathbfit{A}^{+})^T = \mathbfit{A} \mathbfit{A}^{+}$
$(\mathbfit{A}^{+} \mathbfit{A})^T = \mathbfit{A}^{+} \mathbfit{A}$

を満たす行列 $\mathbfit{A}^{+}$ がただ一つ存在して，これを擬似逆行列と呼ぶ。

最小ノルム解

最小二乗推定量 $\hat{\mathbfit{\beta}} = (\mathbfit{X}^T \mathbfit{X})^{-1} \mathbfit{X}^T \mathbfit{y}$ を求める際には，最小二乗推定量 $(\mathbfit{X}^T \mathbfit{X}) \mathbfit{\beta} = \mathbfit{X}^T \mathbfit{y}$ という連立1次方程式を $\mathbfit{\beta}$ について解けばよい。

この際，逆行列が登場するのであるが，逆行列の代わりに擬似逆行列を用いた解を最小ノルム解と呼ぶ。

行列 $\mathbfit{A}$ が正則であるとき，連立方程式

$\begin{align} \mathbfit{A} \mathbfit{x} = \mathbfit{b} \quad (\mathbfit{b} \neq \mathbfit{0}) \tag{1}\\ \\ \end{align}$

は唯一の解を持つ。

一方，行列 $\mathbfit{A}$ が正則でないとき，(1)の連立方程式は1つに定まらず，無数の解を持つ。
(1)式の解を $\tilde{\mathbfit{x}}$ とし， $\mathbfit{A} \mathbfit{x} = \mathbfit{0}$ を満たす $\mathbfit{x}$ ( $\mathrm{Ker}(\mathbfit{A})$ と表す)の要素の1つを $\mathbfit{x}_{null}$ とすると，

$\begin{align} \mathbfit{A} \tilde{\mathbfit{x}} + \mathbfit{A} \mathbfit{x}_{null} &= \mathbfit{A} (\tilde{\mathbfit{x}} + \mathbfit{x}_{null}) \\ \\ &= \mathbfit{b} \\ \\ \end{align}$

なので， $\tilde{\mathbfit{x}} + \mathbfit{x}_{null}$ も(1)式の解となる。

このように(1)式の解は無数に存在数が，解 $\mathbfit{x}$ のなかで，ノルム $\lVert \mathbfit{x} \rVert$ が最も小さくなる解を最小ノルム解と呼ぶ。擬似逆行列 $\mathbfit{A}^{+}$ によって求められる解は最小ノルム解となる。

具体的な擬似逆行列

行列 $\mathbfit{A} \in \mathbb{R}^{N \times M}$ について，行列 $\mathbfit{A}$ の特異値分解 $\mathbfit{A} = \mathbfit{U} \mathbfit{\Sigma} \mathbfit{V}^T$ を用いて，

$\begin{align} &\mathbfit{A}^{+} = \mathbfit{V} \mathbfit{\Sigma}^{+} \mathbfit{U}^T \\ \\ &\mathbfit{U} \in \mathbb{R}^{N \times N}, \quad \mathbfit{V} \in \mathbb{R}^{M \times M} \\ \\ &\mathbfit{\Sigma}^{+} = \begin{bmatrix} \sigma_1^{-1} & & & & & \\ & \ddots & & & & \\ & & \sigma_m^{-1} & & & \\ & & & 0 & & \\ & & & & \ddots & \\ & & & & & 0\\ \end{bmatrix} \in \mathbb{R}^{M \times M} \\ \\ \end{align}$

のように $\mathbfit{A}^{+}$ を定める。

この $\mathbfit{A}^{+}$ に対して，擬似逆行列の4つの条件のうち1つ目の条件 $\mathbfit{A} \mathbfit{A}^{+} \mathbfit{A} = \mathbfit{A}$ を考える。1つ目の条件の左辺を計算すると，

$\begin{align} \mathbfit{A} \mathbfit{A}^{+} \mathbfit{A} &= (\mathbfit{U} \mathbfit{\Sigma} \mathbfit{V}^T) (\mathbfit{V} \mathbfit{\Sigma}^{+} \mathbfit{U}^T) (\mathbfit{U} \mathbfit{\Sigma} \mathbfit{V}^T) \\ \\ &= \mathbfit{U} \mathbfit{\Sigma} \mathbfit{\Sigma}^{+} \mathbfit{\Sigma} \mathbfit{V}^T \\ \\ \end{align}$

となる。

ここで，

$\begin{align} \mathbfit{\Sigma}^{+} \mathbfit{\Sigma} &= \begin{bmatrix} \sigma_1^{-1} & & & & & \\ & \ddots & & & & \\ & & \sigma_m^{-1} & & & \\ & & & 0 & & \\ & & & & \ddots & \\ & & & & & 0\\ \end{bmatrix} \begin{bmatrix} \sigma_1 & & & & & \\ & \ddots & & & & \\ & & \sigma_m & & & \\ & & & \sigma_{m+1} & & \\ & & & & \ddots & \\ & & & & & \sigma_M\\ \end{bmatrix} \\ \\ &= \begin{bmatrix} 1 & & & & & \\ & \ddots & & & & \\ & & 1& & & \\ & & & 0 & & \\ & & & & \ddots & \\ & & & & & 0\\ \end{bmatrix} \end{align}$

となるので，

$\begin{align} \mathbfit{A} \mathbfit{A}^{+} \mathbfit{A} = \mathbfit{U} \begin{bmatrix} \sigma_1 & & & & & \\ & \ddots & & & & \\ & & \sigma_m & & & \\ & & & 0 & & \\ & & & & \ddots & \\ & & & & & 0\\ \end{bmatrix} \mathbfit{V}^T \\ \\ \end{align}$

となる。

一方で， $\mathbfit{A}$ の特異値分解を改めて書くと，

$\begin{align} \mathbfit{A} \mathbfit{A}^{+} \mathbfit{A} = \mathbfit{U} \begin{bmatrix} \sigma_1 & & & & & \\ & \ddots & & & & \\ & & \sigma_m & & & \\ & & & \sigma_{m+1} & & \\ & & & & \ddots & \\ & & & & & \sigma_M\\ \end{bmatrix} \mathbfit{V}^T \\ \\ \end{align}$

なので， $\mathbfit{A} \mathbfit{A}^{+} \mathbfit{A}$ は， $\mathbfit{A}$ の上位特異値を用いた近似になっていることが分かる。

※なお本書P71～72では，「 $\mathbfit{\Sigma}^{+}$ は0でない特異値について，その逆数 $1/\sigma_m$ を対角に並べた行列」とあった。

この定義にしたがうと， $\mathbfit{\Sigma}^{+}$ は $m \times m$ 行列であると読み取れるが，そうすると $\mathbfit{\Sigma}^{+} \mathbfit{\Sigma}$ が定義できなくなる( $\mathbfit{\Sigma}$ は $M \times M$ 行列であるため)。

そのため本記事のように，対角成分に0を補い， $\mathbfit{\Sigma}^{+}$ が $M \times M$ 行列になるように定義したが，その場合でも厳密には $\mathbfit{A} \mathbfit{A}^{+} \mathbfit{A} = \mathbfit{A}$ は成立たない。

まとめと感想

今回は，「第3章回帰分析と最小二乗法」における，擬似逆行列についてまとめた。

製造現場のデータは，測定項目間の相関が強いことが多いため多重共線性やサンプル不足によって $\mathbfit{X}^T \mathbfit{X}$ の逆行列計算が不安定になることが多い。
例えば，センサー数が多いのに稼働試験の回数が限られる初期立ち上げ段階や，似た工程条件で連続的に計測したデータでは，従来の最小二乗推定はノイズの影響を大きく受けやすい。

擬似逆行列の魅力は，「無理やり逆行列を求める」のではなく，「意味のある安定解」を返してくれる点にある。特にSVDベースの計算は数値安定性が高く，さらに特異値を閾値処理して小さい成分を切り捨てることで，現場的には「計測ノイズを含む要因を自然にフィルタリングする」ことになる。

実務的には，