DCGANでテキスト生成を試してみた2

こんにちは、Link-Uの町屋敷です。

前回は、Word2Vec(以下W2V)とDCGANで英語のテキストを生成しているサイトがあったので、日本語でも試しよう!と言う事で実装してみたところなかなか残念な結果だったので、今回はもう少しなんとかならないか続きをします。

生成結果を視覚化する

結果を良くするにはパラメーター調整をしていかなければならないんですが、画像と違って文章生成の場合生成された結果が直感でわかりにくくなっています。

そこで、わかりやすくするためにジェネレーターが生成した結果をWord2Vecで文章に戻すのではなく、それを次元圧縮してプロットします。

    def show_tsne(self, epoch, real, fake, save_pic):
        from sklearn.manifold import TSNE
        col = np.zeros(len(real))
        col = np.append(col, np.ones(len(fake)))
        
        comp = TSNE(n_components=2, random_state=0).fit_transform(np.reshape(np.vstack((real,fake)), (-1, self.img_cols * self.img_rows * self.channels)))
        plt.scatter(comp[:, 0], comp[:, 1], c=col)
        plt.xlim(-30,30)
        plt.ylim(-30,30)
        if save_pic:
            plt.savefig("images/sent_gan_%d.png" % epoch)
        else:
            plt.show()
        plt.close()

本物の文章と偽物の文章を受け取ってTSNEで圧縮してプロット。今回は本物を紫、偽物を黄色でプロットする

いい偽物の文章が生成されていれば、2色のプロットは混ざり合うが、そうでなければ分離される。

ためしに前回の最終結果で試してみるとこんな感じ。

Epoch 100

Epoch 1000

完全に分離されているのがわかります。

MeCabからの品詞の情報を使う

上記の確認方法を使ってDやGのモデル、W2Vのベクトルの長さ、学習率などのパラメーターをいろいろいじくり回しましたが良い結果にはなりませんでした。

出力された結果を見てみると、品詞がガバガバになっていたので、MeCabから品詞の情報を抜き出して、W2Vの後ろに追加したら多少マシになるんじゃないたと思って追加しました。

    def AnalyzeWord(self, word):
        gen = self.mecab.parse(word, as_nodes=True)
        for w in gen:
            return w.posid, w.surface

品詞はposidから取得できます。posidさえあれば単語の活用は単語の原型から変形させればいいので、簡単のために入力する単語をすべて原型にしました。

結果はこんな感じ。

Epoch 100

Epoch 1000

多少マシなはなったがまだひどい。

ちなみに生成された文章はこんな感じ。

掘り起こしぎゅれんずつ欄オアフまた。。

掘り起こしぎゅれんずつプロバスケットボールリーグオアフまた。。

掘り起こしぎゅれんずつシティーハンタースペシャルオアフまた。。

最後には。がつくっていうのくらいしか学習できてない感がある。

品詞の並びを生成してそれを単語に変換する

品詞の並びを生成する

流石にこのままでは終われないので作戦変更(妥協)。

直接雑音から文章を生成するのを諦め、DCGANと2つ使った文章の生成を目指す。

1つめのDCGANで雑音から品詞の列を生成するジェネレーターを生成し、

品詞を単語に変換するニューラルネットを2つめのDCGANで生成することを目指す。

まず雑音から品詞の列を生成する。これが出来ないなら話にならない。

今までW2Vを入力していたところの先程の品詞の列のみを入力してDCGAN。

Epoch 100

Epoch 300

Epoch 500

Epoch 1000

Epoch 5000

どうやら品詞の列ぐらいは生成できるようだ。

学習した重みを保存する。転移学習に使えるよう各レイヤーに名前わつけて、(name = ‘PC1’など)combinedではなけgenerateorの重みを保存する。

品詞の列を生成するニューラルネットをpos_makerと名付ける

    def build_pos_maker(self):
        model = Sequential()
    
        model.add(Dense(256 * int(self.img_rows * self.img_cols), activation="relu",
                          input_dim=self.latent_dim, name = 'PD1'))
        model.add(Reshape((self.img_cols , self.img_rows, 256)))
        model.add(UpSampling2D())
        model.add(Conv2D(128, kernel_size=3, strides=2, padding="same", name = 'PC1'))
        model.add(BatchNormalization(momentum=0.8, name = 'PB1'))
        model.add(Activation("relu"))
        model.add(UpSampling2D())
        model.add(Conv2D(128, kernel_size=3, strides=2, padding="same", name = 'PC2'))
        model.add(BatchNormalization(momentum=0.8, name = 'PB2'))
        model.add(Activation("relu"))
        model.add(UpSampling2D())
        model.add(Conv2D(64, kernel_size=3, strides=2, padding="same", name = 'PC3'))
        model.add(BatchNormalization(momentum=0.8, name = 'PB3'))
        model.add(Activation("relu"))
        model.add(UpSampling2D())
        model.add(Conv2D(32, kernel_size=3, strides=2, padding="same", name = 'PC4'))
        model.add(BatchNormalization(momentum=0.8, name = 'PB4'))
        model.add(Activation("relu"))
        model.add(Conv2D(1, kernel_size=3, strides=1, padding="same", name = 'PC5'))
        model.add(Activation("tanh"))
    
        model.summary()
        
        noise = Input(shape=(self.latent_dim,))
        print(noise)
        pos = model(noise)
    
        return Model(noise, Container(noise, pos)(noise), name='pos_gen')

//略

        self.generator.save_weights('pos_gan_weight.h5')

2つの学習機をつなげる

        # The generator takes noise as input and generates imgs
        z = Input(shape=(self.latent_dim,))
        pos = self.pos_maker(z)
        img = Input(shape=self.txt_shape)
        
        fake = self.generator(pos)
        # For the combined model we will only train the generator
        #Fix pos_maker weights
        self.pos_maker.load_weights('pos_gan_weight.h5', True)
        self.pos_maker.trainable = False

z(雑音)をpos_makerで品詞の列に変換して、generatorで単語にする。pos_makerの重みは更新しない。

2つめの学習機を初期化する

入力に品詞出力に単語のW2V列をセットして学習するこれを初期値とする。

    def train_pos_to_word(self):
        pos_to_w2v = self.generator
        pos_to_w2v.compile(loss='mse',
        optimizer='adam')
        X_train = joblib.load('{0}/pos_labels.pkl'.format(WRITE_JOBLIB_DIR))
        X_train = np.reshape(X_train, (-1,self.img_cols, self.img_rows, 1))
        X_test = joblib.load('{0}/vectorized_words.pkl'.format(WRITE_JOBLIB_DIR))
        X_test = np.reshape(X_test, (-1, self.img_cols, self.img_rows , self.channels))
        pos_to_w2v.fit(X_train, X_test, batch_size=32,
            epochs=12,
            verbose=1)

結果

かかった時間は手元のパソコンとGPUサーバーで2倍ほどの差があった。

Epoch5

Epoch40

Epoch400

長く学習するとモードがなくなって悪化する。このほうほうではすく偏るらしく元のサイトでも課題になっていた。

Epoch40のときの生成分がこちら。

収録オニヅカ原典する通りは単体。

がけっぷち相当し初めてリリース。。。

詳細は記述実際は彼奴。。

すべてリリース必ずた。。。。

巻が第発売全リリース。。

は本来の年リリース。。。

はエピローグ相当する両性それぞれ取消。

が」全巻。。。。

巻から日発売。。。。

スケリグカン・シヌ。。。。。。

第初めてリリースは取消。。。

ミッションイー。。。。。。。

元データが会話文ではないので名詞が多い。

またこの手のやつは大体そうだが、別に意味を理解しているわけではないので、意味不明な文も結構生成されている。

まとめ

DCGANで文章生成は品詞の列は作れることはわかった。

ただ文章生成になると微妙なので、別の手法(別のGANとかRNN系)も試してみて比較したほうが良さそう。

GPUを使って無線LANをクラックする話 / プロファイラの言うことがみんな違うよ編

ついにソメイヨシノも咲いて、「春」って感じですな。

梅も桜も花の形が似ていますが、それもそのはず。も、どちらもバラ科サクラ属なのだそうです。じゃ、他には?という事で調べてみたのですが、モモもそうなんだって。そういえばたしかに、さくらんぼと桃の実ってなんとなく似てますね。

というわけで(?)、今回はモモの花です。

前回までのまとめ

前回までこの「GPUを使って無線LANをクラックする話」では、Pythonの標準プロファイラcProfileを使って測定して議論していましたが、cProfileは実は普通に使う限りではメインスレッドしか測定してくれない事がわかりました(別記事)。マルチスレッドで動いているPythonプログラムのためのプロファイラはいろいろあるけれど(別記事)、C言語のスタックまでプロファイルしてくれる、vmprofというのがどうやら我らがPyritには良さそうです。というわけで、このvmprofを使って再度取り直します。

vmprofを使ってもう一回測定するぜ!→やっぱなんか変だぜ!

というわけでモクモクと測定していきます。

まず可視化のためのサーバを建てる必要があります(デフォルトだとvmprof.comを使うようになっているのですが、このサーバは落ちたまま放置されているようです)。

git clone https://github.com/vmprof/vmprof-server.git
cd vmprof-server

pip install -r requirements/development.txt
python manage.py migrate
python manage.py runserver -v 3

vmprof自体はpipで入るのでもっと簡単です。ただし、さっき自分が建てたWebサーバを指定するのを忘れずに:

pip install vmprof
python -m vmprof --web --web-url "http://127.0.0.1:8000/" pyrit benchmark

で、CPUを1コアだけ使って解くようにした結果がこちら:

全処理時間のうちrun関数が95%使っていて、さらにCPUDevice_solveというC言語の関数が2.5%使っている、などなど。

95%を使っているというrun関数はパスワードを処理するスレッドのエントリーポイントです。CPU/CUDA/OpenCLで共通で、キューから候補のパスワードを取り出し(_gather)、self.solve()を使って計算して別のキューへ投げる(_scatter)という簡単なものです:

    def run(self):
        while not self.shallStop:
            essid, pwlist = self.queue._gather(self.buffersize, timeout=0.5)
            if essid is not None:
                if not self.isTested:
                    self._testComputeFunction(101)
                    self.isTested = True
                t = time.time()
                res = self.solve(essid, pwlist)  # ←ここのメソッド呼び出しでCPU/CUDA/OpenCLが分岐する
                assert len(res) == len(pwlist)
                self.compTime += time.time() - t
                self.resCount += len(res)
                self.callCount += 1
                if self.compTime > 0:
                    # carefully move towards three seconds of execution-time
                    avg = (2 * self.buffersize +
                           (self.resCount / self.compTime * 3)) / 3
                    self.buffersize = int(max(self.minBufferSize,
                                              min(self.maxBufferSize, avg)))
                self.queue._scatter(essid, pwlist, res)

この関数と、そこから先の処理が95%を占めている、というのは理解できます。この「95%」というのは、SIGPROFの基準で95%なので、全スレッドで実際に実行された時間のうち、95%の時間という意味になるようです。言い換えると、sleepしたりしている時間は入らないということです。

しかし、CPUDevice_solveが2.5%とはどういう事か。CPUDevice_solveは、CPUのコアで実行したときに実際にOpenSSLを使ってSHA1を何度も計算する処理です。一番重い、CPUを酷使する処理のはずで、流石にたった2.5%という事はありえないはずです。

…バグかなぁ?

cProfile、ふたたび

cProfileはpython -m cProfile pyrit benchmark と実行するとメインスレッドしかプロファイルしてくれませんが、ソースコードを書き換えれば別のスレッドのプロファイルも取ることができます。これで同じ処理を測定してみます。

さきほどのrun関数の周囲を次のようなコードで囲みます:

    def run(self):
        import cProfile
        prof = cProfile.Profile()
        prof.enable()
        try:
            # ... 元々の処理 ...
        finally:
            prof.disable()
            prof.dump_stats('{}.dat'.format(threading.currentThread().ident))

threading.currentThread().ident はスレッドごとに異なる値で、これをファイル名に使っています。

次のようなコマンドで測定結果を確認してみると(pstatsのドキュメントはこちら):

$ python
Python 2.7.15rc1 (default, Nov 12 2018, 14:31:15) 
[GCC 7.3.0] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import pstats
>>> s = pstats.Stats('140487566550784.dat') # ←ファイル名は毎回異なる
>>> s.sort_stats('cumtime').print_stats()
Fri Mar 29 11:59:24 2019    140487566550784.dat

         3194 function calls in 69.366 seconds

   Ordered by: cumulative time

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
       28   68.735    2.455   68.735    2.455 {method 'solve' of '_cpyrit_cpu.CPUDevice' objects}
       28    0.004    0.000    0.627    0.022 /home/server01/src/Pyrit/cpyrit/cpyrit.py:637(_gather)
       10    0.001    0.000    0.622    0.062 /usr/lib/python2.7/threading.py:309(wait)
       54    0.620    0.011    0.620    0.011 {time.sleep}
        1    0.000    0.000    0.116    0.116 /home/server01/src/Pyrit/cpyrit/cpyrit.py:99(_testComputeFunction)
       27    0.003    0.000    0.004    0.000 /home/server01/src/Pyrit/cpyrit/cpyrit.py:692(_scatter)
      216    0.001    0.000    0.001    0.000 {method 'extend' of 'list' objects}
       55    0.000    0.000    0.000    0.000 /usr/lib/python2.7/threading.py:285(__enter__)
       54    0.000    0.000    0.000    0.000 /usr/lib/python2.7/threading.py:400(notifyAll)
       54    0.000    0.000    0.000    0.000 {sorted}
      139    0.000    0.000    0.000    0.000 {method 'acquire' of 'thread.lock' objects}
... (以下省略) ...

run関数の内側でprof.enable()としているので、上の結果の一覧にrun関数は出てきません。run関数以下の内訳を測定している、ということになります。処理に掛かった全時間が69.366秒、結果の1行目によるとCPUDeviceのsolveメソッドの実行に掛かった時間が68.735秒(約99%)、などなどといったことが分かります。

CPUDeviceのsolveメソッドの実態は、vmprofでCPUDevice_solveとなっていたC言語の関数です。vmprofでは処理の2.5%を占めているというこの関数が、cProfileによるとほぼ99%だそうで…えー、言ってることが、真逆ですねぇ。

まとめ

SIGPROF(アラーム・シグナル)を使った統計的プロファイラであるvmprofを使ってPyritのbenchmarkを測定してみたものの、標準の決定論的プロファイラであるcProfileとは大きく結果が異なる事がわかりました。状況証拠から言って、vmprofのバグじゃないかなぁと思います(しかし残念ながらもう開発も止まってしまった模様…)。前回のような簡単なコードでは結果が一致したからといって信用せず、常に複数のプロファイラを併用して結果を突き合わせたほうが良いのかもしれません。

…マジ、みんなどうやってプロファイル取ってるんだ?と思って調べたらTensorFlowは自分たち専用のプロファイラを作ってた。まっ、そうですよねー。

残念ながら大人の事情により、今回で一端このシリーズは最終回らしいです。機材を貸してくれた人ありがとう。…ちょっとぐらい速くしたかったなぁ。

マルチスレッドで走るPythonプログラムをプロファイルするには?

さて、前の記事では実験とソースコード読解の両方から、python -m cProfile test.pyとして測定した場合、メインスレッドで実行される部分しか計測されていない事がわかりました。

でも、「Deep LearningといえばPython」な昨今、マルチスレッドで動くPythonのプログラムは目の前にいっぱい存在するわけで、それらを測りたいわけで。どげんかせんといかん。

ある意味で本題はここからです。方法をいくつか調べました。

その1:起動したスレッドごとにcProfileを設定する

cProfile.Profileクラスを使うと、python -m cProfile test.pyのようにしてコマンドライン引数を使って行っていたプロファイリングを、コード中からもできます。この方法ならメインスレッド以外にもプロファイラを設定して計測することができます。

こんな感じ:

import threading
import cProfile
import math


class Worker(threading.Thread):
    def __init__(self, n: int):
        threading.Thread.__init__(self)
        self.n = n
        self.sum = 0

    def run(self):
        prof = cProfile.Profile() # プロファイラの作成
        prof.enable() # プロファイラの有効化
        try:
            for i in range(0, self.n):
                self.sum += math.sin(i * i)
        finally:
            prof.disable() # プロファイラを止める
            prof.print_stats('cumtime') # 結果の表示


def main():
    th = Worker(100000000)
    th.start()
    th.join()


if __name__ == '__main__':
    main()

前回と違い、今度こそはmath.sinが出てきます:

% python3 test.py
         100000001 function calls in 18.187 seconds

   Ordered by: cumulative time

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
100000000   18.187    0.000   18.187    0.000 {built-in method math.sin}
        1    0.000    0.000    0.000    0.000 {method 'disable' of '_lsprof.Profiler' objects}

しかしながら、この方法はソースコードを書き換えなければなりません。多くの場合において不可能ではないでしょうが、若干居心地は悪いですし、run関数がC言語で書かれていたりすると、正直かなり大変です。

とはいえ、この方法では確実に枯れているであろう、そしてPython開発チームが公式にメンテナンスしている、cProfileを使うことができます。cProfileに慣れているとか、なんかgithubに上がってる新しいよくわかんないプロダクトは嫌だなぁ、といった人にはこれが良いかもしれません。

その2:vmprof

前回読んだように、cProfileは全ての関数呼び出し(とreturn)を追跡して記録しているのですが、こういうプロファイラを「決定論的プロファイラ」と呼ぶそうです。わざわざそんな名前がつくからには別の種類のプロファイラもあるわけで、vmprof統計的プロファイラです。

具体的に何をどう測るのかというと、(Linuxも含めたUNIXでは)setitimerというOSの提供するインターバル・タイマーを使って、定期的にスレッドを全て止め、libunwindを使ってスタック上のコール・チェーンの状態(今、どの関数がどの関数を呼び出している最中なのか)を記録し、その情報を使って「測定」しているようです。

なので、cProfileと違ってすべての関数呼び出しが把握できるわけではなく、例えば短時間だけ実行されるあまり呼ばれない関数は取りこぼされることがあります。まぁ、そういう関数はプロファイルを取って突き止めたい「ホットな」関数ではないですから、測定できなくてもよい、という判断なのでしょう。

一方で、cProfileはC言語の関数に関しては前回見た通りPythonのインタプリタから呼び出す瞬間に測定していて、そこから先のC言語間の関数呼び出しは関知しない(できない)わけですが、この方法ではC言語の部分まで一緒にプロファイリングできちゃいます。CUDAやC言語の複雑に絡み合うシステムには、こっちの方が有用かもしれません。

結果のビジュアライズも中々わかりやすい。C言語の関数も、シンボル情報がある部分に関しては関数名も表示されます:

run関数でmath.sinの和を計算している処理のうち、integerをdoubleにするPyFloat_AsDoubleがどうやら11%程度使っているらしい、といった情報までわかります。

難点があるとすれば、測定するスレッドが決められない点でしょうか。ソースコードを追いかけた限りでは、必ず全てのスレッドを追いかけて結果を合算しているように読めました。同じクラスから派生するさまざまなクラスが一緒に走る場合(今回高速化しているPyritがまさにそのケースに当てはまる)、ちょっとデータの解釈に困るかもしれません。

C言語のスタックまで計測出来る点は素晴らしいのですが、普通に実行しただけでは必ずしも関数名が見えるとは限らない点はちょっと残念です。python3-dbgパッケージを入れてデバッグ用pythonでも実行してみたのですが、内部の全ての関数シンボルまでは残っていないようでした。VMの中身を深く掘らないといけない時は、pythonを自分でビルドする必要があるでしょう。

このWebインターフェイスなのですが、デフォルトではvmprof.comというウェブサービスにプロファイルした結果をアップロードして表示する(!)というクラウド絶頂期のイマドキらしい中々大胆な作りになっていまして、会社で書いてるプロプライエタリなソースで利用するのにはちょっと憚られるかもしれません。まぁもっとも、ここ最近はずっと落ちているのでそもそもアップロード出来ないんですが(それはそれでどうなんだ)。

もちろん、サーバのソースコードも公開はされていて比較的簡単には使えるので、自分で立てればそのへんは解決可能(というかそうしないと使えない)です。

その3:Pyflame

上のvmprofと同じ統計的プロファイラですが、こちらはLinuxでデバッガを実装するときに使われるシステムコールであるptraceを使っているらしいです(よってLinux専用)。ptraceを使っているおかげで、すでに動いているプロセスに対してアタッチ出来るのが強みだそうな。

Python3.7には未対応な上でIssueには進捗がなく、Python2.7で試してみた限り、Cのスタックも追いかけてくれませんでした(ptraceの力をもってすれば原理的には可能だと思うんですが…)。2019年3月現在は、あんまりおすすめできない選択肢だと思います。

その4:pprofile

こちらは決定論的プロファイラと統計的プロファイラの両方の実装が含まれています。ドキュメントに書いてあるとおり、決定論的・統計的プロファイラをのどちらでもマルチスレッドに対応していて、しかもPythonだけで書かれています。どうなっているのか調べてみたところ、threading.settraceというsys.settraceのマルチスレッド対応版みたいなAPIがあって、それをそのまま利用しているようです。

全てがPythonで書かれている都合上オーバーヘッドが少々大きいのと、C言語のプロファイルは取れません。ただ、著者の主張している通り、動かすのはすごく楽です。測定対象のソースコードを書き換える必要性も一切ありません。

その他:マルチスレッドは(そのままでは)ダメだったやつ

DISるつもりは無いんですが、マルチスレッドはダメそうだったやつも記録しておきます。シングルスレッドのPythonプログラムをデバッグしているときは便利なんじゃないでしょうか。

cprofilev

cProfileの便利なラッパー、だそうな。しかしソースを読んだ限り、便利になった代償として、マルチスレッドでこれを使うのはさらに難しくなっているようです。

line_profiler

cProfileと違って1行ごとに測定できる、というのが売りの決定論的プロファイラ。仕組みもcProfileとだいたい同じで、違いはPyEval_SetProfileの代わりにPyEval_SetTrace使っていること。この関数を使うとPythonの関数の呼び出しの記録だけでなく、Pythonの命令が1行実行されるごとにcallbackが呼ばれるそうな(逆に、C言語の関数のCALL/RETURNは記録されない)。

ここ以外の事情も概ねcProfileと同じ感じでして、コードを修正しなくてもプロファイルを取ってくれるというkernprofコマンドのmain関数はメインスレッドしか測らないものの、手動で測定コードを挿入すればマルチスレッドでも測れるようになっているように読めます(段々めんどくさくなってきて実験はしてません)。

まとめ

プロファイラもまさに十人十色といった感じでした。プロファイラを使う前に、時間を掛けてソースコード読んで、何をどう測っているのか、調べる価値はあると思います。

どれが良いのかは、正直なところ、測定対象のプログラムによると思います。

  • とりあえず動かしたいだけならpprofが一番楽です。
  • pyflameはおすすめできないと書きましたが、今すでに動いているPythonのプログラムをプロファイルできるのはpyflameだけです。
  • Cのコールスタックまでプロファイルを取れるのは、vmprofだけです。
  • cProfileは公式でメンテナンスされていて、枯れているでしょう。

今回ターゲットにしている、CとCUDAとPythonのアマルガムであるPyritの場合は、vmprofが適切かなぁと思います。CaffeやTensorFlowの人たちは何使ってプロファイル取ってるんだろう…。

Pythonの標準プロファイラ、cProfileの中身を追いかける

Python、使ってますか。

わたしは最近、諸事情により、最近PythonとC言語とCUDAを使ったソフトウェア高速化に取り組んでいます(できているとは言っていない)。ちなみに、ディープラーニングじゃないよ。

高速化するときに最も重要な道具は、ベンチマークとプロファイラであることを疑う人は居ないと思います。ベンチマークで高速化の効果を測定し、プロファイラでボトルネックを発見する。これが基本的かつ王道のサイクルですよね。

ところで、プロファイラはどのような仕組みで動いているのでしょうか。追いかけたこと、あります?

何を測っているのか知らずにプロファイラを使う者がおるか!

プロファイラは何を測っていて、何を測っていないのか。

Pythonに最初からついてくるプロファイラ「cProfile」を使ってプロファイリングを行っていたのですが(CUDAとC以前にPythonがボトルネックっぽいので)、わたしはてっきりこのプロファイラはGILを握っている全てのスレッドのそれぞれの動作を測ってるんだと思ってました。だって、「マルチスレッドだけどGILのせいで同時に動くスレッドは1つしかない」というPythonのプログラムの動作をプロファイリングする方法としては、それが一番妥当かつ自然に感じられたからです (でも、今考え直してみると、そうでもない気がしてきた…わからん)。ドキュメントにはマルチスレッドに関しては何にも書いてないし。

しかし、Profilerの出力する結果を見ているとどうにもおかしい気がしてきまして。おおむね予想通りの結果も出てくるけど、どこかおかしい気もする。わたしの予想の方が間違っているのか、プロファイラの測定の誤差や限界なのか。うーん。

なんでそんな考えが浮かんでくるのかといえば、(今回のも含めて)高速化したいプログラムというのは大抵複雑だからです。

というわけで、超単純なサンプルで測ってみました:

import threading
import cProfile
import math

class Worker(threading.Thread):
    def __init__(self, n: int):
        threading.Thread.__init__(self)
        self.n = n
        self.sum = 0

    def run(self):
        for i in range(0, self.n):
            self.sum += math.sin(i * i)

def main():
    th = Worker(100000)
    th.start()
    th.join()

if __name__ == '__main__':
    main()

見たまんまですが、math.sin(i * i)の和を計算するという何の意味もない計算を、メインスレッドとは別のスレッドで行う、無意味だけど時間の掛かるプログラムです。

これをcProfileして出た結果は:

% python -m cProfile -s 'cumtime' test.py
         45 function calls in 21.670 seconds

   Ordered by: cumulative time

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
        1    0.000    0.000   21.670   21.670 {built-in method builtins.exec}
        1    0.000    0.000   21.670   21.670 <string>:1(<module>)
        1    0.000    0.000   21.670   21.670 test.py:17(main)
        5   21.670    4.334   21.670    4.334 {method 'acquire' of '_thread.lock' objects}
        1    0.000    0.000   21.665   21.665 threading.py:1024(join)
        1    0.000    0.000   21.665   21.665 threading.py:1062(_wait_for_tstate_lock)
        1    0.000    0.000    0.005    0.005 threading.py:828(start)
        1    0.000    0.000    0.005    0.005 threading.py:533(wait)
        1    0.000    0.000    0.005    0.005 threading.py:263(wait)
        1    0.000    0.000    0.000    0.000 {built-in method _thread.start_new_thread}
...(略)...

無意味な計算のために20秒ほど掛かっていていますが、その計算に関わっているWorker.run関数やmath.sinは一切リストアップされておらず、スレッドのjoinのみが入っています。

…えー。つまり、cProfileはマルチスレッド対応してないんだねぇ。

…今まで計測した結果、全部無駄だねぇ。…まぁ、そういう事もあるよ。

ソースコードを追いかける

せっかくですから、この機会にソースコードを追いかけて、どのようにプロファイルをしているのか、その仕組みを掴んでおくことにします。事前に予告しておきますが、本当に眺めて見学するだけです。さらに深くソースコードを追いかけたい時のとっかかりにでも使ってください。とはいえ、プロファイラの中身について言及している記事はインターネッツ上にあんまりないので、今回このように記事を書いております。

マルチスレッドのプログラムを手っ取り早くプロファイルしたい人は、読み飛ばして次の「マルチスレッドで走るPythonプログラムをプロファイルするには?」へ進むのも良いかもしれません。

さて。コンソールからpython -m cProfile test.pyを実行すると、cProfile.Profilerの次の関数が最終的には呼ばれます:

# https://github.com/python/cpython/blob/fcd5e84a515e19409840c570730f0728e9fcfc83/Lib/cProfile.py#L97

    def runctx(self, cmd, globals, locals):
        self.enable()
        try:
            exec(cmd, globals, locals) # コマンドライン引数で指定されたファイルの中の文字列を実行する
        finally:
            self.disable()
        return self

exec関数はいわゆるeval関数の親戚で、与えられた文字列(cmd)を、所定の環境の下で実行する組み込み関数です。となると、その周囲を囲うenable/disableが重要そうです。

self.enableの実態を追いかけていくとC言語の関数に行き着きまして:

// https://github.com/python/cpython/blob/62be74290aca26d16f3f55ece7ff6dad14e60e8d/Modules/_lsprof.c#L683

static PyObject*
profiler_enable(ProfilerObject *self, PyObject *args, PyObject *kwds)
{
    int subcalls = -1;
    int builtins = -1;
    static char *kwlist[] = {"subcalls", "builtins", 0};
    if (!PyArg_ParseTupleAndKeywords(args, kwds, "|ii:enable",
                                     kwlist, &subcalls, &builtins))
        return NULL;
    if (setSubcalls(self, subcalls) < 0 || setBuiltins(self, builtins) < 0)
        return NULL;
    PyEval_SetProfile(profiler_callback, (PyObject*)self);
    self->flags |= POF_ENABLED;
    Py_RETURN_NONE;
}

いちばん重要なのは、PyEval_SetProfileの呼び出しだと思います。これはpythonインタプリタの提供するプロファイルのための関数を設定する関数で、cProfileはこの言語処理系の機能をそのまま使っています。ちょっと面白いのは、この関数は外部に公開された、つまり(中の人でなくても)誰でも使えるC言語の関数なことで、やろうと思えば誰でも「ぼくのかんがえたさいきょうのプロファイラ」が作れるということになります。

さて、そのPyEval_SetProfileのソースも見てみましょう:

// https://github.com/python/cpython/blob/8479a3426eb7d1840473f7788e639954363ed37e/Python/ceval.c#L4371

void
PyEval_SetProfile(Py_tracefunc func, PyObject *arg)
{
    PyThreadState *tstate = _PyThreadState_GET();
    PyObject *temp = tstate->c_profileobj;
    Py_XINCREF(arg);
    tstate->c_profilefunc = NULL;
    tstate->c_profileobj = NULL;
    /* Must make sure that tracing is not ignored if 'temp' is freed */
    tstate->use_tracing = tstate->c_tracefunc != NULL;
    Py_XDECREF(temp);
    tstate->c_profilefunc = func;
    tstate->c_profileobj = arg;
    /* Flag that tracing or profiling is turned on */
    tstate->use_tracing = (func != NULL) || (tstate->c_tracefunc != NULL);
}

一行目のPyThreadState_GETはスレッドごとに違うPyThreadStateオブジェクトを返す関数(参考:スレッド局所記憶)です。なので、プロファイラの関数はスレッドごとに設定されます。というわけで、python -m cProfile test.pyで測定した時は、メインスレッド以外の別のスレッドに関してはやっぱり一切計測しません

tstate->c_profilefuncに設定される関数は、C言語とPythonの関数のCALL/RETURNの時に呼ばれるようになっていまして、cProfileでプロファイルを取る時の関数もその情報を使って集計しているのがわかります。

実際にこの関数が呼ばれているところも追いかけておきましょう。

PythonからPythonの関数の呼び出し

Pyrhonの関数のCALL/RETURNに関しては、関数の本体(のバイトコード)を実行する_PyEval_EvalFrameDefaultというインタプリタの核心部みたいな関数の最初と最後で呼ばれます:

# https://github.com/python/cpython/blob/842a2f07f2f08a935ef470bfdaeef40f87490cfc/Python/ceval.c#L940

        if (tstate->c_profilefunc != NULL) {
            /* Similar for c_profilefunc, except it needn't
               return itself and isn't called for "line" events */
            if (call_trace_protected(tstate->c_profilefunc,
                                     tstate->c_profileobj,
                                     tstate, f, PyTrace_CALL, Py_None)) {
                /* Profile function raised an error */
                goto exit_eval_frame;
            }
        }
    }

# https://github.com/python/cpython/blob/842a2f07f2f08a935ef470bfdaeef40f87490cfc/Python/ceval.c#L3575

return_or_yield:
    if (tstate->use_tracing) {
        if (tstate->c_tracefunc) {
            if (call_trace_protected(tstate->c_tracefunc, tstate->c_traceobj,
                                     tstate, f, PyTrace_RETURN, retval)) {
                Py_CLEAR(retval);
            }
        }
        if (tstate->c_profilefunc) {
            if (call_trace_protected(tstate->c_profilefunc, tstate->c_profileobj,
                                     tstate, f, PyTrace_RETURN, retval)) {
                Py_CLEAR(retval);
            }
        }
    }

この関数は、Pythonで書かれた関数を実際に実行するコードです。バイトコードの実行ループがあるので(今回の話とは一切無関係ですが)読んでいて面白い関数ではあります。おすすめです。

PythonからC言語の関数呼び出しは別扱い

C言語の関数の呼び出しは、もうちょっと複雑です。Pythonの内部で使われているVMのオペコードのうち、CALL_METHOD/CALL_FUNCTION/CALL_FUNCTION_KW/CALL_FUNCTION_EX命令が実行された時に呼ばれる、次の関数からC言語の関数の実行が移譲されています:

// https://github.com/python/cpython/blob/842a2f07f2f08a935ef470bfdaeef40f87490cfc/Python/ceval.c#L4648

/* Issue #29227: Inline call_function() into _PyEval_EvalFrameDefault()
   to reduce the stack consumption. */
Py_LOCAL_INLINE(PyObject *) _Py_HOT_FUNCTION
call_function(PyObject ***pp_stack, Py_ssize_t oparg, PyObject *kwnames)
{

    // ... 前処理 ...

    /* Always dispatch PyCFunction first, because these are
       presumed to be the most frequent callable object.
    */
    if (PyCFunction_Check(func)) {
        PyThreadState *tstate = _PyThreadState_GET();
        C_TRACE(x, _PyCFunction_FastCallKeywords(func, stack, nargs, kwnames));
    }
    else if (Py_TYPE(func) == &PyMethodDescr_Type) {
        PyThreadState *tstate = _PyThreadState_GET();
        if (nargs > 0 && tstate->use_tracing) {
            /* We need to create a temporary bound method as argument
               for profiling.
               If nargs == 0, then this cannot work because we have no
               "self". In any case, the call itself would raise
               TypeError (foo needs an argument), so we just skip
               profiling. */
            PyObject *self = stack[0];
            func = Py_TYPE(func)->tp_descr_get(func, self, (PyObject*)Py_TYPE(self));
            if (func != NULL) {
                C_TRACE(x, _PyCFunction_FastCallKeywords(func,
                                                         stack+1, nargs-1,
                                                         kwnames));
                Py_DECREF(func);
            }
            else {
                x = NULL;
            }
        }
        else {
            x = _PyMethodDescr_FastCallKeywords(func, stack, nargs, kwnames);
        }
    }
    else {
        // ... Pythonの関数のときの処理 ...
    }

    // ... 後処理 ...

    return x;
}

この最初の2つのif文がC言語で書かれた関数ないしメソッドの呼び出しに対応しているようです。2番めのPy_TYPE(func) == &PyMethodDescr_Type の方は、調べた限りでは、

>>> po = "po!"
>>> f = po.upper
>>> f()
'PO!'

のような感じで、メソッドにselfだけ束縛した場合に実行される分岐のようです(あんまり自信ない)。

C_TRACEというマクロの定義は次の通りで、callという引数で指定されるC言語関数呼び出しコールの前後でPyTrace_C_CALLとPyTrace_C_RETURNをプロファイラのコールバック関数に対して通知するようになっています。

// https://github.com/python/cpython/blob/842a2f07f2f08a935ef470bfdaeef40f87490cfc/Python/ceval.c#L4617

#define C_TRACE(x, call) 
if (tstate->use_tracing && tstate->c_profilefunc) { 
    if (call_trace(tstate->c_profilefunc, tstate->c_profileobj, 
        tstate, tstate->frame, 
        PyTrace_C_CALL, func)) { 
        x = NULL; 
    } 
    else { 
        x = call; 
        if (tstate->c_profilefunc != NULL) { 
            if (x == NULL) { 
                call_trace_protected(tstate->c_profilefunc, 
                    tstate->c_profileobj, 
                    tstate, tstate->frame, 
                    PyTrace_C_EXCEPTION, func); 
                /* XXX should pass (type, value, tb) */ 
            } else { 
                if (call_trace(tstate->c_profilefunc, 
                    tstate->c_profileobj, 
                    tstate, tstate->frame, 
                    PyTrace_C_RETURN, func)) { 
                    Py_DECREF(x); 
                    x = NULL; 
                } 
            } 
        } 
    } 
} else { 
    x = call; 
}

ここ以外にももっと複雑な関数呼び出しのためのdo_call_core関数もあるのですが、そっちもだいたい似た感じなので省略します。

Pythonの関数のCALL/RETURNは、呼び出された関数の中で記録されるのに対して、C言語の方は呼び出す直前なところがちょっと面白いところでしょうか。もちろん、C言語の関数はPythonのインタプリタから遠く独立していることを踏まえれば、合理的な実装だと思います。

オチはなし:意外と素朴な実装だったなぁ(小学生並みの感想)

…えー、とくにオチはないです。強いて言うなら、意外と素朴な実装でしたね。Linuxカーネルのよくわかんない機能とか使ってるのかなぁという予想も事前にしていたのですが、そんな事はなくて、地道に適当なタイミングで適切に関数を呼び出して記録しているだけでした。

cProfileとprofileの具体的な違いは?

Pythonでプロファイルを取る方法は、cProfileだけではなく、profileというのもあります。こちらは100% Pure Pythonで、cProfileが使えない環境ではこちらを使ってね、とドキュメントには書いてあります。こちらはどのような実装になっているのでしょうか?Pythonのインタプリタの奥深くにプロファイルのための機能が実装されている以上、profileもこれを使っていそうな気はするんですが…。

コードを調べると:

# https://github.com/python/cpython/blob/ad1a25f499362eaf9cbfcafa0b8e2454eb43dcf1/Lib/profile.py#L418

    def runctx(self, cmd, globals, locals):
        self.set_cmd(cmd)
        sys.setprofile(self.dispatcher)
        try:
            exec(cmd, globals, locals)
        finally:
            sys.setprofile(None)
        return self

案の定、sys.setprofileという、PyEval_SetProfileのPython版を使って、C言語の関数ではなくPythonの関数をcallbackとしているだけのように読めます。

1つの引数を取る関数をプロファイルしたい場合、次のようにできます:

import cProfile
import re
cProfile.run('re.compile("foo|bar")')

(お使いのシステムで cProfile が使えないときは代わりに profile を使って下さい)

Python プロファイラ — Python 3.7.3rc1

cProfileのコードはそんなにトリッキーな事をしているようには読めないので、cpythonがコンパイルできる環境ならcProfileも使えると思うんですが、cProfileが使えない環境は実際にあるのだろうか…。

いかがでしたか?

「だからなんだ」と言われればそれでオシマイなんですが、ソースコードを読むのはやっぱり重要だと痛感しました。ドキュメントにはすべては書いてないので。この記事が、いつか、cProfileの挙動に違和感を持った人の参考になりますように。

というわけで、次の記事では、マルチスレッドのPythonプログラムをプロファイルする方法を探していきます。

DCGANでテキスト生成を試してみた

こんにちは、Link-Uの町屋敷です。

今回は生成モデルで有名なGANさんで、Word2Vecを使って英語の文章を生成するものを見つけたので、日本語でもできるのかやってみました。

DCGANのサンプルを動かしてみる

今回はMITライセンスのこのプロジェクトのDCGANをフォークして使っていきます。

まずはそのままのコードを確認していきましょう。

DCGANはGANの一種です。
GAN(Generative Adversarial Networks)は、2つの学習器を互いに争わせることで、学習を進めていきます。2つの学習器とは、
入力された乱数からほしいもの(画像や文章など)を作るジェネレーター(以降G:生成器)と
画像や文章などを入力とし、それがGから生成されたなのかデータセットとして収集してきた本物なのかを見分けるディスクリミネーター(D:識別器)です。

サンプルでは、build_generator()でGの、build_discrimanator内でDの構造が定義されています。

今回はDCGANなので、両方の学習器で畳み込み層を使用しますがプーリング層はありません。

このサンプルでは、コンストラクタで使うネットワークの形と学習方法、損失関数などを定義しています。

Dでは、生成された画像と本物の画像を見分けることを目ざします。生成された画像に0のラベルを、本物の画像に1のラベルをつければ、通常の画像分類と同様なので実装は簡単です。

    # Build and compile the discriminator
        self.discriminator = self.build_discriminator()
        self.discriminator.compile(loss='binary_crossentropy',
            optimizer=optimizer2,
            metrics=['accuracy'])

Gもやっていることは簡単で、ノイズを入力して生成された画像をDに入れて、より高い確率で本物の画像たど判別されるように学習します。

プログラムにするとこうです。

何も考えずに作ってしまうと、Gのモデルの中にDが含まれているので、Gを学習するときにDも学習されてしまいます。

そこでGを学習するときにはDの値を変えないようself.discriminator.trainable = Falseで制限をかけています。

        # The generator takes noise as input and generates imgs
        z = Input(shape=(self.latent_dim,)) #雑音の形の定義
        img = self.generator(z) #生成された画像がimg

         # For the combined model we will only train the generator
        self.discriminator.trainable = False

        # The discriminator takes generated images as input and determines validity
        valid = self.discriminator(img) #Validが判定結果

        # The combined model  (stacked generator and discriminator)
        # Trains the generator to fool the discriminator
        self.combined = Model(z, valid) #入力が雑音で、出力が判定結果のモデル
        self.combined.compile(loss='binary_crossentropy', optimizer=optimizer)

あとは、train関数内でデータを入れたら完成です。

実際に走らせてみましょう。

0Epoch

ただの雑音が

200Epoch

4000Epoch

読める数字になっていっています。

文章生成用に改造する

さてここからが本題です。

先程のサンプルでは28*28の画像を入力にしていましたが、代わりに文字を変換して入力します。

そのままの文字を入力してもどうにもならないので、何らかの変換をしないといけないんですが、それが今回使うWord2Vecです。

正解の文章は前にWikipediaから取ってきた漫画の記事(約5000記事)のデータを使いまわします。

Word2Vecはその名の通り単語をベクトルに変換する技術で、これを使うことによって単語と単語の計算ができるようになります。

        all_sentences = joblib.load('{0}/all_sentences.pkl'.format(WRITE_JOBLIB_DIR))
        sentences = [[word for word in document.lower().split()] for document in all_sentences]
        
        print("Building Word2Vec")
        word_model = Word2Vec(sentences, size=63, min_count=1, window=5) 
        joblib.dump(word_model, '{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))

変換はgensimを使って行います。Word2Vec関数に文章を渡せば辞書が出来ます。。簡単!

先にMecabでの分かち書きを忘れないよう注意。引数のsizeは変換後のベクトルの長さ、windowは。文章中の単語と単語の関係性を計算する距離を表す。この例だと5単語まで。

def VectorizeWord(self):
        all_sentences = joblib.load('{0}/all_sentences.pkl'.format(WRITE_JOBLIB_DIR))
        sentences = [[word for word in document.lower().split()] for document in all_sentences]
        word_model = joblib.load('{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
               
        n_words = 8
        # converted_sentences = [] 
        # converted_sentence = np.zeros(word_model.syn0.shape[0])
        input_data = []
        for s in sentences:
            vectorized_sentence = []
            word_count = 0
            for w in s:
                vector = word_model.wv[w]
                vectorized_sentence.append(vector)
                word_count += 1
                if w == '。':
                    if word_count < n_words:
                        for i in range(n_words - word_count):
                            vector = word_model.wv['。']
                            vectorized_sentence.append(np.append(vector, -1))
                            word_count += 1
                    if word_count == n_words:
                        input_data.append(vectorized_sentence)
                        vectorized_sentence = []
                        word_count = 0
                    elif word_count > n_words:
                        vectorized_sentence = []
                        word_count = 0
        print(np.shape(input_data))
        print(n_words * (word_model.layer1_size + 1))
        input_data = np.reshape(input_data, (-1, n_words * (word_model.layer1_size + 1)))  # 品詞IDの分がたされる
        print(np.shape(input_data))
        print(n_words * len(sentences))
        word_model = joblib.dump(input_data, '{0}/vectorized_words.pkl'.format(WRITE_JOBLIB_DIR))
        print(1)
            

作った辞書を使って単語をベクトルに変換していきます。(vector = word_model.wv[w])の部分。

ついでに1文章に含まれる単語数はまちまちなので、ある一定数に揃えます。一定数に満たない文章は空白で埋めたかったんですが、辞書に登録されてない単語は変換するときにエラーになるので「。」で埋めました。「行き,まし,た,。」 => 「行き,まし,た,。,。,。,。,。」

これで準備完了。先程の画像の代わりに変換した単語を入れていきます。

さてここで問題になるのは、どの形で、単語を渡すかです。

一つは、一単語を行ベクトルを列として(単語数*ベクトル数*1)の画像みたいに渡す方法。

ベクトルを画像で言うRGBみたいなものだと考えて(1*単語数*ベクトル数)もあります。

結論から言うと前者はあまりよろしくなかったです。

生成例がこれなんですけど

epoch0

・ メルチェリーダぷれりゅうどぷれりゅうどガンムマーズクロニクル歌野ゃがんはがちりんにとぶアプトアプトばろそぐいアダルトグッズショップばろケルビム・ゲート・キーパーゃがんはがちりんにとぶキザキ売り込もゃがんはがちりんにとぶゃがんはがちりんにとぶオゲ売り込もキザキ売り込も受け流そ歌野歌野ザ・ビーンズ歌野臨もザ・ビーンズ歌野ザ・ビーンズメルチェリーダ弖虎弖虎弖虎弖虎弖虎メルチェリーダアオノキセキウキウキペディアカードウキウキペディアカード

epoch400

・ 使用二鷹野大遠藤広隆小暮昌広程度桁異なり本当に他紙足取り後ろ他付属ボーナストラックトラック柔弱。。。

epoch0でも兆候が見えてるですが、epoch400ではなんか前の方に名詞固まってるし、トラックトラックみたいな繰り返しがよく起こってます。

とくにこの繰り返しが至るところに出現します。多分ゲネレーターの畳み込み層のせいだと思ったので、入力を後者の方に変更。単語数64も多すぎたので8まで減らした。

結果はこちら

epoch 0

・ 若作り繋っシャンツァカッシーギャグイニチアシブノレ混交デイブレイク

・ 年代デイブレイクビジネスインフラ・リミテッドオフィシャルパンフレット商業ドルドルーヴォドラゴン・サーガジャケットイラストギャラリー

・ 交ぜるスパークオンウェイヴオフィシャルパンフレット描くソガシイナ京田辺イリヤルートドラゴンボールメニュー

・ 福見オールドマン・パーサキュライナーツノートぬおイニチアシブゲームエッセイマンガ木城

・ あかしあ台ノスタルジーゅせんきょうジミナニードルアイフレスポ鬩ホラー

epoch 100

・ 選評再現カレイドスター・ウィンディステージ数。。。。

・ ホワイトナックル・クリムゾンオーブハンター両替キャストバトルオブフェアリーテイル選り抜いゲームパッケージイラスト全て。

・ カードダスステーション終盤ぶった斬る以後後半フルタッチ・アクション。。

・ 温州ニコラ・ケフェウス一般に約前後。。。

・ 時点キディ・ガーランドとおり前後予定同権。。

epoch 200

・ キャラウムカフェウルトラジャンプエッグサイト差し置き。日本一が上回り驚愕

・ ワイディーネ学期エロティック・コメディ。メジャーが上回りから

・ キャラウムカフェ学期エロティック・コメディガールズ・デート年代は上回り共に

・ ドルビーステレオコンサンタルト学期エロティック・コメディガールズ・デートアニヴァーサリーが上回り共に

・ キャラウムカフェ学期エロティック・コメディガールズ・デート以後が上回り奪え

さっきよりは多少マシになったが、まだまだ不自然。

epoch200ではモード崩壊気味でなぜか末尾の「。」が消えた。

ここから改善を始めていくが続きは来月に。

まとめ

一応生成は出来たが精度はまだまだ。

改善点は色々あるので来月はそれをやっていこうと思う。

具体的にはパラメーター調整とかMeCabの情報を使ってConditional GAN化してみるとかそもそもDCGANをやめるとかかなあ。

プログラム全文

'''
MIT License

Copyright (c) 2017 Erik Linder-Norén

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
'''

import re
import json

from timeit import default_timer as timer

import keras.backend as K

import numpy as np
import matplotlib
from matplotlib import pyplot as plt
import joblib

from natto import MeCab
import gensim
from gensim.models.word2vec import Word2Vec

try:
    from keras.engine.topology import Container
except:
    from keras.engine.network import Network as Container
from keras.datasets import mnist
from keras.layers.recurrent import LSTM
from keras.layers.embeddings import Embedding
from keras.models import Model, Sequential
from keras.layers import Dense, Activation, Reshape, UpSampling2D, Conv2D, BatchNormalization, Input, Dropout, ZeroPadding2D, Flatten
from keras.optimizers import Adam
from keras.layers.advanced_activations import LeakyReLU
from keras.losses import binary_crossentropy

TEXT_JSON_DIR = '../WikipediaComic/whole_data'
INFOBOX_JSON_DIR = '.'
INFOBOX_FILE_NAME = 'wiki_infobox_Infobox_animanga_Manga.json'
WRITE_JSON_DIR = '.'
WRITE_JSON_FILE_NAME = 'joined.json' 
WRITE_TEXT_RESULT_DIR = '.'
AUTO_ENCODER_DIR = '.'

WRITE_JOBLIB_DIR = '.'

# Matplotlibの日本語設定
font_path = '/usr/share/fonts/truetype/takao-gothic/TakaoPGothic.ttf'
font_prop = matplotlib.font_manager.FontProperties(fname=font_path)
matplotlib.rcParams['font.family'] = font_prop.get_name()

def cross_entropy_plus_simmilarity(y_true, y_pred, X_pred):
    simi = 0
    for i , v in enumerate(X_pred):
        for j, _ in enumerate(X_pred):
            if i > j:
                simi += K.square(X_pred[i] - X_pred[j])
    n = len(X_pred)
    return simi/(n*(n-1)/2) + binary_crossentropy(y_true, y_pred)

class DCGAN():

    def __init__(self, row, col, additional_featire_count = 0):
        self.additional_feature_count = additional_featire_count
        # Input shape
        self.img_rows = 1
        self.img_cols = col
        self.channels = row + self.additional_feature_count
        
        self.txt_shape = (self.img_cols, 1, self.channels)
        self.latent_dim = 200

        optimizer = Adam(0.0002, 0.5)

        # Build the generator
        self.generator = self.build_generator()
#        self.generator = self.WordGenerator()

        # Build and compile the discriminator
        self.discriminator = self.build_discriminator()
        self.discriminator.compile(loss='binary_crossentropy',
            optimizer=optimizer,
            metrics=['accuracy'])

        # The generator takes noise as input and generates imgs
        z = Input(shape=(self.latent_dim,))
        img = self.generator(z)

        # For the combined model we will only train the generator
        self.discriminator.trainable = False

        # The discriminator takes generated images as input and determines validity
        valid = self.discriminator(img)

        # The combined model  (stacked generator and discriminator)
        # Trains the generator to fool the discriminator
        self.combined = Model(z, valid)
        self.combined.compile(loss='binary_crossentropy', optimizer=optimizer)
        
    def build_generator(self):
    
        model = Sequential()
    
        model.add(Dense(128 * int((self.img_cols) * self.img_rows), activation="relu", input_dim=self.latent_dim))
        model.add(Reshape((self.img_cols , self.img_rows, 128)))
        model.add(Conv2D(128, kernel_size=3, padding="same"))
        model.add(BatchNormalization(momentum=0.8))
        model.add(Activation("relu"))
        model.add(Conv2D(64, kernel_size=3, strides=1, padding="same"))
        model.add(BatchNormalization(momentum=0.8))
        model.add(Activation("relu"))
        model.add(Conv2D(32, kernel_size=3, strides=1, padding="same"))
        model.add(BatchNormalization(momentum=0.8))
        model.add(Activation("relu"))
        model.add(Conv2D(self.channels, kernel_size=3, padding="same"))
        model.add(Activation("tanh"))
    
        model.summary()
    
        noise = Input(shape=(self.latent_dim,))
        img = model(noise)
    
        return Model(noise, Container(noise, img)(noise)) 

    def build_discriminator(self):

        model = Sequential()
        print(self.txt_shape)
        model.add(Conv2D(32, kernel_size=3, strides=2, input_shape=self.txt_shape, padding="same"))
        model.add(LeakyReLU(alpha=0.2))
        model.add(Conv2D(64, kernel_size=3, strides=2, padding="same"))
        model.add(ZeroPadding2D(padding=((0, 1), (0, 1))))
        model.add(BatchNormalization(momentum=0.8))
        model.add(LeakyReLU(alpha=0.2))
        model.add(Dropout(0.5))
        model.add(Conv2D(128, kernel_size=3, strides=2, padding="same"))
        model.add(BatchNormalization(momentum=0.8))
        model.add(LeakyReLU(alpha=0.2))
        model.add(Flatten())
        model.add(Dense(1, activation='sigmoid'))

        model.summary()

        img = Input(shape=self.txt_shape)
        validity = model(img)

        return Model(img, Container(img, validity)(img))

    def train(self, epochs, batch_size=128, save_interval=50):
        word_model = joblib.load('{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
        X_train = joblib.load('{0}/vectorized_words.pkl'.format(WRITE_JOBLIB_DIR))
        print(np.shape(X_train))
        X_train = np.reshape(X_train, (-1, self.img_cols , 1, self.channels))
        print(np.shape(X_train))
        # Adversarial ground truths
        real = np.ones((batch_size, 1))
        fake = np.zeros((batch_size, 1))

        for epoch in range(epochs):

            # ---------------------
            #  Train Discriminator
            # ---------------------

            # Select a random half of images
            idx = np.random.randint(0, X_train.shape[0], batch_size)
            imgs = X_train[idx]

            # Sample noise and generate a batch of new images
            noise = np.random.normal(0, 1, (batch_size, self.latent_dim))
            gen_imgs = self.generator.predict(noise)

            # Train the discriminator (real classified as ones and generated as zeros)
            self.discriminator.trainable = True
            self.generator.non_trainable = False
            d_loss_real = self.discriminator.train_on_batch(imgs, real)
            d_loss_fake = self.discriminator.train_on_batch(gen_imgs, fake)
            d_loss = 0.5 * np.add(d_loss_real, d_loss_fake)

            # ---------------------
            #  Train Generator
            # ---------------------

            # Train the generator (wants discriminator to mistake images as real)
            self.discriminator.trainable = False
            self.generator.non_trainable = True
            g_loss = self.combined.train_on_batch(noise, real)

            # Plot the progress
            print ("%d [D loss: %f, acc.: %.2f%%] [G loss: %f]" % (epoch, d_loss[0], 100 * d_loss[1], g_loss))

            # If at save interval => save generated image samples
            if epoch % save_interval == 0:
                self.show_text(epoch)
    
    def show_text(self, epoch):
        r, c = 5, 5
        noise = np.random.normal(0, 1, (r * c, self.latent_dim))
        gen_text = self.generator.predict(noise)
        word_model = joblib.load('{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
        r, c = 5, 5
        noise = np.random.normal(0, 1, (r * c, self.latent_dim))
        gen_text = self.generator.predict(noise)
        word_model = joblib.load('{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
        with open(WRITE_TEXT_RESULT_DIR + '/generated_epoch{0}.txt'.format(epoch), "w+") as f:
            f.write('Epoch {0}rn'.format(epoch))
            for vector in gen_text:
                word = ''
                s = vector.flatten().reshape(self.img_cols, self.img_rows, self.channels)
                for w in s:   
                    if self.additional_feature_count > 0:
                        word = word + word_model.most_similar(w[:,0:-1*self.additional_feature_count])[0][0]
                    else:
                        word = word + word_model.most_similar(w)[0][0]
                print(word + 'rn')
                f.write(word)
    
    def save_imgs(self, epoch):
        r, c = 5, 5
        noise = np.random.normal(0, 1, (r * c, self.latent_dim))
        gen_imgs = self.generator.predict(noise)

        # Rescale images 0 - 1
        gen_imgs = 0.5 * gen_imgs + 0.5

        fig, axs = plt.subplots(r, c)
        cnt = 0
        for i in range(r):
            for j in range(c):
                axs[i, j].imshow(gen_imgs[cnt, :, :, 0], cmap='gray')
                axs[i, j].axis('off')
                cnt += 1
        fig.savefig("images/mnist_%d.png" % epoch)
        plt.close()


class ProcessWord():

    def __init__(self, word_max, vector_size, do_use_pos = True):
        self.mecab = MeCab()
        self.word_max = word_max
        self.vector_size = vector_size
        self.do_use_pos = do_use_pos
        
    def DcganTrigger(self, epochs=4000, batch_size=64, save_interval=10):
        if self.do_use_pos:
            dg = DCGAN(self.vector_size, self.word_max, 1)
        else:
            dg = DCGAN(self.vector_size, self.word_max, 0)
        dg.train(epochs=epochs, batch_size=batch_size, save_interval=save_interval)

    def ExtractWords(self):
        with open(WRITE_JSON_DIR + '/' + WRITE_JSON_FILE_NAME, 'r', encoding='utf_8') as jw: 
            json_data = json.load(jw)
        
        all_words = [[0]] * len(json_data) 
        with MeCab("-Owakati") as mecab:
            for i, j in enumerate(json_data):
                if i % 100 == 0:
                    print(i)
                text = j['text']
                
                all_words[i] = re.sub(r'[!-~]', "", mecab.parse(text))  # 記号、英語除去
                
        joblib.dump(all_words, '{0}/all_sentences.pkl'.format(WRITE_JOBLIB_DIR), compress=True)
        
    def CreateWord2Vec(self):
        all_sentences = joblib.load('{0}/all_sentences.pkl'.format(WRITE_JOBLIB_DIR))
        sentences = [[self.AnalyzeWord(word)[1] for word in document.lower().split()] for document in all_sentences]
        
        print("Building Word2Vec")
        word_model = Word2Vec(sentences, size=self.vector_size, min_count=1, window=5)
        joblib.dump(word_model, '{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
        # Code tried to prepare LSTM model for word generation
        
    def VectorizeWord(self):
        all_sentences = joblib.load('{0}/all_sentences.pkl'.format(WRITE_JOBLIB_DIR))
        sentences = [[word for word in document.lower().split()] for document in all_sentences]
        word_model = joblib.load('{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
               
        n_words = self.word_max
        # converted_sentences = [] 
        # converted_sentence = np.zeros(word_model.syn0.shape[0])
        input_data = []
        pos_labels = [] 
        for s in sentences:
            vectorized_sentence = []
            pos = []
            word_count = 0
            #print(s)
            for w in s:
                posid, surface = self.AnalyzeWord(w)
                vector = word_model.wv[surface]
                if self.do_use_pos:
                    vector = np.append(vector, float((posid - 34) / 68))  # POSID正規化
                pos.append(float((posid - 34) / 68))
                vectorized_sentence.append(vector)
                word_count += 1
                if w == '。':
                    if word_count < n_words:
                        for i in range(n_words - word_count):
                            vector = word_model.wv['。']
                            if self.do_use_pos:
                                vectorized_sentence.append(np.append(vector, -1))
                            else:
                                vectorized_sentence.append(vector)
                            pos.append(-1)
                            word_count += 1
                    if word_count == n_words:
                        input_data.append(vectorized_sentence)
                        pos_labels.append(pos)
                        vectorized_sentence = []
                        pos =  []
                        word_count = 0
                    elif word_count > n_words:
                        vectorized_sentence = []
                        pos = []
                        word_count = 0
                        
        print(np.shape(input_data))
        print(n_words * (word_model.layer1_size + 1))
        if self.do_use_pos:
            input_data = np.reshape(input_data, (-1, n_words * (word_model.layer1_size + 1)))  # 品詞IDの分がたされる
        else:
            input_data = np.reshape(np.array(input_data), (-1, n_words * (word_model.layer1_size)))
        print(np.shape(input_data))
        print(n_words * len(sentences))
        word_model = joblib.dump(input_data, '{0}/vectorized_words.pkl'.format(WRITE_JOBLIB_DIR))
        joblib.dump(input_data, '{0}/vectorized_words.pkl'.format(WRITE_JOBLIB_DIR))
        joblib.dump(pos_labels, '{0}/pos_labels.pkl'.format(WRITE_JOBLIB_DIR))
        print(1)
        
    def CheckWord2Vec(self, word1, word2, ope):
        word_model = joblib.load('{0}/word2vec.pkl'.format(WRITE_JOBLIB_DIR))
        if ope == '+':
            w = word_model.wv[word1] + word_model.wv[word2]
        if ope == '-':
            w = word_model.wv[word1] - word_model.wv[word2]
        if ope == 'all' or ope == 'ALL':
            self.CheckWord2Vec(word1, word2, '+')
            self.CheckWord2Vec(word1, word2, '-')
            self.CheckWord2Vec(word2, word1, '-')
            return    
        print(word_model.most_similar([w]))
        
    def AnalyzeWord(self, word):
        gen = self.mecab.parse(word, as_nodes=True)
        for w in gen:
            return w.posid, w.surface

if __name__ == '__main__':
    begin = timer()
    pw = ProcessWord(8, 32, False)    
    #pw.ExtractWords()
    #pw.CreateWord2Vec()    
    pw.VectorizeWord()
    pw.DcganTrigger(4000, 256, 10)
    print(timer() - begin)
    pw.CheckWord2Vec('学院', '悪', 'ALL')

GPUを使って無線LANをクラックする話 / モニター、ロック、スケーラビリティ

弊社は最近、それまでの神谷町(東京タワーの近く)から、神田と秋葉原の間くらいに引っ越しをしました。とは言っても、去年の夏の終わりぐらいなので、だいぶ前なのですが。

オフィスの立地で重要なことって、何でしょう。駅からの距離?ビルの高さ?地名のブランド?歓楽街からの近さ?ランチ営業をやってるレストランの多さ?

みなさま、各々重視する点があると思うのですが、「一番近いコンビニが何か」は意外と見過ごされるところではないでしょうか。朝ちょっと寄ってコーヒーを買ったり、お弁当を買ったり。結構お世話になる所だと思います。

というわけで今月の写真は「セブンイレブンのカフェラテ・マシンの『顔』」です。いままでの所は一番近いコンビニが「セブン・イレブン」だったのが「ファミマ!」になってしまって、この眉毛のかわいい顔が見れなくちゃってさ。ちょっとさびしいよね。

…という世間話(?)をしたら、「へ、へぇ…(困惑)」って言われました。

コミュニケーション、難しいね!

前回までのあらすじ

前回分かってしまったことは、仕事を要求するCPUやGPUが増えれば増えるほど、パスワードの解析を依頼する役割を担うPythonのメイン・スレッドは、ロックの取り合いで時間を使い潰すようになってしまうということでした。ここを改善すれば、CPUやGPUがもっと活用できるかも。

このロックは、実際にはどのように使われているのでしょうか?ちょっとここで、Pyritで仕事が分配される様を眺めてみましょう。

Pyritでの「流れ作業」の様子を見学する

Pyritのベンチマークモードでは、メインのPythonスレッドで候補となるパスワードを生成した後に、それらをCPUやGPUに分配しています。で、並列に計算された結果は再度メインのスレッドでCPyritクラスが受け取る、という流れ作業になっています。その流れを見ていきましょう。

仕事の受付け

CPyrit.enqueueに渡されます。このメソッドは、パスワードを生成するPythonのメインスレッドから呼ばれます:

    def enqueue(self, essid, passwords, block=True):
        with self.cv:
            # ... 省略 ...
            passwordlist = list(passwords)
            if len(self.inqueue) > 0 and self.inqueue[-1][0] == essid:
                self.inqueue[-1][1][self.in_idx] = passwordlist
            else:
                self.inqueue.append((essid, {self.in_idx: passwordlist}))
            self.workunits.append(len(passwordlist))
            self.in_idx += len(passwordlist)
            self.cv.notifyAll()

self.cvはthreading.Conditionで、いわゆるモニターと呼ばれている排他制御のためのしくみ、です。「モニタ」って言葉は改めて聞くことはあんまり無い気がするんですが、Javaでおなじみのwait/notify/notifyAllのアレといえば通じる人もいるかもしれません。なにはともあれ、with self.cv: とすることで、このブロックの中のコードが実行されるスレッドが高々1つになるように制御しています。その後の処理は基本的には仕事をself.inqueueその他に適切にセットしているだけ、です。

適切に情報をセットしたら、self.cv.notifyAll() して他のスレッドにnotify(通知)します。具体的に何が起こるのかというと、with self.cv: の中でself.cv.wait() している他のスレッドが全て起こされます。じゃあwith self.cv: しているのは何箇所かというと、これを含めてなんと7箇所もあります。今回はそのうち、たくさん使われていそうなコードパスを眺めていきます。

仕事の受取り

inqueueに入れられた仕事は、CPUやGPUの各スレッドからCPyrit._gatherというメソッドを呼ぶことで受け取られます:

    def _gather(self, desired_size, block=True, timeout=None):
        t = time.time()
        with self.cv:
            passwords = []
            pwslices = []
            cur_essid = None
            restsize = desired_size
            while True:
                self._check_cores()
                for essid, pwdict in self.inqueue:
                   # ... passwordsを埋め、self.inqueueからそのぶん削除する処理 ...
                if len(passwords) > 0:
                    wu = (cur_essid, tuple(passwords))
                    try:
                        self.slices[wu].append(pwslices)
                    except KeyError:
                        self.slices[wu] = [pwslices]
                    self.cv.notifyAll()
                    return wu
                else:
                    if block:
                        if timeout is not None and time.time() - t > timeout:
                            return None, None
                    else:
                        return None, None
                    self.cv.wait(0.1)

self.inqueueからself.slicesに情報を移しつつ(self.slices[wu] = … )、WorkUnitである(ESSID, passwords)のタプルをCPUやGPUのために返します。この情報を使って、CPUはSSEを、GPUはCUDAやOpenCLを使って結果(results )を計算します。今回は情報のやり取りだけに注目するので、実際にどのように計算をするのかは省略する三分クッキング方式で次に進みましょう。

計算した結果を通知するために使われるのが次の_scatter メソッドです:

結果の受付け

CPUやGPUの各スレッドで計算された結果は、CPyrit._scatterメソッドを使って集められます:

    def _scatter(self, essid, passwords, results):
        assert len(results) == len(passwords)
        with self.cv:
            wu = (essid, passwords)
            slices = self.slices[wu].pop(0)
            if len(self.slices[wu]) == 0:
                del self.slices[wu]
            ptr = 0
            for idx, length in slices:
                self.outqueue[idx] = list(results[ptr:ptr + length])
                ptr += length
            for idx in sorted(self.outqueue.iterkeys(), reverse=True)[1:]:
                res = self.outqueue[idx]
                o_idx = idx + len(res)
                if o_idx in self.outqueue:
                    res.extend(self.outqueue[o_idx])
                    del self.outqueue[o_idx]
            self.cv.notifyAll()

self.slicesの情報を消費・削除しつつ、self.outqueue に結果を載せています。難しいところはないですね。

結果の受取り

最後に、Pythonのメインスレッドが、CPyrit.dequeueメソッドを使ってCPUやGPUの各スレッドがパスワードを受取ります:

    def dequeue(self, block=True, timeout=None):
        t = time.time()
        with self.cv:
            if len(self.workunits) == 0:
                return None
            while True:
                wu_length = self.workunits[0]
                if self.out_idx not in self.outqueue 
                 or len(self.outqueue[self.out_idx]) < wu_length:
                    self._check_cores()
            # ... 十分な数の結果が得られない時に待つ処理 ...
                else:
                    reslist = self.outqueue[self.out_idx]
                    del self.outqueue[self.out_idx]
                    results = reslist[:wu_length]
                    self.out_idx += wu_length
                    self.outqueue[self.out_idx] = reslist[wu_length:]
                    self.workunits.pop(0)
                    self.cv.notifyAll()
                    return tuple(results)

self.outqueue から結果を適切に取り出してself.outqueue から削除して、結果をreturnするだけです。

ここから更にCrackerがクラックに成功したかの判定をするわけですが、benchmarkではその処理は省略されているので、今回もここで止めておきます。

ここまでの観察

self.inqueue とself.outqueue の2本のキューがあることがわかりました。inqueueはパスワードを、outqueueにはパスワードから計算した結果をやり取りするために使われています。

しかし、それらを保護するためのthreading.Conditionは2つ…ではなく、1つのself.cv だけです。つまり、実際にはself.inqueueしかいじっていなくても、self.outqueueのための保護も同時に行っていることになります。…なんだか、無駄な感じがしませんか?それぞれのために2つthreading.Conditionを作って管理できないんでしょうか?

しかしながら、これは言うのは簡単なんですが、実際にやるのは難しそうです。まず、self.inqueueとself.outqueueは実は完全に独立しているわけではありません。_scatterや_gatherでいじられるself.slicesを介して微妙に両者は繋がっています。なので、実はCondition変数を2つそのまま使うだけでは正しくモニターを分割できません。self.slices用の排他制御も作りこめば分割できます。が、そういう事をするとデッドロックが起きがちなので、正直避けたい。

Pythonはシングルスレッドでしか動かないので、それぞれのQueueのためにthreading.Conditionを2つにわけられたとしても、Javaの時と違って2スレッドが並列で動くようになったりはしません。が、同じConditionを取り合うことはなくなるので、性能の向上が見込まれる…かもしれません。ロックは一般にスケーラブルではない(待機するスレッドの数が増えれば性能が悪化する)からです。

今回はプロファイリングの結果から言ってあまり実益は無さそうなのですが、もしConditionが2つに分割できればself.cv.wait() が無駄に起こされる回数が減るのは間違いない、ということはコメントさせてください。上のself.cv.wait() をしているコードを丹念に眺めてもらうと分かるのですが、無駄に起こされる事に備えて、特定の条件が満たされているかチェックして、条件が満たされてなかったら「やっぱ駄目じゃん!」と叫びながらもう一度眠りにつくようになっています(これ自体は一般的な使い方です)。さらに余談ながら、この「無駄に起こされるスレッド」が存在するので、self.cv.notify() ではなく、self.cv.notifyAll() を使わねばなりません。もしself.cv.notify()を使ってしまったら、運悪く無関係なスレッドだけが起きて、もう一度self.cv.wait() で眠りにつき、その後は何も起こらなくなり、永遠に応答しなくなってしまう可能性があります。

待機するのは最大何スレッド?

ついでに、with self.cv:self.cv.wait()で待つことになるスレッドの数を概算しましょう。CPUは48スレッド、GPUは4スレッドと仮定します。

  • enqueue:たかだか1スレッド(Pythonのメインスレッド)
  • _gather:最大52スレッド(CPU48コア + GPU4コア)
  • _scatter:最大52スレッド(CPU48コア + GPU4コア)
  • dequeue:たかだか1スレッド(Pythonのメインスレッド)

各スレッドが同時に2つのメソッドで待機することはないので、全体でもself.cvを取り合うのは最大53スレッドと思って間違いないと思います。

ふーむ。ここまで観察したので、次はwith self.cv:self.cv.wait()で待つ時、なにが行われているのか観察しましょう。このthreading.Conditionのソースは、Python本体にあるのでそちらを参照していきます。

threading.Conditionの中を眺める

with self.cv:って具体的に何をしているの?

この↑文字列を書いた時、実際にはself.cv.__enter__という関数が呼ばれます(これはPythonでのお約束です)。で、実装はこんな感じです

    def __enter__(self):
        return self._lock.__enter__()

self._lockはRLock(再入可能ロック)で、この関数はさらにC言語で実装されておりまして

static PyMethodDef rlock_methods[] = {
    {"acquire",      (PyCFunction)(void(*)(void))rlock_acquire,
     METH_VARARGS | METH_KEYWORDS, rlock_acquire_doc},
//....
    {"__enter__",    (PyCFunction)(void(*)(void))rlock_acquire,
     METH_VARARGS | METH_KEYWORDS, rlock_acquire_doc},
// ...
};

__enter__はacquireと同じ関数で、rlock_acquireです。

その先の実装は結構複雑そうですので一端置いておきます。pthreadとか適当にwrapしてるだけかと思ったけどそんなことは無い。

self.cv.wait()は実際には何を待っているのか?

with self.cv:の中では必ずself.cv.wait()を呼ぶことで他のスレッドが処理するのを待っています。この処理は何をしているのかというと

    def wait(self, timeout=None):
        # ... 前処理 ...
        waiter = _allocate_lock()
        waiter.acquire()
        self._waiters.append(waiter)
        saved_state = self._release_save()
        gotit = False
        try:    # restore state no matter what (e.g., KeyboardInterrupt)
            if timeout is None:
                waiter.acquire()
                gotit = True
            else:
                if timeout > 0:
                    gotit = waiter.acquire(True, timeout)
                else:
                    gotit = waiter.acquire(False)
            return gotit
        finally:
            self._acquire_restore(saved_state)
            if not gotit:
                try:
                    self._waiters.remove(waiter)
                except ValueError:
                    pass

新しくlockを作り、このロックに対して二度acquireをしています。このロックは再入可能ではないので、この2度めのacquireでスレッドは必ず中断されて、ほかのスレッドが代わりに起こされます。で、その起こされた別のスレッドがnotifyした時にロックがreleaseされて、二度目のロック以降の中断されたプログラムが実行される、という仕組みになっている、と。

notify/notifyAllするときにはロックに対するacquireは行われません。releaseだけ行います。前回ロックのacquireで時間を使い切っていることがわかりましたが、つまり、時間を食ってるのはこのwith self.cv:の中のacquireか、self.cv.wait()の中のacquireのどちらか、ということになります。

ロックに掛かる時間の内訳を見る

条件lockのacquireをした回数(掛かった時間)with self.cv:した回数(掛かった時間)self.cv.wait()を呼び出した回数(掛かった時間)
1cpu/0gpu1963回(0.514秒)204回(0.002秒)51回(69.438秒)
48cpu/0gpu7355回(3.202秒)1781回(1.988秒)756回(66.133秒)
0cpu/4gpu5689回(7.122秒)5645回(6.685秒)8回(0.762回)
48cpu/4gpu5113回(44.287秒)4676回(43.910秒)309回(3.161 秒)

回数がいまいち釣り合わない気がするんですが(0cpu/4gpuのself.cv.wait()の呼び出しがたった「8回」って本当かいな??)、仕事をとりあえば取り合うようになるほどwith self.cv: を実行する回数が増え、実行時間も掛かるようになっていくことが読み取れます。で、このwith self.cv:は、新しく生成したロックに対して獲得(acquire)しようとするwait()の時と違って、1つのself.cv._lock を獲得しようとしています。

これらの観察結果から、次の仮定を立てます:「with self.cv: したときに起こるself.cv._lock の取り合いがボトルネックの1つである」。

ここ、なんとかならないんでしょうか。次回実験してみましょう。

余談ですが、CPUだけで実行している時、どちらもself.cv.waitに60秒以上掛かっていますが、これらは実際には概ねthread.Sleepの実行に費やされています(つまりスレッドがやることなくて完全に遊んでいる)。

今月のまとめ

  • Pyritのパスワード候補と結果のリスト(inqueue/outqueue)は、threading.Condition(モニター)を使って排他制御されている
  • モニターは複数のロックを使って制御されている:
    • モニター自体のロック(cv._lock)
    • 待機しているスレッドを起こすためのロック(cv.wait()の中で生成されるwaiterという名前のロック)
  • プロファイリングの結果、モニター自体のロックを獲得しようとするスレッドが増えた結果として時間を使い潰していそう
    • このモニターでの排他制御を改善したら全体の性能も改善してくれないかなぁ

明日から東京では花粉が飛散しはじめるそうです。花粉症の方も、まだそうでない方も、みなさまお気をつけて!

GPUを使って無線LANをクラックする話 Pythonのプロファイルを取る回

近所を散歩していたら、もう梅が咲いていました。今日は立春だそうです。早いなー。

わたしは梅のつぼみが大好きなので、これからはしばらく目が離せない日々が続きそうです。

前回のおさらい

前回得た結果は、「CPUとGPUを同時に使うより、GPUだけを使う方が実は倍ぐらい速い!」という、これまた逆説的な結果でした。

未改造のPyritでは「CPUとGPUを同時」に使うか、「CPUだけ」を使うかなので、なまじ「同時」が「CPUだけ」よりは速いが故に、そんな単純な事にも気づかなかった、と。

うーん…なんでこんな事に…。GPUを追加で使うぶんだけ、せめてちょびっとでも性能が向上してほしかったんですが…。

Pythonのコードのプロファイルを取ろう

これまでにCUDAコードのプロファイルを二回ほど取りました(一回目二回目)が、今日はPythonのコードのプロファイルを取りましょう。もちろんPythonはCPUで実行されるわけですが、さらにPythonは高々1つのスレッドでしか同時に動作しないのでした。

今回Pythonのプロファイルを取ろうと画策するのは、次の観察からです:

  • 2CPU、48コアを全部ぶん回すと48倍の性能がでない。htopを観察する限り、かといってCPUを全部使っているわけでもなさそう
    • Pythonの1コアが各コアに仕事を分配しきれていないと考えると説明できます(わかんないけどね)
  • GPUだけを使うにしても、GPUは使い切っていなさそう
    • ここがよくわからない。GPUが消費しきる分のパスワードは生成できそうなので。
  • ダミーパスワードの生成に何割くらい掛けているのかが気になる

というわけでモリモリ取っていきまっしょい。プロファイルのためのツールは公式で提供されていて、cProfileというのを使えば良さそうです。

% python -m cProfile -s tottime <実行ファイル名> <args...>

とやると、関数の内部の実行に掛かった順番でソートされた結果が帰ってきます。tottimeは、関数の中で他の関数を呼び出した時間は除いた、その関数の純粋な実行時間。cumtimeと入れると、その関数から呼び出した関数の時間も含まれます。

例えば、CPUを1コアだけ使った時の表示はこんな感じになります:

(venv) server01@server01:~/src/Pyrit$ time python2.7 -m cProfile -s tottime ./pyrit benchmark
Pyrit 0.5.1 (C) 2008-2011 Lukas Lueg - 2015 John Mora
https://github.com/JPaulMora/Pyrit
This code is distributed under the GNU General Public License v3+

Running benchmark (985.4 PMKs/s)... -  

Computed 985.36 PMKs/s total.
#1: 'CPU-Core (SSE2/AES)': 1045.2 PMKs/s (RTT 2.8)
         79276 function calls (78919 primitive calls) in 69.771 seconds

   Ordered by: internal time

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
     1604   68.890    0.043   68.890    0.043 {time.sleep}
     1963    0.514    0.000    0.514    0.000 {method 'acquire' of 'thread.lock' objects}
        1    0.158    0.158    0.197    0.197 cpyrit.py:29(<module>)
        1    0.094    0.094   69.571   69.571 pyrit_cli.py:1184(benchmark)
       51    0.028    0.001   69.438    1.362 threading.py:309(wait)
    54392    0.011    0.000    0.011    0.000 {method 'random' of '_random.Random' objects}

...(略)...

real	1m9.902s
user	1m9.087s
sys	0m0.302s

ログの一番上の行を見ることで、time.sleepで68.9秒とたくさん時間を使っていること、次に使っているのがacquireメソッドであること、がわかります。

さらに、

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
        1    0.094    0.094   69.571   69.571 pyrit_cli.py:1184(benchmark)
real	1m9.902s

この二行を見ると、このベンチマークを実際に実行するbenchmark関数が69.5秒使っていること、処理時間が69.9秒で、その前後の初期化などに0.4秒ぐらい使われていること(まぁ、予想通りですよね)、などがわかりますよー、と。

この2つの情報を合計すると、CPUを1コアだけ使った時は、69.5秒中68.9秒、ほとんどtime.sleepで時間を潰していること(1コアのCPUが処理するのを待っていること)がなんとなくわかります。

以下前回みたく条件を変えて実行してみましたが、だらだらログを並べても中々わかりにくいので、表にしてみました:

条件1番時間を使っている関数
(全処理に占める割合)
2番時間3番時間4番時間
1cpu/0gpu{time.sleep}
(99.0%)
{method ‘acquire’ of ‘thread.lock’ objects}
(0.7%)
cpyrit.py(<module>)
(0.2%)
pyrit_cli.py(benchmark)
(0.1%)
48cpu/0gpu{time.sleep}
(92.0%)
{method ‘acquire’ of ‘thread.lock’ objects}
(4.5%)
pyrit_cli.py(benchmark)
(2.21%)
cpyrit.py(<module>)
(0.22%)
0cpu/4gpupyrit_cli.py(benchmark)
(61.3%)
{method ‘flush’ of ‘file’ objects}
(14.6%)
{method ‘acquire’ of ‘thread.lock’ objects}
(12.4%)
cpyrit.py(dequeue)
(4.04%)
48cpu/4gpu{method ‘acquire’ of ‘thread.lock’ objects}
(65.2%)
pyrit_cli.py:1184(benchmark)
(16.5%)
cpyrit.py(dequeue)
(5.0%)
{time.sleep}
(3.8%)

だいたい上位に並ぶのはどれも同じ顔ぶればかりなのですが、下のほうに行けばいくほど、つまり、たくさんのCPUやGPUが仕事を要求するようになればなるほど、sleepよりも{method ‘acquire’ of ‘thread.lock’ objects}の処理時間がガンガン増えていく事がわかります。一番下の一番処理を取り合っているところではなんと6割もロックにつぎ込んでいます。これじゃあ、GPUのパスワードのクラックではなく、ロックを取り合うプログラムを実行していたと言っても過言ではありませんな。

そしてdequeue関数の処理内容が増えていくのも気になりますねぇ。この関数は、生成したパスワードをCPUやGPUが受け取るためのもので、threading.Conditionを使ってマルチスレッドの調停を行っております。

この2つから予想されることは、…おそらくこのthreading.Conditionオブジェクトの中にあるロックを取り合ってるんでしょうね…。

あまりにも闇が深そうなので今日はこのぐらいにしておきましょう。

GPUだけを使っているケース(三行目)については追加でコメントさせてください。benchmark関数が一番時間を使っていることから、一生懸命処理するためのパスワードを用意していることはなんとなく察されるのですが、パスワードの生成に使っているはずのrandom関数がリストアップされてこなくて、そのかわりにfile.flushが二番目に食い込んでいます。このflushはCUDA関係っぽい気がしますが、謎です。もうちょっと追いかける価値があると思います。

まとめ

  • 仕事を要求するスレッドが増えれば増えるほど、仕事を用意するスレッドはsleepではなくlockの取り合いで時間を潰すようになる
    • 1つのロックを取り合ってそう
  • GPUだけで処理しているときはパスワードの生成を頑張るみたいだが、実態は要調査
  • 外へ出ろ、梅の花を見ろ!

今月は低電力モードとなっております。寒いからね、仕方ないね(ごめんなさい)。

[Python] 漫画のWikipediaの説明文から発表年を推定する その2 どの単語が分類に有用かを調べる

こんにちはLink-Uの町屋敷です。

前回はWikipediaの漫画の説明文から発表年を推定しました。

そこそこ推定できましたが、そもそも漫画の説明文から発表年を推定してなにがうれしいかって特に生産性は無いんですよね、

しかし、学習器自体に生産性が皆無でも学習器がどんな基準で推定したかがわかれば何か別の有意義な情報が得られる可能性があります。

例えば今回で言うと入力は文章中に出てくる単語なので、時代とともに出現する単語の傾向を知ることが出来ます。

また、学習器の精度を上げていくには、入力する特徴量と答えの関係を色んな角度から見てある仮説を立てて、特徴量を削ったり加工していくんですが、そのあたりとして使うことも出来ます。

そこで前回はアンサンブル系の学習機を使ったときに学習に利用されるfeature_importanceというパラメーターで単語の重要度を見ましたが、

今回は別の方法でやってみます、そこで使うのがLIMEです。

LIMEはどんなアルゴリズムで学習器を作っても入力された特徴量の重要度を出力してくれるスグレモノで、マルチクラスにも対応しています。

LIMEを使う

LIMEは回帰問題でも使えるんですが、マルチクラス分類のほうがわかりやすいのでそちらを使います。

もともと回帰の問題だったのでこれをある期間ごとに分け直します。

def _ConvertLabelForClf(labels):
    for i , v in enumerate(labels):
        if v <= 1990:
            labels[i] = 0  
        elif v <= 2005:
            labels[i] = 1
        elif v <= 2010:
            labels[i] = 2 
        else:
            labels[i] = 3
    return labels

今回は1990年以前、1996-2005,2006-2010,2011-の4クラスに分けました。

前回はGradiantBoostingの回帰の学習器を作りましたが、今回は分類の学習器を作ります。

パワメータ調整などやることはほぼ変わらないです。4クラスの数はすべて同じにはならないので、不均等データのための処理(アンダーサンプリングなど)を行う必要も考えましたが、何もせずにとりあえず回した結果がそんなに悪くなかったのでやってません。

ライブラリはpipで入ります。

def LearnClassifier():
    
    features = joblib.load('{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    labels  = joblib.load('{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
     
    print(np.shape(features))
    print(np.shape(labels))
    
    tf = TfidfTransformer()
    features = tf.fit_transform(features)
    #features = _CompressWithAutoencoder(features)
    labels   = _ConvertLabelForClf(labels)
    
X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.1, random_state=1234)
    
    _, counts = np.unique(y_train, return_counts=True)
    weights = counts/len(y_train)
    weights[0] =  1 - weights[1] - weights[2] - weights[3]
    
    clf = LGBMClassifier(
        learning_rate =0.1, n_estimators=1000,
        max_depth=3,
        objective='multiclass',
    )
    
    print(weights)
    print('Learning Start')
    time = timer()
    clf.fit(X_train,y_train)
    time = timer() - time
    print('Learning Finishn Time: {0}'.format(time))
    
    joblib.dump(clf, '{0}/gradient_boosting_classifier.pkl'.format(WRITE_JOBLIB_DIR))
    y_pred = clf.predict(X_test)
    #for yt, yp in zip (y_test, y_pred):
    #    print((yt,yp))
    
    cm = confusion_matrix(y_test, y_pred, [0,1,2,3])
    print('n')
    print(cm)
    
    f1 = f1_score(y_test, y_pred, [0,1,2,3], average='macro')
    acx = accuracy_score(y_train, clf.predict(X_train), [0,1,2,3])
    acy = accuracy_score(y_test, y_pred, [0,1,2,3])
    
    print('f1 = {0},train_accuracy = {1}, test_accuracy = {2}'.format(f1, acx, acy))

コンフュージョンマトリックスとF1値はこんな感じ

そこそこといったところ。この学習器を使ってLIMEを試してみよう。

def InspectClassifier(dict_param = [5, 0.1]):
    dictionary = corpora.Dictionary.load_from_text('filtered_dic_below{0}_above{1}.txt'.format(dict_param[0], dict_param[1]))
    
    labels  = joblib.load('{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
    features = joblib.load('{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    clf = joblib.load('{0}//gradient_boosting_classifier.pkl'.format(WRITE_JOBLIB_DIR))
    
    
    X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.1, random_state=1234)

まず、データのロード。train_test_splitのシードには先程の学習時に使ったものと同じものを使うこと。

    from lime import lime_tabular
    from sklearn.pipeline import make_pipeline
    c = make_pipeline(clf)
    
    class_names = ['< 1990', '< 2005', '< 2010', '>= 2010']
    
    id2token = np.empty(len(dictionary.token2id), dtype='U64')
    for k, v in dictionary.token2id.items():
        id2token[int(v)] = k
    
    exp = lime_tabular.LimeTabularExplainer(X_train, feature_names=id2token, class_names=class_names)

特徴の列と単語のマッピングにはとgensimで作った辞書はそのまま使えないので加工する。

LimeTextExplainer()があるのにLimeTabularExplainer()を使ってる理由は、テキストのほうに分かち書きされた原文を入力に要求されるため、

英語ならそのままぶち込めばいいが日本がだといちいち分けて原型に戻してやらないといけない。さすがにめんどくさい。あとマルチバイト対応しているのか謎。

あとは公式ドキュメントに沿ってやるだけ。

    idx = 0 
    for idx in range(100):
        lime_result = exp.explain_instance(X_test[idx], c.predict_proba, num_features=100, labels=(0,1,2,3))
        
        lime_result.save_to_file('./exp/exp_{0}_{1}.html'.format(idx, y_train[idx]))
        print('# {0} finished'.format(idx))

公式は結果の表示に

show_in_notebook()

を使ってるけどiPython入れてないと何も起こらないので注意!

また、入力する単語数が多すぎると

こんなのが出て進まない、単語数を減らすために、gensimの辞書生成のパラメーターを調整して単語数を減らす。

試しに単語数を減らす代わりにオートエンコーダーで次元圧縮して突っ込んでみたけど爆死したのでこっちのほうが確実?

結果としてhtmlファイルで出力されるのでブラウザーでみてみよう。

各文章の分類に利用された単語が影響度の高い順に出力されている。

つまり、1990年以下かどうかを表す一眼左の青いグラフが書かれているところに上から順に現像、全集、ビデオとあるが、これは文章中にこれらの単語が一度も出現しないとき、1990年以前でないクラスに文章が分類される可能性が高くなることを示している。条件が0以下つまり文章に無いときなってしまっているのでややこしいが、要は現像とか全集とかの単語がなかったらそんなに昔の作品じゃないんじゃない?ってこと。左上にどの程度の自信があって分類しているかがわかるのも地味に嬉しい。

まとめ

LIMEを使うとどの特徴量がどれだけ分類に貢献しているかがわかりやすい。SVMとかにも使えるのも良い点。

参考サイト

https://github.com/marcotcr/lime

https://qiita.com/fufufukakaka/items/d0081cd38251d22ffebf

[Python] 漫画のWikipediaの説明文から発表年を推定する

こんにちは。Link-Uの町屋敷です。

今回は、テキストデータを解析する一例として、

前回抽出した漫画のWikipediaの文章データを使って、

入力データを説明文、出力データを発表年として、入力データから出力データを推定して行きたいと思います。

また、入力データのどの要素(今回なら単語)がその回帰や分類に効力があるのかを調べる方法も紹介していきたいです。

インフォボックスから発表年のデータを取得する

まず、正解データとして使用する発表年データをインフォボックスから収集しよう。

Wikipediaのインフォボックスはjsonで要素infoboxに文字列として前回保存した。

これを適当なパーサーで処理したら完了! …とはいかない。

まず、前回保存した本文とインフォボックスjsonを結合する。

結合方法は両方のファイルにtitle情報が含まれているので紐づけするだけ。

def LoadTextJsonGenerator(files):
    for jf in files:
        with open(jf) as f:
            json_data = pd.read_json(f, lines=True)
            yield json_data.to_dict()
            
def JoinJsonData(info_json_data):
    text_json_files = sorted(glob(TEXT_JSON_DIR + '/*/*'))
    
    ltjg = LoadTextJsonGenerator(text_json_files)

    text_data = ltjg.__next__()
    for i, ijd in enumerate(info_json_data):    
        while True:
            if not ijd['title']:
                ijd['title'] = 'NULL'
            id = int(ijd['id'])
            if id in text_data['id'].values():
                data_index = list(text_data['id'].values()).index(id)
                info_json_data[i]['text'] = text_data['text'][data_index]
                if ijd['title'] == 'NULL':
                    info_json_data[i]['title'] = text_data['title'][data_index]
                break
            else:
                #print((i,id)) #確認用
                text_data = ltjg.__next__() #StopIterationしたらどっかバグってる

def Preprocess():
    with open(INFOBOX_JSON_DIR + '/' + INFOBOX_FILE_NAME) as infof:
        info_json_data = json.load(infof)
        
    JoinJsonData(info_json_data)
    with open(WRITE_JSON_DIR  + '/' + WRITE_JSON_FILE_NAME, 'w', encoding='utf_8') as jw: 
        json.dump(info_json_data, jw, ensure_ascii=False)

どうやらインフォボックスの要素名は微妙な表記ブレや要素名だけあってデータが入っていないことが割とよくあるので、

表記ブレに対応しつつ、もったいないけどデータが入っていない漫画のデータを捨てる必要がある。

このへんは正規表現で頑張ったらなんとかなる。今回は開始、発表期間、連載期間、発表号に表記ブレしていた。(ソース全文は一番下)

infobox = j['infobox'][0] #複数ある場合でも最初のもののみを使う
publication_year = re.search('| *[開始|発表期間|連載期間|発表号].*?([1|2][8|9|0]dd)[年|.]', infobox, re.MULTILINE | re.DOTALL) #@UndefinedVariable
if not publication_year:
    labels[i] = -1
    vain_count += 1
else:
    labels[i] = publication_year[1]

これで、labelsに年数が入った。データが入っていないときlabelsに-1を入れっておくことで後で対応する文章を消すことができる。

なかなかきれいなポアソン分布。

まあ一定時間で区切ったカウントデータだしね……

テキストを学習機に入れられる形に変形する

次に入力データを加工する

labelsの使えるデータの数を数えると5000ちょいだった。

この数でニューラルネットにに直接ぶち込んでLSTMとかを使って解析するのは厳しいので、BoWに加工する。

具体的には、文章を単語単位に分離して、そのうちの名詞、形容詞、動詞に番号を付けて、文章中のその各単語の出現回数を数える。

with MeCab() as mecab:

文章を単語単位に分離する手法は形態素解析と言われるが、ライブラリで簡単にできる。有名なのはChasenとかMecabだが、今回はMecabとpythonでMecabを使えるようなするNatto-pyを使う。インストール方法は調べたら山程出てくるので省略。

            words = []
            text = j['text']
            for mp in mecab.parse(text, as_nodes=True):
                if not mp.is_eos():
                    feature_splits = mp.feature.split(',')
                    if feature_splits[0] in ['名詞', '動詞', '形容詞']:
                        if feature_splits[1] in ['数']:
                            continue
                        elif feature_splits[2] in ['人名']:
                            continue
                        elif feature_splits[6] in ['*']:
                            continue
                        words.append(feature_splits[6])

if feature_splits[0] in [‘名詞’, ‘動詞’, ‘形容詞’]:の行で品詞を絞って、

その後、出現する単語のうち3とか2001などの数は、答えを書いている可能性があるので除去、

また、人名もその人が活躍する時期はある程度偏ってるはずなので、ほぼ答えになるじゃんと言う事で削除した。

def MakeDict(all_words):
    dictionary = corpora.Dictionary(all_words)
    print(dictionary.token2id)
    for no_below in [5,20,40]:
        for no_above in [0.1,0.3,0.5]:
            dictionary.filter_extremes(no_below=no_below, no_above=no_above)
            dictionary.save_as_text('filtered_dic_below{0}_above{1}.txt'.format(no_below, no_above))

単語への番号付けは、専用の辞書を作って行う。

これもgensimというライブラリがある。Mecabと同様これも情報は大量にあるので省略(例えばここここ

def MakeDict(all_words):
    dictionary = corpora.Dictionary(all_words)
    print(dictionary.token2id)
    for no_below in [5,20,40]:
        for no_above in [0.1,0.3,0.5]:
            dictionary.filter_extremes(no_below=no_below, no_above=no_above)
            dictionary.save_as_text('filtered_dic_below{0}_above{1}.txt'.format(no_below, no_above))

辞書を作る部分がここで、no_belowやno_aboveで作り分けたが今回のデータセットでは違いはほぼなかったので、no_below=5, no_above=0.1を使うことにする。

def MakeFeatures(make_dict = False, dict_param = [5, 0.1]):
    
    all_words = joblib.load('{0}/all_wordss.pkl'.format(WRITE_JOBLIB_DIR))
    labels = joblib.load('{0}/publication_years.pkl'.format(WRITE_JOBLIB_DIR))
    if (make_dict):
        MakeDict(all_words)
    
    dictionary = corpora.Dictionary.load_from_text('filtered_dic_below{0}_above{1}.txt'.format(dict_param[0], dict_param[1]))
    
    dl = len(dictionary)
    features = []
    for w, l in zip(all_words, labels):
        tmp = dictionary.doc2bow(w)
        dense = list(matutils.corpus2dense([tmp], num_terms=len(dictionary)).T[0])
        if not l == -1:
            features.append(dense)
    features = np.array(features)
    features = np.reshape(features, (-1, dl))
    labels = [int(v) for v in labels if v != -1]
    joblib.dump(features, '{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    joblib.dump(labels, '{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))

ここで、先程作った辞書を使ってBoWに変換し、labelが-1のデータを削除する

ここでTFIDFを使ってもいいが、今回はパス。これでデータセットの成形は完了。

集計した単語のうち数が多かったTop50をおいておく。

データの学習

本来はTPOTとかを駆使していろいろな学習機で調査するべきだが、時間がないので決定木グラディアントブースティングマシーン(以下、GBM)を用いた回帰と比較対象用に線形回帰を行う。

グラディアントブースティングを使用できるライブラリはsklearn,XGBoost,lightgbmと大きく3つあるが、機能の多さと実行の速さを考えるとlightgbmが良い。

def TuneXgboostRgr():
    
    features = joblib.load('{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    labels  = joblib.load('{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
    
    X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.1, random_state=1234)
    
    params = { 
        'learning_rate' : [0.3,0.2,0.1]
    }
    
    gs = GridSearchCV(estimator = LGBMRegressor(
                            num_leaves=31,
                            learning_rate =0.1, 
                            n_estimators=1000,
                            max_depth=9,
                            objective='regression',
                            min_sum_hessian_in_leaf=1
                        ),
                        param_grid = params, 
                        cv=5)

    gs.fit(X_train, y_train)
    print('n')
    print(gs.cv_results_)
    print('n')
    print(gs.best_params_)
    print('n')
    print(gs.best_score_)

まず、ハイパーパラメーターのチューニングをする。paramsにチューニングしたい変数と値を配列で入れるとクロスバリデーションもして一番いいパラメーターを調査してくれる。

ただし、paramsに大量の変数を設定すると永遠に終わらなくなるので、1,2種類ずついれて何回もやる。詳しくはグリッドサーチで検索!

def LearnRegressor(clf_name = 'gbm'):
    from sklearn.linear_model import LinearRegression
    features = joblib.load('{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    labels  = joblib.load('{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
    
    X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.1, random_state=1234)
    
    if clf_name == 'gbm':
        rgr = LGBMRegressor( 
            num_leaves=1000, 
            max_depth=9, 
            learning_rate=0.06, 
            n_estimators=100, 
            objective='regression',
        )
    elif clf_name == 'linear':
        rgr = LinearRegression()
    #rgr = XGBRegressor( n_estimators = 100, learning_rate=0.1)
    
    print('Learning Start')
    time = timer()
    rgr.fit(X_train,y_train)
    time = timer() - time
    print('Learning Finishn Time: {0}'.format(time))
    y_pred = rgr.predict(X_test)
    
    joblib.dump(rgr, '{0}/{1}_regressor.pkl'.format(WRITE_JOBLIB_DIR, clf_name))
    for yt, yp in zip (y_test, y_pred):
        print((yt,yp))
    print(mean_squared_error(y_test, y_pred))

ハイパーパラメーターのチューニングが終了したら、テストデータで評価する。

結果がこちら、左が正解、中央がGBM,右が線形回帰

比べてみると線形回帰はたまに未来や19世紀などのありえない数値を予測している。

実際誤差を計算するとGBM 59.1に対して線形回帰は329.5と大差を付けている。

(  正解, GBM , linear)

(2012, 2008, 2011)
(2008, 2011, 2007)
(2004, 2001, 2002)
(2010, 2006, 2025)
(2004, 2006, 2003)
(2006, 2004, 2012)
(2017, 2012, 2041)
(2006, 2007, 2005)
(2008, 2006, 1999)
(1977, 1985, 1987)
(2003, 2006, 2011)
(2008, 2006, 1986)
(2015, 2005, 1998)
(1992, 1998, 2003)
(2017, 2008, 2012)
(2003, 2006, 2000)
(2010, 2006, 2011)
(2010, 2006, 2011)
(2004, 2003, 2002)
(2010, 2006, 1958)
(1985, 2002, 1997)
(1994, 2003, 1988)
(1990, 1996, 2004)
(1999, 2002, 2006)
(1984, 1995, 2010)
(1987, 2000, 1993)
(1962, 2001, 2000)
(1993, 1992, 2005)
(1997, 2006, 2009)
(1988, 1982, 1989)
(1972, 1989, 1959)
(1987, 1999, 1989)
(1968, 1973, 1967)
(1970, 1982, 1970)
(1991, 2000, 2000)
(1998, 2007, 1994)
(1999, 1998, 1996)
(1994, 1993, 2001)
(1968, 1970, 1872)

回帰を行う学習器を生成できました。

これで、入力データのどの要素がその回帰や分類に効力があるのかを調べることが出来ます。

線形回帰は相関係数を計算すればすぐに見れます。

実はGBMのようなブースティングモデルの場合でも同様に簡単に調査できて、

そのまま重要度として保存されています。(参考)

両者の結果を見比べてみましょう。

線形回帰の方は法則がそんなに見えません、いんするとかいう謎単語も含まれてるし……

でも、機械学習で出した結果です!って言われたらこんなもんかとは思いそう。

GBMのほうは、昭和や平成、ビデオ、インターネットや現像といった明らかに時代に関係のあるものが紛れているのがわかります。

やたら社名が入っているのは何故だろう……

でもどちらの信用度が高いかは明らかでしょう。

まとめ

今回は漫画のWikipediaの文章データを使って、発表年を推定する方法を紹介しました。

学習器がどの要素(今回なら単語)を注視しているかも調べましたがある単語が強く関係あるとわかったところであまり使い道がありません。

次回では、別の方法でもっといい調査の方法を紹介します。

プログラム全文

import json
from glob import glob
from functools import reduce
import re
from timeit import default_timer as timer

import numpy as np
import pandas as pd
from natto import MeCab
import joblib
from gensim import corpora, matutils

import matplotlib
from matplotlib import pyplot as plt
from matplotlib.font_manager import FontProperties

from sklearn.ensemble import GradientBoostingRegressor
from sklearn.ensemble import GradientBoostingClassifier

#Xgboostは重いので使わない(5倍くらい違う)
#from xgboost import XGBRegressor
#from xgboost import XGBClassifier

from lightgbm import LGBMRegressor
from lightgbm import LGBMClassifier

from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error

TEXT_JSON_DIR = '../WikipediaComic/whole_data'
INFOBOX_JSON_DIR = '.'
INFOBOX_FILE_NAME = 'wiki_infobox_Infobox_animanga_Manga.json'
WRITE_JSON_DIR = '.'
WRITE_JSON_FILE_NAME = 'joined.json' 

WRITE_JOBLIB_DIR = '.'

#Matplotlibの日本語設定
font_path = '/usr/share/fonts/truetype/takao-gothic/TakaoPGothic.ttf'
font_prop = FontProperties(fname=font_path)
matplotlib.rcParams['font.family'] = font_prop.get_name()

def LoadTextJsonGenerator(files):
    for jf in files:
        with open(jf) as f:
            json_data = pd.read_json(f, lines=True)
            yield json_data.to_dict()
            
def JoinJsonData(info_json_data):
    text_json_files = sorted(glob(TEXT_JSON_DIR + '/*/*'))
    
    ltjg = LoadTextJsonGenerator(text_json_files)

    text_data = ltjg.__next__()
    for i, ijd in enumerate(info_json_data):    
        while True:
            if not ijd['title']:
                ijd['title'] = 'NULL'
            id = int(ijd['id'])
            if id in text_data['id'].values():
                data_index = list(text_data['id'].values()).index(id)
                info_json_data[i]['text'] = text_data['text'][data_index]
                if ijd['title'] == 'NULL':
                    info_json_data[i]['title'] = text_data['title'][data_index]
                break
            else:
                #print((i,id)) #確認用
                text_data = ltjg.__next__() #StopIterationしたらどっかバグってる

def Preprocess():
    with open(INFOBOX_JSON_DIR + '/' + INFOBOX_FILE_NAME) as infof:
        info_json_data = json.load(infof)
        
    JoinJsonData(info_json_data)
    with open(WRITE_JSON_DIR  + '/' + WRITE_JSON_FILE_NAME, 'w', encoding='utf_8') as jw: 
        json.dump(info_json_data, jw, ensure_ascii=False)
    
def ExtractWords():
    with open(WRITE_JSON_DIR  + '/' + WRITE_JSON_FILE_NAME, 'r', encoding='utf_8') as jw: 
        json_data = json.load(jw)
    
    all_words = [[0]] * len(json_data) 
    with MeCab() as mecab:
        for i, j in enumerate(json_data):
            if i % 100 == 0:
                print(i)
            words = []
            text = j['text']
            for mp in mecab.parse(text, as_nodes=True):
                if not mp.is_eos():
                    feature_splits = mp.feature.split(',')
                    if feature_splits[0] in ['名詞', '動詞', '形容詞']:
                        if feature_splits[1] in ['数']:
                            continue
                        elif feature_splits[2] in ['人名']:
                            continue
                        elif feature_splits[6] in ['*']:
                            continue
                        words.append(feature_splits[6])
            all_words[i] = words
            
    joblib.dump(all_words, '{0}/all_wordss.pkl'.format(WRITE_JOBLIB_DIR), compress = True)

def ExtractLabelFromInfobox():
    with open(WRITE_JSON_DIR  + '/' + WRITE_JSON_FILE_NAME, 'r', encoding='utf_8') as jw: 
        json_data = json.load(jw)    
        labels = [0] * len(json_data)
        vain_count = 0
        for i, j in enumerate(json_data):
            if i % 100 == 0:
                print(i)
            infobox = j['infobox'][0] #複数ある場合でも最初のもののみを使う
            publication_year = re.search('| *[開始|発表期間|連載期間|発表号].*?([1|2][8|9|0]dd)[年|.]', infobox, re.MULTILINE | re.DOTALL) #@UndefinedVariable
            if not publication_year:
                labels[i] = -1
                vain_count += 1
            else:
                labels[i] = publication_year[1]
        joblib.dump(labels, '{0}/publication_years.pkl'.format(WRITE_JOBLIB_DIR), compress = True)
    print(vain_count)
     
def MakeDict(all_words):
    dictionary = corpora.Dictionary(all_words)
    print(dictionary.token2id)
    for no_below in [5,20,40]:
        for no_above in [0.1,0.3,0.5]:
            dictionary.filter_extremes(no_below=no_below, no_above=no_above)
            dictionary.save_as_text('filtered_dic_below{0}_above{1}.txt'.format(no_below, no_above))
            
def MakeFeatures(make_dict = False, dict_param = [5, 0.1]):
    
    all_words = joblib.load('{0}/all_wordss.pkl'.format(WRITE_JOBLIB_DIR))
    labels = joblib.load('{0}/publication_years.pkl'.format(WRITE_JOBLIB_DIR))
    if (make_dict):
        MakeDict(all_words)
    
    dictionary = corpora.Dictionary.load_from_text('filtered_dic_below{0}_above{1}.txt'.format(dict_param[0], dict_param[1]))
    
    dl = len(dictionary)
    features = []
    for w, l in zip(all_words, labels):
        tmp = dictionary.doc2bow(w)
        dense = list(matutils.corpus2dense([tmp], num_terms=len(dictionary)).T[0])
        if not l == -1:
            features.append(dense)
    features = np.array(features)
    features = np.reshape(features, (-1, dl))
    labels = [int(v) for v in labels if v != -1]
    joblib.dump(features, '{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    joblib.dump(labels, '{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
            


def LearnRegressor(clf_name = 'gbm'):
    from sklearn.linear_model import LinearRegression
    features = joblib.load('{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    labels  = joblib.load('{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
    
    X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.1, random_state=1234)
    
    if clf_name == 'gbm':
        rgr = LGBMRegressor( 
            num_leaves=1000, 
            max_depth=9, 
            learning_rate=0.06, 
            n_estimators=100, 
            objective='regression',
        )
    elif clf_name == 'linear':
        rgr = LinearRegression()
    #rgr = XGBRegressor( n_estimators = 100, learning_rate=0.1)
    
    print('Learning Start')
    time = timer()
    rgr.fit(X_train,y_train)
    time = timer() - time
    print('Learning Finishn Time: {0}'.format(time))
    y_pred = rgr.predict(X_test)
    
    joblib.dump(rgr, '{0}/{1}_regressor.pkl'.format(WRITE_JOBLIB_DIR, clf_name))
    for yt, yp in zip (y_test, y_pred):
        print((yt,yp))
    print(mean_squared_error(y_test, y_pred))

def TuneXgboostRgr():
    
    features = joblib.load('{0}/features.pkl'.format(WRITE_JOBLIB_DIR))
    labels  = joblib.load('{0}/labels.pkl'.format(WRITE_JOBLIB_DIR))
    
    X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.1, random_state=1234)
    
    params = { 
        'learning_rate' : [0.3,0.2,0.1]
    }
    
    gs = GridSearchCV(estimator = LGBMRegressor(
                            num_leaves=31,
                            learning_rate =0.1, 
                            n_estimators=1000,
                            max_depth=9,
                            objective='regression',
                            min_sum_hessian_in_leaf=1
                        ),
                        param_grid = params, 
                        cv=5)

    gs.fit(X_train, y_train)
    print('n')
    print(gs.cv_results_)
    print('n')
    print(gs.best_params_)
    print('n')
    print(gs.best_score_)

if __name__ == '__main__':
    #Preprocess()
    #ExtractWords()
    #ExtractLabelFromInfobox()
    #MakeFeatures()
    #LearnRegressor()
    InspectRegressor()
    #TuneXgboostRgr()
    
    
    
    
    

参考サイト

https://qiita.com/conta_/items/4b031a44acceb137ec73https://yubais.net/doc/matplotlib/bar.html

https://qiita.com/buruzaemon/items/975027cea6371b2c5ec3https://qiita.com/hoto17296/items/e1f80fef8536a0e5e7dbhttps://qiita.com/yasunori/items/31a23eb259482e4824e2

GPUを使って無線LANをクラックする話:帰ってきたPyrit

次から次へと生まれ続ける「新技術」。今日もまた、「今後ITエンジニアに必要な技術はこれだ!」「キャッチアップ!」「未来が変わる!」と迫ってきます。たしかに、どれもすごい技術なのには間違いなさそうです。

しかしながら、疲弊してきているのも正直なところ本音ではないでしょうか。

…本当に身に付けたい能力ってなんだろう。

このタンポポの隣を通るたびに、いつもわたしは思うのです。

コンクリートの間のわずかな土で、光合成をしながら、冬の寒さにも負けずに、花を咲かせて生きる。

うん、わたしにもこの能力が欲しい!

こう書いてみると、そこらのラノベの主人公のチートな能力よりチートですね、タンポポ。

帰ってきたPyrit

というわけで、前々回山の上で瞑想して考えた「ぼくの考えた最強のベンチマーク・ルール」にしたがってガンガンPyritを書き換えていきましょう。この影響でPyritのリポジトリの構成が諸々変わっているので、使ってる人がいたら気をつけてください(居ないと思うけど)。

現状の問題を煮詰めると、次の2つでした:

  • 以下の処理が、すべてCPUを取り合っている:
    • ランダムなパスワードを生成する処理
    • パスワードをSHA1とかごにょごにょして変形するsolveの処理(GPUにもオフロードされるたぶん一番重い処理)
    • solveした結果とパケットダンプを突き合わせて、パスワードが正解かチェックするCrackerの処理(CPUのみ)
    • GPUのためにデータを詰めて渡す処理
  • CPUとGPUのsolverが、シングルスレッドで生成される「パスワード」の仕事を奪い合っている

というものでした。最初の問題のせいで、CPUやGPUの真の実力がわからず、2つめの問題のせいでマシン全体の実力がわかりません。GPUが1%ぐらいしか使われていなかった衝撃を思い出しておきましょう。

後者のほうはかなり厄介な問題な気がします。というのも、現在のパスワード生成処理は次の極めて簡素なコードで書かれており、小手先の最適化では立ち向かえそうにないからです:

            while time.time() - t < timeout:
                pws = ["barbarbar%s" % random.random() for i in xrange(bsize)]
                cp.enqueue('foo', pws)
                r = cp.dequeue(block=False)
                if r is not None:
                    perfcounter += len(r)
                self.tell("rRunning benchmark (%.1f PMKs/s)... %s" % 
                        (perfcounter.avg, cycler.next()), end=None)

これ1.5倍速くするのは…まぁ…ギリギリ可能かもしれませんが…100倍は無理じゃん?

これらを踏まえつつ、今日はCPUやGPUだけ単体で走らせて測定してみて、何がどうボトルネックになっているのか探っていきましょう。そうすれば、パスワードを今の何倍の速度で生成してやればいいのかのヒントも得られるはずです。

CPUに本気を出してもらう

ベンチマークで測っているのは、一番重い(であろう)solveの処理です。デフォルトではsolveの処理スレッドはCPUのコア数分立ち上がるのですが、CPUは他の処理もしないといけないので、CPUを占有してsolveに専念することができませんし、同じようにCPUで生成されているパスワードも取り合ってしまうでしょう。そこで、立ち上げるスレッドの数を減らします。

.pyrit/configにあるlimit_ncpusでその設定が出来るので、これを減らしていって何が起こるのか調べます:

~/.pyrit$ cat config 
default_storage = file://
limit_ncpus = 0
rpc_announce = true
rpc_announce_broadcast = false
rpc_knownclients = 
rpc_server = false
use_CUDA = false
use_OpenCL = false
workunit_size = 75000

オリジナルのバージョンはlimit_ncpus = 0を指定したときは「デフォルト」扱いになってCPUのコア数分(このマシンでは48個)立ち上がるようになっていたのですが、後でGPUだけで実行するときのために、0を指定したときは本当に0個スレッドを立ち上げる(つまりCPUを一切使わない)ように変更しました(ちなみに、CUDAもCPUも使わない場合、フリーズします)。

測ってみた結果がこちら(PMKというのがいわゆる「パスワード」です):

ncputotal PMK/saverage PMK/s
1989.66989.66
22015.591007.795
43811.58952.895
87264.9908.1125
1613800.76862.5475
3217994.36562.32375
4818600.49387.5102083

で、グラフにするとこんな感じ:

16コアくらいまではだいたいスケールする(16倍とはいかないが、13倍くらいにはなる)んですが、32コアで突然頭打ちになって、コアあたりのパフォーマンスががくっと下がり、だいたい1.8万PMK/secぐらいで張り付きます。というわけで、これが現在のパスワード生成処理の限界なのでは?という仮説が立ち上がります。

GPUに本気を出してもらう

PyritはCUDAをデフォルトでは「全部使う」ようになっているので、これも1台ずつだけ使えるようにしてみました。CPUも一切使わず、GPU1台だけでクラックしていただきます。

で、やってみた結果がこちら:

デバイス名PMK/sec
Tesla V100-PCIE-16GB#1236335
Tesla V100-PCIE-16GB#2232217.35
GeForce GTX 1080 Ti#1222854.14
GeForce GTX 1080 Ti#2222457.37
全部一緒に動かす256639.55
CPUも全部同時に動かす136911.19

桁が多くて目が痛いのですが、単体で走らせた時に22〜23万PMK/sec、全部一緒に動かした時には25万PMK/secぐらい出ております。

これはつまり、あのパスワード生成コードは25万PMK/secは生成できているという事ですね。

…あれ?

ついでにいうと、CPUもGPUも同時に動かすと13万PMK/secという事で半分の速度にまで遅くなってしまいます。CPUとGPUで同時に処理して加速するはずだったのに、足を引っ張りあうのか、減速しております。

…あれ?デフォルトではCPUなしで実行できないし、ひょっとして、これ、誰も気づいてないんでは…。

結局、システム全体のポテンシャルはどれぐらいあるのか?

前に使ったProfilerをもっかい動かして、GPUを何割使ってるのか計測します:

Computed 225057.65 PMKs/s total.
CUDA:
#1: 'CUDA-Device #1 'Tesla V100-PCIE-16GB'': 237348.4 PMKs/s (RTT 0.5)
==41597== Profiling application: venv/bin/python pyrit benchmark
==41597== Profiling result:
            Type  Time(%)      Time     Calls       Avg       Min       Max  Name
 GPU activities:   99.08%  26.9879s       116  232.65ms  27.573ms  270.97ms  cuda_pmk_kernel
                    0.53%  143.46ms       116  1.2367ms  2.1760us  2.9881ms  [CUDA memcpy HtoD]
                    0.40%  108.18ms       116  932.60us  2.3360us  2.7747ms  [CUDA memcpy DtoH]
      API calls:   96.48%  27.1509s       116  234.06ms  27.585ms  272.31ms  cuMemcpyDtoH
                    1.88%  530.46ms         1  530.46ms  530.46ms  530.46ms  cuCtxCreate
                    0.54%  153.24ms       116  1.3211ms  27.996us  3.5144ms  cuMemcpyHtoD
                    0.47%  133.17ms         1  133.17ms  133.17ms  133.17ms  cuCtxDestroy
...

ベンチマークの時間は45秒、そのうちの27秒ほどを使っており、なかなかいい感じです。皮算用をしてみると、理想的な状況下ではGPU一枚あたり30万PMK/secぐらい出てもおかしくないと思います。で、4台あるので全体では120万PMK/sec出てほしい。

CPUは全体で理想的には4.8万PMK/sec出そうでしたから、合計すると125万PMK/secぐらいが理論的な限界値ではないかと皮算用できます。

あくまで皮算用ですけど。

原理的にパスワードはどれくらいの速度で生成できるのか?

ボトルネックだと思っていたパスワード生成はどうやらボトルネックではなさそうです。それを確かめるために、Pyritのコードから抜き出して加工した次のソースで測定します:

import time
import random

cnt = 0
t = time.time()
timeout = 30
while time.time() - t < timeout:
  pws = ["barbarbar%s" % random.random() for i in xrange(50000)]
  cnt += len(pws)

total = time.time()-t
print(cnt / total)

結果は268万PMK/secと出ました。うん、ボトルネックじゃないよここ…「推測はよくない」の例がまた一個増えましたね。…orz

とはいえ、皮算用したボトルネックのだいたい倍ぐらいの速度でしかなく、Pythonとその周辺コードはシングルスレッドでしか走らないことを踏まえると、最終的にボトルネックとして立ちはだかる可能性はまだ十分残されているでしょうな。

今月のまとめと来月の予告

  • CPUだけで走らせると、16コアと32コアの間で速度が頭打ちになる
  • GPUだけで走らせると、CPUとGPUを同時に走らせるより速い
  • 理想的なシステム全体のパフォーマンスは125万PMK/sec程度だと思われるが、パスワードの生成速度はそれをはるかに上回っているのでボトルネックではなさそう

以上の観察を踏まえるに、どうやら問題は最初の問題、つまり「CPUを取り合っている」の問題に尽きそうです。

次回はPythonのプロファイラーを使って、どこで詰まっているのかを調べましょう。

よいお年を〜。たんぽぽには、人間の暦なんか、関係ないけどね。