Sequential Data

Marginals of a chain are chains. Summing (13.2) over \(\mathbf{x}_N\) removes the last factor, because \(\sum_{\mathbf{x}_N}p(\mathbf{x}_N\given\mathbf{x}_{N-1})\)\({}=1\); repeating,

(Sum the innermost variable \(\mathbf{x}_N\) first, then \(\mathbf{x}_{N-1}\), and so on: at each stage the variable being summed appears in exactly one remaining factor, which sums to one.)

The Markov property. By the definition of a conditional distribution and the display above with \(n\) and \(n-1\),

The factor \(p(\mathbf{x}_n\given\mathbf{x}_{n-1})\) appearing here is the model's conditional, and it is also the conditional of the joint: summing \(p(\mathbf{x}_{n-1},\mathbf{x}_n)\)\({}=\sum_{\mathbf{x}_{1:n-2}}p(\mathbf{x}_{1:n})\)\({}=p(\mathbf{x}_{n-1})\,p(\mathbf{x}_n\given\mathbf{x}_{n-1})\) shows that the two agree, so there is no ambiguity in the notation.

Past and future given the present. \(p(\mathbf{x}_{1:N})\) factorises as \(\big[p(\mathbf{x}_1)\prod_{m\le n}p(\mathbf{x}_m\given\mathbf{x}_{m-1})\big]\)\(\big[\prod_{m\gt n}p(\mathbf{x}_m\given\mathbf{x}_{m-1})\big]\)\({}=p(\mathbf{x}_{1:n})\,g(\mathbf{x}_{n:N})\), where \(g\) sums to one over \(\mathbf{x}_{n+1:N}\), so \(g=p(\mathbf{x}_{n+1:N}\given\mathbf{x}_n)\). Dividing by \(p(\mathbf{x}_n)\): \(p(\mathbf{x}_{1:n-1},\mathbf{x}_{n+1:N}\given\mathbf{x}_n)\)\({}=p(\mathbf{x}_{1:n-1}\given\mathbf{x}_n)\,p(\mathbf{x}_{n+1:N}\given\mathbf{x}_n)\).

By d-separation (Section 8.2). In Figure 13.3a every path from \(\mathbf{x}_m\), \(m\lt n-1\), to \(\mathbf{x}_n\) runs along the chain and passes through \(\mathbf{x}_{n-1}\), where the arrows meet head-to-tail. The node \(\mathbf{x}_{n-1}\) is in the conditioning set, so every such path is blocked and \(\mathbf{x}_n\perp\!\!\!\perp\mathbf{x}_{1:n-2}\given\mathbf{x}_{n-1}\), which is (13.3). ∎

Sum and product rules. Exactly as before, summing (13.4) over \(\mathbf{x}_N\), then \(\mathbf{x}_{N-1}\), and so on down to \(\mathbf{x}_{n+1}\), removes one normalised factor at a time (each \(\mathbf{x}_m\) with \(m\gt n\) appears, when it is summed, only in its own factor \(p(\mathbf{x}_m\given\mathbf{x}_{m-1},\mathbf{x}_{m-2})\), because the factors of later variables have already been removed). Hence \(p(\mathbf{x}_{1:n})\)\({}=p(\mathbf{x}_1)p(\mathbf{x}_2\given\mathbf{x}_1)\)\(\prod_{m=3}^{n}p(\mathbf{x}_m\given\mathbf{x}_{m-1},\mathbf{x}_{m-2})\), and the ratio \(p(\mathbf{x}_{1:n})/p(\mathbf{x}_{1:n-1})\) is the single factor \(p(\mathbf{x}_n\given\mathbf{x}_{n-1},\mathbf{x}_{n-2})\). That is (13.122) (which is also (13.123)). For order \(M\) the cancelled factors are \(p(\mathbf{x}_m\given\mathbf{x}_{m-M:m-1})\) and the surviving one is \(p(\mathbf{x}_n\given\mathbf{x}_{n-M:n-1})\).

By d-separation. In Figure 13.3b, take any path from \(\mathbf{x}_m\), \(m\le n-3\), to \(\mathbf{x}_n\), with conditioning set \(\{\mathbf{x}_{n-1},\mathbf{x}_{n-2}\}\). All edges point forward in time. Case 1: the path enters \(\mathbf{x}_n\) from a parent \(\mathbf{v}\in\{\mathbf{x}_{n-1},\mathbf{x}_{n-2}\}\). The edge \(\mathbf{v}\to\mathbf{x}_n\) leaves \(\mathbf{v}\), so the path is head-to-tail or tail-to-tail at \(\mathbf{v}\) (and \(\mathbf{v}\ne\mathbf{x}_m\) since \(m\le n-3\)); \(\mathbf{v}\) is observed, so the path is blocked. Case 2: the path enters \(\mathbf{x}_n\) from a child. Then it visits a node later than \(\mathbf{x}_n\); let \(\mathbf{x}_t\), \(t\gt n\), be the latest node on the path. Both its neighbours on the path are earlier, and edges point forward in time, so both edges point into \(\mathbf{x}_t\): the path is head-to-head there. Every descendant of \(\mathbf{x}_t\) is later than \(n\), hence unobserved, so the path is blocked at \(\mathbf{x}_t\). Either way \(\mathbf{x}_n\perp\!\!\!\perp\mathbf{x}_{1:n-3}\given\mathbf{x}_{n-1},\mathbf{x}_{n-2}\). (This is the directed local Markov property: a node is independent of its non-descendants given its parents.) ∎

An order-\(M\) table assigns a distribution over the \(K\) values of \(\mathbf{x}_n\) to each of the \(K^M\) configurations of \((\mathbf{x}_{n-M},\dots,\mathbf{x}_{n-1})\). A distribution over \(K\) values has \(K-1\) free numbers (they sum to one), so the table has \(K^M(K-1)\). For \(M=1\) this is \(K(K-1)\), the count in the book. In the Gaussian case each \(\mathbf{W}_m\) is \(D\times D\), \(\mathbf{b}\) has \(D\) entries and the symmetric \(\bm\Sigma\) has \(D(D+1)/2\). Tables grow exponentially with \(M\) because they allow any interaction between the \(M\) conditioning values; the autoregressive model allows only additive (linear) effects, which is why it grows linearly. ∎

(a) Take \(m\lt n\) and any conditioning set \(\mathcal{C}\) of observed nodes not containing \(\mathbf{x}_m,\mathbf{x}_n\). The path \(\mathbf{x}_m\leftarrow\mathbf{z}_m\)\({}\to\mathbf{z}_{m+1}\)\({}\to\cdots\)\({}\to\mathbf{z}_n\)\({}\to\mathbf{x}_n\) has intermediate nodes \(\mathbf{z}_m,\dots,\mathbf{z}_n\) only. At \(\mathbf{z}_m\) it is tail-to-tail, at \(\mathbf{z}_{m+1},\dots,\mathbf{z}_n\) head-to-tail, and none of these nodes is in \(\mathcal{C}\) (the latent variables are never observed). A path is blocked only by an observed head-to-tail or tail-to-tail node, or by a head-to-head node with no observed descendant; this path has neither, so it is unblocked and \(\mathbf{x}_m\), \(\mathbf{x}_n\) are not d-separated by \(\mathcal{C}\). A Markov property of order \(M\) would need \(\mathbf{x}_{n}\perp\!\!\!\perp\mathbf{x}_{1:n-M-1}\given\mathbf{x}_{n-M:n-1}\), which the graph does not imply. (Failure of d-separation means the graph allows dependence; part (b) shows it actually occurs.)

(b) With \(\mathbf{A}=\mathbf{I}\) the hidden state is drawn once and kept, so given the state the \(x_m\) are i.i.d. Bernoulli. By Bayes' theorem the posterior odds of state 1 after \(k\) ones in \(n\) trials are

So \(p(\text{state }1\given x_{1:n})=r/(1+r)\), and \(p(x_{n+1}=1\given x_{1:n})\)\({}=0.8\,\frac{r}{1+r}+0.2\,\frac{1}{1+r}\), the stated formula. It depends on the data through \(k=\sum_{m\le n}x_m\), and is strictly increasing in \(k\) (its derivative in \(r\) is \(0.6/(1+r)^2\gt0\)), so flipping any one \(x_m\) changes it. No finite window of recent observations suffices. With \(\mathbf{A}\) close to but not equal to \(\mathbf{I}\) the influence of old observations decays geometrically but never vanishes.

(c) Count the free numbers in \(\bm\pi\), \(\mathbf{A}\) and the emission parameters \(\bm\phi_1,\dots,\bm\phi_K\) of Section 13.2; none of them is indexed by \(n\). ∎

(a) By the sum rule, \(p(z_{n+1,k}=1)\)\({}=\sum_jp(z_{nj}=1)\,A_{jk}\), which is the \(k\)-th entry of \(\mathbf{p}_n\mathbf{A}\). Iterating from \(\mathbf{p}_1=\bm\pi^\T\) gives \(\mathbf{p}_n=\bm\pi^\T\mathbf{A}^{n-1}\).

(b) The set \(\Delta\) of probability vectors is closed and bounded, and \(\mathbf{p}\mapsto\mathbf{p}\mathbf{A}\) maps it into itself (entries stay nonnegative and \(\mathbf{p}\mathbf{A}\mathbf{1}\)\({}=\mathbf{p}\mathbf{1}\)\({}=1\) because the rows of \(\mathbf{A}\) sum to one). Start from any \(\mathbf{p}\in\Delta\) and form the averages \(\bar{\mathbf{p}}_T\)\({}=\frac1T\sum_{t=0}^{T-1}\mathbf{p}\mathbf{A}^t\in\Delta\). Then \(\bar{\mathbf{p}}_T\mathbf{A}-\bar{\mathbf{p}}_T\)\({}=\frac1T(\mathbf{p}\mathbf{A}^T-\mathbf{p})\), whose entries are at most \(2/T\) in absolute value. Since \(\Delta\) is compact, a subsequence \(\bar{\mathbf{p}}_{T_i}\) converges to some \(\bm\pi^{\star\T}\in\Delta\), and passing to the limit gives \(\bm\pi^{\star\T}\mathbf{A}-\bm\pi^{\star\T}\)\({}=\mathbf{0}\).

(c) By (a), \(\mathbf{p}_n\)\({}=\bm\pi^{\star\T}\mathbf{A}^{n-1}\)\({}=\bm\pi^{\star\T}\) for all \(n\). The joint distribution of \((\mathbf{z}_{n:n+s},\mathbf{x}_{n:n+s})\) is obtained from (13.10) by summing out the other variables (the later ones drop out because their factors sum to one, the earlier ones produce the marginal \(\mathbf{p}_n\) of \(\mathbf{z}_n\)):

With homogeneous factors and \(\mathbf{p}_n=\bm\pi^{\star\T}\), the right-hand side does not depend on \(n\).

(d) This is the convergence theorem for finite Markov chains, proved in Result 11.17 (with \(T=\mathbf{A}\)). ∎

Step 1: the log of the joint is linear in the indicators. Take logarithms in (13.10) and use (13.7)–(13.9):

(When some \(A_{jk}=0\) the term \(z_{n-1,j}z_{nk}\ln A_{jk}\) is read as 0 if \(z_{n-1,j}z_{nk}=0\); configurations with \(z_{n-1,j}z_{nk}=1\) then have probability zero and do not contribute to (13.12).)

Step 2: take the posterior expectation. \(Q\) is the expectation of this expression under \(p(\mathbf{Z}\given\mathbf{X},\bm\theta^{\rm old})\). Expectation is linear, and \(\ln\pi_k\), \(\ln A_{jk}\), \(\ln p(\mathbf{x}_n\given\bm\phi_k)\) are constants (they depend on \(\bm\theta\) and the data, not on \(\mathbf{Z}\)), so \(Q\) is the same expression with each indicator replaced by its posterior mean: \(\E[z_{1k}]\), \(\E[z_{n-1,j}z_{nk}]\) and \(\E[z_{nk}]\).

Step 3: those means are \(\gamma\) and \(\xi\). The function \(z_{nk}\) depends on \(\mathbf{Z}\) only through \(\mathbf{z}_n\), so summing \(p(\mathbf{Z}\given\mathbf{X})z_{nk}\) over all other latent variables first gives \(\E[z_{nk}]\)\({}=\sum_{\mathbf{z}_n}p(\mathbf{z}_n\given\mathbf{X})z_{nk}\)\({}=p(z_{nk}=1\given\mathbf{X})\)\({}=\gamma(z_{nk})\), which is (13.15). Likewise \(z_{n-1,j}z_{nk}\) depends only on \((\mathbf{z}_{n-1},\mathbf{z}_n)\) and equals 1 exactly when \(z_{n-1,j}=z_{nk}=1\), giving (13.16). Substituting gives (13.17). ∎

Initial probabilities. Maximise \(F(\bm\pi)=\sum_kg_k\ln\pi_k\) with \(g_k=\gamma(z_{1k})\ge0\), subject to \(\sum_k\pi_k=1\) (Lagrange multipliers). The Lagrangian \(\sum_kg_k\ln\pi_k+\lambda(\sum_k\pi_k-1)\) has derivative \(g_k/\pi_k+\lambda\) in \(\pi_k\); setting it to zero gives \(\pi_k=-g_k/\lambda\). Summing over \(k\) and using the constraint, \(1=-\sum_kg_k/\lambda\), so \(\lambda=-\sum_jg_j\) and \(\pi_k=g_k/\sum_jg_j\), which is (13.18).

Why this is the maximum. \(F\) is a sum of concave functions \(g_k\ln\pi_k\) and the constraint set (the probability simplex) is convex, so a stationary point of the Lagrangian in the interior is the global maximum. A direct check, which also covers \(g_k=0\): with \(\hat\pi_k=g_k/G\), \(G=\sum_jg_j\), Gibbs' inequality (\(\KL(\hat{\bm\pi}\|\bm\pi)\ge0\), toolkit) gives \(\sum_k\hat\pi_k\ln\pi_k\)\({}\le\sum_k\hat\pi_k\ln\hat\pi_k\), and multiplying by \(G\) gives \(F(\bm\pi)\le F(\hat{\bm\pi})\) for every \(\bm\pi\) in the simplex.

Transition probabilities. The second term of (13.17) is \(\sum_j\big[\sum_kc_{jk}\ln A_{jk}\big]\) with \(c_{jk}=\sum_{n=2}^{N}\xi(z_{n-1,j},z_{nk})\). The constraints \(\sum_kA_{jk}=1\) act on each row separately, so each row \(j\) is an instance of the problem just solved, with \(g_k\) replaced by \(c_{jk}\) and one multiplier \(\lambda_j\) per row: \(A_{jk}=c_{jk}/\sum_lc_{jl}\), which is (13.19).

The simplified denominators. \(\gamma(\mathbf{z}_1)\) is a distribution, so \(\sum_j\gamma(z_{1j})=1\). Summing the pairwise marginal over its second argument gives the single marginal: \(\sum_l\xi(z_{n-1,j},z_{nl})\)\({}=p(z_{n-1,j}=1\given\mathbf{X})\)\({}=\gamma(z_{n-1,j})\). ∎

Only the last sum of (13.17) involves \(\bm\phi\), and with separate parameters per state it splits into one weighted Gaussian log likelihood per \(k\). With \(N_k=\sum_n\gamma(z_{nk})\) and dropping \(-\frac D2\ln2\pi\) terms,

Mean. The gradient in \(\bm\mu_k\) (matrix calculus: \(\partial(\mathbf{a}^\T\mathbf{S}\mathbf{a})/\partial\mathbf{a}\)\({}=2\mathbf{S}\mathbf{a}\) for symmetric \(\mathbf{S}\)) is \(\sum_n\gamma(z_{nk})\bm\Sigma_k^{-1}(\mathbf{x}_n-\bm\mu_k)\). Multiplying by \(\bm\Sigma_k\) and setting it to zero gives \(\sum_n\gamma(z_{nk})\mathbf{x}_n=N_k\bm\mu_k\), which is (13.20). \(Q_k\) is a concave quadratic in \(\bm\mu_k\) (its Hessian is \(-N_k\bm\Sigma_k^{-1}\)), and the maximiser does not depend on \(\bm\Sigma_k\), so it is the \(\bm\mu_k\)-part of the joint maximiser.

Covariance. Put \(\bm\mu_k\) at (13.20) and write the quadratic form as a trace (toolkit): with \(\mathbf{S}_k\)\({}=\sum_n\gamma(z_{nk})(\mathbf{x}_n-\bm\mu_k)(\mathbf{x}_n-\bm\mu_k)^\T\) and \(\bm\Lambda=\bm\Sigma_k^{-1}\),

Using \(\partial\ln|\bm\Lambda|/\partial\bm\Lambda\)\({}=\bm\Lambda^{-1}\) and \(\partial\Tr(\bm\Lambda\mathbf{S})/\partial\bm\Lambda\)\({}=\mathbf{S}\) (both for symmetric matrices), the stationarity condition is \(\frac{N_k}{2}\bm\Lambda^{-1}-\frac12\mathbf{S}_k\)\({}=\mathbf{0}\), i.e. \(\bm\Sigma_k=\mathbf{S}_k/N_k\), which is (13.21). Since \(\ln|\bm\Lambda|\) is concave on positive definite matrices and the trace term is linear, \(Q_k\) is concave in \(\bm\Lambda\) and the stationary point is the global maximum (when \(\mathbf{S}_k\) is positive definite; otherwise the likelihood is unbounded, the singularity of Section 9.2.1). ∎

(a) If \(z_k=1\) and \(x_i=1\), all exponents in (13.22) vanish except \(x_iz_k\), so the product is \(\mu_{ik}\), as required. The emission term of (13.17) becomes \(\sum_n\sum_k\gamma(z_{nk})\sum_ix_{ni}\ln\mu_{ik}\)\({}=\sum_k\sum_ic_{ik}\ln\mu_{ik}\), with \(c_{ik}=\sum_n\gamma(z_{nk})x_{ni}\). For each \(k\) this is the problem of Result 13.7 with constraint \(\sum_i\mu_{ik}=1\), so \(\mu_{ik}=c_{ik}/\sum_{i'}c_{i'k}\). Finally \(\sum_{i'}c_{i'k}\)\({}=\sum_n\gamma(z_{nk})\sum_{i'}x_{ni'}\)\({}=\sum_n\gamma(z_{nk})\) because each \(\mathbf{x}_n\) has exactly one entry equal to 1: this is (13.23).

(b) The emission term is \(\sum_k\sum_i\sum_n\gamma(z_{nk})\{x_{ni}\ln\mu_{ik}\)\({}+(1-x_{ni})\ln(1-\mu_{ik})\}\), and each \(\mu_{ik}\) appears in one term without constraint. Its derivative is \(\sum_n\gamma(z_{nk})\{x_{ni}/\mu_{ik}\)\({}-(1-x_{ni})/(1-\mu_{ik})\}\); multiplying by \(\mu_{ik}(1-\mu_{ik})\) and setting to zero gives \(\sum_n\gamma(z_{nk})x_{ni}\)\({}=\mu_{ik}\sum_n\gamma(z_{nk})\), which is (13.23). The function is concave in \(\mu_{ik}\in(0,1)\) (a positive combination of \(\ln\mu\) and \(\ln(1-\mu)\)), so this is the maximum. ∎

By (13.43) below, \(\xi(z_{n-1,j},z_{nk})\)\({}=\alpha(z_{n-1,j})\,A_{jk}\,p(\mathbf{x}_n\given\bm\phi_k)\,\beta(z_{nk})/p(\mathbf{X})\), which is zero for every \(n\) when \(A_{jk}^{\rm old}=0\). So the numerator of (13.19) is zero and \(A_{jk}^{\rm new}=0\). Similarly, by (13.37) and (13.33), \(\gamma(z_{1k})\)\({}=\pi_k\,p(\mathbf{x}_1\given\bm\phi_k)\,\beta(z_{1k})/p(\mathbf{X})\)\({}=0\) when \(\pi_k^{\rm old}=0\), so \(\pi_k^{\rm new}=0\) by (13.18). Induction on the iteration number finishes the proof. (One caveat: if a state \(j\) gets \(\sum_{n\lt N}\gamma(z_{nj})=0\), row \(j\) of (13.19) is \(0/0\); any convention, such as keeping the old row, keeps its zeros.) The converse is not true: a positive entry can be driven towards zero by EM, but a zero can never become positive. ∎

Group the factors of (13.6):

\(h\) sums to one over \(\mathcal{F}_n\) for every \(\mathbf{z}_n\). Sum over \(\mathbf{x}_N\) (only \(p(\mathbf{x}_N\given\mathbf{z}_N)\) contains it, and it sums to one), then over \(\mathbf{z}_N\) (only \(p(\mathbf{z}_N\given\mathbf{z}_{N-1})\) remains with it), then \(\mathbf{x}_{N-1}\), and so on down to \(\mathbf{z}_{n+1}\). The same peeling shows that summing the last pair \((\mathbf{x}_{N},\mathbf{z}_{N})\) out of the length-\(N\) model gives the length-\((N-1)\) model, which is the consistency statement.

\(g\) is the joint of \(\mathcal{P}_n\) and \(\mathbf{z}_n\). Summing the display over \(\mathcal{F}_n\) and \(\mathbf{x}_n\) removes \(h\) and \(p(\mathbf{x}_n\given\mathbf{z}_n)\), so \(p(\mathcal{P}_n,\mathbf{z}_n)\)\({}=g(\mathcal{P}_n,\mathbf{z}_n)\).

Conclusion. Dividing the display by \(p(\mathbf{z}_n)\): \(p(\mathcal{P}_n,\mathbf{x}_n,\mathcal{F}_n\given\mathbf{z}_n)\)\({}=\frac{g(\mathcal{P}_n,\mathbf{z}_n)}{p(\mathbf{z}_n)}\,p(\mathbf{x}_n\given\mathbf{z}_n)\,h(\mathcal{F}_n;\mathbf{z}_n)\). Each factor is a normalised conditional distribution of its own group given \(\mathbf{z}_n\) (the first is \(p(\mathcal{P}_n\given\mathbf{z}_n)\); summing the product over two groups leaves the third, which identifies \(h=p(\mathcal{F}_n\given\mathbf{z}_n)\)), and a joint conditional that is a product of conditionals of the groups is exactly mutual conditional independence. ∎

Two facts are used repeatedly. (i) If \(\mathbf{a}\perp\!\!\!\perp\mathbf{b}\given\mathbf{c}\) then the same holds for any sub-collections of \(\mathbf{a}\) and \(\mathbf{b}\) (sum out the rest). (ii) If \((\mathbf{a},\mathbf{d})\perp\!\!\!\perp\mathbf{b}\given\mathbf{c}\), then \(\mathbf{a}\perp\!\!\!\perp\mathbf{b}\given(\mathbf{c},\mathbf{d})\): indeed \(p(\mathbf{a},\mathbf{b}\given\mathbf{c},\mathbf{d})\)\({}=p(\mathbf{a},\mathbf{d}\given\mathbf{c})p(\mathbf{b}\given\mathbf{c})/p(\mathbf{d}\given\mathbf{c})\)\({}=p(\mathbf{a}\given\mathbf{c},\mathbf{d})p(\mathbf{b}\given\mathbf{c})\), and \(p(\mathbf{b}\given\mathbf{c})\)\({}=p(\mathbf{b}\given\mathbf{c},\mathbf{d})\) by summing over \(\mathbf{a}\). For d-separation recall: a path is blocked if it has an observed node at which the arrows meet head-to-tail or tail-to-tail, or a head-to-head node that is unobserved and has no observed descendant; in Figure 13.5 every \(\mathbf{x}\) node is a leaf, and every path between two nodes runs along the \(\mathbf{z}\)-chain.

(13.24) Result 13.11 at \(n\): \(\mathbf{x}_{1:n}\) (part of \(\mathcal{P}_n\) together with \(\mathbf{x}_n\)) and \(\mathbf{x}_{n+1:N}\) (part of \(\mathcal{F}_n\)) are independent given \(\mathbf{z}_n\), by (i). d-separation: a path from \(\mathbf{x}_m\), \(m\le n\), to \(\mathbf{x}_{m'}\), \(m'\gt n\), runs \(\mathbf{x}_m\leftarrow\mathbf{z}_m\cdots\mathbf{z}_{m'}\)\({}\to\mathbf{x}_{m'}\) and passes through \(\mathbf{z}_n\) head-to-tail (or tail-to-tail if \(m=n\)); \(\mathbf{z}_n\) is observed, so it is blocked.

(13.25) Result 13.11 at \(n\): \(\mathbf{x}_n\perp\!\!\!\perp\mathcal{P}_n\given\mathbf{z}_n\), so \(\mathbf{x}_{1:n-1}\perp\!\!\!\perp\mathbf{x}_n\given\mathbf{z}_n\), which is (13.25). d-separation: every path from \(\mathbf{x}_{1:n-1}\) to \(\mathbf{x}_n\) ends \(\mathbf{z}_{n-1}\to\mathbf{z}_n\to\mathbf{x}_n\), head-to-tail at the observed \(\mathbf{z}_n\).

(13.26) Result 13.11 at \(n-1\): \(\mathcal{P}_{n-1}\cup\{\mathbf{x}_{n-1}\}\) is independent of \(\mathcal{F}_{n-1}\ni\mathbf{z}_n\) given \(\mathbf{z}_{n-1}\); by (i), \(\mathbf{x}_{1:n-1}\perp\!\!\!\perp\mathbf{z}_n\given\mathbf{z}_{n-1}\), which is (13.26). d-separation: paths from \(\mathbf{x}_m\), \(m\le n-1\), to \(\mathbf{z}_n\) pass \(\mathbf{z}_{n-1}\) head-to-tail (tail-to-tail when \(m=n-1\)), and \(\mathbf{z}_{n-1}\) is observed.

(13.27) Result 13.11 at \(n+1\): \(\mathcal{P}_{n+1}\ni\mathbf{z}_n\) is independent of \((\mathbf{x}_{n+1},\mathcal{F}_{n+1})\supset\mathbf{x}_{n+1:N}\) given \(\mathbf{z}_{n+1}\). d-separation: paths from \(\mathbf{z}_n\) to \(\mathbf{x}_{m}\), \(m\gt n\), pass \(\mathbf{z}_{n+1}\) head-to-tail, observed.

(13.28) Result 13.11 at \(n+1\): \(\mathbf{x}_{n+1}\perp\!\!\!\perp\mathcal{F}_{n+1}\given\mathbf{z}_{n+1}\), and \(\mathbf{x}_{n+2:N}\subset\mathcal{F}_{n+1}\). d-separation: paths from \(\mathbf{x}_{n+1}\) to \(\mathbf{x}_m\), \(m\ge n+2\), start \(\mathbf{x}_{n+1}\leftarrow\mathbf{z}_{n+1}\)\({}\to\mathbf{z}_{n+2}\), tail-to-tail at the observed \(\mathbf{z}_{n+1}\).

(13.29) Result 13.11 at \(n\) with (i): \((\mathbf{x}_{1:n-1},\mathbf{z}_{n-1})\), \(\mathbf{x}_n\) and \(\mathbf{x}_{n+1:N}\) are mutually independent given \(\mathbf{z}_n\). Hence

Divide by \(p(\mathbf{z}_{n-1}\given\mathbf{z}_n)\): \(p(\mathbf{X}\given\mathbf{z}_{n-1},\mathbf{z}_n)\)\({}=p(\mathbf{x}_{1:n-1}\given\mathbf{z}_{n-1},\mathbf{z}_n)\)\(\,p(\mathbf{x}_n\given\mathbf{z}_n)\)\(\,p(\mathbf{x}_{n+1:N}\given\mathbf{z}_n)\), and the first factor is \(p(\mathbf{x}_{1:n-1}\given\mathbf{z}_{n-1})\) by (13.26). d-separation: with \(\mathbf{z}_{n-1},\mathbf{z}_n\) observed, paths between the three groups pass \(\mathbf{z}_{n-1}\) or \(\mathbf{z}_n\) head-to-tail or tail-to-tail.

(13.30) In the chain of length \(N+1\) (which marginalises to the given one, by the consistency part of Result 13.11), Result 13.11 at \(n=N+1\) gives \(\mathbf{x}_{N+1}\perp\!\!\!\perp\mathcal{P}_{N+1}\given\mathbf{z}_{N+1}\), and \(\mathbf{X}\subset\mathcal{P}_{N+1}\). d-separation: paths from \(\mathbf{X}\) to \(\mathbf{x}_{N+1}\) end \(\mathbf{z}_N\)\({}\to\mathbf{z}_{N+1}\)\({}\to\mathbf{x}_{N+1}\), head-to-tail at the observed \(\mathbf{z}_{N+1}\).

(13.31) Result 13.11 at \(n=N\) in the chain of length \(N+1\): \((\mathcal{P}_N,\mathbf{x}_N)\) is independent of \(\mathcal{F}_N\)\({}=(\mathbf{z}_{N+1},\mathbf{x}_{N+1})\) given \(\mathbf{z}_N\); so \(\mathbf{z}_{N+1}\perp\!\!\!\perp\mathbf{X}\given\mathbf{z}_N\). d-separation: paths from \(\mathbf{x}_m\), \(m\le N\), to \(\mathbf{z}_{N+1}\) pass \(\mathbf{z}_N\) head-to-tail (tail-to-tail if \(m=N\)), observed.

(13.5) Result 13.11 at \(n\): \(\mathbf{z}_{n-1}\in\mathcal{P}_n\) and \(\mathbf{z}_{n+1}\in\mathcal{F}_n\). ∎

and the last sum contains \(\alpha(\mathbf{z}_{n-1})\) by (13.34), which is (13.36). In the third line, the product rule gives \(p(\mathbf{x}_{1:n}\given\mathbf{z}_n)\)\({}=p(\mathbf{x}_n\given\mathbf{z}_n)\,p(\mathbf{x}_{1:n-1}\given\mathbf{x}_n,\mathbf{z}_n)\), and (13.25) removes \(\mathbf{x}_n\) from the conditioning set. In the seventh, \(p(\mathbf{x}_{1:n-1},\mathbf{z}_n\given\mathbf{z}_{n-1})\)\({}=p(\mathbf{x}_{1:n-1}\given\mathbf{z}_{n-1},\mathbf{z}_n)\,p(\mathbf{z}_n\given\mathbf{z}_{n-1})\) and (13.26) removes \(\mathbf{z}_n\). For (13.37), (13.34) at \(n=1\) is \(p(\mathbf{x}_1,\mathbf{z}_1)\), and the 1-of-\(K\) product picks out \(\pi_kp(\mathbf{x}_1\given\bm\phi_k)\) when \(z_{1k}=1\), as in (13.8)–(13.9). ∎

and the first factor is \(\beta(\mathbf{z}_{n+1})\) by (13.35). In the last line \(p(\mathbf{x}_{n+1:N}\given\mathbf{z}_{n+1})\)\({}=p(\mathbf{x}_{n+1}\given\mathbf{z}_{n+1})\,p(\mathbf{x}_{n+2:N}\given\mathbf{z}_{n+1},\mathbf{x}_{n+1})\) and (13.28) drops \(\mathbf{x}_{n+1}\).

The starting value (13.39). At \(n=N\), (13.33) with the definition of \(\alpha\) reads

and the left-hand side equals \(p(\mathbf{X},\mathbf{z}_N)/p(\mathbf{X})\) by Bayes' theorem. So (13.33) holds at \(n=N\) exactly when \(\beta(\mathbf{z}_N)=1\). This is also what (13.35) means when the set of future observations is empty (the probability of observing nothing is one), and with it (13.38) at \(n=N-1\) gives \(\beta(\mathbf{z}_{N-1})\)\({}=\sum_{\mathbf{z}_N}p(\mathbf{x}_N\given\mathbf{z}_N)p(\mathbf{z}_N\given\mathbf{z}_{N-1})\)\({}=p(\mathbf{x}_N\given\mathbf{z}_{N-1})\), as it should. ∎

(b) Sum (13.33) over \(\mathbf{z}_n\): the left side sums to one because \(\gamma\) is a distribution, so \(p(\mathbf{X})\)\({}=\sum_{\mathbf{z}_n}\alpha(\mathbf{z}_n)\beta(\mathbf{z}_n)\). (Directly: \(\alpha(\mathbf{z}_n)\beta(\mathbf{z}_n)\)\({}=p(\mathbf{x}_{1:n},\mathbf{z}_n)\,p(\mathbf{x}_{n+1:N}\given\mathbf{z}_n)\)\({}=p(\mathbf{X},\mathbf{z}_n)\) by (13.24), and summing over \(\mathbf{z}_n\) gives \(p(\mathbf{X})\).) At \(n=N\), \(\beta(\mathbf{z}_N)=1\) gives (13.42), so the likelihood needs only the forward pass.

(c) By Bayes' theorem and then (13.29),

and \(p(\mathbf{x}_{1:n-1}\given\mathbf{z}_{n-1})p(\mathbf{z}_{n-1})\)\({}=\alpha(\mathbf{z}_{n-1})\), \(p(\mathbf{x}_{n+1:N}\given\mathbf{z}_n)\)\({}=\beta(\mathbf{z}_n)\). Check: summing (13.43) over \(\mathbf{z}_{n-1}\) and using (13.36) gives \(\alpha(\mathbf{z}_n)\beta(\mathbf{z}_n)/p(\mathbf{X})\)\({}=\gamma(\mathbf{z}_n)\), as it must.

(d) Substitute (13.33) into (13.20): the common factor \(1/p(\mathbf{X})\) cancels between numerator and denominator. The same happens in (13.18), (13.19), (13.21) and (13.23). Typo. The book prints the upper limits of the sums in (13.40) as \(n\); they are \(N\). ∎

One step of (13.36) computes, for each of the \(K\) values of \(k\), a sum of \(K\) products \(\alpha(z_{n-1,j})A_{jk}\) and one more multiplication: \(K^2+K\) multiplications and \(K(K-1)\) additions, so \(O(K^2)\) per step and \(O(K^2N)\) for the chain. (13.38) is the same. \(\gamma\) costs \(K\) multiplications per step, and \(\xi\) has \(K^2\) entries per step, each a product of four stored numbers. Direct evaluation of (13.11) sums \(K^N\) terms, each a product of about \(2N\) factors. For example, with \(K=3\) and \(N=100\), \(K^N\approx5\times10^{47}\), against about \(2\cdot9\cdot100\) multiply-adds per pass. The saving comes from the distributive law: in \(\sum_{\mathbf{z}_1}\cdots\sum_{\mathbf{z}_N}\prod(\cdots)\) each \(\mathbf{z}_n\) appears in only two neighbouring transition factors, so the sums can be pushed inside the product one at a time, and \(\alpha(\mathbf{z}_n)\) is exactly the partial sum over all paths that end in a given state at time \(n\). ∎

and \(p(\mathbf{z}_N,\mathbf{X})=\alpha(\mathbf{z}_N)\) by (13.34). The new step is part of the chain of length \(N+1\); by the consistency in Result 13.11 its marginal over the first \(N\) steps is the original model, so \(p(\mathbf{X})\), \(\alpha(\mathbf{z}_N)\) and \(p(\mathbf{z}_N\given\mathbf{X})\) are unchanged by adding it. ∎

The message rules. From Chapter 8, a factor sends \(\mu_{f\to x}(x)\)\({}=\sum f(x,x_1,\dots,x_M)\prod_m\mu_{x_m\to f}(x_m)\) (sum over its other variables, (8.66)), and a variable sends the product of the messages from its other neighbouring factors, \(\mu_{x\to f}(x)\)\({}=\prod_{l\in{\rm ne}(x)\setminus f}\mu_{f_l\to x}(x)\) (8.69). A leaf factor sends itself (8.71) and a leaf variable sends 1 (8.70).

Forward messages. \(\mathbf{z}_{n-1}\) has two neighbours, \(f_{n-1}\) and \(f_n\), so its message to \(f_n\) is the single message from \(f_{n-1}\): (13.47) (variable nodes "perform no computation"). Factor \(f_n\) has neighbours \(\mathbf{z}_{n-1}\), \(\mathbf{z}_n\), so (8.66) gives (13.48), and substituting (13.47) gives (13.49). With (13.46) and (13.50) this is \(\alpha(\mathbf{z}_n)\)\({}=\sum_{\mathbf{z}_{n-1}}p(\mathbf{z}_n\given\mathbf{z}_{n-1})p(\mathbf{x}_n\given\mathbf{z}_n)\alpha(\mathbf{z}_{n-1})\), the recursion (13.36). It starts at the leaf \(h\): \(\mu_{h\to\mathbf{z}_1}\)\({}=h(\mathbf{z}_1)\)\({}=p(\mathbf{z}_1)p(\mathbf{x}_1\given\mathbf{z}_1)\), which is (13.37). Same start and same recursion, so the messages equal \(\alpha\) at every \(n\).

The message is \(\alpha\) by definition, too (Exercise 13.13). By (8.64), the message from \(f_n\) to \(\mathbf{z}_n\) is the sum, over all variables in the part of the tree on the \(f_n\) side of \(\mathbf{z}_n\), of the product of all factors in that part:

where the middle equality is (13.6) for the first \(n\) steps (the marginal of the chain, Result 13.11) and the last is the sum rule. That is (13.34).

Backward messages (Exercise 13.14). Messages from the root \(\mathbf{z}_N\) back to the leaf: \(\mathbf{z}_N\) has the single neighbour \(f_N\), so it sends \(\mu_{\mathbf{z}_N\to f_N}=1\) (8.70), which is \(\beta(\mathbf{z}_N)=1\). Each \(\mathbf{z}_{n+1}\) passes on what it receives from \(f_{n+2}\), and \(f_{n+1}\) applies (8.66): (13.51). With (13.46), \(\beta(\mathbf{z}_n)\)\({}=\sum_{\mathbf{z}_{n+1}}p(\mathbf{z}_{n+1}\given\mathbf{z}_n)\)\(p(\mathbf{x}_{n+1}\given\mathbf{z}_{n+1})\beta(\mathbf{z}_{n+1})\), which is (13.38). By definition (8.67) instead,

by Result 13.11 (the product is \(p(\mathcal{F}_n\given\mathbf{z}_n)\)), which is (13.35).

Marginals. By (8.63) the marginal of a variable node is the product of all incoming messages. The function represented by the factor graph is \(p(\mathbf{X},\mathbf{Z})\), so this "marginal" is \(\sum_{\mathbf{Z}\setminus\mathbf{z}_n}p(\mathbf{X},\mathbf{Z})\)\({}=p(\mathbf{z}_n,\mathbf{X})\): (13.53). Dividing by \(p(\mathbf{X})\) gives (13.54), in agreement with (13.33).

Pairwise marginal (Exercise 13.11). By (8.72), the marginal of the variables of factor \(f_n\) is the factor times the messages arriving at it: \(p(\mathbf{z}_{n-1},\mathbf{z}_n,\mathbf{X})\)\({}=f_n(\mathbf{z}_{n-1},\mathbf{z}_n)\,\mu_{\mathbf{z}_{n-1}\to f_n}(\mathbf{z}_{n-1})\,\mu_{\mathbf{z}_n\to f_n}(\mathbf{z}_n)\). Here \(\mu_{\mathbf{z}_{n-1}\to f_n}\)\({}=\mu_{f_{n-1}\to\mathbf{z}_{n-1}}\)\({}=\alpha(\mathbf{z}_{n-1})\) by (13.47), and \(\mu_{\mathbf{z}_n\to f_n}\)\({}=\mu_{f_{n+1}\to\mathbf{z}_n}\)\({}=\beta(\mathbf{z}_n)\) (the variable passes on the other message). With (13.46) and dividing by \(p(\mathbf{X})\) this is (13.43). ∎

(13.57) is the product rule applied repeatedly: \(p(\mathbf{x}_{1:n})\)\({}=p(\mathbf{x}_1)p(\mathbf{x}_2\given\mathbf{x}_1)\cdots p(\mathbf{x}_n\given\mathbf{x}_{1:n-1})\). (13.58) is (13.55) multiplied out. Substitute (13.58) into (13.36):

and divide by \(\prod_{m=1}^{n-1}c_m\gt0\): (13.59). Summing (13.59) over \(\mathbf{z}_n\), and using \(\sum_{\mathbf{z}_n}\widehat\alpha(\mathbf{z}_n)\)\({}=1\), gives \(c_n\)\({}=\sum_{\mathbf{z}_n}p(\mathbf{x}_n\given\mathbf{z}_n)\sum_{\mathbf{z}_{n-1}}\widehat\alpha(\mathbf{z}_{n-1})p(\mathbf{z}_n\given\mathbf{z}_{n-1})\): the normaliser. Interpretation. The inner sum is the one-step prediction \(p(\mathbf{z}_n\given\mathbf{x}_{1:n-1})\) (sum rule, and \(\mathbf{z}_n\perp\!\!\!\perp\mathbf{x}_{1:n-1}\given\mathbf{z}_{n-1}\) from (13.26)); multiplying by \(p(\mathbf{x}_n\given\mathbf{z}_n)\)\({}=p(\mathbf{x}_n\given\mathbf{z}_n,\mathbf{x}_{1:n-1})\) (from (13.25)) gives \(p(\mathbf{x}_n,\mathbf{z}_n\given\mathbf{x}_{1:n-1})\), whose sum over \(\mathbf{z}_n\) is \(p(\mathbf{x}_n\given\mathbf{x}_{1:n-1})=c_n\). So (13.59) is Bayes' theorem: predict, multiply by the likelihood of the new observation, normalise. ∎

(13.61). By (13.57), \(\prod_{m=n+1}^{N}c_m\)\({}=p(\mathbf{x}_{1:N})/p(\mathbf{x}_{1:n})\)\({}=p(\mathbf{x}_{n+1:N}\given\mathbf{x}_{1:n})\), so (13.60) and (13.35) give (13.61). Unlike \(\beta\), which is a probability of \(N-n\) observations and so tiny, \(\widehat\beta\) is a ratio of two such probabilities and stays of order one. (It is not a distribution over \(\mathbf{z}_n\); what is true is \(\sum_{\mathbf{z}_n}\widehat\alpha(\mathbf{z}_n)\widehat\beta(\mathbf{z}_n)\)\({}=1\), by (13.64).)

(13.62). Substitute (13.60) at \(n\) and \(n+1\) into (13.38):

and divide by \(\prod_{m=n+2}^{N}c_m\). At \(n=N\) the product in (13.60) is empty (equal to 1), so \(\widehat\beta(\mathbf{z}_N)\)\({}=\beta(\mathbf{z}_N)\)\({}=1\). The backward pass reuses the \(c_n\) from the forward pass.

(13.63) is (13.57) at \(n=N\). (13.64). Substitute (13.58) and (13.60) into (13.33): \(\gamma(\mathbf{z}_n)\)\({}=\big(\prod_{m\le n}c_m\big)\)\(\big(\prod_{m\gt n}c_m\big)\)\(\widehat\alpha(\mathbf{z}_n)\widehat\beta(\mathbf{z}_n)\)\(/\prod_{m\le N}c_m\)\({}=\widehat\alpha(\mathbf{z}_n)\widehat\beta(\mathbf{z}_n)\). (13.65). Substitute (13.58) at \(n-1\) and (13.60) at \(n\) into (13.43): the factors \(\prod_{m\le n-1}c_m\) and \(\prod_{m\ge n+1}c_m\) cancel against \(p(\mathbf{X})=\prod_{m\le N}c_m\), leaving exactly one factor \(c_n\) in the denominator. ∎

A conditional independence. Result 13.11 at \(n+1\) says \(\mathcal{P}_{n+1}\)\({}=(\mathbf{z}_{1:n},\mathbf{x}_{1:n})\) is independent of \((\mathbf{x}_{n+1},\mathcal{F}_{n+1})\) given \(\mathbf{z}_{n+1}\). Moving \(\mathbf{x}_{1:n}\) into the conditioning set (fact (ii) in the proof of Result 13.12): \(\mathbf{z}_n\perp\!\!\!\perp\mathbf{x}_{n+1:N}\given(\mathbf{z}_{n+1},\mathbf{x}_{1:n})\). Hence \(p(\mathbf{z}_n\given\mathbf{z}_{n+1},\mathbf{X})\)\({}=p(\mathbf{z}_n\given\mathbf{z}_{n+1},\mathbf{x}_{1:n})\).

Its value. By Bayes' theorem, \(p(\mathbf{z}_n\given\mathbf{z}_{n+1},\mathbf{x}_{1:n})\)\({}\propto p(\mathbf{z}_{n+1}\given\mathbf{z}_n,\mathbf{x}_{1:n})\,p(\mathbf{z}_n\given\mathbf{x}_{1:n})\)\({}=p(\mathbf{z}_{n+1}\given\mathbf{z}_n)\,\widehat\alpha(\mathbf{z}_n)\), using \(\mathbf{z}_{n+1}\perp\!\!\!\perp\mathbf{x}_{1:n}\given\mathbf{z}_n\) (Result 13.11 at \(n\)). Normalising over \(\mathbf{z}_n\) gives the fraction.

The recursion. \(\gamma(\mathbf{z}_n)\)\({}=\sum_{\mathbf{z}_{n+1}}p(\mathbf{z}_n,\mathbf{z}_{n+1}\given\mathbf{X})\)\({}=\sum_{\mathbf{z}_{n+1}}p(\mathbf{z}_n\given\mathbf{z}_{n+1},\mathbf{X})\,\gamma(\mathbf{z}_{n+1})\). At \(n=N\), \(\gamma(\mathbf{z}_N)\)\({}=p(\mathbf{z}_N\given\mathbf{x}_{1:N})\)\({}=\widehat\alpha(\mathbf{z}_N)\). The summand is \(\xi(\mathbf{z}_n,\mathbf{z}_{n+1})\), so the pairwise marginals come out of the same pass. ∎

Two rules for maxima. (R1) If \(a\) does not depend on \(u\), \(\max_u\{a+g(u)\}=a+\max_ug(u)\). (R2) Maxima over several variables can be taken one at a time: \(\max_{u,v}g(u,v)=\max_u\max_vg(u,v)\). Together they let a maximisation be pushed past every term that does not involve the variable, exactly as the distributive law pushes sums past factors (with \((\max,+)\) in place of \((+,\times)\)).

Recursion. Define \(\omega\) by (13.70). For \(n\ge2\), the joint of the first \(n\) steps factorises as \(p(\mathbf{x}_{1:n},\mathbf{z}_{1:n})\)\({}=p(\mathbf{x}_{1:n-1},\mathbf{z}_{1:n-1})\,p(\mathbf{z}_n\given\mathbf{z}_{n-1})\,p(\mathbf{x}_n\given\mathbf{z}_n)\) (from (13.6) for \(n\) steps, Result 13.11). Taking logs and using R2 then R1 twice,

which is (13.68) with \(n\) in place of \(n+1\). At \(n=1\) there is nothing to maximise and \(\omega(\mathbf{z}_1)\)\({}=\ln p(\mathbf{x}_1,\mathbf{z}_1)\)\({}=\ln p(\mathbf{z}_1)\)\({}+\ln p(\mathbf{x}_1\given\mathbf{z}_1)\): (13.69). Finally \(\max_{\mathbf{z}_N}\omega(\mathbf{z}_N)\)\({}=\max_{\mathbf{Z}}\ln p(\mathbf{X},\mathbf{Z})\) by R2, and since \(p(\mathbf{Z}\given\mathbf{X})\)\({}=p(\mathbf{X},\mathbf{Z})/p(\mathbf{X})\) with \(p(\mathbf{X})\) not depending on \(\mathbf{Z}\), and \(\ln\) is increasing, the maximiser is also the most probable sequence a posteriori.

As max-sum. In max-sum (8.93)–(8.94) a factor sends \(\max\) over its other variables of \(\ln f\) plus incoming messages, and a variable sends the sum of messages from its other neighbours. A variable node \(\mathbf{z}_n\) on the chain has one other neighbour, so it passes its message on unchanged: (13.66). Factor \(f_{n+1}\) then sends (13.67). Eliminating (13.66) and using \(\ln f_{n+1}\)\({}=\ln p(\mathbf{z}_{n+1}\given\mathbf{z}_n)\)\({}+\ln p(\mathbf{x}_{n+1}\given\mathbf{z}_{n+1})\) from (13.46), the emission term does not depend on \(\mathbf{z}_n\) and comes out of the max (R1): (13.68). The leaf factor \(h\) sends \(\ln h(\mathbf{z}_1)\) by (8.95)–(8.96) and (13.45), which is (13.69). ∎

Claim: for every \(n\) and \(k\), the path \(\mathbf{z}^{(n,k)}\) obtained by setting \(z_n=k\) and following the pointers back attains the maximum in (13.70) for \(\mathbf{z}_n=k\), i.e. \(\ln p(\mathbf{x}_{1:n},\mathbf{z}^{(n,k)})\)\({}=\omega_{nk}\). For \(n=1\) this is (13.69). If it holds at \(n-1\), let \(j=\psi_n(k)\); the path \(\mathbf{z}^{(n,k)}\) is \(\mathbf{z}^{(n-1,j)}\) followed by \(k\), so by the factorisation of the joint, \(\ln p(\mathbf{x}_{1:n},\mathbf{z}^{(n,k)})\)\({}=\omega_{n-1,j}+\ln A_{jk}\)\({}+\ln p(\mathbf{x}_n\given\bm\phi_k)\)\({}=\omega_{nk}\) by the choice of \(j\) and (13.68). At \(n=N\) and \(k=k_N\), \(\omega_{N,k_N}\)\({}=\max_{\mathbf{Z}}\ln p(\mathbf{X},\mathbf{Z})\). This is the "keep only the best path into each state" argument of p. 631: any path into state \(k\) at time \(n\) that is not the best one into that state can be replaced by the best one without changing the rest of the path, and the total can only increase. ∎

(a) For a fixed guess \(\widehat{\mathbf{Z}}\), the probability that it is entirely right is \(p(\mathbf{Z}=\widehat{\mathbf{Z}}\given\mathbf{X})\), maximised by the posterior mode, which Viterbi finds. The expected number of correct positions is, by linearity, \(\sum_n p(\mathbf{z}_n=\widehat{\mathbf{z}}_n\given\mathbf{X})\)\({}=\sum_n\gamma(\widehat{\mathbf{z}}_n)\); each term involves only its own \(\widehat{\mathbf{z}}_n\), so the sum is maximised term by term by \(\widehat{\mathbf{z}}_n\)\({}=\argmax\gamma(\mathbf{z}_n)\). The two objectives are different (0-1 loss on the sequence against Hamming loss), so their optima differ in general; part (b) is an instance.

(b) The joint is \(p(\mathbf{X},z_1=i,z_2=j)\)\({}=\pi_iB_{ia}A_{ij}B_{jb}\)\({}=\frac13\cdot B_{ia}A_{ij}\cdot\frac12\), because \(B_{jb}=0.5\) for every \(j\). The products \(B_{ia}A_{ij}\) are \(0.4\cdot0.75=0.30\), \(0.4\cdot0.25=0.10\), \(0.35\cdot1=0.35\), \(0.25\cdot1=0.25\), and 0 elsewhere; they sum to 1, so \(p(\mathbf{X})=\frac16\) and the posterior table is exactly these products. Row and column sums give \(\gamma\). Viterbi by hand (in probabilities rather than logs): \(e^{\omega(\mathbf{z}_1)}=\frac13(0.40,0.35,0.25)\). For \(\mathbf{z}_2=1\) the best predecessor is state 1 (\(0.75\cdot0.40/3=0.10\)); for \(\mathbf{z}_2=2\) it is state 2 (\(1\cdot0.35/3\) beats \(1\cdot0.25/3\)); for \(\mathbf{z}_2=3\) it is state 1 (\(0.25\cdot0.40/3\)). Multiplying by \(B_{\cdot b}=0.5\), \(e^{\omega(\mathbf{z}_2)}\)\({}=(0.050,\,0.0583,\,0.0167)\), so \(k_2=2\), \(\psi_2(2)=2\), and the Viterbi path is \((2,2)\) with \(p(\mathbf{X},\mathbf{Z})=0.35/6\), i.e. posterior \(0.35\). (Checked numerically against brute-force enumeration.) The course code shows the same phenomenon in a three-state left-to-right model with observations \(x=(0,0,1)\): the marginal argmaxes give the path \(1,1,3\), which skips a state that cannot be skipped and has probability zero, while Viterbi returns \(1,2,3\). ∎

Given \(\mathbf{z}_n=k\), by the Markov property the next transitions are independent draws from row \(k\) of \(\mathbf{A}\): stay with probability \(A_{kk}\), leave with probability \(1-A_{kk}\). "Stay exactly \(T\) more steps, then leave" has probability \(A_{kk}^T(1-A_{kk})\), and \(A_{kk}^T=\exp(T\ln A_{kk})\). Summing, \(\sum_{T\ge0}A_{kk}^T(1-A_{kk})=1\), and \(\E[T+1]\)\({}=\sum_{T\ge0}(T+1)A_{kk}^T(1-A_{kk})\)\({}=1/(1-A_{kk})\) (differentiate the geometric series). For \(A_{kk}=0.9\) the mean stay is 10 steps, as in Figure 13.8. ∎

Factor graph (Exercise 13.17). The joint of the input-output HMM (conditional on \(\mathbf{U}\)) is \(p(\mathbf{X},\mathbf{Z}\given\mathbf{U})\)\({}=p(\mathbf{z}_1\given\mathbf{u}_1)\)\(\prod_{n\ge2}p(\mathbf{z}_n\given\mathbf{z}_{n-1},\mathbf{u}_n)\)\(\prod_{n\ge1}p(\mathbf{x}_n\given\mathbf{z}_n,\mathbf{u}_n)\). With \(\mathbf{U}\) and \(\mathbf{X}\) fixed, each factor depends on at most two successive latent variables; grouping as in (13.45)–(13.46) gives \(h\cdot\prod_nf_n\), a chain factor graph exactly like Figure 13.15. The same holds for the autoregressive HMM, because \(\mathbf{x}_{n-1},\mathbf{x}_{n-2}\) are fixed numbers in the E step.

Markov property. The posterior \(p(\mathbf{Z}\given\mathbf{X},\mathbf{U})\)\({}\propto h(\mathbf{z}_1)\prod_nf_n(\mathbf{z}_{n-1},\mathbf{z}_n)\) is a chain, so (13.5) holds for it (the same factorisation argument as in Result 13.11). By d-separation in Figure 13.18: \(\mathbf{u}_n\) and \(\mathbf{x}_n\) are connected only to \(\mathbf{z}_n\) and to each other, so removing \(\mathbf{z}_n\) cuts the graph into a past part and a future part; hence every path from \(\mathbf{z}_{n-1}\) to \(\mathbf{z}_{n+1}\) passes through \(\mathbf{z}_n\), and the only one is \(\mathbf{z}_{n-1}\)\({}\to\mathbf{z}_n\)\({}\to\mathbf{z}_{n+1}\), head-to-tail at the observed \(\mathbf{z}_n\). In Figure 13.17 there are also paths through the skip links, for instance \(\mathbf{z}_{n-1}\)\({}\to\mathbf{x}_{n-1}\)\({}\to\mathbf{x}_{n+1}\leftarrow\mathbf{z}_{n+1}\); as the book notes, every such path meets some observed \(\mathbf{x}\) head-to-tail (here \(\mathbf{x}_{n-1}\)), so it is blocked.

Recursions (Exercise 13.18). Apply the sum-product derivation of Result 13.18 to the factors \(h,f_n\): the forward message is \(\alpha(\mathbf{z}_n)\)\({}=\sum_{\mathbf{z}_{n-1}}f_n(\mathbf{z}_{n-1},\mathbf{z}_n)\alpha(\mathbf{z}_{n-1})\) with \(\alpha(\mathbf{z}_1)=h(\mathbf{z}_1)\), the backward message \(\beta(\mathbf{z}_n)\)\({}=\sum_{\mathbf{z}_{n+1}}f_{n+1}(\mathbf{z}_n,\mathbf{z}_{n+1})\beta(\mathbf{z}_{n+1})\) with \(\beta(\mathbf{z}_N)=1\), and the marginals are as stated. By the definition (8.64), \(\alpha(\mathbf{z}_n)\)\({}=\sum_{\mathbf{z}_{1:n-1}}h\prod_{m\le n}f_m\)\({}=p(\mathbf{x}_{1:n},\mathbf{z}_n\given\mathbf{U})\): the factors involve only \(\mathbf{u}\)'s and \(\mathbf{x}\)'s up to time \(n\), and the later factors sum to one. Scaling and Viterbi carry over unchanged. ∎

(a) \(\mathbf{x}_n\) is a common child of \(\mathbf{z}_n^{(1)}\) and \(\mathbf{z}_n^{(2)}\), observed: the path \(\mathbf{z}_n^{(1)}\)\({}\to\mathbf{x}_n\leftarrow\mathbf{z}_n^{(2)}\) is head-to-head at an observed node, so it is unblocked ("explaining away"), and the posterior couples the chains. For the Markov property of chain 1, the green path of Figure 13.20 connects \(\mathbf{z}_{n-1}^{(1)}\) and \(\mathbf{z}_{n+1}^{(1)}\) without passing through \(\mathbf{z}_n^{(1)}\), is head-to-head only at observed nodes and passes unobserved nodes head-to-tail, so it is not blocked by \(\{\mathbf{z}_n^{(1)}\}\cup\mathbf{X}\).

(b) The joint is \(\prod_m\big[p(\mathbf{z}_1^{(m)})\prod_np(\mathbf{z}_n^{(m)}\given\mathbf{z}_{n-1}^{(m)})\big]\)\(\prod_np(\mathbf{x}_n\given\mathbf{s}_n)\). Since the chains evolve independently, \(p(\mathbf{s}_n\given\mathbf{s}_{n-1})\)\({}=\prod_mp(\mathbf{z}_n^{(m)}\given\mathbf{z}_{n-1}^{(m)})\), which in matrix form is the Kronecker product, and the emission depends only on \(\mathbf{s}_n\): this is the HMM (13.10) with \(K^M\) states. Result 13.16 gives \(O(N(K^M)^2)\). With structure, multiply \(\bm\alpha_{n-1}\) (a \(K\times\cdots\times K\) array) by one \(\mathbf{A}^{(m)}\) at a time along its \(m\)-th axis: each multiplication costs \(K^M\cdot K\), and there are \(M\) of them per step, hence \(O(NMK^{M+1})\). Both are exponential in \(M\).

(c) \(\alpha(\mathbf{s}_n)\)\({}=p(\mathbf{x}_{1:n},\mathbf{s}_n)\) is a general function of the \(M\) components: even at \(n=1\), \(\alpha(\mathbf{s}_1)\)\({}=\prod_m\pi^{(m)}(\mathbf{z}_1^{(m)})\,p(\mathbf{x}_1\given\mathbf{z}_1^{(1)},\dots,\mathbf{z}_1^{(M)})\), and the emission term is an arbitrary function of all components (for the linear-Gaussian emission above, its logarithm contains the cross terms \(-\mathbf{z}^{(m)\T}\mathbf{W}^{(m)\T}\bm\Sigma^{-1}\mathbf{W}^{(m')}\mathbf{z}^{(m')}\) for every pair \(m\ne m'\)). A function of \(M\) variables with \(K\) values each and no factorisation needs \(K^M\) numbers. ∎

Unrolling (13.78)–(13.80), \(\mathbf{z}_n\)\({}=\mathbf{A}^{n-1}(\bm\mu_0+\mathbf{u})\)\({}+\sum_{m=2}^{n}\mathbf{A}^{n-m}\mathbf{w}_m\) and \(\mathbf{x}_n=\mathbf{C}\mathbf{z}_n+\mathbf{v}_n\). So the vector of all \(\mathbf{z}\)'s and \(\mathbf{x}\)'s is a constant plus a fixed matrix times the vector \((\mathbf{u},\mathbf{w}_2,\dots,\mathbf{w}_N,\mathbf{v}_1,\dots,\mathbf{v}_N)\), which is Gaussian (independent Gaussian blocks). An affine map of a Gaussian vector is Gaussian, and marginals and conditionals of a Gaussian are Gaussian (toolkit, Result T.14). For (b), take expectations of the unrolled form (the noises have mean zero); \(\mathbf{z}_n\)\({}=\mathbf{A}\mathbf{z}_{n-1}+\mathbf{w}_n\) with \(\mathbf{w}_n\) independent of \(\mathbf{z}_{n-1}\) gives \(\mathbf{S}_n\)\({}=\mathbf{A}\mathbf{S}_{n-1}\mathbf{A}^\T\)\({}+\bm\Gamma\); and for \(n\gt m\), \(\mathbf{z}_n-\E\mathbf{z}_n\)\({}=\mathbf{A}^{n-m}(\mathbf{z}_m-\E\mathbf{z}_m)\)\({}+(\text{noise independent of }\mathbf{z}_m)\), so \(\cov[\mathbf{z}_n,\mathbf{z}_m]\)\({}=\mathbf{A}^{n-m}\mathbf{S}_m\). The \(\mathbf{x}\) moments follow because \(\mathbf{v}_n\) is independent of everything else. ∎

By Result 13.27, \(p(\mathbf{Z}\given\mathbf{X})\) is a Gaussian with some mean \(\mathbf{m}=(\mathbf{m}_1,\dots,\mathbf{m}_N)\). A Gaussian density is largest at its mean (the quadratic form in the exponent is zero there and positive elsewhere), so the joint mode is \(\mathbf{m}\). The marginal of the block \(\mathbf{z}_n\) is Gaussian with mean \(\mathbf{m}_n\), the corresponding block of \(\mathbf{m}\) (the marginal of a partitioned Gaussian keeps its sub-vector of the mean, (2.98)), so its mode is \(\mathbf{m}_n\). Hence the marginal modes, put together, are the joint mode. For the HMM this fails (Result 13.23) because a discrete distribution has no such structure. ∎

From the linear-Gaussian result. Result T.15 (the book's (2.115)): if \(p(\mathbf{x})\)\({}=\N(\mathbf{x}\given\bm\mu,\bm\Lambda^{-1})\) and \(p(\mathbf{y}\given\mathbf{x})\)\({}=\N(\mathbf{y}\given\mathbf{A}\mathbf{x}+\mathbf{b},\mathbf{L}^{-1})\), then \(p(\mathbf{y})\)\({}=\int p(\mathbf{y}\given\mathbf{x})p(\mathbf{x})\d\mathbf{x}\)\({}=\N(\mathbf{y}\given\mathbf{A}\bm\mu+\mathbf{b},\mathbf{L}^{-1}+\mathbf{A}\bm\Lambda^{-1}\mathbf{A}^\T)\). Substitute \(\mathbf{x}\to\mathbf{z}_{n-1}\), \(\bm\mu\to\bm\mu_{n-1}\), \(\bm\Lambda^{-1}\to\mathbf{V}_{n-1}\), \(\mathbf{y}\to\mathbf{z}_n\), \(\mathbf{b}\to\mathbf{0}\), \(\mathbf{L}^{-1}\to\bm\Gamma\).

From the noise equation. If \(\mathbf{z}_{n-1}\)\({}\sim\N(\bm\mu_{n-1},\mathbf{V}_{n-1})\) and \(\mathbf{w}_n\sim\N(\mathbf{0},\bm\Gamma)\) independently, then \(\mathbf{z}_n\)\({}=\mathbf{A}\mathbf{z}_{n-1}+\mathbf{w}_n\) is Gaussian (an affine image of a Gaussian vector) with mean \(\mathbf{A}\bm\mu_{n-1}\) and covariance \(\mathbf{A}\mathbf{V}_{n-1}\mathbf{A}^\T\)\({}+\bm\Gamma\) (covariances of independent terms add). The integral in (13.87) is exactly the density of this \(\mathbf{z}_n\) (the sum rule).

Why it is the prediction. Under the posterior given \(\mathbf{x}_{1:n-1}\), \(\mathbf{z}_{n-1}\)\({}\sim\widehat\alpha(\mathbf{z}_{n-1})\), and \(p(\mathbf{z}_n\given\mathbf{z}_{n-1},\mathbf{x}_{1:n-1})\)\({}=p(\mathbf{z}_n\given\mathbf{z}_{n-1})\) by (13.26) (density version), so the integral is \(p(\mathbf{z}_n\given\mathbf{x}_{1:n-1})\). ∎

Write \(\mathbf{m}=\mathbf{A}\bm\mu_{n-1}\), \(\mathbf{P}=\mathbf{P}_{n-1}\) and \(\mathbf{x}=\mathbf{x}_n\). By (13.86)–(13.87), as functions of \(\mathbf{z}=\mathbf{z}_n\),

Step 1: the exponent. The logarithm of the right-hand side is, up to terms not involving \(\mathbf{z}\),

(expand both quadratics; the cross terms \(-\frac12(-2\mathbf{z}^\T\mathbf{C}^\T\bm\Sigma^{-1}\mathbf{x})\) and \(-\frac12(-2\mathbf{z}^\T\mathbf{P}^{-1}\mathbf{m})\) give the linear term). By completing the square, a density whose log is \(-\frac12\mathbf{z}^\T\mathbf{Q}\mathbf{z}\)\({}+\mathbf{z}^\T\mathbf{b}+\const\) is \(\N(\mathbf{Q}^{-1}\mathbf{b},\mathbf{Q}^{-1})\). Hence

(This is also Result T.15, the book's (2.116), with \(\bm\Lambda=\mathbf{P}^{-1}\), \(\mathbf{L}=\bm\Sigma^{-1}\), \(\mathbf{A}\to\mathbf{C}\), \(\mathbf{b}=\mathbf{0}\).)

Step 2: the covariance in gain form. The Woodbury identity (C.7), \((\mathbf{P}^{-1}+\mathbf{C}^\T\bm\Sigma^{-1}\mathbf{C})^{-1}\)\({}=\mathbf{P}\)\({}-\mathbf{P}\mathbf{C}^\T(\bm\Sigma+\mathbf{C}\mathbf{P}\mathbf{C}^\T)^{-1}\mathbf{C}\mathbf{P}\), gives \(\mathbf{V}_n\)\({}=\mathbf{P}-\mathbf{K}_n\mathbf{C}\mathbf{P}\)\({}=(\mathbf{I}-\mathbf{K}_n\mathbf{C})\mathbf{P}\) with \(\mathbf{K}_n\) as in (13.92): (13.90). Since \(\mathbf{K}_n\mathbf{C}\mathbf{P}\)\({}=\mathbf{P}\mathbf{C}^\T\mathbf{S}_n^{-1}\mathbf{C}\mathbf{P}\)\({}=\mathbf{K}_n\mathbf{S}_n\mathbf{K}_n^\T\), also \(\mathbf{V}_n\)\({}=\mathbf{P}\)\({}-\mathbf{K}_n\mathbf{S}_n\mathbf{K}_n^\T\).

Step 3: the gain identity (C.5). \(\mathbf{V}_n\mathbf{C}^\T\bm\Sigma^{-1}\)\({}=\mathbf{K}_n\), i.e. \((\mathbf{P}^{-1}+\mathbf{C}^\T\bm\Sigma^{-1}\mathbf{C})^{-1}\mathbf{C}^\T\bm\Sigma^{-1}\)\({}=\mathbf{P}\mathbf{C}^\T(\mathbf{C}\mathbf{P}\mathbf{C}^\T+\bm\Sigma)^{-1}\). Multiply on the left by \(\mathbf{P}^{-1}\)\({}+\mathbf{C}^\T\bm\Sigma^{-1}\mathbf{C}\) and on the right by \(\mathbf{C}\mathbf{P}\mathbf{C}^\T+\bm\Sigma\): the claim becomes \(\mathbf{C}^\T\bm\Sigma^{-1}(\mathbf{C}\mathbf{P}\mathbf{C}^\T+\bm\Sigma)\)\({}=(\mathbf{P}^{-1}+\mathbf{C}^\T\bm\Sigma^{-1}\mathbf{C})\mathbf{P}\mathbf{C}^\T\), and both sides equal \(\mathbf{C}^\T\bm\Sigma^{-1}\mathbf{C}\mathbf{P}\mathbf{C}^\T\)\({}+\mathbf{C}^\T\).

Step 4: the mean. From Step 2, \(\mathbf{V}_n\mathbf{P}^{-1}\)\({}=(\mathbf{I}-\mathbf{K}_n\mathbf{C})\mathbf{P}\mathbf{P}^{-1}\)\({}=\mathbf{I}-\mathbf{K}_n\mathbf{C}\); with Step 3,

which is (13.89).

Step 5: the normaliser. Integrate both sides of the first display over \(\mathbf{z}\): the left gives \(c_n\), the right the marginal of \(\mathbf{x}\) in the linear-Gaussian model, \(\N(\mathbf{x}\given\mathbf{C}\mathbf{m},\bm\Sigma+\mathbf{C}\mathbf{P}\mathbf{C}^\T)\) (Result T.15 again): (13.91). As in (13.56), \(c_n=p(\mathbf{x}_n\given\mathbf{x}_{1:n-1})\); since \(\mathbf{e}_n\) is \(\mathbf{x}_n\) minus a function of \(\mathbf{x}_{1:n-1}\), it is \(\N(\mathbf{0},\mathbf{S}_n)\) given the past.

Joseph form. Expand \((\mathbf{I}-\mathbf{K}\mathbf{C})\mathbf{P}(\mathbf{I}-\mathbf{K}\mathbf{C})^\T\)\({}+\mathbf{K}\bm\Sigma\mathbf{K}^\T\)\({}=(\mathbf{I}-\mathbf{K}\mathbf{C})\mathbf{P}\)\({}-(\mathbf{I}-\mathbf{K}\mathbf{C})\mathbf{P}\mathbf{C}^\T\mathbf{K}^\T\)\({}+\mathbf{K}\bm\Sigma\mathbf{K}^\T\); the last two terms are \(-\mathbf{P}\mathbf{C}^\T\mathbf{K}^\T\)\({}+\mathbf{K}(\mathbf{C}\mathbf{P}\mathbf{C}^\T+\bm\Sigma)\mathbf{K}^\T\)\({}=-\mathbf{K}\mathbf{S}\mathbf{K}^\T\)\({}+\mathbf{K}\mathbf{S}\mathbf{K}^\T\)\({}=\mathbf{0}\). (All forms were checked numerically against each other and against conditioning the joint Gaussian.) ∎

(13.93) is the definition \(\widehat\alpha(\mathbf{z}_1)\)\({}=p(\mathbf{z}_1\given\mathbf{x}_1)\)\({}=p(\mathbf{z}_1)p(\mathbf{x}_1\given\mathbf{z}_1)/p(\mathbf{x}_1)\) with \(c_1=p(\mathbf{x}_1)\). Its right-hand side is \(\N(\mathbf{x}_1\given\mathbf{C}\mathbf{z}_1,\bm\Sigma)\N(\mathbf{z}_1\given\bm\mu_0,\mathbf{P}_0)\): the same product as in Result 13.30 with the prediction \((\mathbf{m},\mathbf{P})\) replaced by the prior \((\bm\mu_0,\mathbf{P}_0)\). Steps 1–5 of that proof, with this substitution, give (13.94), (13.95), (13.97) (from (2.116): Exercise 13.23) and (13.96) (from (2.115): Exercise 13.22). Equivalently: start the general recursion from a fictitious \(\bm\mu_{0},\mathbf{V}_0\) with \(\mathbf{A}=\mathbf{I}\), \(\bm\Gamma=\mathbf{0}\) at the first step. ∎

(a) With \(\bm\Sigma=\mathbf{0}\) and \(\mathbf{C}=\mathbf{I}\), (13.92) is \(\mathbf{K}_n\)\({}=\mathbf{P}_{n-1}\mathbf{P}_{n-1}^{-1}\)\({}=\mathbf{I}\) (for \(\mathbf{P}_{n-1}\) invertible, which holds when \(\bm\Gamma\succ0\)). Then (13.89) gives \(\bm\mu_n\)\({}=\mathbf{A}\bm\mu_{n-1}+\mathbf{x}_n\)\({}-\mathbf{A}\bm\mu_{n-1}\)\({}=\mathbf{x}_n\) and (13.90) gives \(\mathbf{V}_n=\mathbf{0}\): the posterior is a point mass at \(\mathbf{x}_n\). (Directly: the emission density in (13.86) collapses to a spike at \(\mathbf{z}_n=\mathbf{x}_n\).)

(b) Induction on \(n\), with \(\mathbf{P}_0=\lambda\mathbf{I}\), \(\lambda\to\infty\). Base: \(\mathbf{K}_1\)\({}=\lambda(\lambda\mathbf{I}+\bm\Sigma)^{-1}\)\({}\to\mathbf{I}\), so \(\bm\mu_1\to\mathbf{x}_1\); \(\mathbf{V}_1\)\({}=(\mathbf{I}-\mathbf{K}_1)\lambda\mathbf{I}\)\({}=\bm\Sigma(\lambda\mathbf{I}+\bm\Sigma)^{-1}\lambda\)\({}\to\bm\Sigma\) (using \(\mathbf{I}\)\({}-\lambda(\lambda\mathbf{I}+\bm\Sigma)^{-1}\)\({}=\bm\Sigma(\lambda\mathbf{I}+\bm\Sigma)^{-1}\)). Step: suppose \(\bm\mu_{n-1}=\bar{\mathbf{x}}_{n-1}\) and \(\mathbf{V}_{n-1}=\bm\Sigma/(n-1)\). With \(\mathbf{A}=\mathbf{I}\), \(\bm\Gamma=\mathbf{0}\), \(\mathbf{P}_{n-1}=\mathbf{V}_{n-1}\), so \(\mathbf{K}_n\)\({}=\frac{\bm\Sigma}{n-1}\big(\frac{\bm\Sigma}{n-1}+\bm\Sigma\big)^{-1}\)\({}=\frac{\bm\Sigma}{n-1}\cdot\frac{n-1}{n}\bm\Sigma^{-1}\)\({}=\frac1n\mathbf{I}\). Then \(\bm\mu_n\)\({}=\bar{\mathbf{x}}_{n-1}\)\({}+\frac1n(\mathbf{x}_n-\bar{\mathbf{x}}_{n-1})\)\({}=\bar{\mathbf{x}}_n\) (the running-mean update) and \(\mathbf{V}_n\)\({}=(1-\frac1n)\frac{\bm\Sigma}{n-1}\)\({}=\frac{\bm\Sigma}{n}\).

(c) Here \(P_{n-1}=V_{n-1}\) and \(K_n=V_{n-1}/(V_{n-1}+\sigma^2)\). From (13.90), \(V_n\)\({}=(1-K_n)V_{n-1}\)\({}=\sigma^2V_{n-1}/(V_{n-1}+\sigma^2)\), i.e. \(1/V_n=1/V_{n-1}+1/\sigma^2\); with \(V_0:=\sigma_0^2\) (so that (13.94)–(13.95) are the case \(n=1\)), induction gives \(1/V_n=1/\sigma_0^2+n/\sigma^2\). For the mean, \(1-K_n=V_n/V_{n-1}\) and \(K_n=V_n/\sigma^2\), so (13.89) reads \(\mu_n/V_n=\mu_{n-1}/V_{n-1}+x_n/\sigma^2\), and induction gives the stated \(\mu_n\). With \(\mu_{\rm ML}=\frac1n\sum x_m\) this is \(\mu_n\)\({}=\frac{\sigma^2}{n\sigma_0^2+\sigma^2}\mu_0\)\({}+\frac{n\sigma_0^2}{n\sigma_0^2+\sigma^2}\mu_{\rm ML}\), which is (2.141), and the variance is (2.142). Sequential Bayesian updating of a static parameter is a Kalman filter. ∎

Step 1: turn (13.99) into a recursion for \(\gamma\). Multiply both sides of (13.99) by \(\widehat\alpha(\mathbf{z}_n)\). On the left, \(c_{n+1}\widehat\alpha(\mathbf{z}_n)\widehat\beta(\mathbf{z}_n)\)\({}=c_{n+1}\gamma(\mathbf{z}_n)\). On the right, write \(\widehat\beta(\mathbf{z}_{n+1})\)\({}=\gamma(\mathbf{z}_{n+1})/\widehat\alpha(\mathbf{z}_{n+1})\), and use the forward recursion (13.85) at \(n+1\), \(c_{n+1}\widehat\alpha(\mathbf{z}_{n+1})\)\({}=p(\mathbf{x}_{n+1}\given\mathbf{z}_{n+1})\,p(\mathbf{z}_{n+1}\given\mathbf{x}_{1:n})\), where \(p(\mathbf{z}_{n+1}\given\mathbf{x}_{1:n})\)\({}=\N(\mathbf{z}_{n+1}\given\mathbf{A}\bm\mu_n,\mathbf{P}_n)\) is the prediction (13.87). The emission factors cancel and \(c_{n+1}\) cancels:

The fraction is \(p(\mathbf{z}_n\given\mathbf{z}_{n+1},\mathbf{x}_{1:n})\) by Bayes' theorem (the numerator is \(p(\mathbf{z}_n\given\mathbf{x}_{1:n})p(\mathbf{z}_{n+1}\given\mathbf{z}_n)\)\({}=p(\mathbf{z}_n,\mathbf{z}_{n+1}\given\mathbf{x}_{1:n})\)), so this is exactly the \(\alpha\)-\(\gamma\) recursion of Result 13.21.

Step 2: the backward conditional. \(p(\mathbf{z}_n\given\mathbf{z}_{n+1},\mathbf{x}_{1:n})\)\({}\propto\N(\mathbf{z}_n\given\bm\mu_n,\mathbf{V}_n)\,\N(\mathbf{z}_{n+1}\given\mathbf{A}\mathbf{z}_n,\bm\Gamma)\). This is the linear-Gaussian model with "\(\mathbf{x}\)" \(=\mathbf{z}_n\) and "\(\mathbf{y}\)" \(=\mathbf{z}_{n+1}\); by Result T.15 (or completing the square as in Result 13.30) the conditional of \(\mathbf{z}_n\) is Gaussian with

The covariance step is Woodbury (C.7) with \(\mathbf{P}_n\)\({}=\bm\Gamma+\mathbf{A}\mathbf{V}_n\mathbf{A}^\T\), and \(\mathbf{J}_n\mathbf{P}_n\mathbf{J}_n^\T\)\({}=\mathbf{V}_n\mathbf{A}^\T\mathbf{P}_n^{-1}\mathbf{P}_n\mathbf{P}_n^{-1}\mathbf{A}\mathbf{V}_n\) is the same matrix; the mean is the same gain argument as Steps 3–4 of Result 13.30 (here the "gain" is \(\mathbf{J}_n\) and the "measurement" is \(\mathbf{z}_{n+1}\)). In other words: given \(\mathbf{z}_{n+1}\), \(\mathbf{z}_n\)\({}=\bm\mu_n\)\({}+\mathbf{J}_n(\mathbf{z}_{n+1}-\mathbf{A}\bm\mu_n)\)\({}+\bm\epsilon_n\) with \(\bm\epsilon_n\)\({}\sim\N(\mathbf{0},\mathbf{V}_n-\mathbf{J}_n\mathbf{P}_n\mathbf{J}_n^\T)\) independent of \(\mathbf{z}_{n+1}\).

Step 3: integrate against \(\gamma(\mathbf{z}_{n+1})\). If \(\mathbf{z}_{n+1}\)\({}\sim\N(\widehat{\bm\mu}_{n+1},\widehat{\mathbf{V}}_{n+1})\), the representation in Step 2 makes \(\mathbf{z}_n\) Gaussian with

which are (13.100) and (13.101). (This is Result T.15's marginal once more.) At \(n=N\), \(\gamma(\mathbf{z}_N)\)\({}=\widehat\alpha(\mathbf{z}_N)\) because \(\widehat\beta(\mathbf{z}_N)=1\), so \(\widehat{\bm\mu}_N=\bm\mu_N\), \(\widehat{\mathbf{V}}_N=\mathbf{V}_N\).

The identity \(\mathbf{A}\mathbf{V}_n\)\({}=\mathbf{P}_n\mathbf{J}_n^\T\). \(\mathbf{J}_n^\T\)\({}=\mathbf{P}_n^{-1}\mathbf{A}\mathbf{V}_n\) because \(\mathbf{V}_n\) and \(\mathbf{P}_n\) are symmetric; multiply by \(\mathbf{P}_n\). ∎

(13.103) (Exercise 13.30). The first expression is (13.65), valid for densities. Substitute (13.84) at \(n-1\), (13.75), (13.76), and \(\widehat\beta(\mathbf{z}_n)\)\({}=\gamma(\mathbf{z}_n)/\widehat\alpha(\mathbf{z}_n)\)\({}=\N(\mathbf{z}_n\given\widehat{\bm\mu}_n,\widehat{\mathbf{V}}_n)/\widehat\alpha(\mathbf{z}_n)\) from (13.98).

Simplify (Exercise 13.31). By (13.85)–(13.87), \(c_n\widehat\alpha(\mathbf{z}_n)\)\({}=\N(\mathbf{x}_n\given\mathbf{C}\mathbf{z}_n,\bm\Sigma)\,\N(\mathbf{z}_n\given\mathbf{A}\bm\mu_{n-1},\mathbf{P}_{n-1})\). The emission factors cancel, leaving

the backward conditional of Step 2 of Result 13.33 (with \(n-1\) in place of \(n\)) times the smoothed marginal of \(\mathbf{z}_n\). So \(\mathbf{z}_{n-1}\)\({}=\bm\mu_{n-1}\)\({}+\mathbf{J}_{n-1}(\mathbf{z}_n-\mathbf{A}\bm\mu_{n-1})\)\({}+\bm\epsilon\) with \(\bm\epsilon\) independent of \(\mathbf{z}_n\)\({}\sim\N(\widehat{\bm\mu}_n,\widehat{\mathbf{V}}_n)\). This pair is jointly Gaussian; the mean of \(\mathbf{z}_{n-1}\) is \(\widehat{\bm\mu}_{n-1}\) by (13.100), its covariance \(\widehat{\mathbf{V}}_{n-1}\) by (13.101), and

(Checked numerically against the cross-covariance block of the full joint Gaussian.) ∎

The terms of \(Q\) containing \(\bm\mu_0,\mathbf{P}_0\) are

with \(\mathbf{R}_0\)\({}=\E[(\mathbf{z}_1-\bm\mu_0)(\mathbf{z}_1-\bm\mu_0)^\T]\)\({}=\E[\mathbf{z}_1\mathbf{z}_1^\T]\)\({}-\E[\mathbf{z}_1]\bm\mu_0^\T\)\({}-\bm\mu_0\E[\mathbf{z}_1]^\T+\bm\mu_0\bm\mu_0^\T\). Mean: \(\partial Q_0/\partial\bm\mu_0\)\({}=\mathbf{P}_0^{-1}(\E[\mathbf{z}_1]-\bm\mu_0)\) (matrix calculus), zero at (13.110); \(Q_0\) is concave in \(\bm\mu_0\) for any \(\mathbf{P}_0\), so this is the maximiser whatever \(\mathbf{P}_0\) is. Covariance: with \(\bm\Lambda=\mathbf{P}_0^{-1}\), \(Q_0\)\({}=\frac12\ln|\bm\Lambda|\)\({}-\frac12\Tr(\bm\Lambda\mathbf{R}_0)\), stationary at \(\bm\Lambda^{-1}=\mathbf{R}_0\) and concave in \(\bm\Lambda\) (as in Result 13.8). At \(\bm\mu_0=\E[\mathbf{z}_1]\), \(\mathbf{R}_0\)\({}=\E[\mathbf{z}_1\mathbf{z}_1^\T]\)\({}-\E[\mathbf{z}_1]\E[\mathbf{z}_1]^\T\), which is (13.111). This is the maximum likelihood fit of a Gaussian to "one observation" \(\mathbf{z}_1\) known only through its posterior. ∎

The objective. (13.112) is the log of (13.75) summed over \(n\ge2\): each term contributes \(-\frac12\ln|\bm\Gamma|\) and the quadratic form. By the trace trick, the expectation is \(\frac12\Tr(\bm\Gamma^{-1}\mathbf{R}(\mathbf{A}))\) with

where \(\mathbf{S}_{11}\)\({}=\sum_{n\ge2}\E[\mathbf{z}_n\mathbf{z}_n^\T]\), \(\mathbf{S}_{10}\)\({}=\sum_{n\ge2}\E[\mathbf{z}_n\mathbf{z}_{n-1}^\T]\), \(\mathbf{S}_{00}\)\({}=\sum_{n\ge2}\E[\mathbf{z}_{n-1}\mathbf{z}_{n-1}^\T]\) (and \(\mathbf{S}_{10}^\T\)\({}=\sum\E[\mathbf{z}_{n-1}\mathbf{z}_n^\T]\)).

Maximise over \(\mathbf{A}\). Using \(\partial\Tr(\mathbf{X}\mathbf{A})/\partial\mathbf{A}\)\({}=\mathbf{X}^\T\), \(\partial\Tr(\mathbf{A}\mathbf{Y})/\partial\mathbf{A}\)\({}=\mathbf{Y}^\T\) and \(\partial\Tr(\mathbf{G}\mathbf{A}\mathbf{S}\mathbf{A}^\T)/\partial\mathbf{A}\)\({}=\mathbf{G}\mathbf{A}\mathbf{S}\)\({}+\mathbf{G}^\T\mathbf{A}\mathbf{S}^\T\)\({}=2\mathbf{G}\mathbf{A}\mathbf{S}\) for symmetric \(\mathbf{G},\mathbf{S}\) (toolkit), with \(\mathbf{G}=\bm\Gamma^{-1}\):

(For the cross terms: \(\Tr(\bm\Gamma^{-1}\mathbf{A}\mathbf{S}_{10}^\T)\)\({}=\Tr(\mathbf{S}_{10}^\T\bm\Gamma^{-1}\mathbf{A})\) has derivative \(\bm\Gamma^{-1}\mathbf{S}_{10}\), and \(\Tr(\bm\Gamma^{-1}\mathbf{S}_{10}\mathbf{A}^\T)\) is the same number.) Setting it to zero, \(\mathbf{A}=\mathbf{S}_{10}\mathbf{S}_{00}^{-1}\): (13.113). The objective is a concave quadratic in \(\mathbf{A}\) (its second-order term is \(-\frac12\Tr(\bm\Gamma^{-1}\mathbf{A}\mathbf{S}_{00}\mathbf{A}^\T)\)\({}\le0\), strictly negative for \(\mathbf{A}\ne\mathbf{0}\) when \(\mathbf{S}_{00}\succ0\)), and the maximiser does not involve \(\bm\Gamma\), which is why \(\mathbf{A}\) can be updated first and then held fixed.

Maximise over \(\bm\Gamma\). With \(\mathbf{A}=\mathbf{A}^{\rm new}\) and \(\bm\Lambda=\bm\Gamma^{-1}\), the objective is \(\frac{N-1}2\ln|\bm\Lambda|\)\({}-\frac12\Tr(\bm\Lambda\mathbf{R})\), concave in \(\bm\Lambda\) with stationary point \(\bm\Gamma=\mathbf{R}(\mathbf{A}^{\rm new})/(N-1)\), which written out is (13.114). Simplification: \(\mathbf{A}^{\rm new}\mathbf{S}_{00}(\mathbf{A}^{\rm new})^\T\)\({}=\mathbf{S}_{10}\mathbf{S}_{00}^{-1}\mathbf{S}_{00}(\mathbf{A}^{\rm new})^\T\)\({}=\mathbf{S}_{10}(\mathbf{A}^{\rm new})^\T\), which cancels the third term, leaving \(\mathbf{R}\)\({}=\mathbf{S}_{11}\)\({}-\mathbf{A}^{\rm new}\mathbf{S}_{10}^\T\). (Both forms, and the fact that small perturbations of \(\mathbf{A}^{\rm new}\) or \(\bm\Gamma^{\rm new}\) decrease \(Q\), were checked numerically.) ∎

Exactly as in Result 13.36, with \(\mathbf{z}_n\) in place of \(\mathbf{z}_{n-1}\), \(\mathbf{x}_n\) in place of \(\mathbf{z}_n\), and \(N\) terms. By the trace trick, the expectation is \(\frac12\Tr(\bm\Sigma^{-1}\mathbf{R}(\mathbf{C}))\) with \(\mathbf{R}(\mathbf{C})\)\({}=\mathbf{T}_{xx}-\mathbf{C}\mathbf{T}_{xz}^\T\)\({}-\mathbf{T}_{xz}\mathbf{C}^\T\)\({}+\mathbf{C}\mathbf{T}_{zz}\mathbf{C}^\T\), \(\mathbf{T}_{xx}=\sum\mathbf{x}_n\mathbf{x}_n^\T\), \(\mathbf{T}_{xz}\)\({}=\sum\mathbf{x}_n\E[\mathbf{z}_n]^\T\), \(\mathbf{T}_{zz}\)\({}=\sum\E[\mathbf{z}_n\mathbf{z}_n^\T]\). (\(\mathbf{x}_n\) is observed, so \(\E[\mathbf{x}_n\mathbf{z}_n^\T]\)\({}=\mathbf{x}_n\E[\mathbf{z}_n]^\T\).) The gradient in \(\mathbf{C}\) is \(\bm\Sigma^{-1}(\mathbf{T}_{xz}-\mathbf{C}\mathbf{T}_{zz})\), zero at \(\mathbf{C}=\mathbf{T}_{xz}\mathbf{T}_{zz}^{-1}\): (13.115), the regression of \(\mathbf{x}\) on \(\mathbf{z}\) with the sufficient statistics replaced by their posterior expectations. Then \(\bm\Sigma=\mathbf{R}(\mathbf{C}^{\rm new})/N\), which is (13.116); and \(\mathbf{C}^{\rm new}\mathbf{T}_{zz}(\mathbf{C}^{\rm new})^\T\)\({}=\mathbf{T}_{xz}(\mathbf{C}^{\rm new})^\T\) cancels the third term. ∎

Replace \(\mathbf{f}\) by its first-order Taylor expansion about the current mean, \(\mathbf{f}(\mathbf{z}_{n-1})\)\({}\approx\mathbf{f}(\bm\mu_{n-1})\)\({}+\mathbf{F}_n(\mathbf{z}_{n-1}-\bm\mu_{n-1})\). The transition becomes linear-Gaussian with an offset, \(\N(\mathbf{z}_n\given\mathbf{F}_n\mathbf{z}_{n-1}+\mathbf{b}_n,\bm\Gamma)\) with \(\mathbf{b}_n\)\({}=\mathbf{f}(\bm\mu_{n-1})\)\({}-\mathbf{F}_n\bm\mu_{n-1}\), and Result 13.29 (with the offset, Result T.15) gives the predicted mean \(\mathbf{F}_n\bm\mu_{n-1}+\mathbf{b}_n\)\({}=\mathbf{f}(\bm\mu_{n-1})\) and covariance \(\mathbf{F}_n\mathbf{V}_{n-1}\mathbf{F}_n^\T\)\({}+\bm\Gamma\). Likewise expand \(\mathbf{h}(\mathbf{z}_n)\)\({}\approx\mathbf{h}(\mathbf{m}_n)\)\({}+\mathbf{H}_n(\mathbf{z}_n-\mathbf{m}_n)\) about the predicted mean, the point where the predictive mass is concentrated; the emission becomes \(\N(\mathbf{x}_n\given\mathbf{H}_n\mathbf{z}_n+\mathbf{d}_n,\bm\Sigma)\) with \(\mathbf{d}_n\)\({}=\mathbf{h}(\mathbf{m}_n)\)\({}-\mathbf{H}_n\mathbf{m}_n\), and Result 13.30 applied to \(\mathbf{x}_n-\mathbf{d}_n\) gives the update, since \(\mathbf{x}_n-\mathbf{d}_n\)\({}-\mathbf{H}_n\mathbf{m}_n\)\({}=\mathbf{x}_n-\mathbf{h}(\mathbf{m}_n)\). The result is exact when \(\mathbf{f},\mathbf{h}\) are affine and accurate when they are nearly linear over the spread of the predictive distribution; it can fail badly (and is overconfident) for strongly nonlinear models, which is where particle filters help. ∎

(13.117). By Bayes' theorem with \(\mathbf{X}_{n-1}\) as background, \(p(\mathbf{z}_n\given\mathbf{X}_n)\)\({}=p(\mathbf{z}_n\given\mathbf{x}_n,\mathbf{X}_{n-1})\)\({}=\frac{p(\mathbf{x}_n\given\mathbf{z}_n,\mathbf{X}_{n-1})\,p(\mathbf{z}_n\given\mathbf{X}_{n-1})}{p(\mathbf{x}_n\given\mathbf{X}_{n-1})}\), and \(p(\mathbf{x}_n\given\mathbf{z}_n,\mathbf{X}_{n-1})\)\({}=p(\mathbf{x}_n\given\mathbf{z}_n)\) by (13.121). The denominator is the integral of the numerator over \(\mathbf{z}_n\). Multiplying by \(f\) and integrating gives the ratio in (13.117). This is a ratio of two expectations under the proposal \(q=p(\mathbf{z}_n\given\mathbf{X}_{n-1})\), \(\E_q[f\,g]/\E_q[g]\) with \(g(\mathbf{z})=p(\mathbf{x}_n\given\mathbf{z})\). Replacing each expectation by the average over the \(L\) samples from \(q\), the factors \(1/L\) cancel and give \(\sum_lf(\mathbf{z}^{(l)})g(\mathbf{z}^{(l)})/\sum_mg(\mathbf{z}^{(m)})\), i.e. (13.117)–(13.118). This is self-normalised importance sampling (Result 11.10, (11.22)–(11.23)) with target \(p(\mathbf{z}_n\given\mathbf{X}_n)\), proposal \(p(\mathbf{z}_n\given\mathbf{X}_{n-1})\) and importance ratio proportional to the likelihood \(p(\mathbf{x}_n\given\mathbf{z}_n)\). By the strong law of large numbers numerator and denominator averages converge to their expectations, and the ratio converges because the denominator's limit \(p(\mathbf{x}_n\given\mathbf{X}_{n-1})\) is positive.

(13.119). By the sum rule and (13.120), \(p(\mathbf{z}_{n+1}\given\mathbf{X}_n)\)\({}=\int p(\mathbf{z}_{n+1}\given\mathbf{z}_n,\mathbf{X}_n)p(\mathbf{z}_n\given\mathbf{X}_n)\d\mathbf{z}_n\)\({}=\int p(\mathbf{z}_{n+1}\given\mathbf{z}_n)p(\mathbf{z}_n\given\mathbf{X}_n)\d\mathbf{z}_n\). Substituting the Bayes expression for \(p(\mathbf{z}_n\given\mathbf{X}_n)\) gives the middle of (13.119), and this is (13.117) with \(f(\mathbf{z}_n)\)\({}=p(\mathbf{z}_{n+1}\given\mathbf{z}_n)\) (for fixed \(\mathbf{z}_{n+1}\)), which gives the mixture. Picking \(l\) with probability \(w^{(l)}\) and then drawing from \(p(\mathbf{z}_{n+1}\given\mathbf{z}_n^{(l)})\) samples this mixture exactly (ancestral sampling of a mixture). ∎

Target recursion. From (13.6) for \(n\) steps and the definition of conditional probability,

Proposal for the trajectory. If the old trajectories have (unnormalised) importance weights \(w_{n-1}\)\({}\propto p(\mathbf{z}_{1:n-1}\given\mathbf{X}_{n-1})/q_{n-1}(\mathbf{z}_{1:n-1})\) for their own proposal \(q_{n-1}\), the extended ones are drawn from \(q_n(\mathbf{z}_{1:n})\)\({}=q_{n-1}(\mathbf{z}_{1:n-1})\,q(\mathbf{z}_n\given\mathbf{z}_{n-1},\mathbf{x}_n)\). Their importance weight is target over proposal:

which is the stated recursion, the last factor being a constant absorbed by the normalisation. Resampling replaces a weighted sample by an unweighted one with (approximately) the same distribution (Result 11.13), so after it every \(w_{n-1}^{(l)}=1/L\).

Optimal proposal. With \(q\)\({}=p(\mathbf{z}_n\given\mathbf{z}_{n-1},\mathbf{x}_n)\)\({}=p(\mathbf{x}_n\given\mathbf{z}_n)p(\mathbf{z}_n\given\mathbf{z}_{n-1})/p(\mathbf{x}_n\given\mathbf{z}_{n-1})\) (Bayes, using (13.121)-type independence), the ratio in the recursion is \(p(\mathbf{x}_n\given\mathbf{z}_{n-1})\), a function of the past only. Given the past the new weight is then deterministic, so its conditional variance is zero, the smallest possible. It is rarely available in closed form, which is why the bootstrap proposal is the default. ∎