SciPyによる線形代数(scipy.linalg)②|SciPy入門 #2
SciPyについて色々と話題になったので公式チュートリアルを元にまとめています。
Linear Algebra (scipy.linalg) — SciPy v1.2.1 Reference Guide
#1では上記の内容をまとめますが、scipy.linalgについてBasic routinesまでの内容を取り扱いました。
#2ではscipy.linalgの後半ということでDecompositionsとMatrix Functionsについてまとめます。
以下目次になります。
1. Decompositions
2. Matrix Functions
3. まとめ
1. Decompositions
・固有値、固有ベクトル(Eigenvalues and eigenvectors)
"The eigenvalue-eigenvector problem is one of the most commonly employed linear algebra operations."とあるように固有値・固有ベクトルの問題は線形代数において最も良く出てくる処理と言っても過言ではないと思われます。scipy.linalgでは固有値・固有ベクトルを求めるにあたってはlinalg.eigを用います。
import numpy as np #すでに行っていれば不要
from scipy import linalg #すでに行っていれば不要
A = np.array([[1, 2], [3, 4]])
la, v = linalg.eig(A)
l1, l2 = laprint(l1, l2) # eigenvalues
print(v[:, 0]) # first eigenvector
print(v[:, 1]) # second eigenvectorprint(np.sum(abs(v**2), axis=0)) # eigenvectors are unitary
v1 = np.array(v[:, 0]).T
print(linalg.norm(A.dot(v1) - l1*v1)) # check the computation
上記のようにlinalg.eigメソッドの引数として行列Aを与えると固有値と固有ベクトルを上記のようにして得ることができます。
・特異値分解(Singular value decomposition)
特異値分解 - Wikipedia
によると特異値分解は「線形代数学における、行列に対する行列分解の一手法である。」とされています。行列の分解に使われる手法で、主成分分析などとも関わりがあります。特異値分解を行うにあたってはlinalg.svdとlinalg.diagsvdを用います。
import numpy as np
from scipy import linalg
A = np.array([[1,2,3],[4,5,6]])
print(A)M,N = A.shape
U,s,Vh = linalg.svd(A)
print(s)
Sig = linalg.diagsvd(s,M,N)print(U)
print(Sig)
print(Vh)
print(U.dot(Sig.dot(Vh)))
分解にあたっての具体的な処理をlinalg.svd、得られた値を元に行列化を行う処理をlinalg.diagsvdで行うことができます。
2. Matrix Functions
SciPyには行列の値を操作する様々なメソッドが実装されています。チュートリアルでは具体的にlinalg.expm、linalg.sinm、linalg.cosm、linalg.tanmなどが実装されています。
from scipy.linalg import expm, sinm, cosm
print(expm(np.zeros*1
print(cosm(a) + 1j*sinm(a))
上記ではオイラーの等式について実際に計算した例です。上記の他にもlinalg.sinhm、linalg.coshmなどが基本的な関数だと実装されています。
また、上図のような理工系でよく用いられる特殊な行列を生成するためのメソッドも実装されています。たとえばlinalg.pascal(4)とすればパスカルの三角形の情報を持った行列を生成することができます。
このようにSciPyでは様々なメソッドが実装されています。
3. まとめ
#1、#2で触れたようにscipy.linalgでは多くの線形代数にまつわるメソッドが実装されています。以前MATLABなどを使用した際に似たような機能を見た記憶があるのですが、最近ではフリーで使用することのできるPythonベースで話が進むことが多いようです。
#3以降ではその他の機能について見ていければと思います。
*1:2,2))))
a = np.array([[1.0, 2.0], [-1.0, 3.0]])
print(expm(1j*a