Ich arbeite derzeit an einem Algorithmus zur Implementierung eines rollierenden Medianfilters (analog zu einem rollierenden Mittelwertfilter) in C. Aus meiner Literaturrecherche geht hervor, dass es zwei einigermaßen effiziente Möglichkeiten gibt, dies zu tun. Das erste besteht darin, das anfängliche Wertefenster zu sortieren und dann eine binäre Suche durchzuführen, um den neuen Wert einzufügen und den vorhandenen bei jeder Iteration zu entfernen.
Die zweite (von Hardle und Steiger, 1995, JRSS-C, Algorithmus 296) baut eine doppelendige Heap-Struktur auf, mit einem Maxheap an einem Ende, einem Minheap am anderen und dem Median in der Mitte. Dies ergibt einen linearen Zeitalgorithmus anstelle eines Algorithmus, der O (n log n) ist.
Hier ist mein Problem: Die Implementierung des ersteren ist machbar, aber ich muss dies auf Millionen von Zeitreihen ausführen, daher ist Effizienz sehr wichtig. Letzteres erweist sich als sehr schwierig umzusetzen. Ich habe Code in der Datei Trunmed.c des Codes für das Statistikpaket von R gefunden, aber er ist ziemlich nicht zu entziffern.
Kennt jemand eine gut geschriebene C-Implementierung für den linearen zeitlich rollierenden Medianalgorithmus?
Bearbeiten: Link zum Trunmed.c-Code http://google.com/codesearch/p?hl=de&sa=N&cd=1&ct=rc#mYw3h_Lb_e0/R-2.2.0/src/library/stats/src/Trunmed.c
Antworten:
Ich habe mir R's
src/library/stats/src/Trunmed.c
ein paar Mal angesehen, da ich auch etwas Ähnliches in einer eigenständigen C ++ - Klasse / C-Subroutine haben wollte. Beachten Sie, dass dies tatsächlich zwei Implementierungen in einer sind (siehesrc/library/stats/man/runmed.Rd
Quelle der Hilfedatei)Es wäre schön zu sehen, dass dies eigenständiger wiederverwendet wird. Machst du Freiwilligenarbeit? Ich kann mit einigen der R-Bits helfen.
Edit 1 : Neben dem Link zur älteren Version von Trunmed.c oben sind hier aktuelle SVN-Kopien von
Srunmed.c
(für die Stuetzle-Version)Trunmed.c
(für die Turlach-Version)runmed.R
für die R-Funktion, die diese aufruftEdit 2 : Ryan Tibshirani hat einen C- und Fortran-Code für schnelles Median-Binning, der ein geeigneter Ausgangspunkt für einen Ansatz mit Fenster sein kann.
quelle
Ich konnte keine moderne Implementierung einer C ++ - Datenstruktur mit Ordnungsstatistik finden, sodass beide Ideen in dem von MAK vorgeschlagenen Link für Top-Codierer implementiert wurden ( Match Editorial : Scrollen Sie nach unten zu FloatingMedian).
Zwei Multisets
Die erste Idee unterteilt die Daten in zwei Datenstrukturen (Heaps, Multisets usw.) mit O (ln N) pro Einfügen / Löschen, sodass das Quantil nicht ohne große Kosten dynamisch geändert werden kann. Das heißt, wir können einen rollierenden Median oder einen rollierenden 75% haben, aber nicht beide gleichzeitig.
Segmentbaum
Die zweite Idee verwendet einen Segmentbaum, der O (ln N) für Einfügen / Löschen / Abfragen ist, aber flexibler ist. Das Beste von allem ist, dass "N" die Größe Ihres Datenbereichs ist. Wenn Ihr rollierender Median ein Fenster von einer Million Elementen hat, Ihre Daten jedoch von 1..65536 abweichen, sind nur 16 Operationen pro Bewegung des rollenden Fensters von 1 Million erforderlich !!
Der c ++ - Code ähnelt dem, was Denis oben gepostet hat ("Hier ist ein einfacher Algorithmus für quantisierte Daten").
Statistische Bäume der GNU-Ordnung
Kurz bevor ich aufgab, stellte ich fest, dass stdlibc ++ Ordnungsstatistikbäume enthält !!!
Diese haben zwei kritische Operationen:
Siehe libstdc ++ Handbuch policy_based_data_structures_test (Suche nach "split and join").
Ich habe den Baum zur Verwendung in einen praktischen Header für Compiler eingeschlossen, die partielle Typedefs im Stil von c ++ 0x / c ++ 11 unterstützen:
quelle
Ich habe eine getan C - Implementierung hier . Einige weitere Details finden Sie in dieser Frage: Rollender Median in der C-Turlach-Implementierung .
Beispielnutzung:
quelle
Ich verwende diesen inkrementellen Medianschätzer:
welches die gleiche Form hat wie der allgemeinere Mittelwertschätzer:
Hier ist eta ein kleiner Lernratenparameter (z. B.
0.001
) undsgn()
die Signumfunktion, die einen von zurückgibt{-1, 0, 1}
. (Verwenden Sie eine Konstanteeta
wie diese, wenn die Daten nicht stationär sind und Sie Änderungen im Laufe der Zeit verfolgen möchten. Andernfalls verwenden Sie für stationäre Quellen eine Arteta = 1 / n
Konvergenz, bei dern
die Anzahl der bisher gesehenen Proben angegeben ist.)Außerdem habe ich den Medianschätzer so geändert, dass er für beliebige Quantile funktioniert. Im Allgemeinen gibt eine Quantilfunktion den Wert an, der die Daten in zwei Brüche unterteilt:
p
und1 - p
. Im Folgenden wird dieser Wert schrittweise geschätzt:Der Wert
p
sollte innerhalb liegen[0, 1]
. Dies verschiebt im Wesentlichen diesgn()
symmetrische Ausgabe der Funktion,{-1, 0, 1}
um sich zu einer Seite zu neigen, wobei die Datenproben in zwei ungleich große Bins aufgeteilt werden (Brüchep
und1 - p
Daten sind kleiner als / größer als die Quantilschätzung). Beachten Sie, dass sichp = 0.5
dies auf den Medianschätzer reduziert.quelle
Hier ist ein einfacher Algorithmus für quantisierte Daten (Monate später):
quelle
Der rollierende Median kann ermittelt werden, indem zwei Partitionen von Zahlen beibehalten werden.
Verwenden Sie zum Verwalten von Partitionen Min Heap und Max Heap.
Max Heap enthält Zahlen, die kleiner als der Median sind.
Min Heap enthält Zahlen, die größer als der Median sind.
Ausgleichsbeschränkung: Wenn die Gesamtzahl der Elemente gerade ist, sollten beide Heaps gleiche Elemente haben.
Wenn die Gesamtzahl der Elemente ungerade ist, hat Max Heap ein Element mehr als Min Heap.
Medianelement: Wenn beide Partitionen die gleiche Anzahl von Elementen haben, ist der Median die Hälfte der Summe aus dem maximalen Element der ersten Partition und dem minimalen Element der zweiten Partition.
Andernfalls ist der Median das maximale Element der ersten Partition.
quelle
Es ist vielleicht erwähnenswert, dass es einen Sonderfall gibt, der eine einfache exakte Lösung hat: Wenn alle Werte im Stream Ganzzahlen innerhalb eines (relativ) kleinen definierten Bereichs sind. Angenommen, sie müssen alle zwischen 0 und 1023 liegen. In diesem Fall definieren Sie einfach ein Array mit 1024 Elementen und eine Anzahl und löschen alle diese Werte. Inkrementieren Sie für jeden Wert im Stream den entsprechenden Bin und die Anzahl. Nachdem der Stream beendet ist, suchen Sie den Behälter, der den höchsten Wert für count / 2 enthält. Dies kann leicht durch Hinzufügen aufeinanderfolgender Behälter ab 0 erreicht werden. Mit derselben Methode kann der Wert einer beliebigen Rangfolge ermittelt werden. (Es gibt eine geringfügige Komplikation, wenn das Erkennen der Behälter-Sättigung und das "Aufrüsten" der Größe der Lagerplätze auf einen größeren Typ während eines Laufs erforderlich ist.)
Dieser Sonderfall mag künstlich erscheinen, ist aber in der Praxis sehr häufig. Es kann auch als Annäherung für reelle Zahlen verwendet werden, wenn sie innerhalb eines Bereichs liegen und eine "gut genug" Genauigkeit bekannt ist. Dies würde für so ziemlich jede Reihe von Messungen an einer Gruppe von "realen" Objekten gelten. Zum Beispiel die Höhen oder Gewichte einer Gruppe von Menschen. Nicht groß genug? Es würde genauso gut für die Längen oder Gewichte aller (einzelnen) Bakterien auf dem Planeten funktionieren - vorausgesetzt, jemand könnte die Daten liefern!
Es sieht so aus, als hätte ich das Original falsch verstanden - es scheint, als würde es einen Schiebefenster-Median anstelle des Medians eines sehr langen Streams wollen. Dieser Ansatz funktioniert immer noch dafür. Laden Sie die ersten N Stream-Werte für das Anfangsfenster, und erhöhen Sie dann für den N + 1. Stream-Wert den entsprechenden Bin, während Sie den Bin entsprechend dem 0. Stream-Wert dekrementieren. In diesem Fall müssen die letzten N Werte beibehalten werden, um die Dekrementierung zu ermöglichen. Dies kann effizient erfolgen, indem ein Array der Größe N zyklisch adressiert wird. Da sich die Position des Medians nur um -2, -1,0,1 ändern kann Bei 2 in jedem Schritt des Schiebefensters müssen nicht alle Bins bis zum Median in jedem Schritt summiert werden. Passen Sie einfach den "Medianzeiger" an, je nachdem, welche Seitenfächer geändert wurden. Zum Beispiel, Wenn sowohl der neue als auch der entfernte Wert unter den aktuellen Median fallen, ändert sich dieser nicht (Offset = 0). Die Methode bricht zusammen, wenn N zu groß wird, um bequem im Speicher gehalten zu werden.
quelle
Wenn Sie in der Lage sind, Werte als Funktion von Zeitpunkten zu referenzieren, können Sie Werte durch Ersetzen abtasten und Bootstrapping anwenden , um einen Bootstrap-Medianwert innerhalb von Konfidenzintervallen zu generieren. Auf diese Weise können Sie einen ungefähren Median mit größerer Effizienz berechnen, als eingehende Werte ständig in eine Datenstruktur zu sortieren.
quelle
Für diejenigen, die einen laufenden Median in Java benötigen ... PriorityQueue ist Ihr Freund. O (log N) einfügen, O (1) aktueller Median und O (N) entfernen. Wenn Sie die Verteilung Ihrer Daten kennen, können Sie viel besser als dies tun.
quelle
}), higher = new PriorityQueue<Integer>();
odernew PriorityQueue<Integer>(10,
. Ich konnte den Code nicht ausführen.Hier ist eine, die verwendet werden kann, wenn die genaue Ausgabe nicht wichtig ist (für Anzeigezwecke usw.). Sie benötigen Totalcount und Lastmedian sowie den neuen Wert.
Erzeugt ziemlich genaue Ergebnisse für Dinge wie page_display_time.
Regeln: Der Eingabestream muss in der Reihenfolge der Seitenanzeigezeit glatt sein, eine große Anzahl (> 30 usw.) aufweisen und einen Median ungleich Null haben.
Beispiel: Ladezeit der Seite, 800 Elemente, 10 ms ... 3000 ms, Durchschnitt 90 ms, realer Median: 11 ms
Nach 30 Eingaben beträgt der Medianfehler im Allgemeinen <= 20% (9 ms..12 ms) und wird immer geringer. Nach 800 Eingaben beträgt der Fehler + -2%.
Ein anderer Denker mit einer ähnlichen Lösung ist hier: Median Filter Super effiziente Implementierung
quelle
Hier ist die Java-Implementierung
quelle
Wenn Sie nur einen geglätteten Durchschnitt benötigen, können Sie schnell / einfach den neuesten Wert mit x und den Durchschnittswert mit (1-x) multiplizieren und dann addieren. Dies wird dann der neue Durchschnitt.
Bearbeiten: Nicht das, wonach der Benutzer gefragt hat und nicht so statistisch gültig, aber gut genug für viele Zwecke.
Ich werde es hier (trotz der Abstimmungen) für die Suche lassen!
quelle