Multivariate Distributions

AdditionalDistributions.jl provides multivariate Gaussian and Student-t distributions with support for rectangular cumulative probabilities.

A rectangular CDF is a probability of the form

The package currently focuses on:

• MvGaussian, which participates directly in the Distributions.AbstractMvNormal interface;

• MvTStudent, a typed wrapper around Distributions.AbstractMvTDist;

• cdf_result support for both the package wrappers and native Distributions.MvNormal / Distributions.MvTDist objects.

AdditionalDistributions.MvGaussian Type

MvGaussian(μ::AbstractVector, Σ::AbstractMatrix)

A multivariate Gaussian distribution backed by Distributions.AbstractMvNormal, with rectangular cumulative probabilities evaluated by the custom Genz-style randomized rank-1 lattice quasi-Monte Carlo backend in AdditionalDistributions.jl.

MvGaussian participates in the Distributions.jl AbstractMvNormal interface, so standard functionality such as pdf, logpdf, rand, insupport, mode, and entropy is inherited from Distributions.jl through a small set of delegated primitives.

MvGaussian(Σ)        # zero-mean version
MvGaussian(μ, Σ)     # with explicit mean and covariance

params(d)            # (μ, Σ)
cdf(d, x)            # P(X₁ ≤ x₁, ..., Xₚ ≤ xₚ)
cdf(d, a, b)         # P(a ≤ X ≤ b)
cdf_result(d, a, b)  # structured numerical result

source

AdditionalDistributions.MvTStudent Type

MvTStudent(ν::Real, μ::AbstractVector, Σ::AbstractMatrix)

A multivariate Student's t distribution backed by a native Distributions.AbstractMvTDist, with rectangular cumulative probabilities evaluated by the custom Genz-style randomized rank-1 lattice QMC backend in AdditionalDistributions.jl.

Σ is the scale/scatter matrix used by Distributions.MvTDist; it is not generally the covariance matrix.

MvTStudent(ν, Σ)
MvTStudent(ν, μ, Σ)

params(d)            # (ν, μ, Σ)
cdf(d, x)            # P(X₁ ≤ x₁, ..., Xₚ ≤ xₚ)
cdf(d, a, b)         # P(a ≤ X ≤ b)
cdf_result(d, a, b)  # structured numerical result

source

AdditionalDistributions.CDFResult Type

CDFResult(value, error, inform, neval, algorithm)

Structured result returned by cdf_result for multivariate rectangular probabilities.

Fields

• value: estimated probability.

• error: estimated absolute integration error.

• inform: convergence/status code.

• neval: requested integration budget.

• algorithm: integration algorithm identifier.

inform codes

• 0: estimated error is within tolerance.

• 1: estimated error is above tolerance for the current budget.

• 2: invalid dimension.

• 3: matrix appears not positive semidefinite during preparation.

CDFResult can be destructured as (value, error, inform) for compatibility with the legacy full=true tuple output.

source

AdditionalDistributions.cdf_result Function

cdf_result(d::Distributions.AbstractMvNormal, a, b;
           m=1000*length(a), abseps=1e-6, releps=1e-6,
           pivot=true, rng=Random.default_rng(),
           antithetic=false, batchsize=0, nshifts=10)

Estimate the rectangular probability P(a ≤ X ≤ b) for any Distributions.AbstractMvNormal and return a CDFResult.

This method is owned by AdditionalDistributions.jl, so it can legally extend the CDF functionality of native Distributions.MvNormal objects without extending Distributions.cdf(::MvNormal, ...).

source

cdf_result(d::Distributions.AbstractMvTDist, a, b;
           m=max(100_000, 10_000*length(a)),
           abseps=1e-6, releps=1e-6, pivot=true,
           rng=Random.default_rng(), antithetic=false,
           batchsize=0, nshifts=8)

Estimate P(a ≤ X ≤ b) for a native Distributions.AbstractMvTDist and return a CDFResult.

The calculation uses the t distribution's location μ and scale/scatter matrix Σ, not its statistical mean and covariance.

source

AdditionalDistributions.mvtcdf Function

mvtcdf(Σ, a, b; ν=0, δ=zeros, maxpts=1000n, abseps=1e-6,
       releps=1e-6, assume_correlation=false, pivot=true,
       antithetic=false, rng=Random.default_rng(), batchsize=0,
       nshifts=nothing)

Rectangular probability for multivariate Gaussian (ν <= 0) and multivariate Student t (ν > 0) distributions using MVSORT plus randomized randomized rank-1 lattice quasi-Monte Carlo with cached CBC lattice construction.

When nshifts is not specified, the core uses 10 randomized shifts for the Gaussian case and 8 randomized shifts for the Student t case.

Returns a named tuple (value, error, inform).

inform codes:

• 0: requested tolerance reached according to the internal error estimate.

• 1: tolerance not reached with maxpts.

• 2: invalid dimension.

• 3: covariance/correlation matrix appears non positive semidefinite.

source

Basic usage

Gaussian rectangular probabilities

using AdditionalDistributions
using LinearAlgebra
using Random

d = 5
Σ = fill(0.5, d, d)
Σ[diagind(Σ)] .= 1.0

mvnormal = MvGaussian(zeros(d), Σ)
lower = fill(-1.0, d)
upper = fill(1.0, d)

res = cdf_result(mvnormal, lower, upper;
    m = 100_000,
    rng = MersenneTwister(1234),
)

res.value
res.error
res.inform

cdf(mvnormal, lower, upper) returns only the probability estimate. Use cdf_result for diagnostics.

Student-t rectangular probabilities

ν = 4.0
mvt = MvTStudent(ν, zeros(d), Σ)

res = cdf_result(mvt, lower, upper;
    m = 1_000_000,
    abseps = 1e-8,
    releps = 1e-8,
    nshifts = 8,
    rng = MersenneTwister(1234),
)

The Student-t CDF is usually harder than the Gaussian CDF because it uses the normal-scale-mixture representation and an additional chi-square coordinate. Small degrees of freedom, high dimension, strong dependence, and tail rectangles may require larger m.

Native Distributions.jl interoperability

The same rectangular-probability backend is available directly for native Distributions.jl multivariate distributions:

using Distributions

dn = MvNormal(zeros(d), Σ)
rn = cdf_result(dn, lower, upper;
    m = 100_000,
    rng = MersenneTwister(1234),
)

dt = MvTDist(ν, zeros(d), Σ)
rt = cdf_result(dt, lower, upper;
    m = 1_000_000,
    nshifts = 8,
    rng = MersenneTwister(1234),
)

This is intentionally provided through cdf_result, a function owned by AdditionalDistributions.jl. The package does not add Distributions.cdf(::MvNormal, ...) or Distributions.cdf(::MvTDist, ...) methods, avoiding type piracy.

For package-owned wrapper types, both one-sided and rectangular scalar CDFs are available:

cdf(mvnormal, upper)         # P(Xᵢ ≤ upperᵢ for all i)
cdf(mvnormal, lower, upper)  # P(lowerᵢ ≤ Xᵢ ≤ upperᵢ for all i)

cdf(mvt, upper)
cdf(mvt, lower, upper)

API summary

cdf(dist, lower, upper; kwargs...)

returns only the probability estimate.

cdf_result(dist, lower, upper; kwargs...)

accepts MvGaussian, MvTStudent, native MvNormal, and native MvTDist objects and returns a structured result:

res.value
res.error
res.inform
res.neval
res.algorithm

Keywords

Keyword	Meaning
m	Integration budget.
abseps	Absolute error tolerance.
releps	Relative error tolerance.
rng	Random number generator for randomized shifts.
nshifts	Number of randomized QMC shifts. Defaults: 10 for Gaussian, 8 for Student-t.
batchsize	Internal batch size. Automatic settings are usually appropriate.
pivot	Enable or disable MVSORT reordering.
antithetic	Optional antithetic reflection when available.

Typical advanced usage:

res = cdf_result(mvt, lower, upper;
    m = 1_000_000,
    abseps = 1e-8,
    releps = 1e-8,
    nshifts = 8,
    rng = MersenneTwister(1234),
)

inform codes

Code	Meaning
0	Estimated error reached the requested tolerance.
1	Estimated error is above tolerance for the current budget.
2	Invalid dimension or integration setup.
3	Matrix appears not positive semidefinite.

inform = 1 is common for difficult high-dimensional or heavy-tailed cases. It means the estimated integration error did not reach the requested tolerance with the current budget. Increase m, adjust tolerances, or compare against an external reference if the result seems suspicious.

Algorithm summary

Gaussian path

The Gaussian rectangular CDF uses: 2. MVSORT variable reordering;

3. Genz-style conditional transformation;

4. a cached CBC rank-1 lattice with randomized shifts and tent-transformed coordinates;

5. batched evaluation to reduce per-point overhead.

Diagonal Gaussian covariance/correlation matrices are handled by an exact product shortcut.

The default number of randomized shifts is:

nshifts = 10

For the default floating-point path, the Gaussian QMC dimension is d - 1. A prime-length CBC rank-1 lattice is selected from the available per-shift budget and its generating vector is cached. Very small budgets and non-floating scalar types retain a Richtmyer fallback.

Student-t path

The Student-t implementation uses the scale-mixture representation. If

then

has a multivariate Student-t distribution.

The algorithm uses the same conditional Gaussian core together with one radial chi-square coordinate. Both the radial coordinate and the conditional Gaussian coordinates use the tent transformation.

The default settings for MvTStudent are:

m = max(100_000, 10_000*d)
nshifts = 8

The Student-t QMC dimension is d: one radial coordinate plus d - 1 conditional Gaussian coordinates.

Important: a diagonal Student-t scale matrix does not imply independent components. The components share a common radial scale, so diagonal Student-t rectangular probabilities are not computed as products of univariate Student-t probabilities.

Because the CBC lattice length is prime-rounded, the actual number of lattice evaluations can be slightly smaller than the nominal budget. CDFResult.neval retains the requested-budget convention for API compatibility.

Reproducibility

Randomized QMC methods depend on random shifts. Use a fixed RNG to obtain reproducible output:

rng = MersenneTwister(1234)
res = cdf_result(mvnormal, lower, upper; rng=rng)

When comparing with another implementation, always report:

• m;

• nshifts;

• abseps;

• releps;

• random seed;

• full CDFResult;

• distribution parameters;

• lower and upper bounds.

Reference comparisons

Gaussian reference tests include Genz-style cases and comparisons with MvNormalCDF.jl.

Student-t reference checks can be performed with R's mvtnorm::pmvt using GenzBretz.

See Benchmarks and Accuracy and Reproducibility for more details.