Birkhoff ergodic theorem for ergodic Markov processes

This question was previously posted on MSE.

This question might be easy but I am really stuck on it.

Let $$M$$ be compact metric space and $$\mathcal B(M)$$ the Borel $$\sigma$$-algebra of M. Consider the discrete-time Markov process, $$\mathbf{X} =\left(\Omega,\{\mathcal F_n\}_{n\in\mathbb N}, \{X_n\}_{n\in\mathbb N}, \{P_n\}_{n\in\mathbb N} , \{\mathbb P_x\}_{x\in M}\right),$$ with state space $$(M,\mathcal B(M)$$ (I am considering that $$0\in\mathbb N$$) i.e.

1. $$(\Omega,\mathcal F_n),$$ is a filtered measurable space,
2. $$X_n:\Omega\to M$$ is $$\mathcal F_n$$ measurable,
3. $$\mathbb P_x [X_0 = x] =1,$$ for every $$x\in M,$$
4. For every $$0 \leq n \leq m \in \mathbb N,$$ $$f:M\to\mathbb R$$ bounded measurable function, and $$x\in M$$ $$\mathbb E_x [f(X_{n+m}) \mid \mathcal F_n] = (P_m f)(X_n), \ \mathbb P_x\ \mathrm{a.s.},$$ where $$P_n$$ is a transition function on $$(M,\mathcal B(M)),$$ i.e. a family of probability maps $$P_n : M\times \mathcal B(M) \to [0,1],$$ such that
• $$P_0(x,\mathrm{d} y) = \delta_x(\mathrm{d}y),$$
• $$P_n(x,\cdot)$$ is a Borel probability measure for every $$x\in M.$$
• For every $$n,m \in \mathbb N$$ and $$A\in\mathcal B(M),$$ $$P_{n+m}(x, A) = \int_{M} P_n(y,A) P_m(x,\mathrm{d} y).$$

Assume that $$\mathbf{X}$$ admits an ergodic stationary measure $$\mu$$ on $$M,$$ i.e. $$\int_{M} P_n(x,A) \mu(\mathrm d x) = \mu(A),\ \forall \ A\in\mathcal B(M),$$ and if $$P_1(x,A) = 1,\ \forall \ x \ \mu\text{-a.s.}\ \in A \Rightarrow \mu(A) = 0\ \text{or }1.$$

Question: I would like to know if under this setup we would have the following ergodic theorem. For every $$f\in L^1(M,\mathcal B(M), \mu),$$ we obtain $$\lim_{n\to\infty} \frac{1}{n} \sum_{i=0}^{n-1}f (X_n(\omega)) = \int_{M} f(x)\mu(\mathrm{d} x),\ \forall\ \omega \text{-}\mathbb P\ \text{a.s.,}$$ where $$\mathbb P(\mathrm{d} y) := \int_{M}\mathbb P_x (\mathrm{d} y) \mu(\mathrm{d} x).$$

Comments regarding my question

Consider $$\mathbf{X}$$ being an ergodic Markov process (using the above notation). For every $$n\in\mathbb N$$ let us consider the projection map \begin{align*} \pi_n : M^{\mathbb N}&\to M\\ (x_m)_{m\in\mathbb N}&\mapsto x_n. \end{align*}

If we define (via Komolgorov extension Theorem) the Borel probability measure $$P_\mu$$ on $$M^\mathbb N$$ as the unique Borel probability, such that given $$A_0,\ldots,A_n \in M,$$ then $$P_{\mu}\left(\{\omega_n\}_{n\in\mathbb N} \in M^{\mathbb N}; x_i\in A_i, \ \forall \ i\in\{0,1,\ldots,n\}\right) = \int_{A_0}\int_{A_1} \ldots \int_{A_{n-1}} P_1(x_{n},A_n) P_1(x_{n-1},\mathrm{d}x_n) \ldots P_1(x_0,\mathrm{d} x_1) \mu(\mathrm{d}x_0).$$

We have that the shift \begin{align*} \theta: (M^{\mathbb N},\mathcal B(M^{\mathbb N}) , P_\pi)&\to (M^{\mathbb N},\mathcal B(M^{\mathbb N}),P_\pi) \\ (x_{n})_{n\in\mathbb N}&\to (x_{n+1})_{n\in\mathbb N}, \end{align*} is an ergodic dynamical system and we have that for every $$f\in L^1(M,\mathcal B(M),\mu)$$ $$\lim_{n\to\infty}\frac{1}{n}\sum_{i=0}^{n-1} f(\pi_i(\omega)) = \int_M f(x) \mu(\mathrm{d} x),\ \forall \ \omega\text{-}P_\pi\ a.s..$$

How do I translate the information of the canonical process (the one above) to the original Markov process $$\mathbf{X}$$? For every $$\omega \in \Omega$$, we have that $$\left(X_n(\omega)\right)_{n\in\mathbb N}\in M^{\mathbb N},$$
and $$\pi_i\left(\left(X_n(\omega)\right)_{n\in\mathbb N}\right) = X_i(\omega).$$

But it is not clear that $$\lim_{n\to\infty}\frac{1}{n}\sum_{i=0}^{n-1} f( X_i(\omega) ) = \int_M f(x) \mu(\mathrm{d} x),\ \forall \ \omega\text{-}\mathbb P\ a.s.,$$ where $$\mathbb P (\mathrm{d} y) = \int_M \mathbb P_x(\mathrm{d} y) \mu(\mathrm{d}x),$$ can anyone help me? Or simply provide a reference for this ergodic Markov chains result without considering the canonical space.

I believe I found an answer. Note that in a similar way that we constructed $$P_\mu,$$ we may define $$P_x$$ as the unique Borel probability on $$M^{\mathbb N}$$, such that given $$A_0,\ldots,A_n \in M,$$ then $$P_{x}\left(\{\omega_n\}_{n\in\mathbb N} \in M^{\mathbb N}; x_i\in A_i, \ \forall \ i\in\{0,1,\ldots,n\}\right) = \int_{A_0}\int_{A_1} \ldots \int_{A_{n-1}} P_1(x_{n},A_n) P_1(x_{n-1},\mathrm{d}x_n) \ldots P_1(x_0,\mathrm{d} x_1) \delta_x(\mathrm{d}x_0),$$ it is clear that $$P_\mu(\mathrm d y) = \int_M P_x(\mathrm d y) \mu(\mathrm x).$$

Consider the measurable inclusion \begin{align} \iota : \Omega &\to M^\mathbb N\\ \omega&\mapsto (X_n(\omega))_{n\in\mathbb N}. \end{align}

We will prove one auxiliary lemma.

Lemma 1: Let $$A\in \mathcal F = \sigma\left(\bigcup_{n\in\mathbb N}\mathcal F_n\right),$$ then $$P_\mu [\iota (A) ] =0 \Rightarrow \mathbb P [A] = 0.$$ Remember that $$\mathbb P = \int_M \mathbb P_x [A] \mu (\mathrm{d}x).$$

Proof. Let $$A$$ be such that $$P_\mu(\iota(A) ) = 0.$$ This means that given $$\varepsilon >0,$$ there exists $$n_0\in\mathbb N,$$ such that, for every $$n>n_0$$. $$P_\mu[\pi_i \in \pi_i(\iota (A)),\ \forall \ i\in\{0,\ldots, n\}] <\varepsilon.$$

Note that $$X_i(A) = \pi_i(\iota(A)), \ \forall \ i\in\mathbb N.$$

Note that for every $$n_0 < n\in\mathbb N,$$ \begin{align*} \varepsilon \geq P_\mu[\pi_i \in \pi_i(A_i),\ i\in\{0,\ldots, n\}] &= \int_M P_x[\pi_i \in X_i(A),\ i\in\{0,\ldots, n\} ] \mu(\mathrm{d}x)\\ &=\int_{M} \mathbb P_x[ X_i \in X_i(A),\ \forall \ i\in\{0,\ldots,n\}] \mu(\mathrm{d}x)\\ &\geq \int_M \int_M \mathbb P_x[A] \mu(\mathrm d x) = \mathbb P[A]. \end{align*} Implying that $$\mathbb P[A]=0.$$

Let $$f\in L^1(M,\mathcal B(M), \mu),$$ consider the set $$B= \left\{y\in M^\mathbb N; \lim_{n\to\infty}\frac{1}{n}\sum_{i=0}^{n-1} f \circ \pi_i (y) \neq \int_M f(x) \mu(\mathrm{d}x)\right\}.$$

We have that $$P_\mu(B) = 0.$$

Note that for every $$n\geq 0,$$ we have that $$f (X_n(\omega)) = f (\pi_n ( \iota(\omega)).$$

We have that $$\lim_{n\to\infty}\frac{1}{n}\sum_{i=0}^{n-1} f \circ X_i (\omega) \neq \int_M f(x) \mu(\mathrm{d}x) \Leftrightarrow \iota(\omega) \in B\Leftrightarrow \omega \in \iota^{-1}(B).$$

By Lemma 1 we have that $$P_\mu (\iota (\iota^{-1}(B))) \leq P_\mu(B) = 0 \Rightarrow \mathbb P [\iota^{-1}(B)] = 0.$$

Therefore, $$\lim_{n\to\infty}\frac{1}{n}\sum_{i=0}^{n-1} f \circ X_i (\omega) =\int_M f(x) \mu(\mathrm{d}x),\ \forall \ \omega \text{-}\mathbb P \ \text{a.s.}. \$$