The Viterbi Algorithm

The Viterbi algorithm is an efficient method for determining the most probable set of states (a single ‘track’ of steps on the time-frequency plane) in a Markov model dependent on data, where the model has a discrete number of states at each step. Rather than computing the probability of every possible track and selecting the most probable, the algorithm maximises this probability after every discrete step. As a result, a partial track which cannot ultimately be the most probable is rejected before the next step is calculated, and only a fraction of all possible tracks need to be computed to find the one that is most probable.

In this work we apply the Viterbi algorithm to a GW strain time-series to find the most probable track of a single variable-frequency signal in the noisy data. We divide the time series into $N$ equal-length and contiguous segments ${\mathbf x}_j$, defining the set $D \equiv \{{\mathbf x}_j\}$. The ‘states’ in the model correspond to the frequencies a signal could have in each segment. A ‘track’ is a list of such frequencies ${\mathbf \nu}\equiv \{\nu_j\}$, where $\nu_j$ is the frequency in the segment ${\mathbf x_j}$.

Our objective is to calculate the most probable track given the data, i.e., the track that maximises $p({\mathbf \nu}\mid D)$. Using Bayes theorem, this posterior probability can be written as

\begin{equation} p({\mathbf \nu} \mid D) = \frac{p({\mathbf \nu})p(D \mid {\mathbf \nu})}{p(D)}, \end{equation}

where $p({\mathbf `:nbsphinx-math:nu`}) $ is the prior probability of the track, $p(D \mid{\mathbf \nu})$ is the likelihood of the track (i.e., the probability of the data given the track) and $p(D)$ is the model evidence (or marginalised likelihood).

The Viterbi algorithm treats the track as the result of a Markovian process, such that the current state depends only on the previous state. It is therefore useful to split the track’s prior into a set of transition probabilities such that

\begin{align} p({\mathbf \nu}) &= p(\nu_{N - 1}, \ldots, \nu_1, \nu_0)\nonumber \\ &= p(\omega_{N} \mid \omega_{N-1}, \ldots, \omega_1,\omega_0)p(\omega_{N-1} \mid \omega_{N-2}, \ldots, \omega_1,\omega_0) \ldots p(\omega_1 \mid \omega_0)p(\omega_0) \\ &= p(\nu_{N - 1} \mid \nu_{N-2})p(\nu_{N-2} \mid \nu_{N-3}) \dots p(\nu_1 \mid \nu_0)p(\nu_0) \nonumber \\ &= p(\nu_0)\prod_{j=1}^{N-1}p(\nu_{j} \mid \nu_{j-1}), \end{align}

where $p(\nu_0)$ is the prior probability that the signal in the first time step has a frequency $\nu_0$ and $p(\nu_{j} \mid \nu_{j-1})$ is the prior ‘transition’ probability for $\nu_j$ given the frequency at the last step was $\nu_{j-1}$.

The noise in each of the segments can be treated as independent, therefore we can write,

\begin{equation} p(D \mid {\mathbf \nu}) = \prod_{j=0}^{N-1}p({\mathbf x_j} \mid \nu_j), \end{equation}

where $p({\mathbf x_j} \mid \nu_j)$ is the likelihood of our signal having a frequency $\nu_j$ in the $j$th segment.

Using the above equations, the posterior probability is then

\begin{equation} p({\mathbf \nu} | D) = \frac{p(\nu_0)p({ \mathbf x_0} | \nu_0) \displaystyle\prod_{j=1}^{N-1}p(\nu_{j} | \nu_{j-1})p({\mathbf x_j} | \nu_j)}{\displaystyle\sum_{S} \left\{p(\nu_0)p({\mathbf x_0} | \nu_0)\displaystyle\prod_{j=1}^{N-1}p(\nu_{j} | \nu_{j-1})p({\mathbf x_j} | \nu_j)\right\}} , \end{equation}

where in the denominator we must sum over all possible tracks $S$. We require the specific track, or set of frequencies, $\hat{\mathbf \nu}$ that maximises the posterior probability, and which therefore maximises the numerator on the right-hand side of the previous equation, i.e.,

\begin{equation} p(\hat{\mathbf \nu} | D) \propto \max_{\mathbf \nu}{\left[p(\nu_0)p({\mathbf x_0} | \nu_0) \prod_{j=1}^{N-1}p(\nu_{j} |\nu_{j-1})p({\mathbf x_j} | \nu_j)\right]}. \end{equation}

This track also maximises the log of the probability,

\begin{equation} \begin{split} \log p(\hat{\mathbf \nu} | D) = \max_{{\mathbf \nu}}{\biggl\{ \log p(\nu_0) + \log p({\mathbf x_0} | \nu_0) } \\ \left. \sum_{j=1}^{N-1} \biggl[ \log p(\nu_{j} | \nu_{j-1}) + \log p({\mathbf x_j} | \nu_j) \biggr] \right\} + \text{const}. \end{split} \end{equation}

The Viterbi algorithm finds the most probable track $\hat{\mathbf \nu}$ by calculating the quantities in the previous equation for each frequency at each time step. In the following sections we explain how this is achieved in practice.

[ ]: