Jak obliczyć oszacowanie gęstości tylnej na podstawie wcześniejszego i prawdopodobnego prawdopodobieństwa?

Próbuję zrozumieć, jak użyć twierdzenia Bayesa do obliczenia tylnej, ale utknąłem w podejściu obliczeniowym, np. W następującym przypadku nie jest dla mnie jasne, jak wziąć iloczyn wcześniejszego i prawdopodobieństwa, a następnie obliczyć tylny:

W tym przykładzie jestem zainteresowany obliczeniem prawdopodobieństwa a posteriori dla i używam standardowej normalnej przed , ale chcę wiedzieć jak obliczyć tylną z wcześniejszego reprezentowanego przez łańcuch MCMC, więc użyję 1000 próbek jako mojego punktu wyjścia.$\mu$$\mu$ $p(\mu)\sim N(\mu = 0, \sigma = 1)$$\mu$

• próbka 1000 z poprzedniego.

  set.seed(0)
  prior.mu      <- 0
  prior.sigma   <- 1
  prior.samples <- sort(rnorm(1000, prior.mu, prior.sigma))
  

• poczyń kilka obserwacji:

  observations <- c(0.4, 0.5, 0.8, 0.1)
  

• i obliczyć prawdopodobieństwo, np $p(y | \mu, \sigma)$:

  likelihood <- prod(dnorm(observations, mean(prior.samplse), sd(prior.samples)))
  

nie do końca rozumiem to:

1. kiedy / jak pomnożyć wcześniej przez prawdopodobieństwo?
2. kiedy / jak znormalizować tylną gęstość?

Uwaga: Jestem zainteresowany ogólnym rozwiązaniem obliczeniowym, które może być problemami uogólnionymi bez rozwiązania analitycznego

Masz kilka pomieszanych rzeczy. Teoria mówi o pomnożeniu wcześniejszej dystrybucji i prawdopodobieństwa, a nie próbek z wcześniejszej dystrybucji. Nie jest też jasne, co masz przed przeorem, czy jest to przeor na coś innego? albo coś innego?

Następnie masz odwrócone prawdopodobieństwo, twoje obserwacje powinny być x z wcześniejszymi losowaniami lub znanymi stałymi stałymi jako średnią i odchylenie standardowe. I nawet wtedy naprawdę byłby to wynik 4 wezwań do zanormowania przy każdej z twoich obserwacji jako x i przy tej samej średniej i standardowym odchyleniu.

Co tak naprawdę nie jest jasne, co próbujesz zrobić. Jakie jest Twoje pytanie? jakie parametry jesteś zainteresowany? jakie masz uprzednie prawa do tych parametrów? czy są inne parametry? czy masz dla nich priorytety lub stałe wartości?

Próba zajmowania się takimi, jakimi jesteś obecnie, tylko bardziej Cię dezorientują, dopóki nie rozwiążesz dokładnie tego, jakie jest twoje pytanie i zaczniesz odtąd działać.

Poniżej dodano ze względu na edycję oryginalnego pytania.

Wciąż brakuje ci niektórych elementów i prawdopodobnie nie rozumiesz wszystkiego, ale możemy zacząć od miejsca, w którym jesteś.

Myślę, że mylisz kilka pojęć. Istnieje prawdopodobieństwo, że pokazuje związek między danymi a parametrami, używasz normalnej, która ma 2 parametry, średnią i odchylenie standardowe (lub wariancję lub precyzję). Potem są wcześniejsze rozkłady parametrów, podałeś normalny przełożony ze średnią 0 i sd 1, ale ta średnia i odchylenie standardowe są całkowicie różne od średniej i odchylenia standardowego prawdopodobieństwa. Aby być kompletnym, musisz albo znać prawdopodobieństwo SD, albo umieścić uprzednio na SD prawdopodobieństwa, dla uproszczenia (ale mniej realnego) założę, że wiemy, że prawdopodobieństwo SD to$\frac12$ (nie ma żadnego ważnego powodu innego niż to działa i różni się od 1).

Możemy więc zacząć podobnie do tego, co zrobiłeś i wygenerować z wcześniejszego:

> obs <- c(0.4, 0.5, 0.8, 0.1)
> pri <- rnorm(10000, 0, 1)

Teraz musimy obliczyć prawdopodobieństwa, jest to oparte na wcześniejszych losowaniach średniej, prawdopodobieństwie z danymi i znanej wartości SD. Funkcja dnorm da nam prawdopodobieństwo jednego punktu, ale musimy pomnożyć razem wartości dla każdej z obserwacji, oto funkcja, która to robi:

> likfun <- function(theta) {
+ sapply( theta, function(t) prod( dnorm(obs, t, 0.5) ) )
+ }

Teraz możemy obliczyć prawdopodobieństwo dla każdego losowania na podstawie wyniku średniego

> tmp <- likfun(pri)

Teraz, aby uzyskać tył, musimy zrobić nowy typ losowania, jednym z podejść podobnych do próbkowania odrzucenia jest pobranie próbek z poprzednich średnich losowań proporcjonalnych do prawdopodobieństwa dla każdego poprzedniego losowania (jest to najbliższy etap pomnożenia, w którym byłeś pytać o):

> post <- sample( pri, 100000, replace=TRUE, prob=tmp )

Teraz możemy spojrzeć na wyniki losowań tylnych:

> mean(post)
[1] 0.4205842
> sd(post)
[1] 0.2421079
> 
> hist(post)
> abline(v=mean(post), col='green')

i porównaj powyższe wyniki z wartościami postaci zamkniętej z teorii

> (1/1^2*mean(pri) + length(obs)/0.5^2 * mean(obs))/( 1/1^2 + length(obs)/0.5^2 )
[1] 0.4233263
> sqrt(1/(1+4*4))
[1] 0.2425356

Nie jest to złe przybliżenie, ale prawdopodobnie lepiej byłoby użyć wbudowanego narzędzia McMC do rysowania od tyłu. Większość z tych narzędzi próbkuje jeden punkt naraz, a nie w partiach takich jak powyżej.

Mówiąc bardziej realistycznie, nie znalibyśmy SD prawdopodobieństwa i potrzebowalibyśmy także wcześniejszego (często wcześniejszy wariant jest $\chi^2$ lub gamma), ale obliczenia są bardziej skomplikowane (przydaje się McMC) i nie ma zamkniętej formy do porównania.

Ogólnym rozwiązaniem jest użycie istniejących narzędzi do wykonywania obliczeń McMC, takich jak WinBugs lub OpenBugs (BRugs in R zapewnia interfejs między R i Bugs) lub pakietów takich LearnBayes w R.