Interacting cluster point process model for epidermal nerve fibers
Introduction
Epidermal nerve fibers (ENFs) are sensory fibers that originate in the dermis and terminate in the epidermis, the outermost layer of the skin. Microscopic analysis of suction-induced skin blisters as a diagnostic tool for neuropathy was introduced by Kennedy et al. (1999), and point process models for this type of data have been studied by Waller et al., 2011, Myllymäki et al., 2012, Olsbo et al., 2013 and Andersson et al. (2016). The data consist of the locations of the origins of nerve fiber bundles as well as the final location of each ENF. The original data are three-dimensional and come in voxels of base 330 × 432 microns, with depth varying between 20 and . Like the studies previously mentioned, we will consider only the projected planar process. One of the medical interests in this type of data lies in studying the effect of diabetes on the spacing of ENFs, where the hypothesis is that nerve fiber loss due to diabetic neuropathy does not result in random nerve removal, but rather in a spatial pattern different from the pattern in healthy patients (Kennedy et al., 1999). In order to evaluate this hypothesis, it is important to have a model for the spatial pattern of ENFs in healthy patients, and this article is a step in this direction.
Our datasets consist of ENF trunks (base points) and their associated terminal nodes of ENFs (end points). Both sets of points are identifiable from images of ENFs where base points are defined as the intersection of ENFs with the dermal–epidermal junction. After entering the epidermis, the ENFs “branch” spreading in different directions ending in terminal nodes which are responsible for sending the sensations of heat and pain to the central nervous system. In this paper, we focus on samples from the thigh from each of four subjects distinguished by labels 171, 224, 230 and 256, and shown in Fig. 1. For example, the panel in Fig. 1 corresponding to patient labeled 171 shows 54 origin points and 95 ENFs in the sample. Although the dataset contains data in three dimensions considering the winding paths of the individual ENFs, or their depth within the skin, our interest lies solely in the spatial pattern of ENF coverage across the skin given by the distribution of base points and corresponding end points in two dimensions. Therefore, we will use the projections of the observed patterns into a box in .
The datasets we analyze have been previously modeled by Olsbo et al. (2013). They propose two spatial point process models for the spatial structure of ENFs. In their first model, both the parent and the offspring are Poisson processes (uniform base and end points) and there is a probabilistic assignment of end points to base points. Their second model (non-orphan cluster) assumes that locations of base points follows a homogeneous Poisson process and then, conditional on the locations of these base points, they model the locations of end points through as the joint distribution of a von Mises distribution for the direction to the end point with respect to the axis and an independent gamma distribution for the distance between the end point and the parent (base) point respectively. Also, they assume independence across base points.
The model we propose is related to the non-orphan cluster model but incorporates a possible repulsion among the base points as well as within the end points in each cluster. The model is the invariant measure of a birth and death process where the birth rate of base points depends on the entire configurations of end points present at the time of the birth whereas the death of a base point (along with all the corresponding end points) is independent of the configuration. Conditionally on the birth of a base point, its corresponding end points follow a repulsive point process.
The manuscript is organized as follows: in Section 2 we provide an overview of interacting point process models, with the proposed interacting cluster point process model being detailed in Section 3. In Section 4 we detail computational approximations necessary for efficient evaluation of the likelihood. Bayesian estimation of the parameters is described in Section 5. In Section 6 we present a simulation study to assess the performance of the estimation procedure, while in Section 7 we apply the tools to the ENF data and verify the adequacy of the model with some goodness-of-fit quantities. Section 8 discusses extensions.
Section snippets
Point processes
Let be a subset of , its Borel sets and consider the set of locally finite counting measures on endowed with the -algebra generated by the mappings . A point process is a -valued random variable. That is, there exists a probability space and a map , is the set of natural numbers including zero, such that
- 1.
For each , is a locally finite counting measure on
- 2.
For each bounded , is an -valued random variable.
There are
An interacting cluster point process as a model for ENFs
The central idea of this paper is to propose a repulsive cluster point process that models the ENFs pattern. The simplest approach is to consider invariant measures of a continuous-time birth and death cluster processes. We assume that our sample represents a snapshot of the process at a fixed time long enough for the process to reach stationarity. There are many birth and death processes that lead to the same stationary distribution.
Initial within-subject temporal replicates of (healthy) ENF
Numerical approximations
To find the maximum likelihood estimator based on (14) poses several problems:
- (i)
the normalizing constant is untractable and cannot be computed;
- (ii)
the integral of , the measure of the shadow, is not easily computed, and
- (iii)
there is an extra randomness due to latent variables that have to be taken into account.
To overcome (i) and (ii) we will use some numerical approximations as described in the next section. To overcome (iii) inference will be performed under the Bayesian
Bayesian approach
Inference can be performed under the Bayesian paradigm and we now discuss the prior specification of the parameters in the model, as well as the resultant posterior distribution and how to perform inference. Let be the parameter vector of the model.
Simulation studies
In the previous discussion we did not consider the radii and as parameters. From a simulation study we propose to estimate the radii from the data as: as twice the minimum distance among all base points, and as the minimum distance among the end points of each cluster. For a justification of these choices, please refer to the Supplementary Materials.
The process’ window can be rescaled to be within a unit window. Therefore, for the
Application to ENF data
In this section, we fit the model separately for each subject, ran the process for 50,000 steps and discarded 2000 as burn-in. The prior distributions over the parameters are the same as in the simulation. Before running the MCMC, we rescaled the window so the parameters can be compared across different patients. We do so by dividing each side of the bounding box by , the length of the largest side of the bounding box. The radii and , in the original box are also divided by , and we
Discussion
We propose a cluster model which is a generalization of the non-orphan cluster model used by Olsbo et al. (2013). The proposed model is a flexible descriptor of the nerve bundles and is richly parametrized to offer tools for diagnostic pathologies. The main novelty of this model is the possibility of incorporating a possible repulsion among the base points as well as within the end points in each cluster. Differently from Olsbo et al. (2013), we perform the inference procedure under the
Acknowledgments
The authors thank Aila Särkkä for fruitful discussions and suggestions that improved the manuscript as well as the associate editor and two anonymous reviewers for their valuable comments. Nancy Garcia was supported by FAPESP, Brazil grants number 2014/26419-0 and 2017/15306-9, and CNPq, Brazil grant 302598/2014-6. Peter Guttorp was supported by NRC, Norway grant number 240838. Guilherme Ludwig was supported by FAPESP, Brazil grants number 2016/09390-4 and 2018/09877-6.
References (33)
- et al.
Perfect simulation for interacting point processes, loss networks and ising models
Stochastic Process. Appl.
(2002) - et al.
The innervation of human epidermis
J. Neurol. Sci.
(1993) - et al.
Development and evaluation of spatial point process models for epidermal nerve fibers
Math. Biosci.
(2013) - et al.
Morphological features of nerves in skin biopsies
J. Neurol. Sci.
(2006) - et al.
Discovering early diabetic neuropathy from epidermal nerve fiber patterns
Stat. Med.
(2016) - et al.
Spatial Point Patterns: Methodology and Applications with R
(2015) - et al.
P-values for composite null models
J. Amer. Statist. Assoc.
(2000) Statistical Decision Theory and Bayesian Analysis
(1985)- et al.
Bayes and Empirical Bayes Methods for Data Analysis, Second Edition
(2000) - et al.
An Introduction to the Theory of Point Processes. Vol. I. Probability and its Applications. New York)
(2003)
Markov Chain Monte Carlo: Stochastic Simulation for Bayesian Inference
Two simple examples for understanding posterior p-values whose distributions are far from uniform
Electron. J. Stat.
Model checking and model improvement
Inference from iterative simulation using multiple sequences
Statist. Sci.
Likelihood inference for spatial point processes
Constrained Monte Carlo maximum likelihood for dependent data
J. R. Stat. Soc. Ser. B Stat. Methodol.
Cited by (4)
Statistical modeling of diabetic neuropathy: Exploring the dynamics of nerve mortality
2023, Statistics in MedicinePairwise interaction Markov model for 3D epidermal nerve fibre endings
2022, Journal of MicroscopyMarked point process analysis of epidermal nerve fibres
2021, Journal of Microscopy