A unified framework on defining depth for point process using function smoothing

The notion of statistical depth has been extensively studied in multivariate and functional data over the past few decades. In contrast, the depth on temporal point process is still under-explored. The problem is challenging because a point process has two types of randomness: 1) the number of events in a process, and 2) the distribution of these events. Recent studies proposed depths in a weighted product of two terms, describing the above two types of randomness, respectively. In this paper, we propose to unify these two randomnesses under one framework by a smoothing procedure. Basically, we transform the point process observations into functions using conventional kernel smoothing methods, and then adopt the well-known functional $h$-depth and its modified, center-based, version to describe the center-outward rank in the original data. To do so, we define a proper metric on the point processes with smoothed functions. We then propose an efficient algorithm to estimated the defined"center". We further explore the mathematical properties of the newly defined depths and study asymptotics. Simulation results show that the proposed depths can properly rank the point process observations. Finally, we demonstrate the new method in a classification task using a real neuronal spike train dataset.


Introduction
The statistical depth is a method to indicate the centrality of a data point with respect to a data cloud. It can provide a center-outward structure that can be used to understand the empirical distribution, similar to the empirical quantiles in univariate data. The concept of depth was firstly proposed by (Tukey 1975) to handle multivariate data. Since then the notion of depth has been widely studied by mathematicians and statisticians. Different forms of depth were proposed and analyzed, based on different usages and criteria. The types of data that those depths were mainly applied to were multivariate data and functional data.
Some well-known multivariate depths include the half-space depth (Tukey 1975), the convex hull peeling depth (Barnett 1976), the Oja depth (Oja 1983), the simplicial depth (Liu et al. 1990), and the Mahalanobis depth (Liu & Singh 1993). A study by Zuo & Serfling (2000) discussed desirable properties of multivariate depth, which include affine invariance, maximality at the center, monotonicity relative to the deepest points, and vanishing at infinity.
As a generalization of multivariate data from finite dimension to infinite dimension, functional data received a lot of attentions recently and functional depths were also well studied.
To name a few, Cuevas et al. (2007) proposed the h-depth, which applied the Gaussian kernel and L 2 norm to form the depth. Cuesta-Albertos & Nieto-Reyes (2008) extended the idea of half-space depth and defined the random Tukey depth. In 2009, the band depth and modified band depth were introduced by López-Pintado & Romo (2009), which were very commonly used in functional data problems. Similar to the work by Zuo & Serfling (2000), Nieto-Reyes et al. (2016) examined the desirable properties of functional depth, including distance invariance, maximality at center, strictly decreasing with respect to the deepest point, upper semi-continuity, and receptivity to convex hull width across the domain and continuity.
As a special type of data, observation from an orderly temporal point process is made up by an ascending sequence of event times. Such observation contains two types of randomness: the number of events and the time locations of these events. If we treat each observation as a vector of the event times, then the dimension of this vector will be a random variable. Given this dimension, this vector will also be random with ascending entries.
We point out that the study on the notion of depth of point process is relatively new. So far the only existing work in this area was done in Liu et al. (2017), where the depth structure was built in two steps: For an observation s, 1) estimating the probability of getting the number of events P (|s|); 2) given the number of events, compute the conditional depth D(s | |s|). The estimation of probability was done through a normalized Poisson mass function and the conditional depth adopted the Mahalanobis depth for multivariate data. The final depth was the multiplication of these two with a weight power r: D(s) = P (|s|) r D(s | |s|). This depth structure involves both types of randomness and satisfies good mathematical properties such as invariant to time-shift and linear transformation, monotone on rays, and upper semi-continuous. However, the combination of the randomness on the number of events and the event time distribution is not natural that the impact of the number of events on the final depth should be adjusted by the hyper-parameter r. The selection of the hyper-parameter could be tricky because an inappropriate value will make one of the two randomnesses dominates the depth value. In addition, the two-step procedure deals with the number of events and the event time independently. A more desirable method should be able to combine both steps in one framework, where the two types of randomness can be measured at the same time.
To deal with these issues, we propose another approach to process the point process observation: a "transformation" through smoothing method (Wand & Jones 1994). The idea is to smooth the observed point process sequence with an appropriate kernel function.
Through this operation, the event time vector will be replaced by a function curve. In Section 2.1.2, we will see that the point process observation and the function curve are oneto-one matched if the kernel function used for smoothing satisfies certain mild conditions. According to this bijective relation, we are able to apply the methods on functional data to the point process observations, such as metrics on functions and functional depths. In this way, both types of randomness are taken into account under one framework. For individual event time, it determines the location of the kernel function curve. Because we take the sum of the kernel functions in the smoothing procedure, a larger number of events will enhance the smoothed curve vertically. In addition, irrespective of the number of entries in the vector, its smoothed curve will be just one function. That is, vectors with different dimensions can be natuarally compared. Therefore, this approach will not suffer the same issues in the framework proposed by Liu et al. (2017) and it has the advantage to utilize existing depth methods for functional data.
The rest of this manuscript is organized as follows. In section 2, we will introduce the details of transforming a vector to a smoothed curve. We will also propose a proper metric for the space of point process observations, followed by a discussion of the properties for this metric. With this metric, we will continue to define new depth methods for the point process observations. In section 3, we will examine the asymptotic theory of the proposed depths. In section 4, we will apply our methods on simulations and experimental datasets to validate the effectiveness of the new method. Finally in section 5, a summary of the paper will be shown, followed by the future work. All mathematical proofs and algorithmic details are given in appendices in the supplementary material.

Equivalent Representation via Function Smoothing
Our goal is to rank point process observations on a finite interval [0, T ] via a kernel smoothing method. We will at first introduce the type of smoothing kernels for this purpose.

Kernel Functions
The basic idea of kernel smoothing on an orderly temporal point process is to assign a kernel to each observed point and then sum over all the assigned kernels to get a smooth function. Often the kernel is a probability density function of a given distribution such as a Gaussian kernel. The kernel function will depend on the time interval [0, T ] for the point process and we denote it as K(·; T ). In general, we propose to use any kernel function which satisfies the following four basic requirements. These requirements are needed in order to achieve good mathematical properties for the new depth function: 1. Continuous and non-negative: 2. Positive at zero: K(0; T ) > 0; 3. Linear independence with shifting: for any n (n ∈ N) different values t 1 < t 2 < · · · < t n , we have: = 0 for any x ∈ [0, T ] ⇐⇒ α 1 = · · · = α n = 0; 4. Scale invariance: for any x ∈ [0, T ] and constant α > 0, K(αx; αT ) ≡ K(x; T ); We define a kernel function to be "proper" if it satisfies the above four conditions. To have an example of a proper function, we can consider the following Gaussian kernel function: where c 1 , c 2 are two positive constants. This result is given in the following lemma, where the proof is given in Appendix A.
Lemma 1. The Gaussian kernel in Eqn. (1) is a proper kernel.

Proper Metric
Before moving on to the depth function, we look for a proper metric to measure the difference between two observed point processes. Since the two point processes, in general, may have different number of events, multivariate measurements cannot be applied directly due to inconsistency in dimensions. As smoothing functions are included in this depth, here we will use the distance on smoothed processes as the distance on point processes. Let the time interval of the point process be [0, T ]. We define the set of observed point process with given dimension l > 0 as: l=0 Ω l is the space of all point processes. Each process x = (x 1 , x 2 , · · · , x l ) ∈ Ω l with l > 0 can be represented using a Dirac delta function in the form: x = l i=1 δ(· − x i ). Let K(·; T ) be the smoothing kernel. Then the smoothed process is a function on [0, T ] in the form: In general, the space of smoothed processes for l events with l > 0 is We point out that the correspondence between a point process and its smoothed version is one-to-one. This is given in the following lemma (the proof is in Appendix B).
Lemma 2. The smoothing procedure is a bijective mapping from Ω to F. Now we are ready to define a metric on the point process space. This definition is based on the classical L p norm on the smoothed processes. This is formally given as follows: Definition 1. For any two point processes s, t ∈ Ω on [0, T ] and the correspondent f s , f t ∈ F given by Eqn.
(2), we define a distance function d K,p between s and t as: Note that the smoothing kernel K(·; T ) and L p norm can influence the distance value.
However, we can prove that the distance d K,p in Eqn.
(3) is a proper metric on Ω for any K and p. This is given in the following theorem (the proof is in Appendix C).
Theorem 1. If the smoothing kernel K satisfies the four conditions to be proper, then the function d K,p in Definition 1 is a proper metric in the point process space Ω. That is, it satisfies non-negativity, identity of indiscernibles, symmetry, and triangle inequality.
Based on Lemma 2 and Theorem 1, we know that the point process space Ω and smoothed process space F are one-to-one, and an L p norm on F can be used to define a proper metric on Ω. This metric helps solve the problem in Ω where a conventional vector metric cannot be directly used as different processes may have different cardinalities. Based on this result, we need to answer two natural questions: 1. For two point processes with the same cardinality, if the event times in the two processes are very close, will the distance d K,p between them be close, too?
2. Conversely, if the distance d K,p between two point processes is close, will they have the same cardinality and their events times are also close?
Question 1 examines if the d K,p distance is continuous with respect to the event times.
We claim that this is true and the conclusion is stated in Proposition 1 as follows. The detailed proof is given in Appendix D.
Proposition 1. For any k ∈ N, suppose y = (y 1 , y 2 , · · · , y k ) is an observed point process k ), n = 1, 2, · · · be a sequence of processes in Ω k . If lim n→∞ x (n) = y, or equivalently, lim n→∞ x (n) i = y i , i = 1, 2, · · · , k, then for any observed point process z = (z 1 , z 2 , · · · , z l ) ∈ Ω, Question 2 examines the inverse continuity. We claim that this inverse continuity is also true and the result is stated in Proposition 2 (see Appendix E for a detailed proof).

Desirable Properties for Depth on Point Process
In general, a depth definition is defined to provide a measure of centrality of a given data point within a data cloud. To measure if such goal is achieved, one often examines desirable mathematical properties corresponding to the centrality measurement. For instance, Zuo & Serfling (2000) proposed the desirable properties for multivariate depth. Later, Nieto-Reyes et al. (2016) provided the desirable properties for functional depth. The commonly studied properties are (1) linear invariance, (2) vanishing at infinity, (3) maximality at the center, and (4) monotonicity.
Motivated by previous studies on desirable properties for multivariate and functional depths, we propose the following five desirable properties for the depth on point process: Suppose D(s; P S ) is the depth function for observed event time vector s with respect to the probability space (Ω, F , P S ) for a random point process S on interval [0, T ], we expect that D satisfies P1 Continuity: D(s; P S ) is continuous with respect to s in Ω.

P2 Linear invariance: For any event time vector
P3 Vanishing at infinity: When the number of events goes to infinity, the depth should go to 0: D(s; P S ) → 0 as |s| → ∞.
P4 Unique maximum at the center: There exists s c ∈ Ω, such that D(s c ; P S ) = max s∈Ω D(s, P S ).
Also, for any x( = s c ) ∈ Ω, D(x; P S ) < D(s c ; P S ). We refer to s c as the "center".
P5 Monotone decreasing from the center: For any s 1 .s 2 ∈ Ω, suppose the center is s c , if We point out that the "center" in P4 was originally defined as a point with symmetry in the data cloud. For example, when dealing with zero-mean multivariate normal samples, we will consider the origin as center, and the corresponding Mahalanobis depth is uniquely maximized at this point. However, for point process data in a finite domain [0, T ], it is difficult to define a geometrically symmetric central point in Ω. Notwithstanding, we expect the proposed depth has a unique maximum point, which results in our notion of "center".
Here we examine if the above five properties are satisfied by the generalized Mahalanobis depth on point process (Liu et al. 2017): 1) If two point processes are close to each other, they should have the same number of events and close event times in the time order. As the classical Mahalanobis depth is continuous with respect to the input vector, the generalized Mahalanobis depth on point process satisfies P1; 2) Linear transformation on the time will not change the number of events, so P (|s|) r will keep the same. Because the classical Mahalanobis depth is linear invariant, P2 is also satisfied; 3) The Poisson mass function will go to 0 as the input goes to infinity, so P (|s|) → 0 as |s| → ∞. Thus, P3 holds, too; 4) D(s | |s|) has maximum value to be 1. However, this maximum in general can be achieved at multiple processes. For example, if the total intensity Λ is an integer, both the population mean with dimension Λ and the population mean with dimension Λ − 1 will achieve the maximum depth. Thus, P4 does not hold as the solution to the maximum is not unique; 5) As a result, P5 does not hold either. In summary, the generalized Mahalanobis depth on point process only satisfies P1 to P3.
Similar to the discussion about the Mahalanobis depth, we will refer to P1 -P5 to explore these properties for other depth functions in the following sections of this paper.

h-depth on Point Process
We have defined a proper metric in the point process space Ω. This metric is based on a smoothing procedure where a point process is equivalently represented by a function in F.
In this section, we will exploit a commonly used functional depth, called h-depth, to define an "h-depth" for point processes in Ω.
At first, we review the notion of h-depth (Cuevas et al. 2007) for functional data. Assume X is a functional random variable in the probability space (Λ, F Λ , P Λ ), where Λ is a subset of L 2 ([0, T ]). Let · denote the conventional L 2 norm. For any function z ∈ Λ, its h-depth is defined as: where G h (·) = exp(− t 2 2h ) is a modified Gaussian kernel with parameter h (> 0) and X is a random function on the probability space (Λ, F Λ , P Λ ). To simplify notation, we use HD(z) to represent HD(z; P Λ ) when the measure P Λ is implicitly known. Given a set of i.i.d.
random functions X 1 , X 2 , · · · , X n ∈ Λ, the sample version of the h-depth of z is given as ).
One important property of the h-depth in Eqn. (4) is its continuity. This is given in the following lemma, where the proof can be found in Appendix F.
Lemma 3. If the functional random variable X satisfies E X < ∞, then HD(·) is a continuous function on Λ.
Based on the notion of h-depth on functional data (in Eqn. (4)) and the proper metric on point process (in Eqn. (3)), we can formally introduce the h-depth to point process observations as follows: Definition 2. Let S be a random point process on [0, T ] in the probability space (Ω, F , P S ).
The h-depth of any s ∈ Ω is defined as: where f s and f S are smoothed curves for s and S (by Eqn. (3)) and G h (t) = exp(− t 2 2h ).
Definition 2 provides the population version of the depth. In practice, we should follow the sample version of h-depth on the given observations, which is given below.
be a sample of event time vectors from a point process on [0, T ]. The empirical h-depth of any s ∈ Ω is defined as: where f s and f S i are smoothed curves for s and With the definitions for both the population and empirical versions of the h-depth on point process, we now examine the mathematical properties based on the discussion in section 2.2. The basic result can be summarized in the following proposition. • P3: D(s; P S ) → 0 when |s| → ∞.
The detailed proof is given in Appendix G. Note that P4 and P5 are not satisfied because h-depth in general may not have a unique maximum point in the data cloud.

Modified h-depth on Point Process
The center has been a critical notion in statistical depths. However, as we have pointed out in Sec. 2.3, the h-depth for a point process may not have a center. In this section, we propose to modify the h-depth by including a "center"-based process in the definition.

Definition and Properties
Many commonly used depth functions have the "center" as the point with maximum depth.
Based on this idea, we propose a center-based new depth function on point process using the smoothing method. The formal definition is given as follows: Definition 4. Let s c be a given "center" point process in Ω on [0, T ]. For any s ∈ Ω, its center-based h-depth is defined to be: where h > 0 is a parameter. f s and f sc are smoothed processes of s and s c , respectively. )], with X to be a functional random variable. If the order of the exponential function and the expectation is switched, then we have ). In this way we obtain the center-based h-depth D(s; s c ).
Remark 2: The L 2 norm in the definition can be generalized to L p , with 1 ≤ p < ∞.
With a center process given in the definition, this modified h-depth for point process is expected to satisfy more desirable mathematical properties than the classical one. Indeed, all the desirable properties in section 2.2 are satisfied and the details are given in the following proposition. • P5: For any s 1 , s 2 ∈ Ω, if d K,2 (s 1 , s c ) < d K,2 (s 2 , s c ), then D(s 1 ; s c ) > D(s 2 ; s c ). The detailed proof for each property is provided in Appendix H. Note that this is a clear advantage over the classical h-depth method, where only properties P1-P3 are satisfied.

Estimation of the Center
In practice, there is no prior knowledge of the center of point process, so a proper estimation of the center will be needed. In this study, we adopt the common notion of the Karcher mean in a metric space (Grove & Karcher 1973). Since we have a well-defined metric between any two point processes in Definition 1, the Karcher mean can be defined as: Definition 5. Let S be a random point process on [0, T ] in the probability space (Ω, F , P S ).
The Karcher mean of S is defined to be where d K,2 (·, ·) is the metric with kernel function K(·; T ) and L 2 norm in Definition 1.
Note that, in general the solution of the minimization may not be unique, so µ K is a set of points. Based on the Karcher mean in Definition 5, its empirical version is given as follows: The empirical Karcher mean also may not be unique, and thereforeS K is a set of solution points, too. For simplicity, we define the following notation: where SSD stands for "Sum of Squared Distance" and t ∈ Ω is the event time vector. With this definition, the empirical Karcher mean is justS (N ) K = arg min t∈Ω SSD(t; S 1 , S 2 , · · · , S N ). For any element inS We will introduce a theoretically-proven convergent algorithm to estimateS K provides an effective searching range. Such algorithm is a combination of simulation annealing and line search. We will at first briefly review the simulation annealing and line search, and then introduce an approach to combine these two methods.
• RJMCMC annealing Simulated annealing (Geman & Geman 1984, Van Laarhoven & Aarts 1987) is a method that makes optimization feasible if we are able to generate samples. To find the global minimum of function f (x), we can build a density function π where T i is a pre-determined decreasing sequence of positive numbers called temperature and lim i→∞ T i = 0 (e.g. T i = C/ log(1 + i) for constant C > 0). Under weak regularity assumptions on f , the density π i will be concentrated on the global minimum points of f . After a large number of iterations of computing π i and sampling from it, we will have samples passing the global minimum point.
Simulation annealing provides a method to find the minimum of a function, but to apply it on the optimization of SSD, a tool that is able to generate vectors with different dimensions will be needed. In this paper, we adopt the reversible-jump Markov Chain Monte Carlo (RJMCMC) method, which was first introduced by Green (1995). Its main idea is the Metropolis Hasting method (Metropolis et al. 1953), but generalizes to allow the candidate to have different dimensions. In this way, it can be used to get samples from densities where the dimension of the random vector is not fixed. One condition to use the RJMCMC method is that the number of possible dimensions should be finite.
According to Theorem 2, the range of dimension is 1, · · · , D (N ) K , so RJMCMC is feasible for the center estimation.
Combining the two methods, we are able to find the Karcher mean through sampling from T i SSD(t; S 1 , S 2 , · · · , S N )]. Because of RJMCMC, the annealing process can search over different dimensions automatically. An example RJMCMC annealing algorithm is provided in Appendix J, where the input to the algorithm can be adjusted before application.
• Line search According to Theorem 2, because the dimensions to search is finite 1, 2 · · · D (N ) K , optimizing within each dimension and comparing the output across dimensions will be another feasible approach, which corresponds to the "line search" idea. An advantage of this method is that the gradient of the SSD can be computed explicitly with the input dimension fixed, as shown in the following proposition.
Proposition 5. If the dimension of the input event time vector t is given and the distance function is d k,2 , then the gradient of SSD with respect to t will be: where t and S i are event time vectors. x is scaler and |t| is the dimension of vector t.
All operations are done element-wise.
If the modified Gaussian kernel K G (t; T ) = c 1 e − c 2 T 2 t 2 is used to smooth the process, then the gradient can be further simplified to: with the function g(·, ·) to be where Φ is the cumulative distribution function for standard Normal distribution N(0, 1); 1 d = [1, 1, · · · , 1] T ∈ R d . All operations are done element-wise.
The proof and computational details are provided in Appendix L. Due to this property, gradient-based methods, such as stochastic gradient descent, can be applied to efficiently do the optimization within a given dimension. Then a comparison of the output will return the solution to the minimum SSD, which is the empirical Karcher mean. An example line search algorithm based on gradient method is provided in Appendix K.
• Combined method Both RJMCMC annealing method and line search method can be used to estimate the empirical Karcher mean, but they have some disadvantages. For the RJMCMC annealing, although it can search over dimensions automatically, it converges slowly to the optimal solution. For line search, it can converge fast to the solution given the dimension, but it needs to search every dimension in the dimension range. As a result, the time cost for the entire line search algorithm is also very large.
Based on their characteristics, if we combine the two methods: do RJMCMC annealing first as a "pre-train" process to locate a range of optimal dimensions and then use line search to do optimization within the narrowed dimension range, the computation efficiency is expected to be significantly improved. Moreover, the output event time vectors from the RJMCMC annealing can be used as the initial values for the line search.
In this way, the number of iterations in the line search will be reduced and it can also help avoid local optimal points. The combined algorithm following this idea is shown in Algorithm 1. The application of this algorithm on simulations and real data will be shown in section 4.

Asymptotic Theory
In this section, we will show that the empirical depth converges to the population depth as the sample size goes to infinity for h-depth and modified h-depth, defined in section 2.
Algorithm 1 Combined method to find empirical Karcher mean Input: the observed event time vectors S 1 , S 2 , · · · , S N .
For RJMCMC annealing: The maximum number of iterations n max ; Initial value x 0 ; Initial dimension k 0 to be the dimension of x 0 ; the upper bound of dimension D For combination: The number of dimensions to be kept in RJMCMC annealing: d r ∈ N.
(1) Pre-train: Use RJMCMC annealing (such as Algorithm 2) with input S 1 , S 2 , · · · , S N , n max , x 0 , k 0 and D is the corresponding sample version with sample size N. Then, for any input s ∈ Ω, This theorem implies the convergence of sample h-depth to the population h-depth, and can be easily proven using the Strong Law of Large Numbers.
For the modified h-depth, the convergence of the depth value will be based on the convergence of the estimated center. As defined in section 2.4.2, the Karcher mean can be treated as the solution in Ω that minimizes the average of squared distances, where the distance is d k,2 with L 2 norm. In fact, we can generalize the Karcher mean definition by using any d k,p distance for p ≥ 1, given in the following form: Definition 7. Let S be a random point process on [0, T ] in the probability space (Ω, F , P S ).
The generalized Karcher mean of S is defined to be where d K,p (·, ·) is the metric with kernel function K(·; T ) and L p norm where p ≥ 1.
It is apparent that Definitions 5 and 6 are special cases with p = 2 in the generalized definitions, respectively. Then the convergence theorem is given as follows: 1. The minimum average of squared distance in {S i } N i=1 will converge to the minimum expected squared distance in (Ω, F , P S ) almost surely, in other words, 2. The empirical generalized Karcher mean converges to the generalized Karcher mean almost surely, in other words,S (N ) The almost-surely convergence from setS is the Kuratowski upper limit (Kuratowski 2014):

Application Results
In this section, we will apply the proposed depth framework, h-depth and modified h-depth, on point process observations in simulations and real experimental data.

Homogeneous Poisson Process
We will at first illustrate the depth methods on observations from HPP(λ) in a finite time interval, where λ is a constant mean value. 100 independent realizations from HPP(0.045) on [0, 100] are generated, where the raster plot of these processes is shown in Figure 1(a).
Basically, the number of events in each process follows a Poisson distribution with mean 4.5 (varying from 1 to 10 in these 100 realizations), and the event time in every observation is uniformly distributed on [0, 100]. We then smooth the realizations using a kernel in Eqn.
Based on the h-depth and modified h-depth in section 2, we compute the sample depth value for each given process where the constant h is set to T = 100. For the modified h-depth, the center is determined in four ways: set by prior information, estimated by RJMCMC annealing, line search and combined method. As the Poisson process is homogeneous, we can intuitively set the center as (20,40,60,80). The output of the center estimation is shown in Table 1. We can see that all the estimated centers have lower SSD than the intuitive center. In particular, the combined method has the superior performance. It has the smallest SSD and also the smallest time cost. We will only use this method for the modified h-depth in comparison with other depth methods.  We then display the corresponding smoothed processes in Figure 3. Note that the generalized Mahalanobis depth is not included as it is not based on smoothing functions. All three plots demonstrate a center-outward decreasing depth value structure that observations whose smoothing curves are closer to the middle of all the curves will have larger depth values. Overall, these simulation results are reasonable and consistent with the basic notion of center-outward ranks.

Inhomogeneous Poisson Process
We will further apply the depth on point process using observations from IPP[λ(t)] to evaluate the perforamnce. The rate function λ(t) is chosen as λ(t) = 3φ(t; 25, 10) + 2φ(t; 75, 10) and the event time interval is still [0, 100], where φ(·; µ, σ) is the density for normal distribution with mean µ and standard deviation σ. As λ(t) contains two peaks, the simulated event times will mainly distribute around 25 and 75. The kernel function is also the modi-  Table 2. It can be seen that the three estimated centers are very similar. The line search method and the combined method both have the lowest SSD, whereas the latter one is much more efficient. That is, the combined method results in the superior performance, and therefore we will only use this method for the modified h-depth in comparison with other depth methods.   (b) By the modified h-depth with center estimated using the combined method.

Method
After computing the depth for all the observations, we obtain the ranking result. Figure   5 shows the top 5 and bottom 5 ranked observations. Similar to the HPP case, we include the result from the generalized Mahalanobis depth with r = 1 for comparison. It can be seen that the outputs from h-depth and modified h-depth are similar -observations with event times around 25 and 75 and number of events around 5 have larger depth. In contrast, for the Mahalanobis depth, we can also observe that processes with number of events around 5 have larger depth. However, rather than concentrating at two peak time locations, the top 5 processes have event time more uniformly distributed. It indicates that h-depth and modified h-depth can better capture the pattern inside the sample. If further taking into account the true λ(t) that the peak at 25 is higher than the peak at 75, we should expect to see more event times concentrated around 25 than around 75. Comparing In summary, based on the result from the HPP and IPP simulations, we can see that the proposed h-depth and modified h-depth both can properly build a center-outward rank on the given point process data.

Experimental data application
We will examine the proposed depth framework using a real spike train recording, where the data were included in the Quantitative Single-Neuron Modeling Competition 2009 (Naud et al. 2009) and were accessible at http://dx.doi.org/10.6080/K0PN93H3. The detailed experiment description can be found in (Carandini et al. 2007).
Briefly, the experiment was performed on rhesus monkeys that retinal input (visual stimulus) was applied and extracellular potentials were recorded for both the retinal (presynaptic) and the geniculate (post-synaptic) simultaneously. There were 10 seconds stimulus, where the first 5 seconds were the same for all trials and the last 5 seconds were unique for each trial. In total, 76 trials were performed. Since the task in the competition was to predict the post-synaptic spikes given the pre-synaptic spikes, only 38 pairs of (post-synaptic, pre-synaptic) were fully given. We will use these spike train data to test the smoothing depth framework through a classification task: In the training set, we have labeled observations for pre-synaptic (abbreviated as "pre-group") and post-synaptic (abbreviated as "post-group"). In the testing set, we compute the depth for each observation in the two groups. The group with a larger depth value will be the predicted label. Before the classification, an exploratory data analysis is performed. The spike times are within the 10 seconds range and the point process plot is shown in Figure 7. 10 example pairs of pre-group, post-group observations are shown. As stated in the experiment description, the stimuli were the same for the first 5 seconds and different for the last 5 seconds. This well explains why the signals have very similar pattern in the first 5 seconds.
Therefore, it is reasonable to separate the trial into two parts, first 5 seconds and last 5 seconds, and we can smooth them in different ways.
After separation, we apply the modified Gaussian kernel to smooth the observations with c 1 = 1, c 2 = 100 for first 5 seconds and c 1 = 1, c 2 = 50 for last 5 seconds. The smoothing curve plot is given in Figure 8. Using the combined method in section 2.4.2, we estimate the Karcher mean for the two groups and the corresponding smoothing curves are shown as thick lines. It can be seen that the Karcher mean curves are able to catch the shape patterns in the first 5 seconds as the observations are nearly aligned, while in the last 5 seconds, the Karcher mean is nearly flat since the observations are very noisy. The classification result is shown in Table 3. All methods are able to provide reasonable classification power on the spike train data. The overall results on the first 5 seconds is better than the last 5 seconds. For the last 5 seconds, it can be seen that the modified h-depth has the best accuracy (76.32%) and F1 scores ( Table 3: The classification result. F1 score (pre) and F1 score (post) are computed according to the recall and precision on pre-group classification and post-group classification.
The boldface indicates best results in each column.
best accuracy (81.58%), while the modified h-depth is better on the post group (F1 score 82.05%) and the other two are better on the pre-group (F1 score 81.58%). In summary, we can see that the smoothing depth framework provides relatively accurate classification in the practical point process observations.

Summary and Future Work
In this paper, we at first proposed a kernel smoothing representation for point process observations. We found that the smoothing procedure builds a bijective mapping between them. We then defined a proper metric distance bewteen the smoothing curves and transformed the problems from point processes to smooth functions. Based on the notion of h-depth on functions, we proposed two methods to form the depth structure on point process observations, i.e. the classical h-depth and the modified h-depth with center given.
We then proposed a center-estimation method for the modified h-depth, with the idea of minimizing the Karcher mean through a combination of the RJMCMC annealing and line search. During the tests on the simulated data, both classical h-depth and modified hdepth resulted in reasonable outcomes: a center-outward decreasing depth structure that observations with smoothing curves closed to the center will have higher depth. In the real neuronal spike train data, we tested the new depth structures with a depth-based classification task and both methods resulted in accurate classifications.
Our investigation is only the starting point of the smoothing curve exploration and there is much to be further studied in the future. We have proposed the modified Gaussian kernel that satisfies the four basic requirements. We will explore other kernel functions that also satisfy the requirements. In the modified h-depth, a combination of RJMCMC annealing and line search has been given as the method to estimate the Karcher mean, which serves as the "center". There should be other ways to define the "center" and other estimation methods. These are all very good topics to further develop depth function on point process. Finally, more functional depth methods, other than h-depth, can be adopted to smoothed point processes, such as the band depth and modified band depth.
Their ranking performance needs thorough investigations.