[SciPy-User] Moving median with O(log n) updates

Wes McKinney wesmckinn@gmail....
Thu Feb 25 15:19:19 CST 2010

```At Pycon this past weekend Raymond Hettinger gave a nice talk about
leveraging Python's built-in data structures to implement various
algorithms in a clear and elegant way.

http://us.pycon.org/2010/conference/schedule/event/86/
video: http://pycon.blip.tv/file/3261295/

He posted a full implementation of a the moving median using a

http://code.activestate.com/recipes/576930-efficient-running-median-using-an-indexable-skipli/

Some time ago for pandas I had implemented moving median is a fairly
naive way (O(n * window)) using the linear-time wirth method
(quickselect would provide equivalent performance). I decided to
translate the implementation above into Cython and tailor it for
ndarray input and float data. In most cases it achieves 5-10x speedup
over the pure Python O(n * log(window)) implementation, and clearly
trouncing other implementations. Here it is:

(the actual rolling_median is in moments.pyx)

and some timings:

In [3]: import pandas

In [4]: arr = randn(10000)

In [5]: %timeit pandas.rolling_median(arr, 1000)
1 loops, best of 3: 401 ms per loop

of course this is roughly linear in the window size:

In [11]: %time pandas.rolling_median(arr, 5000)
CPU times: user 1.07 s, sys: 0.00 s, total: 1.07 s
Wall time: 1.06 s

new version (not yet ready for primetime with NaN-handling etc.)

In [6]: import pandas.lib.tseries as ts

In [7]: %timeit ts.rolling_median(arr, 1000)
10 loops, best of 3: 120 ms per loop

Larger window size in this case has minimal effect:

In [13]: %timeit ts.rolling_median(arr, 5000)
10 loops, best of 3: 120 ms per loop

I'm pretty sure scikits.timeseries is doing something O(n * window^2)

In [8]: import scikits.timeseries.lib as l

In [9]: %time l.mov_median(arr, 1000)
CPU times: user 86.48 s, sys: 0.00 s, total: 86.48 s
Wall time: 86.48 s

I had a fun time with Cython making the implementation fairly good,
but I have the feeling that there's too much Python list access (which
is faster than accessing an object ndarray using the buffer interface
I found). If someone with better C / Cython skills wants to take a
crack at speeding it up even more please be my guest.

- Wes
```