今回は2クラス分類を行うロジスティック回帰というモデルをPythonで実装しながら学びます。分類タスクを機械学習で解く場合の最も基礎的な部分になるので、ぜひ理解しましょう。
2クラス分類とロジスティック回帰
機械学習で行える基本的なタスクの一つが分類問題です。
特に2つのクラスに分類するものを2クラス分類（二値分類）といいます。
例えば
のような感じです。
（ちなみに、3クラス以上に分類する場合は多クラス分類といいます。）
この２クラス分類を解く最も基本的な手法がロジスティック回帰です。
※名前に「回帰」と付いていますが、回帰問題ではなく分類問題を解く手法なので、気を付けてください。
難しく聞こえますが、実は基本的な流れは線形回帰（機械学習の入り口「線形回帰」の実装を Python × NumPy で体験）と変わらないので、線形回帰の流れを抑えていれば大丈夫です。
使用するデータ
今回は、「工場における欠陥製品の数とその工場が正常に稼働しているかどうか」のデータを用います。
・欠陥製品が1～10個：正常稼働
・欠陥製品が11～20個：異常稼働
とする簡易的なデータです。
また、入力（欠陥製品数）に対応する正解データは、2クラス分類問題では0か1かのラベルを用います。
今回は正常稼働の場合をラベル0、異常稼働の場合をラベル1とします。
NumPy配列を使って以下のように作成します。
デフォルトのコード
import numpy as np
import matplotlib.pyplot as plt
# X：20個の入力データ（欠陥製品数 [個]）
X = np.array([1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20])
# Y：Xに対応する正解ラベル（0：正常、1：異常）
Y=np.array([0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1])
plt.yticks([0,0.5,1]) #y軸の目盛り
plt.scatter(X[:10] ,Y[:10], c='b')
plt.scatter(X[10:] ,Y[10:], c='r')
plt.xlabel('Number of defective pridcts')
plt.ylabel('Label')
plt.show()
ロジスティック回帰
今回は上記のデータを使ってロジスティック回帰を説明していきますが、線形回帰と異なるところは「仮説」と「目的関数」だけです。
仮説と目的関数の中身が違いますが、全体的な流れは回帰問題の時と同じで、次の通りです。
1. 仮説
さて、モデルの出力値となる仮説ですが、回帰問題の時のようにデータの分布に沿うような関数が欲しいわけではありません。
今回のような二値分類では、各データに対して正解ラベルが0か1かで与えられているのでした。
そのため、データを与えた時の仮説の値によって、そのデータが0か1かということを決められればいいのです。
そのため、ラベルがちょうど切り替わるところで値が変化し、それ以外のところでは0と1で一定値をとり続ける次のような仮説であれば分類が行えそうです。
実は、「シグモイド関数」という関数がちょうどこのような性質を持っています。
シグモイド関数は、入力値をzとして
で表される関数です。
ロジスティック回帰では、このzに線形回帰での仮説「aX + b」を入れて新たな仮説とします。
仮説をこのように設定することの利点がもう一つあります。それは、出力値を確率として解釈できることです。
シグモイド関数の出力は必ず0～1の間の値なので、例えば「出力値が0.7であれば、クラス1に分類される確率が70%」というような解釈ができます。
つまり、実際に分類する際に「出力が何%以上ならクラス1に分類する」という閾値を決めることができるのです。
以下でシグモイド関数を用いた仮説を実装してみます。
の部分は、Numpyの機能を使って「np.exp(-z)」と書けば計算できます。
デフォルトのコード
import numpy as np
import matplotlib.pyplot as plt
X = np.array([1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20])
Y = np.array([0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1])
##########################################################################
#初期値を0としてパラメータを設定
a=0
b=0
#仮説を定義
z=a*X+b #線形回帰の仮説
h=1/(1+np.exp(-z)) #シグモイド関数
#データと仮説を重ねてプロットします。
plt.yticks([0,1])
plt.scatter(X, Y) #データをプロット
plt.plot(X, h, c='orange') #仮説をプロット
plt.legend((u'Data',u'Regression line')) #凡例
plt.xlabel('Number of defective pridcts')
plt.ylabel('Label')
plt.show()
シグモイド関数を用いて仮説を定義しましたが、パラメータa, bが初期値なので当然先ほどのような形にはなっていません。
ここから、学習によって適切なa, bを見つけていきます。
2. 目的関数
目的関数とは簡単に言うと仮説による予測値と正解値の誤差でしたが、分類問題では回帰問題で用いた二乗和誤差とは違います。
今回の正解値は0か1かの二値ラベルだったので、異常な工場}=1\right))
については})
は1に近い値になり、正常な工場}=0\right))
についてはは0に近い値になってほしいです。
ここから、それぞれの場合について目的関数を次のように定義します。
それぞれグラフにしてみると、次のようになっています。
左のグラフ（}=0)
の時）から、確かに})})
を最小化するとが1に近づくことがわかります。また、反対にが0に近い時は目的関数の値が無限大に大きくなります。
同様に右のグラフ（}=1)
の時）から、})})
を最小化するとが0に近づくことがわかります。また、反対にが1に近い時は目的関数の値が無限大に大きくなります。
よって適切な目的関数を定義できました。
最後に、この二つの場合を次のように一つの式にまとめることができます。
}\log{(h^{(i)})}+(1-Y^{(i)})\log{(1-h^{(i)})}\right))
これは全データに対する目的関数の値を足し合わせたものですが、一つ一つのデータに対しては、ラベル})
の値を代入するときちんと先ほど考えた通りになっていることがわかります。
これが二値分類問題の目的関数であり、「交差エントロピー関数」と呼ばれる関数です。
Pythonでは以下のように書けます。
対数の計算は、Numpyの機能である「np.log」を使用できます。
デフォルトのコード
import numpy as np
X = np.array([1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20])
Y = np.array([0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1])
#仮説のパラメータ初期値
a=0
b=0
#仮説を定義
z=a*X+b
h=1/(1+np.exp(-z))
############################################################
#データ数
m=X.shape[0]
#目的関数（交差エントロピー関数）
J = -(1/m)*((Y*np.log(h)+(1-Y)*np.log(1-h)).sum())
print("初期値の（a,b）での目的関数の値：%f"% J)
初期値での目的関数の値は0.69くらいになったと思います。
それでは、この値を小さくするようにパラメータを更新しましょう。
3. パラメータ更新
ロジスティック回帰でも、線形回帰と同様の考えで勾配降下法によってパラメータを更新します。
そのためパラメータa, bに関する目的関数Jの勾配（微分）を求めなければなりません。
ただし、ロジスティック回帰におけるパラメータの更新式は実は線形回帰とまったく同じになります。
目的関数が違うのに不思議ですが、先ほど定義した目的関数を実際にa, bで微分すれば確かめられます（シグモイド関数σの微分が(1 – σ)σとなることを用いれば簡単です）。
デフォルトのコード
import numpy as np
import matplotlib.pyplot as plt
X = np.array([1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20])
Y = np.array([0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1])
#データ数
m=X.shape[0]
#仮説のパラメータ初期値
a=0
b=0
####################################################################################
#繰り返し回数
iterations = 100000
#学習率
alpha=0.01
#目的関数の値を保存しておくためのリスト
cost=[]
#学習（指定した繰り返し回数だけパラメータを更新）
for iter in range(iterations):
#現在のパラメータで仮説を計算
z=a*X+b
h=1/(1+np.exp(-z))
#パラメータを更新
a = a - (alpha/m)*(((h-Y)*X).sum())
b = b - (alpha/m)*((h-Y).sum())
#更新後のパラメータで仮説と目的関数を計算
z=a*X+b
h=1/(1+np.exp(-z))
J = -((Y*np.log(h)+(1-Y)*np.log(1-h)).sum())/m
#目的関数の値の推移を後で見れるように、現在の目的関数の値をリストに保存します
cost.append(J)
# 学習した結果を表示
print("学習後のa: %f,"% a)
print("学習後のb: %f,"% b)
print("学習後の目的関数の値: %f,"% J)
###################################################################
#目的関数の値の推移と仮説のグラフを並べて表示する
fig = plt.figure()
#1つ目のグラフ
ax1 = fig.add_subplot(1, 2, 1)
ax1.plot(cost) #目的関数の値を保存したリストをプロット
ax1.set_xlabel('Iteration')
ax1.set_ylabel('Cost')
#2つ目のグラフ
ax2 = fig.add_subplot(1, 2, 2) #データと仮説を重ねてプロット
ax2.set_yticks([0,1])
ax2.scatter(X, Y) #データをプロットz=a*X+b
X=np.arange(0,20,0.1)
z=a*X+b #最終的な仮説
h=1/(1+np.exp(-z))
ax2.plot(X, h, c='orange') #仮説をプロット
ax2.legend((u'Data',u'Regression line')) #凡例
plt.xlabel('Number of defective pridcts')
plt.ylabel('Label')
fig.tight_layout()
fig.show()
更新するごとに目的関数の値も小さくなり、仮説もデータによくフィットしています。
未知データの予測
上で得られたパラメータで仮説を定義し、新しくデータxを入れてそれがクラス0とクラス1のどちらに属するかを予測してみます。
デフォルトのコード
import numpy as np
#未知データ
x = 25
#学習済みの仮説
a=1.302372
b=-13.550009
z=a*x+b
h=1/(1+np.exp(-z))
#仮説の値が0.5以上か未満かでクラス分類
if h>=0.5:
print("欠陥製品数が",x,"個の工場は異常稼働している")
else:
print("欠陥製品数が",x,"個の工場は正常稼働している")
まとめ
今回は、「欠陥製品の個数から工場が正常稼働しているかどうかを分類」するタスクを例に、ロジスティック回帰を学びました。線形回帰と流れは同じでしたが、仮説と目的関数（+正解データの与え方）が異なりました。これが機械学習で分類問題を解く方法の基礎になるので、しっかりと要点を抑えておきましょう。
