現在、ローリングメディアンフィルター(ローリング平均フィルターに類似)をC言語で実装するアルゴリズムに取り組んでいる。1つ目は、値の初期ウィンドウをソートし、各反復で新しい値を挿入し、既存の値を削除するためにバイナリサーチを実行することです。
もう1つは(Hardle and Steiger, 1995, JRSS-C, Algorithm 296より)ダブルエンドのヒープ構造を構築し、一方の端にマックスヒープ、もう一方の端にミニヒープ、中央に中央値を置く方法である。これにより、O(n log n)ではなく、線形時間のアルゴリズムが得られる。
前者を実装することは可能だが、何百万もの時系列に対してこれを実行する必要があるため、効率が非常に重要になる。後者を実装するのは非常に難しい。RのstatsパッケージのコードのTrunmed.cファイルにコードを見つけたが、解読不能である。
どなたか、線形時間ローリングメディアンアルゴリズムのよく書かれたC実装を知りませんか?
編集:Trunmed.cコードへのリンク http://google.com/codesearch/p?hl=en&sa=N&cd=1&ct=rc#mYw3h_Lb_e0/R-2.2.0/src/library/stats/src/Trunmed.c
R'のsrc/library/stats/src/Trunmed.c
を何度か見ましたが、スタンドアロンのC++クラス/Cサブルーチンでも同じようなものが欲しかったので。src/library/stats/man/runmed.Rd`(ヘルプファイルのソース)には次のように書かれている。
\details{
Apart from the end values, the result \code{y = runmed(x, k)} simply has
\code{y[j] = median(x[(j-k2):(j+k2)])} (k = 2*k2+1), computed very
efficiently.
The two algorithms are internally entirely different:
\describe{
\item{"Turlach"}{is the Härdle-Steiger
algorithm (see Ref.) as implemented by Berwin Turlach.
A tree algorithm is used, ensuring performance \eqn{O(n \log
k)}{O(n * log(k))} where \code{n <- length(x)} which is
asymptotically optimal.}
\item{"Stuetzle"}{is the (older) Stuetzle-Friedman implementation
which makes use of median \emph{updating} when one observation
enters and one leaves the smoothing window. While this performs as
\eqn{O(n \times k)}{O(n * k)} which is slower asymptotically, it is
considerably faster for small \eqn{k} or \eqn{n}.}
}
}
これをもっと独立した形で再利用してもらえるとうれしいです。ボランティアですか? 私はRの部分を手伝うことができます。
*編集1上記の古いバージョンのTrunmed.cへのリンクの他に、現在のSVNコピーは以下の通りです。
編集2:Ryan Tibshiraniがfast median binningに関するCとFortranのコードを持っている。
平滑化された平均値が必要なだけなら、最新の値にxを掛け、平均値に(1-x)を掛け、それらを足すのが手っ取り早く簡単な方法である。これが新しい平均値となる。
edit:ユーザーが求めたものではなく、統計的に有効でもありませんが、多くの用途には十分です。
しかし、多くの用途には十分である。検索用に(ダウンボートに関わらず)ここに残しておこう!