# RJMCMC

## Bayesian Model Choice

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-Lim_U4nX75SgwH6YoWN%2F-Lim_WA4hIjTk3YVEOoR%2Fdef-1.png?generation=1562069310042792\&alt=media)

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-Lim_U4nX75SgwH6YoWN%2F-Lim_WA6AijfATs8oYGU%2Fdef-2.png?generation=1562069310044904\&alt=media)

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-Lim_U4nX75SgwH6YoWN%2F-Lim_WA8UIatKto9ApE9%2Fdef-3.png?generation=1562069310157102\&alt=media)

## Green's algorithm

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-Lim_U4nX75SgwH6YoWN%2F-Lim_WAAiXdCluMv8vgt%2Falg-green.png?generation=1562069310260330\&alt=media)

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-Lim_U4nX75SgwH6YoWN%2F-Lim_WACHMhz0esMXxki%2Flinear-jacobian.png?generation=1562069309943853\&alt=media)

## Fixed Dimension Reassessment

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-Lim_U4nX75SgwH6YoWN%2F-Lim_WAEtxxDjPhynJ_V%2Ffixed-dims-reassess.png?generation=1562069310037922\&alt=media)

## RJMCMC for change points

The log-likelihood function is

$$
\begin{aligned}
p(y\mid x) &= \sum\_{i=1}^n\log{x(y\_i)}-\int\_0^Lx(t)dt\\
&= \sum\_{j=1}^km\_j\log h\_j-\sum\_{j=0}^kh\_j(s\_{j+1}-s\_j)
\end{aligned}
$$

where $$m\_j=#{y\_i\in\[s\_j,s\_{j+1})}$$.

### H Move

1. choose one of $$h\_0,h\_1,\ldots,h\_k$$ at random, obtaining $$h\_j$$
2. propose a change to $$h\_j'$$ such that $$\log(h\_j'/h\_j)\sim U\[-\frac 12, \frac 12]$$.
3. the acceptance probability is 

   $$
   \begin{aligned}
   \alpha(h\_j, h\_j') &=\frac{p(h\_j'\mid y)}{p(h\_j\mid y)}\times \frac{J(h\_j\mid h\_j')}{J(h\_j'\mid h\_j)}\\
   &= \frac{p(y\mid h\_j')}{p(y\mid h\_j)}\times\frac{\pi(h\_j')}{\pi(h\_j)}\times\frac{J(h\_j\mid h\_j')}{J(h\_j'\mid h\_j)},.
   \end{aligned}
   $$

   Note that 

   $$
   \log\frac{h\_j'}{h\_j}= u ,,
   $$

   it follows that the CDF of $$h\_j'$$ is 

   $$
   F(x) = P(h\_j'\le x) = P(e^uh\le x) = P(u\le \log(x/h)) = \log(x/h) + 1/2
   $$

   and 

   $$
   f(x) = F'(x) = \frac{1}{x},,
   $$

   so 

   $$
   J(h\_j'\mid h\_j) = \frac{1}{h\_j'},.
   $$

   Then we have

   $$
   \begin{aligned}
   \alpha(h\_j,h\_j') &= \frac{p(y\mid h\_j')}{p(y\mid h\_j)}\times \frac{(h\_j')^{\alpha}\exp(-\beta h\_j')}{h\_j^\alpha\exp(-\beta h\_j)}\\
   &=\text{likelihood ratio}\times (h\_j'/h\_j)^\alpha\exp{-\beta(h\_j'-h\_j)}
   \end{aligned}
   $$

### P Move

1. Draw one of $$s\_1,s\_2,\ldots,s\_k$$ at random, obtaining say $$s\_j$$.
2. Propose $$s\_j'\sim U\[s\_{j-1}, s\_{j+1}]$$
3. The acceptance probability is 

   $$
   \begin{aligned}
   \alpha(s\_j,s\_j') &=\frac{p(y\mid s\_j')}{p(y\mid s\_j)}\times \frac{\pi(s\_j')}{\pi(s\_j)}\times \frac{J(s\_j\mid s\_j')}{J(s\_j'\mid s\_j)}\\
   &=\text{likelihood ratio}\times \frac{(s\_{j+1}-s\_j')(s\_j'-s\_{j-1})}{(s\_{j+1}-s\_j)(s\_j-s\_{j-1})}
   \end{aligned}
   $$

   since 

   $$
   \pi(s\_1,s\_2,\ldots,s\_k)=\frac{(2k+1)!}{L^{2k+1}}\prod\_{j=0}^{k+1}(s\_{j+1}-s\_j)
   $$

### Birth Move

Choose a position $$s^\*$$ uniformly distributed on $$\[0,L]$$, which must lie within an existing interval $$(s\_j,s\_{j+1})$$ w\.p 1.

Propose new heights $$h\_j', h\_{j+1}'$$ for the step function on the subintervals $$(s\_j,s^*)$$ and $$(s^*,s\_{j+1})$$. Use a weighted geometric mean for this compromise,

$$
(s^*-s\_j)\log(h\_j') + (s\_{j+1}-s^*)\log(h\_{j+1}')=(s\_{j+1}-s\_j)\log(h\_j)
$$

and define the perturbation to be such that

$$
\frac{h\_{j+1}'}{h\_j'}=\frac{1-u}{u}
$$

with $$u$$ drawn uniformly from $$\[0,1]$$.

The prior ratio, becomes

$$
\frac{p(k+1)}{p(k)}\frac{2(k+1)(2k+3)}{L^2}\frac{(s^*-s\_j)*(s\_{j+1}-s^\*)}{s\_{j+1}-s\_j}\times \frac{\beta^\alpha}{\Gamma(\alpha)}\Big(\frac{h\_j'h\_{j+1}'}{h\_j}\Big)^{\alpha-1}\exp{-\beta(h\_j'+h\_{j+1}'-h)}
$$

the proposal ration becomes

$$
\frac{d\_{k+1}L}{b\_k(k+1)}
$$

and the Jacobian is

$$
\frac{(h\_j'+h\_{j+1}')^2}{h\_j}
$$

### Death Move

If $$s\_{j+1}$$ is removed, the new height over the interval $$(s\_j',s\_{j+1}')=(s\_j,s\_{j+2})$$ is $$h\_j'$$, the weighted geometric mean satisfying

$$
(s\_{j+1}-s\_j)\log(h\_j) + (s\_{j+2}-s\_{j+1})\log(h\_{j+1}) = (s\_{j+1}'-s\_j')\log(h\_j')
$$

The acceptance probability for the corresponding death step has the same form with the appropriate change of labelling of the variables, and the ratio terms inverted.

## Results

The histogram of number of change points is

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-LW1wDlIcpoh2WLMgd80%2F-LW1wFon-Nmly_3eMqaK%2Fhistogram_num_of_change_points.png?generation=1547311321624518\&alt=media)

And we can get the density plot of position (or height) conditional on the number of change points $$k$$, for example, the following plot is for the position when $$k=1,2,3$$

![](https://666993855-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-LJfsESZOIJn_3uGIecs%2F-LW1wDlIcpoh2WLMgd80%2F-LW1wFotVURcxiBLO7Ry%2Fdensity_plot_position.png?generation=1547311321662238\&alt=media)

## References

1. [Green, P. J. (1995). Reversible jump Markov chain Monte Carlo computation and Bayesian model determination. Biometrika, 82(4), 711–732.](https://doi.org/10.1093/biomet/82.4.711)
2. [Sisson, S. A. (2005). Transdimensional Markov Chains: A Decade of Progress and Future Perspectives. Journal of the American Statistical Association, 100(471), 1077–1089.](https://doi.org/10.1198/016214505000000664)
3. [Peter Green's Fortran program AutoRJ](https://people.maths.bris.ac.uk/~mapjg/AutoRJ/)
4. [David Hastie's C program AutoMix](http://www.davidhastie.me.uk/software/automix/)
5. [Ai Jialin, Reversible Jump Markov Chain Monte Carlo Methods. MSc thesis, University of Leeds, Department of Statistics, 2011/12.](http://www1.maths.leeds.ac.uk/~voss/projects/2011-RJMCMC/)