EmissionModels.jl

A Julia package providing emission models for HiddenMarkovModels.jl. It supplies ready-to-use distributions that describe how observations are generated conditioned on the HMM's latent states.

Quick start

using Pkg
Pkg.add(["EmissionModels", "HiddenMarkovModels"])

using EmissionModels
using HiddenMarkovModels

# Create an emission model
dist = PoissonZeroInflated(5.0, 0.3)

# Sample and evaluate densities
x = rand(dist)
logp = logdensityof(dist, x)

# Fit to weighted observations (weights are HMM posterior state probabilities)
obs_seq = [rand(dist) for _ in 1:100]
weight_seq = ones(100)
fit!(dist, obs_seq, weight_seq)

Distribution models

All types implement the HiddenMarkovModels emission interface (rand, logdensityof, fit!).

Count data

Type	Description
PoissonZeroInflated(λ, π)	Zero-inflated Poisson for excess zeros in count data.

Multivariate continuous

Type	Description
MvT(μ, Σ, ν)	Full-covariance multivariate Student's t.
MvTDiag(μ, σ², ν)	Diagonal-covariance multivariate Student's t.

GLM emissions (observation depends on a control vector)

Type	Description
GaussianGLM(β, σ²)	Linear regression with Gaussian noise.
BernoulliGLM(β)	Logistic regression for binary data.
PoissonGLM(β)	Log-linear regression for count data.
MultinomialGLM(B, n)	Softmax regression over categories, for count vectors or choice labels.
MvGaussianGLM(B, Σ)	Multivariate linear regression with shared full covariance.
MvBernoulliGLM(B)	Independent logistic regressions, one per output dimension.
MvPoissonGLM(B)	Independent log-linear regressions, one per output dimension.

GaussianGLM, BernoulliGLM, and PoissonGLM carry a coefficient vector β and emit scalars; the rest carry a coefficient matrix B and emit vectors. All of them support regularization via priors:

β = zeros(3)
glm = GaussianGLM(β, 1.0, RidgePrior(0.5))  # L2 regularization

Each GLM is fit via fit!(glm, y, w; control_seq=X), where control_seq (design matrix X) maps latent states to the regression covariates. Since the GLMs subtype ControlledEmission, a vector of them also works directly as the emissions of a ControlledEmissionHMM.

DDM emissions (two-alternative forced choice)

Type	Description
StimulusCodedDDM(; ν, α, z, τ)	Drift diffusion model whose drift sign follows the per-trial stimulus code.
CoherenceDDM(; k, γ, α, z, τ)	Drift diffusion model whose drift is a power law of signed stimulus coherence.

Both emit (choice, rt) pairs and subtype ControlledEmission, so a vector of them forms a ControlledEmissionHMM — the DDM-HMM, where hidden states are decision-making regimes. The Wiener first-passage-time density and sampler come from SequentialSamplingModels.jl, a weak dependency: the types always construct, but logdensityof, rand, and fit! need using SequentialSamplingModels to activate the extension.

Calcium imaging emissions

Type	Description
CalciumEmission(λ, params)	Autoregressive fluorescence model with the spike count marginalized out.

The model of Keeley, Zoltowski, Charles & Pillow fits an HMM directly to raw calcium traces, with no separate deconvolution step. A latent Poisson spike count drives an AR(p) fluorescence process, and the emission density sums the count out. Only the firing rates λ vary across states; the AR coefficients, spike influx and noise variance live in a shared CalciumParams and are fit with a tied M-step.

using HiddenMarkovModels: ControlledEmissionHMM

λ = [0.15 1.60; 1.70 0.20]                          # neurons × states
dists = calcium_emissions(λ, [0.9 0.85], [1.0, 1.0], [0.02, 0.02])
hmm = ControlledEmissionHMM([0.5, 0.5], [0.97 0.03; 0.03 0.97], dists)

The AR history enters through the control interface, so these are ControlledEmissions like the GLMs and DDMs: each timestep's control is the vector of lagged fluorescence values, built by lagged_controls. Use rand_calcium to simulate (the trace must be rolled out step by step, since each control depends on the previous observation) and init_calcium to initialize a fit from data. For bulk signals with large rates, gaussian=true swaps the truncated Poisson sum for an analytic Gaussian marginal.

ACDC model selection

The Accumulated Cutoff Discrepancy Criterion (ACDC) picks the number of hidden states without penalizing likelihood by parameter count. It inverts each fitted emission through the probability integral transform to recover per-state "stochastic drivers", which are uniform when the model is well specified, and selects the smallest state count whose per-state discrepancies from uniform all fall below a cutoff.

using EmissionModels, HiddenMarkovModels, Distributions

hmm = HMM([0.5, 0.5], [0.95 0.05; 0.05 0.95],
          [Normal(-4.0, 1.0), Normal(4.0, 1.0)])
_, obs_seq = rand(hmm, 3000)

result = component_discrepancies(hmm, obs_seq, KSDiscrepancy())
K = acdc_select([result], 0.05)

It works with any AbstractHMM whose emissions are standard Distributions objects or the types in this package (including the GLMs, via their covariates). Several discrepancy measures are available: KSDiscrepancy, KLDiscrepancy, WassersteinDiscrepancy, MMDDiscrepancy, and SquaredErrorDiscrepancy.

Creating custom emission models

HiddenMarkovModels.jl accepts any type that implements the following interface:

Random.rand(rng::AbstractRNG, dist::MyEmission)
DensityInterface.DensityKind(::MyEmission)      # return HasDensity()
DensityInterface.logdensityof(dist::MyEmission, obs)
StatsAPI.fit!(dist::MyEmission, obs_seq, weight_seq)

See the documentation for details.

Installation

EmissionModels.jl is registered in the General registry:

using Pkg
Pkg.add("EmissionModels")

Requires Julia 1.10 or later.

Contributing

Contributions are welcome. Please follow the Julia Blue Style and add tests for new behavior. Pull requests and issues are appreciated.

Contributors

Thanks goes to these wonderful people (emoji key):

Ryan Senne 🚧 💻 ⚠️ 📖 🤔 👀

Nguyen Nguyen 💻 🤔 👀 ⚠️

dfish 👀

This project follows the all-contributors specification. Contributions of any kind are welcome!

License

EmissionModels.jl is released under the MIT License.