pandasでgroupbyした時に複数の集計関数を同時に適用する

前の記事の続きです。
pandasでデータフレームをgroupbyした時に使える集計関数
ドキュメントのこの記事で参照した部分のすぐ下に、
Applying multiple functions at once
という段落があります。
実はこれ初めて知りました。
今までグルプごとに個数と、平均と、標準偏差を計算したい、みたいな時は、
groupbyして集計を個別に実施して、その結果をmergeするという非常に面倒なことをずっとやっていました。

それが、aggというのを使うと一発でできるようです。


import pandas as pd
from sklearn.datasets import load_iris

# データフレームの準備
iris = load_iris()
df = pd.DataFrame(iris.data, columns=iris.feature_names)
df["target"] = iris.target
df["target_name"] = df.target.apply(lambda x:iris.target_names[x])
del df["target"]

df.groupby("target_name").agg(["count", "mean", "std"])

出力されるのが次です。(ブログのレイアウトの都合上画像で貼り付けます。)

これは便利です。
また、DataFrameのカラム名が2段になっています。
これをみて、indexだけではなく実はcolumnsでも、MultiIndexが使えることを知りました。

pandasでデータフレームをgroupbyした時に使える集計関数

データの集計や分析をpandasで行う時、平均や合計を求めるために、
groupbyを使って集計することがよくあると思います。

非常に手軽に使え流のでなんとなく .sum()や .mean()と書いていたのですが、
そういえば他にどんな関数が使えるのか調べたことがなかったと思ったのでドキュメントを見てみました。
まずここ。
pandas.DataFrame.groupby
平均をとるサンプルコードがありますが求めていた関数の一覧がないですね。

よく読むと、See the user guide for more.とあります。
そのuser guideがこちらです。

Group By: split-apply-combine

この下の方に一覧がありました。

Function

Description

mean()

Compute mean of groups

sum()

Compute sum of group values

size()

Compute group sizes

count()

Compute count of group

std()

Standard deviation of groups

var()

Compute variance of groups

sem()

Standard error of the mean of groups

describe()

Generates descriptive statistics

first()

Compute first of group values

last()

Compute last of group values

nth()

Take nth value, or a subset if n is a list

min()

Compute min of group values

max()

Compute max of group values

グループ化した後に、describe()なんてできたんですね。
少し試してみたのですがこれ便利そうです。
他にもSeriesをスカラーに変換するlambda式なども使えるようです。

kerasのto_categoricalを使ってみる

機械学習の特徴量や正解ラベルをone-hotベクトルにするとき、
自分で実装するか、sklearnのOneHotEncoderを使うことが多いです。
稀に、pandasのget_dummiesを使うこともあります。

ところが、kerasのサンプルコードを読んでいると、to_categoricalというのもよく使われているので確認してみました。
軽く動かしてみると思った通りの挙動をしたので特に必要というわけでもないのですが、
使うライブラリをkerasに統一したいことがあれば利用するかもしれません。

とりあえず使ってみましょう。


import numpy as np
from keras.utils import to_categorical
# テスト用のデータ生成
data = np.random.randint(low=0, high=5, size=10)
print(data)
# One-Hotベクトルに変換
print(to_categorical(data))

# 以下出力
[3 2 1 4 4 1 0 1 0 2]

[[0. 0. 0. 1. 0.]
 [0. 0. 1. 0. 0.]
 [0. 1. 0. 0. 0.]
 [0. 0. 0. 0. 1.]
 [0. 0. 0. 0. 1.]
 [0. 1. 0. 0. 0.]
 [1. 0. 0. 0. 0.]
 [0. 1. 0. 0. 0.]
 [1. 0. 0. 0. 0.]
 [0. 0. 1. 0. 0.]]

to_categorical の2番目の引数(num_classes)として、数値を渡すと、
データの最大値を指定できます。
ただし、データの最大値+1より小さい値を渡すとエラーです。
極端な不均衡データを扱うときなどは念のため指定しておいたほうが安全かも。

試してみたサンプルです。


print(to_categorical(data, 8))

# 出力
[[0. 0. 0. 1. 0. 0. 0. 0.]
 [0. 0. 1. 0. 0. 0. 0. 0.]
 [0. 1. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 1. 0. 0. 0.]
 [0. 0. 0. 0. 1. 0. 0. 0.]
 [0. 1. 0. 0. 0. 0. 0. 0.]
 [1. 0. 0. 0. 0. 0. 0. 0.]
 [0. 1. 0. 0. 0. 0. 0. 0.]
 [1. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 1. 0. 0. 0. 0. 0.]]

Data Pipeline Casual Talk に参加しました

2018年2月13日に、エムスリーさんのオフィスで開催された、Data Pipeline Casual Talkに参加してきました。
実はもともと抽選に漏れていて、補欠だったのですが、開催2時間前くらいに急に補欠から繰り上がり参加になったので慌てて会場に行きました。
一瞬、急すぎるからキャンセルしようかとも思ったのですが、結果的に参加できて非常に幸運でした。

発表一覧

以下感想です。例によって僕の主観が多々入ります。

AI・機械学習チームにおけるデータパイプライン構築

機械学習チームを立ち上げられた時に困ったという各課題については、
自分もまさに経験してきたもので本当に参考になりました。
ログ出力が適切でなかったり、classやtaskの設計が悪かったり、
モデルの再現性の問題やテストの効率の問題など。

Luigi というのを拡張して対応されているということで、
技術力の高さを感じました。
tensorflowやgensimのモデルを同じインタフェースでloadできる仕組みは便利そうですね。
物によって保存や読み込みの方法が違うのでいつも地味に不便な思いをしています。

自分のところではまだそもそも機械学習基盤と呼べるようなものを作れていないので、
luigiも検討対象に加えたいです。

丘サーファーへ「水」を届けるために-これまでとこれから-

発表を聞いていて、金融SEをやっていた時のことを思い出しました。
丘サーファーという表現も面白い。

今の僕の職場では個人情報マスク済みのデータに比較的自由にアクセスでき、
データへのアクセスという面では問題なく業務を進められていますので、やはり恵まれた環境なのでしょう。

ただ、Cloud Composer を活用されているのは参考にしたいです。
生Airflowを使っていて、保守に手が回っていないので。

データ基盤の3分類と進化的データモデリング

論理設計(データモデル)と物理設計(システム構成)を分けて考えられているのが参考になりました。
データパイプラインを設計する時に両端を先に考えて挟み込むように真ん中へ進むのも納得です。
確かに普段の業務でもうまく進んでいるときはこの順番で考えています。

担当者のロカールPCにあるExcelシートが実はデータ基盤の役割を果たしているかもしれない、と聞いて、
即座に具体的なエクセルファイルが思い当たり苦笑いしてしまいました。
各アンチパターンにも思い当たる節が多々あり、今後改善していきたいです。
データ基盤の要素を技術要素と対応させて分けるのもアンチパターンだというのも覚えておこう。

データ分析基盤を「育てる」ための技術

分析作業の主なフローのスライドでまさに自分たち直面している問題が取り上げられていて笑いました。
いろんなところからの依頼が増えてくるとSQLを各作業がどんどん増えて
それで疲弊してしまうのですよね。
良い基盤を作れば解決するというものではなく、
データ基盤を育てていくという考えが大事。

リブセンスのデータ分析基盤とAirflow

Airflowを使ったデータ基盤を構築されています。
僕らの環境とよく似ているので、これも身に覚えがある苦労話に苦笑いする場面が多くありました。
バージョンアップの問題などもまさに。
社員が誰でもSQLをかけるというのは素直にすごいと思います。
特に営業の方たちにまでその文化を広げるのはきっと大変だったのではないかと。
ユーザー数の差があるのはもちろんですが、
それを考慮しても活用具合でずいぶん遅れをとっている気がするので負けないようにしたい。

まとめ

データパイプラインやデータ基盤はその重要性を日々感じているのですが、
専任の担当者もいなくてなかなか手が回らず、いろんな課題意識がある分野でした。
機械学習やデータウェアハウス単体の話に比べて他社の事例もすくなく、
自分たちだけこんなに苦労してるんじゃないかと不安になることもあったので今回のカジュアルトークに参加してよかったです。
だいたいどこも同じような課題に直面されていて、それぞれ工夫して対応されていることがわかりました。
自分が漠然とこんな風にしたいと思っていたことが明文化されていたスライドも多くハッとする場面も多々ありました。
あと、Airflow使ってる会社ってこんなに多かったんですね。
逆に、トレジャーデータは一度も登場しなかった。
今回だけでなく、今後も開催されるそうなので楽しみにしています。

kerasのMNISTデータを読み込んでみる

kerasにはscikit-learnと同じように、いつくかのサンプルデータが付属しています。
その中の一つがMNISTという28*28ピクセルの手書き数字文字のデータです。

scikit-learn にも digits という手書き数字のサンプルデータがありますが、こちらは8*8ピクセルの結構データ量の小さいデータです。
ライブラリの使い方を紹介する上ではこれで問題がないのでよく使っていますが、深層学習をつかずとも十分に判別できてしまうのが短所です。
(このブログの過去の記事でも使ってきたのはこちらです。)

せっかくディープラーニングを試すのであればこちらを使った方がいいと思うので、使い方を紹介します。

ドキュメントはここ。
MNIST database of handwritten digits

早速ですが読み込み方法です。


from keras.datasets import mnist
(x_train, y_train), (x_test, y_test) = mnist.load_data()

インポートしてloadすれば良いというのは scikit-learnのサンプルデータと同じですが、
load_dataの戻り値が少し特徴的です。
長さ2のタプル2個を値に持つタプルを返してきます。
これを受け取るために、上記のような書き方をします。
変数に格納された配列の次元数を見ると次のようになります。


print(x_train.shape) # (60000, 28, 28)
print(y_train.shape) # (60000,)
print(x_test.shape) # (10000, 28, 28)
print(y_test.shape) # (10000,)

x_train などの中身を見れば、ピクセルごとの濃淡として画像データが入っていることがわかるのですが、
せっかくの画像データなので、画像として可視化してみましょう。
matplotlibのimshow という関数を使うと便利です。

matplotlib.pyplot.imshow

最初の16データを可視化したのが次のコードです。
デフォルトの配色はイマイチだったので、文字らしく白背景に黒文字にしました。


import matplotlib.pyplot as plt
fig = plt.figure(figsize=(12, 12))
for i in range(16):
    ax = fig.add_subplot(4, 4, i+1, title=str(y_train[i]))
    ax.imshow(x_train[i], cmap="gray_r")
plt.show()

結果はこちら。
digitsに比べて読みやすい数字画像データが入ってますね。

パイプラインのグリッドサーチ

scilit-learnでグリッドサーチする方法買いた以前の記事の、
次の記事として、パイプラインを使ったモデルのグリッドサーチ方法を書く予定だったのを失念していたので紹介します。
(パイプラインそのもの紹介もまだなのでそのうち記事に起こします。)

公式のドキュメントはこちらです。
sklearn.pipeline.Pipeline
sklearn.model_selection.GridSearchCV

単一のモデルのグリッドサーチとの違いは、サーチ対象のパラメータを指定する時、
“変数名”: [変数のリスト]
で指定したところを、
“ステップ名”__”変数名”: [変数のリスト]
と指定するようにするだけです。

例を見た方がわかりやすいので、irisのデータと、
モデルは PCA + ロジスティック回帰 でやってみます。


# 必要なライブラリのインポート
from sklearn.datasets import load_iris
from sklearn.model_selection import train_test_split
from sklearn.model_selection import GridSearchCV
from sklearn.pipeline import Pipeline
from sklearn.decomposition import PCA
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import classification_report

# データの準備
iris = load_iris()
X_train, X_test, y_train, y_test = train_test_split(
            iris.data,
            iris.target,
            test_size=0.2,
            stratify=iris.target,
)

# モデル(パイプライン)の作成
clf = Pipeline(
    [
        ("pca", PCA()),
        ("lr", LogisticRegression())
    ]
)

# 探索するパラメータの設定
# ここのkeyの指定方法が重要
params = {
    "pca__n_components": [2, 3, 4],
    "lr__penalty": ["l1", "l2"],
    "lr__C": [0.01, 0.1, 1, 10, 100]
}

# グリッドサーチ
gs_clf = GridSearchCV(
    clf,
    params,
    cv=5
)
gs_clf.fit(X_train, y_train)

# 最適なパラメーター
print(gs_clf.best_params_)

# 最適なモデルの評価
best_clf = gs_clf.best_estimator_
print(classification_report(y_test, best_clf.predict(X_test)))

# 以下出力結果

{'lr__C': 100, 'lr__penalty': 'l2', 'pca__n_components': 3}
             precision    recall  f1-score   support

          0       1.00      1.00      1.00        10
          1       1.00      1.00      1.00        10
          2       1.00      1.00      1.00        10

avg / total       1.00      1.00      1.00        30

最終的にテストデータでの評価が正解率100%なのは運が良かっただけです。
train_test_splitの結果次第で、実行するたびにbest_params_も評価も変わります。

次元削減とその後の分類機は、まとめて最適かしたいので非常に便利です。

kerasで学習途中のモデルを保存する

kerasでモデルを学習している時、学習途中のモデルが欲しいことがよくあります。
デフォルトでepockごとに評価スコア出るので、テストデータに限っては過学習する前のモデルが欲しいですし、
前の記事のEarlyStoppingで、patienceを設定しているのであれば、
Stopした最後のバージョンより、そのいくつか前のものの方が評価が高いからです。
このような時、ModelCheckpointというCallBackを使うことで、Epockごとのモデルを保存しておくことができます。

前回同様データの準備やモデルの構築は省略して、ModelCheckpointに必要な部分のコードだけ紹介します。


from keras.callbacks import ModelCheckpoint

# checkpointの設定
checkpoint = ModelCheckpoint(
                    filepath="model-{epoch:02d}-{val_loss:.2f}.h5",
                    monitor='val_loss',
                    save_best_only=True,
                    period=1,
                )
# 学習
history = model.fit(
                    X_train,
                    y_train,
                    epochs=100,
                    batch_size=16,
                    validation_data=[X_test, y_test],
                    callbacks=[early_stopping, checkpoint]
            )

これで、filepathに指定したファイル名でepockごとのモデルが保存できます。
見ての通り、ファイル名にはepock数などを変数で埋め込むことができます。
また、save_best_only=Trueを設定することで、
前のモデルよりも精度が良い場合のみ保存することができます。
ファイル名に変数を含めないように指定しておけば、ベストのモデルだけが保存されている状態になります。

あとは、保存されているモデルをloadして利用すればOKです。

kerasでモデルの学習が進まなくなったら学習を止める方法

機械学習をやっているとき、過学習の抑止や時間の節約のためにモデルの改善が止まった時点で学習を止めたいことがあります。
kerasでは CallBack に EarlyStopping というオブジェクトを設定するおことでそれを実現できます。

モデル本体やデータについてのコードは省略しますので別記事を参照してください、該当部分だけ紹介します。


# インポート
from keras.callbacks import EarlyStopping

# EaelyStoppingの設定
early_stopping =  EarlyStopping(
                            monitor='val_loss',
                            min_delta=0.0,
                            patience=2,
)

# 学習
history = model.fit(
                    X_train,
                    y_train,
                    epochs=100,
                    batch_size=16,
                    validation_data=[X_test, y_test],
                    callbacks=[early_stopping] # CallBacksに設定
            )

これで、 monitor に設定した値が、 patienceの値の回数続けてmin_delta以上改善しないと、
学習がストップします。
monitor には ‘val_loss’ の他、 ‘val_acc’なども設定可能です。

patience の設定が0の時と1の時は挙動が全く同じに見えますね。

特にデメリットも感じられないので、kerasで機械学習を試す時はほぼ確実に使っています。
あまりにも学習が進まないうちに止まってしまう時は、EarlyStopping無しで試したりするのですが、
経験上、EarlyStoppingが悪いことは少なく、モデルの設計が悪かったり、その他どこかにミスがあることが多いです。

kerasで作成したモデルの保存と読み込み

kerasで作成したモデルを保存したり、次回使うときに読み込んだりする方法のメモです。

とりあえず、modelという変数に学習済みのモデルが格納されているとします。
(前回の記事のモデルです。)

 


model.summary()

# 以下出力
_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
dense_1 (Dense)              (None, 16)                48        
_________________________________________________________________
dense_2 (Dense)              (None, 1)                 17        
=================================================================
Total params: 65
Trainable params: 65
Non-trainable params: 0
_________________________________________________________________

pythonでオブジェクトを保存するときは、pickleを使うことが多いですが、
kerasではその方法は推奨しないと明記されています。

How can I save a Keras model?

It is not recommended to use pickle or cPickle to save a Keras model.

その代わりに、saveやload_modelというメソッドが用意されていて、
HDF5形式で保存や読み込みができます。


# 保存
model.save("model.h5")

# 読み込み
from keras.models import load_model
model_0 = load_model("model.h5")

kerasで簡単なモデルを作成してみる

前回の記事でkerasを使えるようにしたので、
動作確認を兼ねて非常に単純なモデルを作ってみました。

線形分離可能なサンプルではつまらないので、半径の異なる2円のToyDataでやってみましょう。
scikit-learnにmake_circlesという関数があるので、これを使います。

早速データを生成して訓練データとテストデータに分け、訓練データの方を可視化してみます。


# 必要なモジュールのインポート
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import make_circles
from sklearn.model_selection import train_test_split
from sklearn.metrics import classification_report
from keras.models import Sequential
from keras.layers.core import Dense

# データの準備
data, target = make_circles(n_samples=2000, noise=0.2, factor=0.3)
X_train, X_test, y_train, y_test = train_test_split(
                                            data,
                                            target,
                                            test_size=0.2,
                                            stratify=target
                                        )
# 訓練データの可視化
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(1, 1, 1, aspect='equal', xlim=(-2, 2), ylim=(-2, 2))
ax.scatter(X_train[y_train == 0, 0], X_train[y_train == 0, 1], marker="o")
ax.scatter(X_train[y_train == 1, 0], X_train[y_train == 1, 1], marker="x")
plt.show()

これが出力された散布図です。
一部ノイズはありますが概ね綺麗に別れていて、動作確認には手頃な問題になっていると思います。

早速ですが、kerasでモデルを作っていきます。
kerasのドキュメントはここです。日本語版があって便利ですね。
特に何も考えず、中間層が1層あるだけのシンプルなニューラルネットワークで作ります。
とりあえず動けば良いので、今回はDropoutもCallbackもなしで。


# モデルの構築
model = Sequential()
model.add(Dense(16, activation='tanh', input_shape=(2,)))
model.add(Dense(1, activation='sigmoid'))
model.summary()
model.compile(
        optimizer='adam',
        loss='binary_crossentropy',
        metrics=['acc']
)

# 以下、出力
_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
dense_1 (Dense)              (None, 16)                48        
_________________________________________________________________
dense_2 (Dense)              (None, 1)                 17        
=================================================================
Total params: 65
Trainable params: 65
Non-trainable params: 0
_________________________________________________________________

モデルができたので、学習させます。
epochsもbatch_sizeも適当です。


# 学習
history = model.fit(
                    X_train,
                    y_train,
                    epochs=100,
                    batch_size=32,
                    validation_data=[X_test, y_test]
            )

正常に学習が進んだことを確認するために、
損失関数と正解率の変動を可視化してみましょう。


# Epoch ごとの正解率と損失関数のプロット
fig = plt.figure(figsize=(12, 12))
ax = fig.add_subplot(2, 1, 1, title="loss")
ax.plot(history.epoch, history.history["loss"], label="train_loss")
ax.plot(history.epoch, history.history["val_loss"], linestyle="-.", label="val_loss")
ax.legend()
ax = fig.add_subplot(2, 1, 2, title="acc")
ax.plot(history.epoch, history.history["acc"], label="train_acc")
ax.plot(history.epoch, history.history["val_acc"], linestyle="-.", label="val_acc")
ax.legend()
plt.show()

できたグラフがこちら。
順調に学習できていますね。

訓練データで評価してみましょう。
もともとやさしい問題なので、なかなかの正解率です。


print(classification_report(y_train, model.predict_classes(X_train)))

# 出力結果
             precision    recall  f1-score   support

          0       0.95      0.96      0.95       800
          1       0.96      0.95      0.95       800

avg / total       0.95      0.95      0.95      1600

最後に、決定境界を可視化してみます。
以前紹介した、等高線のプロットを使います。


# 決定境界の可視化
X_mesh, Y_mesh = np.meshgrid(np.linspace(-2, 2, 401), np.linspace(-2, 2, 401))
Z_mesh = model.predict_classes(np.array([X_mesh.ravel(), Y_mesh.ravel()]).T)
Z_mesh = Z_mesh.reshape(X_mesh.shape)
fig = plt.figure(figsize=(10, 10))
ax = fig.add_subplot(1, 1, 1, aspect='equal', xlim=(-2, 2), ylim=(-2, 2))
ax.contourf(X_mesh, Y_mesh, Z_mesh, alpha=0.2)
plt.scatter(X_test[y_test == 0, 0], X_test[y_test == 0, 1], marker="o")
plt.scatter(X_test[y_test == 1, 0], X_test[y_test == 1, 1], marker="x")
plt.show()

出力がこちら。
真円にはなっていませんが、きちんと円の内側と外側にデータを分離できました。