双スタックスライディングウィンドウ集約

集約とは、あるデータ集合を何らかの形でまとめることだ。総和、要素数、最小値などがそれにあたる。この集約を繰り返し計算したい場面はよくある。たとえばノイズ検出器で「直近30秒間の最大ノイズレベル」を求めたい場合だ。このときデータ上にスライディングウィンドウがあるといい、そのウィンドウ内で集約を行いたい。

集約操作が逆元を持つ二項演算子なら、たとえば整数の加算なら、両端キューを使えば驚くほど簡単に解決できる:

from collections import deque

class SlidingWindowSum:
 def __init__(self):
 self.sum = 0
 self.elems = deque()

 def push(self, x):
 self.sum += x
 self.elems.append(x)

 def pop(self):
 self.sum -= self.elems.popleft()

 def eval(self):
 return self.sum

しかし演算子に逆元がなければどうだろうか? 面白い集約のほとんどがこれに該当する。最小値、分位数、近似カーディナリティ推定(HyperLogLogなど)がそうだ。実のところ、浮動小数点数の総和という単純な操作でさえこの問題を抱えている。浮動小数点の加算は可逆ではないからだ。たとえば入力データに一度でも NaN が現れると、上の素朴なアルゴリズムが返す和は NaN のまま動かなくなる。その値がとっくにウィンドウから外れていても、だ。

6年前、スライディングウィンドウにおける最小値・最大値の維持についてアルゴリズムを考え、cs.stackexchange に投稿した。今ではあのアルゴリズムに意味はないと思っている。幅広い種類の集約に対してこの問題を解く、シンプルで効率的なアルゴリズムが存在するからだ。この記事を書くのは、それを多くの人に知ってもらいたいからだ。もっと広く知られていい手法だと思う。

民間伝承

このアルゴリズムに出会ったのは、より高度な論文を読んでいたときだった。Tangwongsan らによる Low-Latency Sliding-Window Aggregation in Worst-Case Constant Time だ。なぜ「低遅延」なのか。やっていることは以下で説明するものと同じだが、各ステップが O(1) だからだ。とはいえ、その論文は「双スタック」アルゴリズムについても述べている。こちらは償却 O(1) で、しかも比較にならないほどシンプルだ。

面白いのは、論文がこのアルゴリズムを 2011年の Stack Overflow の投稿における「adamax」の功績としている点だ。さらにその投稿は D. Sleator の2001年の講義ノートに触発されたと述べている。だがその講義ノートが扱っているのはスライディングウィンドウ集約ではなく、FIFO キューを実装する古典的な双スタックアルゴリズムとその償却解析だ。結局のところ、これほどシンプルで巧妙なアルゴリズムなら、1970年代の無名の論文にすでに登場していたと分かっても驚かない。

2つのスタック

論文の著者たちと同じく、双スタックアルゴリズムを結合則を満たす任意の集約関数に一般化する。集約を empty()、unit(x)、combine(x, y)、finalize(x) という関数群として抽象化すれば、多くの集約を記述できる。たとえば平均値なら:

empty = lambda: (0, 0)
unit = lambda x: (x, 1)
combine = lambda x, y: (x[0] + y[0], x[1] + y[1])
finalize = lambda x: x[0] / x[1] if x[1] else None

ここで各関数のシグネチャはこうなっている:

fn empty() -> Agg;
fn unit(x: Value) -> Agg;
fn combine(x: Agg, y: Agg) -> Agg;
fn finalize(x: Agg) -> Out;

Value、Agg、Out をわざわざ区別しているのは、整数の総和のような場面では表面的には同じに見えても、文字列の近似カーディナリティ推定では (Value, Agg, Out) = (String, HyperLogLogSketch, u64) が得られるからだ。3つともまったく別の型である。

前置きはこのくらいにして、アルゴリズムはこうだ:

class TwoStackAgg:
 def __init__(self):
 self.values = []
 self.values_agg = empty()
 self.cum_aggs = []

 def push(self, x):
 self.values.append(x)
 self.values_agg = combine(self.values_agg, unit(x))

 def pop(self):
 if not self.cum_aggs:
 cum_agg = empty()
 while self.values:
 cum_agg = combine(unit(self.values.pop()), cum_agg)
 self.cum_aggs.append(cum_agg)
 self.values_agg = empty()
 self.cum_aggs.pop()

 def eval(self):
 return finalize(
 combine(self.cum_aggs[-1], self.values_agg)
 if self.cum_aggs else self.values_agg
 )

以上で、アルゴリズムの全体はこれだけだ。2つのスタック(values と cum_aggs)と、集約値 values_agg が1つ。任意の時点で、values_agg は values の集約結果を保持し、cum_aggs はウィンドウ内で values に含まれないすべての値の累積集約を、逆順で保持している。したがって cum_aggs の末尾の値と values_agg をマージすれば、ウィンドウ全体の集約結果が定数時間で得られる。

巧妙なのはここだ(w をウィンドウサイズとする)。w 回の操作ごとに values をすべて空にし、ローリング集約値を維持しながら部分的な累積集約を cum_aggs にプッシュする。これが償却 O(1) を実現している。w 回の pop ごとに O(w) 回の内部操作が行われ、要素あたりの総作業量は O(1) に抑えられる。個々の操作が定数時間とは限らないにもかかわらず、だ。

図で見るのが一番わかりやすいと思う。[1, 2, ..., 10] の総和を取るとしよう。スライディングウィンドウを4要素に固定すると、eval() を呼ぶたびの状態はこうなる(values_agg はこれらの値の集約結果なので省略する):

 cum_aggs values out
 [] [] = 0
 [] [1] = 1
 [] [1, 2] = 1 + 2
 [] [1, 2, 3] = 1 + 2 + 3
 [] [1, 2, 3, 4] = 1 + 2 + 3 + 4
[4, 3 + 4, 2 + 3 + 4] [5] = 2 + 3 + 4 + 5
 [4, 3 + 4] [5, 6] = 3 + 4 + 5 + 6
 [4] [5, 6, 7] = 4 + 5 + 6 + 7
 [] [5, 6, 7, 8] = 5 + 6 + 7 + 8
[8, 7 + 8, 6 + 7 + 8] [9] = 6 + 7 + 8 + 9
 [8, 7 + 8] [9, 10] = 7 + 8 + 9 + 10
 [8] [9, 10] = 8 + 9 + 10
 [] [9, 10] = 9 + 10
 [10] [] = 10
 [] [] = 0

総メモリ使用量は O(w) で、w はウィンドウの最大サイズだ。解析の都合と例の単純化のため、ウィンドウサイズを w に固定しているが、双スタックアルゴリズム自体はこれを要求しない。eval() の合間に push(x) と pop() を好きなだけ呼び、必要に応じてウィンドウを広げたり狭めたりできる。

浮動小数点は結合則を満たさない

ここまで、集約操作 combine が結合則を満たすことを前提としてきた。つまり:

combine(combine(x, y), z) = combine(x, combine(y, z))

厳密に言えば、浮動小数点の加算はこの性質を満たさない。それでも上のアルゴリズムは十分に有用だ。結果は期待値にかなり近い。補償付き総和アルゴリズム、たとえば Kahan の総和を使えば、さらに近づく。

上のアルゴリズムには、もう一つ有用な性質がある。各集約値は厳密にウィンドウ内の要素だけから構成され、ウィンドウ外の要素は一切関与しない。つまり、ウィンドウ内に NaN や無限大(あるいはその他の外れ値)が現れても、その値が汚染するのはそれを含むウィンドウだけで、残りの計算には影響しない。

NaN や無限大がなくても、この性質は役に立つ。誤差が際限なく伝播しないからだ。たとえばスライディングウィンドウが [1e20, 1] で始まると、素朴なローリング総和はこうなる:

>>> 1e20 + 1 - 1e20 - 1
-1.0

補償付き総和はこうした影響を和らげるが、結果をウィンドウ外の値に依存させないことが、長期的な誤差の蓄積を根本から断ち切る。