Francis Bach
email: francis.bach@inria.fr

Information Theory with Kernel Methods

We consider the analysis of probability distributions through their associated covariance operators from reproducing kernel Hilbert spaces. We show that the von Neumann entropy and relative entropy of these operators are intimately related to the usual notions of Shannon entropy and relative entropy, and share many of their properties. They come together with efficient estimation algorithms from various oracles on the probability distributions. We also consider product spaces and show that for tensor product kernels, we can define notions of mutual information and joint entropies, which can then characterize independence perfectly, but only partially conditional independence. We finally show how these new notions of relative entropy lead to new upper-bounds on log partition functions, that can be used together with convex optimization within variational inference methods, providing a new family of probabilistic inference methods.

Introduction

Characterizing and studying probability distributions through moments has a long history. For distributions supported in a vector space, only an infinite number of moments can characterize the distribution. Beyond traditional polynomial moments, kernel methods based on reproducing kernel Hilbert spaces (RKHS) [START_REF] Schölkopf | Learning with Kernels[END_REF][START_REF] Berlinet | Reproducing Kernel Hilbert Spaces in Probability and Statistics[END_REF] have emerged as a natural tool for studying the interactions between moments and other properties of the underlying distributions (such as independence or conditional independence), as well as providing algorithmic tools to estimate these moments from data.

A natural way of using kernel methods is to consider a feature map ϕ : X → H from the underlying space X and a specific Hilbert space H, and to consider the mean element [START_REF] Bharath | Hilbert space embeddings and metrics on probability measures[END_REF], also referred to as the mean embedding, for a probability distribution p, defined as

µ p = X ϕ(x)dp(x), (1)
which is an element of the Hilbert space H. The key idea behind kernel methods is to study properties of H and its interaction with X only through the kernel function k : X × X → R defined as k(x, y) = ϕ(x), ϕ(y) .

Under classical universality conditions on the kernel function k [START_REF] Bharath | Hilbert space embeddings and metrics on probability measures[END_REF][START_REF] Micchelli | Universal kernels[END_REF], the mapping p → µ p is injective (and thus characterizes the distribution p), and can be estimated from independent and identically distributed (i.i.d.) samples efficiently, with convergence rates proportional to 1/ √ n where n is the number of observations.

Representing a probability distribution p through µ p ∈ H defines explicitly a metric between distributions as (p, q) → µ p -µ q (for the H-norm), which has been extensively used in machine learning and data science [START_REF] Muandet | Kernel mean embedding of distributions: A review and beyond[END_REF], in particular for measuring the dependence between two random variables, leading to algorithms for independent component analysis [START_REF] Cardoso | Dependence, correlation and Gaussianity in independent component analysis[END_REF][START_REF] Bach | Kernel independent component analysis[END_REF] and for feature selection [START_REF] Song | Feature selection via dependence maximization[END_REF], or as a model fitting criterion when estimating parameters of probability distributions [START_REF] Miko Laj Bińkowski | Demystifying MMD GANs[END_REF].

A classical drawback of the underlying geometry is the lack of straightforward connection with classical information theory tools. For example for discrete data, this leads to Euclidean norms between probability mass functions, which is rarely seen as an appropriate geometry for the simplex.

In this paper, we consider the second-order moment, which we call the covariance operator1

Σ p = X ϕ(x)ϕ(x) * dp(x),
which is an operator from H to H, defined through f, Σ p g = f, ϕ(x) g, ϕ(x) dp(x), self-adjoint and positive semi-definite. 2As we show in this paper, it shares similar nice properties with the mean element in terms of universality, estimation from finite samples, and general applicability to all sets X where positive definite kernels can be defined (e.g., structured discrete objects). Furthermore, owing to the many tools from quantum information theory, and in particular the von Neumann entropy tr Σ p log Σ p and its associated divergence (relative entropy), we define these notions and demonstrate their many nice properties and applications, as well as explicit quantitative links with regular Shannon entropies (in fact, we can see the regular notion of entropy as the limit of our kernel entropies when the bandwidth of the kernel goes to zero, and for discrete data, they are equivalent).

The paper is organized as follows:

• We review in Section 2 properties of covariance operators and reproducing kernel Hilbert spaces.

• In Section 3, we review some of the key mathematical results on quantum relative entropies, that we will leverage in further sections.

• In Section 4, we define our notions of kernel entropy and Kullback-Leibler divergence (relative entropy), and give their main properties, which are shared with the classical notions, such as convexity, or additivity for independent variables. In particular our notions of relative entropies will always be lower bounds on the Shannon relative entropy (in this paper, we use interchangeably the terms "relative entropy" and "Kullback-Leibler (KL) divergence").

• In Section 5, we show how we can leverage existing work in kernel methods to estimate these quantities from finite samples, with a convergence rate proportional to 1/ √ n from n observations, and without the need for regularization. We also propose statistically more efficient estimators if integrals of kernel functions under p can be computed.

• In Section 6, we consider product spaces X = X 1 × • • • × X d , and show that for tensor product kernels, we can define notions of mutual information and joint entropies. We can then characterize independence perfectly, but only partially conditional independence, and extend the submodularity property of the regular entropy.

• In Section 7, we show how these new notions of relative entropy lead to new upper-bounds on log partition functions, that can be used together with convex optimization within variational inference methods [START_REF] Wainwright | Graphical Models, Exponential Families, and Variational Inference[END_REF], providing a new family of probabilistic inference methods. Illustrative applications to [0, 1] and {-1, 1} d are explicitly developed.

Covariance operators

In this section, we review the notions of covariance operators, starting with reproducing kernel Hilbert spaces. For more details, see [START_REF] Berlinet | Reproducing Kernel Hilbert Spaces in Probability and Statistics[END_REF][START_REF] Schölkopf | Learning with Kernels[END_REF][START_REF] Shawe | Kernel Methods for Pattern Analysis[END_REF].

Kernels and RKHSs

In this paper we consider a compact set X and probability distributions on X (which we assume to be Borel probability measures), as well as a kernel function k : X × X → R, such that (A1) k is a continuous positive definite kernel on the compact set X, with k(x, x) = 1 for all x ∈ X.

Positive definite kernel functions are functions for which all matrices of pairwise kernel evaluations are positive semi-definite. They are known to define a Hilbert space H (with dot product •, • ) of functions f : X → R, called a reproducing kernel Hilbert space (RKHS), as well as a map ϕ : X → H such that, for all x, y ∈ X, (a) f (x) = f, ϕ(x) , (b) k(•, x) = ϕ(x), and (c) k(x, y) = ϕ(x), ϕ(y) .

A useful representation of functions in H is through linear combinations of kernel functions, that is, for any α ∈ R n and arbitrary elements x 1 , . . . , x n of X, the function f = n i=1 α i k(•, x i ) is in H and its norm is equal to n i,j=1 α i α j k(x i , x j ). It turns out that all functions in H can be generated as appropriate limits of such finite linear combinations.

Note that from k(x, x) = 1, we get ϕ(x) = 1 for all x ∈ X. Note that most of the developments in this paper do not require further structure on the set X, except when characterizing regularity of density functions is needed in Section 4.2, where will assume the existence of a distance.

Universal and characteristic kernels. We note that our covariance operator Σ p happens to be a mean embedding for the feature map x → ϕ(x)ϕ(x) * , from X to the space of operators from H to H. Equipped with the Hilbert-Schmidt norm, we obtain a Hilbert space which happens to be isomorphic 3 to the RKHS obtained from the kernel (x, y) → k(x, y) 2 . Therefore, injectivity of the map p → Σ p can be obtained from sufficient conditions for injectivity of the mean element map, a property for the kernel referred to as characteristic [START_REF] Bharath | Injective Hilbert space embeddings of probability measures[END_REF]. As shown in [START_REF] Fukumizu | Kernel dimension reduction in regression[END_REF]Prop. 5], this is equivalent for the associated RKHS to be dense in L 2 (q) for all probability measures q on X, a sufficient condition being here (because of our compactness assumption) that the RKHS is dense in the set of continuous functions equipped with the uniform norm, a property referred to as universality [START_REF] Steinwart | On the influence of the kernel on the consistency of support vector machines[END_REF]. Note that if (A1) is satisfied, then the universality of k implies the universality of k 2 .

We will see in Section 2.3 below that many kernels are indeed universal, and thus the covariance operator will be characteristic of the corresponding probability distribution.

Covariance operators

Given a probability distribution p on X, we consider the operator Σ p : H → H, defined as

Σ p = X ϕ(x)ϕ(x) * dp(x). (2)
For any f, g ∈ H, we have, by definition: g,

Σ p f = E X∼p f, ϕ(X) g, ϕ(X) = E X∼p f (X)g(X) .
Proposition 1 (Properties of covariance operators) Assume (A1).

(a) For any probability distribution p on X, the operator Σ p : H → H defined in Eq. ( 2), is self-adjoint, positive-semidefinite and has unit trace.

(b) If the kernel k 2 is universal, then the mapping p → Σ p is injective from probability distributions to self-adjoint positive-semidefinite, unit trace operators.

(c) If p has full support in X, then the operator Σ p has a trivial null space.

(d) Given the mean element defined in Eq. ( 1), we have Σ p µ p ⊗ µ p .

Proof (a) The operator is well defined and has unit trace because tr[ϕ(x)ϕ(x

) * ] = ϕ(x) 2 = 1 for all x ∈ X. (b) If Σ p = Σ q , then for all y ∈ X, X k(x, y) 2 [dp(x) -dq(x)] = 0.
The representation of functions in the RKHS associated with k 2 as linear functions of x → k(x, y) 2 , and the universality of k 2 implies that X f (x)[dp(x) -dq(x)] = 0 for all continuous functions f , hence, p = q and the injectivity of the map. Property (c) is a simple consequence of g, Σ p g = E X∼p g(X) 2 , which is valid for all g ∈ H. If this expectation is equal to zero, and since elements of H are continuous functions, we must have g = 0 as soon as p has full support.

Property (d) is a consequence of the non-negativity of the variance of f (X), for f in H.

Since the operators are trace-class, they have a discrete spectrum [START_REF] Reed | Functional analysis[END_REF][START_REF] Brezis | Analyse Fonctionelle[END_REF], composed of a summable sequence of strictly positive eigenvalues tending to zero, and potentially zero (if p does not have full support in X).

The spectral decay is a key quantity for the study of kernel methods (see [21] and references therein), as used in Section 5 to estimate entropies from data.

Throughout this paper, we will consider spectral functions of trace-class self-adjoint operators A, such as log A, exp(A), or more generally f (A), for a function f : R → R. These are defined through the eigenvalue decomposition by keeping eigenvectors and replacing eigenvalues λ by f (λ) (for polynomials, we then recover the usual matrix polynomials).

Uniform distribution and symmetric sets.

In order to define the notion of kernel entropy from kernel relative entropy, we will need to use a specific "base" distribution τ , that will play the role of the uniform distribution when this notion is classically defined (e.g., X is a finite set or a subset of R d , like in all examples in Section 2.3).

For thi base distribution τ on X, we use the notation (instead of Σ τ ):

Σ = X ϕ(x)ϕ(x) * dτ (x).
In this paper, some concepts will be easier for symmetric sets, that is, sets where there exists a set T of bijective transformations t : X → X, such that for all x, y ∈ X, ∃t ∈ T such that t(x) = y. We will then assume that the kernel k and the base distribution τ are compatible with this notion of symmetry, that is, τ (A) = τ (t(A)) for any Borel set A and t ∈ T, and k(x, y) = k(t(x), t(y)) for all transformations t ∈ T.

Then the covariance operator for the base distribution will also exhibit symmetries. More precisely, all eigensubspaces will be invariant by symmetry, that is, if e 1 , . . . , e m ∈ H is an orthonormal basis of an eigensubspace corresponding to an eigenvalue λ, then4 m i=1 e i (t(x)) 2 = m i=1 e i (x) 2 for all transformations t ∈ T and x ∈ X. This implies that the function x → m i=1 e i (t(x)) 2 is constant, equal to its expectation under τ (which is itself equal to mλ).

As a consequence, for a spectral function Σ → f (Σ) based on a real-valued function f : R → R, the dotproduct ϕ(x), f (Σ)ϕ(x) can be expanded as the sum of the contributions coming from all eigensubspaces, with the one associated to the eigenvalue λ leading to f (λ) m i=1 e i (x) 2 being independent of x (and equal to mλf (λ)). This leads notably to ϕ(x), f (Σ)ϕ(x) = tr[Σf (Σ)] for all real-valued functions f defined on R, and to tr[Σf (Σ)] = tr[Σ p f (Σ)] for any probability distributions p. Surjectivity. By definition, every covariance operator is in the closure M of the convex hull of all ϕ(x)ϕ(x) * , for x ∈ X. We will also need the closure A of the span of these vectors. By definition, the mapping p → Σ p is surjective from probability distributions to M. Moreover, if the kernel k is universal,

then 5 M = A ∈ A, A 0, tr A = 1 . (3)
For orthonormal embeddings, where k(x, y) = 0 if x = y, M is (in the appropriate basis) the set of non-negative diagonal operators, while A is the set of diagonal operators (in the same basis).

Note that in most cases A is strictly included in the set of trace class operators from H to H, which is equivalent to the existence of bounded self-adjoint operators M : H → H such that ∀x ∈ X, ϕ(x), M ϕ(x) = 0. Since ϕ(x) 2 = 1 for all x ∈ X, this is equivalent to the existence of a positive bounded operator M so that ∀x ∈ X, ϕ(x), M ϕ(x) = 1, which, with an eigenvalue decomposition leads to a finite or countable family of functions (f i ) i∈I in H such that i∈I f i (x) 2 = 1, that is, a "partition of unity", which are well known to exist in many spaces (see, e.g., [START_REF] Berger | Differential Geometry: Manifolds, Curves, and Surfaces: Manifolds, Curves, and Surfaces[END_REF]Chapter 3] or [29, Section 1.4]).

Examples

Finite set. If X = {x 1 , . . . , x m } is a finite set of cardinality m, then the Gram matrix K ∈ R m×m defined as K ij = ϕ(x i ), ϕ(x j ) completely characterizes the kernel notions up to rotation. If K = I (orthonormal embedding), then there exist m orthonormal vectors u 1 , . . . , u m such that Σ p = m i=1 P(X = i)u i u * i , and thus, all covariance matrices commute and therefore share the same eigenbasis, and we will recover exactly the regular notions of discrete entropy.

We get a symmetric set if K is invariant by permutations, that is, K = αI + (1 -α)11 ⊤ , and as soon as K is invertible, a universal kernel. 6Polynomial kernels. When X is a subset of R d , then we can consider kernels of the form k(x, y) = (α + β x ⊤ y) r , which correspond to ϕ(x) composed of all monomials of order up to r, and the covariance operator is thus composed of traditional moments. Here, we can have k(x, x) = 1 for all x ∈ X only if X is included in a centered sphere. Alternatively, we may only impose k(x, x) 1 for all x ∈ X, which still preserves many of our results in later sections.

Translation-invariant kernels on a compact subset of R d . These are the usual kernels k(x, y) of the form k(x, y) = q(x -y) where q is a function with strictly positive Fourier transform, which are all universal [START_REF] Steinwart | On the influence of the kernel on the consistency of support vector machines[END_REF]. Classical examples are the Gaussian kernels k(x, y) = exp(-x -y 2 2 /σ 2 ), where H is a space of infinitely differentiable functions, or the exponential kernel k(x, y) = exp(-x -y 2 /σ), where H is the Sobolev space of functions with square-integrable derivatives up to order (d + 1)/2. However, they typically do not respect the symmetry of the underlying set (for example for spheres, like below).

Torus. If X = [0, 1], we can consider kernels of the form k(x, y) = k(x -y) where k : R → R is 1-periodic and has a strictly positive Fourier series (note that we use the same notation for the kernel and the underlying function). We then get a universal kernel and a symmetric set (invariant by periodic translation); moreover, the base distribution will be the usual uniform distribution. For example, we will often consider:

k(x) = (1 -e -σ ) 2 1 + e -2σ -2e -σ cos 2πx = 1 1 + sin 2 πx sinh 2 (σ/2)
, for which the Fourier series is k(ω) =

1 0 e -2iπωx k(x)dx = tanh( σ 2 )e -σ|ω| for ω ∈ Z.
This extends naturally to [0, 1] d , by considering tensor products, that is, k(x, y) = d i=1 k(x i -y i ), which for example leads to k(ω) = tanh d ( σ

2 )e -σ ω 1 for ω ∈ Z d in the example above. We will also need an explicit formula for ω∈Z k(ω) log k(ω) = log tanh σ 2 -σ sinh σ .

Spheres. We can also consider X as the unit ℓ 2 -sphere in dimension d, with kernels that are invariant by rotation, leading also to symmetric sets. Kernels are then defined through their expansion in spherical harmonics [START_REF] Smola | Regularization with dot-product kernels[END_REF], with links with neural networks [START_REF] Bach | Breaking the curse of dimensionality with convex neural networks[END_REF].

We could also imagine extensions to Stiefel or Grassman manifolds [START_REF] Edelman | The geometry of algorithms with orthogonality constraints[END_REF], and of course to all sets X where positive definite kernels can be defined (such as graphs, trees, permutations, sets, etc.) [START_REF] Shawe | Kernel Methods for Pattern Analysis[END_REF].

Review of quantum information theoretic quantities

We consider here complex-valued Hilbert spaces, but all results will apply to operators on real Hilbert spaces. In quantum mechanics (see, e.g., [START_REF] Isham | Lectures on Quantum Theory: Mathematical and Structural Foundations[END_REF]), physical states are represented as elements of a Hilbert space H, and a physical system is characterized by a convex combination of orthogonal projections on these states, which is often called a density matrix. Our covariance operators defined as Σ p = X ϕ(x)ϕ(x) * dp(x) can thus be seen as "density operators". It turns out that in quantum information theory (see, e.g., [START_REF] Mark | Quantum Information Theory[END_REF]), information-theoretic quantities are defined based on these density operators. In this section, we review these notions from a mathematical point of view.

Entropy. For a compact positive self-adjoint operator A on the complex Hilbert space H with finite trace, we can define the negative entropy as [START_REF] Araki | Entropy inequalities[END_REF]:

tr A log A = λ∈Λ(A) λ log λ,
where Λ(A) is the set of eigenvalues of A (with the convention that 0 log 0 = 0). It may not be finite if the sequence of eigenvalues is not decreasing sufficiently fast (note that this is unlikely since they are summable). As defined, tr[A log A] is always less than (tr A) log(tr A) and can be equal to -∞.

Kullback-Leibler (KL) divergence (relative entropy). It is defined as

D(A||B) = tr A(log A -log B) ,
for any two positive Hermitian operators with finite trace. It can only be finite if the null space of B is included in the null space of A, that is, if

A = B 1/2 M B
λ i > 0 such that n i=1 λ i = 1, then D n i=1 λ i A i n i=1 λ i B i n i=1 λ i D(A i B i ),
with equality if and only if log B i -log A i does not depend on i. Moreover, if there is equality, there exists an operator M such that for all i, A i = B i M .

(d) Monotonicity of quantum operations: Given operators

C i : H → K and A, B : H → H, i = 1, . . . , n, such that n i=1 C * i C i = I, then D n i=1 C i AC * i n i=1 C i BC * i D(A B).
The last property is classical in quantum information theory as the mapping A → n i=1 C i AC * i is usually referred to as a "quantum operation". Note that quantum operations are usually defined as "completely positive" trace-preserving maps, and that owing to Stinespring representation theorem (see [START_REF] Bhatia | Positive Definite Matrices[END_REF]Theorem 3.1.2]), they are essentially all of that form.

One particular application of monotonicity is to consider n positive self-adjoint operators D i : H → H, i = 1, . . . , n, such that n i=1 D i = I, and consider µ, ν ∈ R n defined through µ i = tr(D i A) and ν i = tr(D i B) for A, B a positive self-adjoint operator with unit trace, then µ and ν are on the simplex, and the monotonicity leads to (see detailed proof in Appendix A.2):

n i=1 µ i log µ i ν i D(A B). (4)
Further results. In Appendix A.1, we provide finer results beyond expectation with respect to a probability measure with finite support, as well as equality cases.

Integral representation. For the KL divergence, we have, by direct integration [START_REF] Ando | Concavity of certain maps on positive definite matrices and applications to Hadamard products[END_REF]:

D(A B) = tr A(log A -log B) = - +∞ 0 tr A(A + λI) -1 -tr A(B + λI) -1 dλ, (5)
which will be used in Section 5.2 for the estimation from a finite sample. As shown in Appendix A, we also have another integral representation formula, that helps in proving additional convexity results.

Kernel entropy and Kullback-Leibler divergence

We define the kernel Kullback-Leibler (KL) divergence, or kernel relative entropy, between two probability distribution p and q as D(

Σ p Σ q ) = tr Σ p (log Σ p -log Σ q ) . (6)
In this section, we list some of its properties, which are similar to the properties of (p, q) → µ p -µ q , but with the added benefit that is has a direct link with Shannon entropy, which is useful in itself when a precise information measure is needed (e.g., for differential privacy [START_REF] Dwork | Calibrating noise to sensitivity in private data analysis[END_REF]), or when used within variational inference (see Section 7).

First, the kernel relative entropy is finite under general conditions.

Proposition 3 (Finiteness of kernel KL divergence) Assume (A1). If p is absolutely continuous with respect to q, with dp dq ∞ α, then both the regular and kernel KL divergences are between zero and log α.

Proof We have then Σ p αΣ q , leading to (note that we must have α 1):

D(Σ p Σ q ) = tr Σ p (log Σ p -log Σ q )
tr Σ p (log(αΣ q ) -log Σ q ) = (log α) tr Σ p = log α.

The positivity follows from Prop. 2, property (b).

This implies that if p has a bounded density with respect to the base measure τ on X, then the relative entropy D(Σ p Σ) is well-defined, and always non-negative (recall that Σ is the covariance operator for the base probability measure τ on X).

We can then define the kernel entropy as

H(Σ p ) = -D(Σ p Σ) -min x∈X ϕ(x), (log Σ)ϕ(x) (7)
= -tr Σ p log Σ p + tr Σ p log Σ -min x∈X ϕ(x), (log Σ)ϕ(x) .

Note that we add terms on top of -tr Σ p log Σ p , one linear term in Σ p , and one constant term, to extend properties of the regular Shannon discrete entropy in all cases. However, everything simplifies in some cases. Indeed, when there are some invariances that are shared by the kernel and the base distribution τ (that is, a symmetric set as defined in Section 2.2), x → ϕ(x), (log Σ)ϕ(x) is constant on X, and thus equal to its expectation under p, which implies min x∈X ϕ(x), (log Σ)ϕ(x) = tr[Σ p log Σ], and then the definition above simplifies to H(Σ p ) = -tr Σ p log Σ p .

The kernel information quantities have natural properties mimicking the traditional quantities, as well as a bound by the regular KL divergence.

Proposition 4 (Properties of kernel entropy and relative entropy) Assume (A1).

(a) The kernel KL divergence defined in Eq. ( 6) is non-negative, and equal to 0 for p = q. If the kernel k 2 is universal, it is zero if and only if p = q.

(b) The kernel entropy defined in Eq. ( 7) is always between 0 and -min x∈X ϕ(x), (log Σ)ϕ(x) (which is equal to -tr Σ log Σ for symmetric sets). If p = τ (base measure), it is equal to the upper-bound, while for Dirac measures with a point mass at any minimizers of x → ϕ(x), (log Σ)ϕ(x) , it is equal to 0 (for symmetric sets, this is all of X). If the kernel k 2 is universal, these conditions are also necessary.

(c) The mapping (p, q) → D(Σ p Σ q ) is convex, and thus the mapping p → H(Σ p ) is concave.

(d) For all probability measures p and q, D(Σ p Σ q ) D(p q). Moreover, the differential entropy with respect to the base measure τ , that is, -X log dp dτ (x)dp(x), is less than H(Σ p )+min x∈X ϕ(x), (log Σ)ϕ(x) .

(e) For all probability measures p and q, D(

Σ p Σ q ) 1 2 Σ p -Σ q 2 * 1 2 Σ p -Σ q 2
HS , where • * denotes the nuclear norm and • HS the Hilbert-Schmidt norm.

Proof Property (a) is the consequence of the non-negativity of the relative entropy when the two arguments have the same trace (here equal to 1), and injectivity of the map p → Σ p when k 2 is universal.

For property (b), the upper-bound results follow similarly. For the non-negativity, we notice that H(Σ p ) defined in Eq. ( 7) is the sum of two non-negative terms -tr[Σ p log Σ p ] (negative because all eigenvalues of Σ p are less than one) and tr Σ p log Σ -min x∈X ϕ(x), (log Σ)ϕ(x) . It is equal to zero if and only the two terms are equal to zero. For the first term, this is equivalent to Σ p being rank one, while for the second one, this is equivalent to p being supported in the set of minimizers. The sufficient condition for being zero thus follows; for the necessary condition, if k 2 is universal, k(x, x ′ ) 2 < 1 if x = x ′ , and thus Σ p being rank one implies that p is a Dirac measure.

Property (c) is a consequence of the convexity of the map (A, B) → tr A(log A-log B) and of the linearity of the mapping p → Σ p . Property (d) is obtained from

D(Σ p Σ q ) = D X ϕ(x)ϕ(x) * dp(x) X dq dp (x)ϕ(x)ϕ(x) * dp(x) X D ϕ(x)ϕ(x) * dq dp (x)ϕ(x)ϕ(x) * dp(x) = X ϕ(x) 2 D 1 dq dp (x) dp(x) = X log dp dq (x) dp(x) = D(p q),
using joint convexity of the quantum relative entropy, and the fact that ϕ(x) = 1 for all x ∈ X. Note that we can only have equality if for all x = y, k(x, y) = 0, that is, we have an orthonormal embedding, which is only possible for discrete X. Note also that the upper-bound remains valid as soon as k(x, x) 1 for all x ∈ X.

Property (e) is simply a consequence of Pinsker inequality for the quantum relative entropy (see [START_REF] Yu | The strong convexity of von Neumann's entropy[END_REF]).

Note that while the relative entropy results do not depend on the choice of the base measure, the notion of kernel entropy depends on this choice (however, for the symmetric sets presented in Section 2.3, the natural choice of the base measure is the uniform measure).

As detailed in Appendix A.4, the properties outlined in Prop. 4, in particular the lower-bound on the relative entropy, is valid for other "f -divergences" [START_REF] Csiszár | Information-type measures of difference of probability distributions and indirect observation[END_REF][START_REF] Topsoe | Some inequalities for information divergence and related measures of discrimination[END_REF], such as the squared Hellinger distance, or the Pearson χ 2 -divergence.

Lower bound on kernel relative entropy

Beyond the upper-bound by the Shannon relative entropy (property (d) in Prop. 4), we can also get a lower bound on the kernel KL divergence, by defining for all y ∈ X, D(y

) = Σ -1/2 ϕ(y)ϕ(y) * Σ -1/2
, which is a positive self-adjoint operator, so that X D(y)dτ (y) = I. We can then apply the monotonicity of the relative entropy with operators D(y) 1/2 , y ∈ X, that is, Eq. ( 4) extended to general expectations. We have

tr[D(y)Σ p ] = X ϕ(x), Σ -1/2 ϕ(y) 2 dp(x) = X h(x, y)dp(x), defining the function h : X × X → R + as h(x, y) = ϕ(x), Σ -1/2 ϕ(y) 2 , which is non negative and such that ∀y ∈ X, X h(x, y)dτ (x) = ϕ(y), Σ -1/2 X ϕ(x)ϕ(x) * dτ (x) Σ -1/2 ϕ(y) = ϕ(y), ϕ(y) = 1.
The function h can thus be seen as a smoothing kernel7 . Thus, by monotonicity of the relative entropy (Appendix A.2), we get:

D(p q) D(Σ p Σ q ),
where p(y) = tr[D(y)Σ p ] = X h(x, y)dp(x) is a smoothed version of p (with the same definition for q). Note that we have directly D(p q) D(p q) from standard Markov chain arguments in information theory8 [18, Section 2.9]. Overall, we have the following sequence of inequalities:

D(p q) D(Σ p Σ q ) D(p q). (8)
These can lead to quantitative bounds between D(Σ p Σ q ) and D(p q), in particular when the smoothing function h is putting most of its mass on pairs (x, y) where x is close to y. This is made explicit below for distributions on a metric space and distributions with Lipschitz-continuous densities.

Small-width asymptotics for metric spaces

We now provide a bound between D(p q) and D(p q) in Eq. ( 8), and thus between D(Σ p Σ q ) and D(p q), when the set X is equipped with a distance d, used to characterize the regularity of density functions.9

This corresponds to bounding the difference in the classical data processing inequality for the Shannon entropy, using regularity properties of the densities of p and q. The following proposition is shown in Appendix B.

Proposition 5 (Joint bound on relative entropy) We assume that p, q have strictly positive Lipschitzcontinuous densities (denoted p and q) with respect to the base measure τ . Using definitions from Section 4.1, we have:

0 D(p q) -D(p q) E(p, q) × sup x∈X X h(x, y)d(x, y) 2 dy, (9)
with

E(p, q) = 2C 2 p/q sup x∈X q(x) p(x) + 4 sup x∈X q(x) p(x) (C p + C q ) 2 sup x∈X p(x) q(x) (1 + C 2 p diam(X) 2 )
, where C p , C q and C p/q are Lipschitz constants satisfying for all x, y ∈ X, q(x) q(y) -1 d(x, y)C q , p(x) p(y) -1 d(x, y)C p , and

p(y) q(y) -p(x) q(x)
d(x, y)C p/q for some distance d on X.

Thus, for Lipschitz-continuous densities, the approximation of D(p q) by D(Σ p Σ q ) is controlled by a kernel dependent quantity, that we now study in special cases.

Discrete orthonormal embeddings. For these kernels where k(x, y) = 1 x=y , then we have Σ = I and thus h(x, y) = 1 x=y ; therefore the quantity sup x∈X X h(x, y)d(x, y) 2 dy is always zero, which is natural here since in this situation, kernel and Shannon entropies are equal.

Torus. We have by definition h(x, y) = ϕ(x), Σ -1/2 ϕ(y) 2 . For a kernel on the torus [0, 1] of the form k(x, y) = k(x -y) where k has Fourier series k, then we have:

h(x, y) = ω∈Z k(ω) 1/2 e 2iπω(x-y) 2 ,
and for the particular choice of k(ω) ∝ e -σ|ω| , we have:

h(x, y) = tanh σ 2 tanh 2 σ 4 1 + sin 2 π(x -y) sinh 2 (σ/4) -2
, and for d(x, y) = | sin π(x -y)|, then, using exact integration, we get:

y+1/2 y-1/2 h(x, y)d(x, y) 2 dx = tanh σ 2 tanh 2 σ 4 sinh 3 (σ/4) 2(1 + sinh 2 (σ/4)) 3/2 σ 2 16
.

For the d-dimensional torus, and, d(x, y) 2 = d j=1 sin 2 π(x j -y j ), we get the bound

dσ 2 16
, and hence a nice scaling in dimension (non-exponential). Overall, in this situation, we thus get an approximation of entropy from Eq. ( 9) up to O(σ 2 ), where σ is the width of the kernel. Thus, if our goal is to estimate the true relative entropy, and not only have a lower bound, we need to have σ going to zero, which corresponds to RKHSs which are bigger and bigger. See discussion in Section 5.3 when estimating quantities from a finite sample.

Note that this extends to all translation invariant kernels on [0, 1] d (and not only for k(ω) ∝ e -σ ω 1 ), and more generally to translation-invariant kernels on compact spaces.

Alternative proof through leverage scores. While we provided a proof of the O(σ 2 ) approximation of the entropy through quantum data processing inequalities, the relationship between the regular entropy and the kernel entropy can also be studied through the integral representation in Eq. ( 5). We indeed have, using Eq. ( 5):

D(Σ p Σ q ) = +∞ 0 tr Σ p (Σ q + λI) -1 -(Σ p + λI) -1 dλ = +∞ 0 X ϕ(x), (Σ q + λI) -1 -(Σ p + λI) -1 ϕ(x) dp(x)dλ = +∞ 0 X ϕ(x), (Σ q + λI) -1 ϕ(x) -ϕ(x), (Σ p + λI) -1 ϕ(x) dp(x)dλ.
We now consider that both p and q have densities with respect to the base measure τ (which we also denote p and q). The main idea, using and extending results from [START_REF] Pauwels | Relating leverage scores and density using regularized Christoffel functions[END_REF], is that for all x and λ, and p smooth enough, the quantity ϕ(x), (Σ p + λI) -1 ϕ(x) , often referred to as a "leverage score" can be approximated as ϕ(x), (p(x)Σ + λI) -1 ϕ(x) .

We thus consider the potential approximation:

+∞ 0 X ϕ(x), (q(x)Σ + λI) -1 ϕ(x) -ϕ(x), ((p(x)Σ + λI) -1 ϕ(x) p(x)dτ (x)dλ = +∞ 0 X tr Σ(q(x)Σ + λI) -1 -tr Σ(p(x)Σ + λI) -1 p(x)dτ (x)dλ = X tr Σ(log p(x)Σ -log q(x)Σ) p(x)dτ (x) = X p(x) log p(x) q(x) dτ (x),
which is exactly the traditional KL divergence D(p q). In order to show a bound similar to Eq. ( 9), a non-asymptotic extension of results from [START_REF] Pauwels | Relating leverage scores and density using regularized Christoffel functions[END_REF] could be carried out.

Estimation from finite sample

Given that the covariance operator is an infinite-dimensional operator, naively we would need eigenvalue decompositions of such operators to compute the kernel information quantities, which is computationally hard. Thus finite-dimensional estimation algorithms are needed. In this section, we provide algorithms for estimating the kernel entropy H(Σ p ) from independent and identically distributed (i.i.d.) samples from p (as required for model fitting [START_REF] Miko Laj Bińkowski | Demystifying MMD GANs[END_REF] or for independent component analysis [START_REF] Muandet | Kernel mean embedding of distributions: A review and beyond[END_REF]), as well as from certain kernel integrals. Estimators from other oracles on p or q, as well as estimators for D(Σ p Σ q ) could derived and analyzed similarly. 10

Estimators

Given Σ p = X ϕ(x)ϕ(x) * dp(x) and x 1 , . . . , x n sampled i.i.d. from p, we consider the natural estimator

Σp = 1 n n i=1 ϕ(x i )ϕ(x i ) * . ( 10
)
The natural estimate for tr Σ p log Σ p is tr Σp log Σp , that can be computed from the kernel matrix K ∈ R n×n defined as

K ij = k(x i , x j ).
Proposition 6 (entropy for empirical covariance estimators) With Σp defined in Eq. ( 10), and the kernel matrix defined above, we have

tr Σp log Σp = tr 1 n K log 1 n K . ( 11
)
Proof The non-zero eigenvectors of Σp belong to the image space of Σp and are thus linear combinations (Kα) i -nλα i ϕ(x i ) = 0, which implies K 2 = nλKα and then Kα = nλα since Kα = 0. Thus, the non-zero eigenvalues of Σp are exactly the ones of 1 n K, and we thus get tr Σp log Σp = tr 1 n K log 1 n K .

f = n j=1 α j ϕ(x j ) for α ∈ R n . Then Σp f = 1 n n i=1 n j=1 α j ϕ(x i )ϕ(x i ) * ϕ(x j ) = 1 n n i=1 (Kα) i ϕ(x i ). Thus, if Kα = nλα, Σp f = λf ,
More generally, for any approximation of Σ p as n i=1 η i ϕ(x i )ϕ(x i ) * , with positive weights η i that sum to one (for example obtained from a finite grid and not sampling), we would get tr Σp log Σp = tr Diag(η)

1/2 K Diag(η) 1/2 log Diag(η) 1/2 K Diag(η) 1/2 .
Running-time complexity. In order to compute the eigenvalue decomposition of K above (needed to compute the matrix logarithm), we need a running time of O(n 3 ). This can be reduced by using column sampling techniques [START_REF] Boutsidis | An improved approximation algorithm for the column subset selection problem[END_REF], also known as "Nyström method" in this context [START_REF] Williams | Using the Nyström method to speed up kernel machines[END_REF]: by projecting the kernel matrix to the span of m of its columns, it leads to a running time in O(m 2 n), with some approximation that could be controlled explicitly [START_REF] Rudi | Less is more: Nyström computational regularization[END_REF].

Analysis

In order to quantify the difference between tr Σp log Σp and tr Σ p log Σ p , we could use regular perturbation arguments for eigenvalues of self-adjoint operators [START_REF] Kato | Perturbation Theory for Linear Operators[END_REF]. However, the function t → t log t has diverging derivatives around zero, which prevents their direct use. Instead, given the integral representation in Eq. ( 5), we see that it will be sufficient to estimate quantities of the form tr (Σ p + λI) -1 Σ p , which are usually referred so as "degrees of freedom" in kernel methods [21]. The natural estimator is tr ( Σp + λI) -1 Σp and its performance has already been thoroughly studied (see [START_REF] Rudi | Generalization properties of learning with random features[END_REF]), and an important quantity there is the maximal leverage score for the base distribution sup x∈X ϕ(x), (Σ + λI) -1 ϕ(x) , which could be replaced by the the average one (which is then equal to the degrees of freedom) if "leverage score sampling" is used (see [21] for more details). For simplicity, we do not consider leverage score sampling, noting that for symmetric sets, it does not change anything.

The next proposition shows that these estimation results for fixed λ can be integrated over λ to get an estimation bound for the entropy (see proof in Appendix C.4), without the need for extra regularization, leading to an estimation in O(1/ √ n).

Proposition 7 Assume (A1) and that p has a density with respect to the base measure which is greater

than α < 1. Assume that c = +∞ 0 sup x∈X ϕ(x), (Σ + λI) -1 ϕ(x) 2 dλ is finite. Given i.i.d. samples x 1 , . . . , x n
from p, and the estimator defined by Eq. ( 11), we have:

E tr Σp log Σp -tr Σ p log Σ p 1 + c(8 log n) 2 nα + 17 √ n 2 √ c + log n .
Note that as a consequence of Jensen's inequality, we always have E tr Σp log Σp tr Σ p log Σ p , as can be seen in simulations (see Figure 1).

Examples

We consider the d-dimensional torus [0, 1] d with a translation-invariant kernel k(x, y) = k(x -y) with k(ω) = k(0)e -σ ω 1 . We have by symmetry ϕ(x),

(Σ + λI) -1 ϕ(x) = tr Σ(Σ + λI) -1 = ω∈Z d k(ω) k(ω) + λ .
As shown in Appendix C.6, we have the bound tr Σ(Σ + λI)

-1 e -σd/2 d! sinh d (σ/2) 1 + log tanh d σ 2 λ d .
Thus, the constant c in Prop. 7 above can be upper bounded by a constant times σ -d (see Appendix C.6 for details). We thus get an overall estimation rate for the kernel entropy proportional to σ -d/2 √ n.

When balancing with the O(σ 2 ) estimation of the entropy for Lipschitz-continuous strictly positive densities in Section 4.2, we can take σ ∝ n -1/(d+4) (lower bandwidth with more observations), leading to a rate of n -2/(d+4) , which is to be compared to the optimal rate equal to n -4/(d+4) for entropy estimation for densities with the same regularity [START_REF] Han | Optimal rates of entropy estimation over Lipschitz balls[END_REF]. Although our primary goal was not to estimate entropies, we get a rate which is the square root of the optimal rate.

Illustrative experiments. We consider X = [0, 1] and p(x) = 4 x-1 2 , with negative differential entropy

1 0 p(x) log p(x)dx = log 2- 1 2 ≈ 0.
1931. We consider estimation from various numbers of i.i.d. observations. See left plot in Figure 1, where we compare two values of σ. Note that the estimation is empirically here in expectation from above, as expected because of Jensen's inequality (as opposed to the right plot). ). Right: from kernel integrals evaluated as uniform samples (method from Section 5.4). As expected we estimate the kernel entropy (limit when n tends to infinity) from above on the left plot, and from below on the right plot, with a faster convergence rate on the right plot.

Estimation from other oracles

If we have access to other oracles on p and q, we can consider other estimators. 11 For example, if we can estimate X k(x, y)k(x, z)dp(x) and X k(x, y)k(x, z)dq(x) for any y, z ∈ X, we can simply sample n points x 1 , . . . , x n from the base measure on X, and, given the projection operator Π on the span of ϕ(x 1 ), . . . , ϕ(x n ), consider Σp = ΠΣ p Π and Σq = ΠΣ q Π, (

and the estimator D( Σp Σq ) which is now a pointwise lower bound 12 on D(Σ p Σ q ) (and thus on D(p q)), which can be beneficial in variational inference, as presented in Section 7.

Noting that ϕ(x i ), Σ p ϕ(x j ) = X k(x i , x)k(x, x j )dp(x) (and similarly for Σ q ), and that K is the matrix of dot-products of ϕ(x 1 ), . . . , ϕ(x n ), the orthogonal projection may be done as

D( Σp Σq ) = D K -1/2 X k(X, x)k(x, X)dp(x)K -1/2 K -1/2 X k(X, x)k(x, X)dq(x)K -1/2 ,
where k(X, x) ∈ R n is the vector with components k(x i , x), i = 1, . . . , n, and K ∈ R n×n the kernel matrix, typically with some added regularization for numetical stability, and incomplete Cholesky decomposition so as to minimize the number of kernel evaluations (see, e.g., [START_REF] Bach | Predictive low-rank decomposition for kernel methods[END_REF] and references therein), with an overall running time complexity of O(n 3 ). A key property is that it always underestimates the kernel relative entropy (and hence the Shannon relative entropy). The following proposition shows that it provides a better estimate (better dependence in n) than with i.i.d. samples (see proof in Appendix C.5).

Proposition 8 Assume (A1) and that p has a density with respect to the base measure which is bounded. Given uniform i.i.d. samples x 1 , . . . , x n from the base distribution, and the estimator defined by Eq. ( 12), we have, if for all λ > 0, sup x∈X ϕ(x), (Σ + λI) -1 ϕ(x)

B log(C/λ) d , for some constants B, C:

0 tr ΠΣ p Π log( ΠΣ p Π) -tr Σ p log Σ p = O log(BC) d ( √ BC + B) exp - 1 4 (Bn) 1/(d+1) .
The proposition above applies to the torus situation in [0, 1] d and the kernel with exponential Fourier series, with now an improved dependence in n and in d, since in that case

B ∼ C ∼ σ d . That is, to reach precision ε, we need n ∝ σ -d log σ d ε d+1 instead of n ∼ σ -d ε -2 .
Other decays for the maximal degrees of freedom are considered in Appendix C.5.

Illustrative experiments. We consider as above p(x) = 4 x -1 2 , and consider estimation from various numbers of i.i.d. observations for the projection. The required integrals are computed by numerical integration. See right plot in Figure 1, where we compare two values of σ. Note that the estimation is as expected here from below, and with a faster rate of convergence than for the estimation method from Section 5.1.

Multivariate probabilistic modelling

Information-theoretic quantities are particularly useful at characterizing statistical dependence between random variables. This is for example the basis of independent component analysis [START_REF] Cardoso | Dependence, correlation and Gaussianity in independent component analysis[END_REF], and can be used for feature selection [START_REF] Song | Feature selection via dependence maximization[END_REF]. We explore how the new notions of entropy can be used in this goal.

Joint entropy, marginal entropy and mutual information

We assume that we have a product set X = X 1 × X 2 , and a feature space H = H 1 ⊗ H 2 , corresponding to the kernel k(x, y) = k 1 (x 1 , y 1 )k 2 (x 2 , y 2 ), and feature map ϕ(x) = ϕ 1 (x 1 ) ⊗ ϕ 2 (x 2 ). We also consider a base distribution on X 1 × X 2 which is the product of the base distributions. The tensor product is here defined as the span of all functions (

x 1 , x 2 ) → f 1 (x 1 )f 2 (x 2 ) for f i ∈ H i , i = 1, 2, with the norm f 1 • f 2 .
Given a joint probability p X 1 X 2 and its marginal probability distributions p X 1 and p X 2 , we can define the covariance operators Σ p X 1 X 2 on H 1 ⊗ H 2 , Σ p X 1 on H 1 , and Σ p X 2 on H 2 , as well as Σ 12 , Σ 1 and Σ 2 the covariance operators of the base distributions on X 1 × X 2 , X 1 and X 2 :

We then define the marginal and joint entropies as:

H(Σ p X 1 ) = -tr Σ p X 1 log Σ p X 1 + tr Σ p X 1 log Σ 1 -min x 1 ∈X 1 ϕ 1 (x 1 ), (log Σ 1 )ϕ 1 (x 1 ) H(Σ p X 2 ) = -tr Σ p X 2 log Σ p X 2 + tr Σ p X 2 log Σ 2 -min x 2 ∈X 2 ϕ 2 (x 2 ), (log Σ 2 )ϕ 2 (x 2 ) H(Σ p X 1 X 2 ) = -tr Σ p X 1 X 2 log Σ p X 1 X 2 + tr Σ p X 1 X 2 log Σ 12 -min x∈X ϕ(x), (log Σ 12 )ϕ(x) .
The product distribution p X 1 p X 2 (with independent components with same marginal distributions as X), which is different in general from the joint distribution p X 1 X 2 , has covariance operator Σ p X 1 p X 2 = Σ p X 1 ⊗ Σ p X 2 , and its entropy is:

H(Σ p X 1 p X 2 ) = -tr Σ p X 1 ⊗ Σ p X 2 log Σ p X 1 ⊗ Σ p X 2 + tr Σ p X 1 ⊗ Σ p X 2 log Σ 12 -min x∈X ϕ(x), (log Σ 12 )ϕ(x) .
We have log Σ p X 1 ⊗Σ p X 2 = log Σ p X 1 ⊗I +log I ⊗Σ p X 2 . Moreover, because the base distribution is separable by assumption, we have ϕ(x), (log Σ 12 )ϕ(x) = ϕ 1 (x 1 ), (log Σ 1 )ϕ 1 (x 1 ) + ϕ 2 (x 2 ), (log Σ 2 )ϕ 2 (x 2 ) , thus, we exactly have:

H(Σ p X 1 p X 2 ) = H(Σ p X 1 ) + H(Σ p X 2 ).
We can therefore define the kernel mutual information (not to be confused with the regular Shannon mutual information) as

I(X 1 , X 2 ) = D(Σ p X 1 X 2 Σ p X 1 ⊗ Σ p X 2 ).
From the developments above, it is equal to

I(X 1 , X 2 ) = H(Σ p X 1 ) + H(Σ p X 2 ) -H(Σ p X 1 X 2 ),
mimicking the traditional equality for regular entropies. Moreover, it is equal to zero if and only if the variables are independent.

Proposition 9 (Characterization of independence) Assume (A1). Then I(X 1 , X 2 ) is a lower bound on the Shannon mutual information, and

I(X 1 , X 2 ) = 0 if and only if Σ p X 1 X 2 = Σ p X 1 ⊗ Σ p X 2 .
If the kernel k 2 is universal, this is equivalent to the random variables X 1 and X 1 being independent.

Proof The first statement is straightforward because of the interpretation as a relative entropy. The second statement is simply the consequence of the fact that Σ p X 1 ⊗ Σ p X 2 is the covariance operator of the product distribution p X 1 p X 2 .

We can make the following observations:

• The extension to mutual information between more than two random variables is straightforward by considering tensor products between all marginals, with the same characterization property of independence.

• Our use of covariance operators is related to the "kernel generalized variance" from [START_REF] Bach | Kernel independent component analysis[END_REF], which considers a joint covariance operator

E ϕ 1 (X 1 )ϕ 1 (X 1 ) * E ϕ 1 (X 1 )ϕ 2 (X 2 ) * E ϕ 2 (X 2 )ϕ 1 (X 1 ) * E ϕ 2 (X 2 )ϕ 2 (X 2 ) * = E ϕ 1 (X 1 ) ϕ 2 (X 2 ) ϕ 1 (X 1 ) ϕ 2 (X 2 ) *
, and compute the mutual information for a Gaussian vector with this covariance operator. The kernel generalized variance can be generalized to more than two random variables, but then the independence characterization property is not satisfied. Moreover, the relationship with the usual mutual information is not as straightforward.

• The extension to conditional entropies is not straightforward. We could define Σ p X 2 |X 1 =x 1 in the natural way, but then we do not have the joint entropy being the sum of the marginal and conditional entropy (otherwise, because of axiomatic characterizations of entropy [START_REF] Csiszár | Axiomatic characterizations of information measures[END_REF], we would get back exactly the Shannon entropy). There is also a potential for exponential families and conditional models akin to generalized models, which remains to be explored.

Data processing inequality

We have a form of data processing inequality, which is classical in quantum information theory (see, e.g., [START_REF] Petz | Monotonicity of quantum relative entropy revisited[END_REF]), which we extend here.

Proposition 10 (Data processing inequality) Assume (A1). We have

D(Σ p X 1 X 2 Σ q X 1 X 2 ) D(Σ p X 1 Σ q X 1 ),
with equality implying (when the kernel k is universal) that the conditional distribution of X 2 given X 1 is the same for p and q, but the reverse implication is not true in general.

Proof The marginal operators Σ p X 1 and Σ q X 1 are obtained by taking "partial traces" which are quantum operations as defined in Section 3, hence the inequality. Following [START_REF] Beth | Inequalities for quantum entropy: A review with conditions for equality[END_REF], we have inequality if and only log

Σ p X 1 X 2 -log Σ p X 1 X 2 = log(Σ p X 1 ⊗ I) -log(Σ q X 1 ⊗ I)
, which is unfortunately not always satisfied when conditional distributions are equal. However, in case of equality, the extra condition (c

) from Prop. 2, becomes Σ p X 1 X 2 (Σ q X 1 ⊗ I) = Σ q X 1 X 2 (Σ p X 1 ⊗ I),
and can be leveraged. Indeed, for all test functions f 1 , g 1 , f 2 , g 2 in H, it leads to:

X 1 ×X 2 X 1 f 1 (x 1 )f 2 (x 2 )g 1 (y 1 )g 2 (y 2 )k 1 (x 1 , y 1 )dp X 1 X 2 (x 1 , x 2 )dq X 1 (y 1 ) = X 1 ×X 2 X 1 f 1 (x 1 )f 2 (x 2 )g 1 (y 1 )g 2 (y 2 )k 1 (x 1 , y 1 )dq X 1 X 2 (x 1 , x 2 )dp 1 (y 1 ).
Considering f 1 and g 1 tending to Dirac functions at x1 (which is possible since the kernel is univeral) we get

X 2 f 2 (x 2 )g 2 (y 2 )dp X 1 X 2 (x 1 , x 2 )dq X 1 (x 1 ) = X 2 f 2 (x 2 )g 2 (y 2 )dq X 1 X 2 (x 1 , x 2 )dp 1 (x 1
), which implies that the conditional distributions of x 2 given x 1 is the same for p and q.

Submodularity and conditional independence

We can now consider three variables X 1 , X 2 , X 3 , and use strong sub-additivity of the quantum relative entropy [START_REF] Beth | Inequalities for quantum entropy: A review with conditions for equality[END_REF][START_REF] Petz | Monotonicity of quantum relative entropy revisited[END_REF].

Proposition 11 (Submodularity) Assume (A1). We have:

H(Σ p X 1 X 2 X 3 ) -H(Σ p X 1 X 2 ) -H(Σ p X 2 X 3 ) + H(Σ p X 2 ) 0,
with equality (for universal kernels) implying that X 1 and X 3 are independent given X 2 , but in general the reverse implication is not true.

Proof As done by [START_REF] Petz | Monotonicity of quantum relative entropy revisited[END_REF], we can simply apply the data processing result (Prop. 10) to X1 = (X 1 , X 2 ), and X2 = X 3 , with the distributions p X1 X2 = p X 1 X 2 X 3 and q X1 X2 = p X 1 p X 2 X 3 .

As a classical consequence, the entropy function is submodular like the regular entropy (see, e.g., [3, Section 6.5]).

Note that we do not obtain a necessary and sufficient characterization of conditional independence, which can be obtained with other tools based on covariance operators [START_REF] Fukumizu | Kernel measures of conditional dependence[END_REF].

Convex duality, log-partition functions and variational inference

Given that our kernel notions of relative entropies are lower bounds on the regular notions, by the traditional convex duality results between maximum entropy and log-partition functions (see, e.g., [START_REF] Wainwright | Graphical Models, Exponential Families, and Variational Inference[END_REF]), we should obtain upper-bounds on log-partition functions. In this section, we show how such bounds can be obtained and computed. We then show that such bounds can also be optimized with respect to the positive definite kernel.

This leads to a new family of variational inference methods, which can then be used in various probabilistic inference tasks [START_REF] Wainwright | Graphical Models, Exponential Families, and Variational Inference[END_REF], and in particular in Bayesian inference [START_REF] Robert | The Bayesian choice: from decision-theoretic foundations to computational implementation[END_REF].

In this section, we will sometimes only assume that for our kernels satisfy (on top of being positive definite) k(x, x) 1 for all x ∈ X, noting that the upper-bound result on the regular Shannon relative entropy is still valid in this case.

Convex duality between operators

We have for self-adjoint operators (regardless of the link with input spaces and probabilities): with the constraints that A 0 and tr A = 1 are automatically satisfied (if tr A = 1, then by replacing M by M + tI, we can make the quantity go to infinity).

We also have (when the constraint for unit trace is removed): with the constraint that A 0 is automatically satisfied (but not tr A = 1), with the optimal M equal to log A -log B.

Bounds on the log-partition function

We can now apply the duality relationships above to covariance operators. Given a bounded function f : X → R and a distribution q on X, our goal is to obtain upper-bounds on the log-partition function log X e f (x) dq(x) which is equal to, by convex duality:

log X e f (x) dq(x) = sup p probability X f (x)dp(x) -D(p q).
In our situation, our operators are subject to belong to the subspace A, equal to the span of all ϕ(x)ϕ(x) * , for x ∈ X, which leads to slight adaptations. We do not necessarily assume that the kernels are universal, so that the hull of all ϕ(x)ϕ(x) * may not be equal to the intersection of A and positive operators with unit trace.

Isotropic kernels. If k(x, x) = 1, for all x ∈ X, then tr Σ p = 1 is equivalent to X dp(x) = 1. We define for a bounded function f : X → R,

a(f ) = sup p probability measure X f (x)dp(x) -D(Σ p Σ q ).
Using convex duality tools from Section 7.1, we have:

a(f ) = sup p non-negative measure inf M X f (x)dp(x) -tr Σ p M + log tr exp(M + log Σ q ) = inf M log tr exp(M + log Σ q ) such that ∀x ∈ X, f (x) ϕ(x), M ϕ(x) .
Note that we can get rid of the constraint X dp(x) = 1 because it is equivalent to tr Σ p = 1.

A relaxation is to relax p not to be a non-negative measure (and use the positivity of Σ p as the only constraint), leading to

b(f ) = sup p measure X f (x)dp(x) -D(Σ q Σ q ). ( 13
)
By construction, we have a(f ) b(f ), and by convex duality,

b(f ) = sup p measure inf M X f (x)dp(x) -tr Σ p M + log tr exp(M + log Σ q ) (14) = inf M log tr exp(M + log Σ q ) such that ∀x ∈ X, f (x) = ϕ(x), M ϕ(x) .
Note that that if the problem above is non-feasible because ϕ is too small (so that f cannot be represented as a quadratic form in ϕ), then we get a value +∞.

If we can write (usually non-uniquely), f (x) = ϕ(x), F ϕ(x) for some self-adjoint operator F , then we have another representation as:

b(f ) = inf N ∈A ⊥ log tr exp(F + N + log Σ q ).
Since D(Σ p Σ q ) D(p q), we have:

b(f ) a(f ) sup p probability X f (x)dp(x) -D(p q) = log X e f (x) dq(x).
As expected, we thus obtain upper bounds on the usual log-partition function.

Note that we could also get lower bounds using tools from Section 4.1, in particular for functions f that can be written as f (x) = X h(x, y)g(y)dτ (y).

Non isotropic kernels. If we only assume that k(x, x) 1 for all x ∈ X, then b(f ) defined in Eq. ( 13) and Eq. ( 14) is not any more an upper-bound but we can define instead

b (c) (f ) = sup p measure, X dp(x)=1, X k(x,x)dp(x) 1 inf M X f (x)dp(x) -tr Σ p M -tr Σ p + tr exp(M + log Σ q ) = sup p measure inf M, c, b 0 X f (x)dp(x) -tr Σ p M -tr Σ p + tr exp(M + log Σ q ) -c X dp(x) -1 -b X k(x, x)dp(x) -1 = inf M, c, b 0 c + b + tr exp(M + log Σ q ) such that ∀x ∈ X, f (x) = c + ϕ(x), (M + I + bI)ϕ(x) = inf M, c c + inf b 0 b + e -b tr exp(M + log Σ q ) such that ∀x ∈ X, f (x) = c + ϕ(x), (M + I)ϕ(x) = inf M, c c + 1 + log tr exp(M + log Σ q ) such that ∀x ∈ X, f (x) = c + ϕ(x), (M + I)ϕ(x) = inf M, c c + 1 + log 1 e tr exp(M + log Σ q ) such that ∀x ∈ X, f (x) = c + ϕ(x), M ϕ(x) ,
with log(a) = log a if a > 1, and a -1 otherwise.

If f can be represented through the operator F , and the constant function through the operator U , then we get

b (c) (f ) = inf M, c c + 1 + log 1 e tr exp(M + log Σ q ) such that F -cU -M ∈ A ⊥ .
We can check that if U = I, then the constraints in the definition of b (c) (f ) are equivalent to tr Σ p , and we exactly recover the expression of b(f ) in Eq. ( 14).

Properties of upper-bounds. The definition of upper-bounds on the log-partition function naturally leads to several questions and applications:

• Can these upper-bounds be tight? In situations like in Section 4.2 where the entropy lower-bound can be made tight (for example when the kernel bandwidth tends to zero), we should also get tight upper-bounds on the log-partition functions, that is, upper-bounds that tend to the classical logpartition function when the bandwidth goes to zero. We leave for future work a detailed study of these approximations.

• Can these upper-bounds be estimated efficiently from simple oracles on the function f ? We provide efficient algorithms in Section 7.4.

• Can these upper-bounds be minimized? We show in Section 7.5 that these upper-bounds happen to be convex in the kernel, so minimizing with respect to the kernel is possible.

Relationship with optimization

When adding a temperature parameter ε > 0, we can extend the traditional link between between optimization and log-partition functions. This corresponds essentially to considering a potential 1 ε f , or defining

b ε (f ) = εb 1 ε f = sup p measure X f (x)dp(x) -εD(Σ q Σ q ).
By convex duality, we get:

b ε (f ) = inf M ε log tr exp 1 ε M + log Σ q such that ∀x ∈ X, f (x) = ϕ(x), M ϕ(x) . When ε tends to zero, then b ε (f ) converges to inf M λ max (M ) such that ∀x ∈ X, f (x) = ϕ(x), M ϕ(x) .
Given that ϕ(x) = 1 for all x ∈ X, by writing M = cI -A for A positive, this is equal to

inf c∈R, A 0 c such that ∀x ∈ X, f (x) = c -ϕ(x), Aϕ(x) ,
which is exactly the optimization formulation of [START_REF] Rudi | Finding global minima via kernel approximations[END_REF] based on "kernel sums-of-squares". Following the traditional relationship between log-partition functions and maxima, we can thus consider our bound on log-partition functions as a smoothed version on the maximum. Therefore, some of the approximation techniques developed for the optimization formulation [START_REF] Rudi | Finding global minima via kernel approximations[END_REF] can be extended to our set-up as well. Moreover, as common in convex optimization, the entropy could be used for smoothing and accelerated optimization algorithms [START_REF] Nesterov | Smooth minimization of non-smooth functions[END_REF].

Computable bounds

In order to compute or approximate b(f ), which depends on the function f : X → R and the feature map ϕ : X → H, we will consider a finite-dimensional approximation φ : X → H where H will be a finite-dimensional subspace of H. All of our approximations will always be upper-bounds on the true log-partition function.

We will always need an efficient representation of à the span of all φ(x) φ(x) * , for x ∈ X, and sometimes of M the hull of all φ(x) φ(x) * for x ∈ X. See examples in Section 7.6 and Section 7.7 for the torus and the hypercube.

In order to estimate the relative entropy D(Σ p Σ q ), we will use D( Σp Σq ) instead. If the kernel (x, y) → k(x, y)-k(x, y) is positive definite, then this is a lower bound, and thus the bounds on the regular quantities are preserved.

In order to compute D( Σp Σq ), we assume we can compute explicitly

X ϕ(x)ϕ(x) * dp(x) and X ϕ(x)ϕ(x) * dq(x).
If φ is obtained from projection on a span of some ϕ(x 1 ), . . . , ϕ(x n ) like in Section 5.4, then this can be done if we can compute expectations of k(x, y)k(x, z) for any y, z ∈ X. We can also build set-specific approximate feature maps (as done in Section 7.6 and Section 7.7).

With an explicit representation of f . If we have an explicit (non-unique) representation of f as f (x) = φ(x), F φ(x) . We then have if φ(x), φ(x) = 1 for all x ∈ X (a similar approximation can be obtained for b (c) (f ) for non-isotropic approximations):

b(f ) = sup C∈A, C 0, tr C=1 tr[C F ] -D(C Σq ) = inf M log tr exp( M + log Σq ) + sup C∈A tr[ C(F -M )] = inf M log tr exp( M + log Σq ) such that M -F ∈ A ⊥ ,
which is a finite-dimensional convex optimization problem. Since the function M → log tr exp( M + log Σq ) is smooth, the optimization problem above can be approximately solved by projected gradient descent (with unit step-size), with iteration:

M ← Π F +A ⊥ M - exp( M + log Σq ) tr exp( M + log Σq )
.

We then obtain the (approximated) optimal C as C = exp( M +log Σq) tr exp( M +log Σq)

. It can be accelerated using extrapolation [START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF]. Note that for any feasible M we obtain a lower bound on b(f ).

Without an explicit representation of f . If an explicit representation of f as f (x) = φ(x), F φ(x) is not available, we can build an approximate representation of f such that f -φ(•), F φ(•) ∞ is as small as possible, and then use the previous algorithm, with this extra approximation factor. On the torus, we could also use the same technique as [START_REF] Woodworth | Non-convex optimization with certificates and fast rates through kernel sums of squares[END_REF] based on approximating the L ∞ -norm with Fourier transforms.

Illustrative experiments. We consider here X = [0, 1], and the kernel with Fourier transforms proportional to e -σ|ω| . We consider the function f (x) = cos(2πx), and approximate log( 1 0 e f (x) dx) ≈ 0.2359. We could either use the finite feature map from Section 5.4 based on random projections (which is a non-isotropic approximation). For simplicity, we consider φ(x) ω = k(ω)e 2iπωx , only for ω ∈ {-r, . . . , r}. We have φ(x) 2 = r ω=-r tanh( σ 2 )e -σ|ω| = tanh( σ 2 ) -1 + 2 1-e -σ(r+1)

1-e -σ = 1 -e -σ(r+1/2) cosh( σ2 ) , and the set A is exactly the set of Toeplitz matrices. See Figure 2 for an illustration with several values of σ.

Kernel learning

The key property we will leverage is the concavity of D(Σ p Σ q ) with respect to the kernel. This allows to maximize efficiently with respect to the kernel to get as close to the Shannon relative entropy D(p q) as possible.

Concavity with respect to the kernel. Given two distributions p and q, then D(Σ p Σ q ) depends on the kernel function k : X × X → R. This happens to be a concave function of k (which is only assumed to be positive definite with no bounds on k(x, x)). Indeed, if k 1 and k 2 are two kernels with feature maps ϕ 1 : X → H 1 and ϕ 2 : X → H 2 , with the associated covariance operators Σ and

η 1 , η 2 0 such that η 1 + η 2 = 1, then with H = H 1 ⊗ H 2 and ϕ(x) = ( √ η 1 ϕ 1 (x), √ η 2 ϕ 2 (x))
, using E 1 = I 0 0 0 and E 2 = 0 0 0 I , using the monotonicity of relative entropy for quantum operations, we have:

η 1 D(Σ (1)
p

Σ (1) q ) + η 2 D(Σ (2) p Σ (2) q ) = D E 1 Σ p E 1 + E 2 Σ p E 2 E 1 Σ q E 1 + E 2 Σ q E 2 D(Σ p Σ q ),
which shows the concavity.

Moreover, the relative entropy is monotonic, that is, if

k 1 -k 2 is positive definite, D(Σ (2) p Σ (2)
q ) D(Σ (1) p Σ (1)
q ). Therefore, when optimizing over the kernel, we only need to look at maximal elements with respect to the positive definite order.

Linear parametrization of kernels. Given a finite-dimensional feature map ψ : X → R d , we look at kernels of the form k(x, y) = ψ(x) ⊤ Λψ(x), where Λ ∈ R d×d is a symmetric positive semidefinite matrix. In order to have the upper-bound on entropy, we consider the set L of matrices Λ such that for all x ∈ X, ψ(x) ⊤ Λψ(x) 1. This is a finite-dimensional convex set, for which we need to know inner approximations (to preserve directions of bounds). We then have that

D(Σ p Σ q ) = D X Λ 1/2 ψ(x)ψ(x) ⊤ Λ 1/2 dp(x) X Λ 1/2 ψ(x)ψ(x) ⊤ Λ 1/2 dq(x)
is convex in Λ (see a direct proof in Appendix A.3). As mentioned earlier, by monotonicity, we only need to consider maximal elements of L.

Relative entropy maximization. Given p and q, from which we can compute

C p = E X∼p ψ(X)ψ(X) ⊤ and C q = E X∼q ψ(X)ψ(X) ⊤ , we maximize D(Λ 1/2 C p Λ 1/2 Λ 1/2 C q Λ 1/2 ) with respect to Λ ∈ L.
We can decompose it as the difference of two convex functions as:

D(Λ 1/2 C p Λ 1/2 Λ 1/2 C q Λ 1/2 ) = tr C 1/2 p ΛC 1/2 p log(C 1/2 p ΛC 1/2 p ) -tr Λ 1/2 C p Λ 1/2 log(Λ 1/2 C q Λ 1/2 ) .
This naturally leads to an optimization algorithm where we lower-bound the first convex function by its tangent and then maximize this affine function minus the other convex functions. This is a parameter-free globally and monotonically convergent algorithm.

Log-partition minimization. Given q and f : X → R, we can try to minimize the lower bound on the log-partition function instead. We assume for simplicity that we know a representation of f as f (x) = ψ(x) ⊤ Gψ(x), that the span B of all ψ(x)ψ(x) ⊤ for x ∈ X is manageable, and that (for simplicity), ψ(x) ⊤ Λψ(x) = 1 for all Λ ∈ L.

We thus need to minimize with respect to Λ ∈ L the following function (which is convex in Λ as the supremum of convex functions):

sup

C∈B, C 0, tr C=1 tr(CG) -D(Λ 1/2 CΛ 1/2 Λ 1/2 C q Λ 1/2 ) = inf M sup C∈B tr C(G -Λ 1/2 M Λ 1/2 ) + log tr exp(M + log(Λ 1/2 C q Λ 1/2 )) = inf M log tr exp(M + log(Λ 1/2 C q Λ 1/2 )) such that G -Λ 1/2 M Λ 1/2 ∈ B ⊥ .
The optimal C and M can be obtained by accelerated projected gradient descent. We can then minimize with respect to Λ by updating Λ like for the entropy maximization, once M (and thus C) have been estimated.

We now provide two illustrating examples, on the torus and the hypercube.

Example 1: Torus [0, 1]

We consider X = [0, 1], and assume that we have m = 2r + 1 frequencies ω ∈ {-r, . . . , +r}, and ψ ω (x) = e 2iπωx . Then we have (C p ) ωω ′ = p(ω -ω ′ ) and (C q ) ωω ′ = q(ω -ω ′ ), and the constraint to be in B is that C p and C q are Toeplitz matrices (note that this requires to know p(ω) and q(ω) for ω ∈ {-2r, . . . , 2r}).

If f has band-limited Fourier series f (zero outside {-2r, . . . , 2r}), then we have f (x) = ψ(x), Gψ(x) , with any matrix G such that ω-ω ′ =δ G ωω ′ = f (δ) for all δ ∈ {-2r, . . . , 2r}. We then have tr

[C p G] = r ω,ω ′ =-r p(ω -ω ′ ) * G ωω ′ = 2r δ=-2r p(δ) * f (δ)
. Moreover, we consider Λ = Diag(η) with η in the simplex, so that, for q(ω) = 1 ω =0 (for the uniform distribution), we have

D(Λ 1/2 C p Λ 1/2 Λ 1/2 C q Λ 1/2 ) = tr C 1/2 p ΛC 1/2 p log(C 1/2 p ΛC 1/2 p ) - r ω=-r η(ω) log η(ω).
Relative entropy maximization. We can maximize with respect to η, using the (here convergent by convexity) the majorization-minimization approach, e.g., by using the affine Taylor expansion of the first term, and minimizing tr

C 1/2 p Diag(η + )C 1/2 p log(C 1/2 p Diag(η + )C 1/2 p ) - r ω=-r η + (ω) log η + (ω) with respect to η + is in the simplex, thus leading to log η + equal to a constant plus diag C 1/2 p log(C 1/2 p Diag(η)C 1/2 p )C 1/2
p . This makes an algorithm in closed form to maximize D(Σ p Σ q ) with respect to η. This is illustrated in Figure 3, where we compare various lower bounds on the negative entropies, either with a fixed kernel, or with kernel learning, showing the benefits of estimating k(ω).

Log-partition minimization. In order to compute (and then minimize) the upper-bound on the logpartition function, we can write

b(f ) = sup p 2r δ=-2r p(δ) * f (δ) -tr C 1/2 p Diag(η)C 1/2 p log(C 1/2 p Diag(η)C 1/2 p ) + r ω=-r η(ω) log η(ω) = inf M log tr exp(M ) + r ω=-r η(ω) log η(ω) s. t. ω-ω ′ =δ η(ω) 1/2 η(ω ′ ) 1/2 M ωω ′ = f (δ), ∀δ{-2r, . . . , 2r}.
Note here that by the classical Caratheodory interpolation theorem (see, e.g., [START_REF] Foias | The Commutant Lifting Approach to Interpolation Problems[END_REF]Chapter 2]), this is exactly a(f ) which is computed.

It can be minimized by (accelerated) projected gradient descent in M , since the function M → log tr exp(M ) is 1-smooth. We simply need to be able to project on the affine subspace, which can be done by solving a linear system. 13 In order to update the kernel η, we use the exact same update as for maximizing the entropy. This is illustrated in Figure 4, where we show benefit of estimating the kernel rather than using a uniform one.

Example 2: Hypercube

We consider X = {-1, 1} d , and ϕ(x) = Diag(η) 1/2 x 1 ∈ R d+1 , for η in the simplex in R d+1 . Here B is exactly the set of all self-adjoint operators with diagonal equal to η. To compute the entropy, we consider the set C of positive matrices with unit diagonal, corresponding to

E x 1 x 1 ⊤ .
13 By Lagrangian duality, the minimizer of To maximize the relative entropy with respect to the uniform distribution, we simply maximize

1 2 M -Z 2 F such that ∀δ ∈ {-r, . . . , r}, tr(M A δ ) = f (δ) obtained as Z + δ λ δ A δ with λ maximizing r δ=-r λ δ ( fδ -tr ZA δ ) -1 2 r δ,δ ′ =-r λ δ λ δ ′ tr A δ A δ ′ .
tr C(p) 1/2 Diag(η)C(p) 1/2 log(C(p) 1/2 Diag(η)C(p) 1/2 ) - d+1 j=1 η j log η j ,
which is a lower bound on x∈{-1,1} d p(x) log p(x) + d log 2. Thus, we get an upperbound on the entropy This is to be compared to the upper-bound of [START_REF] Jordan | Semidefinite relaxations for approximate inference on graphs with cycles[END_REF], which is 1 2 log det C(p) + 1 3

- x∈{-1,1} d p(x) log p(x) d log 2-tr Diag(η) 1/2 C(p) Diag(η) 1/2 log(Diag(η) 1/2 C(p) Diag(η) 1/2 ) + d+1 j=1 η j log η j .
I 0 0 0 + d 2 log πe 2 .
In Figure 5, we compare these bounds on matrices C corresponding to independent components with mean vectors uniformly distributed, showing the potential benefits of the new upperbound on entropy.

Conclusion

In this paper, we presented properties of the von Neumann entropy and relative entropy of covariance operators obtained from reproducing kernel Hilbert spaces. These notions are intimately related to the usual notions of Shannon entropy and relative entropy, share many of their properties, and can be applied to any type of data where positive definite kernels can be defined. They also come with additional computational benefits, with several estimation algorithms with theoretical guarantees. We have also highlighted several properties in terms of multivariate probabilistic modeling or variational inference, but several interesting avenues for future research are worth exploring. First, exploring further approximation guarantees for logpartition function estimation can probably lead to new improved algorithms for variational inference [START_REF] Wainwright | Graphical Models, Exponential Families, and Variational Inference[END_REF] in continuous and discrete domains. Second, we could extend to Rényi entropies [START_REF] Müller-Lennert | On quantum Rényi entropies: A new generalization and some properties[END_REF][START_REF] Frank | Monotonicity of a relative Rényi entropy[END_REF] with similar developments (along the lines of Appendix A.4). Finally, the optimization with respect to the probability measure, which is at the core of log-partition estimation, can be extended to other optimization problems our bound based on quantum entropy with uniform η and the optimized bound). We plot the bounds for independent distributions with random probabilities for each component. We substract the true entropy and compare upper-bounds (with zero meaning that the bounds are tight).

on probability measures, such as already done with similar tools for optimal transport [START_REF] Vacher | A dimension-free computational upper-bound for smooth optimal transport estimation[END_REF] or optimal control [START_REF] Berthier | Infinite-dimensional sums-ofsquares for optimal control[END_REF]. 2 dλ, but there are others, such as, for p ∈ (0, 1]:

is f (t) = -log(t) = 1 -t + ∞ 0 (t -1) 2 λ + t 1 (1 + λ)
-t p = -1 -p(t -1) + sin pπ π +∞ 0 (t -1) 2 λ + t λ p (1 + λ) 2 dλ.
Given A and B two Hermitian operators, we denote by R A the right multiplication by A and L A the left multiplication by A. We have, since L B and R A commute, and log L B = L log B and log R A = R log A :

D(A B) = tr A 1/2 (log A)A 1/2 -tr A 1/2 (log B)A 1/2 = tr (log R A )(A 1/2 )A 1/2 -tr A 1/2 (log L B )(A 1/2 ) = A 1/2 , log R A -log L B A 1/2 = A 1/2 , -log L B R -1 A A 1/2
, where M, N = tr(M * N ) denotes the usual dot-products between operators. Following [START_REF] Beth | Another short and elementary proof of strong subadditivity of quantum entropy[END_REF], we thus have:

D(A B) = A 1/2 , -log L B R -1 A A 1/2 + ∞ 0 (L B R -1 A -I)A 1/2 , (L B R -1 A + λI) -1 (L B R -1 A -I)A 1/2 1 (1 + λ) 2 dλ = tr A -tr B + ∞ 0 R -1/2 A (B -A), (L B R -1 A + tI) -1 R -1/2 A (B -A) 1 (1 + λ) 2 dλ = tr A -tr B + ∞ 0 (B -A), (L B + λR A ) -1 (B -A) 1 (1 + λ) 2 dλ.
Note we can ensure the finiteness of D(A B) as soon as B -1/2 AB -1/2 is bounded. Given that the mapping (u, M ) → u, M -1 u is convex, we obtain the joint convexity of the relative entropy. It can also be shown more directly as shown below.

Proposition 12 Given positive invertible self-adjoint operators A i , B i , then

D n i=1 λ i A i n i=1 λ i B i n i=1 λ i D(A i B i ),
with equality if and only if for all j such that λ j > 0,

log B j -log A j = log n i=1 λ i B i -log n i=1 λ i A i .
Proof Following [START_REF] Beth | Inequalities for quantum entropy: A review with conditions for equality[END_REF]Theorem 7], we have:

D(A i B i )-tr A i (log A-log B) = tr A i (log A i -log A+ log B - log B i ) = tr A i -tr exp(log A -log B + log B i ) + B(A i exp(log A -log B + log B i )), where B(• •) is the Bregman divergence (always non-negative), defined as B(A B) = D(A B) -tr A + tr B. Thus n i=1 λ i D(A i B i ) -D(A B) = n i=1 λ i D(A i B i ) -tr A i (log A -log B) = tr A - n i=1 λ i tr exp(log A -log B + log B i ) + n i=1 λ i B(A i exp(log A -log B + log B i )),
which is greater than zero since the function C → tr exp(K + log C) is concave. There is equality if and only for all i such that λ i > 0, A i = exp(log A -log B + log B i ), because of the Bregman divergences being equal to zero.

More generally, we get for random operators A and B, E D(A B) D(EA EB), with equality if and only if almost surely, log A -log B = log EA -log EB. Note that the concavity of the function B → tr exp(M + log B), which is crucial in deriving concentration inequalities for matrices and operators [START_REF] Tropp | An introduction to matrix concentration inequalities[END_REF], is simply obtained by partial Fenchel dualization of the relative entropy.

We can also provide more refined relationships at equality, following [START_REF] Petz | Sufficient subalgebras and the relative entropy of states of a von Neumann algebra[END_REF].

Proposition 13 Given positive invertible self-adjoint operators A i , B i , and strictly positive real numbers λ i that sum to 1, then

D n i=1 λ i A i n i=1 λ i B i = n i=1 λ i D(A i B i ), implies that for all i, A i B -1 i = AB -1 , where A = n i=1 λ i A i and B = n i=1 λ i B i .
Proof Given the representation above, and the fact that (u, M ) → u, M -1 u is convex, we get the joint convexity of (A, B) → D(A B) is convex. In order to study the equality sign, we can use the lemma below, and for each λ > 0,

L B i + λR A i -1 (B i -A i ) = L B + λR A -1 (B -A)
. By letting λ go to zero, we get the result.

Note that following [START_REF] Petz | Sufficient subalgebras and the relative entropy of states of a von Neumann algebra[END_REF], we can also consider that

(L B i + λR A i ) -1 (B i -A i ) = (L B i R -1 A i + λI) -1 L B i R -1 A i -I I is independent of i,
and is equal to

(L B R -1 A + λI) -1 L B R -1
A -I I. Thus by using representations of a function f : R + → R as

f (t) = (t-1) ∞ 0 dν(λ) t + λ
for a certain signed measure dν, we can deduce that f L B i R -1 A i I is independent of i, and equal to f L B R -1 A I. For example. with f (t) = log t, we get log B i -log A i independent of i. For example. with f (t) = t -1, we get B i A -1 i independent of i. For f (t) = t z , B z i A -z i independent of i, for any complex z.

Lemma 1 If M 1 , . . . , M n are invertible and all λ i are strictly positive and sum to 1, then

n i=1 λ i u ⊤ i M -1 i u i = n i=1 λ i u i n i=1 λ i M i -1 n i=1 λ i u i
if and only if M -1 i u i is independent of i, and equal to M -1 u where M = n i=1 λ i M i and u = n i=1 λ i u i .

Proof We consider the function f (u, M ) = u ⊤ M -1 u, with Taylor expansion

f (u+δ, M +∆) = f (u, M )+2δ ⊤ M -1 u-u ⊤ M -1 ∆M -1 u+2(δ-∆M -1 u) ⊤ M -1 (δ-∆M -1 u)+O( ∆ 3 + δ 3 ).
We now use the following result on the remainder in Jensen's inequality (which itself is a direct consequence of the Taylor expansion with integral remainder taken at X, while expanded at EX):

Ef (X) -f (EX) = 1 0 E (X -EX) ⊤ f ′′ (tX + (1 -t)EX)(X -EX) dt. (15)
In our particular situation, by the invertibility assumption, the matrix M is invertible almost surely, and having equality in Jensen's inequality imposes to have (X -EX) ⊤ f ′′ (tX + (1 -t)EX)(X -EX) = 0 almost surely in t and X. Thus, at t = 0 and X = (u i , M i ), we get:

u i -u -(M i -M )M -1 u = 0, which is exactly M -1 i u i = M -1 u.
A.2 Monotonicity of the relative entropy

Proposition 14 Given operators C i , i = 1, . . . , n, such that n i=1 C * i C i = I, then D n i=1 C i AC * i n i=1 C i BC * i D(A B).
Proof We have, with ω = exp(2iπ/n), using the joint convexity of the relative entropy and the expression

1 n n k=1 ω k(i-j) = 1 i=j : D n i=1 C i AC * i n i=1 C i BC * i n i=1 D(C i AC * i C i BC * i ) = D n i=1 (e i e ⊤ i ) ⊗ E i AE * i n i=1 (e i e ⊤ i ) ⊗ C i BC * i = D 1 n n k=1 ω k(i-j) n i,j=1 (e i e ⊤ j ) ⊗ C i AC * j 1 n n k=1 ω k(i-j) n i,j=1 (e i e ⊤ j ) ⊗ C i BC * j .
Consider the matrices M and N defined by blocks M ij = C i AC * j and N ij = C i BC * j , as well as, the matrix U k defined by blocks (U k ) ij = 1 i=j ω ki I, we get from the equation above:

D n i=1 C i AC * i n i=1 C i BC * i D 1 n n k=1 U k M U * k 1 n n k=1 U k N U * k 1 n n k=1 D U k M U * k U k M U * k = D(M N )
because each U k is unitary. Moreover, we can write M = V AV * and N = V AV * , with a an operator V defined by a column of blocks C i , i = 1, . . . , n, so that

V * V = n i=1 C * i C i = I. Thus D(M N ) = D(A B)
, which shows the desired bound.

Using propositions above, there is is equality if and only: log

(C i BC * i ) -log(C i AC * i ) independent of i and log[U k N U * k ]-log[U k M U * k ] = U k (log N -log M ) U * k independent of k, that is log N -log M block diagonal. Overall, we get equality if log N -log M is block diagonal with equal blocks.
We can also obtain necessary conditions for equality, that is, (C

i BC * i )(C i AC * i ) -1 independent of i, and U k N U * k [U k N U * k ] -1 = U k N M -1 U * k independent of k, that is, N M -1 block-diagonal with equal blocks.
This extends to a random operator C, such that E C * C = I, for which we get:

D E CAC * E CBC * D(A B).
A particular example is obtained from n Hermitian positive operators D 1 , . . . , D n such that n i=1 D i = I, and we get, with µ i = tr(D i A) and ν i = tr(D i B):

n i=1 µ i log µ i ν i D(A B).
This corresponds to, if

D i = k λ ik u ik ⊗ u ik , to C ik = λ 1/2 ik e i u * ik , so that i,k E i AE * i = i,k λ ik e i u * ik M u ik e * i = n i=1
(tr AD i )e i e * i .

A.3 Concavity of relative entropy

In this section, we consider two positive self-adjoint operators, M and N such that M αN for some α.

We consider the function:

f : Λ → D(Λ 1/2 M Λ 1/2 Λ 1/2 N Λ 1/2 ).
It happens to be convex in (M, N ). It is also concave in Λ in the set of positive bounded operators. This is due to the representation, valid for invertible operators:

f (Λ) -tr M Λ + tr N Λ = ∞ 0 Λ 1/2 (N -M )Λ 1/2 , (L Λ 1/2 N Λ 1/2 + λR Λ 1/2 M Λ 1/2 ) -1 Λ 1/2 (N -M )Λ 1/2 1 (1 + λ) 2 dλ = ∞ 0 (N -M ), (L N R -1 Λ + λL -1 Λ R M ) -1 (N -M ) 1 (1 + λ) 2 dλ.
Since the mapping (A, B) → (A -1 + B -1 ) -1 = B -B(A + B) -1 B is matrix-concave [1, Corollary 1.1], we get the desired concavity.

A.4 Extension to other divergences

As mentioned in [START_REF] Beth | Another short and elementary proof of strong subadditivity of quantum entropy[END_REF], both the joint convexity and the monotonicity of the relative entropy D(A B) in Prop. 12 and Prop. 14 can be extended to most operator-convex functions on R * + , and

D f (A B) = A 1/2 , f L B R -1 A A 1/2 ,
with the relative entropy being D -log . Among other interesting cases, we have:

• f (t) = (t -1) 2 , with D f (A B) = tr (B -A) 2 A -1 . • f (t) = (t -1) 2 t , with D f (A B) = tr (B -A) 2 B -1 . • f (t) = (t -1) 2 1 + t , with D f (A B) = A 1/2 , (L B + R A ) -1 A 1/2 . • f (t) = 1 - √ t, with D f (A B) = tr(A) -A 1/2 , [L B R -1 A 1/2 A 1/2 .
Specialized to covariance operators, we have, with A = Σ p and B = Σ q , and following the same reasoning as in the proof of Prop. 4:

D f (Σ p Σ q ) X D f ϕ(x)ϕ(x) * dq dp (x)ϕ(x)ϕ(x) * dp(x) = X f dq dp (x) dp(x) = D f (p q),
which is the reverse f -divergence associated with the function f [START_REF] Csiszár | Information-type measures of difference of probability distributions and indirect observation[END_REF][START_REF] Topsoe | Some inequalities for information divergence and related measures of discrimination[END_REF]. We thus immediately get lower bounds on:

• the square Hellinger distance D f (p q) = X 1- dq dp (x) 1/2 dp(x) = 1 2 X dq dτ (x) 1/2 - dp dτ (x) 1/2 2 dτ (x) for f (t) = 1 - √ t, • the reverse Pearson χ 2 -divergence D f (p q) = X 1- dq dp (x) 2 dp(x) = 1 2 X dq dτ (x)- dp dτ (x) 2 dτ (x) dp/dτ (x) for f (t) = (t -1) 2 , • the Pearson χ 2 -divergence D f (p q) = X 1 - dq dp (x) 2 dp(x) dq/dp(x) = 1 2 X dq dτ (x) - dp dτ (x) 2 dτ (x) dq/dτ (x) for f (t) = (t -1) 2 t , • the Vincze-Le Cam distance D f (p q) = X 1 - dq dp (x) 2 dp(x) 1 + dq/dp(x) equal to 1 2 X dq dτ (x) - dp dτ (x) 2 dτ (x) dq/dτ (x) + dp/dτ (x) for f (t) = (t -1) 2 1 + t .
B Proof of lower bound (Prop. 5)

Proof Starting from p(y) = X h(x, y)dp(x) and q(y) = X h(x, y)dq(x), we can provide an exact expression in the usual data processing inequality for the regular KL divergence [START_REF] Cover | Elements of Information Theory[END_REF], with two joint distributions p(x, y) = p(x)h(x, y) and q(x, y) = q(x)h(x, y), so that, with f (u) = u log u:

D(p q) = D(p q) = X X q(x, y)f p(x, y) q(x, y) dτ (x)dτ (y) X q(y)f X q(x|y) p(x, y) q(x, y) dτ (x) dτ (y) by Jensen's inequality, = X q(y)f p(y) q(y) dτ (y) = D(p q).
Using Eq. ( 15), we thus get

D(p q) -D(p q) = X X q(y)h(x, y) p(x) q(x) - p(y) q(y) 2 1 0 dt t p(x) q(x) + (1 -t) p(y)
q(y) dτ (x)dτ (y).

Since for all x ∈ X, t p(x) q(x) + (1 -t) p(y) q(y) inf x∈X p(x) q(x) , we can then bound as

D(p q) -D(p q) sup x∈X q(x) p(x) X Y q(y)h(x, y) p(x) q(x) - p(y) q(y) 2 dτ (x)dτ (y) = sup x∈X q(x) p(x) X X q(y)h(x, y) p(x) q(x) - p(y) q(y) + p(y) q(y) - p(y) q(y) 2 dτ (x)dτ (y) 2 sup x∈X q(x) p(x) X X q(y)h(x, y) p(x) q(x) - p(y) q(y) 2 dτ (x)dτ (y) +2 sup x∈X q(x) p(x) X q(y) p(y) q(y) - p(y) q(y) 2 dτ (y).
We assume Lipschitz-continuity for p and q in the following form: for all x, y, q(x) q(y) -

1 d(x, y)C q , p (x)
p(y) -1 d(x, y)C p , and p(y) q(y) -p(x) q(x) d(x, y)C p/q for some distance d on X. Then the first term in the bound above is less than 2C 2 p/q sup x∈X q(x) p(x) sup x∈X X h(x, y)d(x, y) 2 dτ (y). For the second term, we write: p(y) q(y) -p(y) q(y) p(y)|q(y) -q(y)| + q(y)|p(y) -p(y)| q(y)q(y) p(y) X |q(y) -q(y)|h(x, y)dτ (x) + q(y) X |p(y) -p(y)|h(x, y)dτ (x) q(y)q(y) C q p(y)q(y) X d(x, y)h(x, y)dτ (x) + C p q(y)p(y) X d(x, y)h(x, y)dτ (x) q(y)q(y)

(C p + C q )p(y) q(y) X d(x, y)h(x, y)dτ (x) (C p + C q )p(y) q(y) X d(x, y) 2 h(x, y)dτ (x) 1/2 . Moreover X p(y) 2 q(y) dτ (y) sup x∈X p(x) q(x) X 2p(y) 2 + 2(p(y) -p(y)) 2 p(y) dτ (y) 2 sup x∈X p(x) q(x) 1 + X (p(y) -p(y)) 2 p(y) dτ (y) 2 sup x∈X p(x) q(x) (1 + C 2 p diam(X) 2 ). Overall, we get D(p q) -D(p q) E(p, q) × sup x∈X h(x, y)d(x, y) 2 dy, with E(p, q) = 2C 2 p/q sup x∈X q(x) p(x) + 4 sup x∈X q(x) p(x) (C p + C q ) 2 sup x∈X p(x) q(x) (1 + C 2 p diam(X) 2 ).

C Convergence of estimation algorithms

In this Appendix, we give convergence proofs for the estimation algorithms from Section 5, starting with concentration inequalities for operators and bounds on degrees of freedom.

Thus,

E|A| P(λ max (B) > u) • max{df(Λ), n} + 1 1 -u E B 2 HS + E λ tr (Σ p + λI) -1/2 B(Σ p + λI) -1/2 .
We have:

P(λ max (B) > u) df(λ) 1 + 3(df max (λ) + u/3) 2 u 4 n 2 exp -nu 2 2df max (λ) + 2u/3 , and
E B 2 HS = 1 n E tr (Σ p + λI) -1 (Σ p -ϕ(x i ) ⊗ ϕ(x i ))(Σ + λI) -1 (Σ p -ϕ(x i ) ⊗ ϕ(x i ))
1 n E tr (Σ p + λI) -1 ϕ(x i ) ⊗ ϕ(x i )(Σ p + λI)

C.3 Degrees of freedom estimation (from projections)

Proposition 16 (Estimation of degrees of freedom (from projections)) Let ϕ : X → H a mapping between a measurable space X and a separable Hilbert space H. Let Σ p = X ϕ(x)ϕ(x) * dp(x), with p a probability distribution and x 1 , . . . , x n ∈ X n i.i.d. samples from q such that ∀x ∈ X, q(x) αp(x), and Π the orthogonal projection on the span of all ϕ(x 1 ), . . . , ϕ(x n ). Let df max (λ) = sup x∈X ϕ(x), (Σ q + λI) -1 ϕ(x) . Then E tr[Σ p (Σ p + λI) -1 ] -tr[ ΠΣ p Π( ΠΣ p Π + λI) -1 ] 16df(αλ)df(µ) exp -9n 32df max (µ) + 8 + 8µ tr (Σ q + αλI) -1 Σ q (Σ q + µI) -1 , if df max (µ) = sup x∈X ϕ(x), (Σ q + µI) -1 ϕ(x) n.

Proof We follow a more direct strategy as the previous proof, following estimation results for random feature expansions and column sampling [START_REF] Rudi | Less is more: Nyström computational regularization[END_REF][START_REF] Bach | On the equivalence between kernel quadrature rules and random feature expansions[END_REF]. 1 1 -u µ(Σ q + µI) -1 , if we have (Σ q + µI) -1/2 (Σ q -Σq )(Σ q + µI) -1/2 uI. This leads to, using Lemma 2, A 2µ 1 -u tr (Σ p + λI) -1 Σ p (Σ q + µI) -1 , with probability greater than 1 -df(µ) 1 + 3(df max (µ) + u/3) 2 u 4 n 2 exp -nu 2 2df max (µ) + 2u/3

. Therefore, we get, with u = 3/4: 0 EA tr(Σ q + αλI) -1 Σ q • df(µ) 1 + 3(df max (µ) + u/3) 2 u 4 n 2 exp -nu 2 2df max (µ) + 2u/3 + 2µ 1 -u tr (Σ q + αλI) -1 Σ q (Σ q + µI) -1 df(αλ)df(µ) 1 + 3(df max (µ) + u/3) 2 u 4 n 2 exp -nu 2 2df max (µ) + 2u/3 + 2µ 1 -u tr (Σ q + αλI) -1 Σ q (Σ q + µI) -1 df(αλ)df(µ) 1 + 3(df max (µ) + 1/4) 2 (3/4) 4 n 2 exp -n(3/4) 2 2df max (µ) + 1/2 +8µ tr (Σ q + αλI) -1 Σ q (Σ q + µI) -1 with u = 3/4.

Using the same reasoning as in the proof of Prop. 15, the probabilistic term is upper-bounded by 16df(αλ)df(µ) exp -9n 32df max (µ) + 8 if df max (µ) n.

C.4 Proof of Proposition 7

Proof We assume that p(x) q(x)/α for all x ∈ X. Given the integral representation of the entropy, we have For λ λ 0 , we have df(λ) df max (λ) = d 1 α C(λ 0 ) = d 0 , which we assume to be less than n. 16df(αλ)df(αελ) exp -9n 32df max (αελ 0 ) + 8 + 8αελ tr (Σ q + αλI) -1 Σ q (Σ q + αελI) 8ελ tr (Σ q + λI) -1 Σ q (Σ q + ελI) 8ελ tr (Σ q + λI) -1 Σ q (Σ q + ελI) -1 dλ, extending integration limits, with df max (αελ 0 ) = n βn . We also have ελ tr (Σ q + λI) -1 Σ q (Σ q + ελI) -1 dλ = σ∈Λ(Σq )

λ ′ 0 σ σ + λ ελ σ + ελ dλ = σ∈Λ(Σq ) εσ 1 -ε λ ′ 0 1 σ + ελ - 1 σ + λ dλ = σ∈Λ(Σq ) εσ 1 -ε 1 ε log σ + ελ ′ σ -log σ + λ ′ σ σ∈Λ(Σq ) σ 1 -ε log (1 + λ ′ εσ -1 ) = 1 1 -ε ελ ′ 0 df(λ)dλ.
Thus, with ε = 1 αλ 0 df -1 (n/β n ),

  and if Σp f = λf with λ = 0 and f = 0 (which implies Kα = 0), then n i=1

Figure 1 :

 1 Figure1: Negative entropy estimation from 20 replications (mean and standard deviation) for two values of the kernel bandwidth σ, as n increases. Left: from i.i.d. samples (method from Section 5.1). Right: from kernel integrals evaluated as uniform samples (method from Section 5.4). As expected we estimate the kernel entropy (limit when n tends to infinity) from above on the left plot, and from below on the right plot, with a faster convergence rate on the right plot.

sup A 0

 0 , tr A=1 tr(AM ) -D(A B) = sup A 0, tr A=1 tr(A(M + log B)) -tr[A log A] = log tr exp(M + log B), with the optimal operator A equal to A = 1 tr exp(M + log B) exp(M + log B). This implies the representation D(A B) = sup M tr AM -log tr exp(M + log B),

sup A 0

 0 tr(AM ) -D(A B) = sup A 0 tr(A(M + log B)) -tr[A log A] = tr exp(M + log B -I), with an optimal A = 1 e exp(M + log B). This implies the representation D(A B) = sup M tr A + tr AM -tr exp(M + log B),

Figure 2 :

 2 Figure 2: Log-partition estimation with different values of σ, as a function of r, the number of considered frequencies.

Figure 3 :

 3 Figure 3: Kernel learning for entropy estimation on the torus [0, 1]. Left: Entropy estimates for increasing order r. Right: optimal kernel for r = 50.

Figure 4 :

 4 Figure 4: Kernel learning for log partition estimation on the torus [0, 1] with increasing order r.

With η constant equal to 1 /

 1 (d+ 1), we get the upper bound: d log 2-1 d + 1 tr C(p) log C(p), which happens to be tight for d = 1.

Figure 5 :

 5 Figure5: Comparison of entropy bounds with increasing dimensions (the log-determinant bound from[START_REF] Jordan | Semidefinite relaxations for approximate inference on graphs with cycles[END_REF], our bound based on quantum entropy with uniform η and the optimized bound). We plot the bounds for independent distributions with random probabilities for each component. We substract the true entropy and compare upper-bounds (with zero meaning that the bounds are tight).

A

  = tr[Σ p (Σ p + λI) -1 ] -tr[ ΠΣ p Π( ΠΣ p Π + λI) -1 ] = λ tr[( ΠΣ p Π + λI) -1 -(Σ p + λI) -1 ] = λ tr[(Σ p + λI) -1 (Σ p -ΠΣ p Π)( ΠΣ p Π + λI) -1 ] using standard expansions of the matrix inverse, = λ tr (Σ p + λI) -1 (Σ p -ΠΣ p Π)(Σ p + λI) -1 +λ tr (Σ p + λI) -1 (Σ p -ΠΣ p Π)( ΠΣ p Π + λI) -1 (Σ p -ΠΣ p Π)(Σ p + λI) -1 .Using Π2 = Π, we can further boundA λ tr (Σ p + λI) -1 (Σ p -ΠΣ p Π)(Σ p + λI) -1 +λ tr (Σ p + λI) -1 (I -Π)Σ p (I -Π) + λ -1 Σ p (I -Π)Σ p (Σ p + λI) -1 using Π2 = Π, = tr (Σ p + λI) -2 (λΣ p -λ ΠΣ p Π + λΣ p -λΣ p Π -λ ΠΣ p + λ ΠΣ p Π + Σ p (I -Π)Σ p ) = 2λ tr (Σ p + λI) -2 Σ p (I -Π) + tr (Σ p + λI) -2 Σ 2 p (I -Π) 2 tr (Σ p + λI) -1 Σ p (I -Π) .We can thus use, with Σq =1 n n i=1 ϕ(x i )ϕ(x i ) * , I -Π = limµ→0 I -Σq ( Σq + µI) -1 = lim µ→0 µ( Σq + µI) -1 µ( Σq + µI) -1

tr) 1 αλ 1 , 1 + λ 0 0 E 1 λ 0 EE

 11010 Σp log Σp -tr Σ p log Σ p = +∞ 0 tr Σ p (Σ p + αλI) -1 -tr Σp ( Σp + αλI) -1 αdλ, it turns out it is possible to truncate large values of λ. Indeed, +∞ λ 1 tr Σ p (Σ p + αλI) -1 -tr Σp ( Σp + αλI) -1 αdλ = tr Σp log( Σp + αλ 1 I) -tr Σ p log(Σ p + αλ 1 I) tr Σp log(I + αλ 1 I) -tr Σ p log(αλ 1 I) = log(1 + αλ 1 ) -log(αλ 1 with a similar bound for the opposite. Thus, we have the bound, with λ 1 λ 0 : EA = tr Σp log Σp -tr Σ p log Σ p 1 αλ df(αλ) -df(αλ) αdλ + λ df(αλ) -df(αλ) αdλ, where df(αλ) = tr Σ p (Σ p + αλI) -1 and df(αλ) = tr Σp ( Σp + αλI) -1 are the usual degrees of freedom for Σ p . We have df(αλ) df max (αλ) sup x∈X ϕ(x), (αΣ + αλI) -1 ϕ(x) = 1 α C(λ),where C(λ) = tr[Σ(Σ + λI)-1 ]. Thus, we get: df(αλ) -df(αλ) αdλ.

C. 5 8 Proof 0 0

 580 Proof of PropositionWe can combine with the proof of Proposition 7, to get:A = tr ΠΣ p Π log ΠΣ p Π -tr Σ p log Σ p = +∞ 0 tr Σ p (Σ p + λI) -1 -tr ΠΣ p Π(ΠΣ p Π + λI) -1 dλ,which is always non-negative and such that, for µ = αελ, if df max (αελ 0 ) n:A λ tr Σ p (Σ p + λI) -1 dλ + λ 1 λ 0 tr Σ p (Σ p + λI) -1 -tr ΠΣ p Π(ΠΣ p Π + λI) -1 dλ +

  1/2 for a certain bounded operator M , and we then have tr A -tr B D(A B) (tr A) log M op .

	The quantity D(A B) -tr A + tr B is always non-negative as the Bregman divergence [15] associated with the convex function A → tr A log A , and may not always be finite. The following properties are classical and are proved from first principles in Appendix A.
	Proposition 2 (Properties of quantum entropy and relative entropy) For Hermitian positive op-
	erators A, B, we have:
	(a) The mappings A → tr[A log A] and (A, B) → D(A||B) are convex.
	(b) We have D(A B) tr A -tr B, with equality if and only if A = B.
	(c) If A i , B i , i = 1, . . . , n, are Hermitian operators and

  -1 ϕ(x i ) ⊗ ϕ(x i ) 1 n df(λ)df max (λ). Moreover E λ tr (Σ p + λI) -1/2 B(Σ p + λI) -1/2 λ E tr (Σ p + λI) -1/2 B(Σ p + λI) -1/2 2Using df(λ) df max (λ) = d, we have the bound:

	This leads to
	E|A|	1 n	df max (λ)df(λ) + + max{df(λ), n}df(λ) 1 + 1 1 n df(λ)df max (λ) 1 -u 3(df max (λ) + u/3) 2 u 4 n 2	exp	-nu 2 2df max (λ) + 2u/3	.
	With u = 3/4, we get:
	E|A|	1 n	df max (λ)df(λ) + 4	1 n	df(λ)df max (λ)
						+ max{df(λ), n}df(λ) 1 +	3(df max (λ) + 1/4) 2 (3/4) 4 n 2	exp	-n(3/4) 2 2df max (λ) + 1/2	.
		E|A|		d √ n	+ 4	d 2 n	+ d max{d, n} 1 +	16(4d + 1) 2 27n 2	exp	-9n 8(4d + 1)	.
	If d n, we further get:
	E|A|	d √ n	+ 4	d 2 n	+ dn 1 +	16(4d + 1) 2 27n 2	exp	-9n 32d + 8	d √ n	+ 4	d 2 n	+ 16dn exp	-9n 32d + 8	.

λ 1 n E ϕ(x), (Σ p + λI) -1 (Σ p + λI) -1 ϕ(x) 2 1 n E ϕ(x), (Σ p + λI) -1 ϕ(x) 2 1 n df max (λ)df(λ).

  We take λ 0 such that d 0 = n log n , with a bound:We take λ 1 = λ 0 + n so that

	Thus						16dn exp		-9n 32d + 8	16dn exp	-9n 32d 0 + 8	,
	and we get the bound, using C(λ) 1 λ :
	A	1 αλ 1	+ nαλ 0 +	0	λ 0	C(λ)dλ +	λ 1 λ 0		d √ n	+ 4	d 2 n	+ 16dn exp	-9n 32d 0 + 8	αdλ
		1 αλ 1	+ nαλ 0 +	0	λ 0	C(λ)dλ +	1 √ n	+ 16n exp	-9n 32d 0 + 8	λ 1 λ 0	C(λ)dλ +	4 nα	λ 0 λ 1	C(λ) 2 dλ.
			A	1 αλ 1	+ nαλ 0 +		0	λ 0	C(λ)dλ +	17 √ n	λ 1 λ 0	C(λ)dλ +	4 nα	λ 0 λ 1	C(λ) 2 dλ
				1 αλ 1	+ nαλ 0 +		0	λ 0	C(λ)dλ +	17 √ n	λ 1 λ 0	C(λ)dλ +	4c nα	.
	The condition on C, implies that C(λ) 2 λ 2				λ λ/2 C(λ ′ ) 2 dλ ′	c, leading to C(λ)	2c λ , leading to
	n log n α	2c λ 0 , and thus λ 0 2c (log n) 2 n 2 α 2 . Then
				A			1 αλ 1	+ 2c	(log n) 2 nα	+ 2c	log n nα	+	17 √ n	λ 1 λ 0	C(λ)dλ +	4c nα	.
	overall, we get:					λ 1 λ 0	C(λ)dλ				0	λ 0 +1	C(λ)dλ + log	λ 0 + n λ 0 + 1	2 √	c 1 + λ 0 + log n. Then,
	A	1 + 8c(log n) 2 nα	+	17 √ n	2 √	c + 2 √	2c	log n nα	+ log n	1 + c(8 log n) 2 nα	+	17 √ n	2 √	c + log n ,
	which is the desired bound.

Note that covariance operators are typically "centered", that is, µpµ * p is subtracted.

In this paper, we use the notation uu * for an element u ∈ H to denote the operator uu * : H → H such that uu * (f ) = f, u u. The usual notation u ⊗ u will be used for tensor products in Section 6.

We can build the simple mapping M → M, ϕ(x)ϕ(x) * obtained from f g * → f g * , ϕ(x)ϕ(x) * = f (x)g(x).

This is a simple consequence of (e1 • t, . . . , em • t) being another orthonormal basis of the same eigensubspace and thus equal to an orthogonal matrix times (e1, . . . , em), so that m i=1 eie * i = m i=1 (ei • t)(ei • t) * , leading to ϕ(x), m i=1 eie * i ϕ(x) constant on X, and thus implying the desired result.

For the non trivial inclusion, we consider A = X ϕ(x)ϕ(x) * dp(x) with X dp(x) = 1, such that A 0. We have for all f ∈ H, X f (x) 2 dp(x) 0. If k is universal, then all positive continuous functions g are limits of functions f 2 , so that X g(x)dp(x) 0, which implies that p is a probability measure.

Note that K invertible implies that the Hadamard (i.e., pointwise) product K • K is invertible, but not vice-versa.

Note here that the smoothing property of h corresponds to the pointwise non-negativity, but that h is also a positive definite kernel.

The distributions p and q are obtained by the same transition kernel, respectively from the distributions p and q.

This section is not crucial to understand the rest of the paper.

With potentially extra terms due to the term tr Σp log Σq .

This section is not crucial to understand the rest of the paper.

Since Π Π * + (I -Π)(I -Π) * = I, the monotonicity of relative entropy (Prop. 14 in Appendix A.2) leads to D(Σp Σq ) D ΠΣp Π + (I -Π)Σp(I -Π) ΠΣq Π + (I -Π)Σq(I -Π) = D ΠΣp Π ΠΣq Π + D (I -Π)Σp(I -Π) (I -Π)Σq(I -Π) , which is greater than D ΠΣp Π ΠΣq Π .

Acknowledgements

The author thanks Alessandro Rudi and Loucas Pillaud-Vivien for fruitful discussions related to this work. This work was funded in part by the French government under management of Agence Nationale de la Recherche as part of the "Investissements d'avenir" program, reference ANR-19-P3IA-0001 (PRAIRIE 3IA Institute). We also acknowledge support the European Research Council (grant SEQUOIA 724063). The comments and suggestions of the anonymous reviewers were greatly appreciated.

A Proof of convexity results on operators

In this section, we provide simple direct proofs of convexity results that we use in this paper. All of these proofs are obtained from the literature, and in particular [START_REF] Lesniewski | Monotone Riemannian metrics and relative entropy on noncommutative probability spaces[END_REF][START_REF] Beth | Another short and elementary proof of strong subadditivity of quantum entropy[END_REF][START_REF] Beth | Inequalities for quantum entropy: A review with conditions for equality[END_REF].

A.1 Joint convexity of the relative entropy

We consider a function f : R * + → R, which is operator-convex, that is, for which the function A → f (A) is convex as an operator-valued function, i.e., f n i=1 λ i A i n i=1 λ i f (A i ). These functions can exactly be represented as:

, where γ, δ > 0 and ν is a positive measure with finite mass on R + [START_REF] Lesniewski | Monotone Riemannian metrics and relative entropy on noncommutative probability spaces[END_REF]. We have f (1) = α. The interesting case for this paper

C.1 Concentration of sums of random self-adjoint operators

We recall in this Appendix classical results on the concentation of covariance operators.

Lemma 2 (Concentration of covariance operators [START_REF] Rudi | Generalization properties of learning with random features[END_REF][START_REF] Tropp | An introduction to matrix concentration inequalities[END_REF]) Let ϕ be mapping ϕ : X → H between a measurable space X and a separable Hilbert space H, such that for all x ∈ X, ϕ(x) R. Given

, assume that λ max (Σ p ) µ. Given x 1 , . . . , x n ∈ X n i.i.d. samples from p, then, for all u > 0,

.

We will apply the lemma above to φ(x) = (Σ p + λI) -1/2 ϕ(x), with Σ p = X ϕ(x) ⊗ ϕ(x)dp(x) and

, and tr Σ p is now tr Σ p (Σ p + λI) -1 = df(λ), and µ = 1, leading to:

C.2 Degrees of freedom estimation

Proposition 15 (Estimation of degrees of freedom) Let ϕ : X → H a mapping between a measurable space X and a separable Hilbert space H. Let Σ p = X ϕ(x)ϕ(x) * dp(x), with p a probability distribution and x 1 , . . . , x n ∈ X n i.i.d. samples from p, and Σp = 1

where df(λ) = tr Σ p (Σ p + λI) -1 , and df max (λ) = sup x∈X ϕ(x), (Σ p + λI) -1 ϕ(x) , which is assumed to satisfy d = df max (λ) n.

Proof We follow closely the proof of [START_REF] Rudi | Generalization properties of learning with random features[END_REF], which starts with:

, as soon as λ max (B) < 1, for

Almost surely, |A| max{df(Λ), n}. Thus, for u < 1,

Polynomial decays. If df(λ) Bλ -1/s , then we have:

, and thus ε

Moreover,

Now turning to asymptotic expansions, we get

We can take β n ∝ log n so that the term exp(-β n /4) is negligible.

We can choose λ 0 such that Bλ

, and then λ (3s-1)/(2s) 1

, finally leading to A = O B s(s+1)/(3s-1) log n n s(s-1)/(3s-1)

.

Exponential decays. If df(λ)

B log(C/λ) d , then we have:

), the first and third terms lead to

With λ 1 ∝ 1 B exp(β n /8), the second term is of order

We can then choose

C.6 Computation of degrees of freedom

We have, when σ tends to zero, with ν

where we have used that, for ν e,

(1 -e -σ )

and