all at once, we can obtain it iteratively using Bayes’ rule. Suppose we have a sequence of outcomes $x_{j}\in\{1,0\}$ where $1$ means heads. The distribution, now parameterized only on $\theta$, is defined recursively as $$\begin{aligned}\psi_{j}(\theta\mid x_{1\ldots j}) & \sim\psi_{j-1}(\theta\mid x_{1\ldots j-1})\Phi(\theta,x_{j})\\
 & =\frac{\psi_{j-1}(\theta\mid x_{1\ldots j-1})\Phi(\theta,x_{j})}{\int_{0}^{1}\psi_{j-1}(\theta\mid x_{1\ldots j-1})\Phi(\theta,x_{j})\D\theta},
\end{aligned}$$ where $$\Phi(\theta,x)=\begin{cases}
\theta & x=1\\
1-\theta & x=0,
\end{cases}$$ and $\psi_{0}(\theta)=1$ is a uniform distribution reflecting no prior knowledge of the bias. The $\sim$ notation means dividing by a constant so that $\int_{0}^{1}\psi_{j}(\theta)\,\D\theta=1$ holds. Our best guess of the bias is ${\mathbb{E}\left[\psi_{j}(\theta)\right]}=\int_{0}^{1}\theta\psi_{j}(\theta)\D\theta$.

Now let us apply this analogy to occupancy prediction. Coin toss outcomes are independent, but occupancy transition probabilities depend on the current state. An any given time there are two possible states, so we need to maintain two distributions per time step. Let $\gamma_{k}\in\{0,1\}$ be the occupancy. The transition probabilities of interest are $$\begin{aligned}p_{k} & ={\mathbb{P}\left(\gamma_{k+1}=1\mid\gamma_{k}=1\right)}={\mathbb{E}\left[f_{k}(p_{k})\right]}\\
q_{k} & ={\mathbb{P}\left(\gamma_{k+1}=1\mid\gamma_{k}=0\right)}={\mathbb{E}\left[g_{k}(q_{k})\right]},
\end{aligned}$$ where the density functions $f_{k}(p_{k})$ and $g_{k}(q_{k})$ are the latest iterations of $f_{k,j}(p_{k})$ and $g_{k,j}(q_{k})$, updated each training instance $j$ using $$\begin{aligned}f_{k,j}(p_{k}\mid\gamma_{1\ldots k+1},\gamma_{k}=1) & \sim f_{k,j-1}(p_{k}\mid\gamma_{1\ldots k})\Phi(p_{k},\gamma_{k+1})\\
f_{k,j}(p_{k}\mid\gamma_{1\ldots k+1},\gamma_{k}=0) & =f_{k,j-1}(p_{k}\mid\gamma_{1\ldots k})\\
g_{k,j}(q_{k}\mid\gamma_{1\ldots k+1},\gamma_{k}=0) & \sim g_{k,j-1}(q_{k}\mid\gamma_{1\ldots k})\Phi(q_{k},\gamma_{k+1})\\
g_{k,j}(q_{k}\mid\gamma_{1\ldots k+1},\gamma_{k}=1) & =g_{k,j-1}(q_{k}\mid\gamma_{1\ldots k}).
\end{aligned}$$ The $\sim$ indicates normalization, and $f_{0}(p_{k})=1$ 