PandasのSeriesの時系列に対して、rollingとapplyを使って数秒ごとのヒストグラムを書きたい

Question

前提・実現したいこと

表題の通りです。時系列データに対して、数秒ごとのヒストグラムを求めようとしているのですが、rolling.applyの戻り値がlistやarray型ではエラーが出るため、うまくできていません。別の方法でも良いので、上記の目的を達成する方法がないでしょうか？

import numpy as np
import pandas as pd

df = pd.DataFrame({"col": np.arange(100)})

def calc_hist_range2(x):
    last_x = x[-1]
    # 最後の値から-10~+10までの範囲のヒストグラムを作成してreturn
    range_ = (last_x - 10, last_x + 10)
    return np.histogram(x, bins=7, range=range_)[0]

df["col"].rolling(3).apply(calc_hist_range2, raw=True)

発生している問題・エラーメッセージ

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-212-acc56cfb208f> in <module>
----> 1 df["col"].rolling(3).apply(calc_hist_range2, raw=True)

~/conda/lib/python3.8/site-packages/pandas/core/window/rolling.py in apply(self, func, raw, engine, engine_kwargs, args, kwargs)
   2057         self, func, raw=False, engine=None, engine_kwargs=None, args=None, kwargs=None,
   2058     ):
-> 2059         return super().apply(
   2060             func,
   2061             raw=raw,

~/conda/lib/python3.8/site-packages/pandas/core/window/rolling.py in apply(self, func, raw, engine, engine_kwargs, args, kwargs)
   1386 
   1387         # name=func & raw=raw for WindowGroupByMixin._apply
-> 1388         return self._apply(
   1389             apply_func,
   1390             center=center,

~/conda/lib/python3.8/site-packages/pandas/core/window/rolling.py in _apply(self, func, center, require_min_periods, floor, is_weighted, name, use_numba_cache, **kwargs)
    586                     result = np.apply_along_axis(calc, self.axis, values)
    587                 else:
--> 588                     result = calc(values)
    589                     result = np.asarray(result)
    590 

~/conda/lib/python3.8/site-packages/pandas/core/window/rolling.py in calc(x)
    574                         closed=self.closed,
    575                     )
--> 576                     return func(x, start, end, min_periods)
    577 
    578             else:

~/conda/lib/python3.8/site-packages/pandas/core/window/rolling.py in apply_func(values, begin, end, min_periods, raw)
   1413             if not raw:
   1414                 values = Series(values, index=self.obj.index)
-> 1415             return window_func(values, begin, end, min_periods)
   1416 
   1417         return apply_func

pandas/_libs/window/aggregations.pyx in pandas._libs.window.aggregations.roll_generic_fixed()

TypeError: only size-1 arrays can be converted to Python scalars

Answer 1

愚直に書くこと自体は難しくはないのでは？

import numpy as np
import pandas as pd

df = pd.DataFrame({"col": np.arange(100)})

def calc_hist_range2(x):
    last_x = x[-1]
    range_ = (last_x - 10, last_x + 10)
    return np.histogram(x, bins=7, range=range_)[0]

data = df["col"].values  # pandasとか要らないよねっ
window = 3
result = []
for i in range(len(data) - window + 1):
    result.append(calc_hist_range2(data[i:i+window]))

print(result)

＃ np.histogramの返す情報はぜんぶ取らないと描画するときとかに困る可能性が大

Comments

[gymel]

ありがとうございます！
生コードをはさむとNaN値処理などのバグで行がずれてしまう可能性を回避するためのテストなどが増えてしまうというコード管理の点と、また行の数が多いときに生コードだと時間がかかりすぎているのでPandasのバックエンドを使って早くしてくれるかなという速度面からPandasでかけないかと質問させていただきました。（自分がPandasに慣れていないだけで上級者であれば良い方法をご存知かと）
もう一週間ほど待って他の回答が来なかった場合にはベストアンサーとさせていただきます！

[hayataka2049]

行ずれとかは、頑張ってください……
速度面に関しては、これでものすごく遅いということはないかと思われます。実データだとどれくらいの大きさになりますか？

[gymel]

100万行のオーダーです。ノートPCで演算しているので生コードだと１回１分程度掛かっちゃっています。
> CPU times: user 1min 53s, sys: 704 ms, total: 1min 54s
> Wall time: 1min 56s
またwindowのサイズを10, 100, 1000と範囲を変えて分布を計算し、別の時系列データに対しても同様の計算を行うため、計３０分くらい掛かっています。

[hayataka2049]

マルチコアを活かせば1/コア数くらいには時間短縮できるかと思われます。
そこから先はかなり大変で、numpy配列のスライス取るのが遅いとかnumpy.histogramが遅いとかです。pandasで処理しても、回答のような直接書いたコードより複雑なことをやってくれるのは確実なので、遅いと思われます。

[gymel]

確かに時間を計測するとnumpy.histogramの遅さがネックですね。。。
ありがとうございます！勉強になりました。

[hayataka2049]

頑張るのであれば、
・結果出力用の配列をあらかじめまとめて固定サイズで確保（100行でbins=7なら100*7とか）。計算のたびにnumpy配列を生成するよりだいぶマシになります。
・その上で、ヒストグラム作成部分のコードをnumbaかcythonあたりで実装する。numba版はexampleが上げられていて、見てわかる通りものすごく複雑ということはないです。これをコピペしてrangeの機能をつけて予めある配列に出力する形にしてついでに型指定付ければいける。
https://numba.pydata.org/numba-examples/examples/density_estimation/histogram/results.html
・numbaでもcythonでもスレッド並列化が効くので、それを使って全コア使い切る。

ですかねえ。テストとか面倒だと思うので、おすすめはしないです。

Answer 2

自己解決しました。
groupbyとrollingのapplyが異なるのは直感的にいまいちですね。。。他に良い方法があれば教えていただきたいです。

df = pd.DataFrame({"col": np.arange(100)})

def calc_hist_range2(x):
    last_x = x[-1]
    range_ = (last_x - 10, last_x + 10)
    hist = np.histogram(x, bins=7, range=range_)
    return pd.Series({k: v for k, v in enumerate(hist[0])})

window = 5
df["col"].groupby(np.arange(df.shape[0]) // window).apply(list).apply(calc_hist_range2)

Comments

[gymel]

すみません。これではwindow置きに計算しているだけでダメでした。