API reference

AbstractStateSpaceModel

Interface for a state-space model. Subtypes implement prior, dyn, and obs to expose their model components. A custom model need not store components in fields with those names. Inference algorithms may require additional structure from the returned components.

DistributionDynamics(f)

Lift a closure f(t, x) -> distribution into a LatentDynamics.

DistributionObservation(f)

Lift a closure f(t, x) -> distribution into an ObservationProcess.

DistributionPrior(dist)

Lift a distribution object into a StatePrior.

LatentDynamics

Transition dynamics of a state-space model. A concrete subtype should implement simulate and/or logdensity, or supply a distribution(dyn, t, x) method from which both are derived. Which methods are required depends on the inference algorithm (e.g. a bootstrap filter needs only simulate, a guided proposal also needs logdensity).

ObservationProcess

Emission process of a state-space model. A concrete subtype should implement logdensity (and optionally simulate for forward simulation), or supply a distribution(obs, t, x) method from which both are derived.

StatePrior

Initial state distribution of a state-space model. A concrete subtype should implement simulate and/or logdensity, or supply a distribution(prior) method from which both are derived. Which methods are required depends on the inference algorithm.

StateSpaceModel(prior, dyn, obs)
StateSpaceModel(model::AbstractStateSpaceModel)

A state-space model composed of an initial state prior, latent dynamics, and an observation process. The one-argument constructor obtains the components through prior, dyn, and obs. It retains those components without copying them; an existing StateSpaceModel is returned unchanged.

distribution(prior::StatePrior)
distribution(dyn::LatentDynamics, t::Integer, x)
distribution(obs::ObservationProcess, t::Integer, x)

Return the distribution associated with a model component. Implementing this method derives simulate and logdensity for free; components may instead implement those directly when no tractable distribution object is available.

Components take no keyword arguments. Dependence on parameters, controls or an outer state is expressed by the component object itself, which downstream packages may build per step.

dyn(model::AbstractStateSpaceModel)

Return the model's latent dynamics, a LatentDynamics.

logdensity(prior::StatePrior, x0)
logdensity(dyn::LatentDynamics, t::Integer, x_prev, x_new)
logdensity(obs::ObservationProcess, t::Integer, x, y)

Evaluate the log-density of a model component, by default that of its distribution.

obs(model::AbstractStateSpaceModel)

Return the model's observation process, an ObservationProcess.

prior(model::AbstractStateSpaceModel)

Return the model's initial state prior, a StatePrior.

simulate(rng::AbstractRNG, prior::StatePrior)
simulate(rng::AbstractRNG, dyn::LatentDynamics, t::Integer, x)
simulate(rng::AbstractRNG, obs::ObservationProcess, t::Integer, x)

Draw from a model component, by default from its distribution.

simulate([rng,] model::AbstractStateSpaceModel, T::Integer)

Simulate a trajectory of length T, returning (x0, xs, ys) where x0 is the initial state, xs the states at times 1:T, and ys the observations at times 1:T.

simulate_from_dist(rng::AbstractRNG, d)

Draw a sample from a distribution object. Defaults to rand(rng, d); specialise it to return static arrays or otherwise control the sample type.

AbstractExecution

How a particle filter evaluates its particle population. Choose it with the execution keyword of ParticleFilter. RBPF, AuxiliaryParticleFilter and ConditionalSMC use the execution of the particle filter they contain.

AbstractLikelihood

Abstract type for backward likelihood representations used in smoothing and ancestor sampling.

Subtypes represent the predictive likelihood p(y | x) in different forms depending on the state space structure (continuous Gaussian vs discrete).

AbstractProposal

Implement distribution(proposal, t, state, observation), or both simulate and logdensity. For an ordinary particle filter state is the previous latent sample. For an RBPF it is the full RBState: proposals can inspect both the outer sample and inner filtering belief, but must draw and score only the new outer sample. The inner filter is advanced by RBPF after this draw.

Per-call wrapper that resolves its component and stamps the activity flags.

AncestorSampling([backward_predictor]) <: AbstractTrajectoryRefreshment

Conditional SMC with ancestor sampling (CSMC-AS / PGAS). At resampling events, the reference particle's ancestor is resampled using backward weights, improving mixing for the full trajectory. The PF resampling trigger is respected: ancestor sampling is performed when ordinary resampling occurs. At skipped steps, every particle retains its ancestor and accumulated weight. For RBPF models, backward predictive likelihoods are computed at the start of each sweep which enable closed form ancestor weights. The optional backward predictor uses the analytical filter's default when omitted. Gaussian defaults apply no jitter; nonzero backward jitter is rejected because it changes the target distribution.

AuxiliaryResampler

A resampling scheme for multistage particle resampling with auxiliary weights

BackwardDiscretePredictor <: AbstractBackwardPredictor

Recursively computes the backward likelihood β_t(i) = p(y_{t:T} | x_t = i) for discrete state-space models, in log-space for numerical stability.

BackwardInformationPredictor(; initial_jitter=nothing, jitter=nothing)

Recursively computes the predictive likelihood p(y_{t:T} | x_t) of a linear-Gaussian model in information form. The jitter fields add numerical-stability shifts to the precision matrix at initialisation and during backward prediction.

Based on https://arxiv.org/pdf/1505.06357.

BackwardSimulation([backward_predictor]) <: AbstractTrajectoryRefreshment

Conditional SMC with backward simulation (CSMC-BS). Runs a full forward filter with particle storage, then samples a trajectory via a backward pass using backward sampling weights. For RBPF models, backward predictive likelihoods are computed on-the-fly during the backward pass.

Note: requires O(N*T) storage since the full particle history must be retained.

CSMCModel{MT, YT} <: AbstractMCMC.AbstractModel

Model wrapper for standalone CSMC sampling via the AbstractMCMC interface.

Fields

• ssm::MT: The state-space model
• observations::YT: Vector of observations

CSMCState{TT}

State of a conditional SMC sampler, containing the current reference trajectory.

The trajectory is a ReferenceTrajectory indexed from 0 (matching the prior at time 0). For RBPF, the trajectory contains outer states only. Inner beliefs are recomputed for the current parameters on every sweep.

CombinedRNG(cpu::AbstractRNG, gpu::AbstractRNG)

Bundle explicitly owned CPU and GPU generators for the positional RNG argument. Scalar draws, distributions and ordinary Julia arrays use cpu; supported GPU array methods and GPU algorithms use gpu. Allocating rand/randn calls return CPU arrays. Use cpu_rng and gpu_rng at backend boundaries.

Both children must be distinct, explicitly owned generators: nested bundles, TaskLocalRNG and identical child objects are rejected. The bundle retains the supplied generators, so drawing through either reference advances the same state. copy copies both children; child generators must support copy for checkpoints.

seed!(rng, seed::Integer) reseeds both children. The CPU receives seed; the GPU receives one UInt64 drawn from a separately seeded copy of the CPU generator, without consuming the live CPU stream. Children must support integer seed!. Reproducibility requires preserving both states and the sequence of operations; concurrent use requires the same care as using the child generators directly.

ConditionalSMC{PF, TR} <: AbstractMCMC.AbstractSampler

Conditional Sequential Monte Carlo sampler with configurable trajectory refreshment.

Fields

• pf::PF: The underlying particle filter (e.g., BF(N), RBPF(BF(N), KF()))
• refreshment::TR: Trajectory refreshment strategy

Examples

ConditionalSMC(BF(100))                                  # Vanilla CSMC (NoRefreshment default)
ConditionalSMC(BF(100), AncestorSampling())              # CSMC with ancestor sampling
ConditionalSMC(RBPF(BF(200), KF()), AncestorSampling())  # Rao-Blackwellised PGAS

CovarianceFactor(F)

Covariance represented by an explicit factor Σ = F*F'. The factor may be rectangular or rank deficient. Square-root algorithms consume F directly, without reconstructing and refactorising Σ. This does not add noise or repair the statistical model.

CovarianceRepair

A positive-definiteness repair strategy applied to a filtered covariance after each update. Repair sits outside the fused Kalman step so its pullback composes by the ordinary chain rule.

DenseParticleContainer{T0, T, WT}

Full-history particle storage: the initial particle states (time 0), plus subsequent particle states, weights, and ancestor indices for times 1, 2, …, T.

T0 is the type of the initial particle states; T is the type of the subsequent particle states. When they match, use the one-argument constructor — pushing subsequent states with a different type will raise an error.

Fields

• initial_states: Time-0 particle states (fixed length N)
• states: Vector of particle populations at times 1..T
• weights: Vector of log-weight arrays at times 1..T
• ancestors: Vector of integer ancestor arrays at times 1..T

Each population retains its CPU or device storage; the outer time vectors live on the CPU.

Constructors

DenseParticleContainer(initial, first_state)
DenseParticleContainer(initial_states)
DenseParticleContainer(initial_states, T)
DenseParticleContainer(initial_states, states_t1, weights_t1, ancestors_t1)

With particle distributions, DenseParticleContainer(initial, first_state) infers both state and weight types from the first completed step. Append later results with push!(history, state). For the raw-vector constructors:

• First form: assumes T == T0.
• Second form: explicit subsequent-state type T (the container starts empty for t≥1). These two raw-vector forms use Float64 weight storage. Use first-step data to infer a different weight type.
• Third form: construct directly from the initial states and the time-1 states, weights, and ancestors — T is inferred from states_t1. Prefer this when the initial and subsequent state types differ so that T is inferred from concrete data.

Ownership

Constructors and append operations copy collection buffers but retain the state objects inside them. Do not mutate stored state objects. Use an explicit deepcopy of an input particle distribution when snapshots of mutable states are needed. Returned trajectories may share their state objects with storage. Device batches copy their batched leaves; trajectory extraction copies selected members so a returned path does not retain full particle populations.

Mutation

Use push!(container, states, weights, ancestors) to append a new time step. Passing a states vector whose eltype does not match T raises an ArgumentError.

DiscreteDynamics(P)

Transition over a finite state space, with P[i, j] = p(x_t = j | x_{t-1} = i).

DiscreteFilter <: AbstractFilter

Forward filtering algorithm for discrete (finite) state-space models. Computes the filtered distribution π_t(i) = p(x_t = i | y_{1:t}).

DiscreteLikelihood <: AbstractLikelihood

A container representing the backward likelihood β_t(i) = p(y | x = i) for discrete state spaces, stored in log-space for numerical stability.

This representation is used in backward filtering algorithms for discrete SSMs (HMMs) and Rao-Blackwellised particle filtering with discrete inner states.

Fields

• log_β::VT: Vector of log backward likelihoods, where log_β[i] = log p(y | x = i)

See also

• BackwardDiscretePredictor: Algorithm that uses this representation

DiscretePrior(α0)

Initial distribution over a finite state space, with α0[i] = p(x_0 = i).

DiscreteSmoother <: AbstractSmoother

Forward-backward smoother for discrete state-space models.

Eigenvalue clipping Σ ↦ V max(Λ, ε) Vᵀ.

GPUExecution()

Initialise a device particle population using the model's GPU initialisation method. Subsequent operations dispatch on that population's storage, while time iteration and scalar decisions and CSMC selected-path likelihoods run on the CPU. The same model must support CPU component resolution as well as batched device methods. This setting does not automatically upload an arbitrary CPU model or make its callbacks GPU-compatible.

Model implementations can extend initialise(::GPUExecution, rng, prior, algo; ref_state) to select their resident device parameters and construct batched particle storage. The BatchedKernels extension provides this method for supported device Gaussian priors.

GaussianPrior(μ0, Σ0)

Gaussian initial-state prior with mean μ0 and covariance Σ0, retaining both array representations. Structured covariance matrices are supported.

GaussianState(μ, Σ)

A multivariate Gaussian with mean μ and covariance Σ, retaining the supplied array representations. Covariances are symmetric by convention; positive-definiteness is not enforced at construction. Kalman initialisation selects full matrix storage independently of model parameters; subsequent calculations retain their natural array representations. Supports rand and logpdf when the covariance meets those operations' requirements.

HierarchicalDynamics(outer, inner)

Dynamics for a Rao-Blackwellisable model. outer is the dynamics for the sampled component; inner is either a constant dynamics atom or a conditioning callable (; t, x_prev, x_new) -> atom.

HierarchicalObservation(inner)

Observation process for a Rao-Blackwellisable model. inner is either a constant observation atom/emission or a conditioning callable (; t, x) -> atom/emission.

HierarchicalPrior(outer, inner)

Prior for a Rao-Blackwellisable model. outer is the prior for the sampled component; inner is either a constant prior atom or a conditioning callable (; x0) -> atom.

HierarchicalState(x, z)

A joint sample from a hierarchical model, with outer component x and a sampled inner component z. This differs from an RBState, whose inner component is a distribution rather than a sample.

InformationLikelihood <: AbstractLikelihood

A container representing an unnormalized Gaussian likelihood p(y | x) in information form, parameterized by natural parameters (λ, Ω).

The unnormalized log-likelihood is given by: log p(y | x) ∝ λ'x - (1/2)x'Ωx

This representation is particularly useful in backward filtering algorithms and Rao-Blackwellised particle filtering, where it represents the predictive likelihood p(y{t:T} | xt) conditioned on future observations.

Fields

• λ::λT: The natural parameter vector (information vector)
• Ω::ΩT: The natural parameter matrix (information/precision matrix)

See also

• natural_params: Extract the natural parameters (λ, Ω)
• BackwardInformationPredictor: Algorithm that uses this representation

Additive diagonal shift Σ + ε I.

KalmanFilter(; repair=NoRepair())

Kalman filter for linear-Gaussian state-space models. repair is a CovarianceRepair strategy applied to the filtered covariance after each update.

LinearGaussianDynamics(A, b, Q)

Linear-Gaussian transition x_t = A x_{t-1} + b + w, w ~ N(0, Q), with array parameters, retaining the representation of Q.

LinearGaussianObservation(H, c, R)

Linear-Gaussian emission y_t = H x_t + c + v, v ~ N(0, R), with array parameters, retaining the representation of R.

NoRefreshment <: AbstractTrajectoryRefreshment

Vanilla conditional SMC with no trajectory refreshment. Uses filter() directly with a reference trajectory pinned to particle 1.

ParameterisedSSM(build, observations)

A parameterised state-space model that maps parameter vectors to concrete SSMs.

Fields

• build: A callable θ -> StateSpaceModel that constructs an SSM from parameters. Fixed model components should be captured via closure.
• observations: The observation sequence y₁:T.

Example

function build_model(θ, fixed)
    b = θ[1:2]
    dyn = LinearGaussianDynamics(fixed.A, b, fixed.Q)
    return StateSpaceModel(fixed.prior, dyn, fixed.obs)
end

pssm = ParameterisedSSM(θ -> build_model(θ, fixed), observations)
model = pssm.build(θ)  # returns a concrete SSM

Particle

A container representing a single particle in a particle filter distribution, composed of a weighted sampled (stored as a log weight) and its ancestor index.

Particle(container::DenseParticleContainer, t::Integer, i::Integer)

Reconstruct the Particle at time t ≥ 1, index i, from the container's stored state, weight, and ancestor index.

ParticleDistribution

A container for particle filters which composes a collection of weighted particles (with their ancestories) into a distibution-like object.

Fields

• particles::VT: Vector of weighted particles
• ll_baseline::WT: Baseline for computing log-likelihood increment. A scalar that caches the unnormalized logsumexp of weights before update (for standard PF/guided filters) or a modified value that includes APF first-stage correction (for auxiliary PF).

ParticleFilter(N, proposal; threshold=1.0, resampler=Systematic(), execution=SerialExecution())

Particle filter with N particles drawn from proposal, an AbstractProposal. Particles are resampled with resampler when the effective sample size falls below threshold * N. execution chooses how the particle population is evaluated; see AbstractExecution. BootstrapFilter(N; kwargs...) proposes from the latent dynamics.

ParticleGibbs{CS, PS, ADT} <: AbstractMCMC.AbstractSampler

Particle Gibbs sampler that alternates between a parameter update (e.g., NUTS) and a trajectory update (conditional SMC).

Fields

• csmc::CS: Conditional SMC sampler for trajectory updates (e.g., ConditionalSMC(RBPF(BF(200), KF()), AncestorSampling()))
• param::PS: Parameter sampler (e.g., AdvancedHMC.NUTS(0.8))
• adtype::ADT: AD backend (ADTypes.AbstractADType). nothing uses AdvancedHMC's default (ForwardDiff). Both regular and hierarchical models support ForwardDiff and Mooncake. Load the package corresponding to the selected AD backend.

Examples

# Regular SSM
ParticleGibbs(ConditionalSMC(BF(100), NoRefreshment()), NUTS(0.8))

# Hierarchical SSM with reverse-mode AD
ParticleGibbs(ConditionalSMC(RBPF(BF(200), KF()), AncestorSampling()), NUTS(0.8); adtype=AutoMooncake())

ParticleGibbsModel{PT, MT} <: AbstractMCMC.AbstractModel

Model for particle Gibbs inference, combining a prior on parameters with a parameterised SSM.

Fields

• prior::PT: Multivariate prior distribution on the parameter vector θ
• param_model::MT: A ParameterisedSSM mapping θ to a concrete SSM

Examples

pssm = ParameterisedSSM(θ -> build_model(θ, fixed), observations)
model = ParticleGibbsModel(MvNormal(zeros(d), 4.0*I), pssm)

ParticleGibbsState{VT, TT, PS}

Internal state of the particle Gibbs sampler.

Fields

• θ: Current parameter vector
• trajectory: Current reference trajectory (ReferenceTrajectory)
• param_state: Parameter sampler state (e.g., AdvancedHMC.HMCState)

ParticleGibbsTransition{VT, NT}

A single transition of the particle Gibbs sampler, containing the parameter values and diagnostics from the parameter sampler.

ParticleTree{T0, T}

Sparse storage of particle ancestry with separate buffers for initial states (time 0) and subsequent states (times 1..T). The initial-state buffer is never expanded or pruned; the subsequent-state buffer grows as needed.

T0 and T are the types of initial and subsequent states respectively.

Fields

• initial_states::Vector{T0}: Time-0 particle states (fixed length N)
• states::Vector{T}: Pool for times 1..T (expandable)
• parents::Vector{Int64}: For each entry in states, the parent's index
• is_penultimate::Vector{Bool}: true if parents[i] indexes initial_states, false if it indexes states
• leaves::Vector{Int64}: Current leaf indices (into initial_states before the first insert, into states thereafter)
• offspring::Vector{Int64}: Offspring counts for states
• free_indices::Stack{Int64}: Unused slots in states
• leaves_in_initial::Ref{Bool}: true until the first insert! call

Constructors

ParticleTree(initial; capacity)
ParticleTree(initial, first_state; capacity)
ParticleTree(initial_states, M)
ParticleTree(initial_states, T, M)
ParticleTree(initial_states, states_t1, ancestors_t1, M)

For particle distributions, ParticleTree(initial) assumes the later state type matches the initial type. ParticleTree(initial, first_state) infers the later type from the first completed step and stores that step, including its links to time zero. Append each later result with push!(tree, state), which validates input, prunes, and inserts it. capacity controls the initial storage allocation and grows automatically. For raw-vector inputs:

• First form: assumes T == T0.
• Second form: explicit subsequent-state type T.
• Third form: construct with time-1 data already inserted — T is inferred from states_t1. Prefer this when initial and subsequent types differ.

Ownership

Collection buffers are copied, but the state objects they contain are shared. Do not mutate stored states, including through a retrieved trajectory. Copy an input particle distribution explicitly with deepcopy when an independent snapshot is needed.

Reference

Jacob, P., Murray L., & Rubenthaler S. (2015). Path storage in the particle filter doi:10.1007/s11222-013-9445-x

RBPF(particle_filter, analytical_filter)

Rao–Blackwellised particle filter for hierarchical models. Particle proposals receive an RBState and return the outer sample; the analytical filter then resolves and updates the conditional inner model. Reference trajectories contain outer samples only.

RBState

A container representing a single state with a Rao-Blackwellised component. This differs from a HierarchicalState which contains a sample of the conditionally analytical state rather than the distribution itself.

Fields

• x::XT: The sampled state component
• z::ZT: The Rao-Blackwellised distribution component

ReferenceTrajectory{T0, T, VT}

A 0-indexed trajectory container where the initial state (x0) may have a different type than the subsequent states (xs). This is common in Rao-Blackwellised settings where the prior's inner distribution differs from the filtered distributions at subsequent time steps.

When T0 == T, the element type collapses to a concrete T (since Union{T,T} == T).

Fields

• x0::T0: State at time 0 (typically from the prior)
• xs::VT: States at times 1, 2, …, T (an AbstractVector{T})

Indices are 0-based: traj[0] === x0, traj[t] === xs[t] for t ≥ 1.

SRKalmanFilter()

Square-root Kalman filter for linear-Gaussian state-space models. Propagates the Cholesky factor of the covariance directly via QR factorisation, keeping the covariance positive semi-definite by construction (no repair is required).

The covariance is represented as Σ = U'U with U upper triangular; the filtering state is a SqrtGaussianState. Explicit CovarianceFactors may supply rank-deficient prior and process-noise roots. Plain covariance arrays are Cholesky-factorized and must be positive definite. Observation noise must be positive definite.

Likelihood differentiation is supported on nondegenerate QR branches. Derivatives through rank changes or zero QR pivots are not supported; singular filtering support does not imply smooth derivatives of a chosen factor representation.

SSMParameterLogDensity(prior, param_model, af, trajectory)
SSMParameterLogDensity(prior, param_model, trajectory)

Log-density for SSM parameters θ conditioned on a fixed trajectory:

log p(θ | trajectory, y) ∝ log p(θ) + log p(trajectory, y | θ)

Implements the LogDensityProblems interface.

Fields

• prior: Multivariate prior distribution on the parameter vector θ
• param_model: A ParameterisedSSM mapping θ to an SSM
• af: Inner analytical filter for HierarchicalSSM (e.g., KalmanFilter()), or nothing for regular SSMs
• trajectory: Current reference trajectory (ReferenceTrajectory indexed from 0)

SSMTrajectory(model, [analytical_filter,] observations)

Turing boundary for a fixed-dimensional continuous trajectory. Its logpdf is the joint trajectory/observation density, not a normalised trajectory prior. Include it exactly once in the model; do not separately add the observations' likelihood. With an analytical filter the trajectory stores outer states only and integrates out the inner states.

The trajectory is flattened for Turing storage. Its state dimensions must be constant over time and parameters. Discrete and variable-dimensional trajectories are not supported by this continuous-vector boundary. Standalone CSMC does not have this restriction.

SerialExecution()

Evaluate particles in order on the calling task, drawing directly from the filter's random number generator. This is the default.

SqrtBackwardInformationPredictor()

Exact Gaussian backward likelihoods evaluated using QR and whitened residual factors. Unlike information-matrix subtraction, this preserves weak and strong information without cancellation. Explicit CovarianceFactors support singular process noise. No covariance or precision repair is applied. Backward messages are used for trajectory updates and are not differentiated by the HMC parameter update.

SqrtGaussianState(μ, U)

A Gaussian in square-root form, where the covariance is Σ = U'U with U upper triangular. Used by the square-root Kalman filter to propagate the covariance factor directly.

SqrtInformationLikelihood(B, r, logscale)

Gaussian likelihood exp(logscale - sum(abs2, B*z-r)/2). B may be rank deficient; this is a likelihood, not a normalized distribution on z. Square-root backward recursions retain its normalization constant, including constants from QR compression.

ThreadedExecution(; blocksize=32, ntasks=nothing)

Evaluate particles concurrently on up to ntasks tasks (Threads.nthreads() when nothing).

Particles are evaluated in consecutive blocks of blocksize. Each population draw takes one value from the filter's random number generator and derives an independent generator for every block from it. Results therefore depend on the generator state, the number of particles and blocksize, but not on ntasks, the number of threads or scheduling. They differ from SerialExecution, which draws from the filter's generator directly. Deterministic particle computations, such as observation updates and ancestor weights, match serial evaluation exactly. Resampling and weight normalisation remain serial.

blocksize is also the smallest unit of work given to a task: prefer larger blocks when each particle is cheap and smaller blocks for expensive particles, such as Rao–Blackwellised ones.

Model code must be thread-safe: draw randomness only from the generator it is given and do not mutate shared state. Every particle must have the same concrete type. Reverse-mode automatic differentiation through threaded evaluation is not supported.

TimeVaryingDynamics(f)

Wrap a conditioning callable f((; t)) -> dynamics so that it occupies a dynamics slot as a dispatchable member of the process hierarchy.

TimeVaryingObservation(f)

Wrap a conditioning callable f((; t)) -> observation so that it occupies an observation slot as a dispatchable member of the process hierarchy.

Pairs a component with its per-field activity flags (a type-domain Boolean tuple).

StateSpaceModel(outer_prior, outer_dyn, inner_prior, inner_dyn, obs)

Shorthand for building a hierarchical model. The explicit component form using HierarchicalPrior/HierarchicalDynamics/HierarchicalObservation is the canonical construction.

_csmc_sample(rng, model, csmc, observations, ref_traj)

Run one conditional SMC sweep, returning (trajectory, log_likelihood).

ref_traj is the reference trajectory from the previous iteration (or nothing for the initial unconditional run). For RBPF, the returned trajectory contains only outer states; legacy inputs containing RBState objects are accepted and stripped. Ancestor sampling respects the underlying resampling trigger and refreshes ancestors only at resampling events. All strategies require a resampler that implements a conditional law (see supports_conditional).

_kalman_adjoints(c, Δμ, ΔΣ, δll, dyn, obs)

Cotangents of the fused Kalman step from its cache c and output cotangents. Per-field activity flags carried by the dynamics/observation components (see WithFlags) skip the adjoints of fields that do not depend on the differentiated parameters; unwrapped atoms are all-active.

Shared reverse recursion propagating the filtered-state cotangent backwards.

_correct_cholesky_sign(R)

Flip rows of an upper triangular matrix so its diagonal is positive, as required for use as a Cholesky factor (QR factorisation may return negative diagonal entries).

_reference_offset(rng, vs, ref_idx, n) -> (u, K)

Draw the reference particle's slot K and the stratum offset u belonging to it, from the index distribution shared by stratified and systematic resampling.

vs is the scaled cumulative weight vector n * F, in which slot m covers (m - 1, m] and particle n covers (vs[n - 1], vs[n]]. Drawing V uniformly on the reference particle's own interval and splitting it into integer and fractional parts therefore realises K ~ Cat((vol(I^m(ref_idx)) / (n * W[ref_idx]))_m) and u ~ Unif(I^K(ref_idx)) jointly, without forming those intervals explicitly.

add_logweight(a, b)

Combine log-weight contributions using ordinary numeric addition. An initial particle weight or an exact bootstrap correction may be an internal zero marker; in that case, return the other contribution unchanged. Use this operation in custom particle methods instead of adding directly to an initial particle's log_w field.

ancestor_weight(particle, dyn, algo, iter, ref_state; kwargs...)

Compute the full (unnormalized) log backward sampling weight for a particle.

This is a convenience function that combines the particle's filtering log-weight with the future conditional density:

\[\log \tilde{w}_{t|T}^{(i)} = \log w_t^{(i)} + \log p(x_{t+1:T}^*, y_{t+1:T} \mid x_{1:t}^{(i)}, y_{1:t})\]

Arguments

• particle::Particle: The candidate particle at time t
• dyn: The latent dynamics model
• algo: The filtering algorithm
• iter::Integer: The time step t
• ref_state: The reference trajectory state at time t+1

Returns

The log backward sampling weight (unnormalized).

See also: future_conditional_density

backward_initialise(algo, obs, ...)

Initialise the backward predictive likelihood at the final time step.

backward_predict(algo, lik, dyn)

Propagate the backward predictive likelihood through the dynamics.

backward_update(algo, lik, obs, ...)

Incorporate an observation into the backward predictive likelihood.

compute_marginal_predictive_likelihood(forward, backward)

Marginal predictive likelihood p(y_{t:T} | y_{1:t-1}) from a one-step predicted Gaussian forward and a backward predictive likelihood backward. Based on Lemma 1 of https://arxiv.org/pdf/1505.06357.

compute_marginal_predictive_likelihood(forward, backward::SqrtInformationLikelihood; include_constant=false)

Gaussian overlap, omitting the suffix-common log normalization by default for accurate relative ancestor weights. Set include_constant=true for a normalized suffix likelihood.

condition_inner(model::HierarchicalSSM, outer_trajectory::AbstractVector)

Return an ordinary StateSpaceModel describing the inner process conditional on the outer trajectory. A ReferenceTrajectory is indexed 0:T; an ordinary one-based vector stores x0 at index 1 and x_t at index t+1. At least the initial state is required.

The prior is resolved once. Dynamics and observations resolve lazily, without materialising per-time matrices, and retain the original component values (including activity flags). The trajectory is borrowed and must remain unchanged while the conditional model is used: rebuild the view when trajectory values or model parameters change. No numerical filtering state or AD workspace is cached. Differentiation through construction can include the trajectory; holding this view fixed during AD instead treats its captured inputs as fixed.

This operation supplies a model, not an algorithm capability: an appropriate analytical filter is still required, and backward prediction/ancestor sampling support is separate.

conditional_sample_ancestors(rng, resampler, weights, ref_idx) -> idxs

Sample ancestor indices from resampler's conditional law given that particle ref_idx survives resampling.

The returned indices satisfy idxs[1] == ref_idx: the reference occupies the first slot, which the particle filters rely on when they propagate the reference trajectory. Only exchangeable schemes may place the reference in a fixed slot directly (Proposition 5 of Finke et al.); others must first be made reindexable by a random cyclical shift, so the remaining indices come back cyclically rotated rather than in the scheme's own slot order. That is immaterial to the filters, which treat all non-reference slots alike.

Drawing unconditionally and overwriting one index afterwards is not a conditional law for any scheme whose ancestors are dependent or not identically distributed. It biases the CSMC kernel and worsens its autocorrelation (Section 5.3 of the reference), which is why this is a separate method rather than a post-hoc fixup.

cpu_rng(rng::AbstractRNG)

Return the CPU child of a CombinedRNG, or rng itself for an ordinary RNG.

create_homogeneous_linear_gaussian_model(μ0, Σ0, A, b, Q, H, c, R)

Convenience constructor for a time-invariant linear-Gaussian state-space model from constant parameter arrays.

execution(algo::AbstractParticleFilter) -> AbstractExecution

The AbstractExecution used to evaluate the particle population. Particle filters without an execution setting evaluate serially.

filter([rng,] model, algo, ys; ref_state=nothing)

Run a filtering algorithm over nonempty observations ys, returning (final_state, total_ll). An empty observation sequence raises an ArgumentError.

future_conditional_density(dyn, algo, iter, state, ref_state; kwargs...)

Compute the log conditional density of the future trajectory given the present state:

\[\log p(x_{t+1:T}, y_{t+1:T} \mid x_{1:t}, y_{1:t})\]

up to an additive constant that does not depend on $x_t$.

This function is the key computational primitive for backward simulation (BS) and ancestor sampling (AS) algorithms. The full backward sampling weight combines this with the filtering weight:

\[\tilde{w}_{t|T}^{(i)} \propto w_t^{(i)} \cdot p(x_{t+1:T}^*, y_{t+1:T} \mid x_{1:t}^{(i)}, y_{1:t})\]

Standard (Markov) Case

For Markovian models, the future conditional density factorizes as:

\[p(x_{t+1:T}, y_{t+1:T} \mid x_{1:t}, y_{1:t}) = f(x_{t+1} \mid x_t) \cdot p(y_{t+1:T} \mid x_{t+1:T})\]

The second factor is constant across candidate ancestors (since $x_{t+1:T}^*$ is fixed), so this function returns only the transition density:

\[\texttt{future\_conditional\_density} = \log f(x_{t+1}^* \mid x_t)\]

Rao-Blackwellised Case

For hierarchical models with Rao-Blackwellisation, the marginal outer state process is non-Markov due to the marginalized inner state. The future conditional density becomes:

\[p(u_{t+1:T}, y_{t+1:T} \mid u_{1:t}, y_{1:t}) \propto f(u_{t+1} \mid u_t) \cdot p(y_{t+1:T} \mid u_{1:T}, y_{1:t})\]

where $u$ denotes the outer (sampled) state. Unlike the Markov case, the second factor depends on the candidate ancestor through $u_{1:t}$. This is computed via the two-filter formula using the forward filtering distribution and backward predictive likelihood.

Arguments

• dyn: The latent dynamics model
• algo: The filtering algorithm
• iter::Integer: The time step t
• state: The candidate state at time t (contains filtering distribution for RB case)
• ref_state: The reference trajectory state at time t+1

Returns

The log future conditional density (up to additive constants independent of $x_t$).

Implementations

• Generic (LatentDynamics, AbstractFilter): Returns logdensity(dyn, iter, state, ref_state)
• Rao-Blackwellised (HierarchicalDynamics, RBPF): Combines outer transition density with marginal predictive likelihood. The ref_state.z must be an AbstractLikelihood (InformationLikelihood or SqrtInformationLikelihood for Gaussian inner states, DiscreteLikelihood for discrete inner states).

See also: compute_marginal_predictive_likelihood, BackwardInformationPredictor, BackwardDiscretePredictor

get_ancestry(tree::ParticleTree)

Return all leaf trajectories as a vector of ReferenceTrajectory, walking from each leaf up to the corresponding initial state.

get_ancestry(container::DenseParticleContainer, i::Integer)

Return the trajectory of particle i as a ReferenceTrajectory, walking the ancestry backwards from the final time to time 0.

gpu_rng(rng::AbstractRNG)

Return the GPU child of a CombinedRNG, or rng itself for an ordinary RNG. This accessor selects a generator; it does not validate backend compatibility.

initialise([rng,] prior, algo; ref_state=nothing)
initialise(execution, rng, prior, particle_filter; ref_state=nothing)

Construct the initial state representation from the prior for a filtering algorithm. Particle filters delegate to the execution-aware overload to choose population storage. Once initialised, subsequent operations dispatch on the state representation.

inner_dynamics(model::HierarchicalSSM, t, x_prev, x_new)
inner_dynamics(dynamics::HierarchicalDynamics, t, x_prev, x_new)
inner_dynamics(component, t, x_prev, x_new)

Resolve the inner transition at time t using the adjacent outer states. Both generative hierarchical operations and conditioned likelihoods use this component-level operation. Longer history dependence requires an augmented outer state for consistent incremental reuse in particle filtering.

inner_loglikelihood(af::AbstractFilter, model::HierarchicalSSM, xs, ys)

Evaluate log p(y₁:T | x₀:T) by conditioning the inner model and calling its ordinary marginal_loglikelihood evaluator. No parameter prior or outer trajectory density is added. The filter must support marginal likelihood evaluation for the resolved inner components. See condition_inner for trajectory indexing and lifetime requirements.

inner_observation(model::HierarchicalSSM, t, x)
inner_observation(observation::HierarchicalObservation, t, x)
inner_observation(component, t, x)

Resolve the inner observation process at time t conditioned on the current outer state. The component form uses resolve(component, (; t, x)).

inner_prior(model::HierarchicalSSM, x0)
inner_prior(prior::HierarchicalPrior, x0)
inner_prior(component, x0)

Resolve the inner prior conditioned on the initial outer state. Define conditioning through callables or resolve methods; the model-level methods are convenience delegates. The component form uses resolve(component, (; x0)); hierarchical components and models delegate to it.

kalman_step_cached(state, d, o, y)

Fused predict-and-update returning (filtered_state, ll_increment, cache). This is the differentiable primitive; the log-likelihood uses the pre-repair innovation covariance and no repair is applied inside it.

kalman_update_cached(state, o, y)

Joseph-form Kalman update returning (filtered_state, ll_increment, cache). The cache exposes the intermediates consumed by the analytic reverse pass.

log_likelihoods(state::DiscreteLikelihood)

Extract the log backward likelihoods from a DiscreteLikelihood.

marginal_loglikelihood(model, ::DiscreteFilter, ys)

Evaluate a finite-state model's deterministic forward likelihood. Emission weights are normalised in log space; impossible observations return -Inf.

marginal_loglikelihood(model, af::KalmanFilter, ys)

Marginal log-likelihood log p(y_{1:T}) of a linear-Gaussian model, computed through the fused Kalman step. Conditional inner models returned by condition_inner use this same evaluator. Observations must be one-based and match the conditional trajectory horizon. Observations must be nonempty. The likelihood total starts with the first increment and retains the scalar type determined by the Kalman calculations.

marginal_loglikelihood(model, ::SRKalmanFilter, ys)

Evaluate the deterministic square-root Kalman likelihood through the ordinary filter.

marginalise!(state::ParticleDistribution)

Compute the log-likelihood increment and normalize particle weights. This function:

1. Computes LSE of current (post-observation) log-weights
2. Calculates llincrement = LSEafter - ll_baseline
3. Normalizes weights using a maximum shift to preserve small log corrections
4. Resets ll_baseline to 0.0

The llbaseline field handles both standard particle filter and auxiliary particle filter cases through a single-scalar caching mechanism. For standard PF, llbaseline equals the LSE before adding observation weights. For APF with resampling, it includes first-stage correction terms computed during the APF resampling step.

natural_params(state::InformationLikelihood)

Extract the natural parameters (λ, Ω) from an InformationLikelihood.

Returns a tuple (λ, Ω) where λ is the information vector and Ω is the information/precision matrix.

outer_logdensity(model::HierarchicalSSM, xs::AbstractVector)

Compute the outer trajectory log-density log p(x₀) + Σₜ log p(xₜ | xₜ₋₁). A ReferenceTrajectory uses indices 0:T; a one-based vector stores x₀ at index 1. Includes all transitions represented by xs, and excludes observations and parameter priors.

predict([rng,] dyn, algo, t, state, y; ref_state=nothing)

Propagate the state representation forward through the dynamics.

resolve(component, ctx)

Resolve a model component against a conditioning context. Constant components (anything in the process hierarchy) resolve to themselves; any other object is treated as a conditioning callable and applied to ctx, a NamedTuple whose fields follow the schema fixed for the component's slot.

A conditioning callable must not subtype the process abstract types; if it does, it is treated as a constant.

rts_backward_step(filtered, d, smoothed_next, predicted=nothing)

Single Rauch-Tung-Striebel backward step. predicted is p(x_{t+1} | y_{1:t}); if omitted it is recomputed from filtered.

smooth([rng,] model, algo, ys; t_smooth=1)

Run a forward-backward smoothing pass, returning (smoothed_state, total_ll) where the smoothed state is p(x_{t_smooth} | y_{1:T}).

step([rng,] model, algo, t, state, y; ref_state=nothing)

Combined predict-and-update step, returning (new_state, ll_increment).

supports_conditional(resampler::AbstractResampler) -> Bool

Whether resampler implements a conditional resampling law, and can therefore be used within ConditionalSMC. Defaults to false; wrappers forward the trait to the scheme they delegate to.

This is unrelated to AbstractConditionalResampler, which describes resamplers that resample only when a trigger fires (such as ESSResampler).

trajectory_logdensity(model::AbstractStateSpaceModel, xs, ys)

Compute the joint log-density of sampled states and observations, including the initial state and every transition. Unlike the four-argument hierarchical method, this method integrates out no states: a hierarchical model therefore takes HierarchicalState values. Trajectory indexing follows condition_inner; observations are one-based.

trajectory_logdensity(model::HierarchicalSSM, af::AbstractFilter, xs, ys)

Compute log p(x₀:T, y₁:T) with the inner states integrated out. This is the sum of outer_logdensity and inner_loglikelihood, using the same model definition. Parameters are supplied by constructing model; their prior density and any parameter transformation Jacobian are deliberately left to the caller (for example, Turing).

Differentiating θ -> trajectory_logdensity(build(θ), af, xs, ys) keeps the outer trajectory fixed while retaining parameter dependence in both the outer and inner processes. This function does not install a particle-Gibbs sampler or manage Turing/AD caches.

two_filter_smooth(filtered, backward_lik)

Combine a forward filtered distribution p(x_t | y_{1:t}) with a backward predictive likelihood p(y_{t+1:T} | x_t) to obtain the smoothed distribution p(x_t | y_{1:T}).

update(obs, algo, t, state, y)

Incorporate observation y into the predicted state, returning (filtered_state, ll_increment).

will_resample(resampler::AbstractResampler, state::ParticleDistribution)

Determine whether a resampler will trigger resampling given the current particle state. For uncondition resamplers, always returns true. For conditional resamplers (e.g., ESSResampler), checks the resampling condition.

with_activity(model, ::Val{flags})

Stamp activity flags (dyn=..., obs=...) into a hierarchical model's inner components so the analytic reverse pass can skip the adjoints of parameter-independent fields. flags must be a type-domain constant built outside the differentiated region. A model not passed through with_activity is all-active.

AlternatingResampler

A resampler wrapper that alternates between resampling and not resampling on each step. This is useful for testing the validity of filters in both cases.

The resampler maintains internal state to track whether it should resample on the next call. By default, it resamples on the first call, then skips the second, then resamples on the third, and so on.

DiscreteGaussianObservation

Observation process that emits Gaussian observations based on discrete inner state. Each discrete state k has associated mean μ[k] and variance σ²[k].

Inner dynamics of the dummy linear Gaussian model.

Linear Gaussian dynamics conditonal on the (previous) outer state (u_t), defined by:
x_{t+1} = A x_t + b + C u_t + w_t

JointDiscreteObservation

Observation process for the joint (product) state space. Joint state index encodes (outer, inner) as: idx = (outer - 1) * K_inner + inner

MixtureObservation

Observation process for discrete-state HMMs that emits Gaussian observations with state-dependent means and unit variance. Used in discrete filter/smoother tests where each discrete state k has mean μs[k].

Optimal importance proposal for a linear Gaussian model.

Transition proposal with covariance inflation, supporting ordinary and RB states.

augmented_kf_drift_posterior(model, observations, drift_indices; kwargs...)

Run a Kalman filter on the augmented model from augment_drift_model and return posterior mean/std for unknown drift components.

central_diff(f, x; h=1e-6)

Central finite-difference gradient of f at x, used as the reference for gradient checks.

check_gradients(f, θ; rtol=1e-6, backend=AutoMooncake(; config=nothing))

Compare the reverse-mode gradient of f at θ against central finite differences. Returns a NamedTuple (; ad, fd, max_abs_error, agrees) so callers can both assert agrees and report the deviation. Requires the AD backend's package to be loaded (Mooncake for the default).

create_dummy_discrete_model(rng, K_outer, K_inner)

Create a dummy hierarchical discrete model and its equivalent joint model.

Returns (joint_model, hier_model) where:

• joint_model: Standard discrete SSM on Kouter * Kinner product space
• hier_model: HierarchicalSSM with discrete outer and inner components

The joint model encodes (outer, inner) pairs as: idx = (outer - 1) * K_inner + inner

make_pd(M)

Project a square matrix to the nearest positive-definite matrix by symmetrising and clipping negative eigenvalues to 1e-8. Returns a symmetric matrix.

with_inner_drift(model::HierarchicalSSM, b)

Return a copy of a dummy linear Gaussian HierarchicalSSM with inner drift replaced by b. The helper preserves the existing inner drift container type (e.g. Vector/SVector).