# re-estimation of emission probabilities in HMM

I am confused about the re-erstimation procedure for emissions in HMMs with Baum-Welch (still). I posted two questions concerning this general topic already and I thought I had cleared up my confusion, but not so. At least not in its entirety. The remaining issue I have is this:

Say I have a HMM with some number of states and a set of 4 emissions. Assume for state $s_k$ I have the following distribution of emissions:

e1    e2    e3    e4
.2    .1    .3    .4

Now suppose the next training sequence is e3 e1

Then the re-estimations for the emissions of state $s_k$ are the following:

$$b_{s_k}(e_1) = \frac{\gamma_1(k)}{\gamma_1(k) + \gamma_2(k)}$$

Let's just say this was equal to .3

$$b_{s_k}(e_2) = \frac{0}{\gamma_1(k) + \gamma_2(k)} = 0$$

$$b_{s_k}(e_3) = \frac{\gamma_2(k)}{\gamma_1(k) + \gamma_2(k)}$$

Let's assume this to be .7 then.

$$b_{s_k}(e_4) = \frac{0}{\gamma_1(k) + \gamma_2(k)} = 0$$

Ok so the re-estimated distribution looks like this then:

e1    e2    e3    e4
.3    .0    .7    .0

After doing this re-estimation for all states with the current sequence I assign the newly estimated ditributions to my HMM by swapping them out for the old distributions and go to the next training sequence.

So, for emissions not observed in the current training sequence I have just wiped out all emission probabilities I had estimated so far. If the next sequence would be for example e2 e4 I would end up wipig the estimations for e1 and e3 again. This cannot be correct, so I must (still) be missing something here. What is that though?

# Edit 1

The way I now went about this is to just replace the non-zero re-estimates and keep the rest and finally re-normalise all emissions. No idea if that actually does the job or messes with the guanrantee of convergence, but I see no other way.

# Edit 2

It does seem to mess with with the rest of the algorithm, because this leads to all transition probabilities per state to converge to 0 with one transition converging to 1.

Why has nobody an answer to this? Is the example confusing? Is my general problem unclear?

# Edit 3

Ok, so just to make sure I got this right: I first build expected values like so: I go over the data set and for each sequence I compute $\sum_{t=1}^{|y|} \xi_t(i, j)$ and accumulate these values for instance in a matirx cell $\xi\text{-acc}[i, j]$

So I end up with

$$\forall i, j : \text{expect-trans}[i, j] \leftarrow \sum_{y \in Y} \sum_{t=1}^{|y|} \xi_t(i, j)$$

These are all the expected values for all transitions $i \to j$.

Then I do the same for the gamma values:

$$\forall i : \text{expect-state}[i] \leftarrow \sum_{y \in Y} \sum_{t=1}^{|y|} \gamma_t(i)$$

These are all the expected values for visiting state $i$.

Then I also still need to do this for the emissions

$$\forall i, e : \text{expect-emiss}[i, e] \leftarrow \sum_{y \in Y} \sum_{t=1}^{|y|} \gamma_t(i) \cdot \chi_{[y_t = e]}$$

These are all teh expected counts for emitting $e$ from $i$.

Once I have all of these expected values I can go over to the maximization step: \begin{align*} &\forall i : \text{HMM.init-prob(i)} & \leftarrow &&\frac{\text{expect-state}[i]}{\sum_{j=1}^N \text{expect-state}[j]}\\ \\ &\forall i, j : \text{HMM.trans-prob(i, j)} & \leftarrow &&\frac{\text{expect-trans}[i, j]}{\text{expect-state}[i]}\\ \\ &\forall i, e : \text{HMM.emiss-prob(i, e)} & \leftarrow &&\frac{\text{expect-emiss}[i, e]}{\sum_{j=1}^N \text{expect-state}[j]}\\ \end{align*}

Have I got this right now? :/

Or is it rather this:

\begin{align*} &\forall y \in Y, i, j, t : \text{expect-trans}[i, j, t] &\leftarrow && \xi_t(i, j) + \text{expect-trans}[i, j, t]\\ \\ &\forall y \in Y i, t : \text{expect-state}[i, t] &\leftarrow && \gamma_t(i) + \text{expect-state}[i, t]\\ \\ &\forall y \in Y, i, e, t : \text{expect-emiss}[i, e, t] &\leftarrow && \gamma_t(i) \cdot \chi_{[y_t = e]} + \text{expect-emiss}[i, e, t] \end{align*}

with the max $t$ depending on the length of the current $y$.

\begin{align*} &\forall i : \text{HMM.init-prob(i)} & \leftarrow &&\frac{\text{expect-state}[i, 1]}{\sum_{j=1}^N \text{expect-state}[j, 1]}\\ \\ &\forall i, j : \text{HMM.trans-prob(i, j)} & \leftarrow &&\frac{\sum_{t=1}^T\text{expect-trans}[i, j, t]}{\sum_{t=1}^T\text{expect-state}[i, t]}\\ \\ &\forall i, e : \text{HMM.emiss-prob(i, e)} & \leftarrow &&\frac{\sum_{t=1}^T \text{expect-emiss}[i, e, t]}{\sum_{t=1}^T \sum_{j=1}^N \text{expect-state}[j, t]}\\ \end{align*}

with $T$ the max index of the rows $\text{expect-trans}[i, j, t]$ and $\text{expect-emiss}[i, e]$ respectively.

You write:

After doing this re-estimation for all states with the current sequence I assign the newly estimated ditributions to my HMM by swapping them out for the old distributions and go to the next training sequence.

in standard Baum-Welch, or any "vanilla" EM algorithm, you gather expected counts from all the training sequences using one set of parameters $\mathbf{b}^{(i)}$ for the current iteration $i$, and only re-estimate the next iteration's $\mathbf{b}^{(i+1)}$ after you have expected counts from the entire dataset.

There are incremental versions of EM that update parameters after examining only one data item or a "minibatch" of several data items. One incremental EM algorithm due to Neal and Hinton (1999) starts with just an initial guess at the expected counts for data item, and then sweeps through the dataset, potentially many times, to optimize the expected counts for each item.

For HMMs, you would start with an initial guess for each $\gamma^{(0)}_{jl}(s_k)$ for each emission type $e_l$ in each sequence $j$, and compute your initial transition $\mathbf{a}^{(0)}$ and emission $\mathbf{b}^{(0)}$ probabilities using those counts. Next, you would visit each sequence, and estimate your new counts $\gamma^{(i+1)}_{jl}(s_k)$ using your current parameters $\mathbf{a}^{(i)}$ and $\mathbf{b}^{(i)}$, and re-estimate the next value for $\mathbf{b}^{(i+1)}$ by replacing your previous expected counts (for the current sequence only!) with the new expected counts: \begin{align} \gamma^{(i+1)}_{jl}( s_k ) & = \gamma^{(i)}_{jl}( s_k ) & \mbox{for }j\neq\mbox{ the current sequence} \\ \gamma^{(i+1)}_{jl}( s_k ) & = \mbox{E}_{\mathbf{b}^{(i)},\mathbf{a}^{(i)}}\left[ \mbox{count}\left( s_k \rightarrow e_l, j \right) \right] & \mbox{for }j = \mbox{ the current sequence} \\ b^{(i+1)}_{sk}( e_l ) & = \frac{ \sum_{j=1}^n \gamma^{(i+1)}_{jl} }{ \sum_{l'=1}^m\sum_{j=1}^n \gamma^{(i+1)}_{jl'} } \end{align} where $m$ is the number of emission types and $n$ is the number of training sequences. Note that you'll also need to re-estimate transition probabilities incrementally in exactly the same fashion.

Neal and Hinton (1999) point out that you can cheaply update the numerator and denominator by storing the sums in summary variables $\tilde{\gamma}^{(i)}_{l}$ and then subtracting the old counts and adding the new counts. For example, here is the numerator: \begin{align} \tilde{\gamma}^{(0)}_{l} & = \sum_{j=1}^n \gamma^{(0)}_{jl} \\ \tilde{\gamma}^{(i+1)}_{l} & = \tilde{\gamma}^{(i)}_l - \gamma^{(i)}_{jl} + \gamma^{(i+1)}_{jl} \end{align}

• I saw in another question of yours that you believe that transition probabilities are updated after every training sequence. This is not true for vanilla EM algorithms like Baum-Welch. For vanilla EM, all model parameters are updated only after gathering expected counts from the entire dataset. For HMMs, you gather expected counts for emissions and transition probabilities from the entire dataset, and only then do you re-estimate $\mathbf{a}$ and $\mathbf{b}$. If you want to update earlier, then you'll need to do something more sophisticated like the incremental EM algorithm I described. – jkpate Apr 23 '17 at 15:42
• Ah ok, that never got stated explictely in all sources I considered :/ However, what do I do with the expected counts that I already have? Do I overwrite them or do I store all of them and in the end compute their respective average? – lo tolmencre Apr 23 '17 at 16:24
• Just sum the expected counts. Your total expected count for $s_k$ emitting $e_l$ is the expected count for sequence one plus the expected count for sequence 2 and so on, all the way through the dataset. – jkpate Apr 23 '17 at 18:03
• Yes, you're right. The example on Wikipedia looks like it's designed with only one sequence in mind. – jkpate Apr 23 '17 at 19:38
• The first one is almost right, except for the emission probability. The denominator should be only the expected count of state $i$: $\forall i, e : \text{HMM.emiss-prob(i, e)} \leftarrow \frac{\text{expect-emiss}[i, e]} {\text{expect-state}[i]}$ – jkpate Apr 24 '17 at 16:12