In this section we sketch our main result by first recalling the construction of the KPCA estimator introduced in Hoffmann [16]. We have at disposal a training set {x1,…,xn}∈D⊂ℝd of n points independently sampled according to some probability distribution P. The input space D is a known compact subset of ℝ^d, but the probability distribution P is unknown and the goal is to estimate the support C of P from the empirical data. We recall that C is the smallest closed subset of D such that ℙ[C] = 1 and we stress that C is in general a proper subset of D, possibly of low dimension.

Classical Principal Component Analysis (PCA) is based on the construction of the vector space V spanned by the first m eigenvectors associated with the largest eigenvalues of the empirical covariance matrix

where x¯=1n∑i=1nxi is the empirical mean. However, if the data do not live on an affine subspace, the set V is not a consistent estimator of the support. In order to take into account non-linear models, following the idea introduced in Schölkopf et al. [17] we consider a feature map Φ from the input space D into the corresponding feature space H, which is assumed to be a Hilbert space, and we replace the empirical covariance matrix with the empirical covariance operator

where μ^n=1n∑i=1nΦ(xi) is the empirical mean in the feature space. As it happens in PCA, we consider the subspace V^m,n of H spanned by the first m- eigenvectors f^1,…,f^m of Tnc^. According to the proposal in Hoffmann [16], we consider the following estimator of the support of the probability distribution P

where τ_n is a suitable threshold depending on the number of examples and

is the projection of an arbitrary point x∈D onto the affine subspace μ^n+V^m,n. We show that, under a suitable assumption on the feature map, called separating property, C^n is a consistent estimator of C with respect to the Hausdorff distance between compact sets, see Theorem 3 of section 4.2.

The separating property was introduced in De Vito et al. [25] and it ensures that the feature space is rich enough to learn any closed subset of D. This assumption plays the same role of the notion of universal kernel [26] in supervised learning.

Moreover, following [25, 27] we extend the KPCA estimator to a class of learning algorithms defined in terms of a low-pass filter function r_m(σ) acting on the spectrum of the covariance matrix and depending on a regularization parameter m ∈ ℕ. The projection of Φ^n(x) onto V^m,n is replaced by the vector

where {f^j}j is the family of eigenvectors of Tnc^ and {σ^j}j is the corresponding family of eigenvalues. The support is then estimated by the set

Note that KPCA corresponds to the choice of the hard-cut off filter

However, other filter functions can be considered, inspired by the theory of regularization for inverse problems [28] and by supervised learning algorithms [29, 30]. In this paper we show that the explicit computation of these spectral estimators reduces to a finite dimensional problem depending only on the kernel K(x, w) = 〈Φ(x), Φ(w)〉 associated with the feature map, as for KPCA. The computational properties of each learning algorithm depend on the choice of the low-pass filter r_m(σ), which can be tuned to out-perform of some specific data set, see the discussion in Rudi et al. [31].

We conclude this section with two considerations. First, in De Vito et al. [25, 27] it is proven a consistency result for a similar estimator, where the subspace V^n,m is computed with respect to the non-centered covariance matrix in the feature space H, instead of the covariance matrix. In this paper we analyze the impact of recentering the data in the feature space H on the support estimation problem, see Theorem 1 below. This point of view is further analyzed in Rudi et al. [32, 33].

Finally note that, our consistency results are based on convergence rates of empirical subspaces to true subspaces of the covariance operator, see Theorem 2 below. The main difference between our result and the result in Blanchard et al. [34], is that we prove the consistency for the case when the dimension m = m_n of the subspace V^m,n goes to infinity slowly enough. On the contrary, in their seminal paper [34] the authors analyze the most specific case when the dimension of the projection space is fixed.

The statistical model is described by a random vector X taking value in D. We denote by P the probability distribution of X, defined on the Borel σ-algebra of D, and by C the support of P.

Since the probability distribution P is unknown, so is its support. We aim to estimate C from a training set of empirical data, which are described by a family X₁, …, X_n of random vectors, which are independent and identically distributed as X. More precisely, we are looking for a closed subset C^n=C^X1,…,Xn⊂D, depending only on X₁, …, X_n, but independent of P, such that

for all probability distributions P. In the context of regression estimate, the above convergence is usually called universal strong consistency [35].

To define the estimator C^n we first map the data into a suitable feature space, so that the support C is represented by a linear subspace.

Assumption 1. Given a Hilbert space H, take Φ:D→H satisfying the following properties:

The space H is called the feature space and the map Φ is called a Mercer feature map.

In the following the norm and scalar product of H are denoted by ||·|| and 〈··〉, respectively.

Assumptions (H1) and (H2) are standard for kernel methods, see Steinwart and Christmann [36]. We now briefly recall some basic consequences. First of all, the map K:D×D→ℝ

is a Mercer kernel and we denote by HK the corresponding (separable) reproducing kernel Hilbert space, whose elements are continuous functions on D. Moreover, each element f∈H defines a function fΦ∈HK by setting f_Φ(x) = 〈f, Φ(x)〉 for all x∈D. Since Φ(D) is total in H, the linear map f ↦ f_Φ is an isometry from H onto HK. In the following, with slight abuse of notation, we write f instead of f_Φ, so that the elements f∈H are viewed as functions on D satisfying the reproducing property

Finally, since D is compact and Φ is continuous, it holds that

Following De Vito et al. [27], we call Φ a separating Mercer feature map if the following the separating property also holds true.

It states that any closed subset C⊂D is mapped by Φ onto the intersection of Φ(D) and the closed subspace span¯{Φ(x)∣x∈C}⊂H. Examples of kernels satisfying the separating property are for D⊂ℝd [27]:

Sobolev kernels with smoothness index s>d2;

As shown in De Vito et al. [25], given a closed set C the equality (2) is equivalent to the condition that for every x0∉C there exists f∈H such that (3)f(x0)≠0  and  f(x)=0 ∀x∈C.

Clearly, an arbitrary Mercer feature map is not able to separate all the closed subsets, but only few of them. To better describe these sets, we introduce the elementary learnable sets, namely

where f∈H. Clearly, Cf is closed and the equality (3) holds true. Furthermore the intersection of an arbitrary family of elementary learnable sets ∩f∈FCf with F⊂H satisfies (3), too. Conversely, if C is a set satisfying (2), select a maximal family {f_j}_j∈J of orthonormal functions in H such that

i.e., a basis of the orthogonal complement of {Φ(x)∣x∈C}, then it is easy to prove that

so that any set which is separating by Φ is the (possibly denumerable) intersection of elementary sets. Assumption (H3) is hence a requirement that the family of the elementary learnable sets, labeled by the elements of H, is rich enough to parameterize all the closed subsets of D by means of (4). In section 4.3 we present some examples.

The Gaussian kernel K(x, w) = e^{−γ|x−w|2} is a popular choice in machine learning, however it is not separating. Indeed, since K is analytic, the elements of the corresponding reproducing kernel Hilbert space are analytic functions, too [36]. It is known that, given an analytic function f ≠ 0, the corresponding elementary learnable set Cf={x∈D∣f(x)=0} is a closed set whose interior is the empty set. Hence also the denumerable intersections have empty interior, so that K can not separate a support with not-empty interior. In Figure 1 we compare the decay behavior of the eigenvalues of the Laplacian and the Gaussian kernels.

The second building block is a low pass filter, we introduce to avoid that the estimator overfits the empirical data. The filter functions were first introduced in the context of inverse problem, see Engl et al. [28] and references therein, and in the context of supervised learning, see Lo Gerfo et al. [29] and Blanchard and Mucke [30].

We now fix some notations. For any f∈H, we denote by f ⊗ f the rank one operator (f ⊗ f)g = 〈gf〉f. We recall that a bounded operator A on H is a Hilbert-Schmidt operator if for some (any) basis {f_j}_j the series ||A||22:=∑j||Afj||2 is finite, ||A||₂ is called the Hilbert-Schmidt norm and ||A||_∞ ≤ ||A||₂, where ||·||_∞ is the spectral norm. We denote by S2 the space of Hilbert-Schmidt operators, which is a separable Hilbert space under the scalar product 〈AB〉2=∑j〈AfjBfj〉.

Assumption 2. A filter function is a sequence of functions r_m:[0, R] → [0, 1], with m ∈ ℕ, satisfying

For fixed m, r_m is a filter cutting the smallest eigenvalues (high frequencies). Indeed, (H4) and (H6) with σ′ = 0 give

On the contrary, if m goes to infinity, by (H5) r_m converges point-wisely to the Heaviside function

Since r_m(σ) converges to Θ(σ), which does not satisfy (5), we have that limm→+∞Lm=+∞.

We fix the interval [0, R] as domain of the filter functions r_m since the eigenvalues of the operators we are interested belong to [0, R], see (23).

Examples of filter functions are

We recall a technical result, which is based on functional calculus for compact operators. If A is a positive Hilbert-Schmidt operator, Hilbert-Schmidt theorem (for compact self-adjoint operators) gives that there exist a basis {f_j}_j of H and a family {σ_j}_j of positive numbers such that

If the spectral norm ||A||_∞ ≤ R, then all the eigenvalues σ_j belong to [0, R] and the spectral calculus defines r_m(A) as the operator on H given by

With this definition each f_j is still an eigenvector of r_m(A), but the corresponding eigenvalue is shrunk to r_m(σ_j). Proposition 1 in the Appendix in Supplementary Material summarizes the main properties of r_m(A).

As anticipated in the introduction, the estimators we propose are a generalization of KPCA suggested by Hoffmann [16] in the context of novelty detection. In our framework this corresponds to the hard cut-off filter, i.e., by labeling the different eigenvalues of A in a decreasing order¹ σ₁ > σ₂ > … > σ_m > σ_m+1 > …, the filter function is

Clearly, r_m satisfies (H4) and (H5), but (H6) does not hold. However, the Lipschitz assumption is needed only to prove the bound (21e) and, for the hard cut-off filter, r_m(A) is simply the orthogonal projector onto the linear space spanned by the eigenvectors whose eigenvalues are bigger than σ_m+1. For such projections [37] proves the following bound

so that (21e) holds true with Lm=2σm+1-σm. Hence, our results also hold for hard cut-off filter at the price to have a Lipschitz constant L_m depending on the eigenvalues of A.

The third building block is made of the eigenvectors of the distribution dependent covariance operator and of its empirical version. The covariance operators are computed by first mapping the data in the feature space H.

As usual, we introduce two random variables Φ(X) and Φ(X) ⊗ Φ(X), taking value in H and in S2, respectively. Since Φ is continuous and X belongs to the compact subset D, both random variables are bounded. We set

where the integrals are in the Bochner sense.

We denote by μ^n and T^n the empirical mean of Φ(X) and Φ(X) ⊗ Φ(X), respectively, and by Tnc^ the empirical covariance operator, respectively. Explicitly,

The main properties of the covariance operator and its empirical version are summarized in Proposition 2 in the Appendix in Supplementary Material.

First of all, we characterize the support of P by means of Q^c, the orthogonal projector onto the null space of the distribution dependence covariance operator T^c. The following theorem will show that the centered feature map

sends the support C onto the intersection of Φc(D) and the closed subspace (I-Qc)H, i.e,

Theorem 1. Assume that Φ is a separating Mercer feature map, then

where Q^c is the orthogonal projector onto the null space of the covariance operator T^c.

Proof. To prove the result we need some technical lemmas, we state and prove in the Appendix in Supplementary Material. Assume first that x∈D is such that Q^c(Φ(x) − μ) = 0. Denoted by Q the orthogonal projection onto the null space of T, by Lemma 2 QQ^c = Q and Qμ = 0, so that

Hence Lemma 1 implies that x ∈ C.

Conversely, if x ∈ C, then as above Q(Φ(x) − μ) = 0. By Lemma 2 we have that Q^c(1 − Q) = ||Q^cμ||⁻²Q^c μ ⊗ Q^cμ. Hence it is enough to prove that

which holds true by Lemma 3.

Our first result is about the convergence of F^n.

Theorem 2. Assume that Φ is a Mercer feature map. Take the sequence {m_n}_n such that

for some constant κ > 0, then

Proof. We first prove Equation (13). Set An=I-rmn(Tnc^). Given x∈D,

where the fourth line is due to (21e), the bound ||An||∞=supσ∈[0,1]|1-rm(σ)|≤1, and the fact that both Φ(x) and μ are bounded by R. By (12b) it follows that

so that, taking into account (24a) and (24c), it holds that

almost surely, provided that lim_n→+∞∥(r_mn(T^c) − (I − Q^c))(Φ(x) − μ)∥ = 0. This last limit is a consequence of (21d) observing that {Φ(x)-μ∣x∈D} is compact since D is compact and Φ is continuous.

We add some comments. Theorem 1 suggests that the consistency depends on the fact that the vector (I-rmn(T^nc))(Φ(x)-μ^n) is a good approximation of Q^c(Φ(x)−μ). By the law of large numbers, T^n and μ^n converge to T and μ, respectively, and Equation (21d) implies that, if m is large enough, (I − r_m(T))(Φ(x)−μ) is closed to Q^c(Φ(x)−μ). Hence, if m_n is large enough, see condition (12a), we expect that rmn(T^nc) is close to rmn(Tc). However, this is true only if m_n goes to infinity slowly enough, see condition (12b). The rate depends on the behavior of the Lipschitz constant L_m, which goes to infinity if m goes to infinity. For example, for Tikhonov filter a sufficient condition is that mn~n12-ϵ with ϵ > 0. With the right choice of m_n, the empirical decision function F^n converges uniformly to the function F(x) = Q^c(Φ(x) − μ), see Equation (13).

If the map Φ is separating, Theorem 1 gives that the zero level set of F is precisely the support C. However, if C is not learnable by Φ, i.e., the equality (2) does not hold, then the zero level set of F is bigger than C. For example, if D is connected, C has not-empty interior and Φ is the feature map associated with the Gaussian kernel, it is possible to prove that F is an analytic function, which is zero on an open set, hence it is zero on the whole space D. We note that, in real applications the difference between Gaussian and Abel kernel, which is separating, is not so big and in our experience the Gaussian kernel provides a reasonable estimator.

From now on we assume that Φ is separating, so that Theorem 1 holds true. However, the uniform convergence of F^n to F does not imply that the zero level sets of F^n converges to C = F⁻¹(0) with respect to the Hausdorff distance. For example, with the Tikhonov filter F^n-1(0) is always the empty set. To overcome the problem, C^n is defined as the τ_n-neighborhood of the zero level set of F^n, where the threshold τ_m goes to zero slowly enough.

Define the data dependent parameter τ^n as

Since F^n∈[0,1], clearly τ^n∈[0,1] and the set estimator becomes

The following result shows that C^n is a universal strongly consistent estimator of the support of the probability distribution P. Note that for KPCA the consistency is not universal since the choice of m_n depends on some a-priori information about the decay of the eigenvalues of the covariance operator T^c, which depends on P.

Theorem 3. Assume that Φ is a separating Mercer feature map. Take the sequence {m_n}_n satisfying (12a)-(12b) and define τ^n by (14). Then

Proof. We first show Equation (15a). Set F(x) = Q^c(Φ(x) − μ) and let E be the event on which F^n converges uniformly to F(x), and F be the event such that X_i ∈ C for all i ≥ 1. Theorem 2 shows that ℙ[E] = 1 and, since C is the support, then ℙ[F] = 1. Take ω ∈ E ∩ F and fix ϵ > 0, then there exists n₀ > 0 (possibly depending on ω and ϵ) such that for all n ≥ n₀ |F^n(x)-F(x)|≤ϵ for all x∈D. By Theorem 1 F(x) = 0 for all x ∈ C and X₁(ω), …, X_n(ω) ∈ C, it follows that |F^n(Xi(ω))|≤ϵ for all 1 ≤ i ≤ n so that 0≤τ^n(ω)≤ϵ, so that the sequence τ^n(ω) goes to zero. Since ℙ[E ∩ F] = 1 Equation (15a) holds true.

We split the proof of Equation (15b) in two steps. We first show that with probability one limn→+∞supx∈C^nd(x,C)=0. On the event E ∩ F, suppose, by contraction, that the sequence {supx∈C^nd(x,C)}n does not converge to zero. Possibly passing to a subsequence, for all n ∈ ℕ there exists xn∈C^n such that d(xn,C^n)≥ϵ0 for some fixed ϵ₀ > 0. Since D is compact, possibly passing to a subsequence, {x_n}_n converges to x0∈D with d(x₀, C) ≥ ϵ₀. We claim that x₀ ∈ C. Indeed,

since xn∈C^n means that ||F^n(xn)||≤τn. If n goes to infinity, since Φ is continuous and by the definition of E and F, the right side of the above inequality goes to zero, so that ||Qc(Φ(x0)-μ)||=0, i.e., by Theorem 1 we get x₀ ∈ C, which is a contraction since by construction d(x₀, C) ≥ ϵ₀ > 0.

We now prove that

For any x∈D, set X_1,n(x) to be a first neighbor of x in the training set {X₁, …, X_n}. It is known that for all x ∈ C,

see for example Lemma 6.1 of Györfi et al. [35].

Choose a denumerable family {_zj}_{j ∈ J} in C such that is dense in C. By Equation (16) there exists an event G with such that ℙ[G] = 1 and, on G, for all j ∈ J

Fix ω ∈ G, we claim that limnsupx∈Cd(x,C^n)=0. Observe that, by definition of τ^n, Xi∈C^n for all 1 ≤ i ≤ n and

so that it is enough to show that limn→+∞supx∈Cd(X1,n(x),x)=0.

Fix ϵ > 0. Since C is compact, there is a finite subset J_ϵ ⊂ J such that {B(z_j, ϵ)}_{j ∈ Jϵ} is a finite covering of C. Furthermore,

Indeed, fix x ∈ C, there exists an index j ∈ J_ϵ such that x ∈ B(z_j, ϵ). By definition of first neighbor, clearly

so that by triangular inequality we get

Taking the supremum over C we get the claim. Since ω ∈ G and J_ϵ is finite,

so that by Equation (17)

Since ϵ is arbitrary, we get limn→+∞supx∈Cd(X1,n(x),x)=0, which implies that

Theorem 3 is an asymptotic result. Up to now, we are not able to provide finite sample bounds on d_H(Ĉ_n, C). It is possible to have finite sample bounds on ||F^n(x)-Qc(Φ(x)-μ)||, as in Theorem 7 of De Vito et al. [25] with the same kind of proof.

The following two examples clarify the notion of the separating condition.

Example 1. Let D be a compact subset of ℝ², H=ℝ6 with the euclidean scalar product, and Φ:D→ℝ6 be the feature map

whose corresponding Mercer kernel is a polynomial kernel of degree two, explicitly given by

Given a vector f=(f1,…,f6)⊤, the corresponding elementary set is the conic

Conversely, all the conics are elementary sets. The family of all the intersections of at most five conics, i.e., the sets whose cartesian equation is a system of the form

where f₁₁, …, f₅₆ ∈ ℝ.

Example 2. The data are the random vectors in ℝ²

where a, c ∈ ℕ, b, d ∈ [0, 2π] and Θ₁, …, Θ_n are independent random variables, each of them uniformly distributed on [0, 2π]. Setting D=[-1,1]2, clearly Xi∈D and the support of their common probability distribution is the Lissajous curve

Figure 2 shows two examples of Lissajous curves. As a filter function r_m, we fix the hard cut-off filter where m is the number of eigenvectors corresponding to the highest eigenvalues we keep. We denote by Ĉnm,τ the corresponding estimator given by (10).

In the first two tests we use the polynomial kernel (18), so that the elementary learnable sets are conics. One can check that the rank of T^c is less or equal than 5. More precisely, if Lis_a,b,c,d is a conic, the rank of T^c is 4 and we need to estimate five parameters, whereas if Lis_a,b,c,d is not a conic, Lis_a,b,c,d is not a learnable set and the rank of T^c is 5.

In the first test the data are sampled from the distribution supported on the circumference Lis1,0,1,π2 (see panel left of Figure 2). In Figure 3 we draw the set Ĉnm,τ for different values of m and τ when n varies. In this toy example n = 5 is enough to learn exactly the support, hence for each n = 2, …, 5 the corresponding values of m_n and τ_n are m_n = 1, 2, 3, 4 and τ_n = 0.01, 0.005, 0.005, 0.002.

In the second test the data are sampled from the distribution supported on the curve Lis_2,0.11,1,0.3, which is not a conic (see panel right of Figure 2). In Figure 4 we draw the set Ĉnm,τ for n = 10, 20, 50, 100, m = 4, and τ = 0.01. Clearly, C^n is not able to estimate Lis_2,0.11,1,0.3.

In the third test, we use the Abel kernel with the data sampled from the distribution supported on the curve Lis_2,0.11,1,0.3 (see panel right of Figure 2). In Figure 5 we show the set Ĉnm,τ for n = 20, 50, 100, 500, m = 5, 20, 30, 50, and τ = 0.4, 0.35, 0.3, 0.2. According the fact that the kernel is separating, C^n is able to estimate Lis_2,0.11,1,0.3 correctly.

We now briefly discuss how to select the parameter m_n and τ_n from the data. The goal of set-learning problem is to recover the support of the probability distribution generating the data by using the given input observations. Since no output is present, set-learning belongs to the category of unsupervised learning problems, for which there is not a general framework accounting for model selection. However there are some possible strategies (whose analysis is out of the scope of this paper). A first approach, we used in our simulations, is based on the monotonicity properties of Ĉnm,τ with respect to m, τ. More precisely, given f ∈ (0, 1), we select (the smallest) m and (the biggest) τ such that at most nf observed points belong to the the estimated set. It is possible to prove that this method is consistent when f tends to 1 as the number of observations increases. Another way to select the parameters consists in transforming the set-learning problem is a supervised one and then performing standard model selection techniques like cross validation. In particular set-learning can be casted in a classification problem by associating the observed example to the class +1 and by defining an auxiliary measure μ (e.g., uniform on a ball of interest in D) associated to −1, from which n i.i.d. points are drawn. It is possible to prove that this last method is consistent when μ(suppρ) = 0.

We now explain the role of the filter function. Given a training set X₁, …, X_n of size n, the separating property (3) applied to the support of the empirical distribution gives that

where μ^n is the empirical mean and I−Qnc^ is the orthogonal projection onto the linear space spanned by the family {Φ(X1)-μ^n,…,Φ(Xn)-μ^n}, which are the centered images of the examples. Hence, given a new point x ∈ D the condition |Qnc^(Φ(x)−μ^n)|| ≤ τ with τ ≪ 1 is satisfied if only if x is close to one of the examples in the training set. Hence the naive estimator {x∈D|‖Qnc^(Φ(x)−μn^)‖≤τ} overfits the data. Hence we would like to replace Qnc^ with an operator, which should be close to the identity on the linear subspace spanned by {Φ(X1)-μ^n,…,Φ(Xn)-μ^n} and it should have a small range. To modulate the two requests, one can consider the following optimization problem

We note that if A is a projection its Hilbert-Schmidt norm ||A||₂ is the square root of the dimension of the range of A. Since

where Tr(A) is the trace, A^⊤ is the transpose and ||A||2=Tr(A⊤A) is the Hilbert-Schmidt norm, then

and the optimal solution is given by

i.e., A_opt is precisely the operator rm(Tnc^) with the Tikhonov filter rm(σ)=σσ+λ and λ=Rm. A different choice of the filter function r_m corresponds to a different regularization of the least-square problem