pandasで数値データを区間ごとに区切って数える

今日偶然知って便利だった関数を紹介します。

まずやりたいこと。
配列でもDatafarameの列でも、数値データを区間(bins)ごとに区切って数えたい場面は多々あると思います。
例えば、
[455,133,666,111,645,236]
みたいなデータを、
100以上,200未満が2個、200以上,300未満が1個、、、といった具合です。
ヒストグラム書くためにやる操作ですね。

これまで僕がこれをやる時100単位でやるのであれば、
各データを100で割り算して小数点以下を切り捨てて、再度100倍するといった面倒なことをやってました。

これを pandas.cut という関数を使うと、非常に容易に行えます。
ドキュメントはこちら。
pandas.cut

早速これを使ってみましょう。
次のコードは、乱数で生成した10個の数値を、100ごとに区切って数えるものです。
結果に出てくる “[300, 400)”といった出力は、文字列に見えますが、Intervalというオブジェクトです。
right=False は指定しないと、n以上m未満ではなく、nより大きいm以下、って区切り方になります。
引数のbinsには区切りたい点のリストを渡していますが、arrange関数で生成して渡しています。


import pandas as pd
import numpy as np  # ダミーデータ生成用

# 300〜999 の整数を10個生成
data = np.random.randint(300, 1000, 10)
# 100 ごとに区切る。
category = pd.cut(data, bins=np.arange(300, 1100, 100), right=False)

print("データとそれが所属する区間")
for d, c in zip(data, category):
    print(d, "・・・", c)

print("\n区間ごとのデータ件数")
print(category.value_counts())

# 以下出力
データとそれが所属する区間
309 ・・・ [300, 400)
305 ・・・ [300, 400)
874 ・・・ [800, 900)
953 ・・・ [900, 1000)
727 ・・・ [700, 800)
950 ・・・ [900, 1000)
384 ・・・ [300, 400)
789 ・・・ [700, 800)
486 ・・・ [400, 500)
501 ・・・ [500, 600)

区間ごとのデータ件数
[300, 400)     3
[400, 500)     1
[500, 600)     1
[600, 700)     0
[700, 800)     2
[800, 900)     1
[900, 1000)    2
dtype: int64

pythonでアルファベットの大文字小文字変換

自然言語処理系の機械学習を行うときやテキストマイニングをやるとき、
前処理としてアルファベットの大文字小文字を統一することがよくあります。
(失われる情報もあるのでしない方が良いという主張も見たことがありますが、僕は大抵の場合小文字に統一します。)

これはpythonの文字列が持っている
str.lower()str.upper()
を使うことで実現できます。
pythonの組み込み型のドキュメントを見ると乗っています。


text = "Hello World!"
print(text.lower())  # hello world!
print(text.upper())  # HELLO WORLD!

これだけだと、記事にしなかったと思うのですが、ドキュメントを読んでいると他にも関数が準備されていることがわかりました。


text = "Hello World!"
# 大文字と小文字を入れ替える
print(text.swapcase())  # hELLO wORLD!

text = "HELLO world!"
# 各単語の1文字目を大文字に、残りを小文字に変換する
print(text.title())  # Hello World!

text = "HELLO world!"
# 最初の文字を大文字に、残りを小文字に変換する
print(text.capitalize())  # Hello world!

正直、使う場面を思いつかないのですが面白いですね。
英語ネイティブな人たちが使うのでしょうか。

pythonで文字と文字コードの相互変換

pythonで文字を文字コードに変えたり、文字コード(整数)をそれが表す文字に変換したりする方法のメモです。
これらはそれぞれ、 ord と chr という組み込み関数で実現できます。

ドキュメントはこちら。
組み込み関数

これらをアスキーコードへの変換やアスキーコードから文字への変換だと説明しているサイトもあるようですが、
ドキュメントに書かれている通りUnicodeに対応しています。
(もしかしたらpython2系の時代はアスキー文字だけだったのかな?)

組み込み関数なので何もインポートせずに利用可能です。
ただし、ordは文字のみを受けつけ、文字列を渡すとエラーになるので注意してください。

単純なサンプル。


print(ord("a"))  # 97
print(ord("あ"))   # 12354
print(chr(97))   # a
print(chr(12354))   # あ

numpyの線形代数モジュールで連立一次方程式を解く

numpyで作成できる配列によく似たオブジェクトのarrayですが、配列同様にネストして多次元のarrayを作成できます。
(日頃から意識せず普通に使ってますね。)
特に2次元のarrayについては、行列として扱うことが可能です。

そしてnumpyには、numpy.linalgという線形代数のモジュールがあり、
行列式や固有値、逆行列の計算などができます。

ドキュメントはこちら。
Linear algebra (numpy.linalg)

今回は基本的な例として次の連立一次方程式を解いてみましょう。
$$
\begin{eqnarray}
&x_0-4&x_1+2&x_2&=7\\
9&x_0-5&x_1+2&x_2&=-2\\
3&x_0-10&x_1+5&x_2&=3
\end{eqnarray}
$$

これは
$$
A = \left[ \begin{matrix}
1&-4&2\\
9&-5&2\\
3&-10&5
\end{matrix}\right]
$$
$$
b = \left[ \begin{matrix}
7\\
-2\\
-3
\end{matrix}\right]
$$
とおいて、 Aに逆行列が存在すれば、(つまりAの行列式が0でなければ)$A^{-1}b$を計算することで解けます。

それではやってみましょう。
行列式はlinalg.det(A)
逆行列はlinalg.inv(A)で算出できます。
行列の式はnp.dot(A,B)A.dot(B)か、A@Bなどの書き方があります。


# 行列の定義
A = np.array(
        [
            [1,  -4, 2],
            [9,  -5, 2],
            [3, -10, 5]
        ]
    )
b = np.array([7, -2, 3]).T

# Aの行列式の確認 (実装の都合でわずかに誤差が出ます。)
print("dat(A)=", np.linalg.det(A))

# 方程式の解
print("[x_0, x_1, x_2]=", np.linalg.inv(A)@b)

# 以下出力
dat(A)= 0.9999999999999893
[x_0, x_1, x_2]= [ -29. -223. -428.]

これで、連立一次方程式が解けました。
Aの行列式が微妙に1にならないのは実装されていているアルゴリズムの都合のようです。

Sympyで方程式を解く

一つ前の自己回帰過程のサンプルを算出する記事のコードの中で、
こっそりとSympyというのを使っていました。

これはpythonで数式処理を行うためのライブラリです。
簡単な使い方などは公式のチュートリアルのイントロを見ていただくとわかります。
SymPy Tutorial » Introduction

これで記事を終えてしまうとあんまりなので、今回は方程式を解く方法紹介します。
シンプルに
$$x^2-5x+6=0$$
を$x$について解いて $x=2,3$という解を得るコードです。


import sympy
# 記号として使いたい文字をsymbols関数で定義
x = sympy.symbols('x')
# 多項式の定義
f = x**2 - 5*x + 6
print(f)  # 出力 : x**2 - 5*x + 6
# solveに渡すことで、 f = 0 という方程式を xについて解いてくれます。
sympy.solve(f, x)  # 出力 [2, 3]

多項式の展開や因数分解、微積分など、結構色々なことができるライブラリなので便利です。
このブログでも順次各機能を紹介したいと思います。

The Zen of Python

Pythonプログラマーならどこかでみたことがあると思われる「The Zen of Python」が、
pep 20 で定義されたものであることを知りました。

PEP 20 — The Zen of Python

ちなみに、読みたくなったら import thisでいつでも読むことができます。
thisって何だろうとずっと思っていたのですが、Easter Eggとして作られたものだったようです。


>>> import this
The Zen of Python, by Tim Peters

Beautiful is better than ugly.
Explicit is better than implicit.
Simple is better than complex.
Complex is better than complicated.
Flat is better than nested.
Sparse is better than dense.
Readability counts.
Special cases aren't special enough to break the rules.
Although practicality beats purity.
Errors should never pass silently.
Unless explicitly silenced.
In the face of ambiguity, refuse the temptation to guess.
There should be one-- and preferably only one --obvious way to do it.
Although that way may not be obvious at first unless you're Dutch.
Now is better than never.
Although never is often better than *right* now.
If the implementation is hard to explain, it's a bad idea.
If the implementation is easy to explain, it may be a good idea.
Namespaces are one honking great idea -- let's do more of those!

とても良いことが書かれていると思いますので、
定期的に読み返したいですね。

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ずつずれているのがちょっと惜しいですね。

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をベイズ最適化のライブラリだって紹介している人を見かけたこともあるのですが、そうではなさそうです。