Algebra utilities
This page documents linear-algebra building blocks used internally by ReactiveMP.jl's built-in message update rules. They are exposed publicly so that custom rules can reuse them without reimplementing common operations.
Matrix constructors
diageye — identity matrix
diageye is a convenience alias for constructing a dense identity matrix of a given size. It is equivalent to Matrix{Float64}(I, n, n) but reads more clearly in rule code:
julia> using ReactiveMP; diageye(3)
3×3 Matrix{Float64}:
1.0 0.0 0.0
0.0 1.0 0.0
0.0 0.0 1.0CompanionMatrix — AR coefficient matrix
CompanionMatrix represents the companion form of an AR(p) coefficient vector θ:
\[C(\theta) = \begin{bmatrix} \theta_1 & \theta_2 & \cdots & \theta_p \\ 1 & 0 & \cdots & 0 \\ 0 & 1 & \cdots & 0 \\ \vdots & & \ddots & \vdots \\ 0 & 0 & \cdots & 0 \end{bmatrix}\]
It is a lazy AbstractMatrix — no allocation is made until an element is accessed or a product is computed. Specialized * methods exploit the sparse structure for efficient matrix-vector and matrix-matrix products.
This matrix appears in the Autoregressive node and Continuous transition node when expressing a higher-order AR process as a first-order state-space model.
PermutationMatrix — structured permutation
PermutationMatrix represents an n×n permutation matrix as a length-n index vector rather than a dense matrix. It supports efficient multiplication with vectors and matrices via specialized mul! dispatch, and its inverse is simply its adjoint.
Permutation matrices appear in normalizing flow layers (PermutationLayer) to shuffle input dimensions between coupling layers.
StandardBasisVector — sparse one-hot vector
StandardBasisVector represents a standard Cartesian basis vector — all zeros except one element — without allocating a dense array. It supports dot products, outer products, and matrix multiplication, all using the sparse structure.
ReactiveMP.diageye — Function
diageye(::Type{T}, n::Int)An alias for the Matrix{T}(I, n, n). Returns a matrix of size n x n with ones (of type T) on the diagonal and zeros everywhere else.
diageye(n::Int)An alias for the Matrix{Float64}(I, n, n). Returns a matrix of size n x n with ones (of type Float64) on the diagonal and zeros everywhere else.
ReactiveMP.CompanionMatrix — Type
CompanionMatrixRepresents a matrix of the following structure:
θ1 θ2 θ3 ... θn-1 θn 1 0 0 ... 0 0 0 1 0 ... 0 0 . . . ... . . . . . ... . . 0 0 0 ... 0 0 0 0 0 ... 1 0
ReactiveMP.PermutationMatrix — Type
PermutationMatrix(ind::Array{T}) creates a permutation matrix with ones at coordinates (k, ind[k]) for k = 1:length(ind).
A permutation matrix represents a matrix containing only zeros and ones, which basically permutes the vector or matrix it is multiplied with. These matrices A are constrained by:
\[ A_{ij} = \{0, 1\}\\ ∑_{i} A_{ij} = 1\\ ∑_{j} A_{ij} = 1\]
An example is the 3-dimensional permutation matrix
\[A = \begin{bmatrix} 1 & 0 & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 0 \end{bmatrix}\]
ReactiveMP.StandardBasisVector — Type
StandardBasisVector{T, len, ind}(scale::T)StandardBasisVector creates a standard Cartesian basis vector of zeros of length len with a single element scale at index ind.
An example is the 3-dimensional standard basis vector for the first dimension
\[e = \begin{bmatrix} 1 \\ 0 \\ 0 \end{bmatrix}\]
Which can be constructed by calling e = StandardBasisVector(3, 1, 1)
In-place scalar and array operations
These functions return alpha * A or -A, reusing the storage of A when the type allows it (i.e. when A is a mutable Array). For immutable arrays or scalars they fall back to a regular allocation. This makes them safe to use in rule code regardless of whether the input is mutable.
ReactiveMP.mul_inplace! — Function
mul_inplace!(alpha, A)Returns alpha * A, modifying and reusing A storage if possible.
See also: negate_inplace!
ReactiveMP.negate_inplace! — Function
Trace and rank-1 update utilities
These functions implement common linear-algebra patterns that appear repeatedly in Gaussian message rules.
mul_trace — allocation-free tr(A·B)
ReactiveMP.mul_trace computes tr(A * B) directly without forming the full product matrix, saving an O(n²) allocation for square matrices.
rank1update — A + x·yᵀ via BLAS
ReactiveMP.rank1update computes A + x * y'. For BlasFloat element types it dispatches to the BLAS ger! routine, which is highly optimized for this pattern.
v_a_vT — v·a·vᵀ
ReactiveMP.v_a_vT computes v * a * v' or v₁ * a * v₂'. When a is a scalar it avoids forming a temporary matrix. Specialized methods exist for StandardBasisVector inputs that exploit the one-hot structure.
ReactiveMP.mul_trace — Function
mul_trace(A, B)Computes tr(A * B) without allocating A * B.
ReactiveMP.rank1update — Function
rank1update(A, x)
rank1update(A, x, y)Helper function for A + x * y'. Uses optimised BLAS version for AbstractFloats and fallbacks to a generic implementation in case of differentiation
ReactiveMP.v_a_vT — Function
v_a_vT(v, a)Computes vav^T efficiently.
Other utilities
ReactiveMP.GammaShapeLikelihood — Type
ν(x) ∝ exp(p*β*x - p*logГ(x)) ≡ exp(γ*x - p*logГ(x))ReactiveMP.ImportanceSamplingApproximation — Type
ImportanceSamplingApproximationThis structure stores all information needed to perform an importance sampling procedure and provides convenient functions to generate samples and weights to approximate expectations
Fields
rng: random number generator objectsnsamples: number of samples generated by default
ReactiveMP.powerset — Function
powerset(iterator)Computes the set of all possible sets of the iterator.
ReactiveMP.besselmod — Function
Modified-bessel function of second kind
mx, vx : mean and variance of the random variable x my, vy : mean and variance of the random variable y rho : correlation coefficient
ReactiveMP.isonehot — Function
isonehot(vec::AbstractVector)Checks if the given vector vec is a one-hot vector, i.e., a vector with exactly one entry approximately equal to one(eltype(vec)) and all other entries approximately equal to zero(eltype(vec)).
Returns true if vec is one-hot, otherwise returns false.