NumPy で線形代数

このページでは、NumPy を用いて線形代数 (Linear Algebra) の計算を解く方法について解説します。

ベクトルのドット積 (点乗積)

ドット積 (a・b) は、np.dot(a, b) で計算できます。
ドット積は、1 次元配列ではベクトルの内積、2 次元配列では行列の積と同じ結果になります。

サンプルコード

import numpy as np

# 整数同士のドット積（通常の掛け算と同じになります）

# 3・4

np.dot(3, 4)

# 行列同士のドット積

# (2, 3)・(2, 3)

np.dot([2, 3], [2, 3])

# 複素数を扱うこともできます

# (2j, 3j)・(2j, 3j)

np.dot([2j, 3j], [2j, 3j])

# 2 次元行列同士のドット積を計算できます。

a = [[1, 0], [0, 1]]

b = [[4, 1], [2, 2]]

np.dot(a, b)

実行結果

>>> np.dot(3, 4)

>>> np.dot([2, 3], [2, 3])

>>> np.dot([2j, 3j], [2j, 3j])

(-13+0j)

>>> a = [[1, 0], [0, 1]]

>>> b = [[4, 1], [2, 2]]

>>> np.dot(a, b)

array([[4, 1],

[2, 2]])

ベクトルの内積

ベクトルの内積 (a・b) は np.inner(a, b) のようにして計算できます。

サンプルコード

# 通常のベクトル同士の内積

# (1, 2, 3)・(0, 1, 0)

a = np.array([1, 2, 3])

b = np.array([0, 1, 0])

np.inner(a, b)

# 多次元ベクトルの内積を計算

a = np.arange(24).reshape((2,3,4))

b = np.arange(4)

np.inner(a, b)

# b がスカラーの場合は以下のように指定します

a = np.eye(2)

np.inner(np.eye(2), 7)

実行結果

>>> a = np.array([1, 2, 3])

>>> b = np.array([0, 1, 0])

>>> np.inner(a, b)

>>> a = np.arange(24).reshape((2,3,4))

>>> a

array([[[ 0, 1, 2, 3],

[ 4, 5, 6, 7],

[ 8, 9, 10, 11]],

[[12, 13, 14, 15],

[16, 17, 18, 19],

[20, 21, 22, 23]]])

>>> b = np.arange(4)

>>> b

array([0, 1, 2, 3])

>>> np.inner(a, b)

array([[ 14, 38, 62],

[ 86, 110, 134]])

>>> a = np.eye(2)

>>> a

array([[ 1., 0.],

[ 0., 1.]])

>>> np.inner(np.eye(2), 7)

array([[ 7., 0.],

[ 0., 7.]])

ベクトルの外積

ベクトルの外積 (a×b) は、np.outer(a, b) のようにして計算できます。

サンプルコード

##マンデルブロット (Mandelbrot)

rl = np.outer(np.ones((5,)), np.linspace(-2, 2, 5))

im = np.outer(1j*np.linspace(2, -2, 5), np.ones((5,)))

grid = rl + im

grid

# 文字列をベクトルとして計算することもできます

x = np.array(['a', 'b', 'c'], dtype=object)

np.outer(x, [1, 2, 3])

実行結果

>>> rl = np.outer(np.ones((5,)), np.linspace(-2, 2, 5))

>>> rl

array([[-2., -1., 0., 1., 2.],

[-2., -1., 0., 1., 2.],

[-2., -1., 0., 1., 2.]])

>>> im = np.outer(1j*np.linspace(2, -2, 5), np.ones((5,)))

>>> im

array([[ 0.+2.j, 0.+2.j, 0.+2.j, 0.+2.j, 0.+2.j],

[ 0.+1.j, 0.+1.j, 0.+1.j, 0.+1.j, 0.+1.j],

[ 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j],

[ 0.-1.j, 0.-1.j, 0.-1.j, 0.-1.j, 0.-1.j],

[ 0.-2.j, 0.-2.j, 0.-2.j, 0.-2.j, 0.-2.j]])

>>> grid = rl + im

>>> grid

array([[-2.+2.j, -1.+2.j, 0.+2.j, 1.+2.j, 2.+2.j],

[-2.+1.j, -1.+1.j, 0.+1.j, 1.+1.j, 2.+1.j],

[-2.+0.j, -1.+0.j, 0.+0.j, 1.+0.j, 2.+0.j],

[-2.-1.j, -1.-1.j, 0.-1.j, 1.-1.j, 2.-1.j],

[-2.-2.j, -1.-2.j, 0.-2.j, 1.-2.j, 2.-2.j]])

>>> x = np.array(['a', 'b', 'c'], dtype=object)

>>> np.outer(x, [1, 2, 3])

array([['a', 'aa', 'aaa'],

['b', 'bb', 'bbb'],

['c', 'cc', 'ccc']], dtype=object)

単位行列 (identity array)

対角成分が1でそれ以外が全て0の行列 (単位行列) は、例えば、N次の単位行列の場合、np.identity(N)で作成できます。

サンプルコード

# 3 次の単位行列

np.identity(3)

# 4 次の単位行列

np.identity(4)

実行結果

>>> np.identity(3)

array([[ 1., 0., 0.],

[ 0., 1., 0.],

[ 0., 0., 1.]])

>>> np.identity(4)

array([[ 1., 0., 0., 0.],

[ 0., 1., 0., 0.],

[ 0., 0., 1., 0.],

[ 0., 0., 0., 1.]])

固有値、固有ベクトル

固有値と固有ベクトルの計算は、numpy.linalg パッケージ (本例では、LAとして読み込み) に含まれる、LA.eig() を用いて計算できます。
なお、固有ベクトルは正規化 (大きさ= 1) された状態で、出力されます。

サンプルコード

# numpy.linalg パッケージを読み込み

from numpy import linalg as LA

# 2 次行列 [[3, 1], [2, 4]] の固有値、固有ベクトル

matrix1 = np.array([[3, 1], [2, 4]])

w, v = LA.eig(matrix1)

w # 固有値

v # 固有ベクトル(正規化済)

# 3 次行列 [[1, -1, 2], [0, 2, -1], [0, 0, 3]] の固有値、固有ベクトル

matrix2 = np.array([[1, -1, 2], [0, 2, -1], [0, 0, 3]])

LA.eig(matrix2)

w # 固有値

v # 固有ベクトル(正規化済)

実行結果

>>> matrix1 = np.array([[3, 1], [2, 4]])

>>> matrix1

array([[3, 1],

[2, 4]])

>>> w, v = LA.eig(matrix1)

>>> w # 固有値

array([ 2., 5.])

>>> v # 固有ベクトル(正規化済)

array([[-0.70710678, -0.4472136 ],

[ 0.70710678, -0.89442719]])

固有値のみを求める場合は、LA.eigvals() を用いて計算できます。

サンプルコード

# 2 次行列 [[3, 1], [2, 4]] の固有値

LA.eigvals(matrix1)

# 3 次行列 [[1, -1, 2], [0, 2, -1], [0, 0, 3]] の固有値

LA.eigvals(matrix2)

実行結果

>>> LA.eigvals(matrix1)

array([ 2., 5.])

>>> LA.eigvals(matrix2)

array([ 1., 2., 3.])

行列式の値 (determinant)

行列式の値は、np.linalg.detを用いて計算できます。内部的には、LU 分解で計算しているようです。

サンプルコード

# 2次の正方行列 [[1, 2], [3, 4]] の行列式の値

a = np.array([[1, 2], [3, 4]])

np.linalg.det(a)

実行結果

>>> a = np.array([[1, 2], [3, 4]])

>>> np.linalg.det(a)

-2.0000000000000004

逆行列 (Inverse Matrix)

逆行列は、np.linalg.inv()で求めることができます。

サンプルコード

# 2 次の行列 A = [[1, 2], [3, 4]] の逆行列 A^-1

from numpy.linalg import inv

a = np.array([[1, 2], [3, 4]])

inv(a)

# 確認 (A・A^-1 = 単位行列になることを確認)

ainv = inv(a)

np.dot(a, ainv)

実行結果

>>> a = np.array([[1, 2], [3, 4]])

>>> inv(a)

array([[-2. , 1. ],

[ 1.5, -0.5]])

>>> ainv = inv(a)

>>> np.round(np.dot(a, ainv))

array([[ 1., 0.],

[ 0., 1.]])

行列式のランク (階数)

Rank (行列の階数) は、np.linalg.matrix_rank 関数にて算出できます。

サンプルコード

# numpy.linalg から matrix_rank 関数をロード

from numpy.linalg import matrix_rank

# 行列 [[1, 2, 3],[3, 2, 1],[2, 4, 6]] のランクを求めます。

matrix1 = np.array([[1, 2, 3], [3, 2, 1], [2, 4, 6]])

matrix1

matrix_rank(matrix1)

# 対角行列の場合は常に次数と同じになります

matrix2 = np.eye(4)

matrix2

matrix_rank(matrix2)

実行結果

>>> matrix1 = np.array([[1, 2, 3], [3, 2, 1], [2, 4, 6]])

>>> matrix1

array([[1, 2, 3],

[3, 2, 1],

[2, 4, 6]])

>>> matrix_rank(matrix1)

>>> matrix2 = np.eye(4)

>>> matrix2

array([[ 1., 0., 0., 0.],

[ 0., 1., 0., 0.],

[ 0., 0., 1., 0.],

[ 0., 0., 0., 1.]])

>>> matrix_rank(matrix2)

零行列

零行列 (全ての成分が 0) の行列は np.zeros() にて作成できます。

サンプルコード

# 長さが 3 の零行列

np.zeros(3)

# 長さが 4 x 4 の零行列

np.zeros([4, 4])

実行結果

>>> np.zeros(3)

array([ 0., 0., 0.])

>>> np.zeros([4, 4])

array([[ 0., 0., 0., 0.],

[ 0., 0., 0., 0.],

[ 0., 0., 0., 0.]])

ノルム

ベクトル a = [1,2,3,4] のノルム ||a||は、LA.norm() を利用して計算できます。

サンプルコード

from numpy import linalg as LA

a = np.array([1, 2, 3, 4])

LA.norm(a)

# 2 次の行列 b = [[1, 2], [3, 4]] のノルム ||b|| は以下のように計算できます。

b = np.array([[1, 2], [3, 4]])

LA.norm(b)

実行結果

>>> a = np.array([1, 2, 3, 4])

>>> LA.norm(a)

5.4772255750516612

>>> b = np.array([[1, 2], [3, 4]])

>>> LA.norm(b)

5.4772255750516612

参考・出典：Linear algebra (numpy.linalg) — NumPy v1.10 Manual

ベクトルのドット積 (点乗積)

ベクトルの内積

ベクトルの外積

単位行列 (identity array)

固有値、固有ベクトル

行列式の値 (determinant)

逆行列 (Inverse Matrix)

行列式のランク (階数)

零行列

ノルム

See also