For $t=1,\dots, n$, let's $r_t\sim\mathcal{N}(0,\,\sigma_t^2)$ and $$\sigma_t^2=\left\{\begin{array}{lcl} \sigma^2 & \text{with probability} & p\\ 1 & \text{with probability} & 1-p \end{array}\right.$$ The sample is $r_t, t=1,\,\dots, n$ and we suppose we know the value of $\sigma$ and our aim is to find the value of $p$ using Gibbs sampling. If the prior of $p$ is $\mathcal{B}eta(1,1)$, the procedure I use to sample $p$ is :
1- Starting from a certain value for $p=p^{0}$.
2- I compute the density of each data $r_t$ in the sample for the two possible values of $\sigma_t^2$ : $$\omega_1 = f(r_t\mid \sigma_t^2=1)\,\, \text{ and } \,\,\omega_{\sigma}=f(r_t\mid \sigma_t^2=\sigma^2)$$ 3- I compute the probability that $\sigma_t^2 = 1$ and $\sigma_t^2 = \sigma^2$ respectively by $$P[\sigma_t^2=1\mid r_t,p] = \frac{p\omega_1}{p\omega_1+(1-p)\omega_{\sigma}}\,\, \text{ and } \,\,P[\sigma_t^2=\sigma^2\mid r_t,p] = 1-P[\sigma_t^2=1\mid r_t,p]$$ 4- Now for each date $t$, I have the probability that $\sigma^2_t=1$ or $\sigma^2_t=\sigma^2$. Then I sample $n_{iter}$ times $\sigma_t^2$ using those probabilities for $t=1,\,\dots,\,n$. At each iteration $j$ ($j=1,\dots,n_{iter}$), the number of time $\sigma_t^2=1$ will be call $n_1^{j}$ and the number of time $\sigma^2_t=\sigma^2$ is $n_{\sigma}^{j}$. Of course $n_{\sigma}^{j}+n_1^{j} = n$ for all $j$.
5- Finally, I sample the posterior of $p$ using $p^{j}\sim \mathcal{B}eta(1+n_{\sigma}^{j},1+n_{1}^{j})$
Problem : This algorithm work fine for $\sigma^2$ far from $1$. But, when $\sigma^2$ is close to $1$, the algorithm doesn't not work.
My understanding : What ever the true value of $p$ is, when $\sigma^2$ is close to $1$, $P[\sigma_t^2=1\mid r_t, p]\approx p$ and does not depend on the data anymore. And then, the posterior is not a posterior as it does not depend on the data. This explanation seems correct because, as $\sigma^2$ is close to $1$ (let's say $\sigma^2=0.9$ or $\sigma^2=1.1$, it is difficult to say if $r_t=0.5$ come from $\mathcal{N}(0,\,1)$ or from $\mathcal{N}(0,\,0.9)$. The algorithm can the suppose all the data come from $\mathcal{N}(0,\,1)$ which lead $p$ to be equal to $0$ or from $\mathcal{N}(0,\,0.9)$ ($p=1$) or even give an equal probability to the two distributions.
My questions :
- Is there another algorithm to better the sampling of the posterior of $p$?
- Theoretically, is it link to the sample size of the number of iteration of my algorithm? In other words, do you think that If I raise the sample size $n$ or the number of iteration $n_{iter}$ I will find the better the sampling of the posterior of $p$?
- What will be the difference between $1$ and $\sigma$ such that the posterior of $p$ can be found properly? Using my empirical trials, $\sigma^2=0.1$, it works, but when $\sigma^2=0.9$ or $\sigma^2=0.8$ it doesn't work.
All kind of contribution will be very helpful. Thanks