Sparse Matrix Implementierung des Kalman Filters?

Ich habe einen Kalman-Filter-basierten Modellierungscode, den ich für eine regionale ionosphärische Kartierungsanwendung in Echtzeit entwickelt habe. Der Code assimiliert Daten von verschiedenen Sensoren mithilfe eines Kalman-Filters in eine Karte (beschrieben durch eine Reihe von Basisfunktionen).

Ich versuche, dies auf einen größeren Bereich und mehr Sensoren zu skalieren, aber der Matrixalgebra-Teil des Kalman-Filters wird aufgrund der großen Matrizen (Tausende von Zeilen / Spalten) sehr langsam. Ich vermute, der beste Weg, um das Laufzeitproblem anzugreifen, besteht darin, die Tatsache zu nutzen, dass diese Matrizen mit 80% oder mehr der gesamten Elemente Null normalerweise sehr spärlich sind. Der Grund dafür ist, dass jeder Sensor einen Vorspannungsparameter hat, der gemeinsam mit den Kartenkoeffizienten geschätzt wird. Dies wird in der Spalte für diesen Sensor in der Kalman H-Matrix als 1 angezeigt, in den Spalten für jeden anderen Sensor- und Kartenkoeffizienten als Null. Es gibt Hunderte von Sensoren, die jeweils 8-10 Beobachtungen in jeder Epoche liefern, daher viele Nullen.

Ich könnte versuchen, die Komponenten des Kalman-Filters mit spärlichen Algorithmen zu implementieren, insbesondere mit Multiplikation und Inversion *, aber ich frage mich, ob es einen noch besseren Ansatz gibt, der den Kalman-Filter in eine andere Form umwandelt, die für Fälle mit Matrizen besser geeignet ist spärlich? Ich weiß, dass ich einen Ensemble-Kalman-Filter oder ähnliches verwenden könnte, aber wenn möglich möchte ich die Optimalität des reinen linearen Kalman-Filters beibehalten. Das gesamte Datenvolumen ist nicht unerschwinglich, nur die großen, spärlichen Matrizen, die sich aus dem linearen Modell ergeben.

In Bezug auf die Implementierung erfolgt dies in IDL, die Kernmatrixalgebra erfolgt jedoch über Aufrufe externer optimierter LA-Bibliotheken (insbesondere ATLAS).

* Ich weiß, dass eine optimale Kalman-Filterimplementierung eine Inversion vermeidet und stattdessen eine UD-Zerlegung verwendet. Ich denke darüber nach, so etwas zu implementieren, damit das die Antwort sein kann, aber ich bin gespannt, ob es angesichts der spärlichen Matrizen eine bessere Lösung gibt.

Bei dünn besetzten Matrizen ist es häufig der Fall, dass dicht ist , obwohl eine Matrix dünn ist. In solchen Fällen ist die Cholesky- oder Faktorisierung von eher spärlich (insbesondere wenn die Zeilen / Spalten von neu angeordnet werden, um das Sparsity-Muster zu verbessern). In den meisten Fällen, wenn Sie die Sparsity ausnutzen möchten und nicht daran interessiert sind Wenn Sie einen iterativen Algorithmus zum Lösen von Systemen mit , ist es besser, eine Faktorisierung der Matrix zu verwenden, als explizit zu berechnen . A - 1 L U A A A A - $A$$A^{-1}$$LU$$A$$A$$A$$A^{-1}$

Insbesondere für die Kalman-Filterung und nicht für die Datenverarbeitung

$S_{k}=(H_{k}P_{k-1,k}H_{k}^{T}+R_{k})^{-1}$

Normalerweise ist es besser, mit einer Faktorisierung von . Da symmetrisch ist und positiv definitiv sein sollte, können Sie dazu eine Cholesky-Faktorisierung oder eine -Faktorisierung verwenden. Sie haben uns gesagt, dass Ihre -Matrix dünn ist, aber Sie haben uns nichts darüber gesagt, ob dünn oder anderweitig strukturiert ist, und natürlich könnte ziemlich dicht sein. S k L D L T H k R k$S_{k}^{-1}$$S_{k}$$LDL^{T}$$H_{k}$$R_{k}$$P$

Ein Grund dafür, dass der Ensemble Kalman Filter (EnKF) und verschiedene Partikelfiltertechniken so beliebt sind, ist, dass für Systeme mit einem extrem großen Zustandsvektor die herkömmliche Kalman-Filterung sehr schwierig wird. EnKF kann für sehr große Zustandsvektoren effizient implementiert werden, wenn diagonal oder nahezu diagonal ist. Diese Fragen wurden von Personen, die auf dem Gebiet der Datenassimilation tätig sind, eingehend behandelt. Ich empfehle daher, Ihre Recherche damit zu beginnen, zu lesen, wie sie mit diesen Problemen umgegangen sind. $R_{k}$

Wir haben einen robusten Algorithmus für den Ensemble Kalman-Filter (und den regulären Kalman-Filter). Es eignet sich gut für spärliche Matrizen und parallele Berechnungen, da es auf orthogonalen Matrizen basiert und sich auf die Quadratwurzel- oder UD-Algorithmen bezieht.

Würde gerne das Papier schicken

Thomas, SJ, J. Hacker und J. Anderson, (2009): Eine robuste Formulierung des Ensembles Kalman Filter Quart J. Royal Met. Soc, Band 135, 507-521,

(PDF vom Verlag ist kostenlos .)

Vor langer Zeit hatte ich die Gelegenheit, an der Reduzierung der Dimensionalität zu arbeiten, die sich mit der Verarbeitung von Daten befasst, die in große Mengen fallen. Die Grundidee dahinter ist, dass die Daten in wenigen Schritten verarbeitet werden, um sie so auszurichten, dass die meisten Informationen daraus berechnet werden können.

Es funktioniert auch ziemlich gut für Matrizen und wird weitgehend verwendet. Das Beste daran ist, dass Sie es nicht einmal programmieren müssen, da bereits Standardbibliotheken dafür verfügbar sind. Wichtige mathematische Tools wie Matlab und Mathematica unterstützen diese Funktionalität ebenfalls unkompliziert.

Es gibt zwei Hauptalgorithmen, die dies erreichen - Hauptkomponentenanalyse und Singularwertzerlegung.

Was diese Algorithmen tatsächlich erreichen, ist das Auffinden der Daten, die Ihren Messwert tatsächlich erheblich beeinflussen. Das Internet ist voll von Informationen zu diesen Algorithmen. Dies zeigt Ihnen, wie Apache es macht.

Wir entschuldigen uns dafür, dass wir die Diskussion für eine Weile nicht erweitert haben - aber ich würde gerne die Zusammenfassung des Papiers veröffentlichen -, gehen aber auch gerne auf die Gründe ein, warum die Berechnungen anders organisiert sind als die Standard-KF von EnKF-Ansätzen