hyperoptで整数値の探索

今回もhyperoptの話題です。
この前の記事では、2変数関数の最小値を探すとき、
それぞれの変数を、一定区間の一様分布の中から探索しました。
そのときは、hp.uniformを使いました。

ただ、場合によっては、変数は実数ではなく整数を取ることもあります。
多くの例があると思いますが、機械学習で言えばDNNの中間層のユニット数などがそうですね。

この場合、 hyperoptで探索する空間を定義するときは、hp.randint(label, upper)というのを使うことで実現できます。
ドキュメントはこちら。
2.1 Parameter Expressions
ただ、これよく見ると最小値の指定ができません。
(ドキュメントに書かれてないだけで実は動くのではと思い試しましたが、本当にできません)

そのような時にどうしようかというのが今回の記事です。
0からではなく、1以上の数で探したいケースや、負の数も扱いたい時などありますからね。

今回使う関数はこちら。
$$f(x, y) = (x-5)^2 + (y+7)^2$$
この関数の$-10 \leq x\leq 10, -10 \leq y\leq 10$ の範囲で最小値を取る点を探します。
正解は明らかに$(x, y) = (5, -7)$ です。

それを実行するコードがこちら。


from hyperopt import hp, fmin, tpe, Trials


# 最小値を探したい関数
def f(x, y):
    return (x-5)**2 + (y+7)**2


# 探索する関数を定義する
def objective(args):
    x = args["x"]
    y = args["y"]

    return f(x, y)


# 探索する空間を定義する
space = {
        'x': -10 + hp.randint('x_', 21),
        'y': -10 + hp.randint('y_', 21),
    }

trials = Trials()

best = fmin(
            fn=objective,
            space=space,
            algo=tpe.suggest,
            max_evals=100,
            trials=trials,
        )

print(best)

# 以下出力
100%|██████████| 100/100 [00:00<00:00, 427.25it/s, best loss: 2.0]
{'x_': 16, 'y_': 2}

ちょっと惜しいですが、 f(x, y)が 2になる点を見つけたようです。

ポイントはここ


space = {
        'x': -10 + hp.randint('x_', 21),
        'y': -10 + hp.randint('y_', 21),
    }

このように定義することで、x_y_ はそれぞれ 0〜20の値を取り、
-10することで、xyは-10から10の値を取ります。
21にはならないので注意。(trials.vals["x_"]の値を確認するとx_=21は探していないことがわかります。)

出力の、{'x_': 16, 'y_': 2}から、
$x=16-10=6, y=2-10=8$ を見つけていることがわかります。

ただ、本当はx,yの値が欲しいのに、x_,y_を戻されると少しだけ不便です。
このようなときは、space_evalという関数が使えるようです。
ドキュメントが見つからないのですが、このページのサンプルコードで使われているのを参考にやってみます。


from hyperopt import space_eval
print(space_eval(space, best))

# 出力
{'x': 6, 'y': -8}

これで必要だった(x, y)の値が取れました。
正解から1ずつずれているのがちょっと惜しいですね。

Data Driven Developer Meetup #5 に参加しました

先日、Data Driven Developer Meetup #5 に参加してきました。
#4にも参加しているので2回目です。

今回の会場はログリーさんのオフィス。職場から近かったので徒歩で行けて便利でした。

メインセッション

データドリブンな小売店舗経営を支えるプロダクト開発の裏側』 株式会社ABEJA 大田黒 (@xecus) さん

まずこの発表資料は ABEJA SIX 2019で使われたものだそうです。
ABEJA SIXは参加したかったのですが業務的に行けなかったので、ラッキーでした。

動画の顔認識や、店内の導線の可視化など技術力が凄まじく高いのはもちろんですが、
「リアル版のグーグルアナリティクス」といった非常にイメージをつかみやすい説明が印象的でした。
これはビジネス的に強そうです。
上のスライド資料だと動画が静止画になってしまいっていましたが、
店舗内入り口のカメラの動画で、性別と年齢が推定されているデモなどは驚きでした。

プロジェクト進行における問題のところは自分も似たような問題に直面したのでよくわかります。
ゴール設定の難しさも、不確実性の多さも、本当にウォーターフォールとは相性が悪いです。

まだ自分の職場では、「研究開発」と言えるほどモデルの開発に専念できる体制も作れておらず、
機械学習のプロジェクトはもっとごちゃごちゃしているのですが、
チームで進める土壌づくりは力を入れないといけないと思いました。

機械学習で動くシステムを実運用に組み込むということ』 株式会社ブレインパッド アナリティクスサービス部 吉田勇太 (@yutatatatata ) さん

これもあるあるのオンパレードで苦笑いする場面の多い発表でした。
画像とテキストという違いはありますが、目視でやっていたことを機械学習で行いたいというニーズは自分も経験しています。
その対象が不均衡データで、
PrecisionとRecallがトレードオフになる厄介さや、Recallが重要だからと言って、Precisionはどこまでも犠牲にしていいわけじゃない苦労など
経験してきたので非常によくわかりました。

また、良いモデルを作るだけではなく業務フローをしっかり考えないといけないこと、
それにあたって、他のメンバーの協力を得ることの重要性など自分の失敗経験もあるので本当に納得です。

このあと、LTはお酒も入って気楽な雰囲気で聞いていました。

LT

学習行動データを蓄積するLearnig Record Storeの開発事例』by @yukinagae さん

まだまだ取り組みが始まったばかりという印象が大きかったのですが、
「Learnig Record Store」というコンセプトをしっかり立てて進められているのが素晴らしいと思いました。
GCPUG覚えておきます。

データドリブンを提供するサービスBrand Officialのアーキテクチャについて』 by @YasutakaYamamoto

久々に Airflowではなくdigdagの方が登場しました。
まだdigdagは使ったことがないのですが一度くらいは調べておきたいですね。
embulkは使っているのですが、よく理解していないので調べなければ。

今回も非常に参考になること、励みになることが多く参加してよかったです。
運営の皆様、登壇者の皆様、ありがとうございました。

hyperoptで多変数関数の最適化

前回前々回とhyperoptを扱ってきましたが、例があまりにもシンプルだったのでもう少しまともなサンプルを出します。
今回はあらかじめ定義された多変数関数が最大値をとる点を探します。
たった2変数の関数ですが、これ以上変数を増やしても本質的にはやることは変わりません。

きちんと最適解の近くを探していることがわかるように最後に可視化も行います。

例にする関数はこちらです。
$$2*e^{-(x-1)^2-(y-2)^2}+e^{-(x+3)^2-(y+4)^2}$$
点$(1, 2)$からちょっとずれたところで最大値を取りそうな関数ですね。
そして、 点$(-3, -4)$付近にもう一つ山があります。

早速やってみましょう。


import numpy as np
from hyperopt import hp,  fmin, tpe, Trials


# 最大値を探したい関数
def f(x, y):
    return 2 * np.exp(-(x-1)**2 - (y-2)**2) + np.exp(-(x+3)**2 - (y+4)**2)


# 探索する関数を定義する
def objective(args):
    x = args["x"]
    y = args["y"]

    # 最大値を検索するため、マイナスをつけて符号を反転
    return - f(x, y)


# 探索する空間を定義する
space = {
        'x': hp.uniform('x', -5, 5),
        'y': hp.uniform('y', -5, 5)
    }

trials = Trials()

best = fmin(
            fn=objective,
            space=space,
            algo=tpe.suggest,
            max_evals=100,
            trials=trials,
        )

print(best)

# 以下出力
100%|██████████| 100/100 [00:00<00:00, 293.11it/s, best loss: -1.9782848207515258]
{'x': 1.0647780122167085, 'y': 1.9180196819617374}

きちんとそれっぽい点を見つけてきましたね。
これまでの例と違うのは、最小値を探索する関数と、探索する空間が、fminの外側で事前に定義されている点です。

あとは trials オブジェクトの中身を確認し、正解付近を重点的に探索していることを確認しましょう。


import matplotlib.pyplot as plt
X = np.linspace(-5, 5, 100)
Y = np.linspace(-5, 5, 100)
xx, yy = np.meshgrid(X, Y)

fig = plt.figure(figsize=(8, 8))
ax = fig.add_subplot(1, 1, 1)
ax.contour(xx, yy, f(xx, yy))
ax.scatter(trials.vals["x"], trials.vals["y"])
plt.show()

できた図がこちらです。

$(-3,-4)$付近はそこそこに、$(1, 2)$付近を重点的に探していることがわかります。

hyperoptで探索したパラメーターの履歴を確認する

前回の記事に続いて、hyperoptの記事です。
前回のサンプルコードでは、100回の探索を経て見つけた最良の値を取得しましたが、
そのほかの結果もみたいということはあると思います。
機械学習のモデルであれば、どのパラメーターが精度に寄与していたとか、
どれはあまり重要でないとか確認したいですからね。

そのような時、hyperopt では Trials というオブジェクトを使います。
ドキュメントのこちらを参考にやってみましょう。
1.3 The Trials Object


from hyperopt import fmin, tpe, hp, Trials
trials = Trials()
best = fmin(
            fn=lambda x: x ** 2,
            space=hp.uniform('x', -10, 10),
            algo=tpe.suggest,
            max_evals=100,
            trials=trials,
        )
print(best)

# 出力
100%|██████████| 100/100 [00:00<00:00, 554.94it/s, best loss: 8.318492974446729e-09]
{'x': 9.120577270352315e-05}

Trials() ってのが途中で出てきた以外は前回の記事と同じコードですね。
この時、trials変数に100回分の履歴が残っています。
試しに一つ見てみましょう。


# trials.trialsに配列型で格納されている結果の一つ目
trials.trials[0]

# 内容
{'state': 2,
 'tid': 0,
 'spec': None,
 'result': {'loss': 96.69903496469523, 'status': 'ok'},
 'misc': {'tid': 0,
  'cmd': ('domain_attachment', 'FMinIter_Domain'),
  'workdir': None,
  'idxs': {'x': [0]},
  'vals': {'x': [-9.833566746847007]}},
 'exp_key': None,
 'owner': None,
 'version': 0,
 'book_time': datetime.datetime(2019, 3, 10, 15, 4, 45, 779000),
 'refresh_time': datetime.datetime(2019, 3, 10, 15, 4, 45, 779000)}

$x=-9.833566746847007$ を試したら 損失関数$x^2$の値は$96.69903496469523$ だったということがわかります。

このほか、次のような変数にそれぞれの値の辞書も格納されています。
trials.results
trials.vals

せっかくなので、本当に正解の$x=0$付近を探索していたのか可視化してみましょう。


import matplotlib.pyplot as plt
fig = plt.figure(figsize=(8, 8))
ax = fig.add_subplot(1, 1, 1)
ax.scatter(trials.vals["x"], trials.losses())
plt.show()

結果がこちら。

きちんと正解付近に点が密集していますね。

hyperoptのインストールと最も簡単な例

昨年参加したpycon2018
実践・競馬データサイエンス というセッションで知った、hyperopt というライブラリのメモです。
とりあえずこの記事では最もシンプルな例で動かすところまで行きます。
ちなみに、pyconの時の資料はこちら。
実践・競馬データサイエンス

なんでも、ハイパーパラメーターをグリッドサーチなどより効率的に探索してくれるとのこと。

公式サイトの通り、 pipでインストールできます。


pip install hyperopt

早速動かしてみましょう。
Basic tutorial

$x^2$ が最小になる$x$(正解は$x=0$です)を探索するコードです。


from hyperopt import fmin, tpe, hp
best = fmin(
            fn=lambda x: x ** 2,
            space=hp.uniform('x', -10, 10),
            algo=tpe.suggest,
            max_evals=100
        )
print(best)

# 以下出力
100%|██████████| 100/100 [00:00<00:00, 570.89it/s, best loss: 0.0009715421527694116]
{'x': 0.031169570942979175}

これで、 $-10<=x<=10$ の範囲を100回探索し、$x=0.031$あたりで最小と見つけてきたようです。 単純すぎていまいちありがたみが感じられませんね。 ただグリッドサーチと比較すると、仮に-10から10を単純に100分割したら区間の幅は0.2になるわけで、 それよりはずっと小さい誤差で最小値付近の値を見つけています。 (とはいえ、グリッドサーチなら $x=0$をピタリと見つけるとは思いますが。) 機械学習のハイパーパラメーターの探索等で役に立つと聞いて調べたものであり、 僕もそのように使うことが多いので、今後の記事でもう少し細かく紹介しようと思います。 一点気になるのは公式ドキュメントのこの記述。

Hyperopt has been designed to accommodate Bayesian optimization algorithms based on Gaussian processes and regression trees, but these are not currently implemented.

[翻訳]
Hyperoptは、ガウス過程と回帰木に基づくベイズ最適化アルゴリズムに対応するように設計されていますが、これらは現在実装されていません。

Qiitaなどで、hpyeroptをベイズ最適化のライブラリだって紹介している人を見かけたこともあるのですが、そうではなさそうです。

DataFrameを特定の列の値によって分割する

python(pandas)で巨大なデータフレームを扱っている時、
ある列の値によって、分割したいことが良くあります。

pandasのgroupbyのユーザーガイドを見ていて、良い方法を見つけたのでそれを紹介します。
Iterating through groups
のところに書かれていますが、groupbyした後そのままsumなどの関数を適用するだけでなく、イテレーターとして使うことで、
各グループに対して順番に処理を行えます。

実際にやってみましょう。
とりあえずサンプルのデータフレームを作ります。


import pandas as pd
import numpy as np
df = pd.DataFrame(
            np.random.randn(10000, 10),
            columns=['x'+str(i) for i in range(10)],
        )
df['y'] = np.random.randint(20, size=10000)

これで、長さが1万行、列が11個あるデータフレーム(df)ができました。
これをy列の値(0〜19)によって、20個のデータフレームに分割します。


df_dict = {}
for name, group in df.groupby('y'):
    df_dict[name] = group

これで、20個のデータフレームに分割できました。

ちなみに、この方法を知る前は次のように書いていました。


df_dict = {}
for name in set(df.y):
    df_dict[name] = df[df.y == name]

コードの量はほとんど同じなのですが、この方法は無駄な比較を何度も行うことになります。

実際、実行時刻を測ってみるとgroupbyを使う方法は、
CPU times: user 7.47 ms, sys: 3.28 ms, total: 10.8 ms
Wall time: 7.71 ms

なのに対して、使わない方法は
CPU times: user 16.4 ms, sys: 1.81 ms, total: 18.2 ms
Wall time: 16.5 ms
となり、倍以上の時間を要しています。

Tableau Publicを使ってみる

職場ではTableau Online と Tableau Desktopを使っています。
非常に便利なツールなので家でも使いたいのですが、これが結構お高い
一方で実はTableauには無料版のTableau Public というのがあります。
前々から存在は知っていたのですが、アカウントを作ってサインインしないと使えないのかなと思い込んでいて、面倒なので使っていませんでした。

ただ、いつまでも先延ばしにするのも変なので、ダウンロードページからダウンロードしてみました。
https://public.tableau.com/ja-jp/s/download
するとなんと、最初にメールアドレスを入力する必要はありますが、
それだけでソフトをダウンロードし、普通にインストールしたら使えました。
めっちゃ手軽。
(逆にいうと、Tableau Publicにパブリッシュするためのアカウントは別途作成する必要があります)

アカウントを作らないとデータを保存できない不便はあるのですが、
私生活におけるデータ分析なら全然問題なく使えそうです。

参考に自分のこれまでの株取引のデータから、
ポジションの保有期間と平均損益(および、保有期間別の取引数)をヒストグラムにしてみました。

損はさっさと切って、利益は伸ばすトレードがちゃんと実践できてるように見えます。
(それならなぜ最近は負けっぱなしなのか、というのは別の話。)

完全にプライベートな話なので公開しませんが、特定の曜日や月に弱い傾向などもサクサクと可視化して見つけられたので、
自分の取引履歴だけでもかなり楽しめそうです。

アカウントの方はそのうち作って、何か公開できるコンテンツができたらあげてみましょうかね。

matplotlibで折れ線グラフ(plot)を書く時に線の書式を指定する

matplotlibで折れ線グラフを書く時に線の書式を指定する方法のメモです。
大抵のケースでは色で見分ければ十分なのですが、たまに点線や破線を使いたいことがあります。
(色が使えない場合もマーカーで見分けてもいいのですけどね。)

ドキュメントはこちらです。
set_linestyle(ls)
matplotlib.pyplot.plot の linestyle のリンクをクリックすると上のページに飛びます。

方法は簡単で、 plot する時に linestyle変数に ‘dashed’ などの文字列や ‘:’ などの記号を指定するだけです。
Noneで、線を消すこともできます。(その場合はマーカーを使いましょう)

試しに4種類の線を引いてみました。
(線の種類の区別を目立たさせるため、色は揃えました。)


import matplotlib.pyplot as plt
import numpy as np

# データ生成
x = np.arange(10)
y0 = np.random.normal(0, 0.5, 10)
y1 = np.random.normal(1, 0.5, 10)
y2 = np.random.normal(2, 0.5, 10)
y3 = np.random.normal(3, 0.5, 10)

plt.rcParams["font.size"] = 14
fig = plt.figure(figsize=(10, 6))
ax = fig.add_subplot(1, 1, 1, title="linestyleの例")
ax.plot(x, y0, linestyle="-", c="b", label="linestyle : solid, '-'")
ax.plot(x, y1, linestyle="--", c="b", label="linestyle : dashed, '--'")
ax.plot(x, y2, linestyle="-.", c="b", label="linestyle : dashdot, '-.''")
ax.plot(x, y3, linestyle=":", c="b", label="linestyle : dotted, ':''")
ax.legend(bbox_to_anchor=(1.05, 1))
plt.show()

出力はこちら。

pythonのdictから値を取得する時にデフォルト値を設定する

今は当たり前のように使っているテクニックですが、初めて知った時は感動したので紹介。

ドキュメントはここ

pythonの辞書オブジェクトから値を取得する時、[ ]を使う方法と、
get()を使う方法があります。


sample_dict = {
    'apple': 'リンゴ',
    'orange': 'オレンジ'
}

print(sample_dict['apple'])  # 'リンゴ'
print(sample_dict.get('orange'))  # 'オレンジ'

タイプ数が少ないので[ ]の方を好むかたが多いと思うのですが、
getの方には、引数がkeyの中に存在しなかった時の初期値を指定できるメリットがあります。

例えば、
sample_dict["grape"] は keyErrorが発生してしまいます。
しかし、
sample_dict.get("grape")はNoneを返しますし、
sample_dict.get("grape", "未定義")とすると、
2個目の引数である、”未定義”を返してくれます。

エラー処理もいらないし、事前に辞書に取得しようとしているkeyが含まれているか確認する必要もなくなり、
非常に便利です。

TreasureDataのTD_PARSE_AGENT関数が便利

前の記事が、若干TreasureData(正確にはPresto)への文句っぽくなったので、
今回はTreasureDataでよく使っている便利関数を一つ紹介します。

それが、TD_PARSE_AGENTです。
ドキュメント: 12. TD_PARSE_AGENT
(これはTreasureDataのUDFなので、他の環境のPrestoでは利用できません)

これは、アクセスログのuser-agentの文字列を解析する関数で、
アクセスデータを分析する時に、OSやブラウザ、PC/スマートフォンの区別、クローラーの除外などの目的で毎日のように使っています。

今、僕がこの記事を書いているブラウザのユーザーagentで試してみましょう。
Mozilla/5.0 (Macintosh; Intel Mac OS X 10_14_3) AppleWebKit/537.36 (KHTML, like Gecko) Chrome/72.0.3626.109 Safari/537.36


SELECT
    TD_PARSE_AGENT('Mozilla/5.0 (Macintosh; Intel Mac OS X 10_14_3) AppleWebKit/537.36 (KHTML, like Gecko) Chrome/72.0.3626.109 Safari/537.36')
-- 出力 
-- {"os":"Mac OSX","vendor":"Google","os_version":"10.14.3","name":"Chrome","category":"pc","version":"72.0.3626.109"}

欲しい情報が一通り揃っていますね。
ちなみに、戻り値はStringではなくMap型です。
そのため、この中から osだけ取りたい時は、下記のようにすると取得できます。
(access_logテーブルのagentという列に、例のUseragentの値が入っていると仮定します)


SELECT
    TD_PARSE_AGENT(agent)['os']
FROM
    access_log
-- 出力 
-- Mac OSX

また、ドキュメントの Exampleの中に、書いてある通り、
category はUserAgentによって7種類の値を取ります。


SELECT TD_PARSE_AGENT(agent)['category'] AS category FROM www_access
> pc // => "pc", "smartphone", "mobilephone", "appliance", "crawler", "misc", "unknown"

clawlerを除外するために使うことが多いです。