MH4522 Spatial Data Science Assignment Questions | NTU

Assignment Brief

The classical kernel estimator1 of the probability density function ϕ(x) of a random variable X is defined by

bϕh(x) := 1⁄nh ∑i=1n φ(x − xi⁄h),
where xi, i = 1, …, n, are n independent samples of X. Here, h > 0 is a positive parameter called the bandwidth, and φ is a bounded probability density function, such that

limx→+∞ x|φ(x)| = 0.
N=1; h=0.1; z =seq(0,1,0.01); kernel=function(z){dnorm(z,0,h/4)};
x=runif(N); kdensity=function(z){sum(as.numeric(lapply(z-x,kernel))/length(x)}
plot(0, xlab = "", ylab = "", type = "l", xlim = c(0,1), col = 0,
ylim=c(0,max(as.numeric(lapply(z,kdensity))),xaxt='n',yaxt='n')
axis(1, at=c(), xlab = "", lwd=2,labels=c(), pos=0,lwd.ticks=2)
axis(2, lwd=2, at = c(1,axTicks(4)), lwd.ticks=2); points(x, rep(0,N), pch=3, lwd = 3, col = "blue")
lines(density(x,width=h),col="purple",lwd=3); lines(z,dunif(z),col="black",lwd=3);
lines(z,as.numeric(lapply(z,kdensity)),col="red",lwd=2,type='l')

This assignment aims to implement a kernel estimation for the intensity of a Poisson point process η on ℝd, d ≥ 1. We assume that the intensity measure µ of η has a C2b density ρ : ℝd → ℝ+ for the Lebesgue measure on (ℝd, B(ℝd)), i.e. µ(dx) = ρ(x)dx, and

IE[η(B)] = µ(B) = ∫B ρ(x)dx, B ∈ B(ℝd).
We also let

∥x∥ = √x12 + ··· + xd2, (x1, …, xd) ∈ ℝd,
denote the Euclidean norm in ℝd, and we denote by φh the Gaussian kernel

φh(u) := 1⁄(2πh2)d/2 e−u2/(2h2), u ∈ ℝ,
with variance h > 0. The following questions are interdependent and should be treated in sequence.

1) Show that for all x ∈ ℝd we have

limh→0 ∫ℝd φh(∥x − y∥)ρ(y)dy1 ··· dyd = ρ(x).
Hint: You may use Taylor’s formula with integral remainder term
ρ(y) = ρ(x) + ∑k=1d (yk − xk)∂ρ⁄∂xk(x) + ∑k,l=1d (yk − xk)(yl − xl) ∫01 (1 − t)∂2ρ⁄∂xk∂xl(x + t(y − x))dt,
x, y ∈ ℝd.

2) Show that the estimator

ρ̂h,0(x) := ∑y∈η φh(∥x − y∥)
of the density ρ(x) is asymptotically unbiased, i.e. we have

limh→0 IE[ρ̂h,0(x)] = ρ(x), x ∈ ℝd.
Hint: Apply Proposition 4.6-a) and the result of Question (1).

3) Show that the asymptotic variance of the estimator ρ̂h,0 satisfies

Var[ρ̂h,0(x)] ≃h→0 ρ(x)⁄(2h)dπd/2, x ∈ ℝd,
i.e.

limh→0 hd Var[ρ̂h,0(x)] = ρ(x)⁄2dπd/2, x ∈ ℝd.
Hint: Apply Proposition 4.6-b) and the result of Question (1).

4) Given a domain A ⊂ ℝd such that 0 < µ(A) < ∞ and f ∈ L1(A, µ), compute the expectation

IE [1{η(A)≥1} 1⁄η(A) ∫A f(x)η(dx)].
Hint: Apply Proposition 4.4 and proceed similarly to the proof of Proposition 4.6-a).

5) Given a domain A ⊂ ℝd such that 0 < µ(A) < ∞ and f ∈ L1(A, µ) ∩ L2(A, µ), compute the variance

Var [1{η(A)≥1} 1⁄η(A) ∫A f(y)η(dy)],
using the quantity

c(A) := IE [ 1⁄η(A) 1{η(A)≥1} ].
Hint: Apply Proposition 4.4 and proceed similarly to the proof of Proposition 4.6-b).

6) Show that for any domain A ⊂ ℝd such that 0 < µ(A) < ∞, we have

c(A) ≤ 2⁄µ(A).
Hint: Write c(A) as a series, and upper bound it term by term.

7) For any h > 0, let Ah ⊂ ℝd denote a domain of finite Lebesgue measure in ℝd, and consider the estimator ρ̂h,1 of the probability density ρ(x)/µ(Ah) defined by

ρ̂h,1(x) := 1{η(Ah)≥1} 1⁄η(Ah) ∑y∈η φh(∥x − y∥).
Show that ρ̂h,1(x) is asymptotically unbiased in the sense that

IE[ρ̂h,1(x)] − ρ(x)⁄µ(Ah) = o(µ(Ah)−1)
as h → 0, i.e.

limh→0 |µ(Ah)IE[ρ̂h,1(x) − ρ(x)⁄µ(Ah)]| = 0, x ∈ ℝd,
provided that µ(Ah) → ∞ as h → 0.

Hint: Apply the results of Question (1) and Question (4).

8) Show that the variance of ρ̂h,1(x) satisfies

limh→0 Var[ρ̂h,1(x)] = 0
provided that µ(Ah)−1 = o(hd).

Hint: Apply the result of Question (5) and use the result of Question (1) as in Question (3), together with the result of Question (6).

9) Show that for any domain A ⊂ ℝd such that 0 < ℓ(A) < ∞ we have

IE[∑y∈η∩A 1⁄ρ(y)] = ℓ(A).

10) Based on a dataset of your choice on a domain A, find the value of h > 0 that minimizes the quantity

IE[(∑y∈η∩A 1⁄ρ̂h,0(y) − ℓ(A))2]
and compare the estimations of the density ρ(x) obtained from ρ̂h,0 and ρ̂h,1 (graphs are welcome).

Examples of datasets include:

Simulated datasets;
The spatstat package; see https://cran.r-project.org/web/packages/spatstat/vignettes/datasets.pdf;
For the scikit-learn package in Python, see https://scikit-learn.org/stable/datasets/real_world.html and this example.
See also:

P. Moraga. Geospatial Health Data – Modeling and Visualization with R-INLA and Shiny. Chapman & Hall/CRC Biostatistics Series. CRC Press, 2020.
P. Moraga. Spatial Statistics for Data Science – Theory and Practice with R. Chapman & Hall/CRC Data Science Series. CRC Press, 2024.