# How to use Metropolis-Hastings algorithm for discrete case

I came across this question that is quite confusing for me. So far, I've dealt or seen many examples of MH and Gibbs sampling for continuous case and haven't really thought about them for discrete case.

Q : Estimate the marginal distribution of X with a Gibbs sampler.

$$X|n,y\sim Bim(n,y)\\y\sim Beta(2,4)\,\,\, , \,\,n\sim Poi(16)$$

(y and n are independent)

I derived joint posterior distribution and each full conditional distribution for each parameter.(I'm not sure about them. There might be some mistakes)

(0)$$P(X,n,y)\propto P(X|n,y)P(n)P(y)\propto \Large\binom{n}{x}\large y^x(1-y)^{n-x}*y(1-y)^3*\Large\frac{16^n}{n! }$$

(1)$$X|n,y\sim Bim(n,y)$$

(2)$$y|X,n\propto \large y^{x+2-1}(1-y)^{n-x+4-1}\sim Beta(x+2, n-x+4)$$

(3)$$n|X,y\propto \Large\binom{n}{x}\large(1-y)^{n-x}\Large\frac{16^n}{n!} \sim ?$$

I don't think I can tell which distributions (3) is. Since I can't think of any distributions, I'm gonna use Metropolis-Hastings algorithm. Then I need proposal distribution which has at least the same support of n.(Since 'n' comes from Poisson distribution, it can take 0,1,2,...)

So my question is...

(1)What kind of proposal distribution should I use? (Binomial? Poisson? or some discrete distributions?)

(2)Does discrete case have the same MH-algorithm steps with continuous case? I mean calculating the ratio, comparing with 1 and deciding whether to accept or not.

Any help would be really appreciated:)

You can sample $$y$$ and $$n$$ since you know their distributions and then sample $$X$$ based on those $$y$$ and $$n$$. That would give you enough to visualise the marginal distribution of $$X$$. For example with R

set.seed(2022)
cases <- 10^5
y <- rbeta(cases,2,4)
n <- rpois(cases,16)
X <- rbinom(cases,n,y)

plot(table(X)/cases)


This is plausible: that simulation gives mean(X) of 5.34868 while $$\mathbb E[X] = \frac{2}{2+4} \times 16 \approx 5.333$$.

• Thank you! Actually I'm a bit confused. Isn't that a conditional distribution of $X|n,y$? What I was expecting is to get a distribution of X by using a Gibbs sampler. I mean sequentially updating each parameters and get $(X_{1},y_{1},n_{1}),(X_{2},y_{2},n_{2}),.....,(X_{k},y_{k},n_{k})$ and then just pull out $X_{1},X_{2},...X_{k}$ to get a marginal distribution of X. Nov 23, 2022 at 10:53
• No - that is a simulation of the marginal distribution of $X$ (it is clearly not binomial) which uses the conditional distribution. Nov 23, 2022 at 11:01
• Thank you! But I'm afraid I have to follow my textbook since it requires me to get each full conditional distribution and then use a Gibbs sampler. Do you have any idea for sampling from (3)? Nov 23, 2022 at 11:24

Since the prompt says use a Gibbs sampler, you can just use Gibbs sampling.

AFAIK the conditional (3) looks correct. Also notice that conditional on $$X$$, $$n$$ can be $$X, X+1, X+2, \dots$$. Eventually, the probabilities become negligible for large values of $$X$$, unless $$y$$ is also extreme. Therefore, you can implement say $$125$$ values of $$n$$ (or $$1250$$ if this doesn't give the grassy oscillation in the time series plots). That is, pretend that $$p(n|X,y)$$ has support $$X, X+1, \dots, X+124$$. To find the normalizing constant (denominator), divide each $$p(n|X,y)$$ by the sum of all the $$\sum_{i=X}^{X+124}p(n_i|X,y)$$. This defines a discrete distribution over 125 numbers. Gibbs can be used to sample from this.

• That's interesting! But I guess I still don't know how to continue it after getting that pmf. By summing and dividing over every support of n, we can get a pmf which is still hard to tell what kind of pmf is. We can't sample from it, can we? If we can't sample, how can we use a Gibbs sampler? Nov 24, 2022 at 4:23
• @Terry It is a pmf, meaning it is a discrete distribution of the generic type. You can sample from it by taking $u$ as $Uniform(0, 1)$ and setting cutoff points to pick which value is sampled from. E.g. if the probability vector is $(p_1,\dots,p_{125})$, choose $X$ if $0\le u < p_1$, $X+1$ if $p_1\le u < p_1+p_2$, etc.
– Vons
Nov 24, 2022 at 4:27
• Ohhhh. That's right. I'll try it and see what I can do:) Thank you! Nov 24, 2022 at 4:55
• Feel free to let me know if it worked :] (I'm curious haha)
– Vons
Nov 24, 2022 at 5:25