Skip to content

Marginalization

After computing a Gaussian approximation to the posterior, INLA provides several methods to marginalize over the latent field variables to obtain univariate marginal distributions. These marginals are essential for inference about individual components of the latent field.

Overview

The marginalize function computes marginal approximations for specified latent variables using different approximation methods:

  • Gaussian Marginalization: Direct marginalization of the Gaussian approximation (fastest, but ignores non-Gaussian structure)

  • Laplace Marginalization: Spline-corrected Gaussian approximation (more accurate for non-Gaussian likelihoods)

Basic Usage

julia
using Latte
using GaussianMarkovRandomFields
using Distributions

# Set up INLA model components
prior_gmrf = GMRF(μ_prior, Q_prior)
obs_model = ExponentialFamily(Bernoulli)  # Non-Gaussian likelihood
θ = NamedTuple()  # Hyperparameters
y = [1, 0, 1, 0]  # Observed data

# Compute Gaussian approximation
obs_lik = obs_model(y; θ...)
ga = gaussian_approximation(prior_gmrf, obs_lik)

# Marginalize selected variables
result = marginalize(ga, obs_lik, 0.0, LaplaceMarginal(), [1, 3]; 
                    prior_gmrf=prior_gmrf)

# Access marginal distributions
marginal_1 = result.marginals[1]  # Marginal for variable 1
marginal_3 = result.marginals[2]  # Marginal for variable 3

# Use standard Distributions.jl interface
μ_1 = mean(marginal_1)
σ_1 = std(marginal_1)
p_1 = pdf(marginal_1, 0.5)
samples = rand(marginal_1, 1000)

Marginalization Methods

Gaussian Marginalization

Latte.GaussianMarginal Type
julia
GaussianMarginal <: MarginalApproximation

Gaussian marginalization: directly marginalize the Gaussian approximation π̃_G. This is the fastest method but ignores non-Gaussian structure.

source

When to use:

  • Gaussian or nearly-Gaussian likelihoods

  • When speed is critical

  • As a baseline comparison

Example:

julia
# Fast Gaussian marginalization
gaussian_result = marginalize(ga, obs_lik, log_prior_θ, GaussianMarginal())

Laplace Marginalization

Latte.LaplaceMarginal Type
julia
LaplaceMarginal <: MarginalApproximation

Laplace marginalization: uses spline-corrected Gaussian approximation. Computes π̃_LA(x_i | θ, y) ≈ π̃_G(x_i | θ, y) * exp(spline(x_i)).

Fields

  • normalize_exactly::Bool: If true, use numerical integration for exact normalization; if false, use Gauss-Hermite approximation (faster)
source

When to use:

  • Non-Gaussian likelihoods (Bernoulli, Poisson, etc.)

  • When accuracy is important

  • For skewed or heavy-tailed marginals

Example:

julia
# Accurate Laplace marginalization
laplace_result = marginalize(ga, obs_lik, log_prior_θ, 
                           LaplaceMarginal(true), [1, 2, 5]; 
                           prior_gmrf=prior_gmrf)

Normalization Options:

  • LaplaceMarginal(true): Exact numerical integration (slower, more accurate)

  • LaplaceMarginal(false): Gauss-Hermite approximation (faster, default)

Main Function

Latte.marginalize Function
julia
marginalize(ga, obs_lik, log_prior_θ, method, indices=1:length(mean(ga));
            prior_gmrf=nothing, augmentation_info=nothing)

Compute marginal approximations for specified latent variables.

Arguments

  • ga: Gaussian approximation (GMRF object)

  • obs_lik: Materialized observation likelihood (contains data and hyperparameters)

  • log_prior_θ::Real: Log-density of hyperparameter prior

  • method::MarginalApproximation: Approximation method

  • indices::Vector{Int}: Variable indices to marginalize (default: all)

  • prior_gmrf: Original prior GMRF (required for Laplace methods, ignored for Gaussian)

  • augmentation_info: Pass the LGM's augmentation_info here when the caller is fitting an AugmentedLatentModel. SimplifiedLaplace uses this to apply a base-coordinate equivalence correction when computing skew for base latents: the skew is evaluated in the original (un-augmented) latent coordinates so that the augmentation introduced for fitting does not alter the reported marginal. Other strategies ignore it. nothing (default) means "treat the model as un-augmented" — appropriate for direct callers and tests that don't go through inla().

  • mean_override: When supplied (a length-length(mean(ga)) vector), VBCMarginal uses it as the corrected latent mean μ* rather than recomputing the per-θ correction; all other methods ignore it. Used by the per-θ INLA hook, which computes μ* once per grid point.

Returns

MarginalResult containing marginal distributions and computation time.

source

Result Structure

Latte.MarginalResult Type
julia
MarginalResult

Container for marginalization results.

Fields

  • indices::Vector{Int}: Indices of marginalized variables

  • marginals::Vector{<:ContinuousUnivariateDistribution}: Marginal distributions

  • method::MarginalApproximation: Approximation method used

  • computation_time::Float64: Computation time in seconds

  • kld_values::Vector{Float64}: Symmetric KLD between Gaussian and corrected marginal per variable

source

The marginal distributions returned are standard Julia Distribution objects that support the full Distributions.jl interface:

  • Gaussian marginals: Normal{Float64} distributions

  • Laplace marginals: SplineAugmentedGaussian{Float64} distributions with lazy computation

SplineAugmentedGaussian Distribution

For Laplace marginalization, the package provides a specialized distribution type:

Latte.SplineAugmentedGaussian Type
julia
SplineAugmentedGaussian{T} <: ContinuousUnivariateDistribution

A distribution representing a Gaussian base with a spline correction factor. This implementation uses on-demand computation with caching for expensive operations like moments (mean, var) and quantiles to ensure high performance in typical use cases (e.g., repeated calls to quantile or mean).

Fields (Internal)

  • base::Normal{T}: The base Gaussian distribution π̃_G.

  • spline: Interpolation object for the log-density correction.

  • normalization_constant::T: The pre-computed log of the normalization constant.

Cached Fields (Internal, Lazy-Loaded)

  • _moments: A cached tuple of (mean, var), computed on first request via Gauss-Hermite quadrature.

  • _cdf_spline: A cached interpolating spline for the CDF.

  • _quantile_spline: A cached interpolating spline for the quantile function (inverse CDF).

source

This distribution implements the full Distributions.jl interface with high-performance lazy computation:

julia
# All standard operations are supported
d = result.marginals[1]  # SplineAugmentedGaussian

# Statistical properties (computed efficiently via Gauss-Hermite quadrature)
μ = mean(d)      # Cached after first computation
σ = std(d)       # Cached after first computation
γ = skewness(d)  # Higher-order moments

# Density evaluation
p = pdf(d, x)
= logpdf(d, x)

# Cumulative distribution (computed via cached interpolation)
F = cdf(d, x)        # Cached spline after first computation
x_p = quantile(d, p) # Inverse CDF via cached spline

# Random sampling (efficient inverse transform method)
samples = rand(d, 1000)  # Uses cached quantile function

Performance Notes

  • First call overhead: Laplace methods compute expensive corrections on first use

  • Subsequent calls: Very fast due to caching (splines, moments, etc.)

  • Memory efficient: Only computes what's requested (moments OR quantiles)

  • Batch processing: Marginalize multiple variables in a single call for efficiency

Mathematical Details

Gaussian Marginalization

Directly extracts marginals from the multivariate Gaussian approximation:

π̃_G(x_i | θ, y) = N(μ_i, Σ_ii)

Laplace Marginalization

Uses spline-corrected Gaussian approximation:

π̃_LA(x_i | θ, y) ≈ π̃_G(x_i | θ, y) × exp(spline(x_i))

where the spline correction accounts for non-Gaussian structure in the likelihood.

Common Patterns

Comparing Methods

julia
# Compare Gaussian vs Laplace for validation
gauss_result = marginalize(ga, obs_lik, log_prior_θ, GaussianMarginal())
laplace_result = marginalize(ga, obs_lik, log_prior_θ, LaplaceMarginal(); 
                           prior_gmrf=prior_gmrf)

# For Gaussian likelihoods, these should be nearly identical
mean_diff = abs(mean(gauss_result.marginals[1]) - mean(laplace_result.marginals[1]))

Posterior Inference

julia
# Extract credible intervals
marginal = result.marginals[1]
ci_lower = quantile(marginal, 0.025)
ci_upper = quantile(marginal, 0.975)

# Posterior probability of positive effect
prob_positive = 1 - cdf(marginal, 0.0)

Custom Variable Selection

julia
# Marginalize specific variables of interest
n = length(mean(ga))
indices = [1, div(n,2), n]  # First, middle, last variables
result = marginalize(ga, obs_lik, log_prior_θ, LaplaceMarginal(), indices; 
                    prior_gmrf=prior_gmrf)

See Also