Spatio-temporal point processes with deep non-stationary kernels
Authors: Zheng Dong
Xiuyuan Cheng
Yao Xie
Venue: ICLR 2023 poster (2023)
Keywords: point process
neural network
non-stationary kernel
low-rank model

Abstract: Point process data are becoming ubiquitous in modern applications, such as social networks, health care, and finance. Despite the powerful expressiveness of the popular recurrent neural network (RNN) models for point process data, they may not successfully capture sophisticated non-stationary dependencies in the data due to their recurrent structures. Another popular type of deep model for point process data is based on representing the influence kernel (rather than the intensity function) by neural networks. We take the latter approach and develop a new deep non-stationary influence kernel that can model non-stationary spatio-temporal point processes. The main idea is to approximate the influence kernel with a novel and general low-rank decomposition, enabling efficient representation through deep neural networks and computational efficiency and better performance. We also take a new approach to maintain the non-negativity constraint of the conditional intensity by introducing a log-barrier penalty. We demonstrate our proposed method's good performance and computational efficiency compared with the state-of-the-art on simulated and real data. 

Published as a conference paper at ICLR 2023
SPATIO-TEMPORAL POINT PROCESSES WITH DEEP NON- STATIONARY KERNELS
Zheng , Xiuyuan Cheng", Yao Xie i H. Milton Stewart School of Industrial and Systems Engineering, Georgia Institute of Technology Department of Mathematics, Duke University yao. x xie@' sye. gatech.edi
ABSTRACT
Point process data are becoming ubiquitous in modern applications, such social networks, health care, and finance. Despite the powerful expressiveness the popular recurrent neural network (RNN) models for point process data, they may not successfully capture sophisticated non-stationary dependencies in the data due to their recurrent structures. Another popular type of deep model for point process data is based on representing the influence kernel (rather than the intensity function) by neural networks. We take the latter approach and develop new deep non-stationary influence kernel that can model non-stationary spatio-temporal point processes. The main idea is to approximate the influence kernel with a novel and general low-rank decomposition, enabling efficient representation through deep neural networks and computational efficiency and better performance. We also take new approach to maintain the non-negativity constraint of the conditional intensity by introducing a log-barrier penalty. We demonstrate our proposed method's s good performance and computational efficiency compared with the maate-of-th-aar art on simulated and real data.
INTRODUCTION
Point process data, consisting of sequential events with timestamps and associated information such as location or category, are ubiquitous in modern scientific fields and real-world applications. The distribution of events is of great scientific and practical interest, both for predicting new events and understanding the events' generative dynamics (Reinhart, 2018). To model such discrete events in continuous time and space, spatio-temporal point processes (STPPs) are widely used in diverse range of domains, including modeling earthquakes (Ogata, 1988; 1998), the spread of infectious discases (Schoenberg et al., 2019; Dong et al., 2021), and wildfire propagation (Hering et al., 2009).
modeling challenge is to accurately capture the underlying generative model event occurrence in general spatio-temmpepppeee temporal point processes (STPP) while maintaining the model efficiency. Specific parametric forms of conditional intensity are proposed in seminal works of Hawkes process (Hawkes, 1971; Ogata, 1988) to tackle the issue of computational complexity in STPPs, which requires evaluating the complex multivariate integral in the likelihood function. They use an exponentially decaying influence kernel measure the influence of a past event over time and assume the influence of all past events is positive and linearly additive. Despite computational simplicity (since the integral of the likelihood function is avoided), such parametric form limits the model's 's practicality in modern applications.
Recent models use neural networks in modeling point processes capture complicated event occurrences. RNN (Du et al., 2016) and LSTM (Mei and Eisner, 2017) have been used by taking advantage of their representation power and capability capturing event temporal dependencies. However, the recurrent structures of RNN-based models cannot capture long-range dependency (Bengio et al., 1994) and attention-bacee structure (Zhang et al., 2020; Zuo et al., 2020) is introduced to address such limitations of RNN. Despite much development, existing models still cannot sufficiently capture aratio-tempora non-stationarity, which are common in real-world data (Graham et al., 2013; Dong al., 2021). Moreover, while RNN-type models may produce strong prediction performance, the models consist of general forms of network layers and the modeling power relies on the hidden states, thus often not easily interpretable.

---

Published as a conference paper at ICLR 2023
200
200
20
t-t'
(a) Kernel matrix of k<t', t) with rank 298.
(b) Kernel matrix of k(t', with rank 7.
Figure 1: An example: equivalent representation of kernel by "displacement", from k(t', t) to k f') can significantly decrease the rank of the kernel function: from 298 to 7, where e' and t represent the historical event time and the current time, respectively. On fitting the one-dimensional synthetic data generated by the model parameterized with (t',t t') outperforms the model parameterized with (l', l). See Section 5 and Appendix C for experiment and formulation details.
A promising approach to overcome the above model restrictions is point process models that combine statistical models with neural network representation, such as Zhu et al. (2022) and Chen et al. (2020), to enjoy both the interpretability and expressive power of neural networks. In particular, the idea is to represent the (possibly non-stationary) influence kernel based on a spectral decomposition and represent the basis functions using neural networks. However, the prior work (Zhu et al., 2022) is not specifically designed for non-stationary kerne and the low-rank representation can be made significantly more efficient, which is the main focus of this paper.
Contribution. this paper, we develop a non-stationary kernel (referred to as DNSK) for (possibly non-stationary) spatio-temporal processes that enjoy efficient low-rank representation, which leads to much improved computational efficiency and predictive performance. The construction is based on an interesting observation that by reparameterize the influence kernel from the original form of k(t', ft. (where is the historical even time, and t is the current time) to an equivalent form k(t', i t') (which thus is parameterized by the displacement t - t' instead), the rank can be reduced significantly, as shown in Figure 1. This observation inspired us to design a much more efficient representation of the non-stationary point processes with much fewer basis functions to represent the same kernel.
In summary, the contributions of our paper include
We introduce an efficient low-rank representation of the influence kernel based on a novel "dis- placement" re-parameterization. Our representation can well-approximate a large class of general non-stationss stationary influence kernels and is generalizable to spatio-tempont kernels (also potentially to data with high- nigh-dimmesional marks). Efficient representation leads to lower computational cost and better prediction power, as demonstrated our experiments.
In model fitting, we introduce log-barrier penalty term in the objective function to ensure the non negative conditional intensity function so the model statistically meaningful, and the problem is numerically stable. This approach also enables the model to learn general influence functions (that can have negative values), which a drastic improvement from existing influence kernel-based methods that require the kernel functions to be non-negative.
Using extensive synthetic and real data experiments, we show the competitive performance of our proposed methods in both model recovery and event prediction compared with the state-of-the-art, such as the RNN-based and transformer based models.
1.1 RELATED WORKS
The original work of A. Hawkes (Hawkes, 1971) provides classic self-exciting point processes for temporal events, which express the conditional intensity function with an influence kernel and base rate. Ogata (1998) proposes parametric form of spatio-temporal influence kernel which enjoys strong model interpretability and efficiency. However, such simple parametric forms own limited expressiveness in characterizing the complex events" dynamic in modern applications (Zhu et al., 2021; Liao et al., 2022).
Neural networks have been widely adopted in point processes (Xiao et al., 2017: Chen al., 2020; Zhu et al., 2021). Du et al. (2016) incorporates recurrent neural networks and Mei and Eisner (2017) continuous-tima invariant of L.TMM to model event influence with exponential decay over time.
2

---

Published as a conference paper at ICLR 2023
These RNN-based models may be unable to capture complicated event dependencies due to the recurrent structure. Zhang et al. (2020); Zuo et al. (2020) introduce self-attentive structures into point processes for their capability memorize long-tern -eg-term influence by dealing with an event sequence as a whole. The main limitation that they assume a dot-product-based score function and assume linearly decaying of event influence. Omi et al. (2019) propose fully-connected neural network to model the cumulative intensity function to go beyond parametric decaying influence. However, the event embeddings are still generated by RNN, and fitting cumulative intensity function by neural networks lacks model interpretability. Note that all the above models tackle temporal events with categorical marks, which are inapplicable in continuous time and location space.
Recent works adopt neural networks in learning the influence kernel function. The kernel introduced in Okawa et al. (2021) uses neural networks to model the latent dynamic of time interval but still assumes an exponentially decaying influence over time. Zhu et al. (2022) proposes kernel representation using spectral decomposition and represents feature functions using deep neural networks to harvest powerful model expressiveness when dealing with marked event data. Our method considers an alternative novel kernel representation that allows the general kernel to be expressed further low-rankly.
BACKGROUND
Spatio-temporal point processes (STPPs) (Reinhart, 2018; Moller and Waagepetersen, 2003) have been widely used to model sequences of random events that happen in continuous time and space. Let H = {44, 6;) } denote the event stream with time I E [0. 7] C R and location 8; E S C Rds of ith event. The event number is also random. Given the observed history H, = ((t,,s) Httt t}
before time t, an STPP is then fully characterized by the conditional intensity function
E [ + A/ x B(s, As)) H4.]
A (1, 8 Hi)
lim Ac40,5a10
(1)
[B(s, As)
where B(s, As) is ball centered at s Rds with radius As, and the counting measure N is defined as the number of events occurring in [t, At] x B(s, As) C Ids++1, Naturally A (f, s(Ht) 2 0 for any arbitrary t and s. the following, we omit the dependency on history Hi and use common shorthand X(t, s). The lig-likelihood of observing 76 on [0, T) x S is given by (Daley et al., 2003)
Jo / A(t. s)dsdt
((H)- = Σ log A (t,, 8;)
Neural point processes parameterize the conditional intensity function by taking advantage of recurrent neural networks (RNNs). In Du et al. (2016), an input vector x, which extracts the information of event t; and the associated event attributes m, (can be event mark or location) is fed into the RNN. A hidden state vector h; is updated by h, = p(h;-1, z;). where p(-) is a mapping fulfilled by recurrent neural network operations. The conditional intensity function on (t,, f4+1 is then defined as X(t) = 5(t, h2)), where is an exponential transformation that guarantees a positive intensity. In Mei and Eisner (2017) the RNN replaced by a continuous-time LSTM module with hidden states (tt) defined on [0. 7'] and a Softplus function 8. Attention based models are introduced in Zuo et al. (2020); Zhang al. (2020) to overcome the inability of RNNs to capture sophisticated event dependencies due to their recurrent structures.
Hawkes process (Hawkes, 1971) is a well-known generalized point process model. Assuming the influences from past events are linearly additive, the conditional intensity function takes the form of
X(t. 4 +
Σ k('',,,'s,s), (M,,'))H.
(3)
where & is an influence kernel function that captures event interactions. Commonly the kernel function assumed to be stationary, that is, l only depends on - l' and s', which limits the model expressivity. In this work, we aim to capture complicated non-statioparity in spatio- tentorll event dependencies by leveraging the strong approximation power of neural networks in kernel fitting.
LOW-RANK DEEP NON-STATIONARY KERNEL
Due to the intricate dependencies between events, it is challenging to choose the form of kernel function that achieves great model expressiveness while enjoying high model efficiency. In this
3

---

Published as a conference paper at ICLR 2023
section, introduce unified model with low-rank deep non-stationary kernel to capture the complex heterogeneity in events' influence over apatio-temporal space.
3.1 KERNEL WITH HISTORY AND SPATIO-TEMPORAN DISPLACEMENT
For the influence kemee function k(t', s', s), , by using the displacements and s as variables, we first re- parameterize the kernel as k(t', t-t', , s',, s s-s'), where the minus in s-s' refers to element-wiss difference between 8 and s' when ds > 1. Then we achieve a finite-rank decomposed representation based on (truncated) singular value decomposition (SVD) for kernel functions (Mollenhauer et al., 2020) (which can be understood as the kernel version matrix SVD, where the cigendecomposition is based on Mercer's Theorem (Mercer, 1909)), and that the decomposed spatial (and temporal) kernel functions can be approximated under shared basis functions (cf. Assumption A.2). The resulting approximate finite-rank representation is written as (details are in Appendix A.1)
B ΣΣ anr (aa((tt--)))uu (s'r 8¹)
k(t', - U,' 80)
1,..., L) are two sets of temporal basis functions that characterize the temporal influence of event and the decaying effect brought by elapsed time t'. Similarly, spatial basis functions {u,, : S 7 R.r = l...., R} capture the spatial influence of event at and the decayed influence after spreading over the displacement of 8 s'. The corresponding
Here [Vi,, : [0, + R,, -
weight our at different spatio-temporal ranks combines each set of basis functions into a weighted summation, leading to the final expression of influence kernel k.
To further enhance the model expressiveness, we use fully-connected neural network to represent each basis function. The history or displacement taken as the input and fed through multiple hidden layers equipped with Softplus non-linear activation function. To allow for inhibiting influence from past events (negative value of influence kernel k), we use a linear output layer for each neural network. For influence kernel with temporal rank L and spatial rank R, we need 2(L + R) independent neural networks for modeling.
The benefits of our proposed kernel framework lies in the following: (i) The kernel parameterization with displacement significantly reduces the rank needed when representing the complicated kernel encountered in practice as shown in Figure 1. (ii) The non-stationarity of original influence of histor- ical events over spatio-temporal space can be conveniently captured by in-homogeneous (1/2 (v [u,],R,, making the model applicable in modeling general STPPs. (iii) The propagating patterns of influence are characterized by {1))) (u, ) which go beyond simple parametrie forms. In particular, when the events' influence has finite range, i.e. there exist Tmax and Omax such that the influence decays to zero if It- L') > Tmax ||s s'll > Omax. we can restrict the parameterization of {p and {e,},² 1 local domain [0, Tmax] x B(0, amax) instead of [0, 7] x S, which fur- ther reduce the model complexity. Details of choosing kernel and neural network architectures are described in Appendix C.
(the class of influence kernel expressed). The proposed deep kernel representation covers a large class of non-stationary kernels generally used in STPPs. In particular, the proposed form of the kernel does not need to be positive semi-definite or even symmetrie (Reinhart, 2018). The low-rank decomposed formulation equation 4 is of SVD-type (cf. Appendix A. While each ypt (and v,) can be viewed as stationary (i.e., shift-invariant), the combination with left modes in the summation enables to model spatio-temporal non-stationarity. The technical assumptions A. and A.2 do not require more than the existence of a low-rank decomposition motivated by kernel SVD. As long as the 2(R L) many functions th,, (P,, and uy, U, sufficiently regular, they can be approximated and learned by a neural network. The universal approximation power of neural networks enables our framework express broad range of general kernel functions, and the low-rank decomposed form reduces the modeling of spatio-tempore kernel to finite many functions on time and space domains (the right modes are on truncated domains), respectively.
Remark 1
EFFICIENT COMPUTATION OF MODEL
We consider model optimization through Maximum likelihood estimation (MLE) (Reinhart, 2018). The resulting conditional intensity function could now be negative by allowing inhibiting historical
4

---

Published as a conference paper at ICLR 2023
influence. A common approach guarantee the non-negativity is adopt nonlinear positive activation function the conditional intensity (Du et al., 2016; Zhu et al., 2022). However, the integral of such nonlinear intensity spatio-temporal space is computationally expensive. To tackle this, we first introduce a log-barrier to the MIE optimization problem to guarantee the non-negativity of conditional intensity function and maintain its linearity. Then we provide computationally efficient strategy that benefits from the linearity of the conditional intensity. The extension of the approach point process data with marks is given in Appendix B.
4.1 MODEL OPTIMIZATION WITH LOG- BARRIER
We re-denote {(H) in equation by E()) in terms of model parameter 0. The constrained MLE optimization problem for model parameter estimation can be formulated as:
min ((0), X(1,)) < 0, VA [0,
Introduce a log-barrier method (Boyd et al., 2004) to ensure the non-negativity of X, and penalize the values of À on a dense enough grid Uhart x Ubar, C [0, x S. The log-barrier is defined as
p(0,b) = Mbbrtt X Ubar,s
Σ
Σ
log(A(te., s b),
(5)
ce
where Cs indicate the index of the gird, and is lower bound of conditional intensity function on the grid guarantee the feasibility of logarithm operation. The MLE optimization problem can be written
- -fa"
JeJ X(t, s)dsdt
min L(0) = -((0)- + p((0,))
Σ log A(ti,
(6)
[ebar [ehrrral
Σ
log(A(to,,,,) -b,, b),
w Mobarartt X Utaar
=1=1
where w is a weight to control the trade-off between log-likelihook and log- barrier; w and can be accordingly during the learning procedure. Details can found in Appendix A.2.
Note that previous works (Du et al., 2016; Mei and Eisner, 2017; Pan et al., 2021; Zuo et al., 2020: Zhu et al., 2022) use scaled positive transformation to guarantee non-negativity conditional intensity function. Compared with them, the log-barrier method preserves the linearity of the conditional intensity function. As shown Table such log-barrier method enables efficient model computation (See more details in Section 4.2) and enhance the model recovery power.
4.2 MODEL COMPUTATION
The log- logllikelihood computation of general STPPs (especially those with general influence function) is often difficult and requires numerical integral and thus time-consuming. Given sequence of events (z, = (44, H> of number n., the complexity of neural network evaluation is of O(n)) for the term of log- sugmmammaiion and of ((K)) (K >> n) when using numerical integration for the double integral term with K sampled points in multi multi-diimnnsionll space. In the following, we circumvent the calculation difficulty by proposing an efficient computation for L(0) with complexity (n)) of neural network evaluation through a domain discretization strategy.
Computation of log-summation.
The first log-summation term in equation 2 can be written as:
n
L
Σ log A(4i, 8:) = Σ log p
On onr (u)) -5))
(7)
CJC Pm 1=1
Note that each th only needs to be evaluated at event time [t and each u, is evaluated at all the event location [m To avoid the redundant evaluations of fo over every pair of (t,, ty), a uniform grid 24, over time horizon [0, Fmax] and evaluate yo on the grid. The value of (t - t.) can be obtained by linear interpolation with values two adjacent grid points of t, 4, By doing
up
5

---

Published as a conference paper at ICLR 2023
so, we only need to evaluate yp: for U4] times on the grids. Note that yp can be simply feed with 0 when f - > Tmax without any neural network evaluation.
Here directly evaluate (,, 83) since numerical interpolation is less accurate in location space. Note that one does not need to evaluate every pair of index (i.j). Instead, we have 1 ;= {(i, ji 0r 8j) will be computed] ((1.5) i S 1 + Tmax J {(1,5) || syll amaa ). We
use O to other pairs of (1,
Computation of integral. A benefit of our approach is that we avoid numerical integration for the conditional intensity function (needed evaluate the likelihood function), since the design of the kernel allows us to decompose the desired integral to integrating basis functions. Specifically, we have
J.' Js All, (l,s))ddl = - +
El 1m1
J. I(1 < <1)k(li, si,
(8)
" -M[ST ΣΣ ur (s))
f-= sP(())tt.
us si)ds Σ (mni((1)))
To compute the integral of (ps, we take the advantage of the pre-computed ii on the grid U. Let Fi(l) == So (pr(T)dr, Then Fi(T (4) can be computed by the linear interpolation of values of Fi at two adjacent grid points (in Uc) of T - lg. In particular, F; evaluated on U4 equals to the cumulative sum of ( divided by the grid width.
The integral of or can be estimated based on grid U in (00, 4max) ) C. Rds since it decays outside the ball. For each 8; fs (s s;) s [ (n
(s)ds, where S 5-5; 5; := {s' s' 8 - 84, 8 E S). Thus the integral is well estimated with the evaluations of e, on grid set 24 ns 8g- Note that in practice we only evaluate v, on Ul once and use subsets of the evaluations for different st. More details about grid-based computation can be found in Appendix A.3.
Computation of log-barrier. The barrier term p(0, b) is calculated in a similar way as equation 7 by replacing (t,, St, (t) with (fe, S b), i.e. we use interpolation calculate pp((e. - .ty) and evaluate u, on subset of (18,, 85)1, C. = 1,..., M
4.3 COMPUTATIONAL COMPLEXITY
Table 1: Comparison of model training time per epoch on a ID and a 3D synthetic data. Time is measured in second.
The evaluation of (u, )RR and (wh over events costs (( + L]m) complexity. The evalu-
ation of (Pt) is of O (TT(()))) complexity since it relies on the grid U4 The evaluation of {u,) ) costs no more than OR ( com- plexity. We note that L, R, Tmax, [2,,, [44,| are all constant that much less than event number n, thus the overall computation complexity will be (n). We compare the model training time per epoch for baseline equipped with a softplus activation function (NSMPP) and our model with log-barrier method (DNSK+Barr e) on 1D synthetic data set and 3D synthetic data The quantitative results in Table 1 demonstrates the efficiency improvement of our model by using log-barrier technique. More details about the computation complexity analysis can be found Appendix A.4.
1D Data I
3D Data set I
Model
Paaninctas
Training time
aPanenctuss
Trxining lime
REMPE
2076
60.040
9415
1720004
DOEE-Barvier
2107
1.299
9225
3.529
5 EXPERIMENT
We use large-scale synthetic and real data sets to demonstrate the superior performance of our model and present the results in this section. Experimental details and results can be found Appendix C. Codes will be released upon publication.
Baselines. We compare our method (DNSK+Barrier) with: (i) RMTPP (RMTPP) (Du et al., 2016); (ii) Neural Hawkes (NH) (Mei and Eisner, 2017); (iii) Transformer Hawkes process (THP) (Zuo et al., 2020); (iv) Parametric Hawkes process (PHP+exp)) with exponentially decaying spatio- temporal kernel; (v) Neual spectral marked point processes (NSMPP) (Zhu et al., 2022); (vi) DNSK without log- barrier but with nen-negative Softplus activation function (DNSK+Soff) Softplus). We note that RMTP?, NE and THP directly model conditional intensity function using neural networks while others learn the influence kernel in the framework of equation 3. In particular, NSMPP designs

---

Published as a conference paper at ICLR 2023
Table Synthetic data results. Testing log-likelihood per event (on the left side of slash, higher the better) and MRE (on the right side of slash, lower the better) are reported in each entry.
Model
ID Data set
ID Data set
ID Data
2D Data set
3D Daza
3D Tatta
-0. -0.467 0.086
mo 25.
1,3 00.212
po 10.155
po /0.117
-2 280 oesy /0.316
33543 oon /0.02
1.315
25.22 174
ce
-0.337 (n oee
2.004 aas,
(. (0.041
PEB +xx
m0:45 o0.00000
-2 2.772% 10.181
1.522k /0.223
-2.737 10.306
-2 2e0 (0.29]
-2-124 nee 10.282
NGME P
-04662 10.078
2.838 /0,
1.478 164
-
22.07 10.156
-2.07
-2 -2.37 pu.37 (0. 193
2.35 poizz (0.250
DNEX- So Elppun
0.455 10.045
om /008T
-1.300 /0.304
-
p. 70.422
-2. me (L.SS
-2.77
10.119
DNSK ernler
(oL45 (1,039
-2.5% 00.16
-1.29% cuoag -0.22 (0/31
-2.585 (caoo:
00,088
-2.49% 00021
-2.25 poceaay 00,822
True Kernel
DNSK Barrier Kerne
NSMPP Kernel
IoCa troking serueccee
Avccaged ever tecting soriinees
1D Data
50
[-t'
30
50
t- !
True Kernel
DNSK+ DNSK+Barrier Barrier Kerne
NSMPP Kernel
for beding aequstice
Avsuaure A aow leod secuintees
1D Data set
20
t-F
20
t-r
Figure 2: Kernel recovery results of 1D Data set 2 and 3. The first three columns show true kernels. kernels learned by DNsK+Barries and kernels learned by NSMPP, respectively. The last two columns show the predicted conditional intensity for testing sequence and the averaged conditional intensity for all testing sequences in each data set, which can be regarded as the event intensity by taking the expectation of the history.
the kernel based on singular value decomposition but parameterizes it without displacement. The model parameters are estimated using the training data via Adam optimization method (Kingma and Ba, 2014). Details of training can be found in Appendix A.2 and C.
5.1 SYNTHETIC DATA EXPERIMENTS
Synthetic data sets. To show the effectiveness of DNSK Barrier, we conduct all the models on three temporal data sets and three
data sets generated by following true kernels: (i) ID exponential kernel (ii) ID non-stationary kernel; (iii) ID infinite rank kernel; (iv) 2D exponential kernel; (v) 3D non-stationary inhibition kerne:; (vi) 3D non-stationary mixture kernel. Data sets are generated using thinning algorithm Daley and Vere-Jones (2008). Each data set composed of 2000 sequences. Details of kernel formulas and data generation can be found in Appendix C.
We consider two performance metrics for testing data evaluation: Mean relative error (MRE) of the predicted intensity and log-likelihoud The true and predicted X (z). X(x) can be calculated using equation 4 with true and learned kernel. The MRE for one test trajectory is defined as
Ix (x) \(z))/A² (x)dx and the averaged MRE over all test trajectories reported. The log- likelihood for observing each testing sequence can be computed according to equation 2, and the average predictive log-likelihood per event is reported. The log-likelihood shows the model's goodness-of-fe fit, and the intensity evaluation further reflects the model's 's ability to recover the underly- ing mechanism of event occurrence and predict the future.
The heat maps in Figure 2 visualize the results of non- non-stationary kernel recovery for DNSK+Barrin and NSMP? on 1D Data set 2 and 3 (The true kernel used in ID Data set is the one in Figure 1). DNSK+Barrick
recovers the true kernel more accurately than NSMPP, indicating the strong repre- sentation power of the low-rank kernel parameterization with displacements. Line charts Figure present the recovered intensities with the true ones (dark grey curves). It demonstrates that our method can accurately capture the temporal dynamics events. In particular, the average conditional intensity À over multiple testing sequences shows the model's ability to recover data non-stationarity over time. While D successfully captures the non-stationarity among data, both RMTPP and NH fail do by showing flat curve of the averaged intensity. Note that THP with
7

---

Published as a conference paper at ICLR 2023
Table 3: Results of real data sets with time and categorical marks. We compare the testing log- likelihood (Testing Ł, higher the better), event time prediction error (lower the better), and event type prediction accuracy (higher the better) of DNSK-Barriez and other baselines on each data set.
Financial
StackOverflow
Model
Testing
Time RMSE
Type Accuracy
Testing
Time RMSE
Type Accuracy
RMTPP
3.890
1.560
0.620
-2.60
9.780
0.459
NE:
33.600
1.560
0.622
-2.550
9.830
0.463
THE -0.938
1.019
0.596
-1.231
11.804
0.436
NSMPP
-3.558
1.276
0.608
-3.182
8.735
0.447
DNSK SoStluus -0,889
0.327
0.621
-2.173
6,416
0.497
DNSK-Barrier -0.709
0.153
0.630
-2.89
4.833
0,508
positional encoding recovers the data non-stationarity (as shown in two figures in the last column). However, our method still outperforms THP which suffers from limited model expressiveness when complicated propagation of event influence is involved (see two figures in the penultimate column).
Tabel summarized the quantitative results of testing log-likelihood and MRE. It shows that DNSK+Barr has superior predictive performance against baselines in characterizing the dynam- ies of data generation in spatio-temporal space. Specifically, with evidently over-pane neee-aaameterizationn for ID Data set generated by stationary exponentially decaying kernel, our model can still approximate the kernel and recover the true conditional intensity without overfitting, which shows the adaptiveness of our model. Moreover, DNSK+ Ba BarK-BBrriee rrierre ei enjoys outstanding performance gain when learning a diverse variety of complicated non-stationary kernels. The comparison between ONSK+Softplus and DNSK+Barries proves that the model with log-barrier achieves better recovery performance by maintaining the linearity of the conditional intensity. THP outperforms RMTPP in non-stationary cases but is still limited due to its pre-assumed parametric form of influence propagation. More results about kernel and intensity recovery can be found in Appendix C.
5.2
REAL DATA RESULTS
Real data sets. We provide comprehensive evaluation of our approach on several real- world data sets: we first use two popular data sets containing time-stamped events with categorical marks to demonstrate the robustness of D Barrier marked STPPs (refer to Appendix B for detailed definition and kernel modeling): (i) Financial Transactions. (Du et al., 2016). This data set contains transaction records of a stock in one day with time unit milliseconds and the action (mark) of each transaction. We partition the events into different sequences by time stamps. (ii) StackOverflow (Leskovec and Krevl, 2014): The data is collected from the website StackOverflow over two years, containing reward records for users who promote engagement in the community. Each user's reward history is treated as a sequence.
Next, we demonstrate the practical versatility of the model using the following spatio-temporal data sets: (i) Southern California earthquake data provided by Southern California Earthquake Data Center (SCEDC) contains time and location information of earthquakes in Southern California. We collect 19,414 records from 1999 to 2019 with magnitude larger than 2.5 and partition the data into multiple sequences by month with average length of 40.2. (ii) Atlanta robbery & burglary data. Atlanta Police Department (APD) provides a proprietary data source for city crime. We extract 3420 reported robberies and 14958 burglaries with time and location from 2013 to 2019. Two crime types preprocessed as separate data sets on a 10-day basis with average lengths of 13.7 and 58.7.
Finally, the model's ability to tackle high-dimensional marks is evaluated with Atlanta textual crime data. The proprietary data set provided by APD records 4644 crime incidents from 2016 to 2017 with time, location, and comprehensive text descriptions. The text information is preprocessed by TF-IDF technique, leading to 2012-dimensional mark for each event. Table summarizes the results of models dealing with categorical marks. Event time and type prediction are evaluated by Root Mean Square Error (RMSE) and accuracy, respectively. We can that DNSK+ MRrSK+Barrin Barrier outperforms the baselines in all prediction tasks by providing less time RMSE and higher type accuracy.
For real-world spatio-tempora data, we report average predictive log- ligeliool per event for the testing since MRE is not applicable. Besides, we perform online prediction for earthquake data to
8

---

Published as a conference paper at ICLR 2023
Rolibery
Banglory
ALL L..00)))))
Inssangle Mfg?
OutedFammple Atla-)
1: I. p.00
In-samplc Mfy.
Out-ol sample Mlge)
Figure 3: Kernel fitting and intensity prediction for ATL robbery & burglary data. First panel cach category shows the learned influence kernel of each crime at fixed geolocation in downtown ATL. Other panels show the predicted conditional intensity for two crimes over space at an in-sample and out-of-sample time, respectively. The dots represent reported events that day.
Table 4: Real data results with spatial and mark information. Testing log-likelihood (higher the better) and prediction mean absolute error of event time and location (lower the better) reported.
South California Earthquake
Atlanta Robbery
Atlanta Burglary
Atlanta Textual Crime
Model
Testing
Time MAE
Location MAE
Testing t
Testing é
Testing !
RMIEP
-1.825(0.00a
6.963
0.602
-2.319
0.216
NH
-18818(0.057
6,880
0.458
-2.44 -2.445go.016)
-0.308
THP
-17784
6.113
0.633
(o.22)
-0.251 (0.0013
PHP +exp
-2,018
8,132
0.287
-2.611
0.382
-10.542
NSMPP
-.11520.187
6,780
0.455
(2.
(.304(0.c0s)
-8. (8.69-(0.350)
Sofppuu
-1.807 (0.00)
1.685
0.192
-2.435 (0.00)
(o.00)
-5.8766
DNSK+Barrier
-1.751 (0.0R0)
(0.080)
1.474
0.431
-2.255,0.008)
0.519
-5.279
1RMPP, NR, and THP are not applicable when dealing with high dimesdiia data.
demonstrate the model predicting ability. The probability density function ftt, s) which represents the conditional probabílity that the next event will occur at (1, 8) given history He can be written
as [(,,s) X(l, s) exp (- Ss Sr, >(T, >(T, v)dTdd) The predicted time and location of the next event can be computed E [la+1] [HHa] = - predict the the time and location of the last event in each sequence. The mean absolute error (MAE) of the predictions is computed. The quantitative results Table 4 show that DNSK+Ba masK+Barrier provides more accurate predictions than other alternatives with higher event logelikel,.ood.
[s L
(l,,s)dsdl E [8n+ [7tc] = Ss f(I, s)dids. We
To demonstrate our model's interpretability and power to capture heterogeneous data characteristics, visualize the learned influence kernels and predicted conditional intensity for crime categories in Figure 3. The first column shows kernel evaluations at fixed geolocation in downtown Atlanta which intuitively reflect the spatial influence of crimes in that neighborhood. The influence of robbery in the downtown area is more intensive but regional, while burglary which is hard to be responded to by police in time would impact a larger neighborhood along major highways of Atlanta. We also provide the predicted conditional intensity over space for two crimes. As we can observe, DNSK BRSSKBBrriee. barrier captures the occurrence of events in regions with a higher crime rate, and crimes of the same category happening in different regions would influence their neighborhoods differently. We note that this example emphasizes the ability of the proposed method to recover data non-stationarity with different sequence lengths, and improve the limited model interpretability of other neural network-basse methods (RMTPP, NH, and THP) in practice.
For Atlanta textual crime data, borrow the idea in Zhu and Xie (2022) by encoding the highly sparse TF-IDF representation into a binary mark vector with dimension d = 50 using Restricted Boltzmann Machine (RBM) (Fischer and Igel, 2012). The average testing log-likelihoods per event for each model are reported in Table 4. The results show that BBrrriirr outperforms PHP texp in Zhu and Xie (2022) and NSMPP by achieving higher testing log likelihooood. We visualize the basis functions of learned influence kernel by
in Figure A.4 in Appendix.
6 CONCLUSION
We propose a deep non-stationary kernel for spatio-temporal point processes using a low-rank parameterization based on displacement, which enables the model to be further low-rank when learning complicated influence kernel and significantly reduces model complexity. The non-negativity of the intensity is guaranteed by lor-barrier method that maintains the linearity of the conditional intensity function. Based on that, we propose a computationally efficient strategy for model estimation. The superior performance of our model is demonstrated using synthetic and real data sets.

---

Published as a conference paper at ICLR 2023
ACKNOWLEDGEMENT
The work is partially supported by NSF DMS-2134037. Z.D. and Y.X. are partially supported by an NSF CAREER CCF-1650913, and NSF DMS-2134022 CMMI-2015787, CMMI-2112533, DMS-1995 1938106, and Foundation.
X.C. is partially supported by NSF and the Alfred P. Sloan
REFERENCES
Alex Reinhart, A review of self-exciting spatio-temporal point processes and their applications. Statistical Science, 3393):29--11,, 2018.
Yosihiko Ogata. Statistical models for earthquake occurrences and residual analysis for point processes. Journal of the American Statistical association, 83(401):9-27, 1988.
Yosihiko Ogata. Space-time point-process models for earthquake occurrences. Annals of the Institute of Statistical Mathematics, 50(2):379-402. 1998.
Frederie Paik Schoenberg, Mare Hoffmann, and Ryan Harrigan. A recursive point process model for infectious diseases. Annals of the Institute of Statistical Mathematics, 71(5):271 2019.
Zheng Dong, Shixiang Zhu, Yao Xie, Jorge Mateu, and Francisco Kodríguez-Cortés. Non- stationary spatio-temporal point process modeling for high-resolution covid-19 data. arXiv preprint arXiv:2109.099 09029, 2021.
Amanda S Hering, Cynthia L Bell, and Mare G Genton. Modeling spatio-temporal wildfire ignition point patterns. Environmental and Ecological Statistics, 16(2):225-250, 2009.
Alan G Hawkes. Spectra of some self-exciting and mutually exciting point processes. Biometrika, 58 (1):83-90, 1971.
Nan Du, Hanjun Dai, Rakshit Trivedi, Utkarsh Upadhyay, Manuel Goomez-Rorigun and Le Song. Recurrent marked temporal point processes: Embedding event history to vector. In Proceedings of the 22nd ACM SIGKDD international conference on knowledge discovery and data mining, pages 1555-1564, 2016.
Hongyuan Mei and Jason M Eisner.
The neural hawkes process: A neurally self-modulating multivariate point process. Advances in neural information processing systems, 30, 2017.
Yoshua Bengio, Patrice Simard, and Paolo Frasconi. Learning long- dependencies with gradient descent is difficult. JEEE transactions on neural networks, 5(2):157-116, 1994.
Qiang Zhang, Aldo Lipani, Omer Kirnap, and Emine Yilmaz. Self-attentive hawkes process. In International conference on machine learning, pages 11183-11193 PMLR, 2020.
Simiao Zuo, Haoming Jiang, Zichong Li, Tuo Zhao, and Hongyuan Zha. Transformer hawkes process. In International conference on machine learning, pages 1 11692-11702. PMLR, 2020.
Michael Graham, Jarno Kiviaho, and Jussi Nikkinen. Short-term and long-term dependencies of the s&p 500 index and commodity prices. Quantitative Finance, 13(4):583-50) 2013.
Shixiang Zhu, Haoyun Wang. Zheng Dong, Xiuyuan Cheng, and Yao Xie. Neural spectral marked point processes. In International Conference Leurning Representations, 2022. URL https: //oopnreviewwwr.m net/ forum ::d:0 robaaoobbbbbbbb..
Ricky TQ Chen, Brandon Amos, and Maximilian Nickel. Neural spatio-temporal point processes. arXiv preprint arXiv:2011.04583 2020.
Shixiang Zhu, Shuang Li, Zhigang Peng, and Yao Xie. Imitation learning of neural spepoolte point processes. IEEE Transactions on Knowledge and Data Engineering. 34(11):5391-544 2021.
10

---

Published as a conference paper at ICLR 2023
Che- Yi Liao,
Garcia, Kamran Paynabar, Zheng Dong, Yao Xie, and Mohammad S Jalali. Tides need stemmed: A locally operating spatio- opatio-temporal mutually exciting point process with dynamic network for improving opioid overdose death prediction. arXiv preprint 27X100221.007570. 2022.
Shuai Xiao, Junchi Yan, Xiaokang Yang, Hongyuan Zha, and Stephen Chu. Modeling the intensity function of point process via recurrent neural networks. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 31, 2017.
Takahiro Omi, Kazuyuki Aihara, et al. Fully neural network based model for general temporal point processes, Advances in neural information processing systems, 32, 2019.
Maya Okawa, Tomoharu Iwata, Yusuke Tanaka, Hiroyuki Toda, Takeshi Kurashima, and Hisashi Kashima. Dynamic hawkes processes for discovering time-evolving communities' states behind diffusion processes. In Proceedings of the 27th ACM SIGKDD Conference on Knowledge Discovery & Data Mining, pages 1276- 1286, 2021.
Jesper Moller and Rasmus Plenge Waagepetersen. Statistical inference and simulation for spatial point processes. CRC press, 2003.
Daryl Daley, David Vere-Jones, et al. An introduction to the theory of point processes: volume I: elementary theory and methods. Springer, 2003.
Mattes Mollenhauer, Ingmar Schuster, Stefan Klus, and Christof Schütte. Singular value decomposi- tion of operators reproducing kernel hilbert spaces. In Advances in Dynamics, Optimization and Computation: A volume dedicated to Michael Dellnitz on the occasion of his 60th birthday, pages 109-131. Springer, 2020.
J Mercer. Functions ofpositive and negativetypeand theircommection with the theory ofintegral equations. Philos. Trinsdictions Rogyal Soc, 209:4-415, 1909.
Stephen Boyd, Stephen P Boyd, and Lieven Vandenberghe. Convex optimization. Cambridge university press, 2004.
Zhimeng Pan, Zheng Wang, Jeff M Phillips, and Shandian Zhe. Self-adaptable point processes with nonparametrie time decays. Advances in Neural Information Processing Systems, 34, 2021.
Diederik P. Kingma and Jimmy Ba. Adam: A method for stochastic optimization, 2014. URL
https: / arw.v.viia.. wwtg//gsg///./oooooo/..../....// abs 111211111. 6980.
Daryl Daley and David Vere-Jones. An Introduction to the Theory of Point Processes. Volume II: General Theory and Structure. Springer, 2008.
Jure Leskovec and Andrej Krevl. Snap datasets: Stanford large network dataset collection, 2014.
Shixiang Zhu and Yao Xie. Spatiotemporra--extuall point processes for crime linkage detection. The Annals of Applied Staristics, 16(2):1151-1170 2022.
Asja Fischer and Christian Igel. An introduction to restricted boltzmann machines. In Iberoamerican congress pattern recognition, pages 14-36. Springer, 2012.
PA W Lewis and Gerald S Shedler. Simulation of nonhomogeneous poisson processes by thinning. Naval research logistics quarterly, 26(3):403-433, 1979.
11

---

Published as a conference paper at ICLR 2023
A ADDITIONAL METHODOLOGY DETAILS
A.1 DERIVATION OF EQUATION 4
We denote T t-tt, s-s', the variables 5' [0, T E [0, Tmax] s' and B(0,, 0max). where the sets S, [30, mmaxx) C R² Viewing the spatial and temporal variables, i.e., (f' and (s's v). as left and right mode variables, respectively, the kernel function SVD (Mollenhauer et al.,
2020; Mercer, 1909) of k gives that
M Σ (unnl(.))ha(33..))
(A.1)
We assume that the SVD can be truncated at lo K with residual of æ for some small 0, and this holds as long as the singular values 0g decay sufficiently fast. To fulfill the approximate finite-rank representation, it suffices to have the scalars on and the functions Jk and h& that the expansion approximates the kernel k even if they are not SVD of the kernel. This leads to the following assumption:
Assumption A.1. There exist coefficients oo and functions (t(t))))) hg(8' v) s.t.
k(L' S
Σ (T))(((ss))))) O(e).
(A.2)
To proceed, one can apply kernel SVD again to lo and hy respectively, and obtain left and right singular functions that potentially differ for different k. Here, we impose that across k k= 1,... KK the singular functions of gk are the same (as shown below, being approximately same suffices) set of basis functions, that is,
ga (1,,T) Σ (()))
As we will truncate I to be up to finite rank again (up to an O(E) residual) we require the (approx- imately) shared singular modes only up to L. Similarly above, technically it suffices to have finite-rank expansion to achieve the ((e) error without requiring them to be SVD, which leads to the following assumption where we assume the same condition for hg:
Assumption A.2. For the 9h and hg in equation A.2, up to an O(e) error,
(i) The K temporal kernel functions (t'(t', r) can be approximated under a same set left and right basis functions, ie., there exist coefficients BEi, and functions (tr('), (P(()) for -,... L, s.t.
gu (1' T) Σ
+ == K.
(A.3)
(ū) The K spatial kernel functions he (h',s' can be approximated under a same set of left and right basis functions, i.e., there exist coefficients "ek, and functions u,,(s',, v, (v) for 1.... R. s.A.
he(s', v) = Σ Tkr (f)c, r 1...., R.
(A.4)
Inserting equation A.3 and equation A.4 into equation A.2 gives the rank-truncated representation of the kernel function. Since K, L., R are fixed numbers, assuming boundedness of all the coefficients and functions, have the representation with the final residual O namely,
A(',T,8,,)) - (h (1), + 0(s).
Defining Ek Ox.Ba oir leads equation 4.
A.2 ALGORITHMS
12

---

Published as a conference paper at ICLR 2023
Algorithm 1 Model parameter estimation
Input: Training set X, batch size M, epoch number E, learning rate 7, constant a 1 update in equation 6.
Initialization: model parameter 0g first epoch 0, 8g-
while < E do
for each batch with size M
1. For ID temporal point process, compute ((0), {A(ber fe.= Mtn.l For spatio-temporal
point process, compute f(0), (X(ter Sc,) o. Maan e lics =I...., Mor.aal
2. Set min (t(ee, Jf.,=1...., -e (or m (( Se, 1.. Rbun lic =1, -,...,,6ar,. c), where is small value to guarantee logarithm feasibility.
3. Compute [(0) (()) au ((, b).
4. Update 0e+1 + 0g 8C
U0.
5, - + 1, w + w end for end while
Algorithm 2 Synthetic data generation
Input: Model X(-), I, S, Upper bound of conditional intensity A. Initialization: H.r t 0, 0
while < T do
Unif(0, 1). 2. t < t un/
1. Sample
3. Sample ~ Unif(S), D Uni
4. À \(tt)
if DĀ S then
n n n 1; tn = sn
Hr + Hr U ((n,, 8n)]
end if
end while
if in >= T' then
return Mr {tnn, sn)}
else
return HT end
A.3 GRID-BASED MODEL COMPUTATION
In this section, we elaborate the details of the grid-based efficient model computation.
In Figure A. we visualize the procedure of computing the integrals of TT (((())t and Ss v.(s & ds in equation respectively. Panel (a) illustrates the calculation of fa TT yo (()dt. As explained in Section 4.2, the evaluations of yor only happens the grid Ui over [0, Tmax] (since (t) = 0 when t Tmax). The value of F(I) = - (P(()) on the grid can be obtained through numerical integration. Then given d; the value of F(T (4) - (pr(t)dt is calculated using linear interpolation of F on two adjacent grid points of t,. Panel (b) shows the computation of Js 0 Sz ds. Given st Is (s s)dss JBtro..m....ee...ee_.ppp - 2y or-ss)ds since v, (s) = when >> amax. Then B(0, amax) is discretized into the grid U. and Is ( sds be calculated based the value of u, on the grid points in U4, n 5; (the deep red dots in Figure A. (bb) using numerical integration.
To evaluate the sensitivity of our model to the chosen grids, we compare the performance of DNSK+Barrier on 3D Data set using grids with different resolutions. The quantitative re- sults of testing log-likelin likelihood and intensity prediction enror reported in Table A.1. We use [24] = 50, [24, = 1500 for the experiments in the main paper. As we can see, the model shows similar performances when higher grid resolution is used and works slightly less accurately but still
13

---

Published as a conference paper at ICLR 2023
9:(c)
Ur POT
B(0 Amaz))
U
0
Tmax
(i(t) 5 (r)dr
0
S -
(b)
Figure A.1: (a) Computation of 1o' (t) computation based on grid 24 Red dots represent grid points. Line segments connecting two light or dark grey dots represent the linear interpolation of yo and Fj. Here 4; is the time of the historical event which is fixed. (b) Computation of [s ur (s s; s;)ds based on grid 2, The background heatmap represents the evaluation of v, over S. Here the fixed 8; is the location of the historical event. The integral is calculated based the values of v, on grid points with dark red color.
Table A.I: Comparison of DNSK + Barrier performance on 3D Data set 2 with different grid resolutions. Testing log-likelihood per event and intensity MRE are reported. The highlighted ones are the results in the main paper.
Spatial resolution: [4,]]
Temporal resolution: [44]]
1000
1500
3000
30
(2.27200.006) /o..ooor//..1022 /0.102
-2.252 002) 00.08
-2.250 /0.088
50
22.257 (0.002)
/aaa//00000 /0.095
-2.251((.o00)) /0.082
-2.249(0.001) 00077
100
-2.25550.001// -2.252(0.001)) 00.0000001//.0811 0.081
-2.250(0.001)00.077
better than other baselines with less number of grid points. It reveals that our choice of grid resolution accurate enough to capture the complex dynamics event occurrences for this non-stationary data, and the model performance is robust different grid resolutions.
In practice, the grids can flexibly chosen to reach the balance of model accuracy and computational efficiency. For instance, the number of uniformly distributed grid points along one dimension be chosen around Omro), where Tig is the average number of events in one observed sequence. Note that [Z4,] or [U,]] would be far less than the total number of observed events because we use thousands of sequences (2000 in our synthetic experiments) for model learning. And the grid size can be even smaller when comes non-Lebesgue-mete measured space.
A.4 DETAILS OF COMPUTATIONAL COMPLEXITY
We provide the detailed analysis of the (nn) computation complexity of L(0) in Section 4.3 as following:
Computation of log-summation. The evaluation of I/ and {attt over events costs O(R + L)n) complexity. The evaluation of { is of O(L(4))) complexity since relies on the 2r, With the assumption that the conditional intensity bounded by a constant C in finite time horizon (Lewis and Shedler, 1979; Daley et al., 2003; Zhu et al., 2022), for each fixed j. the cardinality of set {(i,j) ty t4 S t, Tmax] less than CTmax. which leads (ROOmaax)) complexity of {u, J. evaluation.
Computation of integral. The integration of (1t) only relies on numerical operations of (pr) on grids Ui without extra evaluations of neural networks. The integration { depends on the evaluation on grid U, of O (RR,,,) complexity.
Computation barrier. { H on grid Uhar, is estimated by numerical interpolation of previously computed (p on grid 24, Additional neural network evaluations of {u, } cost no more than O D(RCTmax)) complexity.
14

---

Published as a conference paper at ICLR 2023
B DEEP NON-STATIONARY KERNEL FOR MARKED STPPS
In marked STPPs (Reinhart, 2018), each observed event is associated with additional information describing event attribute, denoted as m M RRdM, Let H {(tt, 8i, m)) denote the event
sequence. Given the observed history H, = {(4,, m,) E H/4 < tj, the conditional intensity function of marked STPPs similarly defined as:
E [N[ At) B(s, As) B(m, Am)) [Ha]
lim s, Ar10,As40,Am.l
[ss, As) (m, Am)]Ar
where (m,, Am) is ball centered at m E RR with radius Am. The log-likelihood observing
H [0, T x S x M is given by
E(H) = Σ log X (f,, -
T So S A(t, m)dmdsdt.
B. KERNEL INCORPORATING MARKS
One of the salient features of our spatio-temporal kernel framework is that can be conveniently adopted in modeling marked STPPs with additional sets of mark basis functions {gg. h4 1 We modify the influence kernel function k accordingly as following:
R
l: 1 1-11 1',,8 m' m)
Virgh1(()) 1' (o) (m)g((mm))
Here m', m E M C RR and {gg. MM networks model the influence of historical mark m'' and current mark m, respectively. Since the mark space M is always categorical and the difference between m'' and m is little practical meaning, use By and hy model m' and separately instead of modeling m'.
represented by independent neural
B.2 LOG- LOG-BARRIER AND MODEL COMPUTATION
The conditional intensity for marked spatio-temporal point processes at (t. 8, m) can be written as:
A(t, m) = + Σ
airy
Σ
ss)) (m4)hy
(E)..,
We need to guarantee the non-negativity of over the space of [0, T x x M. When the total number of unique categorical mark in . is small, the log-barrier be conveniently computed as the summation of À on grids Ubar, x Ubar. x M. In the following we focus on the case that M is high-dimensional with (nn) number of unique marks.
For model simplicity we use non-negative 9g and h2 in this case (which can be done by adding a non-negative activation function the linear output layer in neural networks). We re-write A(t, s, m) and denote as following:
A(1,8, = 1z - Σ
Σ 1r
Oirg
Σ (L..L.m.))E,
(r(1))p(1 - (.),, (m))g(()) hy(m).
(t,0.)
Note that the function in the brackets only with regard to t, We denote it (t, (since is the rth rank of mark). Since hy (m) 2 0, the non-negativity X can be guaranteed by the non-negatyvity of Fy (t, Thus we apply log-barrier method on (t, s). The log-barrier term becomes:
Llve,
1
p(0,b) Qi x Utnrree
log (F,, (f (,,,,)))) b),
q=1
15

---

Published as a conference paper at ICLR 2023
Since model is low-rank, the value of Q will not be large.
For the model computation, the additional evaluations for (g events is of ((Q)) complexity and the evaluations for { only depends on the unique number of marks which at most of O(n). The log- rg--barrer method does not introduce extra evaluation in mark space. Thus the overall computation complexity for DNSK in marked STPPs is still O(n).
C ADDITIONAL EXPERIMENTAL RESULTS
In this section we provide details of data sets and experimental setup, together with additional experimental results.
Synthetic data sets. To show the robustness of our model, we generate three temporal data sets and three spatio-tempon data sets using the following kernels:
ID Data set with exponential kernel: ()) = 1..8-(1-1))
(ii) ID Data set 2 with non-stationary kernel: k(t',t) 0.3(0.5 + 0.5 cos((.2tt)) (
(iii) ID Data set 3 with infinite rank kernel:
> 2- 0.3 + cos(2 + _ ( + 1) k(t',) 0.3 j=1
(iv) 2D Data 1 with exponential kernel: k(l', 1,8 -.55-
(v) 3D Data set 1 with non-stationary inhibition kernel:
cos (10 - 11)
1 2
k(t', t, s', = 0.3(1 0.01t]e 0.01t)e-2((1)) 2((-())
2TOz 2 e00( 1s-8'11
where o. = 0.5, v. = 0.15.
(vi) 3D Data set 2 with non-stationary mixture kernel:
2
M(L' 8,,
Ortu (x)))) 19)
r=1 1=1
1e-e'12
where un (5') = a- +1) U2 (s') = 1-bx(86+1)) (1-888)) 2xo7
20
8') 1
2n1
(1') - 1-all (t') 1 Aill,1p 1')
2x0
e -(t-tt)) -(ttt'') ,2((t
t) (t (t-t-1) I[t-t' 3), and a = 0.3. b, 0.4, a, 0.02, be
0.02, 01 = 0.2, 02 = 0.3, 2, (011, 012 021 a22) = (0.6, 0.15, 0.225, 0.525).
Note that kernel (iii) is the one we illustrated in Figure which is of infinite rank according to the formulas. In Figure the value matrix of k=(1, 1) and k((',, f') are the kernel evaluations same 300 x 300 uniform grid. As we can see, the rank of the value matrix of the same kernel is reduced from 298 after changing the displacement-basee kernel parameterization.
Details of Experimental setup. For RMTPP and NH we test embedding size of [32, 64, 128) and choose 64 for experiments. For THP we take the default experiment setting recommended by Zuo et al. (2020). For NSMPP we use the same model setting in Zhu et al. (2022) with rank 5. Each experiment is implemented by the following procedure: Given the data set, we split 90% of the sequences training and 10% testing set.
We use independent fully-comne connecco neural networks with two-hidden layers for each basis function. Each layer contains 64 hidden nodes. The temporal rank of DNSK+Barries is set to be 1 for synthetic data (i), (ii). (iv), (v), for (vi), and 3 for (iii). The spatial rank for synthetic data (iv), (v) and 2 for (vi). The temporal and spatial rank for real data are both set be 2 through cross
16

---

Published as a conference paper at ICLR 2023
True kC. 1.1)
True k(S, 5.1, -)
True A6, +)
True Average À
o20
0.39
Mark
920
-015
Predict 1,1,111.111 1.1.1)
Predict k(5, 5.1, )
Predict AC,
Predict Average
eas
005
0.3
039
Figure A.2: Kernel recovery results of 2D exponential kernel. The first two columns show true kernel and kernel learned by DNSX+Barrier. The last two columns shows the true and predicted conditional intensity functions of a test sequence. The line charts visualize the conditional intensity average over the ID mark space any given time for the case of presentation. The red dots indicate the time of observed events.
validation. For each real data set, the Tmax chosen be around 7/4 and Smax is 1 for each data set since the location space is normalized before training. The hyper-parament parameter of DNSK+So St ftppus For RMTPP, NH, and THP the batch size is 32 and the learning rate is 10-3, For others, the batch size is 64 and the learning rate is 10-1 The quantitative results are collected by running each experiment for independent times. All experiments are implemented on Google Colaboratory (Pro version) with 25GB RAM and Tesla T4 GPU.
are the same as DNSK+Barries
C.1 SYNTHETIC RESULTS WITH 2D & 3D KERNEL.
In this section we present additional experiment results for the synthetic data sets with 2D exponential and 3D non-stationary mixture kernel. Our proposed model successfully recovers the kernel and event conditional intensity in both case. Note that the recovery of 3D mixture kernel demonstrates the capability of our model handle complex event dependency with mixture patterns by conveniently setting time and mark rank to be more than 1.
C.2 ATLANTA TEXTUAL CRIME DATA WITH HIGH-DIMENSIONSS MARKS
Figure A.4 visualizes the fitting and prediction results DNSK+Barries Our model presents decaying pattern in temporal effect and captures different patterns of spatial influence for incidents in the northeast. Besides, the in-sample and out-of-sample intensity predictions demonstrate the ability of DNSK to characterize the event occurrences by showing different conditional intensities.
17

---

Published as a conference paper at ICLR 2023
X35
Gea
330
E15
611
300
Predict k15,5.3 [-5. 0.53,
Predict k.5, [0.,[O,000 [),
Predict A13.8, -)
Predict
330
3.0
True kts, 5.9, [-0 D.S. 0.55
True kts, 7.(-0.5 0.51. -)
True A13.8,
True A(10.0,
so
605
92
mg
mg
Predict k.5.5.9.-0033 0,51,
Predict kis, (. 0.5),
Predict A(3.8.
Predict A110 0,
8s
100
Figure A.3: Kernel recovery results of 3D non-stationary mixture kernel. The first two columns show snapshots the true kernel and kernel learned by INSK+Barrien The last two columns shows snapshots of the true and predicted conditional intensity functions of test sequence.
23
Oxt-otsmmpp
13
Figure A.4: Model fitting and prediction for high-dimensional real data. First column shows the learned temporal functions. Four panels in the middle shows the learned spatial functions, Deeper color depth indicates higher function value. Last four panels show the predicted conditional intensity over space two in-sample times and two out-of-sample times, respectively. The dots represent event occurrences that day.
18


Reviews:
[{"review_id": "iwJ0iiN7pM", "reviewer": "Reviewer_GLBW", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "The authors of this submission proposed a deep non-stationary kernel (DNSK) modeling Spatio-Temporal Point Processes (STPP) for potentially non-stationary events in continuous time and space. \n\nThe authors focused on Hawkes process assuming that the influences from the past events are linearly additive and in turn modeled the conditional intensity function as $\\lambda(t, s) = \\mu + \\sum_{History} k(t, t', s, s')$, where $k(t, t', s, s')$ is named as the influence kernel function, which captures the spatio-temporal dependency. Specifically, DNSK assumes a factorized kernel along time and space directions, as well as potential event \"marks\" when they are available to handle potential computational challenges. Furthermore, by considering parametrization with \"displacement\", the authors claimed that DNSK achieves low-rank kernels for more efficiency. The log-barrier method is introduced to preserve non-negativity of the conditional intensity function, also maintaining model interpretability and computational efficiency. \n\nExperiments were performed comparing some existing methods: RMTPP, Neural Hawkes (NH), Transformer Hawkes process (THP), Parametric Hawkes process (PHP), Neural spectral marked point processes (NSMPP) with the proposed DNSK without log-barrier but with a non-negative softplus activation function (DNSK+Softplus). The results showed that the proposed method can indeed capture the non-stationarity and achieve good model recovery and event prediction performances.  ", "summary_of_the_review": "The authors proposed DNSK to efficiently model potentially non-stationary events in continuous time and space. Empirical performance comparison results have shown the efficacy of the proposed DNSK methods. Better justification of the adopted kernels, in particular linear combination of MLP-based kernel bases, may be needed. More comprehensive empirical experimental results may need to be provided to confirm the benefits of DNSK compared to the previous similar efforts. \n\n====post-rebuttal====\n\nI truly appreciate the extensive efforts from the authors to address my questions.", "strength_and_weaknesses": "Efficient modeling non-stationary spatio-temporal processes is challenging. The authors followed recent efforts, in particular Zhu et al. 2022, Zou et al. 2020, etc. to develop a \"displacement\"-parametrized neural network based kernels in a Hawkes process model to address the potential computational challenges. The empirical results demonstrated the efficacy of their proposed DNSK. \n\nClearer discussion besides the illustration in Figure 1 may be needed for clearer motivation and insights. For example, it may be necessary to discuss the pros and cons of factorized kernel assumptions and displacement-based parametrization. When introducing Figure 1, the authors may want to briefly mention the computational advantages and difference of the actual event prediction performance before and after using \"displacement\" on the kernels with the same rank. Although the figure did show that the displacement-based kernel matrix with rank 7 also has three peaks as the original kernel with rank 298, it may be necessary to provide numerical kernel performance measure. The authors may also want to discuss more about the reason for building up the kernels using linear combinations of MLP-based basis functions. As MLPs are assumed to have good approximation capability, is it still necessary to have the linear combination of these MLPs to make the kernel more complicated? \n\nThe authors may want to make there empirical performance comparison results more consistent to have all the results reported from all the selected baselines. For example, it may be interesting to explore on the effectiveness of log-likelihood because the potential mismatch of log-likelihood and MSE regarding the THP performance on 1D Data set 2 and 1D Data set 3. For the real-world data, instead of showing the performances of some baselines, the authors may want to present the results from all the applicable baselines. \n\n\n\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "_vfrYi5JBhY", "reviewer": "Reviewer_ZPhD", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This manuscript proposes a method to model non-stationary spatio-temporal events in the framework of hawkes process. The specific method is to construct a more refined kernel function.", "summary_of_the_review": "It is an interesting work on modeling non-stationary spatial-temporal data and the proposed method is solid and demonstrated to be effective.", "strength_and_weaknesses": "Strength:\n1. The idea of using deep non-stationary kernel in the point process to model spatial-temporal data is interesting and somewhat novel.\n2. The method is clear and computational complexity is made.\n3. The experiments are sufficient to demonstrate the effectiveness of the proposed method.\n\n\n\nWeaknesses\n\n1. The kernel k is designed with heuristics, and can it guarantee that k has good properties, such as to be positive definite?\n2. It is said that \" (ii) The non-stationarity of events’ influence over spatial-temporal space can be conveniently captured by non-constant psi_l and u_r\". I wonder if this is like a kind of positional coding of time index? If so, transformer also has the ability to model non-stationary temporal data, and the experimental results show that THP performs comparable as the proposed methods. The authors did not make a detailed analysis of THP, nor did they compare the differences and advantages of THP with the proposed methods.", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "XbvH18m5lE", "reviewer": "Reviewer_D4xW", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "the paper proposes a more general form of kernel function that is typically used in temporal/spatio-temporal point process, by considering an absolute time-dependent component in addition to the relative spatial-time inputs. In addition, authors made another contribution by proposing a more efficient approach  to ensure the positivity of intensity in the form log-barrier to the optimization problem. Empirical results show both accuracy and efficiency gain over baselines. ", "summary_of_the_review": "novel formulation on the new kernel with promising results. ", "strength_and_weaknesses": "Strength:\na new kernel method is proposed by considering absolute time\ncomplexity analysis is provided\nempirical evaluation performance is strong\n\nWeakness:\nadditional empirical study would provide more understanding of the method. \n\nOther comments:\n- A baseline RNN model with the s and t as an input to RNN, which models time-dependent change in function and similar to your base kernel, should be considered. It would highlight the necessity of using kernel. \n- Background does the spatial aspect of the STPP indicate that the change in location always be positive in location? what's an example for such application?\n- does the grid resolution have impact on the learning of the kernel?\n- would a large t value, corresponding to long term dependency, pose any numerical issue ? \n- Eq4: it would be helpful if there is some formal statement on assumptions, the generalization result or its expressivity related to Eq 3.  \n- Table 3: how does DNSK+softplus compare with transformer results? Would using barrier improve THP as well?\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "3e3cmXQ0zXG", "reviewer": "Reviewer_ggHf", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper proposes a deep non-stationary kernel for spatio-temporal point processes using a different parameterization scheme, which reduces the model complexity. The non-negativity of the solution is guaranteed by a log-barrier method which maintains the linearity of the conditional intensity function. In addition, a computationally efficient strategy for model estimation is introduced. Both the synthetic and real data sets are used to validate the superiority of the proposed model. ", "summary_of_the_review": "The major contribution of this paper is to largely reduce the complexity of the existing neural kernel-based approach towards the spatial-temporal point process. The proposed technique is novel and technical sounds. Though no much improvement on performance and the paper is not solving a new problem, the complexity is largely reduced without sacrificing the performance. Thus, I admit the  contribution and novelty of this paper are above-average. There are some minor issue towards statements and descriptions which have been mentioned above. ", "strength_and_weaknesses": "Strength:\n1. The proposed methods significantly reduces the model complexity without sacrificing the performances. The motivation is clear and the contribution is satisfied.\n2. Comprehensive and well designed experiments are conducted.\n3. This paper is easy to follow with good presentation.\n\nWeaknesses:\n1. Some preliminary knowledge needs to be included to make the paper self-contained. For example, it would be better to briefly describe the \"mask point process\".\n2. Some notations need more clarifications. For example, what \"B\" refers to in Equation 2. And what the dimension of \"s\". In Equation 4, there is operation of \"s-s'\", while \"s\" refers to the location which should be a 2-D or 3-D value. How the \"-\" operation is conducted in these values?\n3. More discussions on the motivation of specific techniques can be provided. For example, compared to directly using the neural network, what is the advantage of the kernel in Equation 4. \n4. No clear problem formulation. For better understanding the problem and methodology, it will be great to have a formal problem formulation. \n5. The meaning of \"testing l\" in Table 4 is not explained. \n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "1: The contributions are neither significant nor novel."}, "extra_text": {}}]

Author rebuttal:

Decision: Accept: poster