<!--
Abstract:

CS 461
Introduction to Deep Learning
Lecture 11
-->

# CS 461 - Lecture 13

## Machine Learning Principles

### Applications of HMMs

Bernhard Firner

2025-10-20

---

## Adaboost Regression

* Finally made that decision boundary explanation
  * See recitation for explanation!

---

## Review

* Introduced a new feature to Markov models: hidden states
* Our observations are `emitted` by the hidden states
* The hidden states may `transition` without us knowing

---

## Variables

* Our new models have three variables:
  * $\pi$: the initial state distribution
  * A: the transition matrix
  * B: the emission matrix

---

## EM Algorithm

* Training a HMM is done with the EM algorithm
  * Expectation-Maximization
  * Same as with mixture models

---

## Expectation

* Probabilities are all based upon the hidden states
* From $1:t$:
  * What are the hidden state probabilities given the past observations?
* From $t+1:T$
  * What are the probabilities we were in those hidden states given future observations?
* Assume that our initial state distribution is even:
  * $\pi = [0.5, 0.5]$

---

## Forward

* The expectation is calculated forwards and backwards
  * Forwards:
    * Use the probability of transitioning into each hidden state, $A_j$?
    * And the probability of seeing this token at the state, $B_jO[t]$
* $\alpha = p(x_{t-1})ABO[t] = \alpha[t-1]ABO[t]$
* For $t=0$, $\alpha[0] = p(x_0) \leftarrow \pi BO[0]$

---

## Backwards

* $\beta = A B O[t] p(x_{t+1}|x_{t+1:T}) = A B O[t]\beta[t+1]$
* Treat $\beta[T+1]$ as 1

---

## Dynamic Programming

* By the way
  * The forward and backward pass use dynamic programming
  * There are many pathways through the model 
    * But combine the probabilities as we go and don't need to remember unique paths
* Similar to the assembly line switching problem

---

## State Probability

* $ABO[t]$ shows up in both the forward and backward pass
  * $A$ (the transition matrix) governs how likely we are to be in each state
  * $BO[t]$ is the probability of the observation at $t$ for each hidden state

---

## Expectation Estimation

* Going to normalize values by the chain probability
  * $p(O) = \prod(\alpha)$
  * We know that the observations occurred, so that is the maximum probability we can achieve
* Expected state estimates $\gamma$
* Expected state transition estimates $\xi$

---

## Expectation continued

* $\gamma = \frac{\alpha\beta}{p(O)}$
  * Just the normalized probability of the observations given our current parameters
  * At each time, $t$, so sized $[T,k]$

---

## $\xi$ expectation

* $\xi$ will be size $[T,k,k]$
  * Transition probability matrix at each time
  * Will be built from the forward pass probability, the current transition probability, the probability of the next observation, and the backward probabilities

---

## $\xi$ continued

* $\xi[t] \leftarrow \frac{\alpha[t] A B O[t] \beta[t+1]}{p(O)}$
  * Read this as:
    * the probability of the left chain
    * the probability of seeing the observation given our model emissions
    * the probability of seeing this given the right chain
    * normalized by the entire chain probability

---

## Training and Time-Invariance

* We still assume the time-invariant Markov property
* But now the transitions are between hidden states
* So is the time window over hidden states, or visible states?
  * Could be either, or both!
  * Hidden states could emit n-tuples
* Current we haven't used a window

---

## HMM Uses

* There are generally a few ways to use a HMM:
  * The forward pass, to estimate a sequence likelihood
  * Use $\gamma$, the state expectations, to guess the hidden states of a sequence
  * Learn A, B, $\pi$ for a sequence
  * Use the HMM to generate sequences

---

## Occasionally Unfair Coin

* We trained a HMM on an occasionally unfair coin

```python
# true emission probabilities for a fair coin
fair_emission = [0.5, 0.5]
# true emission probabilities for a loaded die that rolls a 6 50% of the time
unfair_emission = [0.1, 0.9]

actual_emission_matrix = np.array([fair_emission, unfair_emission])

actual_transition_matrix = [ [0.8, 0.2],
                             [0.2, 0.8]]
```

---

## Accuracy

* We may want to know when we're being cheated
  * Just check $\gamma$, from the last expectation

---

## Coin Prediction Accuracy

---

## Coin State Estimates

* At T=5000
  * Slight mistakes, but decent model of the true system

<pre>
$ python3 12_occassionally_unfair_casino.py coin --samples 5000
Converged after 335 steps
Transition probabilities for fixed coin: [0.81646384 0.18353616]
Transition probabilities for fair coin: [0.23310806 0.76689194]
Emission probabilities for fixed coin: [0.11844357 0.88155643]
Emission probabilities for fair coin: [0.51529647 0.48470353]
Log-Likelihood P(O): -2986.7627
Total Expected Time in Fair State (1): 2201.78 / 5000 steps
Total Expected Time in Fixed State (0): 2798.22 / 5000 steps
State prediction accuracy 0.7338
</pre>

---

## Die Prediction Accuracy

---

* At T=5000
  * Better than the coin
  * The die has a more skewed transition matrix

<pre>
$ python3 12_occassionally_unfair_casino.py die --samples 5000
Converged after 151 steps
Transition probabilities for fixed die: [0.95515831 0.04484169]
Transition probabilities for fair die: [0.04374899 0.95625101]
Emission probabilities for fixed die: [0.07943047 0.09693867 0.10160139 0.11055038 0.10491876 0.50656034]
Emission probabilities for fair die: [0.16891674 0.15031408 0.17497145 0.16981986 0.16977329 0.16620459]
Log-Likelihood P(O): -8427.0733
Total Expected Time in Fair State (1): 2535.00 / 5000 steps
Total Expected Time in Fixed State (0): 2465.00 / 5000 steps
State prediction accuracy 0.8456
</pre>

---

## Easier? Harder?

* Difficulty is a combination of factors
  * The transitions
    * $[0.8, 0.2]$ in the coin
    * $[0.9, 0.1]$ in the die
    * Longer die sequences give better statistics
  * Chance of the state oscillating is $0.04$, compared to $0.01$ for the die
* But $0.5$ fixed die change is less than $0.9$ fixed coin

---

## Transition Matrix Modification

* If the coin transition matrix is $[0.99, 0.01]$

<pre>
$ python3 12_occassionally_unfair_casino.py coin --samples 1000
Converged after 100 steps
Transition probabilities for fixed coin: [0.98195615 0.01804385]
Transition probabilities for fair coin: [0.00635834 0.99364166]
Emission probabilities for fixed coin: [0.05223814 0.94776186]
Emission probabilities for fair coin: [0.48337039 0.51662961]
Log-Likelihood P(O): -614.0970
Total Expected Time in Fair State (1): 778.79 / 1000 steps
Total Expected Time in Fixed State (0): 221.21 / 1000 steps
State prediction accuracy 0.98
</pre>

---

## Lessons

* We want hidden state that are distinct
  * Either they emit things rarely found in other states
  * Or they don't transition very frequently
* There is some noise in training
  * Try with multiple random initializations

---

## Discussion

* Unlike some earlier techniques, HMMs are easily applied to both generation and classification.
  * Their structures are also adaptable to specialized applications.
* Similar to mixture models, where we chose the distributions for the experts.
  * In our penguin example, we mixed categorical and gaussian experts.

---

## Other Kinds of HMMs

* We've seen categorical MMs
  * No reason that we couldn't insert gaussian (or other distributions) for continuous tokens.
* Also no reason why we couldn't cluster groups with Markov models
  * Called Mixture Hidden Markov Models (MHMM)
  * Useful when modelling humans navigating a website

-v-

## Look at this paper

---

## HMM Applications

* There are a lot of ways to apply HMMs
* We talked about doing word correction for spam filters
  * Detect when someone is trying to outsmart a Bayesian filter
  * Why? Because our spam filter should be fast
* This is an application for `Profile HMMs`

---

## Profile HMM

* A profile HMM detects a sequence with errors
  * First used in bioinformatics
* Errors (for us) will be insertions or deletions of text

---

## Tuning the Model

* Previously we've never had to tune our models
  * Sure, we've played with the number of clusters
  * But now we're going to change the structure itself
* We are going to add our own expert knowledge into the problem

---

## Formulation

* Let's begin with a word
* Each hidden state of our HMM represents a position in the word

---

## Dealing With Obfuscation

* We could have our hidden states also emit typos with some probability

---

## Messy

* As we think of more convoluted examples, things get messy
  * For example, let's say you are trying to catch the word "violet"
  * Let's say we encounter "|/io1e+"
  * How can we even build an expectation of our hidden state?

---

## The Problem

* If we use a hidden state to represent the true letter
  * let's pick "v"
* If the probability of printing "|/" is also high, then the probability of "v" goes down
  * Soon we aren't searching for a word anymore, just a pile of letters, and we'll match almost anything

---

## The Solution

* Separate the hidden states for obfuscation
* This way the actual word still has the highest probability match

---

## How?

* This helps by separating the high probability "correct" letters from obfuscating noise
* The Obf nodes provide a low probability pathway to get from one real letter to the next
  * But an entire chain that only goes through them will have small probability

---

## Cleanup

* We need an entry point into the string
* And an exit point, to know that we're done (in the event that the last letter is deleted)

---

## Improvement

* Speaking of deleting, deletions should be rare
  * "|/io1e+" deletes 3 characters and inserts 4
  * Can't delete more than 6 letters, but could insert more
* So let's separate insertions from deleting

---

## Current View

---

## Complete Hidden State

* Great(?)
* But what about `|/` for `v`
  * 1 deletion and 2 insertions
* We'll need to allow insertions to remain in the same state

---

## Repeated insertions

---

## Now What?

* That graph describes the transition matrix, A
  * Any path with an arrow is non-zero
* For a word of length $r$, we need $3(r+1)$ hidden states

---

## Transition Matrix

<p><math display="block" xmlns="http://www.w3.org/1998/Math/MathML"><semantics><mrow><mi>A</mi><mo>=</mo><mtable><mtr><mtd columnalign="center" style="text-align: center"><mi>b</mi><mi>b</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>m</mi><mn>1</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>m</mi><mn>2</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>m</mi><mn>3</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>i</mi><mn>1</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>i</mi><mn>2</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>i</mi><mn>3</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>i</mi><mn>4</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>d</mi><mn>1</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>d</mi><mn>2</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><msub><mi>d</mi><mn>3</mn></msub></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mi>e</mi><mi>e</mi><mspace width="1.0em"></mspace></mtd></mtr></mtable><mrow><mo stretchy="true" form="prefix">(</mo><mtable><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.7</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.1</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.7</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.1</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.7</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.1</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.3</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.7</mn></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.75</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.05</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.75</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.05</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.75</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.05</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.2</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.8</mn></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.8</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.05</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.15</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.8</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.05</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.15</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.05</mn></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0.95</mn></mtd></mtr><mtr><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>0</mn><mi>.</mi></mtd><mtd columnalign="center" style="text-align: center"><mn>1</mn><mi>.</mi></mtd></mtr></mtable><mo stretchy="true" form="postfix">)</mo></mrow></mrow><annotation encoding="application/x-tex">
A= 
\begin{matrix}
bb \\
m_1 \\
m_2 \\
m_3 \\
i_1 \\
i_2 \\
i_3 \\
i_4 \\
d_1 \\
d_2 \\
d_3 \\
ee
\quad
\end{matrix}
\begin{pmatrix}
 0. &   0.7  & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.1  & 0.   & 0.   & 0.   \\
 0. &   0.   & 0.7  & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.1  & 0.   & 0.   \\
 0. &   0.   & 0.   & 0.7  & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.1  & 0.   \\
 0. &   0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.3  & 0.   & 0.   & 0.   & 0.7  \\
 0. &   0.75 & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.05 & 0.   & 0.   & 0.   \\
 0. &   0.   & 0.75 & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.05 & 0.   & 0.   \\
 0. &   0.   & 0.   & 0.75 & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.05 & 0.   \\
 0. &   0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.2  & 0.   & 0.   & 0.   & 0.8  \\
 0. &   0.   & 0.8  & 0.   & 0.   & 0.05 & 0.   & 0.   & 0.   & 0.15 & 0.   & 0.   \\
 0. &   0.   & 0.   & 0.8  & 0.   & 0.   & 0.05 & 0.   & 0.   & 0.   & 0.15 & 0.   \\
 0. &   0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.05 & 0.   & 0.   & 0.   & 0.95 \\
 0. &   0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 0.   & 1.   
\end{pmatrix}</annotation></semantics></math></p>

</div>

---

## Forward Pass

* Let's put this into use
* We can take a forward pass with this Markov Model, using it to find misspelled versions of words
* But how do we identify that something is an insertion? Or a deletion?

---

## Forward Pass Modifications

* Only the match states, $m1, m2, m3$, actually match an observable token

---

## Other Matches

* The insertion states match anything they see
* Deletions *could* have occurred if the visible token isn't what a match state expects
* $bb$ is only relevant at step 1
* The probability of $ee$ tells us the probability of a match

---

## Converting Observations

* We'll convert observations into a vector
  * 0 if the current token does not match the hidden state emission
  * 1 if it does
* This will simplify our matrix math

---

## Emissions

* The observation vector, $O$, is used with the emission matrix, B
  * But right now only works for the match states
* $1 - O$ tells us when an insertion is possible
  * And the first insertion state is always possible
    * e.g. *ccccccccat*
* Deletions happen between visible tokens, so they are always possible

---

## Word Vs T

* We may consume more characters than in our word
  * Let's call that size $S$
* We'll give the forward pass sequences of length $S > \lvert Word\rvert$

---

## $\alpha$ Initialization

```python
# Initialization (t=0)
log_alpha = np.full((T+1, 3*(r+1)), -np.inf)
# Begin with the probabilities of seeing the first observation for the emissions of each hidden state
# This will always be the bb state, regardless of any observations
# The first insertion state can occur before we match any character
# Why don't we just say that we begin in the insertion state?
# Because this is how Gordillo and Conde write this
log_alpha[0, r+1] = log_A[0,r+1]
log_alpha[0, 1] = log_A[0,1]
log_alpha[0, 2*(r+1)] = log_A[1, 2*(r+1)]
```

---

## $\alpha$ iteration

```python
for t in range(T):
    # Matching probabilities for any observed tokens
    match_probs = np.exp(log_B[match_slice]) @ obs[t]
    # The first insertion state can always insert. Other insertions only work if there was no match.
    inserts_possible = np.concatenate(([1], 1 - match_probs))
    # Deletes are always possible
    deletes_possible = np.ones((r))
    all_emission_probs = np.log(np.concatenate(([1], match_probs, inserts_possible, deletes_possible, [1])))
    ## Loop over the hidden states
    ## See 17.54 through 17.58 of Murphy's book
    for j in range(3*(r+1)):
        # The emission probabilities are combined with the prior state and transitions
        log_alpha[t+1, j] = np.logaddexp.reduce(log_alpha[t] + log_A[:, j] + all_emission_probs)
```

---

## Full Modification

```python
import numpy as np

def phmm_convert_observations(token, observations):
    """
    Arguments:
        token (str): The token being profiled in the phmm
        observations (list): Observations
    Returns:
        obs (list): Observations converted into a present/absent list of vectors. (0 or 1)
                    There should be len(token) possible values, one for each character, at time t.
    """
    obs = np.full((len(observations), len(token)), 0.0)
    for idx, o_char in enumerate(observations):
        if o_char in token:
            for tidx, char in enumerate(token):
                if char == o_char:
                    obs[idx][tidx] = 1
    return obs

def phmm_forward_pass_log(log_pi, log_A, log_B, obs, token):
    """Computes the forward probabilities (alpha) in log space.
    This is the first part of the equation below:
      p(z_t = j|x_{1:T}) = p(z_t = j|x_{1:t})p(x_{t+1:T}|z_t = j)
    Arguments:
        log_pi (numpy.ndarray): Initial state probabilities
        log_A  (numpy.ndarray): State transition probabilities
        log_B  (numpy.ndarray): Emission probabilities
        obs (list): Observations converted into a present/absent list of vectors. (0 or 1)
                    There should be len(token) possible values, one for each character, at time t.
    """
    T = len(obs)
    r = len(token)
    # The number of hidden states, but split between observed (m and insertion) and unobserved (deletion)
    N = log_A.shape[0]
    # Initialization (t=0)
    log_alpha = np.full((T+1, 3*(r+1)), -np.inf)
    # Begin with the probabilities of seeing the first observation for the emissions of each hidden state
    # This will always be the bb state, regardless of any observations
    #log_alpha[0, :] = log_pi
    # We'll just call that state -1, for simplicity
    match_slice = slice(1, r+1)
    insert_slice = slice(r+1, 2*(r+1))
    del_slice = slice(2*(r+1), -2)
    # The first insertion state can occur before we match any character
    # Why don't we just say that we begin in the insertion state?
    # Because this is how Gordillo and Conde write this
    log_alpha[0, r+1] = log_A[0,r+1]
    log_alpha[0, 1] = log_A[0,1]
    log_alpha[0, 2*(r+1)] = log_A[1, 2*(r+1)]

for t in range(T):
        # Matching probabilities for any observed tokens
        match_probs = np.exp(log_B[match_slice]) @ obs[t]
        # The first insertion state can always insert. Other insertions only work if there was no match.
        inserts_possible = np.concatenate(([1], 1 - match_probs))
        # Deletes are always possible
        deletes_possible = np.ones((r))
        all_emission_probs = np.log(np.concatenate(([1], match_probs, inserts_possible, deletes_possible, [1])))
        ## Loop over the hidden states
        ## See 17.54 through 17.58 of Murphy's book
        for j in range(3*(r+1)):
            # The emission probabilities are combined with the prior state and transitions
            log_alpha[t+1, j] = np.logaddexp.reduce(log_alpha[t] + log_A[:, j] + all_emission_probs)
        ## ee state. Can be transitioned into from the current m, i, and d states
        ## This means that we could change to it and another state at the same time,
        ## which could mess up some assumptions. For now, we won't do this.
        ##log_alpha[t+1,-1] = np.logaddexp.reduce(log_alpha[t+1]+log_A[:,-1])

# If a user wants to know if the entire token was possible, then they should look
    # at log_alpha[t,-1], the ee probability.

# Return the sequence and the likelihood of the entire thing
    # This is less meaningful with pattern matching
    log_likelihood = np.logaddexp.reduce(log_alpha[T-1, :])

return log_alpha, log_likelihood

def initialize_phmm(token):
    """
    Arguments:
        token (list): The token to be profiled
    returns:
        A, B, pi: The transition matrix, emission matrix, and initial state estimation
    """
    r = len(token)
    # The number of hidden states is 3(r+1)
    k = 3*(r+1)
    # Initialize the transition matrix with 0s.
    # It will be sparsely initialized
    transition_matrix = np.zeros((k,k))
    ######### The hidden states are bb, i_0, m_1, d_1, i_1, m_2, ... m_r, d_r, i_r, ee
    # The hidden states are bb, m_1, ..., m_r, i_0, ..., i_r, d_1, ..., d_r, ee
    # Only certain transitions are allowed.
    # The end state is final
    transition_matrix[-1][-1] = 1
    # We'll begin by guessing that the next letter comes at a rate of 0.5 and
    # insertions and deletions are both 0.25
    # Token character states, m, and bb
    for idx in range(0, r):
        # Next insertion
        transition_matrix[idx][idx+r+1] = 0.2
        if idx != r:
            # Next match
            transition_matrix[idx][idx+1] = 0.7
            # Next deletion
            transition_matrix[idx][idx+2*(r+1)] = 0.1
    # The last match state transitions to insertion or ee
    transition_matrix[r][2*r+1] = 0.3
    transition_matrix[r][-1] = 0.7

# Initial insertion states
    for idx in range(r+1, 2*r+1):
        # Stay in insertion: 0.2
        # Move to match 0.75
        # Delete 0.05
        transition_matrix[idx][idx] = 0.2
        transition_matrix[idx][idx-r] = 0.75
        transition_matrix[idx][idx+r+1] = 0.05
    # The ee transition
    transition_matrix[2*r+1][2*r+1] = 0.2
    transition_matrix[2*r+1][-1] = 0.8

# Initial deletion states
    for idx in range(2*(r+1), k-1):
        # Transition to current insertion, next deletion state, or next match
        # Set probabilities at 0.05, 0.15, and 0.8
        transition_matrix[idx][idx-r] = 0.05
        if idx < k-2:
            transition_matrix[idx][idx+1] = 0.15
            transition_matrix[idx][idx-2*r] = 0.8
        else:
            # Transition to ee
            transition_matrix[idx][-1] = 0.95

# On to the emission matrix
    # The grammar is every character in the token
    # Note: If a character shows up twice in the token, it shows up twice in the emission matrix.
    # That means that the characters have their positions encoded in the matrix itself
    # The emission matrix also contains Q, the class of vocabulary not in token, and *, the deletion emission.
    # Note that in "An HMM for detecting spam mail" by Gordilla and Conde, they separated
    # the insertion class into L (all characters not in token) and Q (all non-printing characters)
    # We're going to see if we can do without that, using just Q and *.
    # The order of emission symbols will be [token] + [Q, *]
    # Q and * are really virtual columns.
    # Matching Q just means matching the full row, and matching * means nothing (since the match is unobserved).

emission_matrix = np.zeros((k, r))
    # If the m state has probability 1 for the expected token, and we demand
    # that the tokens have an ordering, then m_i only emits token[i]
    for idx in range(r):
        emission_matrix[1+idx][idx] = 1.0
    # Gordillo and Conde have nonzero probabilities for other states, but I'm
    # just going to simplify things.
    # (Maybe add in some other probabilities for transposition?)

# We must begin in state bb (the initial word state)
    pi = np.full((k), 0.0)
    pi[0] = 1

return transition_matrix, emission_matrix, pi
```

---

## Payoff

* We haven't trained anything, but let's just use it

```python
import argparse
import numpy as np
import pandas
import sys

import hmm_funcs as hmm

# Define the data types for each column in the CSV.
# This assumes a common preprocessed format: 'label', 'subject', 'message_body'.
# 'i4' is a 4-byte integer for the label (0 or 1).
# 'U50' is a Unicode string of up to 50 characters for the subject.
# 'U256' is a Unicode string of up to 256 characters for the message body.
data_type = np.dtype([
    ('label', 'i4'),
    ('subject', 'U50'),
    ('email_to', 'U50'),
    ('email_from', 'U50'),
    ('message', 'U256')
])

def load_trec_csv(filepath):
    """
    Args:
        filepath (str): The path to the CSV file to load.
    Returns:
        numpy.ndarray: Pandas dataset
    """
    print(f"Attempting to load data from: {filepath}")
    try:
        return pandas.read_csv(filepath, quotechar='"')

except IOError as e:
        print(f"Error: Could not open or read file at path '{filepath}'.")
        sys.exit(1)
    except Exception as e:
        print(f"An unexpected error occurred during parsing: {e}")
        sys.exit(1)

if __name__ == "__main__":
    parser = argparse.ArgumentParser()
    parser.add_argument(
        "path",
        type=str,
        help="Path to the Trec2007 csv")
    parser.add_argument(
        "--word",
        required=True,
        type=str,
        help="The word to create a filter for.")

args = parser.parse_args()

# Load the data
    trec_data = load_trec_csv(args.path)
    np_trec = trec_data.to_numpy()

if trec_data.size > 0:
        print(f"Shape of the data: {np.shape}")
        print(f"Data fields: {np.dtype.names}")

# Fields should be label,subject,email_to,email_from,message
        labels = trec_data['label']

print("\nFirst 3 records:")
        # Display the records neatly
        for i in range(min(3, len(trec_data))):
            print(f"  Label: {trec_data['label'][i]}, Subject: '{trec_data['subject'][i]}', Body Snippet: '{trec_data['message'][i][:50]}...'")

# Count the spam and ham
        num_spam = np.sum(labels == 1)
        num_ham = np.sum(labels == 0)

print(f"\nTotal records: {len(labels)}")
        print(f"Spam emails (1): {num_spam}")
        print(f"Ham emails (0): {num_ham}")

phmm_params = [np.log(param) for param in hmm.initialize_phmm(args.word)]

# We could just forward through the corpus, seeing which words are possible matches
    A, B, pi = hmm.initialize_phmm(args.word)
    log_A = np.log(A)
    log_B = np.log(B)
    log_pi = np.log(pi)
    
    num_spam_matches = 0
    num_ham_matches = 0
    for record in range(trec_data.size):
        obs = hmm.phmm_convert_observations(args.word, str(trec_data['message'][record]).lower())
        #print(trec_data['message'][0])
        match_in_record = False

# Loop through chunks of the messages, feeding them into the model
        for idx in range(len(obs) - len(args.word)):
            S = obs[idx:idx+2*len(args.word)]
            log_alpha, log_likelihood = hmm.phmm_forward_pass_log(log_pi, log_A, log_B, S, args.word)
            if np.max(log_alpha[:,-1] > -6):
                match_in_record = True
                print("############################################################")
                print("Match found in:")
                print(trec_data['message'][record][idx:idx+2*len(args.word)])
                print("Probabilities:")
                print(log_alpha[:,-1])
                print("############################################################")
        if match_in_record:
            if labels[record] == 0:
                num_ham_matches += 1
            else:
                num_spam_matches += 1
    print(f"\nTotal records: {len(labels)}")
    print(f"Spam emails (1): {num_spam}")
    print(f"Ham emails (0): {num_ham}")
    print("Found {num_ham_matches} ham matches and {num_spam_matches}.")

sys.exit(0)
```

---

## Results?

<pre>
Match found in:
|/IAGRA- cla
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
        -inf -5.7587489  -5.48896803 -5.44217472 -5.43293618 -5.43107642
 -5.43084088]
</pre>

---

## Good results?

<pre>
Match found in:
V ia gr a
Probabilities:
[        -inf         -inf         -inf         -inf         -inf
         -inf -10.29903327  -8.25430923  -6.47155141  -5.81507491]
</pre>

---

## Bad Results?

<pre>
Match found in:
V ia gr a
Probabilities:
[        -inf         -inf         -inf         -inf         -inf
         -inf -10.29903327  -8.25430923  -6.47155141  -5.81507491]
</pre>

---

## Needs Training

* As with everything, we need to train
* Best method?
  * Imagine two sets of clusters
    * one with regular words
    * one with banned words

---

## Sounds Like Clustering

* It is clustering!
  * Make some HMMs for regular words, some for banned
  * Then train a mixture model and optimize to parameters to best explain the dataset
* To simplify, we could train initial models on a "safe" corpus
  * Published books, for example

---

## Takeaways

* We don't need to do that entire exercise
  * But understand how techniques are used together
* By inserting some knowledge into the model, we can build something robust

---

## Drawback?

* We need to insert our own knowledge into the system
* Some relationship are too copmlicated
  * Pixels in an image, for example
  * Theoretically, we could describe shapes with HMMs, but *could* we?
* That will motivate the next set of ML algorithms

<!--
############################################################
Match found in:
ving announc
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
 -5.45278132 -5.18481769 -5.13768716 -5.12706549 -5.12493309 -5.12448449
 -5.12438986]
############################################################
Match found in:
VI_TRA- exot
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
 -5.27511428 -5.00441677 -4.95731245 -4.94698786 -4.94492939 -4.94450188
 -4.94441285]
############################################################
############################################################
Match found in:
|/IAGRA- cla
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
        -inf -5.7587489  -5.48896803 -5.44217472 -5.43293618 -5.43107642
 -5.43084088]
############################################################
############################################################
Match found in:
vNAbragFvzw2
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
 -5.49336784 -5.22299369 -5.17599125 -5.16560806 -5.16356958 -5.16314929
 -5.16306249]
############################################################
############################################################
Match found in:
Vz+GraGJapha
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
 -5.31181493 -5.04177456 -4.99490237 -4.98447791 -4.98115777 -4.98048352
 -4.98039594]
############################################################
############################################################
Match found in:
via Brazil C
Probabilities:
[       -inf        -inf        -inf        -inf        -inf        -inf
 -6.86504607 -5.27466389 -5.05205104 -5.01223265 -5.00420107 -5.00258657
 -5.00225752]
############################################################

-->

<!--

Review: HMMs

Show curves for accuracy and A and B estimates over sample length with the coins and die.

Discussion:

Extension:

Let's cover a couple of applications.

1. classification of word tokens {SMS, email}
* Just use the current HMM on those datasets

2. detecting and correcting misspelled words
* Create a HMM with all the states described in the paper

Stochastic memoizer for sequence data?

-->