abstract type AbstractCondSampler{d,dC,T,R1,R2} <: ConditionalMapping{d,dC,T}

A conditional mapping equipped with reference maps R1 and R2 from which can be sampled and the distribution be evaluated.

Type Parameters

• d: Dimension of the joint variables.
• dC: Dimension of the conditional variable.
• T: Number type used, e.g. Float64.
• R1: Type of the random variable used for sampling.
• R2: Type of the random variable used for conditioning.

A transport map, applying transport maps to subset of variables.

The most generic implementation of CondSampler and Sampler.

ATTENTION: applies maps in reverse order, i.e., the last map is applied first. T = T1 ◦ T2 ◦ ... ◦ Tn where Ti is the i-th map in the sampler.

abstract type ConditionalMapping{d,dC,T}

Abstract type representing a conditional mapping for variables $(x, y) \in \mathcal{R}^d$ and $y \in \mathcal{R}^{dC}$. Mapping can be either used to map $(z_x, z_y)$ to $(x, y)$, or given $x$ to map $z_y$ to $y$. Such maps are triangular of the form: $\mathccal{T}(z_x, z_y) = (T_x(z_x), T_y(z_x, z_y)) = (x, y)$

• d: The dimension of the input space.
• dC: The dimension of the conditioning space.
• T: Number type used, e.g. Float64.

DownwardClosedTensorizer{d}

A downward closed set is defined by

A feature map of tensorization of polynimials, where $\sigma$ is the function that maps a multiindex to a coefficient index in the tensorization.

InverseMapping{d, dC, T, F}

A reverse mapping is a mapping that is the inverse of another mapping.

const Mapping{d,T} = ConditionalMapping{d,0,T}

A Mapping type that is an alias for ConditionalMapping with dC = 0, hence, non conditional. This type represents a mapping between two sets of measures.

A PSD model where the feature map is made out of orthonormal polynomial functions such as Chebyshev polynomials. There is an own implementation in order to provide orthonormal polynomials, optimal weighted sampling, closed form derivatives and integration, etc. For orthogonal polynomials, the package ApproxFun is used.

Map of type f1(x{1:k-1}) + g(f2(x{1:k-1})) * x_k

The projection mapping is always defined on the uniform distribution on [0,1]^d.

Defined on [0, 1]^d of type \int{0}^{xk} Φ(x{1:k-1}, z)' A Φ(x{1:k-1}, z) dz with A ⪰ 0, tr(A) = 1

abstract type SubsetMapping{d,dC,T,dsub,dCsub} <: ConditionalMapping{d,dC,T}

Abstract type for mappings that work on a subspace of dimension $dsub$ of the input space and $dCsub$ of the conditional space.

Parameters

• d: Dimension of the input space.
• dC: Dimension of the conditional space.
• T: Number type used, e.g. Float64.
• dsub: Dimension of the subspace.
• dCsub: Dimension of the conditional subspace.

TriangularMap{d, dC, T}
A generic implementation of triangular transport maps. In practice, this is
overriten by a specific implementation of a triangular transport map.
The functions defined here can be reused, by defining the `map[k]` and `jacobian[k]`

The map is defned as Q(x)_k = Q(x_k | x_{<k}) and the jacobian as ∂_k Q(x)_k = ∂_k Q(x_k | x_{<k})
where x_{<k} = [x_1, ..., x_{k-1}]

Adaptive self-reinforced maximum likelihood estimation

This function creates a Sampler by adding layers of bridging densities until the residual of the sampler is larger than the smallest residual times (1+ϵsmallest) or larger than the previous residual times (1 + ϵlast).

Arguments

• X_train::PSDDataVector{T}: Training data
• X_val::PSDDataVector{T}: Validation data
• model::PSDModelOrthonormal{dsub, T, S}: Model to be fitted
• β::T: Bridging density parameter
• reference_map::ReferenceMap{d, <:Any, T}: Reference map
• ϵ_last=1e-4: Residual threshold for previous residual
• ϵ_smallest=1e-3: Residual threshold for smallest residual
• return_smallest_residual=false: If true, after conditions is reached the last layers are removed and the sampler with the smallest residual is returned
• residual: Residual function, default is negative log likelihood
• L_max=100: Maximum number of layers
• subsample_data=false: Subsample data in case of large data
• subsample_size=2000: Size of subsample, only used if subsample_data=true
• threading=true: Enable/disable threading
• dC=0: Conditional dimension to estimate π(y | x) where y are the last dC dimensions
• dCsub=0: Conditional dimension of the subspace, at least 1 if dsub < d and dC > 0
• trace_df=nothing: DataFrame to store trace information. If not nothing, the following columns are added at every iteration:
  • :residual: Residual of the sampler
  • :layers: Number of layers
  • :time: Time to add layer

IRLS!(a::PSDModel{T}, X::PSDDataVector{T}, Y::Vector{T}, reweight::Function; λ1=1e-8, trace=false, maxit=5000, tol=1e-6, preeval=true, preevalthresh=5000, ) where {T<:Number}

Minimizes $B^* = \argmin_B L(a_B(x_1), a_B(x_2), ...) + λ_1 tr(B)$ and returns the modified PSDModel with the right matrix B.

Algorithms to create samplers

Maximum likelihood implementation for JuMP solver.

Add penalty to matrices to avoid numerical issues of small eigenvalues.

Stopping rule estimates condition so that P(|X - μ| > δ) < p where X is the sample mean and μ is the true mean

add_support(a::PSDModel{T}, X::PSDdata{T}) where {T<:Number}

Returns a PSD model with added support points, where the model still gives the same results as before (extension of the matrix initialized with zeros).

Remark: This also works for Orthogonal Mapped functions, when the domain endpoints are those of the non mapped functions and the measure is not dependent on x (Lebesgue measure).

TODO: What is for non Lebesgue measure?

Calculates M{σ(i)σ(j)} = \int measure(x) \phi{idim}(x) \phi{j_dim}(x) dx using Gauss-Legendre quadrature.

PDF p(y|x) = p(x, y) / p(x)

conditional_pullback(map, y, x)

Given $(y, x)$ calculate $z_y$ given $x$.

Arguments

• map: A ConditionalMapping object representing the conditional mapping.
• y: A PSDdata object representing the conditioning data.
• x: A PSDdata object representing the input data.

Returns

• y: Conditional pushforward $\mathcal{T}_{\mathcal{Y}}(z_y | x)$.

conditional_pushforward(sampler, z_y, x)

Given $z_y$ and $x$ compute the pushforward of $z_y$ given $x$. If sampler represents a distribution $p(x, y)$ and $z_y \sim \rho_y$, then the pushforward $y$ is distributed as $y \sim p_{Y | X = x}$.

Arguments

• sampler: A ConditionalMapping object representing the conditional mapping.
• u: A PSDdata object representing the conditioning data.
• x: A PSDdata object representing the input data.

Returns

• y: Conditional pushforward $\mathcal{T}_{\mathcal{Y}}(z_y | x)$.

optimizePSDmodel(initial::AbstractMatrix, loss::Function, X::AbstractVector; λ_1::Float64=1e-8, trace::Bool=false, maxit=5000, tol=1e-6, smooth=true, ) where {T<:Number}

Minimizes loss with the constraint of PSD and chooses the right solver depending on the model.

function integral(a::PSDModelPolynomial{d, T, S}, dim::Int; C=nothing)

Integrate the model along a given dimension. The integration constant C gives the x value of where it should start. If C is not given, it is assumed to be the beginning of the interval.

integrate(a::PSDModel{T}, p::Function, χ::Domain; quadraturemethod=gausslegendre, amountquadrature_points=10) where {T<:Number}

returns $\int_χ p(x) a(x) dx$. The idea of the implementation is from proposition 4 in [1]. The integral is approximated by a quadrature rule. The default quadrature rule is Gauss-Legendre.

[1] U. Marteau-Ferey, F. Bach, and A. Rudi, “Non-parametric Models for Non-negative Functions” url: https://arxiv.org/abs/2007.03926

Distribution p(x) = ∫ p(x, y) d y

Marginalize for RBF kernel, where the kernel is defined as k(x, y) = σ^2 * exp(-||x - y||^2 / (2l^2))

Marginalize the model along a given dimension according to the measure to which the polynomials are orthogonal to.

minimize!(a::PSDModel{T}, L::Function, X::PSDDataVector{T}; λ1=1e-8, trace=false, maxit=5000, tol=1e-6, preeval=true, preevalthresh=5000, ) where {T<:Number}

Minimizes $B^* = \argmin_B L(a_B(x_1), a_B(x_2), ...) + λ_1 tr(B)$ and returns the modified PSDModel with the right matrix B.

mul!(a::PSDModel, b::Number)

In place multiplication with a number.

For now, naive implementation of saving and loading samplers.

Conversion of Tensorizer

Φ(a::PSDModelFM, x::PSDdata{T}) where {T<:Number}

Returns the feature map of the PSD model at x.

Φ(a::PSDModelFM, x::PSDdata{T}) where {T<:Number}

Returns the feature map of the PSD model at x.

Φ(a::PSDModelKernel, x::PSDdata{T}) where {T<:Number}

Returns the feature map of the PSD model at x.