原理
求解最佳投影方向，使得同類投影點盡可能的進，異類投影點盡可能的遠
同類投影點距離用同類樣本協方差矩陣表示
\[\omega^T \Sigma_i \omega \quad {第i類樣本協方差} \]
異類投影點距離
\[||\omega^T\mu_0 - \omega^T\mu_1||_2^2 \]
\(\mu_i \, {第i類樣本均值}\)
優化函數
\[J(\omega) = \frac{||\omega^T\mu_0 - \omega^T \mu_1||_2^2}{\omega^T (\Sigma_0 - \Sigma_1) \omega} \]
求上述函數極大值，解出\(\omega\)
定義類間、類內散度矩陣
1.類間散度矩陣
\[S_w = \Sigma_0 + \Sigma_1 = \sum_{x \in X0}(x - \omega)(x - \omega)^T + \sum_{x \in X1}(x - \omega)(x - \omega)^T \]
2.類內散度矩陣
\[S_b = (\mu_0 - \mu_1)(\mu_0 - \mu_1)^T \]
則
\[J(\omega) = \frac{\omega^T S_b \omega}{\omega^T S_w \omega} \]
利用拉格朗日乘數法，可得
\[\omega = S_w^{-1}(\mu_0 - \mu_1) \]
python程序
import numpy as np
import matplotlib.pyplot as plt
M = 2 #屬性個數
N = 50#二分類。每類樣本N個
#隨機生成兩個屬性的N個第一類樣本
feature11 = np.random.randint(0, 7, size = N)
feature12 = np.random.randint(0, 7, size= N)
temp_X1 = np.row_stack((feature11, feature12))
X1 = np.mat(temp_X1)
#隨機生成兩個屬性的N個第二類樣本
feature21 = np.random.randint(5,11, size= N)
feature22 = np.random.randint(7, 14, size= N)
temp_X2 = np.row_stack((feature21, feature22))
X2 = np.mat(temp_X2)
#求投影向量omega
mu1 = np.mat(np.zeros((2,1)))
mu2 = np.mat(np.zeros((2,1)))
X_1t = np.array(X1)
X_2t = np.array(X2)
for i in range(M):
mu1[i, 0] = sum([j for j in X_1t[i,:]])/N
for i in range(M):
mu2[i, 0] = sum([j for j in X_2t[i,:]])/N
#print(mu1, mu2)
s_w1 = np.mat(np.zeros(M))
s_w2 = np.mat(np.zeros(M))
for i in range(N):
s_w1 = s_w1 + (X1[:, i] - mu1)*(X1[:, i] - mu1).T
for i in range(N):
s_w2 = s_w2 + (X2[:, i] - mu2)*(X2[:, i] - mu2).T
s_w = s_w1 + s_w2
Omega = np.linalg.pinv(s_w)*(mu1 - mu2)
#print(Omega)
#畫出散點圖、投影面
fig = plt.figure(1)
plt.scatter(feature11, feature12, marker='+')
plt.scatter(feature21, feature22, marker='*')
xx_1 = np.linspace(0,10,num=50)
yy_1 = Omega[1,0]/Omega[0,0]*xx_1
plt.plot(xx_1,yy_1,color='r')
plt.show()
效果
參考資料
《機器學習》 周志華老師