ExponentialFamilyManifolds

ExponentialFamilyManifolds.jl provides implementations of manifolds for the natural parameters of exponential family distributions, using Manifolds.jl. These manifolds are compatible with ManifoldsBase.jl, enabling optimization of the natural parameters of exponential family distributions using Manopt.jl.

The primary operation in the package is the get_natural_manifold function, which returns the appropriate manifold for the natural parameters of a given exponential family member type, its dimension and (if required) conditioner.

ExponentialFamilyManifolds.get_natural_manifold — Function

get_natural_manifold(::Type{T}, dims, conditioner = nothing)

The function returns a corresponding manifold for the natural parameters of distribution of type T. Optionally accepts the conditioner, which is set to nothing by default. Use empty tuple () for univariate distributions.

julia> using ExponentialFamily, ExponentialFamilyManifolds

julia> ExponentialFamilyManifolds.get_natural_manifold(Beta, ()) isa ExponentialFamilyManifolds.NaturalParametersManifold
true

julia> ExponentialFamilyManifolds.get_natural_manifold(MvNormalMeanCovariance, (3, )) isa ExponentialFamilyManifolds.NaturalParametersManifold
true

The get_natural_manifold function returns NaturalParametersManifold manifold, which is a wrapper around the actual manifold (the base) for the natural parameters that stores extra useful properties and provides the necessary operations for optimization with Manopt.jl.

using ExponentialFamilyManifolds, ExponentialFamily

ExponentialFamilyManifolds.get_natural_manifold(Beta, (), nothing)

ExponentialFamilyManifolds.NaturalParametersManifold{ℝ, Distributions.Beta, Tuple{}, ManifoldsBase.ProductManifold{ℝ, Tuple{ExponentialFamilyManifolds.ShiftedPositiveNumbers{Static.StaticInt{-1}}, ExponentialFamilyManifolds.ShiftedPositiveNumbers{Static.StaticInt{-1}}}}, Nothing}((), ProductManifold(ShiftedPositiveNumbers(static(-1)), ShiftedPositiveNumbers(static(-1))), nothing)

ExponentialFamilyManifolds.NaturalParametersManifold — Type

NaturalParametersManifold(::Type{T}, dims, base, conditioner)

The manifold for the natural parameters of the distribution of type T with dimensions dims. An internal structure, use get_natural_manifold to create an instance of a manifold for the natural parameters of distribution of type T.

Its not advised to use the NaturalParametersManifold to create a manifold, but instead use the get_natural_manifold function.

Natural manifold base

The get_natural_manifold_base function returns the base manifold without the wrapper.

ExponentialFamilyManifolds.get_natural_manifold_base(Beta, (), nothing)

ProductManifold with 2 submanifolds:
 ShiftedPositiveNumbers(static(-1))
 ShiftedPositiveNumbers(static(-1))

The base manifold, however, does not encode the information about the conditioner, hence, it cannot be used for all exponential members. Additionally, it does not encode the type of the underlying exponential family members. For instance, the LogNormal and the univariate Normal distribution share the same base manifold, yet they represent different members of the exponential family of distributions.

ExponentialFamilyManifolds.get_natural_manifold_base — Function

get_natural_manifold_base(M::NaturalParametersManifold)
get_natural_manifold_base(::Type{T}, dims, conditioner = nothing)

Returns base manifold for the distribution of type T of dimension dims. Optionally accepts the conditioner, which is set to nothing by default.

get_natural_manifold_base(::Type{Bernoulli}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Bernoulli distribution.

get_natural_manifold_base(::Type{Beta}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Beta distribution.

get_natural_manifold_base(::Type{Binomial}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Binomial distribution.

get_natural_manifold_base(::Type{Chisq}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Chisq distribution.

get_natural_manifold_base(::Type{Categorical}, dims::Tuple{Int}, conditioner=nothing)

Get the natural manifold base for the Categorical distribution.

get_natural_manifold_base(::Type{Dirichlet}, dims::Tuple{Int}, conditioner=nothing)

Get the natural manifold base for the Dirichlet distribution.

get_natural_manifold_base(::Type{Exponential}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Exponential distribution.

get_natural_manifold_base(::Type{Gamma}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Gamma distribution.

get_natural_manifold_base(::Type{ExponentialFamily.GammaInverse}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the ExponentialFamily.GammaInverse distribution.

get_natural_manifold_base(::Type{Geometric}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Geometric distribution.

get_natural_manifold_base(::Type{Laplace}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Laplace distribution.

get_natural_manifold_base(::Type{LogNormal}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the LogNormal distribution.

get_natural_manifold_base(::Type{NormalMeanVariance}, ::Tuple{}, conditioner = nothing)

Get the natural manifold base for the NormalMeanVariance distribution.

get_natural_manifold_base(::Type{MvNormalMeanCovariance}, dims::Tuple{Int}, conditioner = nothing)

Get the natural manifold base for the MvNormalMeanCovariance distribution.

get_natural_manifold_base(::Type{NegativeBinomial}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the NegativeBinomial distribution.

get_natural_manifold_base(::Type{Pareto}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Pareto distribution.

get_natural_manifold_base(::Type{Poisson}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Poisson distribution.

get_natural_manifold_base(::Type{Rayleigh}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Rayleigh distribution.

get_natural_manifold_base(::Type{Weibull}, ::Tuple{}, conditioner=nothing)

Get the natural manifold base for the Weibull distribution.

get_natural_manifold_base(::Type{WishartFast}, ::Tuple{Int}, conditioner=nothing)

Get the natural manifold base for the WishartFast distribution.

Product manifolds

Some base manifolds are known as Product Manifolds, which consist of several manifolds combined together. For example, the natural parameters of a multivariate Normal distribution form a product of a Euclidean vector manifold and a symmetric negative definite matrix manifold. The partition_point function takes a plain vector of natural parameters and (typically, but not always) returns a partitioned array for each submanifold in the form of an ArrayPartition from RecursiveArrayTools.jl.

M = ExponentialFamilyManifolds.get_natural_manifold(Beta, (), nothing)
p = ExponentialFamilyManifolds.partition_point(M, [ 1.0, 2.0 ])

([1.0], [2.0])

typeof(p)

RecursiveArrayTools.ArrayPartition{Float64, Tuple{SubArray{Float64, 1, Vector{Float64}, Tuple{UnitRange{Int64}}, true}, SubArray{Float64, 1, Vector{Float64}, Tuple{UnitRange{Int64}}, true}}}

The partitioned point functions as a regular vector but encodes the structure of the product manifold, allowing differentiation between the submanifolds within the vector.

ExponentialFamilyManifolds.ManifoldsBase.submanifold_component(p, 1)

1-element view(::Vector{Float64}, 1:1) with eltype Float64:
 1.0

ExponentialFamilyManifolds.ManifoldsBase.submanifold_component(p, 2)

1-element view(::Vector{Float64}, 2:2) with eltype Float64:
 2.0

ExponentialFamilyManifolds.partition_point — Function

partition_point(M::NaturalParametersManifold, p)
partition_point(::Type{T}, dims, point, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold M of type T.

partition_point(::Type{Bernoulli}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Bernoulli.

partition_point(::Type{Beta}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Beta.

partition_point(::Type{Binomial}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Binomial.

partition_point(::Type{Chisq}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Chisq.

partition_point(::Type{Categorical}, dims::Tuple{Int}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Categorical.

partition_point(::Type{Dirichlet}, dims::Tuple{Int}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Dirichlet.

partition_point(::Type{Exponential}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Exponential.

partition_point(::Type{Gamma}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Gamma.

partition_point(::Type{ExponentialFamily.GammaInverse}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type ExponentialFamily.GammaInverse.

partition_point(::Type{Geometric}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Geometric.

partition_point(::Type{Laplace}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Laplace.

partition_point(::Type{LogNormal}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type LogNormal.

partition_point(::Type{NormalMeanVariance}, ::Tuple{}, p, conditioner = nothing)

Converts the point to a compatible representation for the natural manifold of type NormalMeanVariance.

partition_point(::Type{MvNormalMeanCovariance}, dims::Tuple{Int}, p, conditioner = nothing)

Converts the point to a compatible representation for the natural manifold of type MvNormalMeanCovariance.

partition_point(::Type{NegativeBinomial}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type NegativeBinomial.

partition_point(::Type{Pareto}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Pareto.

partition_point(::Type{Poisson}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Poisson.

partition_point(::Type{Rayleigh}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Laplace.

partition_point(::Type{Weibull}, ::Tuple{}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type Weibull.

partition_point(::Type{WishartFast}, ::Tuple{Int}, p, conditioner=nothing)

Converts the point to a compatible representation for the natural manifold of type WishartFast.

Custom generic manifolds

ExponentialFamilyManifolds.jl introduces additional manifolds not included in Manifolds.jl. This is crucial because certain exponential family distributions have natural parameters that require specific manifolds, such as negative definite matrices for the multivariate Gaussian distribution. These manifolds do not implement every operation defined in ManifoldsBase.jl, but they do provide the essential operations needed for optimization with Manopt.jl.

ExponentialFamilyManifolds.ShiftedPositiveNumbers — Type

ShiftedPositiveNumbers(shift)

A manifold representing the positive numbers shifted by shift. The points on this manifold are 1-dimensional vectors with a single element.

ExponentialFamilyManifolds.ShiftedNegativeNumbers — Type

ShiftedNegativeNumbers(shift)

A manifold representing the negative numbers shifted by shift. The points on this manifold are 1-dimensional vectors with a single element.

ExponentialFamilyManifolds.SymmetricNegativeDefinite — Type

SymmetricNegativeDefinite(k)

This manifold represents the set of negative definite matrices of size k × k. Similar to SymmetricPositiveDefinite from Manifolds.jl with the exception that the matrices are negative definite.

ExponentialFamilyManifolds.SinglePointManifold — Type

SinglePointManifold(point)

This manifold represents a set from one point.

Optimization example

Suppose, we have a set of samples from a certain exponential family distribution and we want to estimate the natural parameters of the distribution using the Manopt.jl package.

using ExponentialFamily, Distributions, Plots, StableRNGs

rng  = StableRNG(42)
dist = Beta(24, 6)
data = rand(rng, dist, 200)

histogram(data, xlim = (0, 1), label = "data", normalize=:pdf)

Example block output

using Manopt, ForwardDiff, ExponentialFamilyManifolds

# cost function
function f(M, p)
    ef = convert(ExponentialFamilyDistribution, M, p)
    return -sum((d) -> logpdf(ef, d), data)
end

# gradient function
function g(M, p)
    return ForwardDiff.gradient((p) -> f(M, p), p)
end

M = ExponentialFamilyManifolds.get_natural_manifold(Beta, ())
p = rand(rng, M)
q = gradient_descent(M, f, g, p)

q_ef = convert(ExponentialFamilyDistribution, M, q)
q_η  = getnaturalparameters(q_ef)

([20.396819066270247], [4.357726910129818])

Note that we performed the optimization in the natural parameters space, we can use ExponentialFamily.jl API to convert the vector fo natural parameters to the corresponding mean parameter space:

map(NaturalParametersSpace() => MeanParametersSpace(), Beta, q_η)

2-element Vector{Float64}:
 21.396819066270247
  5.357726910129818

As we can see the result is quite close to the actual distribution, which was used to generate the test data:

params(MeanParametersSpace(), dist)

(24.0, 6.0)

Let's also check the result, by plotting the estimated distribution on top of the data.

histogram(data, xlim = (0, 1), label = "data", normalize=:pdf, fillalpha = 0.3)
plot!(0.0:0.01:1.0, (x) -> pdf(dist, x), label = "actual", fill = 0, fillalpha = 0.2)
plot!(0.0:0.01:1.0, (x) -> pdf(q_ef, x), label = "estimated", fill = 0, fillalpha = 0.5)

Example block output

The difference in KL is quite small as well:

kldivergence(convert(Distribution, q_ef), dist)

0.0035633924185827226

Helpers

ExponentialFamilyManifolds.Negated — Type

Negated(m)

Lazily negates the matrix m, without creating a new matrix. Works by redefining the getindex.

julia> using ExponentialFamilyManifolds

julia> m = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> N = ExponentialFamilyManifolds.Negated(m)
2×2 ExponentialFamilyManifolds.Negated{Int64, Matrix{Int64}}:
 -1  -2
 -3  -4

julia> N[1, 2]
-2

Index

ExponentialFamilyManifolds.NaturalParametersManifold
ExponentialFamilyManifolds.Negated
ExponentialFamilyManifolds.ShiftedNegativeNumbers
ExponentialFamilyManifolds.ShiftedPositiveNumbers
ExponentialFamilyManifolds.SinglePointManifold
ExponentialFamilyManifolds.SymmetricNegativeDefinite
ExponentialFamilyManifolds.get_natural_manifold
ExponentialFamilyManifolds.get_natural_manifold_base
ExponentialFamilyManifolds.partition_point