User API
Documentation for ForneyLab.jl
's user API.
If you want to know how you can extend ForneyLab.jl
(e.g. register new update rules), see Developer API.
Contents
Index
ForneyLab.@RV
ForneyLab.@composite
ForneyLab.@ensureVariables
ForneyLab.@symmetrical
ForneyLab.Addition
ForneyLab.Bernoulli
ForneyLab.Beta
ForneyLab.Categorical
ForneyLab.Clamp
ForneyLab.Contingency
ForneyLab.Dirichlet
ForneyLab.DotProduct
ForneyLab.Equality
ForneyLab.Exponential
ForneyLab.FactorGraph
ForneyLab.Gamma
ForneyLab.GaussianMeanPrecision
ForneyLab.GaussianMeanVariance
ForneyLab.GaussianMixture
ForneyLab.GaussianWeightedMeanPrecision
ForneyLab.InferenceAlgorithm
ForneyLab.LogNormal
ForneyLab.Logit
ForneyLab.MarginalTable
ForneyLab.Message
ForneyLab.Multiplication
ForneyLab.Nonlinear
ForneyLab.PointMass
ForneyLab.Poisson
ForneyLab.PosteriorFactorization
ForneyLab.ProbabilityDistribution
ForneyLab.Probit
ForneyLab.Schedule
ForneyLab.Softmax
ForneyLab.Transition
ForneyLab.Variable
ForneyLab.Wishart
ForneyLab.algorithmSourceCode
ForneyLab.cholinv
ForneyLab.constant
ForneyLab.currentGraph
ForneyLab.currentInferenceAlgorithm
ForneyLab.ensureMatrix
ForneyLab.init
ForneyLab.isApproxEqual
ForneyLab.isRoundedPosDef
ForneyLab.labsbeta
ForneyLab.labsgamma
ForneyLab.leaftypes
ForneyLab.mat
ForneyLab.messagePassingAlgorithm
ForneyLab.placeholder
ForneyLab.step!
ForneyLab.trigammaInverse
Model specification
ForneyLab.@RV
— Macro@RV
provides a convenient way to add Variable
s and FactorNode
s to the graph.
Examples:
# Automatically create new Variable x, try to assign x.id = :x if this id is available
@RV x ~ GaussianMeanVariance(constant(0.0), constant(1.0))
# Explicitly specify the id of the Variable
@RV [id=:my_y] y ~ GaussianMeanVariance(constant(0.0), constant(1.0))
# Automatically assign z.id = :z if this id is not yet taken
@RV z = x + y
# Manual assignment
@RV [id=:my_u] u = x + y
# Just create a variable
@RV x
@RV [id=:my_x] x
ForneyLab.FactorGraph
— TypeA factor graph consisting of factor nodes and edges.
ForneyLab.Variable
— TypeA Variable
encompasses one or more edges in a FactorGraph
.
ForneyLab.currentGraph
— FunctionReturn currently active FactorGraph
. Create one if there is none.
Factor nodes
ForneyLab.Addition
— TypeDescription:
An addition constraint factor node.
f(out,in1,in2) = δ(in1 + in2 - out)
Interfaces:
1. out
2. in1
3. in2
Construction:
Addition(out, in1, in2, id=:some_id)
Base.:-
— Method-(in1::Variable, in2::Variable)
A subtraction constraint based on the addition factor node.
ForneyLab.Bernoulli
— TypeDescription:
Bernoulli factor node
out ∈ {0, 1}
p ∈ [0, 1]
f(out, p) = Ber(out|p) = p^out (1 - p)^{1 - out}
Interfaces:
1. out
2. p
Construction:
Bernoulli(id=:some_id)
ForneyLab.Beta
— TypeDescription:
Beta factor node
Real scalars
a > 0
b > 0
f(out, a, b) = Beta(out|a, b) = Γ(a + b)/(Γ(a) Γ(b)) out^{a - 1} (1 - out)^{b - 1}
Interfaces:
1. out
2. a
3. b
Construction:
Beta(id=:some_id)
ForneyLab.Categorical
— TypeDescription:
Categorical factor node
The categorical node defines a one-dimensional probability
distribution over the normal basis vectors of dimension d
out ∈ {0, 1}^d where Σ_k out_k = 1
p ∈ [0, 1]^d, where Σ_k p_k = 1
f(out, p) = Cat(out | p)
= Π_i p_i^{out_i}
Interfaces:
1. out
2. p
Construction:
Categorical(id=:some_id)
ForneyLab.@ensureVariables
— Macro@ensureVariables(...)
casts all non-Variable
arguments to Variable
through constant(arg).
ForneyLab.Clamp
— TypeDescription:
A factor that clamps a variable to a constant value.
f(out) = δ(out - value)
Interfaces:
1. out
Construction:
Clamp(out, value, id=:some_id)
Clamp(value, id=:some_id)
ForneyLab.constant
— Methodconstant
creates a Variable
which is linked to a new Clamp
, and returns this variable.
y = constant(3.0, id=:y)
ForneyLab.placeholder
— Methodplaceholder(...)
creates a Clamp
node and registers this node as a data placeholder with the current graph.
# Link variable y to buffer with id :y,
# indicate that Clamp will hold Float64 values.
placeholder(y, :y, datatype=Float64)
# Link variable y to index 3 of buffer with id :y.
# Specify the data type by passing a default value for the Clamp.
placeholder(y, :y, index=3, default=0.0)
# Indicate that the Clamp will hold an array of size `dims`,
# with Float64 elements.
placeholder(X, :X, datatype=Float64, dims=(3,2))
ForneyLab.@composite
— MacroThe @composite
macro allows for defining custom (composite) nodes. Composite nodes allow for implementating of custom update rules that may be computationally more efficient or convenient. A composite node can be defined with or without an internal model. For detailed usage instructions we refer to the composite_nodes
demo.
ForneyLab.Contingency
— TypeDescription:
Contingency factor node
The contingency distribution is a multivariate generalization of
the categorical distribution. As a bivariate distribution, the
contingency distribution defines the joint probability
over two unit vectors. The parameter p encodes a contingency matrix
that specifies the probability of co-occurrence.
out1 ∈ {0, 1}^d1 where Σ_j out1_j = 1
out2 ∈ {0, 1}^d2 where Σ_k out2_k = 1
p ∈ [0, 1]^{d1 × d2}, where Σ_jk p_jk = 1
f(out1, out2, p) = Con(out1, out2 | p)
= Π_jk p_jk^{out1_j * out2_k}
A Contingency distribution over more than two variables requires
higher-order tensors as parameters; these are not implemented in ForneyLab.
Interfaces:
1. out1
2. out2
3. p
Construction:
Contingency(id=:some_id)
ForneyLab.Dirichlet
— TypeDescription:
Dirichlet factor node
Multivariate:
f(out, a) = Dir(out|a)
= Γ(Σ_i a_i)/(Π_i Γ(a_i)) Π_i out_i^{a_i}
where 'a' is a vector with every a_i > 0
Matrix variate:
f(out, a) = Π_k Dir(out|a_*k)
where 'a' represents a left-stochastic matrix with every a_jk > 0
Interfaces:
1. out
2. a
Construction:
Dirichlet(id=:some_id)
ForneyLab.DotProduct
— TypeDescription:
out = in1'*in2
in1: d-dimensional vector
in2: d-dimensional vector
out: scalar
in2
|
in1 V out
----->[⋅]----->
f(out, in1, in2) = δ(out - in1'*in2)
Interfaces:
1 i[:out], 2 i[:in1], 3 i[:in2]
Construction:
DotProduct(out, in1, in2, id=:my_node)
ForneyLab.Equality
— TypeDescription:
An equality constraint factor node
f([1],[2],[3]) = δ([1] - [2]) δ([1] - [3])
Interfaces:
1, 2, 3
Construction:
Equality(id=:some_id)
The interfaces of an Equality node have to be connected manually.
ForneyLab.Exponential
— TypeDescription:
Maps a location to a scale parameter by exponentiation
f(out,in1) = δ(out - exp(in1))
Interfaces:
1. out
2. in1
Construction:
Exponential(out, in1, id=:some_id)
ForneyLab.Gamma
— TypeDescription:
A gamma node with shape-rate parameterization:
f(out,a,b) = Gam(out|a,b) = 1/Γ(a) b^a out^{a - 1} exp(-b out)
Interfaces:
1. out
2. a (shape)
3. b (rate)
Construction:
Gamma(out, a, b, id=:some_id)
ForneyLab.GaussianMeanPrecision
— TypeDescription:
A Gaussian with mean-precision parameterization:
f(out,m,w) = 𝒩(out|m,w) = (2π)^{-D/2} |w|^{1/2} exp(-1/2 (out - m)' w (out - m))
Interfaces:
1. out
2. m (mean)
3. w (precision)
Construction:
GaussianMeanPrecision(out, m, w, id=:some_id)
ForneyLab.GaussianMeanVariance
— TypeDescription:
A Gaussian with mean-variance parameterization:
f(out,m,v) = 𝒩(out|m,v) = (2π)^{-D/2} |v|^{-1/2} exp(-1/2 (out - m)' v^{-1} (out - m))
Interfaces:
1. out
2. m (mean)
3. v (covariance)
Construction:
GaussianMeanVariance(out, m, v, id=:some_id)
ForneyLab.GaussianMixture
— TypeDescription:
A Gaussian mixture with mean-precision parameterization:
f(out, z, m1, w1, m2, w2, ...) = 𝒩(out|m1, w1)^z_1 * 𝒩(out|m2, w2)^z_2 * ...
Interfaces:
1. out
2. z (switch)
3. m1 (mean)
4. w1 (precision)
5. m2 (mean)
6. w2 (precision)
...
Construction:
GaussianMixture(out, z, m1, w1, m2, w2, ..., id=:some_id)
ForneyLab.GaussianWeightedMeanPrecision
— TypeDescription:
A Gaussian with weighted-mean-precision parameterization:
f(out,xi,w) = 𝒩(out|xi,w) = (2π)^{-D/2} |w|^{1/2} exp(-1/2 xi' w^{-1} xi) exp(-1/2 xi' w xi + out' xi)
Interfaces:
1. out
2. xi (weighted mean, w*m)
3. w (precision)
Construction:
GaussianWeightedMeanPrecision(out, xi, w, id=:some_id)
ForneyLab.LogNormal
— TypeDescription:
A log-normal node with location-scale parameterization:
f(out,m,s) = logN(out|m, s) = 1/out (2π s)^{-1/2} exp(-1/(2s) (log(out) - m)^2))
Interfaces:
1. out
2. m (location)
3. s (squared scale)
Construction:
LogNormal(out, m, s, id=:some_id)
ForneyLab.Logit
— TypeDescription:
Logit mapping between a real variable in1 ∈ R and a binary variable out ∈ {0, 1}.
f(out, in1, xi) = Ber(out | σ(in1))
>= exp(in1*out) σ(xi) exp[-(in1 + xi)/2 - λ(xi)(in1^2 - xi^2)], where
σ(x) = 1/(1 + exp(-x))
λ(x) = (σ(x) - 1/2)/(2*x)
Interfaces:
1. out (binary)
2. in1 (real)
3. xi (auxiliary variable)
Construction:
Logit(out, in1, xi)
ForneyLab.Multiplication
— TypeDescription:
For continuous random variables, the multiplication node acts
as a (matrix) multiplication constraint, with node function
f(out, in1, a) = δ(out - a*in1)
Interfaces:
1. out
2. in1
3. a
Construction:
Multiplication(out, in1, a, id=:some_id)
ForneyLab.Nonlinear
— TypeDescription:
Nonlinear node modeling a nonlinear relation. Updates for
the nonlinear node are computed through the unscented transform (by default),
importance sampling, or local linear approximation.
For more details see "On Approximate Nonlinear Gaussian Message Passing on
Factor Graphs", Petersen et al. 2018.
f(out, in1) = δ(out - g(in1))
Interfaces:
1. out
2. in1
Construction:
Nonlinear{T}(out, in1; g=g, id=:my_node)
Nonlinear{T}(out, in1; g=g, g_inv=g_inv, id=:my_node)
Nonlinear{T}(out, in1, in2, ...; g=g, id=:my_node)
Nonlinear{T}(out, in1, in2, ...; g=g, g_inv=(g_inv_in1, g_inv_in2, ...), id=:my_node)
Nonlinear{T}(out, in1, in2, ...; g=g, g_inv=(g_inv_in1, nothing, ...), id=:my_node)
where T encodes the approximation method: Unscented, Sampling, or Extended.
ForneyLab.Poisson
— TypeDescription:
Poisson factor node
Real scalars
l > 0 (rate)
f(out, l) = Poisson(out|l) = 1/(x!) * l^x * exp(-l)
Interfaces:
1. out
2. l
Construction:
Poisson(id=:some_id)
ForneyLab.Probit
— TypeDescription:
Constrains a continuous, real-valued variable in1 ∈ R with a binary (boolean) variable out ∈ {0, 1} through a probit link function.
f(out, in1) = Ber(out | Φ(in1))
Interfaces:
1. out (binary)
2. in1 (real)
Construction:
Probit(out, in1, id=:some_id)
ForneyLab.Softmax
— TypeDescription:
Softmax mapping between a real variable in1 ∈ R^d and discrete variable out ∈ {0, 1}^d.
f(out, in1, xi, a) = Cat(out_j | exp(in1_j)/Σ_k exp(in1_k)),
where log(Σ_k exp(x_k)) is upper-bounded according to (Bouchard, 2007).
Interfaces:
1. out (discrete)
2. in1 (real)
3. xi (auxiliary variable)
4. a (auxiliary variable)
Construction:
Softmax(out, in1, xi, a)
ForneyLab.Transition
— TypeDescription:
The transition node models a transition between discrete
random variables, with node function
f(out, in1, a) = Cat(out | a*in1)
Where a is a left-stochastic matrix (columns sum to one)
Interfaces:
1. out
2. in1
3. a
Construction:
Transition(out, in1, a, id=:some_id)
ForneyLab.Wishart
— TypeDescription:
A Wishart node:
f(out,v,nu) = W(out|v, nu) = B(v, nu) |out|^{(nu - D - 1)/2} exp(-1/2 tr(v^{-1} out))
Interfaces:
1. out
2. v (scale matrix)
3. nu (degrees of freedom)
Construction:
Wishart(out, v, nu, id=:some_id)
Scheduling
ForneyLab.MarginalTable
— TypeA MarginalTable
defines the update order for marginal computations.
ForneyLab.PosteriorFactorization
— TypeInitialize an empty PosteriorFactorization
for sequential construction
Initialize a PosteriorFactorization
consisting of a single PosteriorFactor
for the entire graph
Construct a PosteriorFactorization
consisting of one PosteriorFactor
for each argument
ForneyLab.InferenceAlgorithm
— TypeAn InferenceAlgorithm
specifies the computations for the quantities of interest.
ForneyLab.Schedule
— TypeA Schedule
defines the update order for message computations.
ForneyLab.currentInferenceAlgorithm
— FunctionReturn currently active InferenceAlgorithm
. Create one if there is none.
Algorithm assembly
ForneyLab.messagePassingAlgorithm
— FunctionCreate a message passing algorithm to infer marginals over a posterior distribution
Algorithm code generation
ForneyLab.algorithmSourceCode
— FunctionGenerate Julia code for message passing and optional free energy evaluation
Algorithm execution
ForneyLab.Message
— TypeEncodes a message, which is a probability distribution with a scaling factor
ForneyLab.PointMass
— TypePointMass
is an abstract type used to describe point mass distributions. It never occurs in a FactorGraph
, but it is used as a probability distribution type.
ForneyLab.ProbabilityDistribution
— TypeEncodes a probability distribution as a FactorFunction
of type family
with fixed interfaces
Helper
ForneyLab.@symmetrical
— Macro@symmetrical `function_definition`
Duplicate a method definition with the order of the first two arguments swapped. This macro is used to duplicate methods that are symmetrical in their first two input arguments, but require explicit definitions for the different argument orders. Example:
@symmetrical function prod!(x, y, z)
...
end
ForneyLab.cholinv
— MethodMatrix inversion using Cholesky decomposition
ForneyLab.ensureMatrix
— MethodCast input to a Matrix
if necessary
ForneyLab.init
— MethodHelper function to call dynamically generated init
functions
ForneyLab.isApproxEqual
— MethodCheck if arguments are approximately equal
ForneyLab.isRoundedPosDef
— MethodChecks if input matrix is positive definite. We also perform rounding in order to prevent floating point precision problems that isposdef()
` suffers from.
ForneyLab.labsbeta
— MethodWrapper for logabsbeta
function that returns first element of its output
ForneyLab.labsgamma
— MethodWrapper for logabsgamma
function that returns first element of its output
ForneyLab.leaftypes
— Methodleaftypes(datatype)
returns all subtypes of datatype
that are leafs in the type tree.
ForneyLab.mat
— MethodHelper function to construct 1x1 Matrix
ForneyLab.step!
— MethodHelper function to call dynamically generated step!
functions
ForneyLab.trigammaInverse
— MethodSolve trigamma(y) = x
for y
.
Uses Newton's method on the convex function 1/trigramma(y). Iterations converge monotonically. Based on trigammaInverse implementation in R package "limma" by Gordon Smyth: https://github.com/Bioconductor-mirror/limma/blob/master/R/fitFDist.R