分類②・線形判別分析
ノートブックの作成
Colabにアクセスし，新規にノートブックを作成してください． ノートブックのタイトルは chapter8 とします． また，numpy，matplotlib.pyplot，scipy.optimize を導入しておいてください．
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import minimize
データの準備
前回と同じリンゴとメロンの 大きさ（$x_1$） と 重さ（$x_2$） を用います（未知のデータ$u$をリンゴとしてデータに追加）． 各要素はベクトル${\bf x} = (x_1,x_2)$として表すことに注意してください． ここでは，あらかじめデータを標準化しておきます．
[In:]
x1 = np.array([9.5, 11.2, 10.3, 11.5, 13.4, 12.9, 14.2]) #大きさ（cm）
x2 = np.array([278, 298, 382, 443, 1221, 1305, 1512]) #重さ（g）
label = np.array(["red","red","red","red","green","green","green"]) #カテゴリ
x1 = (x1 - x1.mean()) / x1.std()
x2 = (x2 - x2.mean()) / x2.std()
print(x1)
print(x2)
[Out:]
[-1.49090975 -0.41564757 -0.98490402 -0.22589542 0.9758682 0.65961462
1.48187394]
[-0.99454837 -0.95468671 -0.78726775 -0.66568969 0.88492881 1.05234778
1.46491594]
plt.scatter(x1, x2, c=label)
plt.gca().set_aspect("equal", "datalim")
plt.ylim([-1.8, 1.8])
plt.xlabel("x1")
plt.ylabel("x2")
線形判別分析
前回は 分類 における 決定境界 と ベクトル の役割について解説しました． 今回は，決定境界を決めるベクトル${\bf w}=(w_1, w_2)$の導出方法について考えます． 分類のためのアルゴリズムは様々ありますが，今回は 線形判別分析（Linear Discriminant Analysis: LDA） に挑戦してみましょう． 線形判別分析は，データを分類するために，線形（直線） の決定境界を求める方法です．
代表ベクトル
まず最初に，各ラベル（ここでは，リンゴとメロン）の 代表ベクトル を求めます． 代表ベクトルとは重心（平均）のことです． ここで，リンゴの代表ベクトルを${\bf a} = (a_1, a_2)$， メロンの代表ベクトルを${\bf m} = (m_1, m_2)$と表すことにします． 下記のように計算すると，${\bf a} = (-0.78, -0.85)$，${\bf m} = (1.04, 1.13)$となります．
$$ (a_1, a_2) = \left( \frac{-1.49 + -0.42 + -0.98 + -0.23}{4}, \frac{-0.99 + -0.95 + -0.79 + -0.67}{4} \right) $$
$$ (m_1, m_2) = \left( \frac{0.98 + 0.66 + 1.48}{3}, \frac{0.88 + 1.05 + 1.46}{3} \right) $$
[In:]
apple = np.where(label=="red") #リンゴのインデックスを取得
melon = np.where(label=="green") #メロンのインデックスを取得
a = np.array([x1[apple].mean(),x2[apple].mean()])
m = np.array([x1[melon].mean(),x2[melon].mean()])
print(a)
print(m)
[Out:]
[-0.77933919 -0.85054813]
[1.03911892 1.13406417]
plt.scatter(x1, x2, c=label)
plt.scatter(a[0], a[1])
plt.scatter(m[0], m[1])
plt.gca().set_aspect("equal", "datalim")
plt.ylim([-1.8, 1.8])
plt.xlabel("x1")
plt.ylabel("x2")
決定境界
次に，決定境界を決めるベクトル${\bf w}=(w_1, w_2)$について考えます． 前回は${\bf w}=(1, 1)$と勝手に決めてしまいましたが， 線形判別分析では，下記の条件を満たす$w_1$と$w_2$の組み合わせに限定することにします．
$$ |{\bf w}| = w_1^2 + w_2^2 = 1 $$
これは，半径$1$の円周上にベクトル${\bf w}$があることを意味しています． つまりベクトルの 大きさ は固定して，ベクトルの 回転 に注目するということです． 例えば，${\bf w} = \left( \cos(\pi / 6), \sin(\pi / 6) \right) = \left( \frac{\sqrt{3}}{2}, \frac{1}{2} \right)$ は，この条件を満たしています． このとき，決定境界は $x_2 = - \sqrt{3} \cdot x_1$ となります（前回の資料を参照すること）．
[In:]
w1 = np.sqrt(3) / 2
w2 = 1 / 2
w = np.array([w1, w2])
np.linalg.norm(w, ord=2) #ユークリッド距離(ord=2)
[Out:]
0.9999999999999999
目的関数
最適な${\bf w}=(w_1, w_2)$を求めるための目的関数を定めます． ここでも内積の概念を利用します． ベクトル${\bf w}$と，リンゴとメロンの代表値である${\bf r}$と${\bf m}$の内積を求めます． $|{\bf w}|=1$であることから，内積はベクトル${\bf w}$に対する，ベクトル${\bf r}$と${\bf m}$の正射影になります． 正射影された2点間の距離を目的関数$g({\bf w})$とし，これが最大となる${\bf w}=(w_1, w_2)$を求めます． 例えば，${\bf w} =\left( \frac{\sqrt{3}}{2}, \frac{1}{2} \right)$ のときは，$g({\bf w}) = 2.57$となります．
$$ g(w) = |{\bf w} \cdot {\bf a} - {\bf w} \cdot {\bf m}| $$
[In:]
# 正射影された2点間の距離
def gw(w, a, m):
d = np.abs(w.dot(a) - w.dot(m))
return d
gw(w, a, m)
[Out:]
2.5671370685019306
数値解の導出
最後に数値解を求めていきましょう． 回帰で用いた minimize 関数を利用しますが，少し工夫が必要です．
まず，今回の目的は，$g({\bf w})$の 最大化 が目的です． しかし，minimize関数は 最小化 することしか出来ません． そこで，$g({\bf w})$の符号を反転した目的関数を新たに定義します． 負の値を最小化することは，元の正の値を最大化することに等しくなります．
次に，$|{\bf w}| = 1$を表す条件式を定義します． しかし，minimize関数では，eqという比較演算子で条件を定義する必要があり， この比較演算子は「== 0」を意味します． このため，条件式は$|{\bf w}| - 1 == 0$と変形する必要があります．
この制約条件の基で，目的関数を最大化すると， ${\bf w} = (0.68, 0.74)$が解として求まります． 決定境界は$x_2 = - 0.92 \cdot x_1$となり， リンゴとメロンを適切に分類出来ていることがわかります．
[In:]
#目的関数（最小化するため符号を反対に）
def gw2(w, a, m):
d = -1 * gw(w, a, m)
return d
#条件式（ベクトル長から1を引く）
def norm(w):
l = np.linalg.norm(w, ord=2) #ベクトルの大きさ
return l - 1
#制約条件（eqは0と等しいことを表す）
cons = (
{"type": "eq", "fun": norm}
)
w_init = w
result = minimize(gw2, w_init, args=(a, m), constraints=cons, method="SLSQP")
print(result.x)
[Out:]
[0.67569856 0.73717806]
最後に決定境界を求めプロットしてみましょう．
x1_ = np.arange(-3, 3, 0.01)
x2_ = -1 * (result.x[0] / result.x[1]) * x1_
plt.plot(x1_, x2_) #決定境界
plt.scatter(x1, x2, c=label)
plt.scatter(a[0], a[1])
plt.scatter(m[0], m[1])
plt.gca().set_aspect("equal", "datalim")
plt.ylim([-1.8, 1.8])
plt.xlabel("x1")
plt.ylabel("x2")
課題
下記のデータを分類する決定境界を線形判別分析で求めてください．
x1 = np.array([10, 13, 15, 19, 23, 25])
x2 = np.array([17, 18, 18, 16, 17, 16])
label = np.array(["blue","blue","blue","pink","pink","pink"]) #カテゴリ
課題を完成させたら，chapter8.ipynb を保存し， 共有用のリンク と ノートブック（.ipynb） をダウンロードして提出してください． このとき，必ず事前に下記の設定を行ってから提出してください．
- ノートブックの設定で「セルの出力を除外する」のチェックを外す
- ノートブックの変更内容を保存して固定
- 共有設定で「学校法人椙山女学園大学」を「閲覧者」に設定