Gibbs Samplingis a Monte Carlo Markov Chain method for estimating complex joint distributions that draw an instance from the distribution of each variable iteratively based on the current values of the other variables.
The following comment shall be used to determine the posterior (parameters) hyperparameters required for the simulation
Let\({x_{\scriptstyle1\atop\scriptstyle\,}},{x_{\scriptstyle2\atop\scriptstyle\,}},.....{x_{\scriptstylen\atop\scriptstyle\,}}\) be i.i.d. random variables from the normal distribution with unknown mean and unknown variance \({\sigma ^2} > 0.\)
Suppose that for the prior joint distribution, the conditional distribution of μ given\(\tau = {\raise0.7ex\hbox{$1$} \!\mathord{\left/ {\vphantom {1 {{\sigma ^2}}}}\right.\\} \!\lower0.7ex\hbox{${{\sigma ^2}}$}} > 0\) is from the normal distribution with mean \({\mu _{\scriptstyle0\atop\scriptstyle\,}}\)and precision\({\lambda _{\scriptstyle0\atop\scriptstyle\,}}\tau > 0\) and that the marginal distribution of \(\tau \)is the gamma distribution with parameters \({\alpha _{\scriptstyle0\atop\scriptstyle\,}} > 0\,and\,\,{\beta _{\scriptstyle0\atop\scriptstyle\,}} > 0\,\) .
Then for the joint posterior distribution, the conditional distribution of μ given\(\tau = 1/{\sigma ^2} > 0\) is from the normal distribution with mean\({\mu _{\scriptstyle1\atop\scriptstyle\,}}\) and precision\({\lambda _{\scriptstyle1\atop\scriptstyle\,}}\tau > 0\)
The marginal distribution of τ is the gamma distribution with parameters \({\alpha _{\scriptstyle0\atop\scriptstyle\,}} > 0\,and\,\,{\beta _{\scriptstyle0\atop\scriptstyle\,}} > 0\,\) where
\({\beta _{\scriptstyle1\atop\scriptstyle\,}} = {\beta _{\scriptstyle0\atop\scriptstyle\,}} + \frac{1}{2}{s_{\scriptstylen\atop\scriptstyle\,}}^2 + {\frac{{n{\lambda _{\scriptstyle0\atop\scriptstyle\,}}\left( {\overline {{x_{\scriptstylen\atop\scriptstyle\,}}} - {\mu _{\scriptstyle0\atop\scriptstyle\,}}} \right)}}{{2\left( {{\lambda _{\scriptstyle0\atop\scriptstyle\,}} + n} \right)}}^2}.\)
Next, get the data from the two mentioned examples. There is a total of observations. Recall the following product of the likelihood and the prior distribution
\(\xi \left( {\mu ,\tau \mid x} \right)\infty \left\{ {{\tau ^{1/2}}{\rm{ }}exp\left( { - \frac{1}{2}{\lambda _{\scriptstyle1\atop\scriptstyle\,}}\tau {{\left( {\mu - {\mu _{\scriptstyle1\atop\scriptstyle\,}}} \right)}^2}} \right)} \right\}{\tau ^{{\alpha _{\scriptstyle1\atop\scriptstyle\,}}}}\exp \left( { - {\beta _{\scriptstyle1\atop\scriptstyle\,}}{\rm{ }}\tau } \right).\)
This equation may also be written as
\(\xi \left( {\mu ,\tau \mid x} \right)\infty \,\,{\tau ^{{\alpha _{\scriptstyle1 + 1/2 - 1\atop\scriptstyle\,}}}}exp\left\{ { - \tau {\rm{ }}exp\left( { - \frac{1}{2}{\lambda _{\scriptstyle1\atop\scriptstyle\,}}\tau {{\left( {\mu - {\mu _{\scriptstyle1\atop\scriptstyle\,}}} \right)}^2} + {\beta _{\scriptstyle1\atop\scriptstyle\,}}} \right)} \right\}\,\,\,(1)\)
Values \({\alpha _{\scriptstyle1\atop\scriptstyle\,}}\,,\,{\beta _{\scriptstyle1\atop\scriptstyle\,}},{\rm{ }}{\mu _{\scriptstyle1\atop\scriptstyle\,}},{\rm{ }}{\lambda _{\scriptstyle1\atop\scriptstyle\,}}{\rm{ }}\)shall be computed using the \({\rm{n = 30}}\)observations and the comment above.
From Eq\(\left( 1 \right)\). it is true that for fixed τ and μ as arguments of the function, is proportional to a probability density function of the normal distribution with mean\({\mu _{\scriptstyle1\atop\scriptstyle\,}}\) and variance\(1/\left( {\tau {\lambda _{\scriptstyle1\atop\scriptstyle\,}}{\rm{ }}} \right)\)
From Eq\(\left( 1 \right)\). it is true that for fixed μ and τ as arguments of the function, is proportional to a probability density function of the gamma distribution with parameters \({\alpha _{\scriptstyle1\atop\scriptstyle\,}}\, + {\raise0.7ex\hbox{$1$} \!\mathord{\left/ {\vphantom {1 2}}\right.\\} \!\lower0.7ex\hbox{$2$}}\,\) and \(\,{\lambda _{\scriptstyle1\atop\scriptstyle\,}}\tau {\left( {\mu - {\mu _{\scriptstyle1\atop\scriptstyle\,}}} \right)^2} + {\beta _{\scriptstyle1\atop\scriptstyle\,}}\)
The Gibbs Sampling Algorithm:
The steps of the algorithm are
\(\left( {1.} \right)\,\)Pick starting values \({x_2}^{\left( 0 \right)}\) for \(\,{x_2}\) , and let \(\,\,i = 0\,\)
\(\left( {2.} \right)\,\)let be a simulated value from the conditional distribution of \(\,{x_1}\)given that \(\,\,{X_1} = {x_2}^{\left( i \right)}\)
\(\left( {3.} \right)\,\)Let\(\,{x_2}^{\left( {i + 1} \right)\,\,}\,\) be a simulated value from the conditional distribution \(\,{x_2}\) given that \(\,{X_1} = {x_1}^{\left( {i + 1} \right)}\)
\(\left( {4.} \right)\)Repeat steps \(\,2.\,\,and\,3.\) \(\,i\) where\(\,i + 1\)
The number of simulations I , shall be large enough. The given code below gives all required values and the algorithm itself. The initial hyperparameters of the prior distributions are
\({\lambda _{\scriptstyle0\atop\scriptstyle\,}} = 1,\,\,\,{\mu _{\scriptstyle0\atop\scriptstyle\,}} = 1,\,\,\,{\alpha _{\scriptstyle0\atop\scriptstyle\,}}\, = 0.5,\,\,\,\,\,{\beta _{\scriptstyle0\atop\scriptstyle\,}} = 0.5\)
