Explanation regarding Gibbs Sampling

I am new to MCMC and reading a intro paper regarding Gibbs sampling. However, there are two parts in the paper I cannot understand and get stuck. The first part is equation 2.3 in page 168. It says that if $$X^{'}_{j}$$ is drawn from the process of Gibbs sampling, $$X^{'}_{j}\sim f(x|Y^{'}_{j}=y^{'}_{j})$$. It further claims under "general conditions," the distribution of $$X^{'}_{j}$$ converges to $$f(x)$$ as $$n\rightarrow\infty$$. Can anyone help explain what these general conditions are and why/how this statement is true.

The second part is at the beginning of page 169, expression 2.9 says $$\hat{f(x)}=\frac{1}{m}\sum^{m}_{i=1}f(x|y_{i})$$, and the theory behind this equation is $$E[f(x|Y)]=\int f(x|y)f(y)dy$$. This confuses me because my understanding about monte carlo integration is that for example, $$E[X]=\frac{1}{m}\sum^{m}_{i=1}X_{i}$$ is true due to the law of large number that $$\frac{X_{i}}{m}$$ converges to $$X_{i}\times p(X_{i})$$ when $$m\rightarrow\infty$$. However, I don't see how equation 2.9 can be explained by the concept of monote carlo integration. Can anyone explain the intuition or math behind this equation? Thanks!

EDIT 1 (My Answer/Thought for Part 2)

I have thought about this question, and since no one responds yet, I will try my answer first and see whether someone has a better explanation. For the second part, instead of using $$f(x|y_{i})$$, let's focus on the density of a specific value $$x^{a}$$, which is $$f(x^{a}|y_i)$$. For simplicity, now assume $$y_i$$ is randomly drawn from a density $$f(y)$$ instead of $$f(y|x_i)$$ but the result could be generalized to $$f(y|x_i)$$. Both $$f(x|y_i)$$ and $$f(y)$$ are know distribution.

Now we first randomly draw a $$y_i$$ from $$f(y)$$ and once we have $$y_i$$, we further randomly draw a $$x_i$$ from $$f(x|y_i)$$. Suppose we do this process $$n$$ time, and $$n\rightarrow\infty$$. Given a specific $$y^{a}$$ and $$x^{a}$$, we do know the value of $$f(x^{a}|y^{a})$$, which is just the density of $$x^{a}$$ cnditional on $$y^{a}$$. In other words, by Monte Carlo integration, for a given set of $$y_i$$'s and $$x^{a}$$, we can dervive that $$\frac{1}{n}\sum^{n}_{i=1}f(x^{a}|y_i)\approx\int f(x^{a}|y_i)f(y_i)dy_i=\int f(x^{a}|y)f(y)dy=\int \dfrac{f(x^{a},y)}{f(y)}f(y)dy=\int f(x^{a},y)dy=f(x^{a})$$

We can do the same process to dervie $$f(x^{b})$$, $$f(x^{c})$$, and etc. In the end, we will have $$f(x^{a})$$, $$f(x^{b})$$, $$f(x^{c})$$, ..., which all together constrcut $$f(x)$$. Therefore, $$\hat{f(x)}=\frac{1}{m}\sum^{m}_{i=1}f(x|y_{i})$$.

For the convergence of the Gibbs sampler, the central condition is irreducibility: that is, the possibility to reach any part of the support of the target distribution (with density $g$) starting from any part of this support. If the chain is not irreducible, it will remain stuck in a subset of the support. If the chain is irreducible, then it is recurrent (and the Gibbs sampler convergent) since the stationary distribution is the joint target distribution.

Here are two results from our book (written with one of the authors of this "Gibbs for Kids" paper):

Theorem 10.8 For the Gibbs sampler, if the density $g$ satisfies the positivity condition, it is irreducible.

The "positivity condition" was introduced by Hammersley and Clifford (1970), imposing to the support of $g$ to be equal to the product of the conditional distributions used in the Gibbs sampler (for arbitrary conditioning elements). When it holds the joint density $g$ can be expressed in terms of the conditional densities as $$g(y_1,\ldots,y_p) \propto \prod_{j=1}^p \; {g_{\ell_j}(y_{\ell_j}|y_{\ell_1}, \ldots,y_{\ell_{j-1}},y_{\ell_{j+1}}^\prime,\ldots,y_{\ell_p}^\prime) \over g_{\ell_j}(y_{\ell_j}^\prime|y_{\ell_1},\ldots,y_{\ell_{j-1}}, y_{\ell_{j+1}}^\prime,\ldots,y_{\ell_p}^\prime)}$$ for every permutation $\ell$ on $\{1,2,\ldots,p\}$ and every $y' \in \mathcal{Y}$.

Here is a counter-example where the support is made of two convex sets, for which the Gibbs sampler does not explore the entire support:

A second result due to Tierney (1994) is

Lemma 10.9 If the transition kernel associated with the Gibbs sampler is absolutely continuous with respect to the dominating measure, the resulting chain is Harris recurrent.

For the second part of the question, this is a special case of Rao-Blackwellisation. If the $y_i$'s are distributed from the stationary distribution, $f_Y$, then $$\mathbb{E}^Y[f(x|y_i)] =\int_\mathcal{Y} f(x|y_i)f_Y(y_i)\,\text{d}y_i = f_X(x)$$ This is also true when the $Y_i$'s are converging to this stationary distribution by virtue of the ergodic theorem (on averages).