Moment analysis of the Delaunay tessellation field estimator

Moment Anal ysis of the Dela una y Tessella tion Field Estima tor M.N.M. v an Lieshout CWI & Eindh oven University of T e chnolo gy P .O. Bo x 94079, 109 0 GB Amsterdam, Th e Netherlands Abstract The Campb ell–Meck e theorem is u sed to derive explicit expressions for the mean and v ariance of Schaap and V an de W eygaert’s Delaunay tessellation field estimator. Sp ecial attentio n is paid to Poi sson pro cesses. Keywor ds & Phr ases: Campb ell–Mec ke formula, Delaunay tessellation fi eld estimator, gener- alised weigh t function estimator, in t en sit y function, mass preserv ation, Poisso n p oint pro cess, second order factorial moment measure, second order p ro duct densit y . 2000 Mathematics Subje ct Classific ation: 60G55, 62M30. 1 Preliminaries and notation Let ϕ b e a lo cally fin ite p oint p attern in R d arising as a realisation of s imple p oin t pr o cesses Φ on R d [4, 9]. In practice, d ∈ { 1 , 2 , 3 } . W e shall assume that th e p oin ts are in general quadratic p osition [1 1], that is, (a) no d + 2 p oin ts are lo cated on the b ound ary of a sphere, and (b) in the plane no three p oints are co-linear; in higher dimensions, n o k + 1 p oint s lie in a k − 1 dimensional affine subspace for k = 2 , . . . d . These assump tions are satisfied almost surely for r ealisations of a Poisson pr o cess with lo cally fin ite inte nsit y f u nction λ : R d → [0 , ∞ ) or, more generally , for Gibbs p oint pro cesses defined by their probab ility density with resp ect to suc h a Poisso n p ro cess. An y p oint pattern ϕ giv es r ise to t w o interesti ng tessellatio ns. F irst consider the set C ( x i | ϕ ) := { y ∈ R d : || x i − y || ≤ || x j − y || ∀ x j ∈ ϕ } that consists of all p oin ts in R d that are at least as close to x i ∈ ϕ as to an y other p oint of ϕ , which is called the V or ono i c el l of x i . The en s em b le of all V oronoi cells is the V or on oi tessel lation of ϕ [20]. An equiv alen t defi nition is C ( x i | ϕ ) = \ x j 6 = x i ∈ ϕ H ( x i , x j ) , where H ( x i , x j ) is th e closed h alfspace { y ∈ R d : h y − ( x i + x j ) / 2 , x i − x j i ≥ 0 } consisting of p oin ts th at are at least as close to x i as to x j . In R 1 , for x i < x j , H ( x i , x j ) = ( −∞ , ( x i + x j ) / 2]. In the plane, H ( x i , x j ) is the closed halfplane b ounded by the bisecting line L ( x i , x j ) of the segmen t connecting x i and x j that con tains x i . Note that the V oronoi cells are closed and con v ex, but n ot necessarily b ounded. Under our assump tions, in tersections b et w een k = 2 , . . . , d + 1 different V oronoi cells are either empty or of dim en sion d − k + 1 . In p articular, d +1 \ i =1 C ( x i | ϕ ) 6 = ∅ ⇔ b ( x 1 , . . . , x d +1 ) ∩ ϕ = ∅ where b ( x 1 , . . . , x d +1 ) is th e op en ball spann ed by x 1 , . . . x d +1 on its b oun dary , and in that case is a single p oin t, usu ally referr ed to as a ve rtex of the V oronoi diagram. V ertices can b e used to defi n e th e second tessellation of interest to us in this pap er, the Delaunay tessel lation . In deed, supp ose that ϕ cont ains at least d + 1 p oin ts. Eac h V oron oi v ertex arising as the inte rsection of d + 1 cells C ( x i | ϕ ) defines a closed simp lex, th e conv ex h ull of { x 1 , . . . , x d +1 } , whic h is called a Delaunay c el l [5] and denoted by D ( x 1 , . . . , x d +1 ). Note that for d = 1, Delauna y cells are in terv als, whilst in the plane they form triangles. An alternativ e, equiv alen t, edge based construction is to join p oin ts x 1 , x 2 ∈ ϕ that share a common V oronoi b order C ( x 1 | ϕ ) ∩ C ( x 2 | ϕ ) 6 = ∅ in to a Delaunay edge. In this case, x 1 and x 2 are called V or ono i neighb ours . The set of neigh b ours of x 1 in ϕ is denoted by N ( x 1 | ϕ ). Either w a y , the partition of sp ace formed by the Delauna y cells is referred to as the Delaunay tessel lation . The un ion of Delauna y cells con taining x i ∈ ϕ is kno wn as the c ontiguous V or onoi c el l W ( x i | ϕ ) of x i in ϕ . Figure 1: A set of th irt y p oints w ith their V oronoi (d ash ed lines) and Delauna y (solid lines) tessellatio ns. A con tiguous V oronoi cell is indicated by sh ading. F or more details, including an historical accoun t, th e r eader is referr ed to the compr e- hensiv e textb o oks [11, 12]. An illustration is giv en in Figure 1, w hic h was obtained u sing the DELDIR pack ag e [19]. 2 2 Delauna y tessellation field estimator Recen tly , Sc haap and V an de W eyg aert [14, 15] pr op osed to estimate the in tensit y fun ction of a spatial p oint pro cess by the so-called Delaunay tessel lation field estimator (DTFE). The metho d estimates the inte nsit y at p oints in a realisation r ecipro cal to the v olum e of their contig uous V oronoi cell, and distributes these estimated field v alues o v er Delauna y cells b y linear (or other) interpolation. They also consider int erp olation of fi elds x 7→ f ( x ) ∈ R + observ ed at sampling p oin ts. Earlier suggestions to use V oronoi tessellations f or field in terp olation includ e th ose by Ord [13] and Sibs on [16]. Based on extensiv e simula tions, Schaap and V an de W eygaert claim that, in con trast to k ernel esti mators [1], the DTFE preserv es the total mass of the field and fi ne stru ctur al details, app ears to resu lt in smo oth in terp olation, adap ts itself to the lo cal scale and geometry , an d is r elativ ely r obust. T h e limitations of the metho d lie in its sens itivit y to measurement error, b ound ary effects, and triangular artefacts [15]. Our aim in this pap er is a rigorous analysis of this estimator. Throughout this p ap er, let Φ b e a simple p oin t pro cess on R d ha ving realisatio ns in general quadratic p osition for which the exp ected num b er of p oin ts placed in b ounded Borel sets is finite so that its (first order) moment measure exists as a σ -fi nite Bo rel measure. F urtherm ore, assume that the m oment measur e is absolutely contin uous with resp ect to Leb esgue measure with Radon–Nik od ym d eriv ativ e λ : R d → [0 , ∞ ), its intensity function . Definition 1. Consider a p oint pr o c ess Φ observe d in a c onvex b ounde d Bor el subset A of R d . F or x ∈ Φ ∩ A , define d λ ( x ) := d + 1 | W ( x | Φ ∩ A ) | , (1) wher e | · | denotes d -volume. F or any x 0 ∈ A in the interior of some Delaunay c el l, define \ λ ( x 0 ) := 1 d + 1 X x ∈ Φ ∩ D ( x 0 | Φ ∩ A ) d λ ( x ) (2) as the aver age of the e stimate d intensity function values at the d + 1 vertic es x of the D elaunay c el l D ( x 0 | Φ ∩ A ) c ontaining x 0 . A few remarks are in order. S h ould a particular realisation ϕ of Φ happ en to con tain less than d + 1 p oin ts in A , the intensit y function estimate m ay b e set to zero, or (cf. the Lemma b elo w ) to the num b er of p oin ts divided by | A | . On the sides of the Delauna y cells, an y a v eraging ma y b e used – it is a null set. Finally , \ λ ( x 0 ) is set to zero for p oints that do not fall in an y Delauna y cell. Edge effects arise due to the fact that Φ is not observed, only Φ ∩ A , the Delauna y tessellatio n of which partitions th e con v ex h ull of Φ ∩ A ⊆ A . Su c h effects ma y b e dealt with in man y wa ys. F or example, one migh t u se torus corrections, add arb itrary p oin ts on the b oundary of A (the corners for example in th e generic case of a cub e), or dra w lines orthogonal to the edges emanating from p oin ts on the b oundary of th e con v ex h ull, etc. F urth er examples can b e f ou n d in c hapter 6 of [12]. 3 Lemma 1. (Schaap and V an de W e ygaert [14, 15]) L et ϕ b e a r e alisation of the simple p oint pr o c ess Φ c ontaining at le ast d + 1 p oints i n A . Then the estimator of D e finition 1 pr eserves total mass, that is, Z A \ λ ( x 0 ) dx 0 = n ( ϕ ∩ A ) , the numb er of p oints of ϕ i n A . Pro of: W rite D ( ϕ ∩ A ) for the family of Delauna y cells defined by ϕ ∩ A , and note that Z A \ λ ( x 0 ) dx 0 = X D j ∈D ( φ ∩ A ) | D j |   X x ∈ ϕ ∩ D j 1 | W ( x | ϕ ∩ A ) |   = X x ∈ ϕ ∩ A 1 | W ( x | ϕ ∩ A ) |   X D j ∈D ( φ ∩ A ) 1 { x ∈ D j } | D j |   = n ( ϕ ∩ A ) , cf. [15, p. 62 ff.].  3 Mean and v ariance of the Delauna y tessellatio n field esti- mator In this section, we derive the first tw o moment s of the Delauna y tessellation field estimator. Our first r esult concerns the exp ectatio n. Theorem 1. L et Φ b e observe d in a c onvex b ounde d Bor el subset A , and, for a p oint p attern ϕ with n ( ϕ ∩ A ) ≥ d + 1 in gener al quadr atic p osition, set g ( x 0 | x , ϕ ) := P D j ∈D ( φ ∩ A ) 1 { x 0 ∈ D ◦ j ; x ∈ D j } | W ( x | ϕ ∩ A ) | , (3) for x 0 ∈ A \ ϕ , x ∈ ϕ , and let g ( x | x, ϕ ) := ( d + 1) / | W ( x | ϕ ∩ A ) | if x ∈ φ ∩ A . Then the Delaunay tessel lation field estimator define d by (2) and (1) has exp e c tation E h \ λ ( x 0 ) i = Z A E x [ g ( x 0 | x, Φ)] λ ( x ) dx, wher e E x denotes the exp e ctation with r esp e ct to the Palm distribution of Φ at x . F or patterns ϕ with less th an d + 1 p oints falling in A , it is also p ossible to write \ λ ( x 0 ) = P x ∈ ϕ ∩ A g ( x 0 | x , ϕ ) w ith the fu nction g c hosen to suit the particular typ e of edge correction adopted, see Section 2. Pro of: Note that \ λ ( x 0 ) = X x ∈ Φ ∩ A g ( x 0 | x, Φ) . 4 Hence, by the C ampb ell–Mec ke theorem [18], E \ λ ( x 0 ) = Z A E x [ g ( x 0 | x , Φ)] λ ( x ) dx.  Recall that the second order factorial moment me asur e µ (2) is d efined in integ ral terms b y E   6 = X x 1 ,x 2 ∈ Φ f ( x 1 , x 2 )   = Z Z f ( x 1 , x 2 ) dµ (2) ( x 1 , x 2 ) (4) for any non -n egativ e measur able function f . The s um is o v er all p airs of differen t p oin ts. W e shall say that the second order factorial moment measure exists, if it is lo cally finite. If fur thermore µ (2) is absolutely contin uous with r esp ect to the 2-fold p ro duct measure of Leb esgue measure with itself, a Radon–Nik od y m der iv ativ e exists k n o w n as second order pr o duct density and d enoted b y ρ (2) . In this case, (4) reduces to Z Z f ( x 1 , x 2 ) ρ (2) ( x 1 , x 2 ) dx 1 dx 2 . Theorem 2. L et Φ b e observe d in a c onvex b ounde d Bor el sub set A and define the function g by (3). Assume that se c ond or der pr o duct densities exist. Then the D elaunay tessel lation field estimator define d by (2) and (1) has varianc e V ar( \ λ ( x 0 )) = Z A Z A E (2) x,y [ g ( x 0 | x , Φ) g ( x 0 | y , Φ )] ρ (2) ( x, y ) dx dy + Z A E x  g 2 ( x 0 | x, Φ)  λ ( x ) dx −  Z A E x [ g ( x 0 | x, Φ)] λ ( x ) dx  2 , wher e E (2) x,y denotes the two-fold P alm distribution of Φ . Pro of: Remark that E  \ λ ( x 0 ) 2  = E   6 = X x,y ∈ Φ ∩ A g ( x 0 | x , Φ) g ( x 0 | y , Φ )   + E " X x ∈ Φ ∩ A g 2 ( x 0 | x , Φ) # . The cross term on the r igh t hand side is equal to Z A Z A E (2) x,y [ g ( x 0 | x , Φ) g ( x 0 | y , Φ)] ρ (2) ( x, y ) dx dy , see e.g. [4], where E (2) x,y denotes the tw o-fold P alm distribu tion of Φ [8]. Another app eal to the Campb ell–Mec ke theorem yields E " X x ∈ Φ ∩ A g 2 ( x 0 | x, Φ) # = Z A E x  g 2 ( x 0 | x , Φ)  λ ( x ) dx. 5 Finally , the v ariance is obtained u sing Theorem 1.  In general, the integrals inv olve d in Theorems 1 –2 must b e ev aluate d by numerical or sim ulation metho ds. 4 Comparison to a classic estimator The classic estimator of inte nsit y is the ke rnel estimator \ λ B D ( x 0 ) := n (Φ ∩ b ( x 0 , h ) ∩ A ) | b ( x 0 , h ) ∩ A | , x 0 ∈ A. (5) prop osed b y Berman and Diggle [1]. The estimator can b e regarded as a k ernel estimator [17] with k h ( x 0 | x ) = 1 {|| x − x 0 || < h } / | b ( x 0 , h ) ∩ A | , where b ( x 0 , h ) d en otes th e op en b all around x 0 with radius h > 0. The c hoice of b andwidth h determines th e amount of smo othin g. Note that when th e b ound ed observ ation windo w A 6 = ∅ is op en, on e never divides by zero. In fact, a s tronger statemen t can b e made. The function x 7→ | b ( x, h ) ∩ A | is cont in uous and attains its minim um on the closure ¯ A . Since any p oint on the b ound ary ∂ A alw a ys has a neigh b our within distance h in A , inf x ∈ A | b ( x, h ) ∩ A | > 0. F ur ther details may b e foun d e.g. in [3, 6, 18]. Although (5) is a n atur al estimator, it do es not n ecessarily preserv e th e total mass in A [15], nor is it based on a generalised wei gh t function [17]. It is not hard to mo dify the edge correction in (5) to define an estimator [10] that do es pr eserv e total mass and is based on a w eigh t fun ction. Definition 2. Consider a p oint pr o c ess Φ observe d in an op en b ounde d Bor el subset A of R d . F or x 0 ∈ A , define \ λ K ( x 0 ) := X x ∈ Φ ∩ A 1 {|| x − x 0 || < h } | b ( x , h ) ∩ A | . (6) Lemma 2. The estimator of De finition 2 i s a ge ne r alise d weight function e stimator with kernel k h ( x 0 | x ) = 1 {|| x − x 0 || < h } / | b ( x , h ) ∩ A | that pr eserves total mass, that is, Z A \ λ K ( x 0 ) dx 0 = n (Φ ∩ A ) , the numb er of p oints of Φ in A . Pro of: Note that Z A k h ( x 0 | x ) dx 0 = Z A 1 {|| x − x 0 || < h } | b ( x, h ) ∩ A | dx 0 ≡ 1 for all x ∈ A , that is, \ λ K ( · ) is a generalised w eigh t fun ction estimator. F urtherm ore, f or an y realised p oint pattern ϕ , the restriction ϕ ∩ A in A is finite and Z A   X x ∈ ϕ ∩ A 1 {|| x − x 0 || < h } | b ( x , h ) ∩ A |   dx 0 = X x ∈ ϕ ∩ A Z A 1 {|| x − x 0 || < h } | b ( x , h ) ∩ A | dx 0 = n ( ϕ ∩ A ) . 6  Note that the Delaunay tessellation field estimator is based on an adapt ive kernel (3) as it dep ends on the und er lyin g p oin t p attern. Indeed, for ev ery x ∈ A , Z A g ( x 0 | x , φ ) dx 0 = Z A P D j ∈D ( φ ∩ A ) 1 { x 0 ∈ D ◦ j ; x ∈ D j } | W ( x | ϕ ∩ A ) | dx 0 = 1 . A clear adv an tage is that the problem of choosing the bandw idth is a voided. In order to assess the qualit y of th e estimator, we pr o ceed to compute its m ean and v ariance. Theorem 3. L et Φ b e observe d in a b ounde d op en Bor el subset A . Then, the estimator of Definition 2 has exp e ctation E h \ λ K ( x 0 ) i = Z A 1 { x ∈ b ( x 0 , h ) } | b ( x, h ) ∩ A | λ ( x ) dx. Pro of: By the Camp b ell–Mec ke theorem E h \ λ K ( x 0 ) i = E " X x ∈ Φ ∩ A 1 {|| x − x 0 || < h } | b ( x , h ) ∩ A | # = Z A 1 {|| x − x 0 || < h } | b ( x, h ) ∩ A | λ ( x ) dx.  If w e compare Theorem 3 to Theorem 1, the Palm exp ectation E x [ g ( x 0 | x , Φ] is replaced b y k h ( x 0 | x ), as the latter d o es not dep end on th e p oint pro cess Φ . Theorem 4. let Φ b e observe d in a b ounde d op en Bor el subset A and assume that se c ond or der pr o duct densities exist. Then V ar( \ λ K ( x 0 )) = Z Z ( b ( x 0 ,h ) ∩ A ) 2 ρ (2) ( x, y ) − λ ( x ) λ ( y ) | b ( x, h ) ∩ A | | b ( y , h ) ∩ A | dx dy + Z b ( x 0 ,h ) ∩ A λ ( x ) | b ( x, h ) ∩ A | 2 dx. Pro of: Regarding the second momen t, note that E  \ λ K ( x 0 ) 2  = E    6 = X x,y ∈ Φ ∩ A  1 {|| x − x 0 || < h } | b ( x, h ) ∩ A | 1 {|| y − x 0 || < h } | b ( y , h ) ∩ A |     + E ( X x ∈ Φ ∩ A  1 {|| x − x 0 || < h } | b ( x, h ) ∩ A | 2  ) . Then rewrite th e exp ectations as integ rals w ith resp ect to ρ (2) and λ resp ectiv ely to obtain that the v ariance of \ λ K ( x 0 ) is equal to Z b ( x 0 ,h ) ∩ A Z b ( x 0 ,h ) ∩ A 1 | b ( x, h ) ∩ A | | b ( y , h ) ∩ A | ρ (2) ( x, y ) dx dy + Z b ( x 0 ,h ) ∩ A λ ( x ) | b ( x, h ) ∩ A | 2 dx. 7 An app eal to Theorem 3 completes the pr o of.  The result sh ould b e compared to that of Th eorem 2. Similar argument s as those in the pro ofs of Th eorems 3 and 4 applied to the classic Berman–Diggle estimator (5) giv e mean 1 | b ( x 0 , h ) ∩ A | Z b ( x 0 ,h ) ∩ A λ ( x ) dx and v ariance 1 | b ( x 0 , h ) ∩ A | 2 ( Z b ( x 0 ,h ) ∩ A ) λ ( x ) dx + Z ( b ( x 0 ,h ) ∩ A ) 2 h ρ (2) ( x, y ) − λ ( x ) λ ( y ) i dx dy ) . Note that f or x 0 ∈ A ⊖ b (0 , 2 h ) separated by 2 h from the b oundary of A , no edge correction is necessary , and b oth kernel estimators are ident ical. The disadv ant age of k ernel estimators is that they inv olv e a bandwidth p arameter h ; the larger h , the smo other the estimated intensit y function. F or sp ecific mo dels, h ma y b e c h osen b y optimisation of the (in tegrated) mean squared error [6]. In practice, in a p lanar setting, Diggle [6] recommends to choose h prop ortional to n − 1 / 2 , where n is the observed n um b er of p oin ts. F or a fi xed b an d width, neither the Berman–Diggle estimator n or the mo d ification of Definition 2 is unive rsally b etter. F or examples, the reader is referred to [10]. 5 In tensit y estimation for P oisson p oin t pro cesses In general, the in tegrals in v olve d in Theorems 1–4 ha v e to b e ev aluated n umerically . An exception is the case where Φ is a Poisso n p oin t pro cess with a locally finite in tensit y function. Corollary 1. L et Φ b e a Poisson p oint pr o c ess observe d in a c onvex b ounde d Bor el subset A . Then, E h \ λ ( x 0 ) i = Z A E [ g ( x 0 | x, Φ ∪ { x } )] λ ( x ) dx and V ar( \ λ ( x 0 )) = Z A Z A E [ g ( x 0 | x , Φ ∪ { x, y } ) g ( x 0 | y , Φ ∪ { x, y } )] λ ( x ) λ ( y ) dx dy + Z A E  g 2 ( x 0 | x, Φ ∪ { x } )  λ ( x ) dx −  Z A E [ g ( x 0 | x , Φ ∪ { x } )] λ ( x ) dx  2 . Pro of: F or a P oisson pro cess, the Pa lm d istribution at x is equal to the sup erp osition of its distribution P with an extra p oin t at x , the t w o-fold Palm d istribution P (2) x,y is the su p erp o- sition of P with x and y . F urthermore, ρ (2) ( x, y ) = λ ( x ) λ ( y ) is a pro duct d en sit y . Plugging these results into the expr essions of Theorems 1–2 completes the pr o of.  8 Corollary 2. let Φ b e a Poisson p oint pr o c ess observe d in a b ounde d op en Bor el subset A and assume that se c ond or der pr o duct densities e xist. Then, V ar( \ λ K ( x 0 )) = Z b ( x 0 ,h ) ∩ A λ ( x ) | b ( x, h ) ∩ A | 2 dx. Pro of: Use that ρ (2) ( x, y ) = λ ( x ) λ ( y ) and apply Theorem 4 .  The v ariance of the Berman–Diggle estimator is R b ( x 0 ,h ) ∩ A λ ( x ) dx/ | b ( x 0 , h ) ∩ A | 2 . F or stationary P oisson pro cesses, even more can b e said. In the remainder of this section, define g as in (3) with A = R d . Theorem 5. L et Φ b e a stationary P oisson p oint pr o c ess in R d with intensity λ > 0 . Then, the Delaunay tessel lation field estimator d λ (0) is asymptotic al ly unbiase d. Pro of: Let b ( x, y 1 , . . . , y d ) b e the op en ball spanned by the p oint s x , y 1 , . . . , y d on its top ologica l b ound ary , and let D ◦ ( x, y 1 , . . . , y d ) b e the op en simplex that is th e in terior of the con v ex h ull of { x, y 1 , . . . , y d } . Recall that the p oints x, y 1 , . . . , y d define a V oronoi v ertex, or, equiv alen tly , a Delauna y cell if and only if there are no p oints in b ( x, y 1 , . . . , y d ). By Corollary 1, asymptotically E h d λ (0) i = λ Z R d E [ g (0 | x, Φ ∪ { x } ) ] dx = λ Z E   6 = X { y 1 ,...,y d }⊂ Φ 1 { 0 ∈ D ◦ ( x, y 1 , . . . , y d ); b ( x, y 1 , . . . , y d ) ∩ (Φ ∪ { x } ) = ∅} | W ( x | Φ ∪ { x } ) |   dx = λ Z E   6 = X { y 1 ,...,y d }⊂ Φ 1 { 0 ∈ D ◦ ( x, y 1 , . . . , y d ); b ( x, y 1 , . . . , y d ) ∩ Φ = ∅} | W ( x | Φ ∪ { x } ) |   dx = λ Z E   6 = X { z 1 ,...,z d }⊂ Φ − x 1 {− x ∈ D ◦ (0 , z 1 , . . . , z d ); b (0 , z 1 , . . . , z d ) ∩ Φ − x = ∅} | W (0 | Φ − x ∪ { 0 } ) |   dx = λ Z E   6 = X { z 1 ,...,z d }⊂ Φ 1 {− x ∈ D ◦ (0 , z 1 , . . . , z d ); b (0 , z 1 , . . . , z d ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 } ) |   dx b y stationarit y . Hence, b y F ub ini’s theorem, E h d λ (0) i = λ E   P 6 = { z 1 ,...,z d }⊂ Φ | D ◦ (0 , z 1 , . . . , z d ) | 1 { b (0 , z 1 , . . . , z d ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 } ) |   = λ E  | W (0 | Φ ∪ { 0 } ) | | W (0 | Φ ∪ { 0 } ) |  = λ.  9 The asymptotic v ariance of the Delauna y tessellation field estimator increases quadrati- cally w ith λ with a constant multiplier that dep en ds on the dimension. Th e pro of rests on the follo wing t w o lemmata. Lemma 3. L et Φ b e a stationar y Poisson p oint pr o c ess in R d with intensity λ > 0 . Then, C ( λ, d ) := Z Z E [ g (0 | x, Φ ∪ { x, y } ) g (0 | y , Φ ∪ { x, y } )] λ ( x ) λ ( y ) dx dy = λ 2 Z E 1  | W (0 | Φ ∪ { 0 , x } ) ∩ W ( x | Φ ∪ { 0 , x } ) | | W (0 | Φ ∪ { 0 , x } ) | | W ( x | Φ ∪ { 0 , x } ) |  dx, wher e E 1 denotes exp e ctation with r esp e ct to a u ni t intensity Poisson p oint pr o c ess. By the Nguyen–Zessin form ula [9], alternativ ely C ( λ, d ) = λ E   1 | W (0 | Φ ∪ { 0 } ) | X y ∈N (0 | Φ ∪{ 0 } ) | W (0 | Φ ∪ { 0 } ) ∩ W ( y | Φ ∪ { 0 } ) | | W ( y | Φ ∪ { 0 } ) |   = λ 2 E 1   1 | W (0 | Φ ∪ { 0 } ) | X y ∈N (0 | Φ ∪{ 0 } ) | W (0 | Φ ∪ { 0 } ) ∩ W ( y | Φ ∪ { 0 } ) | | W ( y | Φ ∪ { 0 } ) |   . Pro of: W rite Φ d − 1 for sets of d − 1 distinct p oin ts in Φ. Then, as λ ( x ) ≡ λ is constan t, and g (0 | x, Φ ∪ { x, y } ) g (0 | y , Φ ∪ { x, y } ) v anishes when x and y do not b elong to th e same Delauna y cell contai ning 0 in its in terior, C ( λ, d ) = λ 2 Z Z E   X z ∈ Φ d − 1 1 { 0 ∈ D ◦ ( x, y , z ); b ( x, y , z ) ∩ Φ = ∅} | W ( x | Φ ∪ { x, y } ) | | W ( y | Φ ∪ { x, y } ) |   dx dy = λ 2 Z Z E   X z ∈ Φ − x ; d − 1 1 {− x ∈ D ◦ (0 , y − x, z ); b (0 , y − x, z ) ∩ Φ − x = ∅} | W (0 | Φ − x ∪ { 0 , y − x } ) | | W ( y − x | Φ − x ∪ { 0 , y − x } ) |   dx dy . Because of stationarit y , C ( λ, d ) = λ 2 Z Z E   X z ∈ Φ d − 1 1 {− x ∈ D ◦ (0 , y − x, z ); b (0 , y − x, z ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 , y − x } ) | | W ( y − x | Φ ∪ { 0 , y − x } ) |   dx dy = λ 2 Z Z E " P z ∈ Φ d − 1 1 {− x ∈ D ◦ (0 , y , z ); b (0 , y , z ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 , y } ) | | W ( y | Φ ∪ { 0 , y } ) | # dx dy . Scaling by λ 1 /d yields that λ − 2 C ( λ, d ) is equal to Z Z E " P z ∈ Φ d − 1 1 {− λ 1 /d x ∈ D ◦ (0 , λ 1 /d y , λ 1 /d z ); b (0 , λ 1 /d y , λ 1 /d z ) ∩ λ 1 /d Φ = ∅} λ − 1 | W (0 | λ 1 /d Φ ∪ { 0 , λ 1 /d y } ) | λ − 1 | W ( λ 1 /d y | λ 1 /d Φ ∪ { 0 , λ 1 /d y } ) | # dx dy . 10 Since λ 1 /d Φ is a unit int ensit y Poi sson p oin t pro cess, we obtain λ − 2 C ( λ, d ) = Z Z E 1 " P z ∈ Φ d − 1 1 {− x ∈ D ◦ (0 , y , z ); b (0 , y , z ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 , y } ) | | W ( y | Φ ∪ { 0 , y } ) | # dx dy . An app eal to F ub ini’s theorem to in tegrat e out o v er x completes the pr o of.  Lemma 4. L et Φ b e a stationar y Poisson p oint pr o c ess in R d with intensity λ > 0 . Then, C ′ ( λ, d ) := Z E  g 2 (0 | x, Φ ∪ { x } )  λ ( x ) dx = λ 2 E 1  1 | W (0 | Φ ∪ { 0 } ) |  . wher e E 1 denotes exp e ctation with r esp e ct to a u ni t intensity Poisson p oint pr o c ess. Pro of: Using λ ( x ) = λ and argueing as in the pro of of Theorem 5, we get C ′ ( λ, d ) = λ Z E     6 = X { z 1 ,...,z d }⊂ Φ 1 {− x ∈ D ◦ (0 , z 1 , . . . , z d ); b (0 , z 1 , . . . , z d ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 } ) |   2   dx = λ Z E   6 = X { z 1 ,...,z d }⊂ Φ 1 {− x ∈ D ◦ (0 , z 1 , . . . , z d ); b (0 , z 1 , . . . , z d ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 } ) | 2   dx (7) as − x b elongs to a s ingle Delauna y in terior. W rite Φ d for sets of d distinct p oints in Φ and scale eac h p oint in (7) by λ 1 /d to obtain that C ′ ( λ, d ) is equal to λ Z E   X z ∈ Φ d 1 {− λ 1 /d x ∈ D ◦ (0 , λ 1 /d z ); b (0 , λ 1 /d z ) ∩ λ 1 /d Φ = ∅} λ − 2 | W (0 | λ 1 /d Φ ∪ { 0 } ) | 2   dx whic h, since λ 1 /d Φ is a u nit rate Po isson pr o cess reduces to = λ 2 Z E 1 " P { z 1 ,...,z d }⊂ Φ 1 {− x ∈ D ◦ (0 , z 1 , . . . , z d ); b (0 , z 1 , . . . , z d ) ∩ Φ = ∅} | W (0 | Φ ∪ { 0 } ) | 2 # dx. The sum of d -v olumes of Delauna y cells in v olving 0 is th at of its con tiguous V oronoi cell, and w e conclude th at C ′ ( λ, d ) = λ 2 E 1  1 | W (0 | Φ ∪ { 0 } ) |  .  The ab ov e results can b e s ummarised as follo ws. Theorem 6. L et Φ b e a stationary P oisson p oint pr o c ess in R d with intensity λ > 0 . Then, the Delaunay tessel lation field estimator d λ (0) has asymptotic varianc e c d λ 2 with c d = E 1   1 | W (0 | Φ ∪ { 0 } ) |    1 + X y ∈N (0 | Φ ∪{ 0 } ) | W (0 | Φ ∪ { 0 } ) ∩ W ( y | Φ ∪ { 0 } ) | | W ( y | Φ ∪ { 0 } ) |      − 1 . 11 Note that the classic Berman–Diggle estimator (5) is asymptotically unbiased with v ari- ance λ ω − 1 d h − d , where ω d is the volume of the unit ball in R d . In words, the Berman–Diggle estimator is more efficien t wh en ev er the av erage num b er of p oin ts p er test ball exceeds 1 /c d . 6 P oisson pro cesses on the line F or one-dimensional Po isson pro cesses, the distribu tion of th e con tiguous V oronoi cell can b e calculate d exp licitly for arbitrary in tens ity functions. F or simplicit y , assume that A = [ − w, w ] is an interv al of radiu s w > 0 either side of the origin. The follo wing lemma is w ell-kno wn . Lemma 5. L et Φ b e a P oisson p oint pr o c ess on [ − w, w ] with finite intensity function λ and write Λ( a, b ) = R b a λ ( x ) dx for the moment me asur e of ( a, b ) for any − w ≤ a ≤ b ≤ w . F or x ∈ ( − w , w ) , define the r andom variables Φ − ( x ) := max { y ∈ {− w } ∪ (Φ ∩ [ − w, x )) } ; Φ + ( x ) := min { y ∈ { w } ∪ (Φ ∩ ( x, w ]) } . Then, their distribution functions ar e given by F − ( t ) = exp [ − Λ( t, x )] for t ∈ ( − w , x ) , with an atom of mass P (Φ − ( x ) = − w ) = exp [ − Λ( − w , x )] at − w , r esp e c tively F + ( s ) = 1 − exp [ − Λ( x, s )] for s ∈ ( x, w ) with an atom at w of mass P (Φ + ( x ) = w ) = exp [ − Λ( x, w )] . M or e over, for fixe d x , Φ + ( x ) and Φ − ( x ) ar e indep endent r andom variables. 6.1 Exp ectation of the DTFE for Poisson pro c esses on the line Note that on th e real lin e, the con tiguous V oronoi cell W ( x | (Φ ∪ { x } ) ∩ [ − w , w ]) is the in terv al [Φ − ( x ) , Φ + ( x )]. Thus, Lemma 5 can b e used to calculate the moment s of the Delauna y tessellatio n field estimator. I n this section, we sh all deal with edge effects by p lacing t w o ghost p oints at th e b ord ers − w and w . Theorem 7. L et Φ b e a Poisson p oint pr o c ess observe d in A = [ − w , w ] for some w > 0 with lo c al ly finite intensity function λ : R → [0 , ∞ ) . Then, for x 0 ∈ A , E h \ λ ( x 0 ) i = Z x 0 − w Z w x 0 Λ( t, s ) λ ( s ) λ ( t ) s − t e − Λ( t,s ) dt ds + Λ( − w, w ) e − Λ( − w ,w ) 2 w + Z w x 0 Λ( − w, s ) λ ( s ) w + s e − Λ( − w ,s ) ds + Z x 0 − w Λ( t, w ) λ ( t ) w − t e − Λ( t,w ) dt. (8) 12 Pro of: Fix x 0 6 = x ∈ ( − w , w ), and let ϕ b e a realisation of Φ, which w e augment by − w and w in order to obtain b ound ed Delauna y cells. Since almost sur ely , x 6∈ Φ and x 0 6∈ Φ, assume x 0 , x 6∈ ϕ , and consider g ( x 0 | x, ϕ ∪ { x } ) as defined in (3). Note that x 0 b elongs to a single Delaunay cell inte rior. If x is no endp oint of this cell, g ( x 0 | x, ϕ ∪ { x } ) = 0. Otherw ise, g ( x 0 | x , ϕ ∪ { x } ) = 1 / ( ϕ + ( x 0 ) − ϕ − ( x 0 )), cf. Lemm a 5 . First, assume x < x 0 . By Lemma 5 ap p lied to th e p oin t x , E [ g ( x 0 | x , Φ ∪ { x } )] = Z x − w Z w x 0 dF − ( t ) dF + ( s ) s − t = e − Λ( − w ,w ) 2 w + + Z w x 0 λ ( s ) w + s e − Λ( − w ,s ) ds + Z x − w λ ( t ) w − t e − Λ( t,w ) dt + Z x − w Z w x 0 λ ( s ) λ ( t ) s − t e − Λ( t,s ) dt ds. Similarly , for x 0 < x , E [ g ( x 0 | x, { x } ∪ Φ)] = Z x 0 − w Z w x dF − ( t ) dF + ( s ) s − t = e − Λ( − w ,w ) 2 w + Z w x λ ( s ) w + s e − Λ( − w ,s ) ds + Z x 0 − w λ ( t ) w − t e − Λ( t,w ) dt + Z x 0 − w Z w x λ ( s ) λ ( t ) s − t e − Λ( t,s ) dt ds. By Theorem 1, the exp ectation of th e Delauna y tessellation fi eld estimator is as s tated for x 0 ∈ ( − w , w ). It remains to consider x 0 = − w or w . In the fi rst case, ϕ − ( x 0 ) is replaced by − w ; for x 0 = w , ϕ + ( x 0 ) b ecomes w in the ev aluation of g ( x 0 | x , ϕ ∪ { x } ). Th us, for example, E [ g ( − w | x, Φ ∪ { x } )] = e − Λ( − w ,x ) Z w x dF + ( s ) w + s = e − Λ( − w ,w ) 2 w + Z w x λ ( s ) w + s e − Λ( − w ,s ) ds, with a similar expression for x 0 = w . Up on in tegratio n, (8) is obtained, under the conv ent ion that inte grals ov er int erv als of zero length v anish .  In general, (8) m ust b e ev aluated n umerically . F or the homogeneous Poisson pro cess, analytic ev aluation is p ossible. I n fact, it can b e sho wn that the estimator is unbiased even near the b orders of the observ ation inte rv al. Corollary 3. L et Φ b e a stationary Poisson p oint pr o c ess observe d in A = [ − w , w ] for some w > 0 with intensity λ > 0 . Then, the Delaunay tessel lation field estimator \ λ ( x 0 ) is unbiase d for al l x 0 ∈ A . Pro of: F or a stationary Poi sson p oin t pro cess, the d ouble in tegral in (8) reduces to λ  e λx 0 − e − λw  ×  e − λx 0 − e − λw  and in particular v anish es for x 0 = − w or w . T he th r ee b order correction terms are equ al to λe − 2 λw , to λe − λw ( e − λx 0 − e − λw ), and to λe − λw ( e λx 0 − e − λw ), resp ectiv ely . Th e sum of all four terms is λ , so th e estimator is unbiased.  Note that the Berman–Diggle estimator is unbiased as w ell, but that this ma y not b e true for (6) d ue to edge correction near the b order. 13 6.2 V ariance of the DTFE for Poisson pro cesses on the line In this section, we derive the asymptotic v ariance of the Delauna y tessellation field estimator for a stationary P oisson pro cess on the line. The r esult can b e used to approximat e the v ariance when the underlyin g intensit y function is smo othly v arying. Theorem 8. L et Φ b e a stationar y Poisson p oint pr o c e ss observe d in A = [ − w , w ] for some w > 0 with intensity λ > 0 . Then, as w → ∞ , the Delaunay tessel lation field estimator d λ (0) has asymptotic varianc e 2 λ 2 (2 − π 2 / 6) ≈ 0 . 7 λ 2 . The result should b e compared to λ/ (2 h ) for the Berman–Diggle k ernel estimator [1 ], see also [10]. If 2 λh > 1 . 4, that is the a v erage n um b er of p oin ts p er b in at least 1 . 4, kernel estimation is th e b etter choic e. Naturally , in order to compute \ λ ( x 0 ), t w o p oint s of the underlying pro cess are used. In order to giv e the pro of, some sp ecial fun ction th eory is needed. Let x > 0. Recall that the exp on ential integral is defined as E 1 ( x ) = Z ∞ 1 e − tx t dt = Z ∞ x e − u u du. Its in tegral satisfies E 2 ( x ) = Z ∞ x E 1 ( s ) ds = e − x − xE 1 ( x ) . In the limit, E 1 (0) = ∞ and E 2 (0) = 1. F u rthermore, Z ∞ 0 u e u E 1 ( u ) 2 du = 2 − π 2 6 . See for example [7] for fu rther details. W e shall also need the equation Z c 0 e ax E 1 ( ax ) dx = γ + log ( ac ) + e ac E 1 ( ac ) a where a and c are strictly p ositiv e constant s, and γ ≈ 0 . 577 is the Euler-Masc h eroni constan t. Pro of: By Theorem 3, asymptotically E h d λ (0) i = λ . F or th e v ariance, b y Theorem 2, w e need to ev aluate t w o fu r ther integrals. No w, argueing as in the pro of of Theorem 7 , Z A E  g 2 ( x 0 | x, Φ ∪ { x } )  λ ( x ) dx = Z x 0 − w Z w x 0 Λ( t, s ) λ ( s ) λ ( t ) ( s − t ) 2 e − Λ( t,s ) dt ds + Λ( − w, w ) e − Λ( − w ,w ) 4 w 2 + Z w x 0 Λ( − w, s ) λ ( s ) ( w + s ) 2 e − Λ( − w ,s ) ds + Z x 0 − w Λ( t, w ) λ ( t ) ( w − t ) 2 e − Λ( t,w ) dt. (9) Since the intensit y fun ction is constan t and w e to ok x 0 = 0, (9) reduces to λe − 2 λw 2 w + λe − λw Z w 0 λe − λs w + s ds + λe − λw Z 0 − w λe λt w − t dt + Z 0 − w Z w 0 λ 3 e λt e − λs s − t dt ds. 14 Clearly , the fi r st term ab o ve con v erges to 0 as w → ∞ . Du e to symmetry , the t w o midd le terms are equal. Note that 2 λ Z w 0 λe − λ ( s + w ) s + w ds = 2 λ 2 Z 2 λw λw e − u u du = 2 λ 2 [ E 1 ( λw ) − E 1 (2 λw )] , whic h con v erges to zero as w → ∞ . Moreo ve r, λ 3 Z 0 −∞ Z ∞ 0 e λt e − λs s − t dt ds = λ 3 Z 0 −∞ E 1 ( − λt ) dt = λ 2 E 2 (0) = λ 2 . T o calculate the double in tegral in Theorem 2, let x 6 = y b e p oin ts of ( − w, w ), fix x 0 6∈ { x, y , − w , w } , and let ϕ b e a realisation of Φ, w hic h w e augmen t by − w and w in order to obtain b oun ded Delauna y cells. Since almost su rely none of x , y or x 0 lie in Φ, assume x 0 , x, y 6∈ ϕ , and consider g ( x 0 | x, ϕ ∪ { x, y } ) as defin ed in (3). Note that x 0 b elongs to a single Delauna y cell int erior. If x and y are n ot b oth endp oints of this cell, g ( x 0 | x, ϕ ∪ { x, y } ) g ( x 0 | y , ϕ ∪ { x, y } ) = 0. Otherw ise, without loss of generalit y , x < x 0 < y , and g ( x 0 | x , ϕ ∪ { x, y } ) = 1 / ( y − ϕ − ( x 0 )) and g ( x 0 | y , ϕ ∪ { x, y } ) = 1 / ( ϕ + ( x 0 ) − x ). Th us, for x < x 0 and y > x 0 , let F − and F + b e the cum ulativ e distribu tion functions of Φ − ( x 0 ) and Φ + ( x 0 ). By Lemma 5 , E [ g ( x 0 | x, Φ ∪ { x, y } ) g ( x 0 | y , Φ ∪ { x, y } )] = Z x − w Z w y dF − ( t ) dF + ( s ) ( y − t ) ( s − x ) = Z w y λ ( s ) ( w + y ) ( s − x ) e − Λ( − w ,s ) ds + Z x − w λ ( t ) ( y − t ) ( w − x ) e − Λ( t,w ) dt + e − Λ( − w ,w ) ( w + y ) ( w − x ) + Z x − w Z w y λ ( s ) λ ( t ) ( y − t ) ( s − x ) e − Λ( t,s ) dt ds. By symmetry , Z A Z A E [ g ( x 0 | x , Φ ∪ { x, y } ) g ( x 0 | y , Φ ∪ { x, y } )] λ ( x ) λ ( y ) dx dy = 2 e − Λ( − w ,w ) Z x 0 − w λ ( x ) w − x dx Z w x 0 λ ( y ) w + y dy + 2 Z w x 0 λ ( s ) e − Λ( − w ,s )  Z x 0 − w λ ( x ) s − x dx Z s x 0 λ ( y ) w + y dy  ds + 2 Z x 0 − w λ ( t ) e − Λ( t,w )  Z x 0 t λ ( x ) w − x dx Z w x 0 λ ( y ) y − t dy  dt + 2 Z x 0 − w Z w x 0 λ ( t ) λ ( s ) e − Λ( t,s )  Z x 0 t λ ( x ) s − x dx Z s x 0 λ ( y ) y − t dy  dt ds. (10) F or x 0 ∈ {− w , w } , form ula (10) holds true un d er the con v entio n that integral s ov er inte rv als of zero length v anish , as in this case x 0 cannot b elong to any Delauna y cell with endp oints x < x 0 < y . 15 Next, w e plug in x 0 = 0 and λ ( · ) ≡ λ , and consid er eac h integ ral in (10) in tur n. Th e main term is the four fold in tegral Z 0 − w Z w 0 Z x − w Z w y 2 λ 4 e λt e − λs ( y − t ) ( s − x ) dx dy dt ds. Its limit as w → ∞ is 2 λ 4 Z 0 −∞ Z ∞ 0 e λ ( y − x ) " Z x −∞ e − λ ( y − t ) y − t dt Z ∞ y e − λ ( s − x ) s − x ds # dx dy = 2 λ 4 Z 0 −∞ Z ∞ 0 e λ ( y − x ) E 1 ( λ ( y − x )) 2 dx dy = 2 λ 3 Z 0 −∞ Z ∞ − λx e u E 1 ( u ) 2 dx du = 2 λ 2 Z ∞ 0 Z ∞ y e u E 1 ( u ) 2 dy du = 2 λ 2 Z ∞ 0 u e u E 1 ( u ) 2 du = 2 λ 2 (2 − π 2 / 6) , up on a change of int egration order. The first term in (10) reduces to 2 e − 2 λw ( λ log 2) 2 for a homogeneous Po isson p ro cess, whic h tends to zero as w → ∞ . It remains to consider the sum of the tw o three f old in tegrals in (10) Z 0 − w Z w 0 Z w y 4 λ 3 e − λ ( s + w ) ( s − x ) ( y + w ) dx dy ds whic h can b e written as 4 λ 3 Z w 0 Z w 0  Z s 0 dy y + w  e − λ ( s + w ) s + x dx ds ≤ 4 λ 3 log 2 Z w 0 e − λw + λx Z w 0 e − λ ( s + x ) s + x ds ! dx = 4 λ 3 log 2 Z w 0 e − λw + λx [ E 1 ( λx ) − E 1 ( λx + λw )] dx = 4 λ 2 h ( λ, w ) log 2 , where h ( λ, w ) = e − λw Z λw 0 e u [ E 1 ( u ) − E 1 ( u + λw )] du =  e − λw + e − 2 λw  Z λw 0 e u E 1 ( u ) du − e − 2 λw Z 2 λw 0 e u E 1 ( u ) du = e − λw γ + ( e − λw + e − 2 λw ) log ( λw ) − e − 2 λw log(2 λw ) + E 1 ( λw )(1 + e − λw ) − E 1 (2 λw ) tends to zero as w → ∞ . The pro of is fin ished u p on collectio n of all terms.  As a corollary , the pro of give s an expression for the second moment of the Delauna y tessellatio n field estimator of the inte nsit y function for P oisson pro cesses with not necessarily constan t lo cally fin ite intensit y fun ction on interv als of the f orm [ − w, w ] by com bin ing (9)– (10). A slightl y simpler pro of can b e obtained b y an app eal to Theorem 6 , but such a pro of cannot b e generalised to non-homogeneous P oisson pro cesses. 16 7 Discussion In th is pap er, w e analysed Sc haap and V an de W eyga ert’s Delauna y tessel lation field estimator [14, 15] for the in tensit y fu nction of a p oin t p r o cess. W e expressed its mean an d v ariance in term s of the first and second order factorial m omen t measur es of the und erlyin g p oin t pro cess, and placed the estimator in the con text of ad ap tive k ernel estimation. W e then fo cussed on Po isson p oin t p ro cesses, and s h o wed that f or stationary P oisson pro cesses, the DTFE is asymptotically un biased with a v ariance that is prop ortional to the squared in tensity . The p rop ortionalit y constan t dep ends on the dimension. F or d = 1, explicit calculat ion is p ossible. F or d = 2, w e used th e DELDIR pac k age [19] to obtain C ( λ, 2) ≈ 0 . 8 λ 2 and C ′ ( λ, 2) ≈ 0 . 6 λ 2 , see Lemma 3 and 4. Note that in the plane it is p ossible to write mean and v ariance as rep eated in tegral s in th e spirit of Calk a [2], b ut explicit ev aluation seems difficult. Sim u lations for the case d = 3 of most inte rest to cosmologists can b e f ound in Sc haap’s Ph.D. thesis [15]. Ac kno wledgemen t The author is grateful to Dr. N.M. T emme for access to [7]. References [1] M. Berman and P .J. Diggle (1989). Estimating we igh ted in teg rals of the second-order in tensit y of a spatial p oin t p ro cess. Journal of the R oyal Statistic al So ci e ty Series B , 51:81– 92. [2] P . Calk a (2003). An explicit expression for the distribu tion of the n um b er of s id es of th e t ypical P oisson–V oronoi cell. A dvanc es in Applie d Pr ob ability (SGSA) , 35:863–87 0. [3] N.A.C. Cressie (1991, revised edition 1993). Statistics for sp atial data. John Wiley and Sons, New Y ork. [4] D.J. Daley and D. V ere–Jones (1988, revised edition 2003, 2008). An intr o duction to the the ory of p oint pr o c esses. S pringer, New Y ork. [5] B.N. Delaunay (1934). Sur la sph` ere vide. Iszvestia Akademii Nauk SSSR, Otdelenie Matematicheskikh i Estestvennykh N auk , 7:793–8 00. [6] P .J. Diggle (1983 , revised edition 2003). Statistic al analysis of sp atial p oint p atterns. Academic Press, London. [7] M. Geller and E.W. Ng (1969 ). A table of inte grals of th e exp onentia l inte gral. Journal of R ese ar ch of the National Bu r e au of Standar ds , 73B:191–2 10. [8] K.-H. Hanisc h (1982). On inv ersion f orm ulae for n -fold Palm distrib utions of p oin t p ro- cesses in LC S -spaces. M athematisch e N achrichten , 106:171– 179. 17 [9] M.N.M. v an Lieshout (2000 ). Markov p oint pr o c esses and their applic ations. Im p erial College Press, Lond on. [10] M.N.M. v an L ieshout (2007). Edge corrected non-parametric in tensit y function estima- tors for h eterogeneo us P oisson p oint pr o cesses. EURANDOM Rep ort 2007-04 2. [11] J. Møller (1994). L e ctur es on r andom V or onoi tessel lations. Lecture Notes in Statistics 87. Sprin ger-V erlag, New Y ork. [12] A. Ok ab e, B. Bo ots, and K . S ugihara (1992). Sp atial tessel lations. Conc e pts and applic a- tions of V or onoi diagr ams. With a for ewor d by D.G. Kendal l. Wiley S eries in Probability and Mathematical S tatistics. John Wiley and Sons, C hic h ester. [13] J.K. Ord (1978). Ho w many trees in a forest? Mathematic al Scie ntist , 3:23–33 . [14] W.E. Sc haap and R. v an de W eygaert (2000). Letter to the Editor. C on tin uous fi elds and discrete samples: reconstruction th rough Delauna y tessellations. Astr onomy and Astr ophysics , 363:L2 9–L32. [15] W.E. Sc haap (2007) . DTFE. The Delaunay T essel lation Field Estimator. Ph.D. Thesis, Univ ersit y of Gronin gen. [16] R. S ibson (1981 ). A brief d escription of n atural neighbour interp olatio n. In : V. Barnett (Ed.) Interpr eting multivariate data. John Wiley and S on s , New Y ork, 21–36. [17] B.W. Silv erman (1986). Density estimation for statistics and data analysis. Chapman and Hall, London . [18] D. Sto y an, W.S. Kend all, and J. Mec ke (1987, revised edition 1995). Sto chastic g e ometry and its applic ations. John Wiley and Sons, Chichester. [19] R. T urner. DELDIR, An R pac k age to construct the Delaunay triangulation and the Diric h let (V oronoi) tessellation of a planar p oint set. htt p://www.m ath.unb.ca / rolf . [20] G. V oronoi (1908 ). Nouv elles app licatio ns des p aram ` etres con tinues ` a la th´ eorie des formes quadratiques. R ei ne und Angewandte mathematik 134:198–28 7. 18