※ この記事を書く際に利用したコードに関連するライブラリのバージョンは以下になります．
Python: 3.9.7
numpy: 1.20.3
matplotlib: 3.4.2
この記事で説明すること
一般逆行列の基本的なイメージと計算方法について説明します（主に参考文献[1]の4.2節に基づいています）．
行列を使った連立1次方程式ぐらいの線形代数の知識があれば理解できると思います．
具体的には以下について，pythonのコードを示しながら説明します．
- 一般逆行列とはなんなのか
- 一般逆行列がどのように計算されるかの確認
-
一般逆行列を使って解が無数にある連立1次方程式や，解がない連立1次方程式の それっぽい1つの解 を求める→ 長くなりそうなので，一般逆行列の初歩の初歩（その2: 連立1次方程式への応用）に回しました．
一般逆行列とはなんなのか
一般逆行列は，いろんな行列に対して逆行列っぽいものを考えられるように一般化したもの，というイメージです．
「疑似逆行列」や「一般化逆行列」とも呼ばれます．
つまり，逆行列が定義できない非正則な行列や正方でない行列に関しても，逆行列っぽいものを定めたものです．
※ 一般逆行列には様々なものがあるらしいですが，今回は ムーア・ペンローズの一般逆行列 のことを指します．
例えば，以下の4つの行列を考えた場合，逆行列が定義できるのは $A$ だけです．
- $B$ は行列式が0であり，正則でないので，逆行列を定義できません．
- $C,D$ は共に正方行列ではないので，逆行列を定義できません．
A =
\begin{bmatrix}
2 & 1 \\
1 & 3
\end{bmatrix}
,~
B =
\begin{bmatrix}
6 & 2 \\
3 & 1
\end{bmatrix}
,~
C =
\begin{bmatrix}
2 & 1 & 0 \\
-1 & 4 & 1
\end{bmatrix}
,~
D =
\begin{bmatrix}
2 & 1 \\
-1 & -2 \\
0 & 8
\end{bmatrix}
しかし，一般逆行列は $A,B,C,D$ 全てに対して定義できます．
また，$A$ の一般逆行列は $A$ の逆行列と一致します（一般逆行列は逆行列の一般化と考えられるため）．
一般逆行列の計算方法
特異値分解を利用して計算できる
一般逆行列の計算には 行列の 特異値分解 を用います．(特異値分解 - Wikipedia)
$A = USV^T$ と特異値分解できるとき，$A$ の一般逆行列 $A^-$ は$A^- = VS^{-1}U^T$ と分解できます．
つまり，$A$ を特異値分解して得られた $U, S, V$ を用いれば，$A^-$ を求めることができるということです．
$S$ は特異値が対角成分に並んだ対角行列となるので，$S^{-1}$ は非ゼロの特異値($S$ の非ゼロ要素)を逆数にした行列となります．
なぜこのように特異値分解を利用して求められるか？などの導出についてはこの記事では扱いません．
とりあえず，特異値分解を利用すれば簡単に一般逆行列が計算できることが伝われば幸いです．
ライブラリを用いて計算してみる
以下で，先に述べた非正則行列 $B$ の一般逆行列を求めます．
numpy および，線形代数関連の計算に特化した numpy.linalg を利用します．
import numpy as np
nl = np.linalg
では，以下の行列 $B$ の一般逆行列を計算してみます．
B = \begin{bmatrix}
6 & 2 \\
3 & 1
\end{bmatrix}
numpy.linalg.pinv を使って簡単に求められる
逆行列を求める関数として numpy.linalg.inv がありますが，一般逆行列を求める関数は numpy.linalg.pinv です．これを使います．
B = np.array([[6,2],[3,1]])
print('B^- = \n {0}\n'.format(nl.pinv(B)))
# 出力
# B^- =
# [[0.04 0.12]
# [0.02 0.06]]
結果から，$B$ の一般逆行列を，
B^- = \begin{bmatrix}
0.04 & 0.12 \\
0.02 & 0.06
\end{bmatrix}
と簡単に求めることができました．
特異値分解から計算した結果とあってるか確認
ここでは，numpy.linalg.pinv を使って求めた結果が，特異値分解から構成した一般逆行列と一致しているか確認してみます．
特異値分解は numpy.linalg.svd でできます．返り値は，$B=USV^T$ と分解した時の行列 $U,S,V$ です．$S$ は対角行列なので，対角成分が1次元のリストで返ってきます．
# SVDを行う
u, s, v = nl.svd(B)
print('Bの特異値 = {0}'.format(s))
# 出力
# Bの特異値 = [7.07106781e+00 4.77629022e-16]
$B$ の一般逆行列 $B^- = VS^{-1}U^T$ と分解できるので，上で求めた $U,S,V$ から $B^-$ を構成します．
$S^{-1}$ の計算の際には，特異値( $S$ の対角成分)の逆数を計算しなければいけません．
この時，0とみなす特異値の閾値を eps として定義しておいて，それより大きい特異値のみ逆数をとることにします．
eps = 1e-10 # 0とみなす特異値の閾値
# t <= eps なる特異値 t は0にし，そうでないものは逆数1/tをとる
s = [1/t * int(t > eps) for t in s]
# sを2次元ベクトルから2*2の対角行列にする
sigma = np.zeros((2,2))
for i in range(sigma.ndim):
sigma[i][i] = s[i] # 対角成分のみ設定
print('S^-1 = \n{0}'.format(sigma))
# 特異値分解の結果から一般逆行列を構成する
Bp = np.dot(np.dot(v.T, sigma), u.T)
print('B^- = \n {0}\n'.format(Bp))
# 出力
# S^-1 =
# [[0.14142136 0. ]
# [0. 0. ]]
#
# B^- =
# [[0.04 0.12]
# [0.02 0.06]]
上の結果から，numpy.linalg.pinv を用いて計算した $B$ の一般逆行列と一致することが確認できました．
ついでに，最初に述べた行列 $C, D$ に関しても一般逆行列を求めると以下のようになります．
一般に $m \times n$ 型の行列の一般逆行列は $n \times m$ 型になることがわかると思います．
C = np.array([[2, 1, 0], [-1, 4, 1]])
print('C^- = \n{0}'.format(nl.pinv(C)))
D = np.array([[2, 1], [-1, -2], [0, 8]])
print('D^- = \n{0}'.format(nl.pinv(D)))
# 出力
# C^- =
# [[ 0.44186047 -0.10465116]
# [ 0.11627907 0.20930233]
# [-0.02325581 0.05813953]]
#
# D^- =
# [[ 0.40729483 -0.18541033 -0.09726444]
# [-0.00911854 -0.01823708 0.12158055]]
一般逆行列の初歩の初歩（その2: 連立1次方程式への応用）の方では，この記事の続きとして，一般逆行列を使って解が無数にある場合や，解が存在しない場合でも，連立方程式のそれっぽい解を求めることについて説明します．
参考
[1] 金谷健一, これなら分かる最適化数学 ―基礎原理から計算手法まで― (共立出版), 2005