Consistent Estimator

11.2.3 Varianzreduzierung

Eine flüchtige Untersuchung der Spektralschätzungen für die Drahtmessreihen in den Abbildungen 5 und 6 zeigt eine erhebliche Variabilität über die Frequenzen hinweg, so dass es schwierig ist, die Gesamtstruktur in den Spektralschätzungen ohne eine angemessene Menge an Studien zu erkennen. Alle direkten Spektralschätzer leiden unter dieser inhärenten Choppiness, die durch Berücksichtigung der Verteilungseigenschaften von S ^ X (d) (f) erklärt werden kann. Erstens, wenn f nicht zu nahe an 0 oder f (N) ist und wenn SW (N) eine milde Regelmäßigkeitsbedingung erfüllt, dann ist 2S ^ w (d) (f) / Sw (f) = dx22; dh der rv 2S ^ w (d) (f) / Sw (f) ist in der Verteilung ungefähr gleich einem Chi-Quadrat rv mit 2 Freiheitsgraden. Wenn keine Verjüngung verwendet wird, gilt f als “nicht zu nahe” an 0 oder f (N), wenn 1 / (n-p) Δt< f< f (N) -1 / (n-p) Δt; wenn eine Verjüngung verwendet wird, müssen wir 1 / (n-p) Δt durch einen größeren Term ersetzen, der die vergrößerte Breite der Zentralkeule des Spektralfensters widerspiegelt (zum Beispiel ist der Term für die Hanning-Datenverjüngung ungefähr 2 / (n-p) Δt, so dass f “nicht zu nahe” ist, wenn 2 / (n-p) Δt < f < f (N) -2 / (n-p) Δt).

Da ein Chi-Quadrat rv xv2 mit v Freiheitsgraden eine Varianz von 2u hat, haben wir die Approximation V=Sw2(f). Dieses Ergebnis ist unabhängig von der Anzahl der Wt, wir haben: im Gegensatz zu Statistiken wie dem Stichprobenmittelwert unabhängiger und identisch verteilter Gaußscher rvs nimmt die Varianz von S ^ W (d) (f) nicht auf 0 ab, wenn die Stichprobengröße n-p größer wird (außer im uninteressanten Fall Sw(f) = 0). Dieses Ergebnis erklärt die choppiness der direkten spektralen Schätzungen in den Abbildungen 5 und 6. In der statistischen Terminologie ist S ^ W (d) (f) ein inkonsistenter Schätzer von Sw (f).

wobei Wm(⋅) ein Glättungsfenster ist, dessen Glättungseigenschaften durch den Glättungsparameter m gesteuert werden. Mit anderen Worten, der Schätzer S ^ W (lw) (⋅) wird erhalten, indem ein Glättungsfenster mit dem direkten Spektralschätzer S ^ w (d) (⋅) gefaltet wird. Ein typisches Glättungsfenster hat das gleiche Aussehen wie ein Spektralfenster. Es gibt einen zentralen Lappen mit einer Breite, die durch den Glättungsparameter m eingestellt werden kann: Je breiter dieser zentrale Lappen ist, desto glatter wird S ^ W(w) (⋅). Es kann auch eine Reihe von lästigen Seitenkeulen geben, die ein Auslaufen des Fensters verursachen. Das Vorhandensein von Glättungsfensterlecks wird leicht erkannt, indem Diagramme von S ^ W (lw) (⋅) und S ^ W (d) (⋅) überlagert und nach Frequenzbereichen gesucht werden, in denen ersteres keine geglättete Version des letzteren zu sein scheint.

Die statistischen Eigenschaften von S ^ W(lw)(.) sind wegen des folgenden großen Probenergebnisses handhabbar. Wenn S ^ W (d) (⋅) tatsächlich das Periodogramm ist (dh wir haben die Werte von Wt nicht angegeben), ist die Menge der rvs S ^ W (d) (j / (n−p) Δt), j = 1,2, …, j,, ungefähr paarweise unkorreliert, wobei jeder rv proportional zu a χ22, rv ist (hier ist J die größte ganze Zahl, so dass J / (n-p) < 1/2). Wenn wir die Verjüngung verwendet haben, um Sw (d) (⋅) zu bilden, gilt eine ähnliche Aussage für einen kleineren Satz von rvs, die auf einem gröberen Gitter mit gleich beabstandeten Frequenzen definiert sind — wenn der Grad der Verjüngung zunimmt, nimmt die Anzahl der ungefähr unkorrelierten rvs ab. Unter der Annahme, dass die sdf Sw (⋅) langsam über Frequenzen variiert (Vorweißung hilft, dies wahr zu machen) und dass die zentrale Keule des Glättungsfensters im Vergleich zu den Variationen in Sw (⋅) ausreichend klein ist, folgt daraus, dass S ^ W (d) (f) in Gl. (11.15) kann durch eine lineare Kombination von unkorrelierten χ22 rvs angenähert werden. Ein Standardargument “äquivalente Freiheitsgrade” kann dann verwendet werden, um die Verteilung von S ^ W (lw) (f) zu approximieren. (siehe Gl. (11.17) später).

wobei C^τ.W (d), ist der acvs-Schätzer in Gl. (11.9) entsprechend S ^ W(d)(.) und {wt.m} ist ein Verzögerungsfenster (dies kann als inverse Fourier-Transformation des Glättungsfensters Wm(⋅) angesehen werden). In der Tat, weil S ^ W(d)(.) ein trigonometrisches Polynom ist, können alle diskreten Windungen der Form ΣkWm(f−fk’)S^W(d)(fk’) auch über Gl. (11.16) mit einer geeigneten Auswahl von wt,m-Werten (Einzelheiten siehe Abschnitt 6.7). Unsere beiden praktischen Möglichkeiten der Berechnung S ^ W(l, w)(.) ergeben somit äquivalente Schätzer. Es sei denn, die diskrete Faltung ist ausreichend kurz, Gl. (11.16) ist rechnerisch schneller zu verwenden.

in guter Näherung,wobei v als äquivalente Freiheitsgrade für S ^ W(lw)(f) bezeichnet wird und durch v=2 (n−p)BwΔt / Ch gegeben ist. Dabei ist Bw ein Maß für die Bandbreite des Glättungsfensters Wm(⋅)a)n d kann über BW=1/Δt∑T=−(n−p−1)n−p−1wT,m2 berechnet werden;;auf der anderen Seite hängt Ch nur von der Verjüngung ab, die auf die Werte von angewendet wird Wtund kann über berechnet werden Ch=(n−p)∑1= p+ 1nht4beachten Sie, dass, wenn wir uns nicht explizit verjüngen, dann ht= 1 / n−pund daher Ch> 1; Für eine typische Datenverjüngung sagt uns die Cauchy-Ungleichung, dass Ch> 1 (zum Beispiel Ch≈1,94 für die Hanning-Datenverjüngung). Die äquivalenten Freiheitsgrade für S ^ W (lw) (f) nehmen somit zu, wenn wir die Glättungsfensterbandbreite erhöhen, und nehmen ab, wenn wir den Verjüngungsgrad erhöhen. Gleichung (11.17) sagt uns, dass E≈SW (f) und dass V≈SW2 (f) / v, also nimmt mit zunehmendem v V ab.

Die Prozentpunkte ην(α) sind in zahlreichen Lehrbüchern tabelliert oder können mit einem von Best und Roberts angegebenen Algorithmus berechnet werden

welches eine Breite hat, die unabhängig von S ^W(lw)(.). Dies ist der Grund für das Zeichnen von SDF-Schätzungen auf einer Dezibel- (oder logarithmischen) Skala.

wobei m als positive ganze Zahl angenommen wird und τ= τ / m. Dieses Verzögerungsfenster ist einfach zu berechnen und weist Seitenkeulen auf, deren Hüllkurve als f-4 abklingt, so dass Glättungsfensterlecks selten ein Problem darstellen. In guter Näherung ist die Glättungsfensterbandbreite für das Parzen-Lag-Fenster gegeben durch Bw=1,85/(mΔt). Wenn m zunimmt, nimmt die Bandbreite des Glättungsfensters ab, und der resultierende Verzögerungsfensterschätzer wird weniger glatt. Die zugehörigen äquivalenten Freiheitsgrade sind näherungsweise durch v=3,71(n-p)/(mCh) gegeben. Das Parzen-Lag-Fenster für m = 32 und das zugehörige Glättungsfenster sind in Abbildung 7 dargestellt.

Abb.7. Parzen Verzögerungsfenster (a) und das entsprechende Glättungsfenster (b) für m = 32. Die Glättungsfensterbandbreite beträgt N = 0,058.

Als Beispiel zeigt Abbildung 8 (a) einen nachgefärbten Verzögerungsfensterschätzer für die Drahtwellenmessdaten (die durchgezogene Kurve) zusammen mit dem entsprechenden nachgefärbten direkten Spektralschätzer (die Punkte, diese zeigen dieselbe Schätzung wie in Abbildung 6 (b)). Hier wurde das Parzen-Lag-Fenster mit einem Wert von m= 237 für den Glättungsfensterparameter verwendet (die entsprechenden äquivalenten Freiheitsgrade v sind 64). Dieser Wert wurde nach einigen Experimenten gewählt und scheint einen Verzögerungsfensterschätzer zu erzeugen, der alle wichtigen spektralen Merkmale erfasst, die durch den direkten Spektralschätzer für Frequenzen zwischen 0,4 und 4,0 Hz angezeigt werden (beachten Sie jedoch, dass dieser Schätzer die Spitze zwischen 0,0 und 0,4 Hz ziemlich schlecht ausschmiert). Wir haben auch ein Kreuzmuster gezeichnet, dessen vertikale Höhe die Länge eines 95% igen Konfidenzintervalls für 10 ⋅ log10 (SX (f)) (basierend auf dem postcolored lag window Estimator) und dessen horizontale Breite die Glättungsfensterbandbreite darstellt.

Abb.8. Postcolored Parzen lag Fenster Spektralschätzung-durchgezogene Kurve auf Grundstück (ein)—und WOSA Spektralschätzung-durchgezogene Kurve auf (b)—für Drahtwellenmessgerät Zeitreihe. Der Glättungsfensterparameter für das Parzen-Lag-Fenster war m = 237, was v = 64 äquivalente Freiheitsgrade ergab. Die WOSA-Spektralschätzung wurde unter Verwendung eines Hanning-Datenmodells an Blöcken mit 256 Datenpunkten gebildet, wobei sich benachbarte Blöcke um 50% überlappen. Die äquivalenten Freiheitsgrade für diese Schätzung sind v = 59.

WOSA Spektralschätzer. Betrachten wir nun den zweiten gemeinsamen Ansatz zur Varianzreduzierung, nämlich die überlappende Segmentmittelung von Welch (Welch; Carter und Referenzen darin). Die Grundidee besteht darin, eine Zeitreihe in mehrere Blöcke zu unterteilen (z., Segmente), berechnen Sie eine direkte Spektralschätzung für jeden Block und erzeugen Sie dann die WOSA-Spektralschätzung, indem Sie diese Spektralschätzungen zusammen mitteln. Im Allgemeinen dürfen sich die Blöcke überlappen, wobei der Grad der Überlappung durch den Grad der Verjüngung bestimmt wird — je schwerer der Grad der Verjüngung, desto mehr sollten sich die Blöcke überlappen (Thomson ). Mit Ausnahme des Anfangs und des Endes der Zeitreihe werden Datenwerte, die sich in einem Block stark verjüngen, in einem anderen Block leicht verjüngt, so dass wir intuitiv “Informationen” zurückerobern, die aufgrund der Verjüngung in einem Block aus Blöcken, die ihn überlappen, verloren gehen. Da es rechnerisch effizient implementiert werden kann (unter Verwendung des Fast-Fourier-Transformationsalgorithmus) und weil es sehr lange Zeitreihen (oder Zeitreihen mit einem zeitlich variierenden spectmm) verarbeiten kann, ist das WOSA-Schätzschema die Grundlage für viele der kommerziellen Spektrumanalysatoren auf dem Markt.

wobei nn die Gesamtzahl der Blöcke ist und s ein ganzzahliger Verschiebungsfaktor ist, der 0< s≤ns und s (nB-1) = n-ns erfüllt (beachten Sie, dass der Block für j = 0 Datenwerte χ1 ,…,Xns, während der Block für j=nB-1, Xn-ns+1,…,Xe).

( hier ht=0 per Definition für alle t>ns). Wenn wir uns auf den Fall einer 50% igen Blocküberlappung (d. H. s = ns / 2) mit einer Hanning-Datenverjüngung beziehen (eine gängige Empfehlung in der technischen Literatur), kann dies durch die einfache Formel v≈36nB21 (19nB-1) angenähert werden. Wenn also die Anzahl der Blöcke NB zunimmt, nehmen auch die äquivalenten Freiheitsgrade zu, was einen Spektralschätzer mit reduzierter Varianz ergibt. Sofern SX (⋅) keine relativ strukturlose sdf aufweist, können wir nB jedoch nicht beliebig klein machen, ohne eine starke Verzerrung der einzelnen direkten Spektralschätzer zu verursachen, die hauptsächlich auf Auflösungsverlust zurückzuführen ist. (Einzelheiten zu den obigen Ergebnissen finden Sie in Abschnitt 6.17.)

Abbildung 8(b) zeigt einen WOSA-Spektralschätzer für die Drahtwellenmessdaten (die durchgezogene Kurve). Diese Serie hat n = 4096 Datenwerte. Einige Experimente zeigten, dass eine Blockgröße von ns = 256 und die Hanning-Datenverjüngung eine vernünftige Wahl für die Schätzung der sdf zwischen 0,4 und 4,0 Hz unter Verwendung von WOSA sind. Bei einer Blocküberlappung von 50% beträgt der Verschiebungsfaktor s = ns / 2 = 128; die Gesamtzahl der Blöcke beträgt nB = 1_δ (n−ns) + 1 = 31; und v, die äquivalenten Freiheitsgrade, beträgt ungefähr 59. Die 31 einzelnen direkten Spektralschätzungen, die zusammen gemittelt wurden, um die WOSA-Schätzung zu bilden, sind in Abbildung 8 (b) als Punkte dargestellt.

Wir haben auch ein “Bandbreite / Konfidenzintervall” kreuz und quer ähnlich dem in Abbildung 8 (a) aufgetragen, aber jetzt ist die “Bandbreite” (dh die horizontale Breite) der Abstand in der Frequenz zwischen ungefähr unkorrelierten spektralen Schätzungen. Das Maß für die Bandbreite ist eine Funktion der Blockgröße ns und der in WOSA verwendeten Datenmenge. Für den Hanning-Kegel beträgt die Bandbreite ungefähr 1,94 / (nsΔt). Die Kreuzungen in den Abbildungen 8(a) und 8 (b) sind ziemlich ähnlich, was darauf hinweist, dass die statistischen Eigenschaften des nachgefärbten Parzen-Lag-Fensters und der WOSA-Spektralschätzungen vergleichbar sind: tatsächlich stimmen die tatsächlichen Schätzungen eng überein, wobei die WOSA-Schätzung etwas glatter aussieht.

Multitaper Spektralschätzer. Eine interessante Alternative zur Lag-Fenster- oder WOSA-Spektralschätzung ist der Multitaper-Ansatz von Thomson . Die Multitaper-Spektralschätzung kann als eine Möglichkeit angesehen werden, einen direkten Spektralschätzer mit mehr als nur zwei äquivalenten Freiheitsgraden (typische Werte sind 4 bis 16) zu erzeugen. Als solches unterscheidet sich das Multitaper-Verfahren im Geiste von den anderen beiden Schätzern dadurch, dass es nicht versucht, stark geglättete Spektren zu erzeugen. Eine Erhöhung der Freiheitsgrade von 2 auf nur 10 reicht jedoch aus, um die Breite eines 95% -Konfidenzintervalls für die sdf um mehr als eine Größenordnung zu verkleinern und damit die Variabilität in der spektralen Schätzung so weit zu reduzieren, dass das menschliche Auge die Gesamtstruktur leicht erkennen kann. Detaillierte Diskussionen zum Multitaper-Ansatz finden Sie in und Kapitel 7 von . Hier skizzieren wir lediglich die wichtigsten Ideen.

( Thomson plädiert dafür, die S ^ k, X (mt) (f) adaptiv zu gewichten, anstatt sie einfach zusammen zu mitteln. Ein Vergleich dieser Definition für S ^ k, X(mt)(⋅) mit Gl. (118) zeigt, dass S ^ k, X (mt) (⋅) tatsächlich nur ein direkter Spektralschätzer ist, so dass der Multitaper-Schätzer nur ein Durchschnitt von direkten Spektralschätzern ist, die einen orthonormalen Satz von Verjüngungen verwenden. Unter bestimmten milden Bedingungen übersetzt sich die Orthonormalität der Verjüngungen in den Frequenzbereich als ungefähre Unabhängigkeit jedes einzelnen S ^ k, X (mt) (f); d. H. S ^ j.X (mt) (f). Ungefähre Unabhängigkeit impliziert wiederum, dass 2KS ^ k, X (mt) (f) / SX (f) = dX22k ungefähr, so dass die äquivalenten Freiheitsgrade für S ^ X (mt) (f) gleich der doppelten Anzahl der verwendeten Datenverjüngungen sind.

Mit anderen Worten, vorbehaltlich der Einschränkung, orthogonal zu allen dpss-Verjüngungen niedrigerer Ordnung zu sein, ist die dpss-Verjüngung k-ter Ordnung in dem eingeschränkten Sinne “optimal”, dass die Seitenkeulen ihres Spektralfensters so weit wie möglich unterdrückt werden, gemessen durch das Konzentrationsverhältnis. Methoden zur Berechnung der DPSS-Datenverjüngungen werden in Kapitel 8 erläutert.

In einer Reihe von Artikeln hat Slepian (und darin enthaltene Referenzen) die Natur von dpss eingehend untersucht. Eine wichtige Tatsache, die er diskutiert, ist, dass das Konzentrationsverhältnis λk (n, W) streng abnimmt, wenn k in einer Weise zunimmt, so dass λk (n, W) nahe der Einheit für k < 2nW Δt ist, wonach es sich schnell nähert 0 mit zunehmendem k (der Wert 2nWΔt wird manchmal als Shannon-Zahl bezeichnet). Da λk (n, W) nahe der Einheit sein muss, damit {ht, k} eine anständige Datenverjüngung ist, ist die Multitaper-Spektralschätzung auf die Verwendung von höchstens — und in der Praxis normalerweise weniger – 2nWΔt orthonormalen DPSS-Verjüngungen beschränkt.

Ein Beispiel für die Multitaper-Spektralschätzung ist in Abbildung 9 dargestellt. Die linke Spalte der Diagramme zeigt die DPSS-Daten k-ter Ordnung für n = 4096, nW = 4 / Δt und k im Bereich von 0 (oberes Diagramm) bis K-1 = 5 (unteres Diagramm). Die dünnen horizontalen Linien in jedem dieser Diagramme zeigen den Nullpegel an, so dass, während der dpss nullter Ordnung überall streng positiv ist (aber ziemlich nahe bei 0 in der Nähe von t = 1 und t = n), die Verjüngungen höherer Ordnung sowohl positive als auch negative Werte annehmen. Beachten Sie auch, dass die Verjüngung nulliger Ordnung Werte der Zeitreihen nahe t = 1 und t = n stark herabsetzt, dass diese Werte jedoch sukzessive durch die Verjüngungen höherer Ordnung stärker gewichtet werden (eine Interpretation von Multitapering ist, dass die Verjüngungen höherer Ordnung Informationen “verlieren”, wenn nur eine einzige Datenverjüngung verwendet wird). Die durchgezogene Kurve in Abbildung 9 (b) zeigt eine Multitaper-Spektralschätzung S^ X(mt) (⋅) für die Drahtwellenmessdaten basierend auf diesen 6 DPSS-Verjüngungen, während die Punkte die sechs einzelnen direkten Spektralschätzungen S ^K.X(mt) (⋅). Beachten Sie, dass die Anzahl der von uns verwendeten Verjüngungen unter der Shannon-Zahl 2nWΔt = 8 liegt und dass v, die äquivalenten Freiheitsgrade, hier 2K = 12 ist. Die Multitaper-Spektralschätzung sieht viel zerhackter aus als entweder die Lag-Fenster-Spektralschätzung von Abbildung 8 (a) oder die WOSA-Schätzung von Abbildung 8 (b), die beide eine deutlich höhere Anzahl äquivalenter Freiheitsgrade aufweisen (v = 64 bzw. v = 59). Nichtsdestotrotz ist die Variabilität in der Multitaper-Spektralschätzung klein genug, so dass das Auge die Gesamtstruktur leicht erkennen kann (vgl. S ^ X (mt) (⋅) mit den beiden Spektralschätzungen in Abbildung 5), und weil es nicht stark geglättet ist, kann die Multitaper-Schätzung die Spektralstruktur in der Nähe von f = 0 deutlich besser erfassen.

Abb.9. Multitaper Spektralschätzung

Basierend auf Leistungsgrenzen argumentiert Bronez [16], dass der Multitaper-Spektralschätzer statistische Eigenschaften aufweist, die WOSA für sdfs mit sehr hohen Dynamikbereichen überlegen sind (mehr Forschung ist jedoch erforderlich, um zu überprüfen, ob diese Grenzen in der Praxis zu einem tatsächlichen Vorteil führen). Im Vergleich zum Prewhitening ist Multitapering in Situationen nützlich, in denen Leckage ein Concem ist, aber es ist nicht praktikabel, Prewhitening-Filter sorgfältig zu entwerfen (dies tritt beispielsweise in der Explorationsgeophysik aufgrund des enormen Volumens an routinemäßig gesammelten Zeitreihen auf). Schließlich stellen wir fest, dass Thomson und Chave [17 ein ansprechendes Schema beschreiben, bei dem Multitapering in Verbindung mit WOSA verwendet wird.