Hidden Markov Model (HMM)

A Hidden Markov Model (HMM) is a probabilistic model for sequential data where an unobserved (hidden) state follows a Markov chain, and each state generates an observed output with state-dependent probabilities. It is important in AI/ML for inferring latent structure, smoothing, and decoding sequences under uncertainty. Example: speech recognition, where hidden phoneme states emit acoustic feature vectors.

Imagine you’re watching someone behind a curtain. You can’t see what they’re doing, but you can hear clues—like footsteps, humming, or silence. A Hidden Markov Model (HMM) works like that: it assumes there are “hidden” situations (called states) you can’t directly observe, and you only see noisy hints (called observations).

It also assumes the hidden state changes over time in a step-by-step way, where the next state mostly depends on the current one. In AI, HMMs are used for things like speech recognition, where the words are hidden but the sounds are observed.

import numpy as np  # numerical arrays for time-series features
import pandas as pd  # tabular handling for timestamps and observations
from hmmlearn.hmm import GaussianHMM  # HMM with Gaussian emissions (common for real-valued signals)

# --- Create a realistic "wearable sensor" time series (e.g., activity intensity) ---
rng = np.random.default_rng(7)  # reproducible randomness for a stable demo
n = 240  # e.g., 240 minutes of data

# Simulate 3 latent regimes: resting (low), walking (medium), running (high)
true_states = np.repeat([0, 1, 2], [90, 90, 60])  # piecewise regimes over time
x = np.empty(n)  # observed sensor feature (real-valued)
x[true_states == 0] = rng.normal(loc=0.2, scale=0.10, size=(true_states == 0).sum())  # resting
x[true_states == 1] = rng.normal(loc=1.0, scale=0.20, size=(true_states == 1).sum())  # walking
x[true_states == 2] = rng.normal(loc=2.2, scale=0.25, size=(true_states == 2).sum())  # running

df = pd.DataFrame({"t": pd.date_range("2026-01-01", periods=n, freq="min"), "intensity": x})  # timestamped series

# --- Fit an HMM to recover latent regimes from the observed time series ---
X = df[["intensity"]].to_numpy()  # hmmlearn expects shape (n_samples, n_features)
hmm = GaussianHMM(n_components=3, covariance_type="diag", n_iter=200, random_state=7)  # 3 hidden states
hmm.fit(X)  # learn transition probs + emission parameters from data

z_hat = hmm.predict(X)  # *** KEY LINE *** infer the most likely hidden state sequence for each time step

# Summarize what each inferred state "looks like" in observation space
state_means = hmm.means_.ravel()  # per-state mean intensity learned by the model
order = np.argsort(state_means)  # sort states low->high to map to intuitive regimes
df["state_hat"] = pd.Categorical(order.searchsorted(z_hat), ordered=True)  # relabel states by mean level

print("Learned state means (sorted low->high):", np.sort(state_means))
print("Learned transition matrix (rows sum to 1):\n", hmm.transmat_)
print(df.head(8))  # show first few timestamps with inferred state labels

• The inferred sequence z_hat is the per-timestep latent regime assignment produced by GaussianHMM.predict, typically used for segmentation or regime detection before downstream modeling.
• GaussianHMM.fit estimates both the transition matrix (transmat_) and emission parameters (means_, covars_), which can be used as features (e.g., regime durations, switching rates) in ML pipelines.
• State labels are arbitrary; sorting by means_ (as done here) avoids misinterpreting “state 0/1/2” as ordered regimes without checking.
• Common mistake: passing a 1D array into hmmlearn; it requires a 2D feature matrix of shape (n_samples, n_features).

library(depmixS4)   # HMMs via dependent mixture models (widely used in R)
library(dplyr)      # convenient data manipulation

set.seed(7)  # reproducibility for the simulated time series
n <- 240     # e.g., 240 minutes of wearable sensor data

# Simulate 3 latent regimes: resting (low), walking (medium), running (high)
true_states <- c(rep(1, 90), rep(2, 90), rep(3, 60))  # piecewise regimes over time
x <- numeric(n)  # observed sensor feature (real-valued)
x[true_states == 1] <- rnorm(sum(true_states == 1), mean = 0.2, sd = 0.10)  # resting
x[true_states == 2] <- rnorm(sum(true_states == 2), mean = 1.0, sd = 0.20)  # walking
x[true_states == 3] <- rnorm(sum(true_states == 3), mean = 2.2, sd = 0.25)  # running

df <- tibble(
  t = seq.POSIXt(from = as.POSIXct("2026-01-01 00:00:00"), by = "min", length.out = n),
  intensity = x
)

# Fit a 3-state Gaussian-emission HMM to the observed intensity series
mod <- depmix(intensity ~ 1, data = df, nstates = 3, family = gaussian())  # intercept-only emissions per state
fit_mod <- fit(mod, verbose = FALSE)  # estimate transitions + emission parameters

post <- posterior(fit_mod)  # *** KEY LINE *** compute per-time posterior and most-likely state (Viterbi-like label)

# Extract and order states by their fitted emission means for interpretability
pars <- getpars(fit_mod)  # packed parameter vector (includes emission means and transition logits)
state_means <- sapply(1:3, function(s) {
  # depmix stores Gaussian response parameters per state; "response" pulls them in structured form
  response(fit_mod)[[s]][[1]]@parameters$mean
})
order <- order(state_means)  # low->high mean intensity
df_out <- df %>%
  mutate(state_hat = match(post$state, order))  # relabel inferred states by mean level

print(sort(state_means))     # learned emission means (sorted low->high)
print(head(df_out, 8))       # first few timestamps with inferred state labels

• post$state from posterior is the per-timestep inferred state path, commonly used to tag each time point with a regime for feature engineering or monitoring.
• depmix with intensity ~ 1 fits state-specific intercepts (means) for a Gaussian emission; adding covariates to the formula is a typical way to make emissions depend on additional signals.
• State indices are not ordered; sorting by fitted emission means (via response) helps prevent mislabeling regimes when comparing runs or models.
• Common mistake: assuming posterior() returns only hard labels—inspect posterior probabilities too when you need uncertainty-aware decisions (e.g., filtering low-confidence segments).