1. 程式人生 > >講一下numpy的矩陣特征值分解

講一下numpy的矩陣特征值分解

6.2 大於 存在 code 假設 個性 4.3 向量 3.2

主要還是調包:

from numpy.linalg import eig

特征值分解: A = P*B*PT 當然也可以寫成 A = PT*B*P 其中B為對角元為A的特征值的對角矩陣。

首先A得正定,然後才能在實數域上分解,

>>> A = np.random.randint(-10,10,(4,4))
>>> A
array([[  6,   9, -10,  -1],
       [  5,   9,   5,  -5],
       [ -8,   7,  -4,   4],
       [ -1,  -9,   0,   6]])

>>> C = np.dot(A.T, A) >>> C array([[126, 52, -3, -69], [ 52, 292, -73, -80], [ -3, -73, 141, -31], [-69, -80, -31, 78]]) >>> vals, vecs = eig(C) >>> vals array([357.33597086, 174.10172008, 8.84429957, 96.71800949]) >>> vecs array([[-0.28738314, -0.51589436, -0.38221983, -0.71075449], [
-0.87487263, 0.12873861, -0.24968051, 0.39456798], [ 0.2572149 , -0.69304313, -0.33950158, 0.58161018], [ 0.29300052, 0.48679627, -0.82237845, -0.02955945]])

故使用時應先將特征值轉換為矩陣:

>>> Lambda = np.diag(vals)
>>> Lambda
array([[357.33597086,   0.        ,   0.        ,   0.        ],
       [  0.        , 
174.10172008, 0. , 0. ], [ 0. , 0. , 8.84429957, 0. ], [ 0. , 0. , 0. , 96.71800949]]) >>> np.dot(np.dot(vecs, Lambda), vecs.T) # 與C=A.T*A相等 array([[126., 52., -3., -69.], [ 52., 292., -73., -80.], [ -3., -73., 141., -31.], [-69., -80., -31., 78.]]) >>> np.dot(np.dot(vecs.T, Lambda), vecs) array([[171.65817919, 45.58778569, 53.20435074, 13.37512137], [ 45.58778569, 125.15670964, 28.22684299, 134.91290105], [ 53.20435074, 28.22684299, 129.48789571, 80.5284382 ], [ 13.37512137, 134.91290105, 80.5284382 , 210.69721545]])

故驗證了使用np中的eig分解為A=P*B*PT 而不是A=PT*B*P,其中P=vecs

C = vecs * np.diag(vals) * vecs.T # 這裏簡寫*為矩陣乘法

然後再來看使用np中的eig分解出來的vec中行向量是特征向量還是列向量是特征向量,只需驗證:A*vecs[0] = vals[0]*vecs[0]

>>> np.dot(C, vecs[0])
array([-12.84806258, -80.82266859,   6.66283128,  17.51094927])
>>> vals[0]*vecs[0]
array([-102.69233303, -184.34761071, -136.58089252, -253.97814676])

>>> np.dot(C, vecs[:,0])
array([-102.69233303, -312.62346098,   91.91213634,  104.69962583])
>>> vals[0]*vecs[:, 0]
array([-102.69233303, -312.62346098,   91.91213634,  104.69962583])

後者兩個是相等的,故使用np中的eig分解出的vecs的列向量是特征向量。

然後我們可以驗證P是單位正交矩陣:

>>> np.dot(vecs.T, vecs)
array([[ 1.00000000e+00, -7.13175042e-17, -2.45525952e-18,
         2.75965773e-16],
       [-7.13175042e-17,  1.00000000e+00,  2.49530948e-17,
        -5.58839097e-16],
       [-2.45525952e-18,  2.49530948e-17,  1.00000000e+00,
        -7.85564967e-17],
       [ 2.75965773e-16, -5.58839097e-16, -7.85564967e-17,
         1.00000000e+00]])

>>> np.dot(vecs, vecs.T)
array([[ 1.00000000e+00,  2.97888865e-16, -2.68317972e-16,
         1.69020590e-16],
       [ 2.97888865e-16,  1.00000000e+00, -4.40952204e-18,
        -6.24188690e-17],
       [-2.68317972e-16, -4.40952204e-18,  1.00000000e+00,
        -1.13726775e-17],
       [ 1.69020590e-16, -6.24188690e-17, -1.13726775e-17,
         1.00000000e+00]])

# 可以看到除對角元外其他都是非常小的數

即PT*P = P*PT = E , PT=P-1。事實上,在求解P的過程中就使用了施密特正交化過程。

另一方面,我們從數學角度來看:

首先補充一些數學知識:

首先AB相似:P-1*A*P=B,AB合同:CT*A*C=B,

二次型:系數在K中的一個n元二次多項式。由其生成的矩陣稱為二次型的矩陣,二次型的矩陣一定是對稱矩陣!

正定矩陣:實二次型xT*A*x > 0, x為列向量。

性質:假設A為正定矩陣

1、正定矩陣特征值全大於0

2、行列式 |A| >0

3、A合同於單位陣E,即存在可逆方陣C, s.t. CT*E*C = A = CT*C, 顯然可得A為對稱正定

正交矩陣:A*AT=AT*A=E ,

性質:

1、A的各行/列是單位向量且兩兩正交

2、AT=A-1

3、|A|=1

4、(Ax,Ay)=(x,y)x,y∈R

酉矩陣:A*AH=AH*A=E 顯然為正交矩陣在復數域上的推廣。其中H為共軛轉置。

性質:

1、A的各行/列是單位向量且兩兩正交

2、AH=A-1

3、|A|=1

(這裏補充一個厄米特矩陣:AH = A)

正規矩陣:A*AH=AH*A (以上的矩陣均有這個性質,故正規矩陣最為廣泛)

正規矩陣的充要條件是:存在酉矩陣U,使得A酉相似於對角矩陣B,即UH*A*U=U-1*A*U=B。

A = P*B*P ,其中B為對角元素為A的特征值的對角陣,P的列向量為特征值對應的特征向量(因為B每行乘以P每列)

講一下numpy的矩陣特征值分解