I want to find sample from posterior distribution but it is complex, so i want to sample the posterior using MCMC. The posterior distribution contains integration (because the reliability function is not in closed form). My question is,How can i sample from this posterior. The full conditional posterior is not well known,so i canot to use Gibbs sampling.I will use Metropolis-Hastings algorithm but what i can do with the integration in posterior distribution. I want to use R programm. enter image description here

  • $\begingroup$ Depending on your specific case, you can use software like Stan for this kind of thing. It implements a probabilistic programming language that runs a form of Hamiltonian Monte-Carlo sampling $\endgroup$ – shadowtalker Jun 20 '17 at 3:27

Without further details about your model, your question remains unclear. The integral can usually be bypassed by introducing an auxiliary variable corresponding to the integrand of the integral. This means that, if the target density writes down as $$g(x)\int_\mathfrak{Y} f(x,y)\,\text{d}y$$ one can define a secondary target$$h(x,y)=g(x)f(x,y)$$and run an MCMC algorithm on that new target. (You can e.g. check the keywords demarginalisation and latent variable in our book.)

For instance, replace enter image description here

with a joint density on $(t_i,\phi_i)$ \begin{align*}(\cos\phi_i)^{\beta_1/\delta-1}&(\sin\phi_i)^{\beta_2/\delta-1}\\ &\times\left[ (t_i\cos\phi_i/\theta_1)^{\beta_1/\delta}+(t_i\sin\phi_i/\theta_2)^{\beta_2/\delta}\right]^{\delta-2}\\ &\times\left[ (t_i\cos\phi_i/\theta_1)^{\beta_1/\delta}+(t_i\sin\phi_i/\theta_2)^{\beta_2/\delta}+\frac{1}{\delta}-1\right] \end{align*} that need be integrated within the target of the MCMC algorithm.

| cite | improve this answer | |
  • $\begingroup$ Dear, I attach the full conditional posterior for one parameter $\endgroup$ – Manal Jun 20 '17 at 12:11
  • $\begingroup$ For each integral in your expression, introduce a latent variable and remove the integral. This means using (n-k) latent variables when an integral is set to the power (n-k). $\endgroup$ – Xi'an Jun 21 '17 at 6:21
  • $\begingroup$ Thanx for respond,do you have reference (l hope be easy). I did not study it before. $\endgroup$ – Manal Jun 21 '17 at 14:31
  • $\begingroup$ Do you have a good and easy to understand reference? $\endgroup$ – Manal Jun 22 '17 at 10:23
  • $\begingroup$ George Casella and I wrote an iintroduction to Monte Carlo methods a few years ago. There is no "easy" reference in that you have to understand the notion of demarginalisation and latent variables to operate. $\endgroup$ – Xi'an Jun 22 '17 at 10:27

If you can, it may be easiest to go with the solution pointed out by Xi’an. If that is not an option, then here are some alternatives.

Some programming languages have an integration routine in-buildt and accessible from within their Bayesian procedures (e.g. SAS allows you to use adaptive quadrature via the CALL QUAD statement from within PROC MCMC for a user defined function defined in a preceding PROC FCMP).

If that is not the case, then, if your programming language allows you to, you can program your own integration routine. Adaptive quadrature may be an obvious candidate, because it is not very challenging to program. In some MCMC sampling schemes that may not be attractive (e.g. for Hamiltonian Monte Carlo, I believe adaptive quadrature would cause issues within derivatives due to the parameter value dependent number of terms). If you can get bounds on the error terms, then non-adaptive quadrature may do the job, but that is often challenging. For Hamiltonian Monte Carlo we may well be talking about Stan, and in that case we could also use the ordinary differential equation solver it provides to do our integration.

| cite | improve this answer | |
  • $\begingroup$ Dear, I attach the full conditional posterior for one parameter $\endgroup$ – Manal Jun 20 '17 at 12:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.