Wie kann man die Eigenschaften der Kovarianzmatrix sicherstellen, wenn man ein multivariates normales Modell mit maximaler Wahrscheinlichkeit anpasst?

22

Angenommen, ich habe das folgende Modell

yi=f(xi,θ)+εi

Dabei ist , ein Vektor erklärender Variablen, die Parameter der nichtlinearen Funktion und , wobei natürlich Matrix.x i θ f iN ( 0 , ) K × KyiRKxiθfεiN(0,Σ)ΣK×K

Das Ziel ist die übliche Schätzung von und . Die naheliegende Wahl ist die Maximum-Likelihood-Methode. Log-Wahrscheinlichkeit für dieses Modell (vorausgesetzt, wir haben eine Stichprobe ) sieht so aus& Sigma; ( y i , x i ) , i = 1 , . . . , nθΣ(yi,xi),i=1,...,n

l(θ,Σ)=n2log(2π)n2logdetΣi=1n(yif(xi,θ))Σ1(yf(xi,θ)))

Nun scheint dies einfach zu sein, die Log-Wahrscheinlichkeit wird angegeben, Daten eingegeben und ein Algorithmus für die nichtlineare Optimierung verwendet. Das Problem ist, wie sichergestellt werden kann, dass definitiv positiv ist. Die Verwendung von beispielsweise in R (oder einem anderen nichtlinearen Optimierungsalgorithmus) garantiert mir nicht, dass eindeutig positiv ist.ΣΣoptimΣ

Die Frage ist also, wie sichergestellt werden kann, dass eindeutig positiv bleibt. Ich sehe zwei mögliche Lösungen:Σ

  1. als umparametrieren wobei eine obere Dreiecksmatrix oder eine symmetrische Matrix ist. Dann ist immer positiv und kann frei sein.R R ' R Σ RΣRRRΣR

  2. Verwenden Sie die Profilwahrscheinlichkeit. Leiten Sie die Formeln für θ^(Σ) und Σ^(θ) . Beginnen Sie mit einigem θ0 und iterieren Sie Σ^j=Σ^(θ^j1) , θ^j=θ^(Σ^j1) bis zur Konvergenz.

Gibt es einen anderen Weg und was ist mit diesen beiden Ansätzen, funktionieren sie, sind sie Standard? Dies scheint ein ziemlich normales Problem zu sein, aber die schnelle Suche hat mir keine Hinweise gegeben. Ich weiß, dass auch eine Bayes'sche Schätzung möglich wäre, aber im Moment würde ich mich nicht darauf einlassen wollen.

mpiktas
quelle
Ich habe das gleiche Problem in einem Kalman-Algorithmus, aber das Problem ist viel komplizierter und nicht so einfach, den Hamilton-Trick anzuwenden. Ich frage mich dann, ob es einfacher wäre, einfach . Auf diese Weise erzwinge ich, dass der Code keinen Fehler ausgibt und ändere die Lösung nicht. Dies hat auch den Vorteil, dass dieser Term das gleiche Vorzeichen wie der letzte Teil der Wahrscheinlichkeit haben muss. Irgendwelche Ideen? log(detΣ+1)
econ_pipo

Antworten:

6

Unter der Annahme, dass Sie sich beim Erstellen der Kovarianzmatrix automatisch um das Symmetrieproblem kümmern, ist Ihre Log-Wahrscheinlichkeit wenn aufgrund des Terms in nicht positiv bestimmt ist Modell richtig? Um einen numerischen Fehler zu vermeiden, wenn würde ich vorberechnen und, falls dies nicht positiv ist, die Log-Wahrscheinlichkeit gleich -Inf machen, andernfalls fortfahren. Sie müssen die Determinante sowieso berechnen, so dass dies Sie keine zusätzliche Berechnung kostet. Σ log d e t Σ d e t Σ < 0 d e t ΣΣlogdet Σdet Σ<0det Σ

Makro
quelle
5

Wie sich herausstellt, können Sie die maximale Profilwahrscheinlichkeit verwenden, um die erforderlichen Eigenschaften sicherzustellen. Sie können für gegebene beweisen , dass , durch maximiert l( θ ,Σ)θ^l(θ^,Σ)

Σ^=1nich=1nε^ichε^ich,

woher

ε^ich=yich-f(xich,θ^)

Dann ist es möglich, das zu zeigen

ich=1n(yich-f(xich,θ^))Σ^-1(y-f(xich,θ^)))=cOnst,

Daher müssen wir nur maximieren

lR(θ,Σ)=-n2LogdetΣ^.

In diesem Fall erfüllt natürlich alle erforderlichen Eigenschaften. Die Beweise sind identisch für den Fall, dass linear ist, was in der Zeitreihenanalyse von JD Hamilton auf Seite 295 zu finden ist, daher habe ich sie weggelassen.fΣf

mpiktas
quelle
3

Eine alternative Parametrisierung für die Kovarianzmatrix sind die Eigenwerte und "Givens" -Winkel . p ( p - 1 ) / 2 θ i jλ1,...,λpp(p-1)/2θichj

Das heißt, wir können schreiben

Σ=GTΛG

wobei orthonormal ist undG

Λ=dicheinG(λ1,...,λp)

mit .λ1...λp0

In der Zwischenzeit kann in Form von Winkeln, , mit und eindeutig parametrisiert werden . [1]p ( p - 1 ) / 2 θ i j i = 1 , 2 , . . . , P - 1 j = i , . . . , p - 1Gp(p-1)/2θichjich=1,2,...,p-1j=ich,...,p-1

(Details hinzugefügt werden)

[1]: Hoffman, Raffenetti, Ruedenberg. "Verallgemeinerung von Eulerwinkeln zu N-dimensionalen orthogonalen Matrizen". J. Math. Phys. 13, 528 (1972)

charles.y.zheng
quelle
Die Matrix ist eigentlich orthogonal, weil eine symmetrische Matrix ist. Dies ist der Ansatz, den ich empfehlen wollte - Im Grunde genommen bedeutet dies, den Vektor und die Modellfunktion so zu drehen , dass die Fehler unabhängig sind, und dann OLS auf jede der gedrehten Komponenten anzuwenden (glaube ich). Σ y i f ( x i , θ )GΣyichf(xich,θ)
Wahrscheinlichkeitslogik
2

In Anlehnung an die Lösung von charles.y.zheng möchten Sie möglicherweise modellieren , wobei eine Diagonalmatrix und eine Cholesky-Faktorisierung einer Rangaktualisierung von . Sie müssen nur dann die Diagonale von positiv halten, um positiv eindeutig zu halten . Das heißt, Sie sollten die Diagonale von und die Elemente von schätzen, anstatt schätzen . Λ C Λ Λ Σ Λ C ΣΣ=Λ+CCΛCΛΛΣΛCΣ

shabbychef
quelle
Darf ich unter diagonalen Elementen in diesen Einstellungen irgendetwas haben, solange die Diagonale positiv ist? Wenn Matrizen auf diese Weise in Anzahl simuliert werden, sind nicht alle positiv bestimmt.
11.
Λ ist eine Diagonalmatrix.
Shabbychef