Multinomial symmetric positive definite matrices

MultinomialSymmetricPositiveDefinite <: AbstractMultinomialDoublyStochastic

The symmetric positive definite multinomial matrices manifold consists of all symmetric $n×n$ matrices with positive eigenvalues, and positive entries such that each column sums to one, i.e.

\[\begin{aligned} \mathcal{SP}^+(n) \coloneqq \bigl\{ p ∈ ℝ^{n×n}\ \big|\ &p_{i,j} > 0 \text{ for all } i=1,…,n, j=1,…,n,\\ & p^\mathrm{T} = p,\\ & p\mathbf{1}_n = \mathbf{1}_n,\\ & a^\mathrm{T}pa > 0 \text{ for all } a ∈ ℝ^{n}\backslash\{\mathbf{0}_n\} \bigr\}, \end{aligned}\]

where $\mathbf{1}_n$ and $\mathbf{0}_n$ are the vectors of length $n$ containing ones and zeros, respectively. More details about this manifold can be found in [DH19].

Constructor

MultinomialSymmetricPositiveDefinite(n::Int; parameter::Symbol=:type)

Generate the manifold of matrices $\mathbb R^{n×n}$ that are symmetric, positive definite, and doubly stochastic.

is_flat(M::MultinomialSymmetricPositiveDefinite)

Return whether the MultinomialSymmetricPositiveDefinite M is flat, which is the case if and only if its dimension is one.

manifold_dimension(M::MultinomialSymmetricPositiveDefinite)

Return the dimension of the MultinomialSymmetricPositiveDefinite manifold,

\[\operatorname{dim}_{\mathcal{SP}^+(n)} = \frac{n(n-1)}{2},\]

the dimension of the polytope of symmetric doubly stochastic matrices, see Section 1 of [Dav15], of which the manifold is an open subset, see pp. 36 and 41 of [Dou20].

project(M::MultinomialSymmetricPositiveDefinite, p, X)

Project X onto the tangent space at p on the MultinomialSymmetricPositiveDefinite M. The manifold is an open subset of the MultinomialSymmetric matrices and has their geometry, see Section 3.3, p. 39 of [Dou20], so with the symmetric part $X_{\mathrm{s}} = \frac{1}{2}(X+X^{\mathrm{T}})$ the projection reads

\[ \operatorname{proj}_p(X) = X_{\mathrm{s}} - (α\mathbf{1}_n^{\mathrm{T}} + \mathbf{1}_n α^{\mathrm{T}}) ⊙ p, \qquad (I_n+p)α = X_{\mathrm{s}}\mathbf{1}_n,\]

where $⊙$ denotes the elementwise product and $\mathbf{1}_n$ the vector of $n$ ones.

retract(M::MultinomialSymmetricPositiveDefinite, p, X, ::ProjectionRetraction)

Compute the projection retraction of the MultinomialSymmetric matrices. It projects $p⊙\exp(X⨸p)$ onto the doubly stochastic matrices, where $⊙,⨸$ are elementwise multiplication and division and $\exp$ is the elementwise exponential. The result is positive definite only for tangent vectors X small enough, see Section 3.3, p. 39 of [Dou20].

Random.rand!(
    rng::AbstractRNG, M::MultinomialSymmetricPositiveDefinite, p::AbstractMatrix;
    vector_at = nothing, kwargs...,
)

Generate a random point on MultinomialSymmetricPositiveDefinite manifold. If vector_at is given, a tangent vector is drawn with rand! of the embedding MultinomialSymmetric instead, and kwargs... are passed on. The steps are as follows:

1. Generate a random totally positive matrix a. Construct a vector L of n random positive increasing real numbers. b. Construct the Vandermonde matrix V based on the sequence L. c. Perform LU factorization of V in such way that both L and U components have positive elements. d. Convert the LU factorization into LDU factorization by taking the diagonal of U and dividing U by it, V=LDU. e. Construct a new matrix R = UDL which is totally positive.
2. Project the totally positive matrix R onto the manifold of MultinomialDoubleStochastic matrices.
3. Symmetrize the projected matrix and return the result. If the projection did not converge or the result is not positive definite, start over from step 1.

This method roughly follows the procedure described in https://math.stackexchange.com/questions/2773460/how-to-generate-a-totally-positive-matrix-randomly-using-software-like-maple