単純パーセプトロンをPythonで組んでみる

いきなり自分でハードル上げてみました（笑）。ちなみに何故単純パーセプトロンを最初に持ってきたのか？というと、id:echizen_tmさんのブログ記事でも触れておられる通り

機械学習には大きく分けて「識別関数」「識別モデル」「生成モデル」の3つの種類がある。このなかで識別関数は確率を使わないので初心者が入門するのに最適。
識別関数で有名なのはSVM(Support Vector Machine、サポートベクターマシン)。名前を聞いたことがある人も多いと思う。そこで早速SVMを学ぼうとすると敷居が高くて挫折しがち。

実はSVMは(大雑把に言うと)パーセプトロンという基礎的な識別関数に「マージン最大化」と「カーネル関数」という考え方を導入したもの。なので機械学習入門者は最初にパーセプトロンを学ぶのが良いと思われる。

それゆえ、実際に僕も以前Matlabで糞コード書きながら勉強してた時はやはり単純パーセプトロンから入ったからです。ということで、僕のブログでも単純パーセプトロンから入ろうと思います。

単純パーセプトロンとは何ぞや

以前社内の勉強会で話した時のスライドをSlideShareに上げておきます（ソース部分はMatlabコードだったため割愛）*1。

Simple perceptron by TJO from Takashi J OZAKI

さぁこれで後はスライドに書いた定義をそのままPythonでベタっと書くだけ！・・・ということで、スライド一通りお読みください。今回はたまたまスライドがあったので、手抜きになってしまってすみません*2。

実際にPythonで糞コードを書いてみる

注）ここから先は、最低でも上記slideshareを参考にして一度はPythonで組んでみた後で、読むようにして下さい！さもないと、機械学習のアルゴリズムを知った上で実装するというスキルは身に付かないです！

・・・というわけで僕が書いてみたものは以下の通り。恥を忍んで公開します。糞コードなのは僕のコーディングスキルが糞だからです（泣）。ソースはGitHubにあります。

# -*- coding: utf-8 -*-

# 演算用にNumPyを、プロット用にmatplotlibをimport
import numpy as np
import matplotlib.pyplot as plt

# 識別関数の本体：y=w'xを計算してるだけ
def predict(wvec,xvec):
    out=np.dot(wvec,xvec)
    if out>=0:
        res=1
    else:
        res=-1
    return [res,out]

# 学習部：識別関数に学習データを順繰りに入れて、
# 重みベクトルを更新する
def train(wvec,xvec,label):
    [res,out]=predict(wvec,xvec)
    c=0.5 # 学習係数。大き過ぎてもまともに学習してくれないので1未満ぐらいで
    if out*label<0:
        wtmp=wvec+c*label*xvec
        return wtmp
    else:
        return wvec

# 以下本体
if __name__=='__main__':
    
    item_num=100 # 学習データは100個
    loop=1000    # 収束判定は一切してないけどとりあえず1000回ループ
    init_wvec=[1,-1,1] # 重みベクトルの初期値、適当
    
    # 学習データはxy平面の第1象限と第3象限に50個ずつ
    # np.random.randomでばらつかせてみた
    x1_1=np.ones(int(item_num/2))+10*np.random.random(int(item_num/2))
    x1_2=np.ones(int(item_num/2))+10*np.random.random(int(item_num/2))
    x2_1=-np.ones(int(item_num/2))-10*np.random.random(int(item_num/2))
    x2_2=-np.ones(int(item_num/2))-10*np.random.random(int(item_num/2))
    z=np.ones(int(item_num/2)) # こいつはバイアス項、今回は無視して1に固定

    # 学習データを1つのマトリクスにまとめる
    x1=np.c_[x1_1,x1_2,z]
    x2=np.c_[x2_1,x2_2,z]
    class_x=np.array(np.r_[x1,x2])

    # 教師ラベルを1 or -1で振って1つのベクトルにまとめる
    label1=np.ones(int(item_num/2))
    label2=-1*np.ones(int(item_num/2))
    label_x=np.array(np.r_[label1,label2])

    # NumPy.ndarrayはappend()メソッドが使えないので
    # 糞コードらしく初期arrayを作ってそこに垂直に足していく策を取る
    wvec=np.vstack((init_wvec,init_wvec))
    
    # ループ回数の分だけ繰り返しつつ、重みベクトルを学習させる
    for j in range(loop):
        for i in range(item_num):
            wvec_new=train(wvec[-1],class_x[i,:],label_x[i])
            wvec=np.vstack((wvec,wvec_new))
    w=wvec[-1] # 重みベクトルを垂直に足していった最後のものを採用する
    print w

    # 分離直線を引く
    x_fig=range(-15,16)
    y_fig=[-(w[1]/w[0])*xi-(w[2]/w[1]) for xi in x_fig]

    # 漫然とプロットする
    plt.scatter(x1[:,0],x1[:,1],marker='o',color='g',s=100)
    plt.scatter(x2[:,0],x2[:,1],marker='s',color='b',s=100)
    plt.plot(x_fig,y_fig)
    plt.show()

ということで、実行すると大体こんな感じになります。

そうそう、収束判定一切かけてないので、デタラメな分離直線になることもあります（真横とか）。収束判定かけたい人は、例えば2回目以降のループで重みベクトルwvecの更新量に一定の打ち切り基準を設けるとかすれば良いかなーと思います*3。

ポイントは、__main__以外には2つしか関数が要らないところ。しかもその関数は超カンタン。ぶっちゃけ、このコード書く時に苦労したのはPythonの制御構造とNumPy arrayの使い方のところぐらいです*4。

機械学習は「書いてみれば意外とカンタン」

というわけで、Pythonで書いてみましたがものすごーく簡単でした*5。

パーセプトロンに限らず、機械学習の多くは数学的アルゴリズムの定義こそ厳めしくて難解に見えますが、実は式自体の形がコードで書くのに適してる*6ことが多いので、コードにしてしまえばカンタンということが往々にしてあります。

なので、機械学習のアルゴリズムに興味を持ったら、まずは自分で手を動かしてコードを組んでしまうことをぜひお薦めします。概要もわかるし、挙動もわかるし、パラメータを色々試せばそのアルゴリズムの特徴もわかるので*7、やって損はないですよ！

*1:余談だけどslideshareをどうやってembedするか分からなかったので調べてしまった→ http://staff.hatenablog.com/entry/2012/02/27/174925

*2:次回以降はちゃんと定義式の簡単な説明とかブログ上でやろうと思ってます。。。いや本気ですよ！？

*3:糞コードしか書けない人間の言葉なので説得力ありませんが

*4:でも糞コードで汚ねえぞコラ、とか突っ込まないでー

*5:Python自体が難しいと思ったあなたは残念ながら僕と同レベルです笑

*6:ΣとかΠとか内積とか

*7:SVMの汎化性能なんかは絶対にコード書いてパラメータ変えながら分離超平面を眺めてみないとわからない