Let $X$ be a measurable space, and let $\Phi$ be a stationary, discrete-time Markov process with state space $X$ and stationary measure $p$ on $X$. Let $f:X\to R$ be a bounded harmonic function for $\Phi$, meaning that for $p$-almost all $x\in X$, $f(x)$ is equal to its expectation after one transition, $E_x(f)$.
Call now a measurable subset $S$ of $X$ invariant if for $p$-almost all $x\in X$, we have that the probability of transition from $x$ to $S$ is 1 if $x\in S$ and 0 otherwise. Invariant measurable subsets of $X$ form a $\sigma$-algebra.
Is $f$ measurable, in general, for this $\sigma$-algebra?
A reference would also be appreciated.
Yes.
Let me white the process $\Phi$ as $(Y_n)_{n\geq 0}$, and $\mathcal F_n = \sigma(Y_0, \dots, Y_n)$ the canonical filtration. We denote the indicator function of $A$ by $[A]$, and the sigma-algebra generated by the invariant sets by $\mathcal S$.
Since $f$ is a bounded harmonic function, $(f(Y_n))_{n\geq 0}$ is a bounded $(\mathcal F_n)$-martingale, thus by the $L^q$ martingale convergence theorem converges almost surely and in any $L^q$, $q\in [1,\infty)$, towards some r.v. $F_\infty$.
Claim 1: with probability one, $f(Y_0) = f(Y_1)$.
proof. For every $\epsilon>0$ we can find $n$ such that $\| f(Y_{n+1}) - f(Y_n) \|_1 \leq \epsilon$. By the stationarity of $(Y_n)_{n\geq 0}$, for every $k\geq 0$ and $\epsilon >0$ $$\|f(Y_{k+1}) - f(Y_k)\|_1 = \| f(Y_{n+1}) - f(Y_n) \|_1 \leq \epsilon \implies \mathbb \|f(Y_{k+1}) - f(Y_k) \|_1 = 0, $$ meaning $f(Y_{k+1}) = f(Y_k)$ with probability one. The claim follows. $\square$
With this we can finally answer OP’s question:
Claim 2: every harmonic and bounded $f$ is $\mathcal S$-measurable.
proof. Let $B$ be a measurable set. Then \begin{align} &\mathbb P(f(Y_1) = f(Y_0))=1 \\ &\implies 0 = \mathbb E\left([Y_0 \in f^{-1}(B)] [Y_1 \notin f^{-1}(B)] \right) = \int_{f^{-1}(B)} \mathbb P_y(f(Y_1) \notin B) dp(y) \end{align} meaning that for $p$-a.a. y in $f^{-1}(B)$, $$ \mathbb P_y(Y_1 \in f^{-1}(B)) = 1 . $$ The same applies to the complement of $B$: for $p$-a.a. y in $X \setminus f^{-1}(B)$, $$ \mathbb P_y(Y_1 \notin f^{-1}(B)) = 1 . $$ This means precisely that $f^{-1}(B)$ is an invariant set, thus is in $\mathcal S$. Since our reasoning applies for every measurable $B$, we have just shown that $f$ is $\mathcal S$-measurable. $\square$
Additional note: The proof of Claim 1 actually shows that $\mathbb P( (f(Y_n))_{n\geq 0}$ is a constant sequence$)=1$. Any bounded function $g$ such that $(g(Y_n))_{n\geq 0}$ is a a.s. constant sequence is harmonic. Such a function is $[S]$ when $S$ is invariant.
Additional note 2: the proof also works if one assumes $f(Y_0) \in L^q$ for some $q>1$ instead of $f$ bounded.