StatsFragments

Python, R, Rust, 統計, 機械学習とか

Python pandas プロット機能を使いこなす

pandas は可視化のための API を提供しており、折れ線グラフ、棒グラフといった基本的なプロットを簡易な API で利用することができる。一般的な使い方は公式ドキュメントに記載がある。

これらの機能は matplotlib に対する 薄い wrapper によって提供されている。ここでは pandas 側で一処理を加えることによって、ドキュメントに記載されているプロットより少し凝った出力を得る方法を書きたい。

補足 サンプルデータに対する見せ方として不適切なものがあるが、プロットの例ということでご容赦ください。

パッケージのインポート

import matplotlib.pyplot as plt
plt.style.use('ggplot')

import matplotlib as mpl
mpl.__version__
# '1.5.0'

import numpy as np
np.__version__
# '1.10.1'

import pandas as pd
pd.__version__
# u'0.17.0'

折れ線グラフ、棒グラフ系

サンプルデータとして 気象庁から 東京、札幌、福岡の 2014 年の日別平均気温データをダウンロードし、pd.read_csv で読み込んだ。

df = pd.read_csv('data.csv', index_col=0, header=[0, 1, 2], skiprows=[0], encoding='shift-jis')
df = df.iloc[:, [0, 3, 6]]
df.columns = [u'東京', u'札幌', u'福岡']
df.index = pd.to_datetime(df.index)
df.head()

f:id:sinhrks:20151115214421p:plain

まずは単純な折れ線グラフを描く。

df.plot()

f:id:sinhrks:20151115214525p:plain

次に、適当な単位 (ここでは月次) でグルーピングして棒グラフを書きたい。集計には pd.Grouper を利用する。

参考 Python pandas アクセサ / Grouperで少し高度なグルーピング/集計 - StatsFragments

monthly = df.groupby(pd.Grouper(level=0, freq='M')).mean()
monthly

f:id:sinhrks:20151115215001p:plain

このとき、index が 時系列 ( DatetimeIndex ) の場合、そのままでは x 軸がタイムスタンプとして表示されてしまう。

monthly.plot.bar(color=['#348ABD', '#7A68A6', '#A60628'])

f:id:sinhrks:20151115215123p:plain

x軸 を簡単にフォーマットしたい場合、Index.format() メソッドを利用して index を文字列表現に変換するとよい。

monthly.index.format()
# ['2014-01-31', '2014-02-28', '2014-03-31', '2014-04-30', '2014-05-31', '2014-06-30',
#  '2014-07-31', '2014-08-31', '2014-09-30', '2014-10-31', '2014-11-30', '2014-12-31']
monthly.index = monthly.index.format()
monthly.plot.bar(color=['#348ABD', '#7A68A6', '#A60628'])

f:id:sinhrks:20151115214908p:plain

また、各棒の塗り分けには colormap を指定することもできる。

monthly.plot.bar(cmap='rainbow')

f:id:sinhrks:20151115215206p:plain

折れ線グラフ / 棒グラフを一つのプロットとして描画する場合は以下のようにする。.plot メソッドmatplotlib.axes.Axes インスタンスを返すため、続くプロットの描画先として その Axes を指定すればよい。

ax = monthly[u'札幌'].plot(legend=True)
monthly[[u'東京', u'福岡']].plot.bar(ax=ax, rot=30)

f:id:sinhrks:20151115215325p:plain

また、matplotlib 側で 極座標Axes を作成しておけば、レーダーチャートのような描画もできる。

補足 一部 線が切れているのは負値のため。

fig = plt.figure()
ax = fig.add_axes([0.1, 0.1, 0.9, 0.9], polar=True)
indexer = np.arange(1, 13) * 2 * np.pi / 12
monthly.index = indexer
monthly.append(monthly.iloc[0]).plot(ax=ax)
ax.set_xticks(indexer)
ax.set_xticklabels(np.arange(1, 13));

f:id:sinhrks:20151115215334p:plain

分布描画系

各都市の気温の分布が見たい、といった場合にはヒストグラムを描画する。

df.plot.hist(bins=100, alpha=0.5)

f:id:sinhrks:20151115215608p:plain

より細かい単位でグループ分けしてグラフを描きたい場合は、まず plt.subplots で必要な Axes を生成する。続けて、pandas でグループ化したデータを各 Axes に対して描画すればよい。

ここでは月ごとに サブプロットにして カーネル密度推定したグラフを描画してみる。

fig, axes = plt.subplots(4, 3, figsize=(8, 6))
plt.subplots_adjust(wspace=0.5, hspace=0.5)
for (ax, (key, group)) in zip(axes.flatten(), df.groupby(pd.Grouper(level=0, freq='M'))):
    ax = group.plot.kde(ax=ax, legend=False, fontsize=8)
    ax.set_ylabel('')
    ax.set_title(key, fontsize=8)

f:id:sinhrks:20151115215728p:plain

同様に箱ヒゲ図も描ける。

fig, axes = plt.subplots(4, 3, figsize=(8, 6))
plt.subplots_adjust(wspace=0.5, hspace=0.5)
for (ax, (key, group)) in zip(axes.flatten(), df.groupby(pd.Grouper(level=0, freq='M'))):
    ax = group.plot.box(ax=ax)
    ax.set_ylabel('')
    ax.set_title(key, fontsize=8)

f:id:sinhrks:20151115215752p:plain

散布図系

通常の散布図は以下のようにして描ける。

df.plot(kind='scatter', x=u'札幌', y=u'福岡')

f:id:sinhrks:20151115220121p:plain

各点を適当に色分けしたい場合、c キーワードで各点の色を指定する。また、colormap も合わせて指定できる。ここで指定した値は連続値として扱われるため colorbar が表示されている (非表示にすることもできる)。

df.plot(kind='scatter', x=u'札幌', y=u'福岡', c=df.index.month, cmap='winter')

f:id:sinhrks:20151115220154p:plain

各点をグループ別に塗り分ける ( 離散値として扱う ) 場合は、単一の Axes に対して グループ化したデータを順番に描画していけばよい。

cmap = plt.get_cmap('rainbow')
colors = [cmap(c / 12.0) for c in np.arange(1, 13)]
colors
# [(0.33529411764705885, 0.25584277759443558, 0.99164469551074275, 1.0),
#  (0.17058823529411765, 0.49465584339977881, 0.96671840426918743, 1.0),
#  ...
#  (1.0, 0.25584277759443586, 0.12899921653020341, 1.0),
#  (1.0, 1.2246467991473532e-16, 6.123233995736766e-17, 1.0)]

fig, ax = plt.subplots(1, 1)
for i, (key, group) in enumerate(df.groupby(pd.Grouper(level=0, freq='M')), start=1):
    group.plot(kind='scatter', x=u'札幌', y=u'福岡', color=cmap(i / 12.0), ax=ax, label=i)

f:id:sinhrks:20151115220859p:plain

また、ダミーのカラムを作成することで月次の気温の分布を散布図としてみることもできる。

df['month'] = df.index.month
df.plot(kind='scatter', x='month', y=u'札幌')

f:id:sinhrks:20151115220521p:plain

x 軸方向に散らす。

df['month2'] = df.index.month + np.random.normal(loc=0, scale=0.05, size=len(df.index.month))
df.plot(kind='scatter', x='month2', y=u'札幌')

f:id:sinhrks:20151115220531p:plain

バブルチャートも描ける。前処理として札幌の気温を離散化し、月次のデータ数を集計する。

df2 = df.copy()
df2[u'札幌2'] = df2[u'札幌'] // 10 * 10 
df2['count'] = 1
df2 = df2.groupby(['month', u'札幌2'])['count'].count().reset_index()
df2.head()

f:id:sinhrks:20151115220541p:plain

バブルの大きさは s キーワードで指定する。

df2.plot.scatter(x='month', y=u'札幌2', s=df2['count'] * 10)
df.plot(kind='scatter', x='month2', y=u'札幌')

f:id:sinhrks:20151115220549p:plain

さらに変形して plt.imshow でヒートマップを描画する。imshowpandasAPI にはないため、DataFrame.pipe でデータを渡す。

piv = pd.pivot_table(df, index=u'札幌2', columns=df.index.month, values=u'札幌', aggfunc='count')
piv = piv.fillna(0)
piv

f:id:sinhrks:20151115220601p:plain

piv.pipe(plt.imshow, cmap='winter')
ax = plt.gca()
ax.invert_yaxis()
ax.set_yticks(np.arange(4))
ax.set_yticklabels(piv.index);

f:id:sinhrks:20151115220608p:plain

まとめ

pandas での前処理 + 可視化機能の組み合わせを利用して、より柔軟にプロットを行う方法を記載した。pandas の裏側は ndarray のため、最後の例のように pandas 側に API がないプロットも簡単に描ける。

Pythonによるデータ分析入門 ―NumPy、pandasを使ったデータ処理

Pythonによるデータ分析入門 ―NumPy、pandasを使ったデータ処理

Rust で重回帰

簡単なアルゴリズムを実装しつつ Rust を学びたい。環境はこちら。

  • Rust(Nightly) rustc 1.6.0-nightly (d5fde83ae 2015-11-12)
  • rust-csv 0.14.3 : CSV の読み込みに利用
  • nalgebra 0.3.1 : 行列/ベクトルの処理に利用

Python や R と比較すると、型のために Vecslice が透過的に扱えず戸惑う。また、 Generics を使った場合に所有権まわりのエラーが出て解決できなかったため、浮動小数f64 のみの定義とした。もう少し理解できたら書き直したい。

extern crate csv;
extern crate nalgebra;

use std::vec::Vec;
use nalgebra::{DVec, DMat, Inv, Mean, ColSlice, Iterable};

fn main() {

    // R の "trees" データを使用
    let data = "Girth,Height,Volume
8.3,70,10.3
8.6,65,10.3
8.8,63,10.2
10.5,72,16.4
10.7,81,18.8
10.8,83,19.7
11.0,66,15.6
11.0,75,18.2
11.1,80,22.6
11.2,75,19.9
11.3,79,24.2
11.4,76,21.0
11.4,76,21.4
11.7,69,21.3
12.0,75,19.1
12.9,74,22.2
12.9,85,33.8
13.3,86,27.4
13.7,71,25.7
13.8,64,24.9
14.0,78,34.5
14.2,80,31.7
14.5,74,36.3
16.0,72,38.3
16.3,77,42.6
17.3,81,55.4
17.5,82,55.7
17.9,80,58.3
18.0,80,51.5
18.0,80,51.0
20.6,87,77.0";

    let mut ncols = 0;

    // http://burntsushi.net/rustdoc/csv/
    let mut reader = csv::Reader::from_string(data).has_headers(true);

    let mut x: Vec<f64> = vec![];
    for row in reader.decode() {
        let vals: Vec<f64> = row.unwrap();
        // 1 列目を目的変数, 残りは説明変数とする
        for i in 0 .. vals.len() {
            x.push(vals[i]);
        }
        // 説明変数の数 (各行は必ず同数のデータを含むと仮定)
        ncols = vals.len();
    }

    // レコード数
    let nrows: usize = x.len() / ncols;

    // ベクトルから nrows x ncols の DMat (動的サイズの行列) を作成
    // http://nalgebra.org/doc/nalgebra/struct.DMat.html
    let dx = DMat::from_row_vec(nrows, ncols, &x);

    // 重回帰
    let coefs = lm_fit(&dx);
    println!("結果 {:?}", &coefs);
}

fn lm_fit(data: &DMat<f64>) -> Vec<f64> {
    // 重回帰

    // 列ごとの平均値を計算
    let means = data.mean();

    let nrows = data.nrows();
    let nfeatures = data.ncols() - 1;

    // 偏差平方和積和行列
    // DMat::from_fn で、行番号 i, 列番号 j を引数とする関数から各要素の値を生成できる
    let smx = DMat::from_fn(nfeatures, nfeatures,
                            |i, j| sum_square(&data.col_slice(i + 1, 0, nrows),
                                              &data.col_slice(j + 1, 0, nrows),
                                              means[i+1], means[j+1]));
    // 偏差積和行列
    // DMat と Vec では演算ができないため、こちらも一列の DMat として生成
    let smy = DMat::from_fn(nfeatures, 1,
                            |i, j| sum_square(&data.col_slice(i + 1, 0, nrows),
                                              &data.col_slice(0, 0, nrows),
                                              means[i], means[0]));
    // 偏回帰係数を計算し、Vec に変換
    let mut res = (smx.inv().unwrap() * smy).to_vec();

    // 切片を計算し、0 番目の要素として挿入
    let intercept = (1..means.len()).fold(means[0], |m, i| m - res[i-1]*means[i]);
    res.insert(0, intercept);
    return res
}

fn sum_square(vec1: &DVec<f64>, vec2: &DVec<f64>, m1: f64, m2: f64) -> f64 {
    // 平方和
    let mut val: f64 = 0.;
    for (v1, v2) in vec1.iter().zip(vec2.iter()) {
        val += (v1 - m1) * (v2 - m2);
    }
    return val;
}

実行結果。左から切片と各説明変数の偏回帰係数。

# 結果 [10.81637050444831, -0.04548349100904309, 0.19517974893553702]

R の結果。

lm(Girth ~ Height + Volume, data = trees)
# 
# Call:
# lm(formula = Girth ~ Height + Volume, data = trees)
# 
# Coefficients:
# (Intercept)       Height       Volume  
#    10.81637     -0.04548      0.19518

2016/12/26追記 上では逆行列を使ったが、利用するパッケージに線形方程式を解く関数 / メソッドが用意されていればそちらを使った方が数値的に安定 / 高速。

github.com

岩波データサイエンス Vol.1

ご恵贈いただきました。 ありがとうございます! あわせてタスクもいただきました (下部)。

書籍のコンテンツ

各章ごとの内容は id:sfchaos さんが詳しく紹介されています。

d.hatena.ne.jp

まだ すべて読めていないのですが、以下 3 点がよいポイントだと思います。

  • 理論 と サンプルプログラム 両方の記載がある
  • BUGS, Stan, PyMC3 と主要なパッケージが網羅されている
  • サンプルは単純な回帰だけでなく 時系列 / 空間ベイズを含む

補足 書籍には コラム "Pythonとは" という データ分析視点での Python 紹介があるのですが、中身は結構な pandas 推しでした。著者の方、いったい何者なんだ...。

Stan 入門

依頼により、著者の松浦さんが作成した RStan サンプルの PyStan 版を作成させていただきました。

以下リポジトリの "Example..." フォルダ中に含まれる Jupyter Notebook ( .ipynb ) を開くと、 GitHub 上でプログラムと出力を確認することができます。

github.com

RStanPyStanAPI 差異については公式ドキュメントに記載がありますがアップデートがされていません。下表にサンプル中で利用したものを簡単にまとめます。

RStan PyStan 概要
stan(...) pystan.stan(...) モデルのコンパイルとサンプリングの実行
stan_model(...) pystan.StanModel(...) モデルのコンパイル
sampling(...) StanModel.sampling(...) サンプリングの実行
extract(...) StanFit4model.extract() サンプリング結果を取得
traceplot(...) StanFit4model.plot(...) サンプリングの経過をプロット

書籍や Rstan のサンプルと合わせてご参照ください。

岩波データサイエンス Vol.1

岩波データサイエンス Vol.1

pandas 0.17.0 の主要な変更点

先日 10/9 に pandas 0.17.0 がリリースされた。直近のバージョンアップの中では かなり機能追加が多いリリースとなった。

重要な変更は リリースノート にハイライトとして列挙しているのだが、これらはある程度 pandas を使いこなしている方向けの記載となっている。

そのため、ここでは よりライトなユーザ向けに重要と思われる変更を書く。特に、ユーザ側のプログラムに影響がある以下の3点について記載する。

  1. ソート API の統合 ( sort_values / sort_index )
  2. 重複削除 API の改善 ( drop_duplicates / duplicated )
  3. .plot アクセサの追加

準備

import numpy as np
import pandas as pd

np.__version__
# '1.10.1'

pd.__version__
# u'0.17.0'

1. ソート API の統合 ( sort_values / sort_index )

これまで一貫性がなかったソート API が以下のように統合された。以前のメソッド / 引数は deprecate されており、将来のバージョン ( 通例では 0.19.0 以降 ) で削除される。

値によるソート ( sort_values として統合):

0.16.2以前 0.17.0以降
Series.order() Series.sort_values()
Series.sort() Series.sort_values(inplace=True)
DataFrame.sort(columns=...) DataFrame.sort_values(by=...)

ラベル ( Index ) によるソート ( sort_index として統合):

0.16.2以前 0.17.0以降
Series.sort_index() Series.sort_index()
Series.sortlevel(level=...) Series.sort_index(level=...)
DataFrame.sort_index() DataFrame.sort_index()
DataFrame.sortlevel(level=...) DataFrame.sort_index(level=...)
DataFrame.sort() DataFrame.sort_index()

それぞれ、引数として以下のオプションが利用できる。

  • inplace: 破壊的にソートするかどうかを指定する。既定は False (非破壊)
  • axis: ソートの方向(軸)を指定する
  • ascending: ソート順序 (昇順/降順) を指定する。既定は True (昇順)
s = pd.Series(list('afcedb'), index=[1, 3, 2, 4, 0, 5])
s
# 1    a
# 3    f
# 2    c
# 4    e
# 0    d
# 5    b
# dtype: object

# 値によるソート
s.sort_values()
# 1    a
# 5    b
# 2    c
# 0    d
# 4    e
# 3    f
# dtype: object

# ラベル (index) によるソート
s.sort_index()
# 0    d
# 1    a
# 2    c
# 3    f
# 4    e
# 5    b
# dtype: object

補足 もともとの仕様は NumPy (というか Python) の .sort メソッドの挙動 ( 破壊的なソート) にあわせたものだった。今回の変更により NumPy とは差異が発生する。

a = np.array([3, 2, 1])
a
# array([3, 2, 1])

a.sort()
a
# array([1, 2, 3])

2. 重複削除 API の改善 ( drop_duplicates / duplicated )

重複した値の削除を行う drop_duplicates で、重複値の全削除ができるようになった。

0.16.2以前 0.17.0以降
.drop_duplicates() .drop_duplicates()
.drop_duplicates(take_last=True) .drop_duplicates(keep='last')
- .drop_duplicates(keep=False)

削除方法の指定は これまでの take_last から keep オプション に変更された。

s = pd.Series(list('abebdcd'))
s
# 0    a
# 1    b
# 2    e
# 3    b
# 4    d
# 5    c
# 6    d
# dtype: object

# 重複がある場合、最初に出現した値を残す
s.drop_duplicates()
# 0    a
# 1    b
# 2    e
# 4    d
# 5    c
# dtype: object

# 重複がある場合、最後に出現した値を残す
s.drop_duplicates(keep='last')
# 0    a
# 2    e
# 3    b
# 5    c
# 6    d
# dtype: object

# 重複した値を残さない
s.drop_duplicates(keep=False)
# 0    a
# 2    e
# 5    c
# dtype: object

また、値が重複しているかどうかを調べる .duplicated でも keep オプションが利用できる。

# 重複がある場合、最初に出現した値以外を True 
s.duplicated()
# 0    False
# 1    False
# 2    False
# 3     True
# 4    False
# 5    False
# 6     True
# dtype: bool

# 重複がある場合、最後に出現した値以外を True
s.duplicated(keep='last')
# 0    False
# 1     True
# 2    False
# 3    False
# 4     True
# 5    False
# 6    False
# dtype: bool

# 重複 全てを True
s.duplicated(keep=False)
# 0    False
# 1     True
# 2    False
# 3     True
# 4     True
# 5    False
# 6     True
# dtype: bool

Index の値から重複を全削除する場合は Index.duplicated を以下のように利用すればよい。

s = pd.Series(np.arange(6), index=list('abcbda'))
s
# a    0
# b    1
# c    2
# b    3
# d    4
# a    5
# dtype: int64

s[~s.index.duplicated(keep=False)]
# c    2
# d    4
# dtype: int64

3. .plot アクセサの追加

DataFrameSeries をプロットする際のグラフの種類が、 .plot.bar().plot.hist() のように .plot をアクセサとして指定できるようになった。

これまでのメソッド呼び出し ( .plot(kind='bar') ) は引き続き利用できる。

df = pd.DataFrame(np.random.randn(5, 10))

# df.plot(kind='box') と同じ
df.plot.box()

f:id:sinhrks:20151017203318p:plain

その他

また、上記ほどではないが重要な内容を列挙する。気になるものがあれば詳細はリンク先で。

  • GIL 解放 ( Dask 利用時のパフォーマンス向上) (詳細)
  • pandas.io.data を deprecate し、別パッケージ pandas-datareader として分離 (詳細)
  • pd.to_datetime のエラー処理の変更、日時パース時の挙動統一 (詳細)
  • Index 同士の比較演算時の一部挙動の変更 (詳細)
  • ターミナル上での日本語データ表示時の位置補正オプションの追加 (詳細)

まとめ

pandas 0.17.0 での特に重要な変更点 3 点を記載した。

  1. ソート API の統合 ( sort_values / sort_index )
  2. 重複削除 API の改善 ( drop_duplicates / duplicated )
  3. .plot アクセサの追加

これら以外の変更点については リリースノート を一読お願いします。

10/18 編集 コメントご指摘により誤記を修正

Pythonによるデータ分析入門 ―NumPy、pandasを使ったデータ処理

Pythonによるデータ分析入門 ―NumPy、pandasを使ったデータ処理

PyConJP 2015: pandas/Daskについてお話させていただきました

10日、11日と PyCon JP に参加させていただきました。ご参加いただいた皆様、スタッフの皆様ありがとうございました。資料はこちらになります。

pandas internals

パフォーマンス向上のための pandas 内部実装の説明といくつかの TIPS について。そのうち翻訳するかもしれません。

speakerdeck.com

Dask: 軽量並列分散フレームワーク (LT)

speakerdeck.com

元ネタ

以下のエントリをベースに、それぞれ新しい内容を追加しています。

sinhrks.hatenablog.com

sinhrks.hatenablog.com

Python XGBoost + pandas 連携の改善

一部 こちらの続き。その後 いくつかプルリクを送りXGBoostpandas を連携させて使えるようになってきたため、その内容を書きたい。

sinhrks.hatenablog.com

できるようになったことは 以下 3 点。

  1. DMatrix でのラベルと型の指定
  2. pd.DataFrame からの DMatrix の作成
  3. xgb.cv の結果を pd.DataFrame として取得

補足 XGBoost では PyPI の更新をスクリプトで不定期にやっているようで、同一バージョンに見えても枝番が振られていたりして見分けにくい。記載は本日時点のこのコミットの情報。

%matplotlib inline
import numpy as np
import xgboost as xgb
from sklearn import datasets

import matplotlib.pyplot as plt
plt.style.use('ggplot')

xgb.__version__
# '0.4'

1. DMatrix でのラベルと型の指定

これは pandas 関係ない。DMatrix に以下 2 つのプロパティが追加され、任意のラベル/型が指定できるようにした (指定しない場合はこれまでと同じ挙動)。

  • feature_names: 変数のラベルを指定。ラベルは英数字のみ。
  • feature_types: 変数の型を指定。指定できるのは、
    • q: 量的変数 (既定)
    • i: ダミー変数
    • int: 整数型
    • float: 浮動小数点型。

まずは feature_names のみを指定。

iris = datasets.load_iris()

features = ['SepalLength', 'SepalWidth', 'PetalLength', 'PetalWidth']
dm = xgb.DMatrix(iris.data, label=iris.target, feature_names=features)
dm.feature_names
# ['SepalLength', 'SepalWidth', 'PetalLength', 'PetalWidth']

DMatrix で指定した feature_namesxgb.Booster へも引き継がれ、.get_dump 時 やプロットに利用できる。

params={'objective': 'multi:softprob',
        'eval_metric': 'mlogloss',
        'eta': 0.3,
        'num_class': 3}

np.random.seed(1)
bst = xgb.train(params, dm, num_boost_round=18)
xgb.plot_importance(bst)

f:id:sinhrks:20151003080826p:plain

次に feature_types。これは内部的な学習には全く関係がなく、.get_dumpxgb.plot_tree をする際の ラベルのフォーマットを決めるだけ。例えば以下のようなデータがあったとする。

data = np.array([[1, 0.1, 0], [2, 0.2, 1], [2, 0.2, 0]])
data
# array([[ 1. ,  0.1,  0. ],
#        [ 2. ,  0.2,  1. ],
#        [ 2. ,  0.2,  0. ]])

このとき、1 列目を int, 2 列目を float, 3 列目をダミー変数 i として扱いたければ以下のように指定すると、 plot_tree したときの見た目が若干変わる。全変数 既定 (q: 量的変数) でもとくに不自由ないので、強いこだわりのある方以外は指定する必要はない。

dm = xgb.DMatrix(data, feature_names=['A', 'B', 'C'],
                 feature_types=['int', 'float', 'i'])
dm.feature_types
# ['int', 'float', 'i']

補足 これまで、変数のラベル/型の指定は "featmap.txt" のような設定ファイルから別途指定する必要があった。この指定が Python で直接できるようになった。

2 pd.DataFrame からの DMatrix の作成

pandas.DataFrame から DMatrix が作成できるようにした。このとき、feature_namesfeature_typesDataFrame の定義から適当に設定される。

補足 DMatrixlabel には pd.Series が渡せる (前からできた)。

import pandas as pd

# 表示行数を指定
pd.set_option('display.max_rows', 8)
# 一行に表示する文字数を指定
pd.set_option('display.width', 120)

# pd.DataFrame を作成
train = pd.DataFrame(iris.data, columns=['SepalLength', 'SepalWidth', 'PetalLength', 'PetalWidth'])
train
#      SepalLength  SepalWidth  PetalLength  PetalWidth
# 0            5.1         3.5          1.4         0.2
# 1            4.9         3.0          1.4         0.2
# 2            4.7         3.2          1.3         0.2
# 3            4.6         3.1          1.5         0.2
# ..           ...         ...          ...         ...
# 146          6.3         2.5          5.0         1.9
# 147          6.5         3.0          5.2         2.0
# 148          6.2         3.4          5.4         2.3
# 149          5.9         3.0          5.1         1.8
# 
# [150 rows x 4 columns]

# DMatrix を作成
dm = xgb.DMatrix(train, label=pd.Series(iris.target))

dm.feature_names
# ['SepalLength', 'SepalWidth', 'PetalLength', 'PetalWidth']

dm.feature_types
# ['q', 'q', 'q', 'q']

DMatrix に変換できるのは、列が int64, float64, boolDataFrame のみ。他の型を自動でダミー変数に変換したりはしない。

df = pd.DataFrame([[1, 0.1, False], [2, 0.2, True], [2, 0.2, False]],
                  columns=['A', 'B', 'C'])
df
#    A    B      C
# 0  1  0.1  False
# 1  2  0.2   True
# 2  2  0.2  False

df.dtypes
# A      int64
# B    float64
# C       bool
# dtype: object

dm = xgb.DMatrix(df)

dm.feature_names
# ['A', 'B', 'C']

dm.feature_types
# ['int', 'q', 'i']

# object 型を含むためエラー
df = pd.DataFrame([[1, 0.1, 'x'], [2, 0.2, 'y'], [2, 0.2, 'z']],
                  columns=['A', 'B', 'C'])

df
#    A    B  C
# 0  1  0.1  x
# 1  2  0.2  y
# 2  2  0.2  z

df.dtypes
# A      int64
# B    float64
# C     object
# dtype: object

xgb.DMatrix(df)
# ValueError: DataFrame.dtypes must be int, float or bool

3. xgb.cv の結果を pd.DataFrame として取得

クロスバリデーションを行う xgb.cv はこれまで実行結果を文字列で返していたため、その結果をパースもしくは目視で確認する必要があった。この返り値を pd.DataFrame もしくは np.ndarray とし、プログラムで処理がしやすいようにした。

train = pd.DataFrame(iris.data, columns=['SepalLength', 'SepalWidth', 'PetalLength', 'PetalWidth'])
dm = xgb.DMatrix(train, label=pd.Series(iris.target))

cv = xgb.cv(params, dm, num_boost_round=50, nfold=10)
cv
#     test-mlogloss-mean  test-mlogloss-std  train-mlogloss-mean  train-mlogloss-std
# 0             0.753459           0.027033             0.737631            0.003818
# 1             0.552303           0.048738             0.526929            0.005102
# 2             0.423481           0.066469             0.390115            0.005873
# 3             0.339942           0.082163             0.295637            0.006148
# ..                 ...                ...                  ...                 ...
# 46            0.208001           0.259161             0.018070            0.001759
# 47            0.208355           0.261166             0.017898            0.001724
# 48            0.208468           0.261520             0.017755            0.001703
# 49            0.208566           0.260967             0.017617            0.001686
# 
# [50 rows x 4 columns]

cv.sort(columns='test-mlogloss-mean')
#     test-mlogloss-mean  test-mlogloss-std  train-mlogloss-mean  train-mlogloss-std
# 11            0.175506           0.177823             0.058496            0.004715
# 13            0.175514           0.194641             0.045222            0.004117
# 14            0.176150           0.203265             0.040493            0.004050
# 12            0.176277           0.186902             0.051074            0.004423
# ..                 ...                ...                  ...                 ...
# 3             0.339942           0.082163             0.295637            0.006148
# 2             0.423481           0.066469             0.390115            0.005873
# 1             0.552303           0.048738             0.526929            0.005102
# 0             0.753459           0.027033             0.737631            0.003818
# 
# [50 rows x 4 columns]

pandas がインストールされていない、もしくは as_pandas=False を指定した場合、返り値は np.ndarray となる。

cv = xgb.cv(params, dm, num_boost_round=50, nfold=10, as_pandas=False)
cv
# array([[ 0.7534586 ,  0.02703308,  0.7376308 ,  0.00381775],
#        ...
#        [ 0.2085664 ,  0.26096721,  0.0176169 ,  0.00168606]])

np.argmin(cv, axis=0)
# array([11,  0, 49, 49])

まとめ

XGBoost と pandas をより便利に使うため、以下 3 点の修正を行った。

  1. DMatrix でのラベルと型の指定
  2. pd.DataFrame からの DMatrix の作成
  3. xgb.cv の結果を pd.DataFrame として取得

これ使って Kaggle で大勝利したい...と思ったのですが、他の日本勢に使っていただいたほうが勝率高いはずなのでシェアします。

Python Dask で 並列 DataFrame 処理

はじめに

先日のエントリで少し記載した Dask について、その使い方を書く。Dask を使うと、NumPypandasAPI を利用して並列計算/分散処理を行うことができる。また、Dask は Out-Of-Core (データ量が多くメモリに乗らない場合) の処理も考慮した実装になっている。

sinhrks.hatenablog.com

上にも書いたが、DaskNumPypandas を置き換えるものではない。数値計算のためのバックエンドとして NumPypandas を利用するため、むしろこれらのパッケージが必須である。

DaskNumPypandasAPI を完全にはサポートしていないため、並列 / Out-Of-Core 処理が必要な場面では Dask を、他では NumPy / pandas を使うのがよいと思う。pandasDask のデータはそれぞれ簡単に相互変換できる。

補足 とはいえ都度の変換は手間なので、pandas の処理実行時に Dask を利用するオプションをつける という検討はされている。

インストール

pip もしくは conda で。

pip install dask

準備

まずは必要なパッケージを import する。

import numpy as np
import pandas as pd
import dask.dataframe as dd

np.__version__
# '1.9.3'

pd.__version__
# '0.16.2'

# バージョン表示のためにインポート。
import dask
dask.__version__
# '0.7.1'

pandas から Dask への変換

サンプルデータは すでにメモリ上にある pd.DataFrame とする。

df = pd.DataFrame({'X': np.arange(10), 
                   'Y': np.arange(10, 20),
                   'Z': np.arange(20, 30)},
                  index=list('abcdefghij'))
df
#    X   Y   Z
# a  0  10  20
# b  1  11  21
# c  2  12  22
# d  3  13  23
# e  4  14  24
# f  5  15  25
# g  6  16  26
# h  7  17  27
# i  8  18  28
# j  9  19  29

pandas のデータ構造から Dask に変換するには dd.from_pandas。2つめの引数で データをいくつのパーティションに分割するかを指定している。結果は dask.dataframe.DataFrame ( dd.DataFrame ) となる。

divisions はデータがどこで分割されたかを示す。表示から、1つ目のパーティションには .index が "a" 〜 "e" までのデータが、2つ目のには "f" 〜 "j" までのデータが含まれていることがわかる。

重要 dd.DataFrame の処理全般について、行の順序は保証されない。各パーティションには divisions で示される .index を持つ行が含まれるが、パーティション内が常にソートされているとは限らない。

ddf = dd.from_pandas(df, 2)
ddf
# dd.DataFrame<from_pandas-b08addf72f0693a9fa1bb6c21d91f3d4, divisions=('a', 'f', 'j')>

# DataFrame の列名
ddf.columns
# ('X', 'Y', 'Z')

# DataFrame の index
ddf.index
# dd.Index<from_pandas-b08addf72f0693a9fa1bb6c21d91f3d4-index, divisions=('a', 'f', 'j')>

# DataFrame の divisions (パーティションの分割箇所)
ddf.divisions
# ('a', 'f', 'j')
 
# DataFrame のパーティション数
ddf.npartitions
# 2

dd.DataFrame からデータを取得する (計算処理を実行する) には .compute()。結果、元の pd.DataFrame が得られる。

ddf.compute()
#    X   Y   Z
# a  0  10  20
# b  1  11  21
# c  2  12  22
# d  3  13  23
# e  4  14  24
# f  5  15  25
# g  6  16  26
# h  7  17  27
# i  8  18  28
# j  9  19  29

内部処理

ここから、dd.DataFrame でどういった処理ができるのか、内部動作とあわせて記載する。といっても難しいことは全くやっていない。

まず、データ全体について 1 加算する処理を考える。これは 各パーティションごとに 1 加算して連結するのと同じ。

ddf + 1
# dd.DataFrame<elemwise-5b9ae0407254158903343113fac41421, divisions=('a', 'f', 'j')>

(ddf + 1).compute()
# 略

f:id:sinhrks:20150924210655p:plain

次に、列ごとの合計をとる処理。これは、各パーティションごとに列の合計をとって連結し、もう一度 合計をとる処理と同じ。

列の合計をとるため、結果は Series になる。ddf.sum() の時点では .compute() が呼ばれていないため実際の計算は行われていないが、Dask が結果の型 や divisions を推測し正しく表示している。

ddf.sum()
# dd.Series<dataframe-sum--second-7ba12c9d58c17f61406b84b6c30d7a26, divisions=(None, None)>

ddf.sum().compute()
# X     45
# Y    145
# Z    245
# dtype: int64

f:id:sinhrks:20150924210729p:plain

Dask ではこのような形で、計算処理をパーティションごとに並列 / Out-Of-Core 実行できる形に読み替えている。これらの処理は内部的には Computational Graph ( Dask Graph ) として表現され、.compute() によって実行される。

各処理の Dask Graph は、.visualize() メソッドを利用して確認できる。Graph 上で縦につながっていない処理同士は並列で実行できる。

ddf.sum().visualize()

f:id:sinhrks:20150924212007p:plain

各列の平均をとる場合、内部的には各列の .sum()と 各列の .count() をそれぞれ計算して除算。

ddf.mean().compute()
# X     4.5
# Y    14.5
# Z    24.5
# dtype: float64

ddf.mean().visualize()

f:id:sinhrks:20150924212112p:plain

DataFrame 同士の演算や、演算をチェインすることもできる。互いのパーティションが異なる場合はそれらが一致するよう調整が行われる。

((ddf - (ddf * 2)) == - ddf).visualize() 

f:id:sinhrks:20150924212136p:plain

また、累積関数 ( cumxxx ) や ウィンドウ関数 ( rolling_xxx ) なども利用できる。

ddf.cumsum().compute()
#     X    Y    Z
# a   0   10   20
# b   1   21   41
# c   3   33   63
# d   6   46   86
# e  10   60  110
# f  15   75  135
# g  21   91  161
# h  28  108  188
# i  36  126  216
# j  45  145  245

ddf.cumsum().visualize()

f:id:sinhrks:20150924213152p:plain

concat, join などの 連結 / 結合もできる。通常の演算と同じく、dd.DataFrame 同士のパーティションは適当に調整される。

df2 = pd.DataFrame({'A': np.arange(5), 
                    'B': np.arange(10, 15)},
                   index=list('adefg'))
df2
#    A   B
# a  0  10
# d  1  11
# e  2  12
# f  3  13
# g  4  14

ddf2 = dd.from_pandas(df2, 2)
ddf2
# dd.DataFrame<from_pandas-667963fc37e22688843f02da80df5963, divisions=('a', 'f', 'g')>

ddf.join(ddf2).compute()
#    X   Y   Z   A   B
# a  0  10  20   0  10
# b  1  11  21 NaN NaN
# c  2  12  22 NaN NaN
# d  3  13  23   1  11
# e  4  14  24   2  12
# f  5  15  25   3  13
# g  6  16  26   4  14
# h  7  17  27 NaN NaN
# i  8  18  28 NaN NaN
# j  9  19  29 NaN NaN

ddf.join(ddf2).visualize()

f:id:sinhrks:20150924214828p:plain

サポートされている処理の一覧は以下のAPIドキュメントを。一部利用できない引数が明記されていないが、次バージョンにて改訂

9/26 追記 処理結果については、行の順序以外は pandas の処理と一致するはず。例外は quantile のような percentile をとる処理。これらは Out-Of-Core 処理のための近似アルゴリズムを使っており、正確な値とずれることがある。

実データでの利用例

こちらが良エントリ (英語)。

まとめ

Dask を利用して DataFrame を並列処理する方法を記載した。手順は、

  • dd.from_pandas を利用して pd.DataFramedd.DataFrame へ変換。
  • 実行したいメソッド / 演算を dd.DataFrame に対して適用。
  • .compute() で計算を実行し、結果を取得する。計算処理は Dask にて自動的に並列化される。

最後、pandas 0.16.2 時点では並列処理による速度向上は大きくはない。これは Python の GIL (Global Interpreter Lock ) により並列実行できる処理が限定されているため。今月中にリリース予定の pandas 0.17.0 では いくつかの処理で Cython から明示的に GIL 解放するよう実装を変更しており、並列化による速度向上は大きくなる。

参考 Pandas Releasing the GIL

Python 次世代の多次元配列パッケージ群

このところ、たびたび NumPy 後継が...とか 並列処理が...という話を聞くので、この秋 注目の多次元配列パッケージをまとめたい。

バックエンド系

NumPy のように数値計算処理を自前で実装しているパッケージ。

DyND

Blaze プロジェクトのひとつ。C++ 実装 + Python バインディングGitHub にいくつか Example があがっているが、複合型やカテゴリカル型、GroupBy 操作がサポートされていて熱い。ラベルデータも NumPy より簡単に実装できそうだ。

speakerdeck.com

並列分散系

自身では直接 数値計算処理を行わず、バックエンド ( 主に NumPy )を利用して並列/分散処理を行うパッケージ。1 物理PC/複数コアでの並列計算を主用途とし、NumPy, pandas では少し苦しいが PySpark などを使うほどじゃない...という場合に利用するもの。

Dask

Blaze プロジェクトのひとつ。Pure-Python 実装。主用途は単一物理PCでの複数コア計算だが、複数PCに処理を分散させることもできる。

NumPyAPI で並列計算を行う dask.array のほか、toolz 相当の操作を行う dask.bagpandas 相当の処理を行う dask.dataframe など一連のデータ構造が揃っている。

開発者はベースパッケージである NumPytoolzpandas のコミッタとの兼任が多い。自分もコミット権限をいただいており、dask.dataframe への API 追加を行なっている。

speakerdeck.com

DistArray

NumPyAPI で並列計算を行う DistArray のみをサポート、ほかのデータ構造はない。Pure-Python 実装。見た感じ、並列処理の基本的な考え方は Dask と同じようだ。Enthought が好きな方はこちらを使えばよいかと。

docs.enthought.com

bolt

NumPyAPI でローカル計算/分散処理を行う。Pure-Python実装。ローカル計算は NumPy をそのまま使い (並列処理しない)、分散処理は Spark で行う。海外の方が期待の新星っぽい扱いをしていたので気になったのだが、現時点で実装されている API は多くはない。

ほかにもいくつかプロジェクトがあるが、自分としては Dask を流行らせたいので、その記事を書きます。

9/24 追記 書きました。

sinhrks.hatenablog.com

NumPy でつくる俺々データ構造

はじめに

Python での数値計算の基盤をなす NumPy 、直感的なスライスやブロードキャスト、関数のベクトル適用など大変便利だ。

import numpy as np
np.__version__
# '1.9.2'

np.array([1, 2, 3])
# array([1, 2, 3])

np.array([1, 2, 3])[:2]
# array([1, 2])

np.array([1, 2, 3]) + 1
# array([2, 3, 4])

が、用途によっては NumPy 標準ではその機能を実現できない場合がある。例えば、

  • 配列とメタデータをひとつのクラスで扱いたい
  • 配列への入力や型を制約/検証したい
  • 自作クラスを NumPy の Universal Functions (ufunc) に対応させたい
  • 新しい型 ( dtype ) を作りたい

こういったとき、NumPy の動作を継承した自作のクラスが作れるとうれしい。方法として、大きく二つがある。

  1. NumPy ndarray を継承する
  2. NumPy Array Interface を実装する

どちらを使うべきかは目的による。サブクラス化しないほうが柔軟な処理ができるが、互換をとるための実装は多くなる。pandas の主要データ構造では 柔軟性を重視して ndarray の継承をやめ、現在は 2 の方法を取っている。

ndarray のサブクラスを作成する方法は以下のドキュメントに記載されている。きちんとやりたい場合はこちらを読むのがよい。

このエントリでは 2. Array Interface を使って最低限の動作をさせる方法を書く。

クラスの定義

サンプルとして、名前をもつ配列クラス NamedArray を考える。__init__での初期化時に、名前をあらわす name と配列の値 values を渡すことにする。

class NamedArray(object):
    
    def __init__(self, name, values):
        self.name = name
        self.values = np.array(values)
        
    def __repr__(self):
        return 'NamedArray: {0}: {1}'.format(self.name, self.values)
    
    def __len__(self):
        return len(self.values)
        
n = NamedArray('x', [1, 2, 3])
n
# NamedArray: x: [1 2 3]

len(n)
# 3

このクラスに ndarray と互換の挙動をさせたい。まず、NamedArray インスタンスから ndarray が直接作成できるとうれしい。が、NamedArraynp.array に渡すと shape のない ndarray (0-dimmentional array) が作成されてしまう。

np.array(n)
# array(NamedArray: x: [1 2 3], dtype=object)

np.array(n).shape
# ()

理由は np.arrayNamedArrayスカラーとして扱うため。期待の動作をさせるためには、np.arrayNamedArray クラスが配列であることを教える必要がある。

具体的には NamedArray に Array Interface ( __array__ ) を定義すればよい。以下のように __array__ メソッドを追加したクラス NPNamedArray を作成すると、期待通り values の値から ndarray が作成されていることがわかる。

class NPNamedArray(NamedArray):
    
    def __array__(self):
        return self.values
    
n = NPNamedArray('x', [1, 2, 3])
np.array(n)
# array([1, 2, 3])

np.array(n).shape
# (3,)

ufunc 処理のオーバーライド

また、Array Interface をもつインスタンスには ufunc が適用できるようになる。引数は ufunc 内で ndarray に変換されるため、返り値は ndarray になる。

np.sin(n)
# array([ 0.84147098,  0.90929743,  0.14112001])

このときの内部の動きをみてみる。ufunc 適用後の処理は __array_wrap__ メソッドによって制御することができる。

__array_wrap__ メソッドには out_arrcontext 2 つの引数を定義する必要がある。これらに何が渡されているのか、print して確認する。

class NPNamedArray(NamedArray):
    
    def __array__(self):
        return self.values
    
    def __array_wrap__(self, out_arr, context=None):
        print('out_arr: ', out_arr)
        print('context: ', context)
        return out_arr
    
n = NPNamedArray('x', [1, 2, 3])
np.sin(n)
# ('out_arr: ', array([ 0.84147098,  0.90929743,  0.14112001]))
# ('context: ', (<ufunc 'sin'>, (NamedArray: x: [1 2 3],), 0))
# array([ 0.84147098,  0.90929743,  0.14112001])

out_arr には ufunc が適用された結果の ndarray が、 context には ufunc を呼び出した際のシグネチャ ( ufunc 自体, ufunc の引数, ufunc のドメイン ) が tuple で渡ってきている。__array_wrap__ の返り値を NPNamedArray に書き換えてやれば、期待する動作ができそうだ。

補足 np.isrealnp.ones_like など、Universal Functions に記載されていても __array_wrap__ を通らないものもある。

補足 ufunc のドメインは NumPy の場合 (おそらく) 全部 0。

class NPNamedArray(NamedArray):
    
    def __array__(self):
        return self.values
    
    def __array_wrap__(self, out_arr, context=None):
        return NPNamedArray(self.name, out_arr)
    
n = NPNamedArray('x', [1, 2, 3])
np.sin(n) 
# NamedArray: x: [ 0.84147098  0.90929743  0.14112001]

ここまでの時点で、NPNamedArray には __add__ メソッドを定義していないため、通常の加算はエラーになる。

n + 1
# TypeError: unsupported operand type(s) for +: 'NPNamedArray' and 'int'

が、ufunc である np.addndarray を含む演算は利用できるようになる。これは先日のエントリに記載した通り ndarray 側に Array Interface をもつクラスとの演算の定義があるため。

np.add(n, 1)
# NamedArray: x: [2 3 4]

n + np.array([5, 5, 5])
# NamedArray: x: [6 7 8]

また、__array_wrap__context 引数の値によって、特定の ufunc の処理をオーバーライドすることができる。例えば np.sin を適用したときに任意の値を返すこともできる。

class NPNamedArray(NamedArray):
    
    def __array__(self):
        return self.values
    
    def __array_wrap__(self, out_arr, context=None):
        if context[0] is np.sin:
            return u'返したい値'
        return NPNamedArray(self.name, out_arr)
    
n = NPNamedArray('x', [1, 2, 3])
np.sin(n)
# u'ひまわり'

とはいえ、上のように個々の関数について条件分岐を作成するのは手間だ。こんなとき、 ufunc 自体のプロパティを参照することで より簡単に書ける場合がある。

例えば、特定の型への処理をサポートする ufunc のみを利用する場合は .types プロパティを参照すればよい。このプロパティは ufunc がどのような引数/返り値の型を取るかを示す Array-protocol Type String を返す。

参考

np.sin.types
# ['e->e', 'f->f', 'd->d', 'g->g', 'F->F', 'D->D', 'G->G', 'O->O']

np.cos.types
# ['e->e', 'f->f', 'd->d', 'g->g', 'F->F', 'D->D', 'G->G', 'O->O']

np.add.types
# ['??->?', 'bb->b', 'BB->B', 'hh->h', 'HH->H', 'ii->i', 'II->I', 'll->l', 'LL->L',
#  'qq->q', 'QQ->Q', 'ee->e', 'ff->f', 'dd->d', 'gg->g', 'FF->F', 'DD->D', 'GG->G',
#  'Mm->M', 'mm->m', 'mM->M', 'OO->O']

まとめ

NumPy 互換の動作を実現する方法のうち、Array Interface を利用した方法を記載した。自分で一から実装するよりは簡単だと思う。

Pythonによるデータ分析入門 ―NumPy、pandasを使ったデータ処理

Pythonによるデータ分析入門 ―NumPy、pandasを使ったデータ処理

Python でパイプ演算子を使いたい


ネタ記事です。/ This is a joke post which makes no practical sense.


はじめに

Python pandas では主要な操作を以下のようにメソッドチェインの形で書くことができる。

# Python (pandas)
df.assign(x=df['y'] + df['z']).groupby('x').sum()

pandas v0.16.2 で DataFrameSeries.pipe というメソッドが追加され、このチェインを外部の関数/メソッドに対して連結できるようになった。利用例は以下のリンクを。

補足 matplotlib でも v1.5.0 で ラベルデータ対応 が追加され、各関数が .pipe から利用できるようになる予定。

このメソッドチェインによる "処理の連結" を、R ( {magrittr} ) のようにパイプ演算子 %>% を使って統一的に書きたい。

# R (magrittr + dplyr)
df %>% mutate(x = y + z) %>% group_by(x) %>% summarize_each(funs(sum))

R と異なり、 Python では自作の演算子を定義することはできない。Python演算子の中では 右方向ビットシフト演算子 >> がパイプ演算子っぽく見えなくもない。これを使って R のような動作をさせたい。

演算子のオーバーライド

まずはサンプルデータを作成する。

import numpy as np
import pandas as pd
pd.__version__
# '0.16.2'

df = pd.DataFrame({'X': [1, 2, 3, 4, 5],
                   'Y': ['a', 'b', 'a', 'b', 'a']})
df
#    X  Y
# 0  1  a
# 1  2  b
# 2  3  a
# 3  4  b
# 4  5  a

普通に .pipe を使うとこのような感じになる。

func = lambda x: x
df.pipe(func)
#    X  Y
# 0  1  a
# 1  2  b
# 2  3  a
# 3  4  b
# 4  5  a

この処理を >> 演算子を使って書きたい。そのままでは、演算方法が未定義のため TypeError となる。

df >> func
# TypeError: unsupported operand type(s) for >>: 'DataFrame' and 'function'

Python では各演算子に対応する処理はマジックメソッドとして実装する。そのため、>> に対応する __rshift__.pipe と同じ処理をさせればよい。また、R のように 右辺には関数を渡したいので、pandas のメソッドを適当にラップする関数 pipify を定義する。

pd.DataFrame.__rshift__ = pd.DataFrame.pipe

def pipify(f):
    def wrapped(*args, **kwargs):
        return lambda self: f(self, *args, **kwargs)
    return wrapped

head = pipify(pd.DataFrame.head)

できた。

df >> head(n=2)
#    X  Y
# 0  1  a
# 1  2  b

groupby も同じように書ける。集約関数は DataFrameDataFrameGroupBy 双方から使うため以下のように定義する。

groupby = pipify(pd.DataFrame.groupby)

def sum(*args, **kwargs):
    return lambda self: self.sum(*args, **kwargs)

また、groupby した結果からも >> で演算を連結したい。そのためにDataFrameGroupBy.__rshift__ を定義する。

def _pipe(self, func):
    return func(self)

pd.core.groupby.DataFrameGroupBy.__rshift__ = _pipe

できた。

df >> groupby('Y') >> sum()
#    X
# Y   
# a  9
# b  6

df >> sum()
# X       15
# Y    ababa
# dtype: object

右辺側クラスによるオーバーライド

DataFrame が左辺にある場合、パイプ演算は上記のように定義できた。加えて、DataFrame が右辺にある場合も同じようなことがしたい。

以下のような記法で、左辺の辞書から DataFrame を作成したい。これも (当然) そのままではできない。

data = {'X': [1, 2, 3, 4, 5],
        'Y': ['a', 'b', 'a', 'b', 'a']} 

data >> pd.DataFrame()
# TypeError: unsupported operand type(s) for >>: 'dict' and 'DataFrame'

Python では、左辺側に 右辺を処理できるマジックメソッドがない場合は、右辺側のマジックメソッド ( DataFarme.__rrshift__ ) が利用される。これを適当に定義すればよい。

def _get(a, b):
    return pd.DataFrame(b)

pd.DataFrame.__rrshift__ = _get

できた。

data = {'X': [1, 2, 3, 4, 5],
        'Y': ['a', 'b', 'a', 'b', 'a']} 

data >> pd.DataFrame()
#    X  Y
# 0  1  a
# 1  2  b
# 2  3  a
# 3  4  b
# 4  5  a

NumPy に対するオーバーライド

DataFramenp.ndarray からも作成することができる。が、>> ではうまく動かない。

data = np.array([[1, 2], [3, 4]])

pd.DataFrame(data)
#    0  1
# 0  1  2
# 1  3  4

data >> pd.DataFrame()
# TypeError: ufunc 'right_shift' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''

ndarray のマジックメソッドは スカラーndarray ライクな引数を汎用的に処理できるように定義されている。具体的には、引数が ndarray のサブクラス、もしくは Array Interface ( __array__ ) をもつとき、ndarray によって処理が行われる。

DataFramendarray のサブクラスではないが Array Interface をもつ。右辺が不正なため通常の例外が送出されてしまい、右辺側のマジックメソッドに処理が移らない。

参考 Subclassing ndarray — NumPy v1.11 Manual

issubclass(pd.DataFrame, np.ndarray)
# False

df.__array__()
# array([[1, 'a'],
#        [2, 'b'],
#        [3, 'a'],
#        [4, 'b'],
#        [5, 'a']], dtype=object)

ここで DataFrame 側のマジックメソッドを利用したい場合、上記リンクに記載されている __array_priority__ プロパティを設定すればよい。DataFrame 側の優先度を上げると期待した処理が行われることがわかる。

pd.DataFrame.__array_priority__ = 1000

data >> pd.DataFrame()
#    0  1
# 0  1  2
# 1  3  4

まとめ

R 使ったほうがいい。Use R!

補足 こんなパッケージもあります。

2016/10/18 追記

続き(?)です。

sinhrks.hatenablog.com