State Space Models - Teil 3
MARSS
und die echten Daten
Alles was wir bis jetzt mit MARSS
gemacht haben, war eher eine Fingerübung, ein vorsichtiges rantasten. Jetzt wollen wir die Estimation-Bestie entfesseln und auf echte Daten loslassen. Bisher haben wir nur “Terra” geformt, also das was am naheliegensten ist: selbst simulierte Daten. Wer kennt es nicht: Man kann sich ja andere Daten gerne anschauen, aber simuliert wird immer noch zuhause. Im nächsten Schritt wollen wir uns aber hinauswagen – etwas neues probieren. Wir wollen “State” und “Space” formen als wäre es ein Spiel. Als spielten wir “STATE N’ SPACE FORMING MARSS
”!
Welche Daten entsprechen unserem Geschmack? Wir wählen in Österreich am 29. September den Nationalrat. Zur Wahl stehen offensichtlich 9 Parteien (ÖVP, SPÖ, FPÖ, GRÜNE, NEOS, KPÖ, BIER, Keine, Madeleine Petrovic). Jeder möchte schon vorher wissen wie die Wahl ausgeht, denn je präziser die Vorhersage, desto unwichtiger scheint die eigene Stimme zu sein und desto weniger das schlechte Gewissen seine Stimme nicht abgegeben zu haben. “Trendprognosen”, wie zuletzt bei der EU-Wahl, sind derzeit hoch im Kurs und MARSS
öffnet uns die Türe diese Daten zu bändigen und ordentliche Forecasts zu machen.
Allerdings stellen uns diese Daten vor einige Herausforderungen:
- Wir haben nun eine multivariate Zeitreihe (es gibt ja den Verlauf mehrerer Parteien und nicht bloß eine Zeitreihe wie in den bisherigen 2 Blogpsts – Danke Demokratie)
- Wir haben es hier mit Anteilen zu tun. Eine Partei kann nicht über 100 %, oder unter 0 % der Stimmen haben. Auch die Summe der Anteile aller Parteien kann nicht über 100 % liegen. Wir müssen die Daten also irgendwie geschickt transformieren.
- Die Daten der Umfragen sind nicht äquidistant. In den bisherigen 2 Blogposts sind wir implizit davon ausgegangen, dass ein Sprung von Zeitpunkt 1 zu 2 genauso groß ist wie von 2 zu 3. Alle Abstände wurden als gleich weit voneinander entfernt angenommen. Umfragedaten sind mal dichter gestaffelt, mal ist der Abstand zwischen den Umfragen größer. Das macht uns keine Freude aber wir müssen mit diesem Umstand irgendwie konstruktiv umgehen.
Ok, die Fronten sind geklärt. STATE N’ SPACE FORMING MARSS
kann beginnen! Wir versuchen aus unserer Statistik-Kartenhand die passenden und natürlich hochpotentesten Karten zum richtigen Zeitpunkt zu spielen!
Multivariate Zeitreihe
Unsere Daten sind multivariat – man kann es drehen und wenden wie man will. Niemand mag multivariate Daten, niemand mag Methoden zu multivariaten Daten, jeder wünscht sich sein y rank und schlank wie einen Vektor. Das MARSS
package kann glücklicherweise out of the box mit multivariaten Zeitreihen umgehen – die Dinge werden trotzdem komplizierter. y wird also vom Zeitreihenvektor zur Zeitreihenmatrix. Damit ist es aber nicht getan. Das ganze Modell wird deutlich komplexer, d. h. wir müssen uns auch mit (fast) allen Komponenten des Modells auseinandersetzen wie sie in Gleichung 1 dargestellt sind. Warum? Wir müssen MARSS
direkt die jeweiligen Matrizen übergeben und so unser maßgeschneidertes Modell zu definieren.
Zeile 1 von Gleichung 1 zeigt die Modellierung der “True Values” – also der wahren Werte des Modells, quasi die state equation. In Zeile 2 kommen die Beobachtungen ins Spiel! Zeile 3 sagt nur etwas über den Startzustand der ersten Beobachtungen aus. Was steckt sonst noch drinnen?
- \(x_{m \times T}\) \(\,\,\rightarrow\) hier stecken die states drinnen, also unsere True Values.
- \(B_{m \times m}\) \(\,\,\rightarrow\) Diese Matrix kombiniert die states vom Zeitpunkt t-1 und erklärt damit den aktuellen, latenten state vom Zeitpunkt t!
- \(Q_{m \times m}\) \(\,\,\rightarrow\) Kovarianzmatrix der einzelnen Zustände. In unseren bisherigen Modellen, wäre diese Matrix einfach eine Diagnoalmatrix gewesen. In einem multivarianten Modell ist die Vorstellung, dass die Zeitreihen sich nicht unanbhängig voneinander entwickeln, ja naheliegend. Dementsprechend sollten wir die Q-Matrix im Auge behalten.
- \(R_{n \times n}\) \(\,\,\rightarrow\) Diese Matrix ist Teil der Beobachtungsgleichung und stellt den “klassischen” Fehlerterm dar, also das “weiße Rauschen”, das dem systematischen Verlauf der True Values hinzugefügt wird. Hier ist es auch möglich Kovarianzen zu formulieren.
- \(y_{n \times T}\) \(\,\,\rightarrow\) Das ist unsere Zeitreihe – also die Daten die uns vorliegen und die wir fitten wollen!
- \(Z_{m \times n}\) \(\,\,\rightarrow\) Definiert die Beziehung der states zu den Beobachtungen.
- \(u, a\) sind typischerweise 0, aber das wissen wir ohnehin schon von Teil 1 und 2, und stellen die Erwartungswerte der multivariaten Normalverteilungen dar, die quasi den Möglichkeitsrahmen des Measurement Errors vorgeben.
Hier ist eine detaillierte Darstellung des MARSS
Modells samt Beispielen.
\[ \begin{aligned} \begin{bmatrix}x_1\\ x_2\end{bmatrix}_t &= \begin{bmatrix}b_{11}&b_{12}\\b_{21}&b_{22}\end{bmatrix} \begin{bmatrix}x_1\\x_2\end{bmatrix}_{t-1} + \begin{bmatrix}w_1\\ w_2\end{bmatrix}_t,\quad \begin{bmatrix}w_1\\ w_2\end{bmatrix}_t \sim \textup{\textrm{MVN}}\begin{pmatrix}\begin{bmatrix}u_1\\u_2\end{bmatrix},\begin{bmatrix}q_{11}&q_{12}\\q_{21}&q_{22}\end{bmatrix} \end{pmatrix} \\ \\ \begin{bmatrix}y_1\\ y_2\\ y_3\end{bmatrix}_t &= \begin{bmatrix}z_{11}&z_{12}\\ z_{21}&z_{22}\\ z_{31}&z_{32}\end{bmatrix} \begin{bmatrix}x_1\\ x_2\end{bmatrix}_t + \begin{bmatrix}v_1\\ v_2\\ v_3\end{bmatrix}_t, \begin{bmatrix}v_1\\ v_2\\ v_3\end{bmatrix}_t \sim \textup{\textrm{MVN}}\begin{pmatrix}\begin{bmatrix}a_1\\ a_2\\ a_3\end{bmatrix}, \begin{bmatrix}r_{11}&r_{12}&r_{13}\\r_{21}&r_{22}&r_{23}\\r_{31}&r_{32}&r_{33}\end{bmatrix} \end{pmatrix} \\ \\ \begin{bmatrix}x_1\\ x_2\end{bmatrix}_0 &\sim \textup{\textrm{MVN}}\begin{pmatrix}\begin{bmatrix}\pi_1\\ \pi_2\end{bmatrix},\begin{bmatrix}\nu_{11}&\nu_{12}\\ \nu_{21}&\nu_{22}\end{bmatrix} \end{pmatrix} \quad or\quad \begin{bmatrix}x_1\\ x_2\end{bmatrix}_1 \sim \textup{\textrm{MVN}}\begin{pmatrix}\begin{bmatrix}\pi_1\\ \pi_2\end{bmatrix},\begin{bmatrix}\nu_{11}&\nu_{12}\\ \nu_{21}&\nu_{22}\end{bmatrix} \end{pmatrix} \end{aligned} \tag{1}\]
Aber um zu wissen wie wir die Modelle aufbauen, muss man zuerst wissen wie die Daten nun tatsächlich aussehen.
Daten
Um eine Vorhersage für die Nationalratswahlen zu machen brauchen wir österreichische Umfragedaten, die im wesentlichen abfragen: “Wen würden Sie wählen wenn kommenden Sonntag Nationalratswahl wäre?”. Wir finden hier Umfragedaten aus unterschiedlichsten Ländern. Im folgenden werden die österreichischen Daten (at.csv) runtergeladen und direkt weiterverarbeitet und umfassend aufbereitet. Eine grobe Übersicht bietet Abbildung 1. Hier wird der Verlauf der Anteile jeder Partei nach Umfrageinstitut getrennt dargestellt.
Was bedeutet aber “aufbereiten”? Wie soll das Endprodukt aussehen? Unsere Daten brauchen:
- Datum: Wann fand die Erhebung statt? Umfragen finden typischerweise nicht an einem einzelnen Tag statt, sondern in einem Zeitraum von mehreren Tagen. Um jeder Umfrage ein konkretes Datum zuzuordnen wurde das Datum des Umfragestarts gewählt. Unterscheidet man allerings nicht zwischen den Umfrageinstituten im Modell, so kann es vorkommen, dass Umfragen zur selben Zeit stattfinden – also das selbe Start-Datum aufweisen.
MARSS
kann mit dieser etwaigen Gleichzeitigkeit nicht umgehen. Daher wurde in solchen Fällen ein passendes Datum zwischen Start und Enddatum der Umfrage gewählt, so dass es zu keinen Überschneidungen kommt. Damit sind alle Zeitpunkte unique! - Anteile je Partei: Wie hoch ist der Anteil jener Leute, die diese Partei wählen würden? Diese müssen sich auf 1 summieren (Achtung Rundungsfehler).
Eine sehr übersichtliche Anzahl an Merkmalen. Allerdings müssen hier noch 2 weiter Schritte durchgeführt werden. Es muss Äquidistanz hergestellt, und die Daten transformiert werden.
Äquidistanz
In den meisten State-Space Model packages wird davon ausgegangen, dass die Zeitabstände zwischen den Beobachtungen gleich sind. Bei unseren Daten ist das nicht der Fall, da die einzelnen Umfragen bei Bedarf durchgeführt werden und üblicherweise nicht in fixen Abständen. MARSS
kann zwar nicht mit unterschiedlichen Zeitabständen umgehen, dafür machen dem package Missing Values keine Kopfzerbrechen. Diese Stärke können wir für uns nutzen. Wir erstellen eine Datenmatrix die alle Tage zwischen Anfang und Ende der Zeitreihe beinhaltet. Da wir an den meisten Tagen keine Beobachtung haben, ergibt das eine Menge Missing Values. Das löst unser Problem mehr oder weniger elegant, denn nun sind alle Zeitpunkte gleichabständig – auch wenn nicht alle Zeitpunkte Beobachtungen ausweisen! Die meisten Zeitpunkte weisen tatsächlich Missing Values auf!
Transformation
Im Fall von Umfragedaten wie unseren, müssen wir mit Anteile umgehen und der Tatsache zurechtkommen, dass unsere Variablen auf Zahlen zwischen 0 und 1 festgelegt sind, es sich also nicht überall zwischen \((-\infty,\infty)\) gemütlich machen können. Das betrachten wir natürlich mit großer Sorge. Es gibt 3 Möglichkeiten mit dieser Unannehmlichkeit umzugehen:
- Die falsche Methode: Wir modellieren mittels linearen Modells ohne Transformation. Problem: Es können (und werden) Werte > 1 und < 0 resultieren. “No, no, no, no, nooo! Don’t do that to me!” Ist das Originalzitat von gleichartigen Daten die mit solchen Methoden traktiert wurden.
- Die richtige Methode 1: Wir modellieren die Anteile direkt z. B. mittels einer Dirichlet Verteilung und modellieren dessen Parameter über die Zeit mittels eines State-Space Models. Spoiler: Ich habe dazu kein R-Package gefunden.
- Die richtige Methode 2: Transformieren! Wir transformieren bei k Parteien diese k Zeitreihen so, dass k-1 Zeitreihen resultieren deren mögliche Ausprägungen zwischen \(-\infty\) und \(\infty\) liegen. Dafür müssen wir uns nichts selbst ausdenken oder herleiten , sondern dafür gibt es packages und Lösungen die uns mit diesem Problem helfen, wie z. B. das compositions package.
Unser Modell
Um einen groben Überblick über das Modell zu bekommen: Was sind die wesentlichsten Komponenten? Insbesondere jene, die wegen des multivariaten Charakters der Daten anders aussehen als in Teil 1 und Teil 2 dargestellt:
1. Q-Matrix
2. B-Matrix
3. Z-Matrix
Q Matrix
Die Matrix ist so aufgebaut, dass sie sowohl die Varianzen als auch Kovarianzen zwischen verschiedenen Elementen (Levels und Slopes) der (latenten) Zustandsvariablen beinhaltet. Wenn ein Parametername übergeben wird, werden diese Parameter geschätzt, wenn eine Zelle 0 enthält wird kein entsprechender Parameter geschätzt. Die einfachste Variante wäre eine Diagonalmatrix – aber wir powern den Q-Matrix Enhancer und iterieren uns bis zur Konvergenz (hoffentlich).
Die Diagonalelemente der Matrix, wie etwa Q.var_lvl1
oder Q.var_slp1
, stehen für die Varianz der jeweiligen Zustandsvariablen. Die Nicht-Diagonalelemente wie z. B. Q.cov_lvl1_lvl2
oder Q.cov_slp1_slp2
, repräsentieren die Kovarianzen zwischen den verschiedenen Variablen. Wie in der Matrix zu sehen ist werden nur Kovarianzen zwischen allen Tupels von Slopes und Levels angenommen. Intuitiv könnte man sagen: Wenn eine Partei Zugewinne verzeichnet, der Slope also positiv ist, muss zumindest eine andere Partei verlieren. Tatsächlich können auch 2 oder alle anderen Parteien verlieren. Diese Kovarianzparameter versuchen diese Information aus den Daten zu extrahieren. Man nimmt also an, dass im Falle eines steigenden Slopes nicht alle anderen beliebig sinken. Daher werden Abhängigkeiten zwischen den Slopes, und ebenso zwischen den Levels geschätzt, aber eben nicht zwischen Levels und Slopes. Dies ist sowohl eine inhaltliche Entscheidung (was wäre die Interpretation?) als auch eine Frage der Sparsamkeit. Abbildung Abbildung 4 zeigt später die Struktur dieser Matrix nochmal auf eine andere Weise.
\[\begin{gather*} \small \begin{bmatrix} Q.var\_lvl1&0&Q.cov\_lvl1\_lvl2&0& \ldots\\ 0&Q.var\_slp1&0&Q.cov\_slp1\_slp2& \ldots\\ Q.cov\_lvl1\_lvl2&0&Q.var\_lvl2&0& \ldots\\ 0&Q.cov\_slp1\_slp2&0&Q.var\_slp2& \ldots\\ \vdots & \vdots & \vdots & \vdots \end{bmatrix} \end{gather*}\]
B - Matrix
Die B-Matrix enthält ebenfalls selbst keine Parameter. Sie ist im Wesentlichen aus Blöcken von identen \(2 \times 2\) Matrizen aufgebaut ist. Und zwar genau aus jenen \(2 \times 2\) Matrizen, die die Parameter (x) zu den beiden aus Teil 1 und 2 bekannten Zustandsgleichungen kombinieren und aus ihnen unser bekanntes Local-Linear Trend Model machen. Die x der Vergangenheit (t-1) werden zu den x der Zukunft zu kombiniert.
\[\begin{gather*} \begin{bmatrix} x_{X1}\\ x_{X2}\\ x_{X3}\\ x_{X4}\\ x_{X5}\\ x_{X6}\\ x_{X7}\\ x_{X8}\\ x_{X9}\\ x_{X10} \end{bmatrix}_{t}= \begin{bmatrix} 1&1&0&0&0&0&0&0&0&0\\ 0&1&0&0&0&0&0&0&0&0\\ 0&0&1&1&0&0&0&0&0&0\\ 0&0&0&1&0&0&0&0&0&0\\ 0&0&0&0&1&1&0&0&0&0\\ 0&0&0&0&0&1&0&0&0&0\\ 0&0&0&0&0&0&1&1&0&0\\ 0&0&0&0&0&0&0&1&0&0\\ 0&0&0&0&0&0&0&0&1&1\\ 0&0&0&0&0&0&0&0&0&1 \end{bmatrix} \begin{bmatrix} x_{X1}\\ x_{X2}\\ x_{X3}\\ x_{X4}\\ x_{X5}\\ x_{X6}\\ x_{X7}\\ x_{X8}\\ x_{X9}\\ x_{X10} \end{bmatrix}_{t-1} \end{gather*}\]
Z-Matrix
Diese Matrix ist dafür verantwortlich, die Beziehung zwischen den beobachteten Daten und den zugrunde liegenden latenten Zustandsvariablen herzustellen, d. h. die Zustandsvariablen werden gruppiert. Jede beobachtete Variable (y_transf1
bis y_transf5
) wird als lineare Kombination von Paaren von Zustandsvariablen (x_X1
bis x_X10
) dargestellt. Das bedeutet, dass jede beobachtete Zeitreihe von zwei Zustandsvariablen beschrieben wird – diese sollen zusammen die beobachtete Dynamik beschreiben.
\[\begin{gather*} \begin{bmatrix} y_{transf1}\\ y_{transf2}\\ y_{transf3}\\ y_{transf4}\\ y_{transf5} \end{bmatrix}_{t}= \begin{bmatrix} 1&1&0&0&0&0&0&0&0&0\\ 0&0&1&1&0&0&0&0&0&0\\ 0&0&0&0&1&1&0&0&0&0\\ 0&0&0&0&0&0&1&1&0&0\\ 0&0&0&0&0&0&0&0&1&1 \end{bmatrix} \begin{bmatrix} x_{X1}\\ x_{X2}\\ x_{X3}\\ x_{X4}\\ x_{X5}\\ x_{X6}\\ x_{X7}\\ x_{X8}\\ x_{X9}\\ x_{X10} \end{bmatrix}_{t} \end{gather*}\]
Ergebnisse und Forecasts
Nachdem wir nun lange auf Matrizen geschaut haben (oder auch nicht) wird es endlich noch spannender, denn wir wollen wissen wie die jeweiligen Parteien abschneiden. Für alle Leser:innen die gespannt auf die Ergebnisse von “Liste Petrovic” und “Keine” warten: Diese beiden Parteien gehen in der Kategorie “Other” auf. Als vorsichtige Prognose gleich vorweg: Hoch werden die beiden Parteien das Match nicht gewinnen.
Abbildung 2 zeigt einen gemittelten Forecast aus 3 unterschiedlichen Modellen, unser Model Ensamble (ME):
- Das Modell mit der “Enhanced” Q-Matrix, also mit der Extraportion Parameter (unser Best Fitting Model BFM)
- Ein Modell das als Q-Matrix eine Diagonalmatrix annimmt, aber für jede Partei wird eine eigene Varianz geschätzt (sehr einfaches Modell)
- Ein noch sparsameres Modell: Gleiche Varianz für alle Parteien in der Q-Matrix (noch viel einfacheres Modell)
Abbildung 3 zeigt unser Best Fitting Modell (BFM) mit den extra steilen und kurvigen Slopes für den Sensation-Seeker unter den Politik-Insider.
Warum ist es prinzipiell klug Forecast Ergebnisse zu mitteln? Weil die Annahmen von Modellen typischerweise nie ganz richtig sind (Achtung: Euphemismus) und sich durch die Kombination von unterschiedlichen Modellen die typischerweise unterschiedliche Annahmen haben, im günstigsten Fall alle diese kleinen Sonderbarkeiten, Vereinfachungen oder Überkomplexitäten der einzelnen Modelle rausmitteln würden. Oft wird bei einem derartigen Ensamble beim Mitteln anhand von AIC oder BIC gewichtet. In unserem Fall ist das etwas witzlos, denn das Modell mit der Enhanced Q-Matrix laut AIC und BIC deutlich besser performt also die anderen beiden Modelle – daher natürlich auch: Best Fitting Modell (BFM). Damit wäre die Gewichtung zu 100 % auf dem BFM, was das Ensamble ad absurdum führt. Also gehen in diesem Fall alle Modelle gleichermaßen stark in das Ensamble ein.
Warum?
Weil die Verläufe des BFM einfach sehr sehr unrealistisch aussehen. Jeder der schon einmal etwas mit Statistik zu tun hat weiß, dass das wirklich zu viel des Guten ist. Das Modell overacted und überinterpretiert wohl die gefundenen Effekte.
Partei | Last Polls | Forecast | Trend |
---|---|---|---|
FPÖ | 27 % | 26 % | |
ÖVP | 23 % | 24 % | |
SPÖ | 22 % | 21.2 % | |
NEOS | 10 % | 9.8 % | |
GRÜNE | 8 % | 8.3 % | |
Other | 3 % | 3.7 % | |
BIER | 5 % | 3.7 % | |
KPÖ | 2 % | 3.3 % |
Partei | Last Polls | Forecast | Trend |
---|---|---|---|
FPÖ | 27 % | 26.1 % | |
ÖVP | 23 % | 23.3 % | |
SPÖ | 22 % | 21.6 % | |
NEOS | 10 % | 8.7 % | |
GRÜNE | 8 % | 7.9 % | |
Other | 3 % | 6.7 % | |
KPÖ | 2 % | 5.6 % | |
BIER | 5 % | 0.1 % |
Unsere Modelle unterscheiden sich primär anhand ihrer Q-Matrizen. Also betrachten wir die Q-Matrix unseres BFM um vielleicht so diese rasanten predicts einigermaßen erklären zu können. Dazu müssen wir die Kovarianzen die uns MARSS
zurückgibt in Korrelationen umwandeln. Das Ergebnis ist in Abbildung 4 zu sehen. Achtung! Die Slopes und Levels sind keiner einzelnen Partei zuzuordnen, da es sich um die transformierten Werte handelt1! Zu sehen sind einige kleine und moderate Korrelationen, allerdings auch einige sehr hohe/niedrige Korrelationen – Slope 1 und Slope 3 sind zu 0.99 miteinander korreliert, lvl 2 und 3 zu -0.93! Diese Werte sind wohl ähnlich unrealistisch hoch wie die daraus resultierenden Predicts. Da aber die Annahme dass die Slopes und Levels völlig unkorreliert sind auch unrealistisch scheinen, ist eine Enesamle von Modellen vermutlich eine gerechte und praktikable Lösung.
Das war der dritte und letzte Teil zu den State Space Models. Hier sind alle Daten des dritten Teils zum Download bereit, wer Lust hat selbst noch weiter Modelle zu rechnen oder Grafiken zu erstellen.
Fußnoten
Werden k Variablen transformiert, resultieren nach der Transformation k-1 Variablen!↩︎