Aperçu des sections

  • Bayesian parameter inference

  • Bayesian model choice

  • The case of linear regression models

  • Simulation of random variables

    • Exponential Distribution Simulation

      The code simulates samples from an exponential distribution with a rate parameter of 2 using inverse transform sampling. It then plots a histogram of the generated data and overlays the theoretical exponential density using curve().

      Normal Distribution Simulation

      The code generates samples from a standard normal distribution using the Box-Muller transform. It creates a histogram and overlays the theoretical normal density.

      Multivariate Normal Simulation

      The code simulates a bivariate normal distribution with a specified mean and covariance matrix using the Cholesky decomposition method. The resulting samples are then plotted.

      Gamma Distribution using Acceptance-Rejection

      This part simulates samples from a Gamma(2,1) distribution using the acceptance-rejection method. The optimal instrumental distribution used is exponential with rate 1/2. The code calculates and displays the efficiency of the algorithm (i.e., the number of rejections).

      The acceptance-rejection method is repeated using a different instrumental distribution, Exp(1/10), and the efficiency is recalculated.

      Target Distribution (Modified Gamma)

      The final part targets a more complex distribution: x * exp(-x) * (1 + cos^2(x)). The acceptance-rejection method is used again, and the code plots the target distribution along with the instrumental exponential distribution for comparison.


    • # solution

      ptilde <- function(x){0.7*exp(-x^2/2)+0.3*exp(-x^2/2+7*(x-7/2))}
      ratio <- function(x){ptilde(x)/dnorm(x,2.1,3)}
      curve(ratio,-10,20,ylim=c(0,12))
      abline(a=11,b=0)


    • This code implements Bayesian inference using the rejection sampling method to simulate from the posterior distribution. The problem involves estimating the posterior distribution for a model with a normal likelihood and a Cauchy prior

      Likelihood xθN(θ,1)x \mid \theta \sim N(\theta, 1) and x=0.8x = 0.8

      Prior θCauchy(0,1)\theta \sim \text{Cauchy} which is a heavy-tailed distribution with no finite mean or variance

      Target The posterior distribution is proportional to the product of the likelihood and the prior, which is intractable

      We approximate the posterior distribution by generating samples using rejection sampling and compute the posterior expectation (Bayesian estimate) using a Monte Carlo approximation

      Target Function target.tilde(theta):

      Represents the unnormalized posterior: likelihood exp(xθθ2/2)\exp(x \cdot \theta - \theta^2/2) times the prior 11+θ2\frac{1}{1 + \theta^2} (Cauchy)

      Proposal Distribution proposal.tilde(theta)

      The proposal distribution is the Cauchy distribution 11+θ2\frac{1}{1 + \theta^2} which is a good candidate due to the prior also being Cauchy

      Rejection Sampling Procedure

      A sample yy is drawn from the Cauchy distribution using rt(1, 1)

      A uniform random variable u[0,1]u \in [0, 1] is generated

      The acceptance probability ρ\rho is computed as the ratio of the target to the proposal, scaled by M.tildeM.tilde

      If uρu \leq \rho, the sample is accepted and stored in simu[i]; otherwise, the loop continues until a valid sample is accepted

      nsimu: Keeps track of the total number of simulations attempted (including rejected samples)

      Monte Carlo Estimate

      After generating the samples, the posterior expectation (or any other statistic of interest) can be approximated using the empirical mean of the simulated values stored in simu

      This method leverages rejection sampling and the Monte Carlo approach to perform Bayesian inference for models where the posterior distribution cannot be computed exactly

  • Monte Carlo and Monte Carlo Markov Chains methods

  • Bayesian essentials with R