Polynomial volume estimation and its applications

Polynomial volume estimation and its applications Antonio Cuevas, Beatriz Pateiro-López Version: Accepted Manuscript HOW TO CITE Antonio Cuevas, Beatriz Pateiro-López (2018) Polynomial volume estimation and its applications, Journal of Statistical Planning and Inference, Volume 196, pp 174-184, DOI: 10.1016/j.jspi.2017.11.005. FUNDING Research has been funded by project INNPAR2D (reference code MTM2016-76969- P), from the Ministry of Economy and Competitiveness and the European Regional Development Fund (ERDF). © 2017 Elsevier B.V. All Rights reserved. This manuscript version is made available under the CC-BY- NC-ND 4.0 license (http://creativecommons.org/licenses/by-ncnd/4.0/) Polynomial volume estimation and its applications Antonio Cuevasaand Beatriz Pateiro-L´opezb aDepartamento de Matem´aticas, Universidad Aut´onoma de Madrid, Spain bDepartamento de Estad´ıstica, An´alisis Matem´atico y Optimizaci´on, Universidad de Santiago de Compostela, Spain Abstract Given a compact set S⊂Rdwe consider the problem of estimating, from a random sample of points, the Lebesgue measure of S,µ(S), and its boundary measure, L(S) (as defined by the Minkowski content of ∂S). This topic has received some attention, especially in the two-dimensional case d= 2, motivated by applications in image analysis. A new method to simultaneously estimate µ(S) and L(S) from a sample of points inside Sis proposed. The basic idea is to assume that Shas a polynomial volume, that is, that V(r) := µ{x: d(x, S)≤r}is a polynomial in rof degree d, for all rin some interval [0, R). We develop a minimum distance approach to estimate the coefficients of V(r) and, in particular µ(S) and L(S), which correspond, respectively, to the independent term and the first degree coefficient of V(r). The strong consistency of the proposed estimators is proved. Some numerical illustrations are given. Keywords: Set estimation; volume estimation; boundary length estimation 1 Introduction The general background. Set estimation The theory of set estimation is closely linked to nonparametric statistics and stochastic geometry. The general goal of this theory is to estimate a compact set S⊂Rdfrom a random sample of points; see, e.g., Cuevas (2009) for a short overview. Some relevant applications appear in different areas, including ecology [estimation of the habitat of a species or the home range of typical individuals; see Getz and Wilmers (2004), Kie et al. (2010)], econometrics [estimation of the efficient boundary in productivity analysis; see Simar and Wilson (2000)], image analysis (Willett and Novak, 2007; Jang, 2006), nonparametric quality control (Ba´ıllo and Cuevas, 2006), and clustering (Rinaldo and Wasserman, 2010). 1 In this setting, a natural aim is the estimation of some functionals of S, in particular the volume and boundary measure of S. Some notations and basic definitions Let us consider the d-dimensional Euclidean space Rd, equipped with the usual inner product h·,·i and the corresponding norm k·k. Given a set A⊂Rd, we denote by Ac, int(A) and ∂A the complement, interior and boundary of A, respectively. We denote by B(x, r) the closed ball with centre xand radius r. Also, for any compact set A⊂Rdwe will denote (with a slight abuse of notation) by B(A, ε) the closed ε-neighbourhood, or ε-parallel set, of A,B(A, ε) = {x∈Rd:d(x, A)≤ε}, where d(a, C) = inf {ka−ck:c∈C}. If Aand C are non-empty compact subsets of Rd, the Hausdorff distance between Aand Cis defined by dH(A, C) = inf{ε > 0 : A⊂B(C, ε), C ⊂B(A, ε)}.Denote by µ(S) the d-dimensional Lebesgue measure of S. Thus µ(S) is the volume of Sfor d= 3 and the area for d= 2. When no confusion is possible we will sometimes use these terms even for the general case S⊂Rd. The boundary measure of S(i.e., the “perimeter” or the “surface area” of S) is often defined in terms of the Minkowski content, L0(S), or its one-sided version, L(S), given by L0(S) = lim ε→0 µ(B(∂S, ε)) 2ε,or L(S) = lim ε→0 µ(B(S, ε)\S) ε,(1) provided that these limits do exist and are finite. Typically, the values L0(S) and L(S) coincide for regular enough sets; see Ambrosio, Colesanti and Villa (2008) for details and additional references on the Minkowski content. In what follows, we will mostly use L(S). Let νbe a Borel measure on Rd. Let A, C be Borel sets with finite ν-measure. The (pseudo) distance in measure between Aand Cis given by dν(A, C) = ν(A∆C),where ∆ denotes the symmetric difference between Aand C, that is, A∆C= (A\C)∪(C\A). We often use either the Lebesgue measure µor a probability measure in the role of ν. Statement of the problem Suppose we have a random sample Xn={X1, . . . , Xn}of iid observations from a random variable Xwith absolutely continuous probability distribution PX≡Pand compact support S⊂Rd. We want to estimate the volume (Lebesgue measure) of S,µ(S), and the surface measure, L(S), as defined in (1). 2 A brief overview of boundary measure and volume estimation To gain some perspective, and for comparison purposes with our proposal here, we list here some (mostly recent) contributions on the topic of volume and boundary measure estimation. These contributions can be organized according to different criteria depending on i) the assumed model to generate the sample observations, ii) the functional (volume or surface area) to be estimated, iii) the type of estimator; in some cases the estimation is undertaken in a plug-in fashion, as a by-product of a set estimator ˆ Sof Sso that µ(S) and L(S) are estimated by µ(ˆ S) and L(ˆ S). In other cases direct estimators of the volume and boundary measure are constructed, not relying on any previous estimation of S. Whatever the chosen approach some restrictions must be imposed on the set S. It is clear that the family of all compact sets with a finite boundary measure is huge and the task of estimating both µ(S) and L(S) from a finite sample of points seems hopeless unless some additional shape conditions are imposed. We now summarize some contributions on the topic according to the criteria i)-iii) listed above. First, let us consider the references dealing with estimation under convexity-type restrictions on the basis of the “inside model”, i.e., all the sample points X1, . . . , Xnare taken inside the target set Saccording to a distribution with support S. If Sis assumed to be convex, then the natural estimator for µ(S) is µ(Sn), where Sndenotes the convex hull of the sample points. Some deep results concerning convergence rates and asymptotic distribution for this estimator can be found in Brunel (2016) and Pardon (2011). Of course the estimator µ(Sn) is typically biased since in general Sn(S. The unbiased estimation of µ(S) when Sis convex is addressed in Baldin and Reiss (2016). Although the assumption of convexity for Sis natural and appealing from different points of view, it is also quite restrictive for many practical applications. This has motivated the use of several extensions of the notion of convex set. One of them is α-convexity: a closed set S is said to be α-convex when it can be expressed as the intersection of the complements of a family of open balls of radius α. This definition is clearly inspired in the characterization of 3 a closed convex set as an intersection of closed half-spaces; see Cuevas et al. (2012) for some background and references. Assuming that Sis α-convex the volume µ(Sn,α) of the α-convex hull of the sample provides a natural, biased, estimator for µ(S); see Rodr´ıguez-Casal (2007). An improved (bias corrected) version of this estimator, has been proposed by Arias-Castro et al. (2016). It achieves the minimax convergence rate under the regularity assumption that both Sand Scare α-convex. Regarding the estimation of the perimeter L(S), let us mention the case in which S⊂ R2satisfies the above mentioned α-convexity condition for a given α > 0. Now, a natural estimator of L(S) from an inside sample, is the corresponding perimeter of the α-convex hull of the sample, L(Sn,α). In Cuevas et al. (2012, Th. 6) it is proved, under mild conditions, that in this α-convex bivariate case we have L(Sn,α)→L(S), almost surely (a.s.), as n→ ∞. The non-trivial computational aspects can be dealt with using the R-package alphahull by Pateiro-L´opez and Rodr´ıguez-Casal (2010). Other related interesting ideas on the estimation of the perimeter, relying on the use of the so-called α-shape, are analyzed in Arias-Castro and Rodr´ıguez-Casal (2016). Now, let us focus on the references dealing with the estimation of the boundary measure under the following “inside-outside” model: assume that the target set Sfulfils S⊂(0,1)d. Under the “inside-outside” model we have independent identically distributed (iid) observations (X1,IS(X1)),...,(Xn,IS(Xn)) of a random variable (X, IS(X)) where Xis uniformly distributed on [0,1]dand ISstands for the indicator function of S. Thus, under this model we also have sample data outside Sand we assume that for each Xiwe know IS(Xi), that is, we are able to decide whether Xi∈Sor Xi∈Sc. A plug-in type estimator of the boundary Minkowski content, see Cuevas et al. (2007); Armend´ariz et al. (2009), is: Ln=µ(Tn(εn)) 2εn ,with Tn={z∈[0,1]d:∃Xi∈B(z, εn)∩S, and Xj∈B(z, εn)∩Sc}. Ak-NN version of this idea can be found in Cuevas et al. (2013). A different estimator based on Delaunay triangulations has been proposed by Jim´enez and Yukich (2011). In fact, the technique proposed by these authors allows for the estimation 4 of more general surface integrals under quite general shape restrictions. Still, the considered sampling model requires to have data points inside and outside the target set S. Practical motivations The perimeter L(S) and the area µ(S) are, obviously, basic functionals of primary interest in the analysis of a set S⊂R2. The so-called “compactness index” CI(S) = Perimeter(S)2 Area(S)=L(S)2 µ(S) has been used in shape analysis (with applications in medicine). See Montero and Bribiesca (2009) for a survey. Roughly speaking, CI(S) measures how “irregular” the shape of Sis. The square root of CI(S) has been sometimes called “contour index”; see Canzonieri and Carbone (1998) for an application in oncology. Even if Sis completely known, we might want to have good estimators of L(S) and µ(S) based on Monte Carlo samples. Thus, the estimators of L(S) and µ(S) provide a sort of stochastic algorithms to approximate these quantities as well as the compactness index CI(S). 2 The assumption of polynomial volume: its geometric meaning As we have indicated in the introduction, when only an inside sample is available, most usual estimators of µ(S) and L(S) use a geometric shape condition on S(such as convexity or αconvexity) which is incorporated to the estimator Sn, via the “hull-principle” (that is, Snis defined as the minimal set fulfilling the imposed condition). Then µ(S) and L(S) are estimated by µ(Sn) and L(Sn), which entails additional problems for the practical evaluation of these quantities. We will follow here a different strategy: first we will assume a condition on Sexpressed in “analytic” or “algebraic” terms: it simply consists on imposing that the volume of the parallel set V(r) = µ(B(S, r)) is a polynomial on some domain [0, R). Second, we will see that the coefficients of such polynomial have a direct interpretation in terms of the parameters of interest, µ(S) and L(S). Finally, the target parameters will be obtained by minimizing the distance between the function V(r) and a consistent estimator Vn(r) of this function. So we start by defining our crucial assumption. 5 Definition 1. A compact set S⊂Rdis said to fulfil the polynomial volume property if there exist constants θ0, . . . , θd∈Rand R > 0such that µ(B(S, r)) = θ0+θ1r+. . . +θdrd,for all r∈[0, R).(2) This condition has been recently employed in a statistical context by Berrendero et al. (2014). However these authors use a sample model quite different to that considered here. Geometric aspects of the polynomial volume assumption First, it is clear that under condition (2), the Lebesgue measure of Sis µ(S) = θ0and the one-sided Minkowski content (1) of Sis L(S) = θ1. Second, we should mention that assumption (2) is closely related to the concept of convexity. The classical Steiner’s theorem in convex geometry establishes that any compact convex set satisfies condition (2) for all r∈[0,∞); see e.g. Morvan (2008, Ch. 16). In Heveling et al. (2004) it is proved that, for the two-dimensional case d= 2, the validity of (2) in [0,∞) is in fact equivalent to the convexity of S. These authors also give counterexamples to prove that such equivalence does not hold for d= 3. Third, and foremost, there is a class of sets, much wider than that of convex sets, which satisfies condition (2). This is the class of sets with positive reach. This notion, first introduced by Federer (1959), has a clear intuitive meaning. It is formally defined as follows. Let Unp(S) be the set of points x∈Rdwith a unique metric projection, denoted by ξS(x), on S. This means that, for x∈Unp(S), ξS(x) is the unique point fulfilling d(x, S) = kx−ξS(x)k. For x∈S, let reach(S, x) = sup{r > 0 : int(B(x, r)) ⊂Unp(S). Then, the reach of Sis defined by r0:= reach(S) = inf reach(S, x) : x∈S. We will say that Sis a set with positive reach if reach(S)>0. The condition that Shas a positive reach is clearly a sort of smoothness assumption on S, expressed in purely geometric terms, with no direct use of differentiability properties. In simple, informal terms, the assumption r0>0 rules out the existence of sharp inwards peaks in the boundary of S. Positive reach is a well-known extension of the crucial notion of convexity. Indeed it is not difficult to show that for a closed set S, reach(S) = ∞if and only if Sis convex. In the pioneering paper by Federer (1959), the author proved the following result, establishing the relation between positive reach and polynomial volume. 6 Theorem (Federer, 1959, Th. 5.6 and Th. 5.19). Let S⊂Rdbe a compact set with r0= reach(S)>0. Then there exist unique values Φ0(S),...,Φd(S)such that µ(B(S, r)) = d X i=0 rd−iωd−iΦi(S),for 0≤r < r0,(3) where ω0= 1 and, for j≥1,ωjdenotes the j-dimensional measure of the unit ball in Rj. Furthermore, Φ0(S)coincides with the Euler characteristic of S. As a consequence, if Sis a compact set with positive reach, then Φd(S) = µ(S) and the one-sided Minkowski content, L(S), in (1) always exists and corresponds to the coefficient of the first-degree term in (3). Federer’s result is in fact much deeper than stated here, since in particular the Φi(S) have an interpretation in terms of curvature measures. For our purposes, besides the independent term (= µ(S)) and the coefficient of the first order term (= L(S)), it is important that the term Φ0(S) in (3) equals the Euler characteristic of S. This is an integer valued topological invariant (i.e., it is preserved by homeomorphisms). For example, if Sis a closed ball, Φ0(S) = 1, and the same holds for any other compact set Shomeomorphic to the closed ball. In practice, this means that the highest order coefficient can be sometimes assumed to be known in (3). This amounts to impose a further, not too restrictive, geometric condition. In particular, the additional restriction θd=µ(B(0,1)) in (2) is meaningful and useful in many cases. As a further appealing feature of the polynomial volume assumption, let us mention that many simple interesting sets not-fulfilling the positive reach property do in fact satisfy (2). This is the case of the “pac-man” figure obtained in R2by erasing in the unit ball all the points in the first quadrant (see Figure 1), or two tangent spheres in R3, and many others; see Heveling et al. (2004) for details. 3 Estimation of µ(S)and L(S)under the polynomial volume assumption In this section we establish the theoretical basis of our method. In particular, it is shown how to consistently estimate the coefficients θ0, . . . , θdin (2). 7 Figure 1: The set in gray has not positive reach but satisfies the polynomial volume property (2). 3.1 An auxiliary result The following elementary result is mentioned here just as an auxiliary tool in the proof of Theorem 1. It establishes that the Hausdorff convergence of Snto Sentails the convergence in measure for the corresponding parallel sets, that is, dµ(B(Sn, r), B(S, r)) →0 for r > 0. No further assumption is needed. In intuitive terms, this means that the parallel set is always “regular enough” to ensure the convergence in measure. In particular, note that if we take the sample itself, i.e., Sn=Xn, the assumption dH(Sn, S)→0 is obviously true, with probability one, so that the estimation in measure of B(S, r) can be done even using the simplest possible estimator, B(Xn, r) = ∪n i=1B(Xi, r); see Devroye and Wise (1980), Kim and Korostelev (2000), Cuevas and Rodr´ıguez-Casal (2004). Lemma 1. Let Sbe a compact non-empty set in Rd. Let Snbe a sequence of compact sets such that dH(Sn, S)→0. Then, for r > 0, dµ(B(Sn, r), B(S, r)) →0.(4) Proof. We have to prove that µ(B(Sn, r)∆B(S, r)) →0. Given ε∈(0, r/2] we have, for n large enough, µ(B(Sn, r)\B(S, r)) ≤µ(B(Sn, r)\B(S, r −ε)) ≤µ(B(S, r +ε)\B(S, r −ε)),(5) since by the Hausdorff convergence, for nlarge enough, Sn⊂B(S, ε) and, therefore, B(Sn, r)⊂ B(S, r +ε). Second, we have also eventually µ(B(S, r)\B(Sn, r)) ≤µ(B(S, r)\B(Sn, r −ε)) ≤µ(B(S, r)\B(S, r −2ε)),(6) since, again by the Hausdorff convergence, for nlarge enough, S⊂B(Sn, ε) and, therefore, B(S, r −2ε)⊂B(Sn, r −ε). Now, the result follows directly from (5) and (6), together with 8 |ˆ θ0−θ0|/θ0|ˆ θ1−θ1|/θ1|ˆ θ2−θ2|/θ2 Sample size Mean SD Mean SD Mean SD n= 500 0.06673 0.00825 0.01544 0.00336 0.00203 0.00076 n= 1000 0.04242 0.00494 0.00972 0.00194 0.00129 0.00046 n= 5000 0.01460 0.00135 0.00333 0.00052 0.00044 0.00015 Table 6: Mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for S=B(0,1), i= 0,1,2. For each sample, the estimations ˆ θiare obtained using the interpolation method from Vn(r1), Vn(r2) and Vn(r3), where r1, r2and r3are randomly selected in the interval [1,2]. |ˆ θ0−θ0|/θ0|ˆ θ1−θ1|/θ1 Sample size Mean SD Mean SD n= 500 0.15705 0.01577 0.00411 0.00206 n= 1000 0.09960 0.00903 0.00262 0.00122 n= 5000 0.03385 0.00257 0.00086 0.00040 Table 7: Mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for S=B(0,5) \int(B(0,4)), i= 0,1 (assumming θ2= 0 known). For each sample, the estimations ˆ θiare obtained using the interpolation method from Vn(r1) and Vn(r2), where r1and r2are randomly selected in the interval [2,3.5]. 4.2 Some three-dimensional examples: the ball and the torus We now show the practical performance of our method in the three-dimensional Euclidean space. We have considered two different sets S⊂R3, that satisfy the polynomial volume assumption (2). The estimations of the coefficients θi,i= 0,...,3, were obtained using the interpolation method described in (11) since, as commented before, it is computationally faster and the results obtained in the bidimensional case with this methodology were comparable to those obtained with the methodology of Theorem 1. Thus, let us first consider the unit ball S=B(0,1) in R3. Note that Ssatisfies (2) for all r∈R. In Table 8 we show, for different sample sizes, the mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 samples from the uniform distribution on S,i= 0,...,3. The 15 estimations ˆ θiare obtained by solving the system of equations in (11) where rj,j= 1,...,4, are randomly selected in the interval [10,15]. We can compare the results in Table 8 with those obtained if we estimate θ0and θ1with ˆ θ0=µ(Sn) and ˆ θ1=L(Sn), where Snis the convex hull of the sample and µ(Sn) and L(Sn) its volume and surface area, respectively, see Figure 2. Results are summarized in Table 9. Figure 2: Left, S=B(0,1) in R3. Middle, uniform sample of size n= 100 in Sand convex hull of the sample in green. Right, in red, Sn(r) for a given value of r. |ˆ θ0−θ0|/θ0|ˆ θ1−θ1|/θ1|ˆ θ2−θ2|/θ2|ˆ θ3−θ3|/θ3 Sample size Mean SD Mean SD Mean SD Mean SD n= 2000 0.13073 0.00505 0.05617 0.00264 0.02290 0.00119 0.00000 0.00000 n= 5000 0.08381 0.00264 0.03578 0.00133 0.01453 0.00060 0.00000 0.00000 n= 10000 0.05974 0.00153 0.02537 0.00078 0.01027 0.00034 0.00000 0.00000 Table 8: Mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for S=B(0,1) in R3, i= 0,...,3. For each sample, the estimations ˆ θiare obtained using the interpolation method from Vn(rj), with j= 1,...,4, where rjare randomly selected in the interval [10,15]. As an example of non-convex set in R3, we consider a torus Swith major radius 5 and minor radius 1, see Figure 3. In Table 10, we show the mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for S(for simplicity we only show the results corresponding 16 |ˆ θ0−θ0|/θ0|ˆ θ1−θ1|/θ1 Sample size Mean SD Mean SD n= 2000 0.09688 0.00426 0.05642 0.00264 n= 5000 0.06213 0.00218 0.03594 0.00133 n= 10000 0.04422 0.00129 0.02549 0.00078 Table 9: Mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for S=B(0,1) in R3,i= 0,1. Now ˆ θ0=µ(Sn) and ˆ θ1=L(Sn), where Snstands for the convex hull of the sample and µ(Sn) and L(Sn) denote its volume and surface area, respectively. to ˆ θ0and ˆ θ1). Again, the estimations ˆ θiare obtained by solving the system of equations in (11) where now rjare randomly selected in the interval [3,4]. It is important to note that, unlike the two-dimensional non-convex case, we cannot compare our results in Table 10 to any other procedure based on a more specific information on the set of interest (we are not aware of any implementation that supports the computation of the α-convex hull of the sample for d= 3). Figure 3: Left, torus Swith major radius 5 and minor radius 1. Middle, uniform sample of size n= 100 in S. Right, in red, Sn(r) for a given valur of r. 4.3 Possible extensions of the method Lemma 1 establishes that the Hausdorff convergence of Snto Sentails the convergence in measure for the corresponding parallel sets. Since Sn=Xnis dH-consistent, our proposal is 17 |ˆ θ0−θ0|/θ0|ˆ θ1−θ1|/θ1 Sample size Mean SD Mean SD n= 2000 0.48674 3.49394 0.14970 1.41482 n= 5000 0.21870 0.19926 0.06006 0.09791 n= 10000 0.14743 0.03044 0.03932 0.01578 Table 10: Mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for a torus Swith major radius 5 and minor radius 1, i= 0,1. For each sample, the estimations ˆ θiare obtained using the interpolation method from Vn(rj), with j= 1,...,4, where rjare randomly selected in the interval [3,4]. to estimate V(r) using this very simple estimator. But of course, the volume of the parallel set might be also approximated using more sophisticated estimators Sn. For example, if we assume in advance that Sis convex (or α-convex), this information can be incorporated by estimating the “true” polynomial function V(r) of the target set Sin terms of the volume functions of the sample convex hull (or the sample α-convex hull). These estimators of V(r) should be, in principle, more accurate, as they carry additional information on the set S. However, the computation of the parallel volumes for such sets is far from simple, even in dimension d= 2, and especially for the α-convex hull, whose structure is considerably involved. So, in those cases the direct plug-in estimators considered in Tables 2 and 9 look as a quite natural alternative. If no shape restriction is assumed, one could also think of estimating Sby using Sn= B(Xn, εn), where εnis a sequence of smoothing parameters, see Devroye and Wise (1980). In that case, V(r) might be estimated from Vn(r) = µ(B(Xn, εn+r)). To illustrate these ideas, we consider a regular polytope. Let S= [0,1]2be the unit square in R2. For different sample sizes, we generate B= 500 samples from the uniform distribution on S. For each sample, we applied the “interpolation method” described in Subsection 4.1 with Vn(rj) computed as µ(B(Sn, rj)), where Snis (a) the sample Xn, (b) the convex hull of the sample and (c) the Devroye-Wise estimator B(Xn, εn) for different values of the smoothing parameter εn. Results are summarized in Table 11. The best results are obtained using the Devroye-Wise estimator, for which the corresponding calculations can be made quite efficiently (at least for dimensions 2 and 3). Note that, however, the results in case (c) depend largely on an adequate choice of the smoothing parameter εn. The dilated convex hull (b) slightly outperforms our 18 proposal (as a drawback, it is harder to compute in higher dimensions). The overall conclusions suggest a good performance of our estimator, in terms of both, economy of assumptions and computational simplicity. The results of the comparison of (a) and (c) in the unit cube S= [0,1]3⊂R3can be found in the Supplementary material document. We also include in the Supplementary material a brief simulation study with non uniform samples. |ˆ θ0−θ0|/θ0|ˆ θ1−θ1|/θ1|ˆ θ2−θ2|/θ2 Sample size Mean SD Mean SD Mean SD n= 500 (a) 0.05758 0.00782 0.04363 0.01136 0.00085 0.00032 (b) 0.03197 0.00872 0.04704 0.01112 0.00016 0.00042 (c) εn= 0.005 0.03833 0.00794 0.03580 0.01136 0.00085 0.00031 εn= 0.01 0.01892 0.00805 0.02799 0.01132 0.00084 0.00031 εn= 0.04 0.10083 0.00889 0.01949 0.01050 0.00080 0.00029 n= 1000 (a) 0.03612 0.00485 0.03123 0.00786 0.00057 0.00019 (b) 0.01842 0.00669 0.03356 0.00811 0.00020 0.00081 (c) εn= 0.005 0.01664 0.00493 0.02339 0.00786 0.00057 0.00019 εn= 0.01 0.00472 0.00344 0.01563 0.00771 0.00056 0.00019 εn= 0.04 0.12414 0.00558 0.03148 0.00786 0.00053 0.00018 n= 5000 (a) 0.01219 0.00133 0.01436 0.00347 0.00021 0.00007 (b) 0.00506 0.00449 0.01539 0.00380 0.00018 0.00006 (c) εn= 0.005 0.00761 0.00136 0.00655 0.00341 0.00021 0.00006 εn= 0.01 0.02757 0.00139 0.00306 0.00212 0.00021 0.00006 εn= 0.04 0.15062 0.00165 0.04843 0.00347 0.00020 0.00006 Table 11: Mean value and standard deviation of |ˆ θi−θi|/θiover B= 500 repeats for S= [0,1]2,i= 0,1,2. For each sample, the estimations ˆ θiare obtained using the interpolation method from Vn(r1), Vn(r2) and Vn(r3), where r1, r2and r3are randomly selected in the interval [1,2]. Vn(rj) is computed as µ(B(Sn, rj)), where Sn is (a) the sample Xn, (b) the convex hull of the sample and (c) the Devroye-Wise estimator B(Xn, εn), for different values of εn. 19 4.4 Some technical details All computations were carried out in R, (R Core Team (2017)). In the computation of the volume of the union of a family of 3D balls, we have used the Structural Bioinformatics Library (SBL), a C++/Python API by Cazals and Dreyfus (2016). The convex hull of the sample in R3, as well as its volume and surface area is computed using the R-package geometry by Habel et al. (2015). For the dilated convex hull in R2we used the R-package polyclip by Johnson and Baddeley (2017). 5 Discussion The assumption of polynomial volume, as stated in Definition 1 is reasonably general for practical purposes, including applications in image analysis. As discussed above, condition (2) not only applies to the broad class of sets with positive reach but also covers other sets with inward non-smooth peaks in their boundaries. Unlike other methods for boundary measure estimation (see Cuevas et al. (2007), Jim´enez and Yukich (2011), among others) the polynomial volume (PV) method outlined here only requires a sample inside the set S⊂Rd: no external sample points are needed. Also, no auxiliary set estimator Snof Sis required, except for the sample itself. Likewise, the PV method can perform the simultaneous estimation of µ(S) and L(S) and we do not require uniformity over Sfor the distribution of the sample points. The role of the parameter Rin our PV method is worth of some comments. Such parameter concerns the set S, that is, the target of the estimation. Hence it should be seen as a regularity parameter (or a smoothness parameter) very much in the same way as the number of derivatives or the H¨older exponent that we assume in an underlying density function to be estimated. Therefore Rshould not be confused with the tuning or smoothing parameters commonly found in nonparametrics (e.g. in density estimation). The crucial difference lies in the fact that such smoothing or tuning parameters appear only in the expressions of the estimators: they usually measure the extent at which such estimators are “close” to the data and are typically chosen to asymptotically achieve efficiency in the estimation. Note that in our case we must also fix the interval [a, b]⊂(0, R) where the minimum distance will be performed but aand bare not properly smoothing parameters to be chosen in an optimal way. 20 R Figure 4: For the set in gray in the left, (2) holds for [0, R). For the set in gray in the right (2) only holds for a small interval [0, R/6), where R/6 is radius of the small indentations in the boundary. If Ris large enough, any choice of bclose to Rand a“not too close” to 0 will do the job. However, it should be noted that the value of Rin assumption (2) is relevant. If Ris too small, the PV method will not work in practice, unless extremely large sample sizes are considered. This is quite intuitive, see Figure 4: if (2) only holds for a small interval [0, R) then, Swill be relatively irregular and, consequently, µ(S) and L(S) will be harder to estimate. The PV method is flexible enough to incorporate additional geometric information on the target set S, whenever such information can be translated in terms of the coefficients of V(r). This is the case, of the Euler Poincar´e characteristic of Swhich, as explained in Section 2, is directly related to the higher order coefficient θd. In the numerical experiments of Section 4 we have shown the effect of incorporating the knowledge of this coefficient to the PV estimation process. In general, the empirical results shown in Section 4 suggest that the method is competitive, even if we compare it with other procedures based on a more specific information on the set S. Again, we should stress the ease of implementation even for case S⊂R3which is very difficult to address with other approaches. The main theoretical open problem connected to the PV method is the study of the asymptotic distribution and convergence rates of the estimates ˆ θj. This is a non-trivial issue which could be perhaps tackled with techniques similar to those employed in the study of minimum distance estimators. 21 Supplementary material. Supplementary material file includes some additional numerical experiments with data in R3and non-uniform samples. Acknowledgements. This work has been partially supported by Spanish Grants MTM2016- 78751-P (A. Cuevas) and MTM2016-76969-P (B. Pateiro-L´opez). References Ambrosio, L., Colesanti, A. and Villa, E. (2008). Outer Minkowski content for some classes of closed sets. Math. Ann. 342, 727–748. Amenta, N., Choi, S., Dey, T.K., Leekha, N. (2002). A simple algorithm for homeomorphic surface reconstruction. Internat. J. Comput. Geom. Appl. 12, 125-141. Arias Castro, E., Rodr´ıguez-Casal, A. (2016). On estimating the perimeter using the alphashape. Ann. Inst. H. Poincar´e Probab. Statist. 53, 1051-1068. Arias-Castro, E., Pateiro-L´opez, B. and Rodr´ıguez-Casal, A. (2016). Minimax estimation of the volume of a set with smooth boundary. Manuscript arXiv:1605.01333v1. Armend´ariz, I., Cuevas, A. and Fraiman, R. (2009). Nonparametric estimation of boundary measures and related functionals: asymptotic results. Adv. in Appl. Probab. 41, 311–322. Ba´ıllo, A. and Cuevas, A. (2006). Parametric versus nonparametric tolerance regions in detection problems. Computational Statistics 21, 523–536. Baldin, N. and Reiss, M. (2016). Unbiased estimation of the volume of a convex body. Manuscript,Stoch. Proc. Appl. 126, 3716–3732. Berrendero, J.R., Cholaquidis, A., Cuevas, A. and Fraiman, R. (2014). A geometrically motivated parametric model in manifold estimation. Statistics 48, 983-1004. Brunel, V. E. (2016). Adaptive estimation of convex and polytopal density support. Probability Theory and Related Fields 164, 1-16. 22 Canzonieri, V. and Carbone, A. (1998). Clinical and biological applications of image analysis in non-Hodgkin’s lymphomas. Hematol. Oncol. 16, 15–28. Cazals, F. and Dreyfus, T. (2016) The Structural Bioinformatics Library: modeling in biomolecular science and beyond. Research Report RR-8957, Inria. http://sbl.inria.fr Cuevas, A. and Rodr´ıguez-Casal, A. (2004). On boundary estimation. Adv. in Appl. Probab. 36, 340–354. Cuevas, A. (2009). Set estimation: Another bridge between statistics and geometry. BEIO 25, 71–85. Cuevas, A., Fraiman, R. and Pateiro-L´opez, B. (2012). On statistical properties of sets fulfilling rolling-type conditions. Adv. in Appl. Probab. 44, 311–329. Cuevas, A., Fraiman, R. and Rodr´ıguez-Casal, A. (2007). A nonparametric approach to the estimation of lengths and surface areas. Ann. Stat. 35, 1031–1051. Cuevas, A., Fraiman, R. and Gy¨orfi, L. (2013) Towards a universally consistent estimator of the Minkowski content. ESAIM Probab. Stat. 17, 359–369. DeVore, R. A. (1986). Approximation of functions. In Approximation Theory, C. de Boor, ed., pp. 1-20. Proceedings of Symposia in Applied Mathematics, 36. AMS Short Course Lecture Notes. American Mathematical Society, Providence, RI. Devroye, L. and Wise, G. (1980). Detection of abnormal behavior via nonparametric, estimation of the support. SIAM J. Appl. Math. 3, 480–488. Federer, H. (1959). Curvature measures. Trans. Amer. Math. Soc. 93, 418–491. Getz, W. M. and Wilmers, C. C. (2004). A local nearest-neighbor convex-hull construction of home ranges and utilization distributions. Ecography 27, 489–505. Habel, K., Grasman, R., Gramacy, R. B., Stahel, A. and Sterratt, D. C. (2015) geometry: Mesh Generation and Surface Tesselation. R package version 0.3-6. https://CRAN.R- project.org/package=geometry 23 Heveling, M., Hug, D. and Last, G. (2004). Does polynomial parallel volume imply convexity? Math. Ann 328, 469–479. Jang, W. (2006). Nonparametric density estimation and clustering in astronomical sky survey. Comp. Statist. & Data Anal. 50, 760–774. Johnson, A. and Baddeley, A. (2017). polyclip: Polygon Clipping. R package version 1.6-1. https://CRAN.R-project.org/package=polyclip Jim´enez, R. and Yukich, J. E. (2011). Nonparametric estimation of surface integrals. Ann. Statist. 39, 232–260. Kie, J.G., Matthiopoulos, J., Fieberg, J. Powell, R. A., Cagnacci, F., Mitchell, M. S., Gaillard, J.-M. and Moorcroft, P. R. (2010). The home-range concept: are traditional estimators still relevant with modern telemetry technology? Phil. Trans. R. Soc. B 365, 2221–2231. Kim, J. C. and Korostel¨ev, A. P. (2000). Estimation of smooth functionals in image models. Math. Methods Statist. 9, 140–159. Korostel¨ev, A. P. and A. B. Tsybakov (1993). Minimax theory of image reconstruction, Volume 82 of Lecture Notes in Statistics. Springer-Verlag, New York. Laha, R. G. and Rohatgi, V. K. (1979). Probability Theory. Wiley, New York. Montero, R.S. and Bribiesca, E. (2009). State of the art of compactness and circularity measures. International Mathematical Forum 4, 1305–1335. Moore, M. (1984). On the estimation of a convex set. Ann. Statist. 12, 1090-1099. Morvan, J. M. (2008). Generalized Curvatures. Springer-Verlag, Berlin Pardon, J. (2011). Central limit theorems for random polygons in an arbitrary convex set. Ann. Probab. 39, 881-903. Pateiro-L´opez, B. and Rodr´ıguez-Casal, A. (2010) Generalizing the convex hull of a sample: The R package alphahull. J. Statist. Softw. 5, 1–28. 24