Kroki, aby ustalić rozkład boczny, kiedy może być wystarczająco prosty, aby mieć postać analityczną?

Zapytano o to również w Computational Science.

Próbuję obliczyć bayesowskie oszacowanie niektórych współczynników dla autoregresji, z 11 próbkami danych: gdzie jest Gaussa ze średnią 0 i wariancją Wcześniejszy rozkład na wektorze jest Gaussa ze średnią i ukośną macierzą kowariancji z wpisy ukośne równe .

Y_{i} = μ + α \cdot Y_{i - 1} + ϵ_{i}

$Y_{i} = \mu + \alpha\cdot{}Y_{i-1} + \epsilon_{i}$

ϵ_{i}

$\epsilon_{i}$

σ_{e}^{2}

$\sigma_{e}^{2}$

(μ, α)^{t}

$(\mu, \alpha)^{t}$

(0, 0)

$(0,0)$

σ_{p}^{2}

$\sigma_{p}^{2}$

W oparciu o formułę autoregresji oznacza to, że rozkład punktów danych ( $Y_{i}$ ) jest normalny ze średnią $\mu + \alpha\cdot{}Y_{i-1}$ i wariancją $\sigma_{e}^{2}$ . Tak więc gęstość wszystkich punktów danych $(Y)$ łącznie (przy założeniu niezależności, co jest w porządku dla programu, który piszę), będzie wynosić:

p (Y | (μ, α)^{t}) = \prod_{i = 2}^{11} \frac{1}{\sqrt{2 π σ_{e}^{2}}} \exp \frac{- (Y_{i} - μ - α \cdot Y_{i - 1})^{2}}{2 σ_{e}^{2}} .

$p(Y \quad | (\mu, \alpha)^{t}) = \prod_{i=2}^{11}\frac{1}{\sqrt{2\pi\sigma_{e}^{2}}}\exp{\frac{-(Y_{i} - \mu - \alpha\cdot{}Y_{i-1})^{2}}{2\sigma_{e}^{2}}}.$

Twierdzeniem Bayesa możemy wziąć iloczyn powyższej gęstości z wcześniejszą gęstością, a wtedy potrzebujemy tylko stałej normalizującej. Mam przeczucie, że powinno to działać jako rozkład Gaussa, więc możemy się martwić o stałą normalizującą na końcu, zamiast jawnie obliczać ją za pomocą całek względem $\mu$ i $\alpha$ .

Z tą częścią mam problem. Jak obliczyć mnożenie wcześniejszej gęstości (która jest wielowymiarowa) i iloczynu gęstości danych jednowymiarowych? Tylny musi mieć czystą gęstość $\mu$ i $\alpha$ , ale nie widzę, jak można to uzyskać z takiego produktu.

Wszelkie wskazówki są naprawdę pomocne, nawet jeśli po prostu skierujesz mnie we właściwym kierunku, a następnie muszę iść i zrobić bałaganiarską algebrę (co próbowałem już kilka razy).

Na początek jest postać licznika z reguły Bayesa:

\frac{1}{(2 π σ_{e}^{2})^{5} \cdot 2 π σ_{p}^{2}} \exp [\frac{1}{2 σ_{e}^{2}} \sum_{i = 2}^{11} (Y_{i} - μ - α \cdot Y_{i - 1})^{2} - \frac{μ^{2}}{2 σ_{p}^{2}} - \frac{α^{2}}{2 σ_{p}^{2}}] .

$\frac{1}{(2\pi\sigma_{e}^{2})^{5}\cdot{}2\pi\sigma_{p}^{2}} \exp{\biggl [ \frac{1}{2\sigma_{e}^{2}}\sum_{i=2}^{11}(Y_{i} - \mu - \alpha\cdot{}Y_{i-1})^{2} - \frac{\mu^{2}}{2\sigma_{p}^{2}} - \frac{\alpha^{2}}{2\sigma_{p}^{2}} \biggr ] }.$

Problem polega na tym, jak zobaczyć, że zmniejsza się to do gęstości Gaussa . $(\mu, \alpha)^{t}$

Dodany

Ostatecznie sprowadza się to do następującego ogólnego problemu. Jeśli otrzymasz jakieś wyrażenie kwadratowe, takie jak jak to ująć w formę kwadratową dla niektórych macierzy 2x2 ? W prostych przypadkach jest to dość proste, ale jakiego procesu używasz, aby uzyskać średnie oszacowania, i ?

A μ^{2} + B μ α + C α^{2} + J μ + K α + L

$A\mu^{2} + B\mu\alpha + C\alpha^{2} + J\mu + K\alpha + L$

(μ - \hat{μ}, α - \hat{α}) Q (μ - \hat{μ}, α - \hat{α})^{t}

$(\mu-\hat{\mu},\alpha-\hat{\alpha})Q(\mu-\hat{\mu},\alpha-\hat{\alpha})^{t}$

Q

$Q$

\hat{μ}

$\hat{\mu}$

\hat{α}

$\hat{\alpha}$

Uwaga: wypróbowałem prostą opcję rozszerzenia formuły macierzowej, a następnie próbowałem zrównać współczynniki jak powyżej. Problem w moim przypadku polega na tym, że stała wynosi zero, a następnie otrzymuję trzy równania w dwóch niewiadomych, więc niedokładne jest jedynie dopasowanie współczynników (nawet jeśli założę symetryczną macierz kwadratową). $L$

bayesian mathematical-statistics posterior

— Ely
źródło

Moja odpowiedź na [to pytanie] ( stats.stackexchange.com/questions/22852/… ) może być pomocna. Zauważ, że potrzebujesz pierwszej do pierwszej obserwacji - iteracje na tym się kończą.

— probabilityislogic

Nie rozumiem, dlaczego potrzebuję tego w tym przypadku. Przedziały czasowe mam traktować tak, jakby były uwarunkowane niezależnie od obserwacji. Zauważ, że iloczyn gęstości złącza wynosi właśnie od . Nie sądzę, żebym miał tu pobierać sekwencyjnie aktualizowaną formułę, tylko jedną formułę dla tylnego .

i = 2..11

$i=2..11$

p ((μ, α)^{t} | Y)

$p((\mu,\alpha)^{t}\quad |Y)$

— ely

„Wielowariantowa” we wcześniejszym nie jest sprzeczna z „jednowymiarową” w gęstości danych, ponieważ są to gęstości w .

p (α, μ)

$p(\alpha,\mu)$

y_{i}

$y_i$

— Xi'an

Wskazówka, która znajdowała się w mojej odpowiedzi na poprzednią odpowiedź, to przyjrzeć się temu, jak zintegrowałem parametry - ponieważ zrobisz tutaj dokładnie te same całki. Pytanie zakłada, że parametry wariancji są znane, więc są one stałymi. Wystarczy spojrzeć na zależność od licznika. Aby to zobaczyć, zauważ, że możemy napisać: $\alpha,\mu$

p (μ, α | Y) = \frac{p (μ, α) p (Y | μ, α)}{\int \int p (μ, α) p (Y | μ, α) d μ d α}

$p(\mu,\alpha|Y)=\frac{p(\mu,\alpha)p(Y|\mu,\alpha)}{\int\int p(\mu,\alpha)p(Y|\mu,\alpha)d\mu d\alpha}$

= \frac{\frac{1}{(2 π σ_{e}^{2})^{5} \cdot 2 π σ_{p}^{2}} \exp [- \frac{1}{2 σ_{e}^{2}} \sum_{i = 2}^{11} (Y_{i} - μ - α \cdot Y_{i - 1})^{2} - \frac{μ^{2}}{2 σ_{p}^{2}} - \frac{α^{2}}{2 σ_{p}^{2}}]}{\int \int \frac{1}{(2 π σ_{e}^{2})^{5} \cdot 2 π σ_{p}^{2}} \exp [- \frac{1}{2 σ_{e}^{2}} \sum_{i = 2}^{11} (Y_{i} - μ - α \cdot Y_{i - 1})^{2} - \frac{μ^{2}}{2 σ_{p}^{2}} - \frac{α^{2}}{2 σ_{p}^{2}}] d μ d α}

$=\frac{\frac{1}{(2\pi\sigma_{e}^{2})^{5}\cdot{}2\pi\sigma_{p}^{2}} \exp{\biggl [ -\frac{1}{2\sigma_{e}^{2}}\sum_{i=2}^{11}(Y_{i} - \mu - \alpha\cdot{}Y_{i-1})^{2} - \frac{\mu^{2}}{2\sigma_{p}^{2}} - \frac{\alpha^{2}}{2\sigma_{p}^{2}} \biggr ] }}{\int\int \frac{1}{(2\pi\sigma_{e}^{2})^{5}\cdot{}2\pi\sigma_{p}^{2}} \exp{\biggl [ -\frac{1}{2\sigma_{e}^{2}}\sum_{i=2}^{11}(Y_{i} - \mu - \alpha\cdot{}Y_{i-1})^{2} - \frac{\mu^{2}}{2\sigma_{p}^{2}} - \frac{\alpha^{2}}{2\sigma_{p}^{2}} \biggr ] }d\mu d\alpha}$

Zauważ, można przesunąć pierwszy współczynnik się podwójnej całki na mianowniku, i kasuje się z licznikiem. Możemy również pobrać sumę kwadratów a także anuluje. Całka, z jaką mamy teraz, jest (po rozszerzeniu kwadratu): $\frac{1}{(2\pi\sigma_{e}^{2})^{5}\cdot{}2\pi\sigma_{p}^{2}}$ $\exp{\biggl [ -\frac{1}{2\sigma_{e}^{2}}\sum_{i=2}^{11}Y_{i}^{2} \biggr ]}$

= \frac{\exp [- \frac{10 μ^{2} + α^{2} \sum_{i = 1}^{10} Y_{i}^{2} - 2 μ \sum_{i = 2}^{11} Y_{i} - 2 α \sum_{i = 2}^{11} Y_{i} Y_{i - 1} + 2 μ α \sum_{i = 1}^{10} Y_{i}}{2 σ_{e}^{2}} - \frac{μ^{2}}{2 σ_{p}^{2}} - \frac{α^{2}}{2 σ_{p}^{2}}]}{\int \int \exp [- \frac{10 μ^{2} + α^{2} \sum_{i = 1}^{10} Y_{i}^{2} - 2 μ \sum_{i = 2}^{11} Y_{i} - 2 α \sum_{i = 2}^{11} Y_{i} Y_{i - 1} + 2 μ α \sum_{i = 1}^{10} Y_{i}}{2 σ_{e}^{2}} - \frac{μ^{2}}{2 σ_{p}^{2}} - \frac{α^{2}}{2 σ_{p}^{2}}] d μ d α}

$=\frac{\exp{\biggl [ -\frac{10\mu^2+\alpha^2\sum_{i=1}^{10}Y_{i}^{2}-2\mu\sum_{i=2}^{11}Y_i-2\alpha\sum_{i=2}^{11}Y_{i}Y_{i-1}+2\mu\alpha\sum_{i=1}^{10}Y_i}{2\sigma_{e}^{2}} - \frac{\mu^{2}}{2\sigma_{p}^{2}} - \frac{\alpha^{2}}{2\sigma_{p}^{2}} \biggr ] }}{\int\int \exp{\biggl [ -\frac{10\mu^2+\alpha^2\sum_{i=1}^{10}Y_{i}^{2}-2\mu\sum_{i=2}^{11}Y_i-2\alpha\sum_{i=2}^{11}Y_{i}Y_{i-1}+2\mu\alpha\sum_{i=1}^{10}Y_i}{2\sigma_{e}^{2}} - \frac{\mu^{2}}{2\sigma_{p}^{2}} - \frac{\alpha^{2}}{2\sigma_{p}^{2}} \biggr ] }d\mu d\alpha}$

Teraz możemy użyć ogólnego wyniku z normalnego pliku pdf.

\int \exp (- a z^{2} + b z - c) d z = \sqrt{\frac{π}{a}} \exp (\frac{b^{2}}{4 a} - c)

$\int \exp\left(-az^2+bz-c\right)dz=\sqrt{\frac{\pi}{a}}\exp\left(\frac{b^2}{4a}-c\right)$ Wynika to z wypełnienia kwadratu na i zauważenia, że nie zależy od . Zauważ, że wewnętrzna całka nad ma tę postać z i i . Po wykonaniu tej całki okaże się, że pozostała całka na

- a z^{2} + b z

$-az^2+bz$

c

$c$

z

$z$

μ

$\mu$

a = \frac{10}{2 σ_{e}^{2}} + \frac{1}{2 σ_{p}^{2}}

$a=\frac{10}{2\sigma^2_e}+\frac{1}{2\sigma^2_p}$

b = \frac{\sum_{i = 2}^{11} Y_{i} - α \sum_{i = 1}^{10} Y_{i}}{σ_{e}^{2}}

$b=\frac{\sum_{i=2}^{11}Y_i-\alpha\sum_{i=1}^{10}Y_i}{\sigma_{e}^{2}}$

c = \frac{α^{2} \sum_{i = 1}^{10} Y_{i}^{2} - 2 α \sum_{i = 2}^{11} Y_{i} Y_{i - 1}}{2 σ_{e}^{2}} + \frac{α^{2}}{2 σ_{p}^{2}}

$c=\frac{\alpha^2\sum_{i=1}^{10}Y_{i}^{2}-2\alpha\sum_{i=2}^{11}Y_{i}Y_{i-1}}{2\sigma_{e}^{2}}+ \frac{\alpha^{2}}{2\sigma_{p}^{2}}$

α

$\alpha$ ma również tę postać, więc możesz ponownie użyć tej formuły, używając innej litery . Powinieneś być w stanie napisać swój późniejszy w postaci gdzie jest macierzą

a, b, c

$a,b,c$

\frac{1}{2 π | V |^{\frac{1}{2}}} \exp [- \frac{1}{2} (μ - \hat{μ}, α - \hat{α}) V^{- 1} (μ - \hat{μ}, α - \hat{α})^{T}]

$\frac{1}{2\pi|V|^{\frac{1}{2}}}\exp\left[-\frac{1}{2}(\mu-\hat{\mu},\alpha-\hat{\alpha})V^{-1}(\mu-\hat{\mu},\alpha-\hat{\alpha})^T\right]$

V

$V$

2 \times 2

$2\times 2$

Daj mi znać, jeśli potrzebujesz więcej wskazówek.

aktualizacja

(uwaga: poprawna formuła powinna wynosić zamiast ) $10\mu^2$ $\mu^2$

jeśli spojrzymy na kwadratową formę, którą napisaliście w aktualizacji, zauważymy, że istnieje współczynników ( ma znaczenia dla a posteriori, ponieważ zawsze możemy dodać dowolną stałą, która anuluje w mianowniku). Mamy też niewiadomych . Jest to zatem „dobrze postawiony” problem, o ile równania są liniowo niezależne. Jeśli rozwiniemy kwadratowy otrzymujemy: $5$ $L$ $5$ $\hat{\mu},\hat{\alpha},Q_{11},Q_{12}=Q_{21},Q_{22}$ $(\mu-\hat{\mu},\alpha-\hat{\alpha})Q(\mu-\hat{\mu},\alpha-\hat{\alpha})^{t}$

Q_{11} (μ - \hat{μ})^{2} + Q_{22} (α - \hat{α})^{2} + 2 Q_{12} (μ - \hat{μ}) (α - \hat{α})

$Q_{11}(\mu-\hat{\mu})^2+Q_{22}(\alpha-\hat{\alpha})^2+2Q_{12}(\mu-\hat{\mu})(\alpha-\hat{\alpha})$

= Q_{11} μ^{2} + 2 Q_{21} μ α + Q_{22} α^{2} - (2 Q_{11} \hat{μ} + 2 Q_{12} \hat{α}) μ - (2 Q_{22} \hat{α} + 2 Q_{12} \hat{μ}) α +

$=Q_{11}\mu^{2} + 2Q_{21}\mu\alpha + Q_{22}\alpha^{2} - (2Q_{11}\hat{\mu}+2Q_{12}\hat{\alpha})\mu - (2Q_{22}\hat{\alpha}+2Q_{12}\hat{\mu})\alpha +$

+ Q_{11} {\hat{μ}}^{2} + Q_{22} {\hat{α}}^{2} + 2 Q_{12} \hat{μ} \hat{α}

$+Q_{11}\hat{\mu}^2+Q_{22}\hat{\alpha}^2+2Q_{12}\hat{\mu}\hat{\alpha}$

Porównując współczynnik drugiego rzędu, otrzymujemy co mówi nam, jak wygląda (odwrotna) macierz kowariancji. Również dwie nieznacznie bardziej skomplikowane równania dla po podstawieniu do . Można je zapisać w postaci macierzy jako: $A=Q_{11},B=2Q_{12},C=Q_{22}$ $\hat{\alpha},\hat{\mu}$ $Q$

- (\begin{matrix} 2 A & B \\ B & 2 C \end{matrix}) (\begin{matrix} \hat{μ} \\ \hat{α} \end{matrix}) = (\begin{matrix} J \\ K \end{matrix})

$-\begin{pmatrix}2A & B \\ B & 2C\end{pmatrix} \begin{pmatrix}\hat{\mu} \\ \hat{\alpha}\end{pmatrix} = \begin{pmatrix}J \\ K\end{pmatrix}$

Tak więc szacunki są podane przez:

(\begin{matrix} \hat{μ} \\ \hat{α} \end{matrix}) = - {(\begin{matrix} 2 A & B \\ B & 2 C \end{matrix})}^{- 1} (\begin{matrix} J \\ K \end{matrix}) = \frac{1}{4 A C - B^{2}} (\begin{matrix} B K - 2 J C \\ B J - 2 K A \end{matrix})

$\begin{pmatrix}\hat{\mu} \\ \hat{\alpha}\end{pmatrix} = -\begin{pmatrix}2A & B \\ B & 2C\end{pmatrix}^{-1}\begin{pmatrix}J \\ K\end{pmatrix}=\frac{1}{4AC-B^2}\begin{pmatrix}BK-2JC \\ BJ-2KA\end{pmatrix}$

Pokazuje, że nie mamy unikalnych szacunków, chyba że . Teraz mamy: $4AC\neq B^2$

\begin{array}{cc} A = \frac{10}{2 σ_{e}^{2}} + \frac{1}{2 σ_{p}^{2}} & B = \frac{\sum_{i = 1}^{10} Y_{i}}{σ_{e}^{2}} & C = \frac{\sum_{i = 1}^{10} Y_{i}^{2}}{2 σ_{e}^{2}} + \frac{1}{2 σ_{p}^{2}} \\ J = - \frac{\sum_{i = 2}^{11} Y_{i}}{σ_{e}^{2}} & K = - \frac{\sum_{i = 2}^{11} Y_{i} Y_{i - 1}}{σ_{e}^{2}} \end{array}

$\begin{array}{c c} A=\frac{10}{2\sigma^2_e}+\frac{1}{2\sigma^2_p} & B=\frac{\sum_{i=1}^{10}Y_i}{\sigma_{e}^{2}} & C=\frac{\sum_{i=1}^{10}Y_{i}^{2}}{2\sigma^2_e}+\frac{1}{2\sigma^2_p} \\ J=-\frac{\sum_{i=2}^{11}Y_i}{\sigma_{e}^{2}} & K=-\frac{\sum_{i=2}^{11}Y_{i}Y_{i-1}}{\sigma_{e}^{2}} \end{array}$

Zauważ, że jeśli zdefiniujemy dla i przyjmiemy limit wówczas szacunki dla są podane przez zwykłe najmniejsze kwadraty oszacowanie i gdzie i . Tak więc szacunki tylne są średnią ważoną między szacunkami OLS a wcześniejszymi szacunkami . $X_i=Y_{i-1}$ $i=2,\dots,11$ $\sigma^2_p\to\infty$ $\mu,\alpha$ $\hat{\alpha}=\frac{\sum_{i=2}^{11}(Y_{i}-\overline{Y})(X_{i}-\overline{X})}{\sum_{i=2}^{11}(X_{i}-\overline{X})^2}$ $\hat{\mu}=\overline{Y}-\hat{\alpha}\overline{X}$ $\overline{Y}=\frac{1}{10}\sum_{i=2}^{11}Y_i$ $\overline{X}=\frac{1}{10}\sum_{i=2}^{11}X_i=\frac{1}{10}\sum_{i=1}^{10}Y_i$ $(0,0)$

— prawdopodobieństwo prawdopodobieństwa
źródło

Nie jest to szczególnie pomocne, ponieważ wspomniałem konkretnie, że nie chodzi tu o mianownik. Mianownik jest tylko stałą normalizującą, co stanie się oczywiste po zredukowaniu licznika do postaci Gaussa. Więc sztuczki do obliczania całek w mianowniku są matematycznie naprawdę fajne, ale po prostu nie są potrzebne w mojej aplikacji. Jedynym problemem, z którym potrzebuję rozwiązania, jest manipulowanie licznikiem.

— ely

Ta odpowiedź daje zarówno licznik, jak i mianownik. Licznik wykazuje prawidłowy wielomian drugiego stopnia w który prowadzi do normalnej postaci kwadratowej, co jest podkreślone przez prawdopodobieństwo logiczne.

(α, μ)

$(\alpha,\mu)$

— Xi'an

@ems - obliczając stałą normalizującą zbudujesz wymaganą formę kwadratową. będzie zawierać warunki potrzebne do wypełnienia kwadratu

— prawdopodobieństwo

Nie rozumiem, jak to daje ci kwadratową formę. Obliczyłem dwie całki w mianowniku, używając opublikowanej przez ciebie tożsamości całkowej Gaussa. W końcu mam po prostu ogromną, niechlujną stałą. Wydaje się, że nie ma wyraźnego sposobu na przyjęcie tej stałej i przekształcenie jej w coś, co jest wyznacznikiem mocy 1/2, itd. Nie wspominając, że nie rozumiem, jak to wszystko wyjaśnia, jak obliczyć nową „ mean vector ' .. O to prosiłem o pomoc w pierwotnym pytaniu.

(\hat{μ}, \hat{α})^{t}

$(\hat{\mu},\hat{\alpha})^{t}$

— ely

Ogromne dzięki za szczegółowe dodanie. Popełniłem głupie błędy, próbując wykonać algebrę, aby obliczyć formę kwadratową. Twoje komentarze na temat relacji do estymatora OLS są bardzo interesujące i docenione. Myślę, że to przyspieszy mój kod, ponieważ będę mógł pobierać próbki z formy analitycznej, która ma wbudowane, zoptymalizowane metody. Mój pierwotny plan polegał na użyciu Metropolis-Hastings do pobierania próbek, ale było to bardzo powolne. Dzięki!

— ely