# Joint Modeling of Event Sequence and Time Series with Attentional Twin Recurrent Neural Networks

###### Abstract

A variety of real-world processes (over networks) produce sequences of data whose complex temporal dynamics need to be studied. More especially, the event timestamps can carry important information about the underlying network dynamics, which otherwise are not available from the time-series evenly sampled from continuous signals. Moreover, in most complex processes, event sequences and evenly-sampled times series data can interact with each other, which renders joint modeling of those two sources of data necessary. To tackle the above problems, in this paper, we utilize the rich framework of (temporal) point processes to model event data and timely update its intensity function by the synergic twin Recurrent Neural Networks (RNNs). In the proposed architecture, the intensity function is synergistically modulated by one RNN with asynchronous events as input and another RNN with time series as input. Furthermore, to enhance the interpretability of the model, the attention mechanism for the neural point process is introduced. The whole model with event type and timestamp prediction output layers can be trained end-to-end and allows a black-box treatment for modeling the intensity. We substantiate the superiority of our model in synthetic data and three real-world benchmark datasets.

## 1 Introduction

Event sequences are becoming increasingly available in a variety of applications. Such event sequences, which are asynchronously generated with random timestamps, are ubiquitous in areas such as e-commerce, social networks, electronic health data, and equipment failures. The event data can carry rich information not only about the event attribute (e.g., type, participator) but also the timestamp indicating when the event takes place. A major line of research [1] has been devoted to studying event sequence, especially exploring the timestamp information to model the underlying dynamics of the system, whereby point process has been a powerful and elegant framework in this direction.

Being treated as a random variable when the event is stochastically generated in an asynchronous manner, the timestamp makes the event sequence of point processes fundamentally different from the time series [2] evenly-sampled from continuous signals because the asynchronous timestamps reflect the network dynamic while the time for time-series is deterministic..

However these time series data, when available, provide timely updates of background environment where events occur in the temporal point process, such as temperature for computing servers or blood pressure for patients. Many complex systems posses such time series data regularly recorded along with the point processes data.

While there have been many recent works on modeling continuous-time point processes [3, 4, 5, 6, 7] and time series [8, 9, 10], most of them treat these two processes independently and separately, ignoring the influence one may have on the other over time. To better understand the dynamics of point processes, there is an urgent need for joint models of the two processes, which are largely inexistent to date. There are related efforts in linking the time series and event sequence to each other. In fact, one popular way to convert a time series to an event sequence is by detecting multiple events (e.g., based on thresholding the stock price series [11]) from the series data. On the other hand, statistical aggregation (e.g., total number of counts) is often performed on each time interval with equal length to extract aligned time series data from the event sequences. However such a coarse treatment can lead to the key information loss about the actual behavior of the process, or at least in a too early stage.

Recent progresses on modeling point process includes mathematical reformulations and optimization techniques [12, 13, 14, 5] and novel parametric forms [15, 16, 17, 18] as carefully designed by human prior knowledge to capture the characters of the dataset in their study. One major limitation of those model is that the specified form of point process limits its capability to capture the dynamic of data. Moreover, it may suffer from misspecification, for which the model is not suitable for the data. Recent works e.g., [13] start to turn to non-parametric form to fit the structure of a point process, but their method is under the Hawkes process formulation, which runs the risk of unknown model complexity and can be inappropriate for point processes that disobey the self or mutual-exciting rule assumed by the Hawkes model [19]. In another recent work [6], the authors proposed a semi-parametric pesudo point process model, which assumes a time-decaying influence between events and the background of intensity is constant. Besides, the event type is regarded as the mark associated with a univariate process.

In this paper, we view the conditional intensity of a point process as a nonlinear mapping from the joint embedding of time series and past event data to the predicted transient occurrence intensity of future events with different types. Such a nonlinear mapping is expected to be complex and flexible enough to model various characters of real event data for its application utility, e.g., failure prediction, social network analysis, and disease network mining as will be empirically studied in the experiment part of this paper. To overcome the disadvantages associated with the explicit parametric form of intensity function we bypass direct modeling of the intensity function and directly model the next event time and dimension. Neural networks are our choice to model this nonparametric mapping.

We utilize the state-of-the-art deep learning techniques to efficiently and flexibly model the intensity function of the point processes. Instead of predefining the form of point process, we turn to the synergic multi-RNNs (specifically twin-RNNs in this paper) as a natural way to encode such nonlinear and dynamic mapping, in an effort for modeling an end-to-end nonlinear intensity mapping without any prior knowledge. To further improve its interpretability, we infuse the attention model to improve its capability for both prediction and relational mining among event dimensions.

Key idea and highlights. Our model interprets the conditional intensity function of a point process as a nonlinear mapping, which is synergetically established by a composite neural network with two RNNs as its building blocks. As illustrated in Fig. 1, time series (top row) and event sequence (bottom row) are distinct to each other. The underlying rationale is that time series is more suitable to carry the synchronously and regularly updated (i.e. in a fixed pace) or constant profile features. In contrast, the event sequence can compactly catch event driven, more abrupt information, which can affect the conditional intensity function over longer period of time.

More specifically, We first argue that many conditional intensity functions can be viewed as an integration of two effects: i) spontaneous background component inherently affected by the internal (time-varying) attributes of the individual; ii) effects from history events. Meanwhile, most information in the real world can also be covered by continuously updated features like age, temperature, and asynchronous event data such as clinical records, failures. This motivates us to devise a general approach. Then, we use one RNN whose units are aligned with the time points of a time series, and another RNN whose units are aligned with events. The time series RNN can timely update the intensity function while the event sequence RNN is used to efficiently capture the long-range dependency over history. They can interact with each other through synergic non-linear mapping. This allows fitting arbitrary dynamics of point process which otherwise will be difficult or often impossible to be specified by a parameterized model restricted to certain assumptions.

As an extension to the conference version [20]^{1}^{1}1The main extensions include: i) in contrast to [20] where the so-called pseudo multi-dimensional point process [21] is adopted which involves only a single intensity function for all types of event sequence, here we separately model the intensity functions for each event type leading to the so-called genuine multi-dimensional point process via recurrent neural networks; ii) based on the resulting multi-dimensional point process model, we incorporate a new attention mechanism to improve the interpretability of the prediction model. This expands the dependency model in RNN from the recent hidden variable to a set of recent ones which further improves its modeling capability; iii) a more thorough analysis and verification via devised simulation experiments; iv) performing more experiments from both social network data and healthcare data. Note that the added attention model enables the relation discovery capability as also verified in both simulation based and real-world data. While the conference paper [20] only deals with event prediction rather than relation mining., the overall highlights of this paper are:

i) To the best of our knowledge, this is the first work to jointly interpret and instantiate the conditional intensity function with fused time series and event sequence RNNs. This opens up the room for connecting the neural network techniques to traditional point process that emphasizes more on specific model driven by domain knowledge. The introduction of a full RNN treatment lessen the efforts for the design of (semi-)parametric point process model and its complex learning algorithms which often call for special tricks e.g. [22] that prohibiting the wide use for practitioners. In contrast, neural networks and specifically RNNs, are becoming off-the-shelf tools and getting widely used recently.

ii) We model the genuine multi-dimensional point process through recurrent neural networks. Previous work [6] use the RNN to model so-called pseudo multi-dimensional point process [21]. Actually, they regard the event sequence as a univariate point process and treat the dimension as the mark associated with events. Consequently, there exists only one intensity function for all the processes instead of one per each dimension. On the contrary, in our work the process is the result of the superposition of sub-processes for each dimension. In this way, we can separate the parameters of each dimension as well as capture their interactions. This leads to more effective simulation and learning algorithm.

iii) To improve the interpretability, it is also perhaps the first time, to our knowledge, an attention based RNN model for point process is proposed. For multi-dimensional point process, our proposed attention mechanism allows each dimensional has its own attention parameters. One typical resulting utility involves decision support and causality analysis [23].

iv) Our model is simple and general and can be end-to-end trained. We target three empirical application domains to demonstrate the superiority of the proposed method, namely, predictive maintenance, social network analysis and disease relation mining. The state-of-the-art performance on relational mining, event type and timestamp prediction corroborates its suitability to real-world applications.

## 2 Related Work and Motivation

We review the related concepts and work in this section, which is mainly focused on Recurrent Neural Networks (RNNs) and their applications in time series and sequence data, respectively. Then we discuss existing point process methods and their connection to RNNs. All these observations motivate the work of this paper.

Recurrent neural network. The building block of our model is the Recurrent Neural Networks (RNNs) [24, 25] and its modern variants e.g., Long Short-Term Memory (LSTM) units [26, 27] and Gated Recurrent Units (GRU) [28]. RNNs are dynamical systems whose next state and output depend on the present network state and input, which are more general models than the feed-forward networks. RNNs have long been explored in perceptual applications for many decades, however it can be very difficult for training RNNs to learn long-range dynamics in part due to the vanishing and exploding gradients problem. LSTMs provide a solution by incorporating memory units that allow the network to learn when to forget previous hidden states and when to update hidden states given new information. Recently, RNNs and LSTMs have been successfully applied in large-scale vision [29], speech [30] and language [31] problems.

RNNs for series and event data. From application perspective, we consider two main scenarios in this paper: i) RNNs for synchronized series with evenly spaced interval e.g., time series or indexed sequence with pure order information e.g., language; ii) asynchronous sequence with timestamp e.g., event data.

i) Synchronized series: RNNs have been a long time a natural tool for standard time series modeling and prediction [32, 33], whereby the indexed series data point is fed as input to an (unfold) RNN. In a broader sense, video frames can also be treated as time series and RNN are widely used in recent visual analytics works [34] and so for speech [30]. RNNs are also intensively adopted for sequence modeling tasks [28] when only order information is considered.

ii) Asynchronous event: In contrast, event sequence with timestamp about their occurrence, which are asynchronously and randomly distributed over the continuous time space, is another typical input type for RNNs [6, 35] (despite its title for ’time series’). One key differentiation against the first scenario is that the timestamp or time duration between events (together with other features) is taken as input to the RNNs. By doing so, (long-range) event dependency can be effectively encoded.

Interpretability and attention model. Prediction accuracy and model interpretability are two goals of many successful predictive methods. Existing works often have to suffer the tradeoff between the two by either picking complex black box models such as deep neural network or relying on traditional models with better interpretation such as Logistic Regression often with less accuracy compared with state-of-the-art deep neural network models. Despite the promising gain in accuracy, RNNs are relatively difficult to interpret. There have been several attempts to interpret RNNs [35, 36, 37]. However, they either compute the attention score by the same function regardless of the affected point’s dimension [35], or only consider the hidden state of the decoder for sequence prediction [36, 37]. As for multi-dimensional point process, past events shall influence the intensity function differently for each dimension. As a result, we explicitly assign different attention function for each dimension which is modeled by respective intensity functions, thus leading to an infectivity matrix based attention mechanism which will be detailed later in this paper.

Point processes. Point process is a mathematically rich and principled framework for modeling event data [1]. It is a random process whose realization consists of a list of discrete events localized in time. The dynamics of the point process can be well captured by its conditional intensity function whose definition is briefly reviewed here: for a short time window , represents the rate for the occurrence of a new event conditioned on the history :

where is the expectation of the number of events happened in the interval given the historical observations . The conditional intensity function has played a central role in point processes and many popular processes vary on how it is parameterized. Some typical examples include:

1) Poisson process [38]: the homogeneous Poisson process has a very simple form for its intensity function: . Poisson process and its time-varying generalization are both assumed to be independent of the history.

2) Reinforced Poisson processes [39, 15]: the model captures the ‘rich-get-richer’ mechanism characterized by a compact intensity function, which is recently used for popularity prediction [15].

3) Hawkes process [19]: Recently, Hawkes process has received a wide attention in network cascades modeling [14, 5], community structure [40], viral diffusion and activity shaping[7], criminology [41], optimization and intervention in social networks [42], recommendation systems [43], and verification of crowd generated data [44]. As an illustration example intensively used in this paper, we particularly write out its intensity function is:

where is the infectivity matrix, indicating the directional influence strength from dimension to . It explicitly uses a triggering term to model the excitation effect from history events where the parameter denotes the decaying bandwidth. The model is originally motivated to analyze the earthquake and its aftershocks[45].

4) Reactive point process [16]: it can be regarded as a generalization of the Hawkes process by adding a self-inhibiting term to account for the inhibiting effects from history events.

5) Self-correcting process [46]: its background part increases steadily, while it is decreased by a constant every time a new event appears.

We summarize the above forms in Table I. It tries to separate the spontaneous background component and history event effect explicitly. This also motivates us to design an RNN model that can flexibly model various point process forms without model specification.

Model | Background | History event effect |
---|---|---|

Poisson process | ||

Reinforced poisson process | ||

Hawkes process | ||

Reactive point process | ||

Self-correcting process |

Note: is Dirac function, is the time-decaying kernel and can be constant or time-varying function.

## 3 Network Structure and Learning

In this section, we will present the proposed network structure along with the learning algorithm for modeling the behavior of dynamic events.

### 3.1 Brief on RNN as building block

Taking a sequence as input, the RNN generates the hidden states , also known high-level representation of inputs [24, 25]:

where is the profile associated with each event, and is a non-linear function, and are parameters to be learned. One common choice for non-linear function is Sigmoid or tanh, who suffers from vanishing-gradients problem [47] and poor long-range dependency modeling capability. In contrast, we implement our RNN with Long Short Term Memory (LSTM) [26, 27] for its popularity and capability for capturing long-range dependency. In fact, other RNN variants e.g. Gated Recurrent Units (GRU) [28] can also be alternative choices, while the analysis of the consequence of this particular choice is not the focus of our paper. To make the presented paper self-contained, we reiterate the formulation of LSTM as follows:

where denotes element-wise multiplication and the recurrent activation is the Logistic Sigmoid function. Unlike standard RNNs, the Long Short Term Memory (LSTM) architecture uses memory cells to store and output information, allowing it to better discover long-range temporal relationships. Specifically, , , , and are respectively the input gate, forget gate, output gate, and cell activation vectors^{2}^{2}2In subsection 3.1 we slightly abuse the notations and in fact their effects are only valid in this subsection and have no relation to other notations used in the other parts of this paper.. By default, the value stored in the LSTM cell is maintained unless it is added to by the input gate or diminished by the forget gate . The output gate controls the emission of the memory from the LSTM cell.
Compactly, we represent the LSTM system via the following equation:

### 3.2 Infectivity matrix based attention mechanism

Now we give further formal notational definitions used in this paper. Time series data is denoted by , e.g. a patient’s temperature and blood pressure recorded in different dimensions of the vector . Event data is represented as , where is the dimension representing a categorical information (e.g a mark or type or an agent), where is the finite set of all event types, and is occurrence time. The former can timely affect the transient occurrence intensity of events and the latter can often abruptly cause jump and transition of states of agents and capture long-range event dependency [48].

As shown in Fig. 3, for our proposed network, these two sources of data are fed separately to two RNNs (we call it a twin RNN structure in this paper) whose outputs are further combined to serve as the inputs of subsequent layers. For event sequence with length , we can generate a hidden variable sequence as the high-level representation of input sequence in RNN. To predict the dimension (e.g. event type) and time for the -th event, the history , prior to that event should be utilized. The most recent is often regarded as a compressed representation of .

One may argue the necessity for involving a twin RNN structure, since the instantaneous time series data can be sampled and fed into a single RNN along with the event data when an event occurs. However, there are particular advantages for adopting such a twin-RNN structure. Empirically, it has been shown in some recent study [34] that using separate RNN for each time series data e.g., video and time series sensors can lead to better prediction accuracy than a single RNN fed with the combination of the multiple time series data. More importantly, the events can occur in arbitrary timestamp i.e., they are asynchronous while the time series data is often sampled with equal time interval being a synchronous sequence. This inherent difference inspires us to model these two different sequence data via separate RNN as their dynamics can be rather varying.

However, there are still two limitations for the above approach: i) The prediction capability or model expressiveness is limited: in fact only the recently updated hidden variable is used for prediction regardless of the length of input sequence. ii) The interpretability is limited. As we compress all information into a fixed vector , it is hard to infer which event contributes most to the prediction. For example in the problem of multi-dimensional Hawkes process learning, one important goal is to uncover the hidden network structure, infectivity matrix , from real-world event sequences, such as the influence strength between users in social network [14, 5], or progression relationship between event types [17]. Uncovering the hidden structure is also stressed in causal analysis [23], which gives much evidence for prediction result. This calls for particular mechanisms to improve its flexibility and interpretability.

In this work, we devise a temporal attention mechanism to enable interpretable prediction models for point processes. Events from a certain dimension may have higher influence on some dimensions. We exploit this observation to make the trained neural network model more interpretable and expressive. To achieve this, we first expand the representation of to be a set of vectors instead of only the most recent , referred as context vectors. Each of them is localized to its respective preceding event of interest from inputs . Inspired by the Hawkes process, the influence strength from to is introduced and it is modeled by:

(1) |

where is the feature vector to be learned for the particular prediction dimension and is the score function which gives the influence strength from to . Once the influence strength are computed, we can generate the representation vector for the next layer:

(2) |

where is the attention function, which computes the final representation of . Here we choose the widely used soft attention mechanism [36], whose influence from former events is in an additive form [5]:

(3) |

Note that hard attention [37] only assigns 0 or 1 to the influence strength , which is too rough to capture fine-grained influence. Since the cost is differentiable with respect to the parameters, we can easily train the network end-to-end using backpropagation.

After the model is trained i.e. the parameter is fixed, for each testing event sequence and its computed by Eq.1, we define the infectivity matrix to reflect the mutual influence among dimensions as , where represents the average of all divided by .

The attention mechanism is depicted in Fig. 2. This attention mechanism can allow the network to refer back to the preceding hidden variables , instead of forcing it to encode all information into the recent . It can retrieve from internal memory and choose what to attend to.

Finally we point out that in the context of point process, our attention mechanism is different from the existing work [35, 36, 37] in that they only consider one-way effect over the sequence. In another word, their approaches are current state agnostic. This leads to a vector representation (similar to the role of the vector used in Eq.1) for the weight variables instead of a two-way infectivity matrix . Moreover, to make the model tractable, we use a parameterized form by Eq.1 to represent the two-way weights.

### 3.3 Network structure

Now we give the full description of our network whose overview is depicted in Fig. 3. The dashed part in the right of figure is illustrated in detail in Fig. 2. For time series data e.g., temperature, blood pressure, they are sampled evenly over time. We use to indicate the dense feature vector sampled at different timestamps. Those signals are expected to reflect the states of each dimension and drive the occurrence intensity of events. Hence we have:

For event sequence , we can generate hidden states through LSTM, which can capture long-range dependency of events. First we project the dimension to a low-dimensional embedding vector space. Then the embedding vector combined with timestamps is fed to LSTM. We use the following equation to represent the process:

where denotes the embedding vector of input and the embedding matrix to learn.

For dimension , its final representation of is , which is obtained through the attention mechanism introduced in Eq. 3. For the score function of Eq. 1, we specify it by:

(4) |

When the context vector is similar to , the attention function produces a large score. Otherwise a small one. In an extreme case, the score is zero when the context vector is orthogonal to feature vector of dimension . To promote sparsity of infectivity matrix, we threshold the score with a minus operation. The threshold can control the degree of sparsity which is set to throughout this paper. Note the form of Eq. 4 is also used in [36, 37] to model the one-way attention weight . From this, it is clear that our model for the attention is two-way between dimension to as shown in Eq. 1.

To jointly model event sequence and time series, we combine them into one synergic layer as illustrated in Fig. 3:

(5) |

where is the concatenation of the two vectors. The synergic layer can be any function, coupling two data sources together. Here we use the Sigmoid function. As a result, we obtain a representation for the output dimension . We can use this representation to compute the intensity for each dimension and then simulate its next occurrence time and its dimension jointly. Here we take a more efficient approach by firstly predicting the next event’s dimension and then further predicting the occurrence timestamp based on the predicted event dimension. Note that the intensity function is modeled implicitly within the neural network architecture and we directly model the timing and dimension of the events. In this way, we overcome the expensive computation cost from explicit parametric form of intensity function.

To predict the next event’s dimension , we apply the Softmax operation to those representations where is the number of event dimensions.

(6) |

where are model parameters to learn. The optimal dimension (as the prediction result) is computed by selecting the corresponding maximum element in :

(7) |

After we obtain the optimal dimension , we use the representation to derive occurrence time following:

(8) |

where are model parameters for learning.

### 3.4 End-to-end learning

The likelihood of observing a sequence along with time series signals can be expressed as follows:

(9) |

where the weight parameters are set as the inverse of the sample count in that dimension against the total size of samples, to weight more on those dimensions with fewer training samples. This is in line with the importance weighting policy for skewed data in machine learning [49].

For the second term, the underlying rationale is that we not only encourage correct prediction of the coming event dimension, but also require the corresponding timestamp of the event to be close to the ground truth. We adopt a Gaussian penalty function:

As shown in Fig. 3, the output from the timestamp prediction layer is fed to the classification loss layer to compute the above penalty given the actual timestamp . We adopt RMSprop gradients [50] which have been shown to work well on training deep networks to learn these parameters.

By directly optimizing the loss function, we learn the prediction model in an end-to-end manner without the need for sophisticated or carefully designed algorithms (e.g., Majorization-Minimization techniques [51, 5]) used in generative Point process models. Moreover, as pointed out by recent work [18], another limitation for the generative point process model is that they are aimed to maximize the joint probability of all observed events via a maximum likelihood estimator, which is not tailored to the prediction task.

## 4 Experiments and Discussion

We evaluate the proposed approach on both synthetic and real-world datasets, from which three popular application scenarios are covered: social network analysis, electronic health records (EHR) mining, and proactive machine maintenance. The first two scenarios involve public available benchmark dataset: MemeTracker and MIMIC, while the latter involves a private ATM maintenance dataset from a commercial bank headquartered in North America.

The code is based on Theano running on a Linux server with 32G memory, 2 CPUs with 6 cores for each: Intel(R) Xeon(R) CPU E5-2603 v3@1.60GHz. We also use 4 GPU:GeForce GTX TITAN X with 12G memory backed by CUDA and MKL for acceleration.

### 4.1 Baselines and evaluation metrics

We compare the proposed method to the following algorithms and state-of-the-art methods:

1) Logistic model: We use Logistic regression for event timestamp prediction and an independent Logistic classification model for event type prediction. To make sure the proposed method and the logistic model use the same amount of information, the predictor features in the regression are comprised of the concatenation of feature vectors for sub-windows of all active time series RNN.

2) Hawkes Process: To enable multi-type event prediction, we use a Multi-dimensional Hawkes process [14, 7]. The Full, Sparse, LowRankSparse indicate the different types of Hawkes Process model. The full model has no regularization on infectivity matrix while the Sparse and LowRankSparse ones have sparse and lowrank-sparse regularization, respectively. In the following, Hawkes process indicates LowRankSparse model if not explicitly mentioned. The inputs are comprised of event sequences with dimensions and timestamps.

3) Recurrent Marked Temporal Point Processes (RMTPP): [6] uses a neural network to model the event dependency flexibly. The inputs are event sequences with continuous signals sampled when they happen. The method can only sample features of transient time series when the events happen and use partially parametric form for the base intensity and a time-decaying influence from former event to the current one. Another difference is that it assumes an independent distribution between time and marks and predicts the dimension and time independently given the learned representation of the history.

4) TRPP: This method uses time series and event sequences to collaboratively model the intensity function of point process. The model uses RNN to model the non-linear mapping from history to the predicted marker and time, and treats RNN as a block box without much interpretability. Here we rename it as Twin Recurrent Point Processes (TRPP). This method is the one presented in the conference version of this paper [20].

5) ERPP: We also include TRPP’s degraded version by only keeping the event sequence as input and removing the time series RNN, which is termed by Event Recurrent Point Processes (ERPP). Including the term ‘event’ also helps distinguish it from the existing term RPP: Reinforced Poisson Processes [39, 15].

6) Markov Chain (MC): The Markov chain refers to the sequence of random variables, with the Markov property that future state only depends on the current state and is conditionally independent of the history. The order of Markov chain indicates how many recent states on which the future state depends. As this model can only learn the transition probability of dimensions, we use it to predict the dimensions without taking the time into account. The optimal order of Markov chain is determined by the performance on separate validation dataset.

7) Continuous Time Markov Chain (CTMC): The CTMC is a special type of semi-Markov model, which models the continuous transition among dimensions as a Markov process. It predicts the next dimension with the earliest transition time, therefore, it can jointly predict time and dimension.

8) Homogeneous Poisson Process: This method implements the most basic point process model in which the intensity function is constant and events occur independently. It can estimate interval-event gaps.

9) Self-correcting Process: When the occurrence of an event decrease the probability of other events, we are facing a variant of point process called self-correcting processes. Its intensity function is shown in Table I and it can estimate the inter-event time.

Inline with the above TRPP and ERPP methods, we term our model Attentional Twin Recurrent Point Processes (ATRPP). To further study the effect of the time series RNN, we also evaluate the baseline version without using this channel, which is termed as Attentional Event Recurrent Point Processes (AERPP).

Evaluation metrics. We use several common metrics for performance evaluation. For the next event dimension prediction, we adopt Precision, Recall, F1 Score and Confusion matrix. For event time prediction, we use the Mean Absolute Error (MAE) which measures the absolute difference between the predicted time point and the actual one. For the infectivity matrix, we use RankCorr [14] to measure whether the relative order of the estimated influence strength is correctly recovered, when the true infectivity matrix is available. The RankCorr is defined as the averaged Kendall rank correlation coefficient^{3}^{3}3https://en.wikipedia.org/wiki/Kendall_rank_correlation_coefficient between each row of ground-truth and estimated infectivity matrix. RelErr measures the relative deviation between estimated and and ground-truth , which is defined as the average of .

### 4.2 Experiments on synthetic data

The test on synthetic data provides a quantitative way for evaluating the performance of compared methods.

Synthetic data generation. We use simulated data with known ground-truth to quantitatively verify our model by aforementioned metrics. Specifically, we generate cascades from multi-dimensional hawkes process via the Thinning algorithm [52]. We choose for the number of event dimensions. The background intensity term is set uniformly at random: . Mutual influence is set similarly to . Half of the elements in the infectivity matrix are set to 0 by random in order to mimic the sparsity of influence between dimensions in many real-world problems. For the stability of the simulation process, the matrix is scaled such that its spectral radius is no larger than one. The decaying parameter in the Hawkes process is set to . The detailed description of these parameters can be found in Sec. 1 for the Hawkes process. We simulate 5000 independent cascades, 3000 for training, 1000 for validating, and 1000 for testing, respectively. To generate time series signals, , with some explanation capability to the background intensity, we sample from for all dimensions , where, is a Gaussian noise, .

Experimental results. The performance on synthetic data is shown in Fig. 4. The relative error of infectivity matrix, RankCorr is demonstrated to verify the capability of uncovering hidden network structure among those dimensions i.e., the nodes in the network. The accuracy of time and dimension prediction is shown in order to compare the predictive performance. Our model achieves a better prediction performance, and meanwhile uncovers the infectivity matrix better than the alternatives. The self-correcting process suffers from model misspecification and performs worse than other point process models as the events are self-exciting not self-correcting. Our non-parametric model can learn from data and generalize well without prior knowledge of data.

### 4.3 Predictive machine maintenance

Predictive maintenance is a sound testbed for our model. It involves equipment risk prediction to allow for proactive scheduling of corrective maintenance. Such an early identification of potential concerns helps deploy limited resources more efficiently and cost effectively, reduce operations costs and maximize equipment uptime. Predictive maintenance is adopted in a wide variety of applications such as fire inspection, data center and electrical grid management e.g. [16]. For its practical importance in different scenarios and relative rich event data for modeling, we target our model to a real-world dataset of more than 1,000 automated teller machines (ATMs) from a global bank headquartered in North America.

We have no prior knowledge on the dynamics of the complex system and the task can involve arbitrarily working schedules and heterogeneous mix of conditions. It takes much cost or even impractical to devise specialized models.

The studied dataset is comprised of the event logs involving error reporting and failure tickets, which is originally collected from 1,554 ATMs. The event log of error records includes device identity, timestamp, message content, priority, code, and action. A ticket (TIKT) means that maintenance will be conducted. Statistics of the data is presented in Table II. The error type indicates which component encounters an error: 1) printer (PRT), 2) cash dispenser module (CNG), 3) Internet data center (IDC), 4) communication part (COMM), 5) printer monitor (LMTP), 6) miscellaneous e.g., hip card module, usb (MISC). The time series here consists of features: i) the inventory information: ATM type, age, operations, temperatures; ii) event frequency for each event in the recent an hour interval. The event types and their occurrence time from the ATMs are an event sequences. Therefore, there are 1554 sequences in total, which are randomly divided into training (50%), validating (20%) and testing (30%) portions.

type | total | max | min | mean | std |
---|---|---|---|---|---|

TIKT | 2226(–) | 10(137.04) | 0(1.21) | 2.09(31.70) | 1.85(25.14) |

PRT | 9204(–) | 88(210.13) | 0(0.10) | 8.64(12.12) | 11.37(21.41) |

CNG | 7767(–) | 50(200.07) | 0(0.10) | 7.29(15.49) | 6.59(23.87) |

IDC | 4082(–) | 116(206.61) | 0(0.10) | 3.83(23.85) | 5.84(30.71) |

COMM | 3371(–) | 47(202.79) | 0(0.10) | 3.16(22.35) | 3.90(29.36) |

LMTP | 2525(–) | 81(207.93) | 0(0.10) | 2.37(22.86) | 4.41(34.56) |

MISC | 1485(–) | 32(204.41) | 0(0.10) | 1.39(24.27) | 2.54(34.38) |

model | precision | recall | F1 score | MAE |
---|---|---|---|---|

Poisson | —– | —– | —– | 4.76 |

SelfCorrecting | —– | —– | —– | 4.65 |

Markov Chain | 0.530 | 0.591 | 0.545 | —– |

CTMC | 0.516 | 0.554 | 0.503 | 5.16 |

Logistic | 0.428 | 0.375 | 0.367 | 4.51 |

Hawkes | 0.459 | 0.514 | 0.495 | 5.43 |

RMTPP | 0.587 | 0.640 | 0.607 | 4.31 |

TRPP | 0.607 | 0.661 | 0.626 | 4.18 |

ERPP | 0.559 | 0.639 | 0.599 | 4.37 |

ATRPP | ||||

AERPP | 0.599 | 0.672 | 0.617 | 3.98 |

Table III shows the averaged performance of the proposed method compared to the alternatives. The confusion matrix for the seven event types are shown in Fig. 5 by all methods. Not surprisingly, for both event type and timestamp prediction, our main approach, i.e., ATRPP outperforms by a notable margin. ATRPP report 0.634 F1 score and 3.92 MAE while AERPP reaches 0.617 F1 score and 3.98 MAE. Obviously, this verifies that synergically modeling event sequence and time series can boost the performance of predicting future events and time.Interestingly, all point process based models obtain better results on this task which suggests they are more promising compared to classical classification models. Indeed, our methodology provides an end-to-end learning mechanism without any pre-assumption in modeling point process. All these empirical results on real-world tasks suggest the efficacy of our approach, especially in capturing the temporal dynamics of events data.

Visualization of influence pattern. We visualize the infectivity matrix of ATRPP as in Fig. 6. Each node denotes one dimension which here represents one type of events. The directed edge means the influence strength from source node to destination node. The size of nodes and depth of color is proportional to the weighted degree of nodes, which indicate the total influence of the node has on others. The width of edges is is proportional to the strength of influence. Self-loop edges are located at the right of nodes without an arrow. Note this setting applies to the subsequent two experiments. As is shown, it’s obvious that TIKT (maintenance) have a strong influence over all types of errors (breakdown) as maintenance can greatly decrease the probability of breakdown of machines. Also, self-loop edge of TIKT node is too small to see which indicates that maintenance has low correlation itself. All types of errors have self-loop, indicating a recurrent pattern of errors. The breakdown of communication module (COMM) often leads to disfunction of cash dispenser module (CNG), printer (PRT) and internet data center (IDC). Besides, internet data center (IDC) problems influence cash dispenser module (CNG) and printer (PRT) weakly.

### 4.4 Social network analysis

In line with the previous works for information diffusion tracking [42, 14, 53], the public dataset MemeTracker^{4}^{4}4http://memetracker.org is used in this paper, which contains more than 172 million news articles or blog posts from various online media. The information, such as ideas, products and user behaviors, propagates over sites in a subtle way. For example, when Mark Zuckerberg posted ”I just killed a pig and a goat”, the meme appeared on the theguardian, Fortune, Business Insider one after the other and it became viral. This cascade can be regarded as a realization of an information diffusion over the network. By looking at the observed diffusion history, we want to know through which site and when a particular meme is likely to spread. Besides, we want to uncover the hidden diffusion structure from those meme cascades, which is useful in other social network analysis applications, e.g., marketing by information maximization [54]. From online articles, we collect more than 1 million meme cascades over different websites. For each meme cascade, we have the timestamp when sites mention a specific meme. For the experiments here, we use the top 500 media sites with the largest number of documents and select the meme cascades diffuse over them as done in previous works [42, 14]. As a result, we obtain around 31 million meme cascades, which are randomly split into training (50%), validating (%20) and testing (%30) parts. The event sequences are the meme cascades, which contain the website (dimension) and the timestamp. We count the times that a meme is mentioned during an hour over all the websites and use it as the time-series to reflect the hotness of the meme and the activity of websites.. The time interval of meme is shown in Fig. 9(a).

As the ground truth of network is unknown, we proceed by following the protocol as designed and adopted in [55, 56, 14]. We create a graph , for which each node is a website. If a post on site has a hyperlink pointed to site , then we create a directed edge with weight 1. If multiple hyperlinks exist between two sites, then the weight is accumulated. We use the graph as the ground truth and compare it with the inferred infectivity matrix from meme cascades. The prediction performance is evaluated by Accuracy@k, which evaluates whether the true label is within the top predicted dimensions.

model | accuracy@10 | accuracy@5 | MAE |
---|---|---|---|

Poisson | —– | —– | 1.63 |

SelfCorrecting | —– | —– | 1.70 |

Markov Chain | 0.563 | 0.472 | —– |

CTMC | 0.513 | 0.453 | 1.69 |

Logistic | 0.463 | 0.416 | 1.72 |

Hawkes | 0.623 | 0.563 | 1.68 |

RMTPP | 0.679 | 0.589 | 1.55 |

TRPP | 0.681 | 0.592 | 1.52 |

ERPP | 0.673 | 0.586 | 1.56 |

ATRPP | |||

AERPP | 0.678 | 0.589 | 1.45 |

The prediction performance is shown in Table IV from which one can observe our model outperforms the alternatives. The rank correlation is shown in Fig. 8. Our models ATRPP, AERPP can better uncover the infectivity matrix than the competitive methods in terms of the correlation rank metric.

In order to visualize the learned network, we use community detection algorithm [57] with resolution 0.9 [58] over learned directed network of ATRPP, which renders 18 communities. The resolution parameter controls the resolution of detected communities. It is set to lower values to get more communities (smaller ones), and is set higher to get fewer communities (bigger ones). Therefore, the communities can vary from the macroscale in which all nodes belong to the same community, to the microscale in which every node forms its own community. Fig. 7 shows some examples of those communities. Some media domains, like liverpooldailypost.com dominate in the cluster and what they publish usually spread to others and get viral.

### 4.5 Electronic health records mining

Uncovering the disease progression relation is important in healthcare analytics, which helps take preventive measures before fatal diseases happen. MIMIC-III (Medical Information Mart for Intensive Care III) is a large, publicly available dataset^{5}^{5}5https://mimic.physionet.org, which contains de-identified health-related data during 2001 to 2012 for more than 40,000 patients. It includes information such as demographics, vital sign measurements, diagnoses and procedures. For each visit, a patient receives multiple diagnoses, with one as the primary one. We filter out 937 patients and 75 diseases. The age, weight, heart rate and blood pressure are used as time series signals of patients. The distribution of time intervals between every two visits is shown in Fig. 9(b). We have used
the sequences of 600 patients to train, 100 to evaluate and the rest for test.

model | accuracy | MAE |
---|---|---|

Poisson | —— | 0.562 |

SelfCorrecting | —— | 0.579 |

Markov Chain | 77.53% | —— |

CTMC | 73.62% | 0.583 |

Logistic | 69.36% | 0.643 |

Hawkes | 78.37% | 0.517 |

RMTPP | 82.52% | 0.546 |

TRPP | 82.26% | 0.513 |

ERPP | 78.23% | 0.521 |

ATRPP | ||

AERPP | 83.96% | 0.503 |

Similar to MemeTracker, community detection algorithm [57] with resolution 0.9 is applied on the learned directed network from ATRPP. The results show cohesion within communities, which demonstrates the effectiveness of our attention mechanism. Note that some edges are too small to be visible. Specifically, Fig. 10(a) is about liver diseases. The node of alcohol cirrhosis has a large self-loop edge, which means this disease generally repeats many times. Besides, it has a large edge towards complications of transplanted kidney and hepatic encephalopathy, which means alcohol cirrhosis has a high probability of developing into the these two diseases. Fig. 10(b) is about respiratory diseases. Similarly, mechanical complication of tracheostomy and pseudomonal pneumonia have large self-loop edges, indicating they often relapse. Mechanical complication of tracheostomy has an edge towards pseudomonal pneumonia while the reverse, interestingly, does not exist. Fig. 10(c) shows the graph for alcohol-related diseases. Obsessed in alcohol has impact on reaction to indwelling urinary catheter regarding urinary system and hemorrhage of gastrointestinal tract regarding gastrointestinal system. Fig. 10(d) is about heart and blood related diseases. Three different parts of body form a progression line, consisting of diseases of trachea and bronchus (respiratory disease), paroxysmal ventricular tachycardia (heart rate disturbance) and acute vascular insufficiency of intestine (intestinal disease). The other line is septicemia leads to acute myocardial infarction of other anterior wall (heart attack), which results in paroxysmal ventricular tachycardia (heart rate disturbance). Fig. 10(e) is about blood metabolic diseases. Hemorrhage complicating a procedure leads to pulmonary embolism and infarction and complications due to renal dialysis device, implant, and graft. Diabetes with ketoacidosis results in urinary tract infection. For the observed strong association as stated above, we conjecture it might be due to either causal or correlation relationship, which can provide supporting evidence and implication for clinical staff and is subject to further analysis by health practitioners.

Table V reports the predictive performance of various models. ATRPP outperforms alternatives in both disease types and time prediction. Here first-order Markov Chain outperforms higher-order models, the reason might be due to the fact that visits of patients are sparse and there is not enough data to train the higher order Markov chains.

## 5 Conclusion

We conclude this paper with Fig. 11 and identify our proposed method as a recurrent point process model. To elaborate, Hawkes process uses a full explicit parametric model and RMTPP misses the dense time series features to model time-varying base intensity, assumes a partially parametric form for it and model the pseudo multi-dimensional point process. We make a further step by proposing an interpretable model which is simple and general and can be trained end-to-end. Most importantly, our model can uncover the subtle network structure and provide interpretable evidence for predicting result. The extensive experiments in this paper have clearly suggested its superior performance in synthetic and real-world data, even when we have no domain knowledge on the problem at hand. This is in contrast to existing point process models where an assumption about the dynamics is often needed to be specified beforehand.

## Acknowledgments

The authors thank Robert Chen with Emory University School of Medicine, and for helpful discussions and suggestions on the study of the computational experimental results on the MIMIC dataset. We are also thankful to Changsheng Li who shared us the ATM log data from IBM to allow us to perform the predictive maintenance study on real-world data.

## References

- [1] O. Aalen, O. Borgan, and H. Gjessing, Survival and event history analysis: a process point of view. Springer Science & Business Media, 2008.
- [2] D. C. Montgomery, C. L. Jennings, and M. Kulahci, Introduction to time series analysis and forecasting. John Wiley & Sons, 2015.
- [3] L. Li, H. Deng, A. Dong, Y. Chang, and H. Zha, “Identifying and labeling search tasks via query-based hawkes processes,” in KDD, 2014.
- [4] L. Yu, P. Cui, F. Wang, C. Song, and S. Yang, “From micro to macro: Uncovering and predicting information cascading process with behavioral dynamics,” 2015.
- [5] M. Farajtabar, Y. Wang, M. G. Rodriguez, S. Li, H. Zha, and L. Song, “Coevolve: A joint point process model for information diffusion and network co-evolution,” in Advances in Neural Information Processing Systems, 2015, pp. 1954–1962.
- [6] N. Du, H. Dai, R. Trivedi, U. Upadhyay, M. Gomez-Rodriguez, and L. Song, “Recurrent marked temporal point processes: Embedding event history to vectore,” in KDD, 2016.
- [7] M. Farajtabar, N. Du, M. G. Rodriguez, I. Valera, H. Zha, and L. Song, “Shaping social activity by incentivizing users,” in Advances in neural information processing systems, 2014, pp. 2474–2482.
- [8] G. E. Box, G. M. Jenkins, G. C. Reinsel, and G. M. Ljung, Time series analysis: forecasting and control. John Wiley & Sons, 2015.
- [9] C. Chatfield, The analysis of time series: an introduction. CRC press, 2016.
- [10] V. Guralnik and J. Srivastava, “Event detection from time series data,” in Proceedings of the fifth ACM SIGKDD international conference on Knowledge discovery and data mining. ACM, 1999, pp. 33–42.
- [11] E. Bacry, I. Mastromatteo, and J.-F. Muzy, “Hawkes processes in finance,” Market Microstructure and Liquidity, vol. 1, no. 01, p. 1550005, 2015.
- [12] E. Lewis and E. Mohler, “A nonparametric em algorithm for multiscale hawkes processes,” Journal of Nonparametric Statistics, 2011.
- [13] K. Zhou, H. Zha, and L. Song, “Learning triggering kernels for multi-dimensional hawkes processes.” in ICML (3), 2013, pp. 1301–1309.
- [14] ——, “Learning social infectivity in sparse low-rank networks using multi-dimensional hawkes processes.” in AISTATS, 2013.
- [15] H. Shen, D. Wang, C. Song, and A. Barabási, “Modeling and predicting popularity dynamics via reinforced poisson processes,” in AAAI, 2014.
- [16] S. Ertekin, C. Rudin, and T. H. McCormick, “Reactive point processes: A new approach to predicting power failures in underground electrical systems,” The Annals of Applied Statistics, vol. 9, no. 1, pp. 122–144, 2015.
- [17] E. Choi, N. Du, R. Chen, L. Song, and J. Sun, “Constructing disease network and temporal progression model via context-sensitive hawkes process,” in ICDM. IEEE, 2015.
- [18] H. Xu, W. Wu, S. Nemati, and H. Zha, “Patient flow prediction via discriminative learning of mutually-correcting processes,” TKDE, 2016.
- [19] A. G. Hawkes, “Spectra of some self-exciting and mutually exciting point processes,” Biometrika, 1971.
- [20] X. Shuai, J. Yan, X. Yang, H. Zha, and S. Chu, “Modeling the intensity function of point process via recurrent neural networks,” in AAAI, 2017.
- [21] T. J. Liniger, “Multivariate hawkes processes,” PhD thesis, Swiss Federal Institute Of Technology, Zurich, 2009.
- [22] J. Yan, Y. Wang, K. Zhou, J. Huang, C. H. Tian, H. Y. Zha, and W. S. Dong, “Towards effective prioritizing water pipe replacement and rehabilitation,” in IJCAI, 2013.
- [23] H. Xu, M. Farajtabar, and H. Zha, “Learning granger causality for hawkes processes,” in ICML, 2016.
- [24] J. L. Elman, “Finding structure in time,” Cognitive Science, vol. 14, pp. 179–211, 1990.
- [25] R. Pascanu, T. Mikolov, and Y. Bengio, “On the difficulty of training recurrent neural networks,” in ICML, 2013.
- [26] S. Hochreiter and J. Schmidhuber, “Long short-term memory,” Neural computation, vol. 9, no. 8, pp. 1735–1780, 1997.
- [27] A. Graves, “Generating sequences with recurrent neural networks,” in arXiv:1308.0850, 2013.
- [28] J. Chung, C. Gulcehre, K. Cho, and Y. Bengio, “Empirical evaluation of gated recurrent neural networks on sequence modeling,” in arXiv:1412.3555, 2014.
- [29] K. Gregor, I. Danihelka, A. Graves, D. Rezende, and D. Wierstra, “Draw: A recurrent neural network for image generation,” in ICML, 2015.
- [30] A. Graves, A. rahman Mohamed, and G. Hinton, “Towards end-to-end speech recognition with recurrent neural networks,” in ICML, 2014.
- [31] I. Sutskever, O. Vinyals, and Q. V. Le., “Sequence to sequence learning with neural networks,” in NIPS, 2014.
- [32] J. T. Connor, R. D. Martin, and L. E. Atlas, “Recurrent neural networks and robust time series prediction,” IEEE Transactions on Neural Networks, vol. 5, no. 2, pp. 240–254, 1994.
- [33] R. Chandra and M. Zhang, “Cooperative coevolution of elman recurrent neural networks for chaotic time series prediction,” Neurocomputing, vol. 86, pp. 116–123, 2012.
- [34] A. Jain, A. Singh, H. S. Koppula, S. Soh, and A. Saxena, “Recurrent neural networks for driver activity anticipation via sensory-fusion architecture,” in ICRA, 2016.
- [35] E. Choi, M. T. Bahadori, J. Sun, J. Kulas, A. Schuetz, and W. Stewart, “Retain: An interpretable predictive model for healthcare using reverse time attention mechanism,” in NIPS, 2016.
- [36] K. Xu, J. Ba, R. Kiros, K. Cho, and A. Courville, “Show, attend and tell: Neural image caption generation with visual attention,” in ICML, 2015.
- [37] K. Cho, A. Courville, and Y. Bengio, “Describing multimedia content using attention-based encoder-decoder networks,” IEEE Transactions on Multimedia, vol. 17, no. 11, pp. 1875–1886, 2015.
- [38] J. F. C. Kingman, Poisson processes. volume 3. Oxford university press, 1992.
- [39] R. Pemantle, “A survey of random processes with reinforcement,” Probability Survey, vol. 4, no. 0, pp. 1–79, 2007.
- [40] L. Tran, M. Farajtabar, L. Song, and H. Zha, “Netcodec: Community detection from individual activities,” in Proceedings of the 2015 SIAM International Conference on Data Mining. SIAM, 2015, pp. 91–99.
- [41] E. Lewis, G. Mohler, P. J. Brantingham, and A. Bertozzi, “Self-exciting point process models of insurgency in iraq,” UCLA CAM Reports 10, vol. 38, 2010.
- [42] M. Farajtabar, X. Ye, S. Harati, L. Song, and H. Zha, “Multistage campaigning in social networks,” in Advances in Neural Information Processing Systems, 2016, pp. 4718–4726.
- [43] S. A. Hosseini, K. Alizadeh, A. Khodadadi, A. Arabzadeh, M. Farajtabar, H. Zha, and H. R. Rabiee, “Recurrent poisson factorization for temporal recommendation,” arXiv preprint arXiv:1703.01442, 2017.
- [44] B. Tabibian, I. Valera, M. Farajtabar, L. Song, B. Schölkopf, and M. Gomez-Rodriguez, “Distilling information reliability and source trustworthiness from digital traces,” arXiv preprint arXiv:1610.07472, 2016.
- [45] Y. Ogata, “Statistical models for earthquake occurrences and residual analysis for point processes,” J. Amer. Statist. Assoc., vol. 83, no. 401, pp. 9–27, 1988.
- [46] V. Isham and M. Westcott, “A self-correcting pint process,” Advances in Applied Probability, vol. 37, pp. 629–646, 1979.
- [47] R. Pascanu, T. Mikolov, and Y. Bengio, “On the difficulty of training recurrent neural networks.” ICML (3), vol. 28, pp. 1310–1318, 2013.
- [48] A. G. Hawkes and D. Oakes, “A cluster process representation of a self-exciting process,” Journal of Applied Probability, 1974.
- [49] A. Rosenberg, “Classifying skewed data: Importance weighting to optimize average recall.” in INTERSPEECH, 2012, pp. 2242–2245.
- [50] Y. N. Dauphin, H. de Vries, J. Chung, and Y. Bengio, “Rmsprop and equilibrated adaptive learning rates for non-convex optimization,” in arXiv:1502.04390, 2015.
- [51] J. Yan, Y. Wang, K. Zhou, J. Huang, C. H. Tian, H. Y. Zha, and W. S. Dong, “Towards effective prioritizing water pipe replacement and rehabilitation,” in IJCAI, 2013.
- [52] Y. Ogata, “On lewis’ simulation method for point processes,” IEEE Transactions on Information Theory, vol. 27, no. 1, pp. 23–31, 1981.
- [53] M. G. Rodriguez, J. Leskovec, D. Balduzzi, and B. Schölkopf, “Uncovering the structure and temporal dynamics of information propagation,” Network Science, vol. 2, no. 01, pp. 26–65, 2014.
- [54] H. Zhuang, Y. Sun, J. Tang, J. Zhang, and X. Sun, “Influence maximization in dynamic social networks,” 2013.
- [55] M. Gomez Rodriguez, J. Leskovec, and A. Krause, “Inferring networks of diffusion and influence,” in Proceedings of the 16th ACM SIGKDD international conference on Knowledge discovery and data mining. ACM, 2010, pp. 1019–1028.
- [56] M. Gomez Rodriguez, D. Balduzzi, B. Schölkopf, G. T. Scheffer et al., “Uncovering the temporal dynamics of diffusion networks,” in 28th International Conference on Machine Learning (ICML 2011). International Machine Learning Society, 2011, pp. 561–568.
- [57] V. D. Blondel, J.-L. Guillaume, R. Lambiotte, and E. Lefebvre, “Fast unfolding of communities in large networks,” Journal of statistical mechanics: theory and experiment, vol. 2008, no. 10, p. P10008, 2008.
- [58] R. Lambiotte, J.-C. Delvenne, and M. Barahona, “Laplacian dynamics and multiscale modular structure in networks,” arXiv preprint arXiv:0812.1770, 2008.