Wasserstein Learning of Deep Generative Point Process Models

Shuai Xiao, Mehrdad Farajtabar, Xiaojing Ye, Junchi Yan, Le Song, Hongyuan Zha

Introduction

Event sequences are ubiquitous in areas such as e-commerce, social networks, and health informatics. For example, events in e-commerce are the times a customer purchases a product from an online vendor such as Amazon. In social networks, event sequences are the times a user signs on or generates posts, clicks, and likes. In health informatics, events can be the times when a patient exhibits symptoms or receives treatments. Bidding and asking orders also comprise events in the stock market. In all of these applications, understanding and predicting user behaviors exhibited by the event dynamics are of great practical, economic, and societal interest.

Temporal point processes is an effective mathematical tool for modeling events data. It has been applied to sequences arising from social networks , electronic health records , e-commerce , and finance . A temporal point process is a random process whose realization consists of a list of discrete events localized in (continuous) time. The point process representation of sequence data is fundamentally different from the discrete time representation typically used in time series analysis. It directly models the time period between events as random variables, and allows temporal events to be modeled accurately, without requiring the choice of a time window to aggregate events, which may cause discretization errors. Moreover, it has a remarkably extensive theoretical foundation .

However, conventional point process models often make strong unrealistic assumptions about the generative processes of the event sequences. In fact, a point process is characterized by its conditional intensity function – a stochastic model for the time of the next event given all the times of previous events. The functional form of the intensity is often designed to capture the phenomena of interests. Some examples are homogeneous and non-homogeneous Poisson processes , self-exciting point processes , self-correcting point process models , and survival processes . Unfortunately, they make various parametric assumptions about the latent dynamics governing the generation of the observed point patterns. As a consequence, model misspecification can cause significantly degraded performance using point process models, which is also shown by our experimental results later.

To address the aforementioned problem, the authors in propose to learn a general representation of the underlying dynamics from the event history without assuming a fixed parametric form in advance. The intensity function of the temporal point process is viewed as a nonlinear function of the history of the process and is parameterized using a recurrent neural network. Apparently this line of work still relies on explicit modeling of the intensity function. However, in many tasks such as data generation or event prediction, knowledge of the whole intensity function is unnecessary. On the other hand, sampling sequences from intensity-based models is usually performed via a thinning algorithm , which is computationally expensive; many sample events might be rejected because of the rejection step, especially when the intensity exhibits high variation. More importantly, most of the methods based on intensity function are trained by maximizing log likelihood or a lower bound on it. They are asymptotically equivalent to minimizing the Kullback-Leibler (KL) divergence between the data and model distributions, which suffers serious issues such as mode dropping . Alternatively, Generative Adversarial Networks (GAN) have proven to be a promising alternative to traditional maximum likelihood approaches .

In this paper, for the first time, we bypass the intensity-based modeling and likelihood-based estimation of temporal point processes and propose a neural network-based model with a generative adversarial learning scheme for point processes. In GANs, two models are used to solve a minimax game: a generator which samples synthetic data from the model, and a discriminator which classifies the data as real or synthetic. Theoretically speaking, these models are capable of modeling an arbitrarily complex probability distribution, including distributions over discrete events. They achieve state-of-the-art results on a variety of generative tasks such as image generation, image super-resolution, 3D object generation, and video prediction .

The original GAN in minimizes the Jensen-Shannon (JS) and is regarded as highly unstable and prone to miss modes. Recently, Wasserstein GAN (WGAN) is proposed to use the Earth Moving distance (EM) as an objective for training GANs. Furthermore it is shown that the EM objective has many advantages as the loss function correlates with the quality of the generated samples and reduces mode dropping . Moreover, it leverages the geometry of the space of event sequences in terms of their distance, which is not the case for an MLE-based approach. In this paper we extend the notion of WGAN for temporal point processes and adopt a Recurrent Neural Network (RNN) for training. Importantly, we are able to demonstrate that Wasserstein distance training of RNN point process models outperforms the same architecture trained using MLE.

In a nutshell, the contributions of the paper are: i) We propose the first intensity-free generative model for point processes and introduce the first (to our best knowledge) likelihood-free corresponding learning methods; ii)-3mm We extend WGAN for point processes with Recurrent Neural Network architecture for sequence generation learning; iii) In contrast to the usual subjective measures of evaluating GANs we use a statistical and a quantitative measure to compare the performance of the model to the conventional ones. iv) Extensive experiments involving various types of point processes on both synthetic and real datasets show the promising performance of our approach.

Proposed Framework

In this section, we define Point Processes in a way that is suitable to be combined with the WGANs.

For inhomogeneous Poisson process, M(B)=∫Bλ(x)dxM(B)=\int_{B}\lambda(x)dx, where the intensity function λ(x)\lambda(x) yields a positive measurable function on SS. Intuitively speaking, λ(x)dx\lambda(x)dx is the expected number of events in the infinitesimal dxdx. For the most common type of point process, a Homogeneous Poisson process, λ(x)=λ\lambda(x)=\lambda and M(B)=λ∣B∣M(B)=\lambda|B|, where ∣⋅∣|\cdot| is the Lebesgue measure on (S,B)(S,\mathcal{B}). More generally, in Cox point processes, λ(x)\lambda(x) can be a random density possibly depending on the history of the process. For any point process, given λ(⋅)\lambda(\cdot), N(B)∼Poisson(∫Bλ(x)dx)N(B)\sim\text{Poisson}(\int_{B}\lambda(x)dx). In addition, if B1,…,Bk∈BB_{1},\ldots,B_{k}\in\mathcal{B} are disjoint, then N(B1),…,N(Bk)N(B_{1}),\ldots,N(B_{k}) are independent conditioning on λ(⋅)\lambda(\cdot).

For the ease of exposition, we will present the framework for the case where the events are happening in the real half-line of time. But the framework is easily extensible to the general space.

2 Temporal Point Processes

A particularly interesting case of point processes is given when SS is the time interval [0,T)[0,T), which we will call a temporal point process. Here, a realization is simply a set of time points: ξ=∑i=1nδti\xi=\sum_{i=1}^{n}\delta_{t_{i}}. With a slight notation abuse we will write ξ={t1,…,tn}\xi=\{t_{1},\ldots,t_{n}\} where each tit_{i} is a random time before TT. Using a conditional intensity (rate) function is the usual way to characterize point processes.

For Inhomogeneous Poisson process (IP), the intensity λ(t)\lambda(t) is a fixed non-negative function supported in [0,T)[0,T). For example, it can be a multi-modal function comprised of kk Gaussian kernels: λ(t)=∑i=1kαi(2πσi2)−1/2exp⁡(−(t−ci)2/σi2),\lambda(t)=\sum_{i=1}^{k}\alpha_{i}(2\pi\sigma_{i}^{2})^{-1/2}\exp\left(-(t-c_{i})^{2}/\sigma_{i}^{2}\right), for t∈[0,T)t\in[0,T), where cic_{i} and σi\sigma_{i} are fixed center and standard deviation, respectively, and αi\alpha_{i} is the weight (or importance) for kernel ii.

A self-exciting (Hawkes) process (SE) is a cox process where the intensity is determined by previous (random) events in a special parametric form: λ(t)=μ+β∑ti<tg(t−ti),\lambda(t)=\mu+\beta\sum_{t_{i}<t}g(t-t_{i}), where gg is a nonnegative kernel function, e.g., g(t)=exp⁡(−ωt)g(t)=\exp(-\omega t) for some ω>0\omega>0. This process has an implication that the occurrence of an event will increase the probability of near future events and its influence will (usually) decrease over time, as captured by (the usually) decaying fixed kernel gg. μ\mu is the exogenous rate of firing events and α\alpha is the coefficient for the endogenous rate.

In contrast, in self-correcting processes (SC), an event will decrease the probability of an event: λ(t)=exp⁡(ηt−∑ti<tγ).\lambda(t)=\exp(\eta t-\sum_{t_{i}<t}\gamma). The exp⁡\exp ensures that the intensity is positive, while η\eta and γ\gamma are exogenous and endogenous rates.

We can utilize more flexible ways to model the intensity, e.g., by a Recurrent Neural Network (RNN): λ(t)=gw(t,hti),\lambda(t)=g_{w}(t,h_{t_{i}}), where htih_{t_{i}} is the feedback loop capturing the influence of previous events (last updated at the latest event) and is updated by hti=hv(ti,hti−1)h_{t_{i}}=h_{v}(t_{i},h_{t_{i-1}}). Here ww, vv are network weights.

3 Wasserstein-Distance for Temporal Point Processes

where ss is a fixed limiting point in border of the compact space SS and the minimum is over all permutations of 1…m1\ldots m. The second term penalizes unmatched points in a very special way which will be clarified later. Appendix B proves that it is indeed a valid distance measure.

Interestingly, in the case of temporal point process in [0,T)[0,T) the distance between ξ={t1,…,tn}\xi=\left\{t_{1},\ldots,t_{n}\right\} and ρ={τ1,…,τm}\rho=\left\{\tau_{1},\ldots,\tau_{m}\right\} is reduced to

where the time points are ordered increasingly, s=Ts=T is chosen as the anchor point, and ∣⋅∣|\cdot| is the Lebesgue measure in the real line. A proof is given in Appendix C. This choice of distance is significant in two senses. First, it is computationally efficient and no excessive computation is involved. Secondly, in terms of point processes, it is interpreted as the volume by which the two counting measures differ. Figure 1-b demonstrates this intuition and justifies our choice of metric in Ξ\Xi and Appendix D contains the proof.

Equation (1) is computationally highly intractable and its dual form is usually utilized :

However, solving the dual form is still highly nontrivial. Enumerating all Lipschitz functions over point process realizations is impossible. Instead, we choose a parametric family of functions to approximate the search space fwf_{w} and consider solving the problem

where w∈Ww\in\mathcal{W} is the parameter. The more flexible fwf_{w}, the more accurate will be the approximation.

It is notable that W-distance leverages the geometry of the space of event sequences in terms of their distance, which is not the case for MLE-based approach. It in turn requires functions of event sequences f(x1,x2,...)f(x_{1},x_{2},...), rather than functions of the time stamps f(xi)f(x_{i}). Furthermore, Stein’s method to approximate Poisson processes is also relevant as they are defining distances between a Poisson process and an arbitrary point process.

4 WGAN for Temporal Point Processes

5 Ingredients of WGANTPP

The generator transforms a given sequence to another sequence. Similar to we use Recurrent Neural Networks (RNN) to model the generator. If the input and output sequences are ζ={z1,…,zn}\zeta=\left\{z_{1},\ldots,z_{n}\right\} and ρ={t1,…,tn}\rho=\left\{t_{1},\ldots,t_{n}\right\} then the generator gθ(ζ)=ρg_{\theta}(\zeta)=\rho works according to

Here hih_{i} is the kk-dimensional history embedding vector and ϕgh\phi^{h}_{g} and ϕgx\phi^{x}_{g} are the activation functions. The parameter set of the generator is θ={(Agh)k×1,(Bgh)k×k,(bgh)k×1,(Bgx)1×k,(bgx)1×1}.\theta=\left\{\left(A^{h}_{g}\right)_{k\times 1},\left(B^{h}_{g}\right)_{k\times k},\left(b^{h}_{g}\right)_{k\times 1},\left(B^{x}_{g}\right)_{1\times k},\left(b^{x}_{g}\right)_{1\times 1}\right\}. Similarly, we define the discriminator function who assigns a scalar value fw(ρ)=∑i=1naif_{w}(\rho)=\sum_{i=1}^{n}a_{i} to the sequence ρ={t1,…,tn}\rho=\left\{t_{1},\ldots,t_{n}\right\} according to

where the parameter set is comprised of w={(Adh)k×1,(Bdh)k×k,(bdh)k×1,(Bda)1×k,(bda)1×1}.w=\left\{\left(A^{h}_{d}\right)_{k\times 1},\left(B^{h}_{d}\right)_{k\times k},\left(b^{h}_{d}\right)_{k\times 1},\left(B^{a}_{d}\right)_{1\times k},\left(b^{a}_{d}\right)_{1\times 1}\right\}. Note that both generator and discriminator RNNs are causal networks. Each event is only influenced by the previous events. To enforce the Lipschitz constraints the original WGAN paper adopts weight clipping. However, our initial experiments shows an inferior performance by using weight clipping. This is also reported by the same authors in their follow-up paper to the original work. The poor performance of weight clipping for enforcing 1-Lipschitz can be seen theoretically as well: just consider a simple neural network with one input, one neuron, and one output: f(x)=σ(wx+b)f(x)=\sigma(wx+b) and the weight clipping w<cw<c. Then,

It is clear that when 1/∣σ′(wx+b)∣<c1/|\sigma^{\prime}(wx+b)|<c, which is quite likely to happen, the Lipschitz constraint is not necessarily satisfied. In our work, we use a novel approach for enforcing the Lipschitz constraints, avoiding the computation of the gradient which can be costly and difficult for point processes. We add the Lipschitz constraint as a regularization term to the empirical loss of RNN.

We can take each of the (2L2)2L\choose 2 pairs of real and generator sequences, and regularize based on them; however, we have seen that only a small portion of pairs (O(L)O(L)), randomly selected, is sufficient. The procedure of WGANTPP learning is given in Appendix E

Remark The significance of Lipschitz constraint and regularization (or more generally any capacity control) is more apparent when we consider the connection of W-distance and optimal transport problem . Basically, minimizing the W-distance between the empirical distribution and the model distribution is equivalent to a semidiscrete optimal transport . Without capacity control for the generator and discriminator, the optimal solution simply maps a partition of the sample space to the set of data points, in effect, memorizing the data points.

Experiments

Synthetic datasets. We simulate 20,000 sequences over time [0,T)[0,T) where T=15T=15, for inhomogeneous process (IP), self-exciting (SE), and self-correcting process (SC), recurrent neural point process (NN). We also create another 4 (=C43=C_{4}^{3}) datasets from the above 4 synthetic data by a uniform mixture from the triplets. The new datasets IP+SE+SC, IP+SE+NN, IP+SC+NN, SE+SC+NN are created to testify the mode dropping problem of learning a generative model. The parameter setting follows:

i) Inhomogeneous process. The intensity function is independent from history and given in Sec. 2.2, where k=3,α=,c=,σ=k=3,\alpha=,c=,\sigma=. ii) Self-exciting process. The past events increase the rate of future events. The conditional intensity function is given in Sec. 2.2 where μ=1.0,β=0.8\mu=1.0,\beta=0.8 and the decaying kernel g(t−ti)=e−(t−ti)g(t-t_{i})=e^{-(t-t_{i})}. iii) Self-correcting process. The conditional intensity function is defined in Sec. 2.2. It increases with time and decreases by events occurrence. We set η=1.0,γ=0.2\eta=1.0,\gamma=0.2. iv) Recurrent Neural Network process. The conditional intensity is given in Sec. 2.2, where the neural network’s parameters are set randomly and fixed. We first feed random variable from uniform distribution, and then iteratively sample events from the intensity and feed the output into the RNN to get the new intensity for the next step.

Real datasets. We collect sequences separately from four public available datasets, namely, health-care MIMIC-III, public media MemeTracker, NYSE stock exchanges, and publications citations. The time scale for all real data are scaled to , and the details are as follows:

i) MIMIC. MIMIC-III (Medical Information Mart for Intensive Care III) is a large, publicly available dataset, which contains de-identified health-related data during 2001 to 2012 for more than 40,000 patients. We worked with patients who appear at least 3 times, which renders 2246 patients. Their visiting timestamps are collected as the sequences. ii) Meme. MemeTracker tracks the meme diffusion over public media, which contains more than 172 million news articles or blog posts. The memes are sentences, such as ideas, proverbs, and the time is recorded when it spreads to certain websites. We randomly sample 22,000 cascades. iii) MAS. Microsoft Academic Search provides access to its data, including publication venues, time, citations, etc. We collect citations records for 50,000 papers. iv) NYSE. We use 0.7 million high-frequency transaction records from NYSE for a stock in one day. The transactions are evenly divided into 3,200 sequences with equal durations.

2 Experimental Setup

Details. The WGAN generator and discriminator are both RNNs with 64 hidden neurons without batch normalization, and the activation function is tanh. The Adam optimization method with learning rate 1e-4, β1=0.5,β2=0.9\beta_{1}=0.5,\beta_{2}=0.9, is applied and the batch size is 256.

Baselines. We compare the proposed method of learning point processes (i.e., minimizing sample distance) with maximum likelihood based methods for point process. To use MLE inference for point process, we have to specify its parametric model. The used parametric model are inhomogeneous Poisson process (mixture of Gaussian), self-exciting process, self-correcting process, and RNN. For each data, we use all the above solvers to learn the model and generate new sequences, and then we compare the generated sequences with real ones.

Evaluation metrics. Although our model is an intensity-free approach we will evaluate the performance by metrics that are computed via intensity. For all models, we work with the empirical intensity instead. Note that our objective measures are in sharp contrast with the best practices in GANs in which the performance is usually evaluated subjectively, e.g., by visual quality assessment. We evaluate the performance of different methods to learn the underlying processes via two measures: 1) The first one is the well-known QQ plot of sequences generated from learned model. The quantile-quantile (q-q) plot is the graphical representation of the quantiles of the first data set against the quantiles of the second data set. From the time change property of point processes, if the sequences come from the point process λ(t)\lambda(t) , then the integral Λ=∫titt+1λ(s)ds\Lambda=\int_{t_{i}}^{t_{t+1}}\lambda(s)ds between consecutive events should be exponential distribution with parameter 1. Therefore, the QQ plot of Λ\Lambda against exponential distribution with rate 1 should fall approximately along a 45-degree reference line. The evaluation procedure is as follows: i) The ground-truth data is generated from a model, say IP; ii) All 5 methods are used to learn the unknown process using the ground-truth data; iii) The learned model is used to generate a sequence; iv) The sequence is used against the theoretical quantiles from the model to see if the sequence is really coming from the ground-truth generator or not; v) The deviation from slope 1 is visualized or reported as a performance measure. 2) The second metric is the deviation between empirical intensity from the learned model and the ground truth intensity. We can estimate empirical intensity λ′(t)=E(N(t+δt)−N(t))/δt\lambda^{\prime}(t)=E(N(t+\delta t)-N(t))/\delta t from sufficient number of realizations of point process through counting the average number of events during [t,t+δt][t,t+\delta t], where N(t)N(t) is the count process for λ(t)\lambda(t). The L1L_{1} distance between the ground-truth empirical intensity and the learned empirical intensity is reported as a performance measure.

3 Results and Discussion

Synthetic data. Figure 2 presents the learning ability of WGANTPP when the ground-truth data is generated via different types of point process. We first compare the QQ plots in the top row from the micro perspective view, where QQ plot describes the dependency between events. Red dots legend-ed with Real are the optimal QQ distribution, where the intensity function generates the sequences are known. We can observe that even though WGANTPP has no prior information about the ground-truth point process, it can estimate the model better except for the estimator that knows the parametric form of data. This is quite expected: When we are training a model and we know the parametric form of the generating model we can find it better. However, whenever the model is misspecified (i.e., we don’t know the parametric from a priori) WGANTPP outperforms other parametric forms and RNN approach. The middle row of figure 2 compares the empirical intensity. The Real line is the optimal empirical intensity estimated from the real data. The estimator can recover the empirical intensity well in the case that we know the parametric form where the data comes from. Otherwise, estimated intensity degrades considerably when the model is misspecified. We can observe our WGANTPP produces the empirical intensity better and performs robustly across different types of point process data. The fact that the empirical intensity estimated from MLE-IP method are good and QQ plots are very bad indicates the inhomogeneous Poisson process can capture the average intensity (Macro dynamics) accurately but incapable of capturing the dependency between events (Micro dynamics). To testify that WGANTPP can cope with mode dropping, we generate mixtures of data from three different point processes and use this data to train different models. Models with specified form can handle limited types of data and fail to learn from diverse data sources. The last row of figure 2 shows the learned intensity from mixtures of data. WGANTPP produces better empirical intensity than alternatives, which fail to capture the heterogeneity in data. To verify the robustness of WGANTPP, we randomly initialize the generator parameters and run 10 rounds to get the mean and std of deviations for both empirical intensity and QQ plot from ground truth. For empirical intensity, we compute the integral of difference of learned intensity and ground-truth intensity. Table 1 reports the mean and std of deviations for intensity deviation. For each estimators, we obtain the slope from the regression line for its QQ plot. Table 1 reports the mean and std of deviations for slope of the QQ plot. Compared to the MLE-estimators, WGANTPP consistently outperforms even without prior knowledge about the parametric form of the true underlying generative point process. Note that for mixture models QQ-plot is not feasible.

Real-world data. We evaluate WGANTPP on a diverse real-world data process from health-care, public media, scientific activities and stock exchange. For those real world data, the underlying generative process is unknown, previous works usually assume that they are certain types of point process from their domain knowledge. Figure 3 shows the intensity learned from different models, where Real is estimated from the real-world data itself. Table 2 reports the intensity deviation. When all models have no prior knowledge about the true generative process, WGANTPP recovers intensity better than all the other models across the data sets.

Analysis. We have observed that when the generating model is misspecified WGANTPP outperforms the other methods without leveraging the a priori knowledge of the parametric form. However, when the exact parametric form of data is known and when it is utilized to learn the parameters, MLE with this full knowledge performs better. However, this is generally a strong assumption. As we have observed from the real-world experiments WGANTPP is superior in terms of performance. Somewhat surprising is the observation that WGANTPP tends to outperform the MLE-NN approach which basically uses the same RNN architecture but trained using MLE. The superior performance of our approach compared to MLE-NN is another witness of the the benefits of using W-distance in finding a generator that fits the observed sequences well. Even though the expressive power of the estimators is the same for WGANTPP and MLE-NN, MLE-NN may suffer from mode dropping or get stuck in an inferior local minimum since maximizing likelihood is asymptotically equivalent to minimizing the Kullback-Leibler (KL) divergence between the data and model distribution. The inherent weakness of KL divergence renders MLE-NN perform unstably, and the large variances of deviations empirically demonstrate this point.

Conclusion and Future Work

We have presented a novel approach for Wasserstein learning of deep generative point processes which requires no prior knowledge about the underlying true process and can estimate it accurately across a wide scope of theoretical and real-world processes. For the future work, we would like to explore the connection of the WGAN with the optimal transport problem. We will also explore other possible distance metrics over the realizations of point processes, and more sophisticated transforms of point processes, particularly those that are causal. Extending the current work to marked point processes and processes over structured spaces is another interesting venue for future work.

References

Appendix A Data flow of Wassterstein learning for point process

Figure 4 illustrates the data flow for WGANTPP.

Appendix B Proof that ∥⋅∥⋆\|\cdot\|_{\star} is a norm

It is obvious that ∥⋅∥⋆\|\cdot\|_{\star} is nonnegative and symmetric. If ∥ξ−ρ∥⋆=0\|\xi-\rho\|_{\star}=0, then m=nm=n and there is a assignment σ\sigma such that xi=yσ(i)x_{i}=y_{\sigma(i)} for all i=1,…,ni=1,\dots,n.

Now we prove that ∥⋅∥⋆\|\cdot\|_{\star} has triangle inequality. WLOG, assume that ξ={x1,…,xn}\xi=\{x_{1},\dots,x_{n}\}, ρ={y1,…,yk}\rho=\{y_{1},\dots,y_{k}\} and ζ={z1,…,zm}\zeta=\{z_{1},\dots,z_{m}\} where n≤k≤mn\leq k\leq m. Define the permutation σ^\hat{\sigma} on {1,…,k}\{1,\dots,k\} by

where the last equality is due to the fact that the minimization is taken over all permutations σ\sigma of {1,…,m}\{1,\dots,m\}, and σ^\hat{\sigma} is a fixed permutation of {1,…,k}\{1,\dots,k\} where k≤mk\leq m. This completes the proof.

Appendix C Proposed ∥⋅∥⋆\|\cdot\|_{\star} Distance on the Real Line

In this section, we prove that finding the distance between sequences ξ\xi and ρ\rho,

in the case of temporal point process in [0,T)[0,T), i.e., ξ={t1<t2<…<tn}\xi=\left\{t_{1}<t_{2}<\ldots<t_{n}\right\} and ρ={τ1<τ2<…<τm}\rho=\left\{\tau_{1}<\tau_{2}<\ldots<\tau_{m}\right\}, reduces to

Here, without loss of generality n≤mn\leq m is assumed. The choice of s=Ts=T is basically padding the shorter sequences with TT. Given, the sequences have the same length now, we claim that the identity permutation i.e., σ(i)=i\sigma(i)=i is the minimizer in (14). We proceed by a proof by contradiction. Assume that the minimizer is NOT the identity permutation. Then, find the first ii such that σ(i)≠i\sigma(i)\neq i. Then, Σ(i)=j\Sigma(i)=j where j>ij>i. Therefore, there should be a k>ik>i such that σ(k)=i\sigma(k)=i. Then, if you change the permutation according to σ(i)=i\sigma(i)=i and σ(k)=j\sigma(k)=j the cost will change by

Given i<ji<j and i<ki<k, it is easy to see that Δ>0\Delta>0. This means that we’ve found a better permutation which contradicts our assumption. Therefore, the optimal permutation will match the event points in an increasing order one by one.

Appendix D Equivalence of the ∥⋅∥⋆\|\cdot\|_{\star} Distance and Difference in Count Measures

The count measure of a temporal point process is a special case of the one defined for point processes in general space in Section 2.1. For a Borel subset B⊂S=[0,T)B\subset S=[0,T) we have N(B)=∫t∈Bξ(t)dtN(B)=\int_{t\in B}\xi(t)dt. With a little abuse of notation we write N(t):=N([0,t))=∫0tξ(t)dtN(t):=N([0,t))=\int_{0}^{t}\xi(t)dt. Figure 1 is a good guidance through this paragraph. Starting from time 0 the first gap in count measure starts from min⁡(t1,τ1)\min(t_{1},\tau_{1}) and ends in max⁡(t1,τ1)\max(t_{1},\tau_{1}). Therefore, there is difference equal to s1=max⁡(t1,τ1)−min⁡(t1,τ1)=∣t1−τ1∣s_{1}=\max(t_{1},\tau_{1})-\min(t_{1},\tau_{1})=|t_{1}-\tau_{1}| in the count measure. Similarly, the second block of difference has volume of s2=∣t2−τ2∣s_{2}=|t_{2}-\tau_{2}|, and so on. Finally, for m>nm>n the (n+i)(n+i)-th block make a difference of sn+i=T−τn+is_{n+i}=T-\tau_{n+i}. Therefore, the area (L1L_{1} distance) between the two sequences is a equal to S=∑i=1msiS=\sum_{i=1}^{m}s_{i}. On the other hand by looking (15) we observe that ∥ξ−ρ∥⋆=∑i=1msi\|\xi-\rho\|_{\star}=\sum_{i=1}^{m}s_{i}. Therefore, by choice of s=Ts=T as an anchor point, the distance we have is exactly the area between the two count measures.

Appendix E WGANTPP algorithm

The procedure of WGANTPP learning is given in Algorithm 1