1

在抛硬币中，我们想计算 p(theta|Data)，其中 theta 是基础参数。

  • 先验遵循带有参数 a 和 b 的 beta 分布。
  • 可能性遵循伯努利分布，这给了我们出现正面的概率。

下面是代码实现：

    a = 1  # a and b are the beta distribution's parameters
    b= 1
    num = 1e5 #Number of candidate theta values
    z= 17220 #Number of heads
    N= 143293 #Total number of flips

    Theta = seq(0.07,0.12, length.out= num)
    prior = dbeta(Theta, a,b) #Compute the prior at each value 

    likelihood = Theta^z *(1-Theta)^(N-z)

    pData = likelihood * prior /sum(likelihood * prior) #Compute evidence
    posterior = likelihood*prior / pData

我想验证后验是否等于解析解 beta(a+z, N-z+b)。但是，由于由于θ值很小，似然性等于0，因此证据的概率是Nan，后验概率也是如此。

我已经尝试计算对数似然度，但它给了我一个很大的负数，在取指数时等于 0。

 Theta = seq(0.07,0.12, by= num_steps)
 lprior = log(dbeta(Theta, a,b)) #Compute the log prior at each value 

 llikelihood = log(Theta)*z + log(1-Theta)*(N-z) #log likelihood

 lpData = llikelihood + lprior - sum(llikelihood + lprior) #compute evidence

  lposterior = llikelihood+lprior - lpData
  posterior = exp(lposterior)
  plot(Theta, posterior, type="l")
  lines(Theta, exp(llikelihood), type="l")
  lines(Theta, exp(lprior), type="l")

如果我的最终目标是有一个漂亮的图表来显示后验、可能性和先验，就像这样

这个

我应该如何计算每个值？

4

1 回答 1

0

这个答案是由@JosephClarkMcIntyre 通过评论部分提供的。

这是一个摘要：

  • 在伯努利试验中，当 N - 试验的总数 - 和 z - 成功的总数很大并且基础参数 theta 很小时，最好只在对数空间中操作，而不要取指数。
  • 此外，由于对数函数在增加，比较两个分布的对数后验类似于比较后验。
  • 上述实现是错误的，因为计算证据的公式不正确。p(证据) = sum(likelihood*prior), p(log_evidence)= sum(log_likelihood +log_prior)

这是最终代码，其中先验、可能性和证据都在日志空间中：

  a = 1  # a and b are the beta distribution's paramteres
  b= 1
  num_steps = 1e5
  z= 17220 #Number of heads
  N= 143293 #Total number of flips

  Theta = seq(from=0.07,to=0.12, length.out= num_steps)
  lprior = dbeta(Theta, a,b,log=TRUE) #Compute the log prior at each value 

  llikelihood = log(Theta)*z + log1p(-Theta)*(N-z) #log likelihood

  lpData = sum(llikelihood + lprior) #compute log of the evidence

  lposterior = llikelihood+lprior - lpData
  plot(Theta,log(dbeta(Theta,a+z,N-z+b)))
  plot(Theta, lposterior, type="l")

但是，分析和计算的对数后验与图中所示不同。

如果您认为此答案存在缺陷，请随时发表评论或解释为什么分析和计算的对数后验不一样。^^

于 2019-01-17T08:57:43.300 回答