Как я могу вычислить апостериорную оценку плотности из априорной и вероятностной?

Я пытаюсь понять, как использовать теорему Байеса для вычисления апостериорного значения, но я застреваю с вычислительным подходом, например, в следующем случае мне не ясно, как взять произведение предыдущего и вероятности, а затем вычислить задний:

Для этого примера меня интересует вычисление апостериорной вероятности и я использую стандартную норму перед µ p ( µ ) ∼ N ( µ = 0 , σ = 1 ) , но я хочу знать, как рассчитать апостериорный из до μ, который представлен цепью MCMC, поэтому я буду использовать 1000 образцов в качестве отправной точки.$\mu$$\mu$ $p(\mu)\sim N(\mu = 0, \sigma = 1)$$\mu$

• образец 1000 от предшествующего.

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

• сделать некоторые наблюдения:

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

• и вычислить вероятность, например, :$p(y | \mu, \sigma)$

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

что я не совсем понимаю, это:

1. когда / как умножить предшествующее на вероятность?
2. когда / как нормализовать заднюю плотность?

пожалуйста, обратите внимание: меня интересует общее вычислительное решение, которое может быть обобщаемыми задачами без аналитического решения

Вы перепутали несколько вещей. Теория говорит о умножении предыдущего распределения и вероятности, а не выборок из предыдущего распределения. Также не ясно, что у вас есть априор, означает ли это априор? или что-то другое?

Тогда у вас все изменится по вероятности, ваши наблюдения должны быть x с предшествующими ничьями или известными фиксированными константами как средним и стандартным отклонением. И даже тогда это действительно будет продуктом 4 вызовов dnorm с каждым из ваших наблюдений как x и тем же средним и стандартным отклонением.

Что на самом деле не ясно, что вы пытаетесь сделать. Какой у Вас вопрос? какие параметры вас интересуют? что у вас есть по этим параметрам? Есть ли другие параметры? у вас есть приоры или фиксированные значения для них?

Попытки действовать так, как вы сейчас делаете, будут еще больше сбивать вас с толку, пока вы не решите точно, что ваш вопрос, и оттуда не поработаете.

Ниже добавляется должное после редактирования исходного вопроса.

Вам все еще не хватает некоторых частей, и, вероятно, не все понимают, но мы можем начать с того, где вы находитесь.

$\frac12$

Таким образом, мы можем начать аналогично тому, что вы сделали, и сгенерировать из предыдущего:

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

Теперь нам нужно вычислить вероятности, это основано на предыдущих извлечениях среднего значения, вероятности с данными и известной ценности SD. Функция dnorm даст нам вероятность отдельной точки, но нам нужно умножить вместе значения для каждого из наблюдений, вот функция для этого:

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

Теперь мы можем вычислить вероятность для каждого розыгрыша от предыдущего для среднего

> tmp <- likfun(pri)

Теперь, чтобы получить апостериор, нам нужно сделать новый тип розыгрыша. Один из подходов, подобный выборке отбраковки, состоит в том, чтобы отбирать образцы из предыдущих средних розыгрышей, пропорциональных вероятности для каждого предыдущего розыгрыша (это самый близкий к шагу умножения, которым вы были спрашивать о):

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

Теперь мы можем посмотреть на результаты задних розыгрышей:

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

и сравнить приведенные выше результаты со значениями в замкнутой форме из теории

> (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

Неплохое приближение, но, вероятно, будет лучше использовать встроенный инструмент McMC для рисования сзади. Большинство из этих инструментов выбирают одну точку за раз, а не партиями, как указано выше.

$\chi^2$

Общее решение состоит в том, чтобы использовать существующие инструменты для выполнения вычислений McMC, такие как WinBugs или OpenBugs (BRugs в R дает интерфейс между R и Bugs) или упаковывать такие LearnBayes в R.