From Spatial Counts to a Likelihood
A spatial point pattern records locations only within a specified observation window. For likelihood inference, both the observed points and the area in which they could have been observed matter: the points contribute intensity at their locations, while the unoccupied parts of the window contribute evidence through the total intensity there. Treating an unobserved region as if it had been surveyed changes the likelihood and can change the estimate.
Let \(W\subset\mathbb{R}^d\) be a bounded Lebesgue-measurable observation window with \(0<|W|<\infty\), where \(|W|\) denotes its \(d\)-dimensional volume. Work on a probability space \((\Omega,\mathcal{F},\mathbb{P})\). Let \(N(A)\) count observed points in each measurable \(A\subseteq W\). Suppose these counts come from a Poisson point process with intensity density \(\lambda_\theta(x)\ge0\), relative to volume, for parameter \(\theta\). Assume \(\Lambda_\theta(W)=\int_W\lambda_\theta(x)\,dx<\infty\). The observed pattern consists of \(n=N(W)\) locations \(x_1,\ldots,x_n\in W\). Since the intensity is a density, locations on any volume-zero set, including a regular window boundary of zero volume, occur with probability zero.
Here the likelihood uses the observed locations and the observed count. To see its origin, put \(\Lambda=\int_W\lambda_\theta(x)\,dx\). When \(\Lambda>0\), the count \(N(W)\) has a Poisson distribution with mean \(\Lambda\); conditional on \(N(W)=n\), the locations have density \(\lambda_\theta(x)/\Lambda\) on \(W\). The conditional-location fact follows from the independent Poisson counts: for a finite partition \(A_1,\ldots,A_m\) of \(W\), conditioning independent Poisson counts on their sum \(n\) gives multinomial cell probabilities \(\Lambda_\theta(A_j)/\Lambda\), exactly those of \(n\) independent locations with that density. Combining the count probability and the conditional density for an ordered list gives $$ \frac{e^{-\Lambda}\Lambda^n}{n!}\prod_{r=1}^n\frac{\lambda_\theta(x_r)}{\Lambda} =\frac{e^{-\Lambda}}{n!}\prod_{r=1}^n\lambda_\theta(x_r). $$ The factor \(1/n!\) depends on the count but not on \(\theta\); dropping it gives the stated likelihood. If \(\Lambda=0\), there are no points almost surely, and the empty-pattern likelihood is \(1\), consistent with the formula.
Constant Intensity
For a constant intensity \(\lambda_\theta(x)=\theta\) with \(\theta\ge0\), the likelihood depends on the window through its volume, not just through the observed coordinates. Write \(V=|W|\). For \(n>0\) and \(\theta>0\), the log-likelihood is \(n\log\theta-V\theta\). Its derivative is \(n/\theta-V\) and its second derivative is \(-n/\theta^2<0\). The likelihood is zero at \(\theta=0\) when \(n>0\), so the unique maximum is \(\widehat{\theta}=n/V\). If \(n=0\), the likelihood is \(e^{-V\theta}\), maximized at \(\theta=0\) when that value is allowed. If the parameter space is instead \(\theta>0\), the empty-pattern likelihood has supremum \(1\) as \(\theta\) decreases to zero, but no maximum within that parameter space.
Worked Example: Estimating a Constant Intensity
Let \(W=[0,2]\times[0,3]\), so its area is \(2\cdot3=6\). Suppose the observed locations are \((\tfrac12,\tfrac12)\), \((1,1)\), and \((\tfrac32,2)\). Under constant intensity \(\theta\ge0\), there are \(n=3\) points. The likelihood, omitting only a parameter-independent factor, is $$ L(\theta)=e^{-6\theta}\theta^3,\qquad \theta>0. $$ The log-likelihood is \(\ell(\theta)=3\log\theta-6\theta\). Its derivative is \(3/\theta-6\), which vanishes at \(\widehat{\theta}=1/2\). Its second derivative is \(-3/\theta^2<0\), and \(L(0)=0\); thus this is the maximum over \(\theta\ge0\).
Check: Substituting \(\theta=1/2\) into the score gives \(3/(1/2)-6=6-6=0\). The log-likelihood there is \(3\log(1/2)-6(1/2)=-3\log2-3\), and the likelihood is \(e^{-3}(1/2)^3=e^{-3}/8\). The estimate also equals count divided by observed area: \(3/6=1/2\).
Covariates and Log-Linear Intensities
A common spatial model represents intensity through covariates \(z(x)\in\mathbb{R}^p\): $$ \lambda_\beta(x)=\exp\!\big(\beta^\mathsf{T}z(x)\big),\qquad \beta\in\mathbb{R}^p. $$ The exponential ensures positive intensity. Its log-likelihood is $$ \ell(\beta)=\sum_{r=1}^n\beta^\mathsf{T}z(x_r)-\int_W\exp\!\big(\beta^\mathsf{T}z(x)\big)\,dx. $$ The integral is the model's expected number of points in the observed window. In general it depends on \(\beta\), and maximizing the likelihood requires accounting for it. A finite likelihood integral at one parameter value alone does not guarantee that derivatives can be taken under the integral.
Proof. For concavity, take any two parameter vectors \(\beta_0,\beta_1\) in the domain and \(a\in[0,1]\). At each \(x\), the exponential function is convex, so $$ e^{(a\beta_0+(1-a)\beta_1)^\mathsf{T}z(x)} \le a e^{\beta_0^\mathsf{T}z(x)}+(1-a)e^{\beta_1^\mathsf{T}z(x)}. $$ Integrating preserves the inequality. Therefore the negative integral in \(\ell\) is concave as a function of \(\beta\). The data term \(\sum_r\beta^\mathsf{T}z(x_r)\) is linear, so adding it preserves concavity. Where differentiation under the integral is justified, differentiating the data term gives \(\sum_r z(x_r)\), and differentiating the integral once gives \(\int_W z(x)e^{\beta^\mathsf{T}z(x)}dx\). Differentiating again gives the matrix integral of \(z(x)z(x)^\mathsf{T}e^{\beta^\mathsf{T}z(x)}\), with a minus sign in the log-likelihood. For any vector \(v\in\mathbb{R}^p\), the quadratic form of this Hessian is $$ v^\mathsf{T}\nabla^2\ell(\beta)v=-\int_W (v^\mathsf{T}z(x))^2e^{\beta^\mathsf{T}z(x)}\,dx\le0. $$ This also verifies concavity wherever the Hessian formula applies. At an interior differentiable maximum, each directional derivative must vanish, giving the score equation. \(\square\)
One convenient sufficient regularity condition is that \(z\) is bounded almost everywhere on the finite-volume window. On any compact neighborhood of a fixed \(\beta\), the exponential and its first two derivatives are then bounded by integrable constants, so differentiation under the integral is valid. More generally, suitable integrable bounds for the derivatives on a neighborhood suffice. Finiteness of the likelihood integral and score integral at a single parameter is not enough to assert the Hessian formula.
Worked Example: A Covariate Model with an Interior Maximum
On \(W=[0,1]\times[0,1]\), consider \(\lambda_\beta(x,y)=e^{\beta x}\), where \(\beta\in\mathbb{R}\). The two observed locations are \((\tfrac14,\tfrac14)\) and \((\tfrac14,\tfrac34)\), so their \(x\)-coordinates sum to \(1/4+1/4=1/2\). Since \(z(x,y)=x\) is bounded, the score and Hessian formulas apply. The log-likelihood is $$ \ell(\beta)=\frac{\beta}{2}-\int_0^1\int_0^1e^{\beta x}\,dy\,dx. $$ At \(\beta=0\), the integral is \(\int_0^1\int_0^1 1\,dy\,dx=1\), so \(\ell(0)=-1\). The score there is \(1/2-\int_0^1\int_0^1x\,dy\,dx=1/2-1/2=0\). The Hessian at zero is \(-\int_0^1\int_0^1x^2\,dy\,dx=-1/3\).
For every finite \(\beta\), the Hessian is \(-\int_0^1x^2e^{\beta x}dx<0\), because the integrand is positive for \(x>0\). Thus the log-likelihood is strictly concave, and its zero score at \(\beta=0\) makes \(0\) the unique global maximizer. Check: The data contribution to the score is \(1/4+1/4=1/2\), while the fitted integral contribution at zero is \(\int_0^1 x\,dx=1/2\); their difference is exactly zero.
What the Observation Boundary Changes
The integral in the likelihood is over \(W\), the region where points were observed. It is not automatically an integral over a larger scientific region, nor over a neighborhood of the window. For a homogeneous process, the relevant exposure is the area or volume of \(W\). For a spatially varying intensity, the relevant exposure is the integral of that intensity over the portion actually observed. If a covariate is defined on a larger region, values outside \(W\) do not add observed exposure under this likelihood.
Geometric boundaries matter particularly when a window cuts through a region used to describe the model. A cell or neighborhood that is only partly inside \(W\) contributes only its in-window portion to \(\int_W\lambda_\theta(x)dx\). There is no general rule to “correct” a point-location likelihood by expanding the observation window: doing so treats unobserved area as though its events had been recorded. A different sampling design, such as conditioning on the total count or observing only a selected subset of points, may have a different likelihood and must be modeled as such.
Worked Example: A Window Clipped from a Larger Region
Suppose the underlying model on \(D=[0,2]\times[0,1]\) has intensity \(\lambda_\theta(x,y)=\theta(1+x)\), with \(\theta>0\), but observations are available only in \(W=[0,1]\times[0,1]\). There are two observed points, at \(x\)-coordinates \(1/4\) and \(3/4\). The point-intensity product is $$ \theta^2(1+\tfrac14)(1+\tfrac34)=\theta^2\cdot\tfrac54\cdot\tfrac74=\tfrac{35}{16}\theta^2. $$ The correct integrated intensity is over \(W\): $$ \int_0^1\int_0^1\theta(1+x)\,dy\,dx =\theta\int_0^1(1+x)\,dx =\theta(1+\tfrac12)=\tfrac32\theta. $$ Consequently \(L_W(\theta)=\tfrac{35}{16}\theta^2e^{-3\theta/2}\), up to a parameter-independent factor, and \(\ell_W'(\theta)=2/\theta-3/2\). The score vanishes at \(\widehat{\theta}=4/3\); since \(\ell_W''(\theta)=-2/\theta^2<0\), this is the unique maximum.
Check: At \(4/3\), the score is \(2/(4/3)-3/2=3/2-3/2=0\). If one incorrectly integrates over all of \(D\), the integral is \(\theta\int_0^2(1+x)dx=\theta(2+2)=4\theta\). That incorrect likelihood gives score \(2/\theta-4\) and estimate \(1/2\), not \(4/3\). The discrepancy is caused by including the unobserved strip \(1<x\le2\) as exposure.
Common Mistakes and Practical Checks
A frequent error is to retain the product of intensities at observed points but omit the negative integrated-intensity term. For the constant-intensity example with three points in area \(6\), this wrong calculation would maximize \(\theta^3\) by increasing \(\theta\) without bound. The correction is to include the empty-space contribution \(e^{-6\theta}\): the log-likelihood \(3\log\theta-6\theta\) has score \(3/\theta-6\), and its maximum is \(1/2\). Another error is to use the area of a nominal study region rather than the region actually observed. Before maximizing, verify the window in the integral, the intensity at every observed location, and the parameter domain. For covariate models, verify separately that the integral is finite and that any claimed derivatives are justified.
The concavity result is useful because any interior solution to the score equations is a global maximizer on a convex parameter domain, provided the likelihood is defined there. Strict concavity can give uniqueness, but only when the data and covariates make curvature strictly negative in every nonzero parameter direction. If there is a direction \(v\ne0\) with \(v^\mathsf{T}z(x)=0\) almost everywhere and \(v^\mathsf{T}z(x_r)=0\) for every observed point, moving the parameter along \(v\) does not change the likelihood. Hence identifiability and uniqueness should not be assumed merely from writing down a log-linear model.
Check Your Understanding
For each question, specify the window and parameter assumptions before calculating.
- A homogeneous Poisson process has four observed points in a window of area \(8\). Find the maximum-likelihood estimate of its intensity and state the log-likelihood score equation.
- For a log-linear model \(\lambda_\beta(x)=e^{\beta z(x)}\), write the score equation when \(z\) is bounded on a finite-area window.
- Explain why finiteness of the integrated intensity at a single \(\beta\) does not, by itself, justify the Hessian formula.
- In the clipped-window example, which region belongs in the integrated-intensity term, and why?
- What parameter-domain issue arises for a constant-intensity model with no observed points if the allowed parameter values are strictly positive?