Der effizienteste Weg, um eine lange Zeitreihe Python zu filtern
Ich habe eine große Zeitreihe, sagen wir 1e10, die sich aus der Aufzeichnung neuronaler Aktivität, d. h. Spannungen, ergibt. Bevor ich weitere Analysen mache, möchte ich diese Daten zwischen 300 Hz und 7000 Hz filtern. Im Folgenden poste ich den Code für den von mir entworfenen Butterworth-Filter. Wie mache ich diesen Filter schneller? Es dauert zu lange, um ausgeführt zu werden.
Aktualisieren: Beispieldaten
Hier ist ein Link zu Daten wie in jeder Zeile von data .
Bei der Formatierung steht jede Zeile für eine andere Aufnahmequelle und jede Spalte für einen Zeitpunkt. Die Daten werden mit 20.000 Hz abgetastet.
%Vor%Aktualisieren
Wie die meisten NumPy können auch SciPy-Funktionen lfilter eine mehrdimensionale Eingabe annehmen und so erzeugt map unnötigen Overhead. Das heißt, man kann
als
%Vor% Standardmäßig arbeitet lfilter auf der letzten Nicht-Singleton-Achse. Für eine 2D-Matrix bedeutet dies, dass die Funktion auf jede Zeile angewendet wird, was genau das ist, was ich brauche. Leider ist es immer noch zu langsam.
1 Antwort
Erstens, Ihr Datensample ist in einem proprietären Format, habe ich recht? Selbst mit der Biosig-Toolbox für Python kann dieses Format nicht gelesen werden. Vielleicht irre ich mich, aber es ist mir nicht gelungen, es zu lesen.
Also werde ich meine Antwort auf künstliche Daten stützen, die von einem Rössler-Oszillator erzeugt werden. Es ist ein chaotischer 3D-Oszillator, der oft im Bereich der nichtlinearen Zeitreihenanalyse verwendet wird.
%Vor%Kommen Sie nun zu Ihren Funktionsdefinitionen:
%Vor%Ich habe sie unverändert gelassen. Ich erzeuge einige Daten mit meinem Oszillator, aber ich nehme nur die dritte Komponente davon. Ich füge etwas Gaußsches Rauschen hinzu, um etwas herausfiltern zu können.
%Vor% Kommen wir nun zu der Geschwindigkeitsfrage. Ich verwende das timeit -Modul, um die Ausführungszeiten zu überprüfen. Diese Anweisungen führen die Filterung 100 Mal aus und messen die Gesamtzeit. Ich mache diese Messung für Ordnung 2 und Ordnung 8 (ja, Sie wollen schärfere spektrale Kanten, ich weiß, aber warten Sie)
Die Ausgabe dieser beiden print-Anweisungen lautet:
%Vor%Das macht einen Faktor von 20! Die Verwendung von Auftrag 8 für einen Butterworth-Filter ist definitiv keine gute Idee. Ich kann mir keine Situation vorstellen, in der das Sinn machen würde. Um die anderen Probleme zu untersuchen, die bei der Verwendung eines solchen Filters auftreten, betrachten wir die Wirkung dieser Filter. Wir wenden Bandpassfilter auf unsere Daten an, einmal mit der Reihenfolge 8 und einmal mit der Reihenfolge 2:
%Vor%Lasst uns jetzt einige Plots machen. Zuerst einfache Linien (mir war die x-Achse egal)
%Vor%Das gibt uns:
Oh, was ist mit der grünen Linie passiert? Seltsam, nicht wahr? Der Grund ist, dass Butterworth-Filter der Ordnung 8 zu einer ziemlich instabilen Sache werden. Schon mal von Resonanz-Katastrophe / Katastrophe gehört? So sieht es aus.
Die spektralen Leistungsdichten dieser Signale können wie folgt dargestellt werden:
%Vor% 
Hier sehen Sie, dass die grüne Linie schärfere Kanten hat, aber um welchen Preis? Künstliche Spitze bei ca. 300 Hz. Das Signal ist vollständig verzerrt.
Was sollen Sie tun?
- Niemals Verwenden Sie einen Butterworth-Filter der Ordnung 8
- Verwenden Sie niedrigere Reihenfolge, wenn es ausreicht.
- Wenn nicht, erstellen Sie einen FIR-Filter mit den Parks-McGlellan- oder Remez-Exchange-Algorithmen. Es gibt zum Beispiel scipy.signal.remez.
Noch ein Hinweis: Wenn Sie sich für die Phase Ihres Signals interessieren, sollten Sie vor und zurück in der Zeit filtern . lfilter verschiebt die Phasen andernfalls. Eine Implementierung eines solchen Algorithmus, der gemeinhin als filtfilt bezeichnet wird, finden Sie unter my Github-Repository .
Noch ein Programmiertipp:
Wenn Sie die Situation von Pass-Through-Parametern haben (vier Parameter von Dies reduziert die Code-Redundanz und macht Ihren Code weniger fehleranfällig. Schließlich ist hier das vollständige Skript, um es einfach kopieren zu können, um es auszuprobieren. butter_bandpass_filter werden nur an butter_bandpass übergeben, können Sie *args und **kwargs .
Tags und Links python python-2.7 numpy scientific-computing signal-processing