Up

固有値と固有ベクトル

固有分解/スペクトル分解

 

実対称行列Aの固有分解(eigen decomposition, characteristic decomposition;スペクトル分解, spectral decomposition)は、numpy.linalg.eigh(A)で簡単に求めることができる。

例えば、実対称行列に対してスクリプト

w, v = numpy.linalg.eigh(A)

を実行すると、wに固有値、vに固有ベクトルが返される。このとき、配列wには、固有値が小さいものから順に並べられていて、行列vには固有値の並びに対応する順に固有ベクトルが列ベクトルとして並べられている。これを、固有値の絶対値の大きいものからの順に並べ直し、対応して並べ直した固有ベクトルは行ベクトルとして返すようにラッパーを掛けたものを用意した。

ここで、固有値あるいは固有ベクトルとは、以下のように定義されるものである(岡本安晴「統計学を学ぶための数学入門[下]培風館、第4.1節」)。

行列とベクトルおよび数値(スカラー)が次式

を満たすとき、数値を行列の固有値(eigenvaluecharacteristic value)と呼び、ベクトルを固有値に対する固有ベクトル(eigenvectorcharacteristic vector)と呼ぶ。

次の実対称行列は個の固有値をもつ。それらを、・・・、と書き、固有値に対する固有ベクトルを長さが1であるように調整したものをとおくと、次式(1)が成り立つ。

ここで、

である。

いま、

、 

とおくと、式(1)は

と書ける。

このとき、次式(2)が成り立つ。

式(2)から、固有値の大きいものが先にある方が、近似の精度の観点からは使い易いと思われる。また、固有ベクトルも、行ベクトルとして用意されている方が、取り出しやすい。このための関数eigen_symを以下のように用意した。ファイル名はeigen4me.pyである。続いて説明する使用例のスクリプトファイルとともに圧縮ファイルFilesEigen.zipにまとめた。

==============  eigen4me.py  ========================

import numpy as np

import math

 

def eigen_sym(A):

    """

            Re-arranged results from numpy.linalg.eigh

            The eigen values and vectors are re-ordered

            according to descending order of absolute values

            of eigen values

            Returned values:

                w : eigen values

                vrow : eigen vectors, each vector as a row vector

           

    """

    p = A.shape[0]

    w = None

    vrow = None

    if p == 1:

        w = np.full((1,),A[0][0])

        vrow = np.full((1,1), 1)

    else:

        w, v = np.linalg.eigh(A)

        s_w = np.empty(p)

        for i in range(p):

            s_w[i] = -math.fabs(w[i])

        idx = np.argsort(s_w)

        w = w[idx]

        vrow = v.T[idx]

 

        return w, vrow

==============================================

実対称行列Aを引数として関数eigen_sym

eigen_sym(A)

と呼び出せば、返り値はreturn

return w, vrow

によって、wに固有値が絶対値の大きいものから小さいものへの順に、vrowには対応する固有ベクトルが行ベクトルとして設定されたものが返される。wshape(p,)であり、vrowshape(p, p)である。関数eigen_symの使用例を次のスクリプトに示す。

=============  main.py  ==================

import numpy as np

from eigen4me import eigen_sym

 

u1 = np.array([[1, 1, 1, 1]])

w1 = -1.0

u2 = np.array([[1, 1, -1, -1]])

w2 = 0.5

u3 = np.array([[1, -1, 1, -1]])

w3 = 2.0

A = w1 * u1.transpose() @ u1 + w2 * u2.transpose() @ u2 +\

        w3 * u3.transpose() @ u3

print(A)

 

w, vrow = eigen_sym(A)

print('w.shape = ', w.shape)

print('vrow.shape = ', vrow.shape)

 

p = A.shape[0]

 

for i in range(p):

    print('w = ', w[i])

    print('vrow = ', vrow[i])

 

B = np.full((p, p), 0.0)

for i in range(p):

    trow = np.array([vrow[i]])      #   Change shape from (p,) to (1, p)

    B += w[i] * trow.transpose() @ trow

print('B =\n', B)

 

D = np.diag(w)

print('D =\n', D)

ReA = vrow.transpose() @ D @ vrow

print('ReA =\n', ReA)

=====================================

 

実行結果は以下のようになる。

[[ 1.5 -2.5  0.5 -3.5]

 [-2.5  1.5 -3.5  0.5]

 [ 0.5 -3.5  1.5 -2.5]

 [-3.5  0.5 -2.5  1.5]]

w.shape =  (4,)

vrow.shape =  (4, 4)

w =  8.0

vrow =  [-0.5  0.5 -0.5  0.5]

w =  -4.0

vrow =  [-0.5 -0.5 -0.5 -0.5]

w =  2.000000000000001

vrow =  [-0.5 -0.5  0.5  0.5]

w =  -1.2511631413972606e-15

vrow =  [ 0.5 -0.5 -0.5  0.5]

B =

 [[ 1.5 -2.5  0.5 -3.5]

 [-2.5  1.5 -3.5  0.5]

 [ 0.5 -3.5  1.5 -2.5]

 [-3.5  0.5 -2.5  1.5]]

D =

 [[ 8.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00]

 [ 0.00000000e+00 -4.00000000e+00  0.00000000e+00  0.00000000e+00]

 [ 0.00000000e+00  0.00000000e+00  2.00000000e+00  0.00000000e+00]

 [ 0.00000000e+00  0.00000000e+00  0.00000000e+00 -1.25116314e-15]]

ReA =

 [[ 1.5 -2.5  0.5 -3.5]

 [-2.5  1.5 -3.5  0.5]

 [ 0.5 -3.5  1.5 -2.5]

 [-3.5  0.5 -2.5  1.5]]

 

実対称行列の固有分解が確認できる。

 

固有分解の説明は、<岡本安晴「ないPythonデータ分析――多変量解析、ベイズ統計分析(PyStan、PyMC――丸善出版>でも行っている。

 

Up