2012-06-26 14 views
5

Die Verwendung von Rprof zeigte, dass der dct im dtt-Paket der Haupttäter in einem Stück R-Code war, der ziemlich langsam lief. Wenn ich es im stats-Paket für fft austausche (was nicht die gleiche Transformation ist, aber die gleiche Zeit für die Berechnung benötigt), hat sich meine Laufzeit dramatisch verbessert. In der Tat waren 2/3 meiner Rprof-Leitungen vorher dct-Anrufe und nur 3 Leitungen von ungefähr 600 waren fft-Anrufe, nachdem sie den Schalter gewechselt hatten.Wie führt man eine * Fast * DCT (diskrete Kosinustransformation) in R durch?

Wird die dct-Implementierung im dtt-Paket nicht mit der schnellen diskreten Fourier-Transformation durchgeführt? Wenn ja, gibt es ein Paket, das eines hat? (Ich weiß, man könnte die Daten verdoppeln und dann Koeffizienten für den dct aus diesen fft-Koeffizienten extrahieren, aber ein geradliniger schneller dct wäre sicherlich netter, und da müsste eigentlich irgendwo stehen).

+0

Versuchen Sie, das Negative zu vermeiden (und so scheint, der Anwalt des Teufels/Kämpfer-Poster zu sein) in Titeln ;-) –

+0

'Bibliothek (sos); findFn ("schnelle Fourier-Cosinus-Transformation") "offenbart eine" rfft "-Funktion in dem" wavethresh ". Es sieht so aus, als ob es durch einen Pack-/Entpackmechanismus funktioniert - basierend auf einer rostigen Erinnerung an * Numerical Recipes * könnte dies wahrscheinlich noch effizienter gemacht werden, aber vielleicht ist es nützlich? –

+0

Ich verstehe, welche Motivation es für dct sein könnte, nicht das schnelle dct zu verwenden. Einer, der schnelle dct ist viel, viel schwerer zu schreiben. Auch ein Vektor der Größe N hängt von der Primfaktor - Faktorisierung von N ab. Einige dct - Versionen packen also den Vektor auf die nächsthöhere Potenz von 2 für die Performance, aber dann sind die zurückgegebenen Koeffizienten nicht wirklich "die "DCT-Koeffizienten dieses Vektors, wie sie von der Auffüllung betroffen sind. Zur gleichen Zeit, in Situationen, in denen N auf eine Potenz von 2 gesetzt werden kann, führt die schnelle Version zu einer wahnsinnigen Leistungssteigerung. –

Antwort

7

Es scheint keine schnelle dct zu geben, aber es gibt eine fft (schnelle Fourier-Transformation) im stats-Paket, also hier ist, wie Sie die schnelle dct mit fft bekommen könnte.

Verwenden Sie dies auf eigene Gefahr. Ich habe es nicht ernsthaft überprüft. Ich überprüfte es auf ein paar Vektoren verschiedener Größen und es ergab die gleichen Ergebnisse wie Funktion dct im Paket dtt auf denen. Wenn jemand mich überprüfen möchte, indem Sie es mit der Ausgabe von dct vergleichen, dann zögern Sie nicht, dies zu tun und Ihre Ergebnisse zu posten.

Nehmen Sie Ihren Vektor und erweitern Sie ihn zu einem Vektor, der doppelt so lang ist wie folgt: Wenn Ihr Vektor v = (1,2,3) ist, dann verdoppeln Sie die Einträge zu w = (1,2,3,3,2, 1). Beachten Sie die Bestellung. Wenn Ihr Vektor v = (1,2,4,9) ist, dann verdoppeln Sie die Einträge zu w = (1,2,4,9,9,4,2,1)

Sei N die Länge von Ihr ORIGINAL-Vektor (bevor Sie seine Länge verdoppelt haben).

dann werden die ersten N Koeffizienten der .5 * FFT (w)/exp (Komplex (imaginär = pi/2/N) * (seq (2 * N) -1)) sollte mit computing dct einverstanden (v) , außer es sollte in fast allen Fällen dramatisch schneller sein.

Überlegungen zur Geschwindigkeit. Wenn Sie Faktor N initiieren, dann ist die Zeit, die es benötigt, um den schnellen dct zu berechnen, wie die Zeit, die benötigt wird, um einen langsamen dct für jeden dieser Primfaktoren durchzuführen. Also, wenn N 2^K ist, ist es so, als würde man einen K verschiedenen langsamen dct-Transformationen auf einem Vektor der Länge zwei machen, so dass sein wirklich schnell ist. Wenn N prim ist (Worst-Case-Szenario), dann gibt es überhaupt keine Beschleunigung. Die größte Beschleunigung liegt bei Vektoren, die eine Zweierpotenz haben.

Hinweis: Der obige R-Code sieht unglaublich unfreundlich aus, also lassen Sie mich sagen, was vor sich geht. Nach der doppelten Länge auf die richtige Weise sind die ersten N Koeffizienten des fft, die Sie bekommen, fast das Richtige. Die Koeffizienten müssen jedoch etwas optimiert werden. Sei P für e^(pi * i/2/N). Lassen Sie den ersten Koeffizienten in Ruhe. Teile den zweiten Koeffizienten durch P, dividiere den dritten durch P^2, dividiere den vierten durch P^3, etc ... Teile dann alle die Koeffizienten durch 2 (einschließlich des ersten), um mit der Normalisierung übereinzustimmen, für die R verwendet wird der dct.

Diese sollte geben das gleiche wie die Verwendung der dct-Funktion in Paket dtt, aber in fast allen Fällen dramatisch schneller sein.