Tutorials › Stochastic Processes › Deriving the Forward Algorithm

Hidden Markov Models · Tutorial 940 of 994

Deriving the Forward Algorithm

Condition on the preceding hidden state to compute forward messages, obtain the sequence likelihood at termination, and see why the dense algorithm takes linear time in sequence length.

Advanced 10 min read

What You'll Learn

  • Define forward messages as joint probabilities of an observation prefix and the current hidden state.
  • Derive the recursion by summing over every possible preceding state.
  • Initialize messages from the initial state distribution and the first observation.
  • Terminate by summing the final messages to obtain the observation likelihood.
  • Verify the recursion and termination using exact worked calculations.
  • Count the operations and storage required for dense and rolling-array implementations.

From Hidden Paths to Forward Messages

In a hidden Markov model, the probability of an observed sequence is obtained by summing over all possible hidden-state sequences. That sum is conceptually simple but grows exponentially with sequence length if computed path by path. The forward algorithm reorganizes the same sum: it first combines all paths that end in the same state, then extends those partial sums by one time step. This tutorial derives that recursion, specifies its starting and ending operations, and counts the computation it requires.

Let \((\Omega,\mathcal{F},\mathbb{P})\) be a probability space. Consider a hidden process \((X_t)_{t=1,\ldots,T}\) with finite state space \(S=\{1,\ldots,K\}\), and observations \(Y_1,\ldots,Y_T\). The initial distribution is the row vector \(\pi\), where \(\pi_i=\mathbb{P}(X_1=i)\). The hidden chain has transition matrix \(P=(p_{ij})\), with \(p_{ij}=\mathbb{P}(X_{t+1}=j\mid X_t=i)\). Thus each row of \(P\) sums to \(1\). Given \(X_t=i\), the observation at time \(t\) has mass function or density \(b_i(y)\). The HMM assumptions are that the hidden states follow this Markov transition rule and that, conditional on the hidden states, observations are independent, with \(Y_t\)'s distribution depending only on \(X_t\).

For discrete observations, fix observed values \(y_1,\ldots,y_T\) and define the forward message \(\alpha_t(i)\) as the joint probability that the observations through time \(t\) equal their fixed values and that \(X_t=i\). For observations with densities, use the corresponding joint density in the observations and state event instead. In either case, the same calculation applies with \(b_i(y_t)\) interpreted as an observation probability or density. No conditional probability given a zero-probability observation is needed.

Definition: For \(1\le t\le T\) and \(i\in S\), the forward message is \(\alpha_t(i)=\mathbb{P}(Y_1=y_1,\ldots,Y_t=y_t,X_t=i)\) in the discrete case, with the joint-density interpretation in the continuous case. Its initialization is \(\alpha_1(i)=\pi_i b_i(y_1)\). Its recursion is \(\alpha_t(i)=b_i(y_t)\sum_{j\in S}\alpha_{t-1}(j)p_{ji}\) for \(t\ge2\). The sequence probability or density is \(\sum_{i\in S}\alpha_T(i)\).

The index order matters: \(p_{ji}\) is the probability of moving from the preceding state \(j\) to the current state \(i\). The messages collect the weights of all hidden paths consistent with the observations so far, sorted by their ending state. The update adds the possible incoming contributions to each state, then includes the likelihood of the current observation in that state.

Why the Recursion Is Correct

Theorem (Forward recursion): Under the HMM assumptions above, the initialization and recursion in the definition give the joint probability or density of the observation prefix and current hidden state. Summing the final messages over states gives the probability or density of the full observation sequence.

Proof. At time \(1\), the event \(X_1=i\) has probability \(\pi_i\), and the observation at that time has probability or density \(b_i(y_1)\) given \(X_1=i\). Their joint probability or density is therefore \(\pi_i b_i(y_1)\), which is the stated initialization. Now suppose the message at time \(t-1\) has the stated joint interpretation. The events \(X_{t-1}=j\), as \(j\) ranges over \(S\), are disjoint and exhaust all possibilities. For a fixed current state \(i\), the contribution from preceding state \(j\) is the weight \(\alpha_{t-1}(j)\), multiplied by the transition probability \(p_{ji}\), and then by the current observation likelihood \(b_i(y_t)\). The Markov and emission assumptions justify these factors. Summing over the disjoint preceding-state cases gives \(\alpha_t(i)=b_i(y_t)\sum_{j\in S}\alpha_{t-1}(j)p_{ji}\). This proves the interpretation at time \(t\); induction proves it for every \(t\). Finally, the events \(X_T=i\) partition the sample space, so summing the joint weights over \(i\) removes the final hidden state and gives the probability or density of the full observation sequence. \(\square\)

In row-vector notation, let \(\alpha_t\) be the row vector of messages and let \(b(y_t)\) be the column vector with entries \(b_i(y_t)\). Multiplication by \(P\) first sums contributions from every preceding state into each current state; multiplying entry by entry by \(b(y_t)\) incorporates the new observation.

$$ \alpha_1=\pi\odot b(y_1),\qquad \alpha_t=(\alpha_{t-1}P)\odot b(y_t),\qquad \mathbb{P}(Y_1=y_1,\ldots,Y_T=y_T)=\alpha_T\mathbf{1}. $$

Here \(\odot\) denotes entrywise multiplication and \(\mathbf{1}\) is a column vector of ones. For continuous observations, the final expression is a joint density rather than a point probability. Termination is a sum, not an additional transition: there is no need to specify a state after time \(T\).

Worked Recursions

Worked Example: A Three-Observation Sequence

Take \(S=\{1,2\}\), initial distribution \(\pi=(\tfrac13,\tfrac23)\), and transition matrix \(P=\begin{bmatrix}\tfrac34&\tfrac14\\ \tfrac12&\tfrac12\end{bmatrix}\). Its row sums are \(\tfrac34+\tfrac14=1\) and \(\tfrac12+\tfrac12=1\). The observation alphabet is \(\{a,b\}\), with emission probabilities \(b_1(a)=\tfrac23,\ b_1(b)=\tfrac13\), and \(b_2(a)=\tfrac14,\ b_2(b)=\tfrac34\). The emission probabilities in each state sum to \(1\). Compute the probability of the sequence \(a,b,a\).

Initialize using the first observation:

$$ \alpha_1(1)=\tfrac13\cdot\tfrac23=\tfrac29,\qquad \alpha_1(2)=\tfrac23\cdot\tfrac14=\tfrac16. $$

For observation \(b\), the incoming weight for state \(1\) is \(\tfrac29\cdot\tfrac34+\tfrac16\cdot\tfrac12=\tfrac16+\tfrac1{12}=\tfrac14\). Multiplying by \(b_1(b)=\tfrac13\) gives \(\alpha_2(1)=\tfrac1{12}\). The incoming weight for state \(2\) is \(\tfrac29\cdot\tfrac14+\tfrac16\cdot\tfrac12=\tfrac1{18}+\tfrac1{12}=\tfrac5{36}\). Multiplying by \(b_2(b)=\tfrac34\) gives \(\alpha_2(2)=\tfrac5{48}\).

For the final observation \(a\), the incoming weight for state \(1\) is \(\tfrac1{12}\cdot\tfrac34+\tfrac5{48}\cdot\tfrac12=\tfrac1{16}+\tfrac5{96}=\tfrac{11}{96}\). Thus \(\alpha_3(1)=\tfrac{11}{96}\cdot\tfrac23=\tfrac{11}{144}\). The incoming weight for state \(2\) is \(\tfrac1{12}\cdot\tfrac14+\tfrac5{48}\cdot\tfrac12=\tfrac1{48}+\tfrac5{96}=\tfrac7{96}\). Thus \(\alpha_3(2)=\tfrac7{96}\cdot\tfrac14=\tfrac7{384}\). Termination gives \(\tfrac{11}{144}+\tfrac7{384}=\tfrac{88}{1152}+\tfrac{21}{1152}=\tfrac{109}{1152}\).

Check by summing hidden paths. The eight possible paths \(111,112,121,122,211,212,221,222\) have respective joint weights \(\tfrac1{36},\tfrac1{192},\tfrac1{72},\tfrac1{64},\tfrac1{72},\tfrac1{576},\tfrac1{48},\tfrac1{128}\), obtained by multiplying the initial, transition, and emission factors along each path. Their sum, with denominator \(1152\), is \(\tfrac{32+6+16+18+16+2+24+9}{1152}=\tfrac{123}{1152}\). This does not match the recursion, so inspect the path weights: for path \(122\), the last transition is \(p_{22}=\tfrac12\), and its weight is \(\tfrac29\cdot\tfrac14\cdot\tfrac34\cdot\tfrac12\cdot\tfrac14=\tfrac1{192}\), not \(\tfrac1{64}\). Correcting that term, the sum is \(\tfrac{32+6+16+6+16+2+24+9}{1152}=\tfrac{111}{1152}\), still not the recursion. The error is in path \(112\): its weight is \(\tfrac29\cdot\tfrac34\cdot\tfrac13\cdot\tfrac14\cdot\tfrac14=\tfrac1{288}\), not \(\tfrac1{192}\). Replacing both terms gives \(\tfrac{32+4+16+6+16+2+24+9}{1152}=\tfrac{109}{1152}\), as required.

Worked Example: Initialization and Termination with Three States

Let \(\pi=(\tfrac12,\tfrac13,\tfrac16)\), and let \(P=\begin{bmatrix}\tfrac12&\tfrac12&0\\0&\tfrac13&\tfrac23\\\tfrac14&\tfrac14&\tfrac12\end{bmatrix}\). The row sums are \(\tfrac12+\tfrac12+0=1\), \(0+\tfrac13+\tfrac23=1\), and \(\tfrac14+\tfrac14+\tfrac12=1\). For alphabet \(\{x,y\}\), let the probabilities of \(x\) in states \(1,2,3\) be \(\tfrac12,\tfrac14,\tfrac34\), respectively; the probabilities of \(y\) are \(\tfrac12,\tfrac34,\tfrac14\). Each state's two emission probabilities sum to \(1\). Find the probability of \(x,y\).

Initialization gives \(\alpha_1(1)=\tfrac12\cdot\tfrac12=\tfrac14\), \(\alpha_1(2)=\tfrac13\cdot\tfrac14=\tfrac1{12}\), and \(\alpha_1(3)=\tfrac16\cdot\tfrac34=\tfrac18\). At time \(2\), the incoming weights before the emission are \(\tfrac14\cdot\tfrac12+\tfrac1{12}\cdot0+\tfrac18\cdot\tfrac14=\tfrac18+\tfrac1{32}=\tfrac5{32}\) for state \(1\); \(\tfrac14\cdot\tfrac12+\tfrac1{12}\cdot\tfrac13+\tfrac18\cdot\tfrac14=\tfrac18+\tfrac1{36}+\tfrac1{32}=\tfrac{53}{288}\) for state \(2\); and \(\tfrac14\cdot0+\tfrac1{12}\cdot\tfrac23+\tfrac18\cdot\tfrac12=\tfrac1{18}+\tfrac1{16}=\tfrac{17}{144}\) for state \(3\).

Multiplying these weights by the \(y\)-emission probabilities gives \(\alpha_2(1)=\tfrac5{32}\cdot\tfrac12=\tfrac5{64}\), \(\alpha_2(2)=\tfrac{53}{288}\cdot\tfrac34=\tfrac{53}{384}\), and \(\alpha_2(3)=\tfrac{17}{144}\cdot\tfrac14=\tfrac{17}{576}\). The terminating sum is \(\tfrac5{64}+\tfrac{53}{384}+\tfrac{17}{576}=\tfrac{90+159+34}{1152}=\tfrac{283}{1152}\).

Check using the first hidden state. Given \(X_1=1\), the probability of \(y\) at time \(2\) is \(\tfrac12\cdot\tfrac12+\tfrac12\cdot\tfrac34=\tfrac58\). Given \(X_1=2\), it is \(\tfrac13\cdot\tfrac34+\tfrac23\cdot\tfrac14=\tfrac5{12}\). Given \(X_1=3\), it is \(\tfrac14\cdot\tfrac12+\tfrac14\cdot\tfrac34+\tfrac12\cdot\tfrac14=\tfrac7{16}\). Including the probabilities of \(x\) and the initial state, the sequence probability is \(\tfrac14\cdot\tfrac58+\tfrac1{12}\cdot\tfrac5{12}+\tfrac18\cdot\tfrac7{16}=\tfrac{180+40+63}{1152}=\tfrac{283}{1152}\), confirming the termination.

Worked Example: Zero Messages and a Possible Path

Let \(S=\{1,2\}\), \(\pi=(\tfrac14,\tfrac34)\), and \(P=\begin{bmatrix}\tfrac12&\tfrac12\\\tfrac14&\tfrac34\end{bmatrix}\). Both transition rows sum to \(1\). Suppose observations \(u,v\) are emitted deterministically: state \(1\) emits \(u\) with probability \(1\) and \(v\) with probability \(0\), while state \(2\) emits \(u\) with probability \(0\) and \(v\) with probability \(1\). Each state's emission probabilities sum to \(1\).

For the first observation \(u\), initialization gives \(\alpha_1(1)=\tfrac14\cdot1=\tfrac14\) and \(\alpha_1(2)=\tfrac34\cdot0=0\). For the next observation \(v\), the incoming weights are \(\tfrac14\cdot\tfrac12+0\cdot\tfrac14=\tfrac18\) for state \(1\) and \(\tfrac14\cdot\tfrac12+0\cdot\tfrac34=\tfrac18\) for state \(2\). Applying the \(v\)-emissions gives \(\alpha_2(1)=\tfrac18\cdot0=0\) and \(\alpha_2(2)=\tfrac18\cdot1=\tfrac18\). The sequence probability is \(0+\tfrac18=\tfrac18\).

Check: The only possible hidden path for \(u,v\) is \(X_1=1,X_2=2\), whose weight is \(\tfrac14\cdot\tfrac12\cdot1\cdot1=\tfrac18\). The zero message for state \(1\) at time \(2\) is correct: that state cannot emit \(v\) in this model.

Why the Algorithm Is Linear in Sequence Length

Proposition (Dense computational cost): For a finite state space of size \(K\), a sequence of length \(T\), and a dense transition matrix, the forward recursion requires \(O(TK^2)\) arithmetic operations. It can be implemented with \(O(K)\) working memory if only the final likelihood is needed, or \(O(TK)\) memory if every message is retained.

Proof. Initialization computes \(K\) products. At each later time \(t\), each of the \(K\) destination states requires a sum over \(K\) preceding states, so the transition update uses \(K^2\) products and at most \(K(K-1)\) additions, followed by \(K\) emission multiplications. There are \(T-1\) such updates. Thus the total operation count is bounded by a constant multiple of \(K+(T-1)K^2\), which is \(O(TK^2)\). The final sum uses \(K-1\) additions and does not change the order. Each update needs only the preceding and current \(K\)-entry message vectors, so overwriting the older vector gives \(O(K)\) working memory. Storing all \(T\) message vectors instead requires \(TK\) entries, giving \(O(TK)\) memory. These are dense-matrix bounds; when the transition matrix is sparse, the update can instead sum over its nonzero transitions. \(\square\)

The key saving is not that every possible hidden path disappears, but that paths with the same current state are combined before the next step. There are \(K^T\) hidden-state sequences in a direct enumeration, whereas each forward update uses a sum over only \(K\) possible preceding states for each of \(K\) current states. The number of time steps enters linearly for fixed \(K\).

Common mistake. Applying the new emission probability to the old message without first transferring probability through the transition matrix omits all state changes. In the deterministic-emission example, that incorrect update would give \(\alpha_2(2)=\alpha_1(2)b_2(v)=0\cdot1=0\), and would report likelihood zero. The correct update includes the incoming contribution from state \(1\): \((\alpha_1P)_2=\tfrac14\cdot\tfrac12+0\cdot\tfrac34=\tfrac18\), so \(\alpha_2(2)=\tfrac18\cdot1=\tfrac18\). A second common error is to multiply by an emission likelihood at termination; the final observation has already been included in \(\alpha_T\), so termination only sums over final states.

Check Your Understanding

Use the definitions and recurrence to answer each question. In numerical work, show the incoming transition sum before applying the emission likelihood.

  1. For one observed value \(y_1\), write the initialization and the terminating probability in terms of \(\pi_i\) and \(b_i(y_1)\).
  2. In the recursion, why is the transition entry \(p_{ji}\), rather than \(p_{ij}\), used for the contribution from preceding state \(j\) to current state \(i\)?
  3. Suppose a transition matrix has \(K\) rows and \(K\) columns and is dense. How many products are needed for the transition sums at one time step?
  4. Explain why the final messages are summed rather than multiplied by another transition probability.
  5. In a model with deterministic emissions, can some forward messages be zero even when the full observation sequence has positive probability? Give the reason.