1. Scalable pattern counting in temporal network data
Pattern counting is one of the fundamental problems in data mining (Chakrabarti2006curriculum, ; Han2011book, )
. A particularly important case is counting patterns in graph data, which is used within a variety of network analysis tasks such as anomaly detection
(Noble2003anomaly, ; Sun2007GraphScope, ), role discovery (Henderson2012RolX, ; Rossi2015role, ), and clustering (Rohe2013blessing, ; Benson2016hoo, ; Tsourakakis2017scalable, ). These methods typically make use of features derived from the frequencies of small graph patterns—usually called motifs (Milo2002motifs, ) or graphlets (Przulj2004graphlet, ) (we adopt the “motif” terminology in this paper)—and are used across a range of disciplines, including social network analysis (Leskovec2010signed, ; Ugander2013subgraphs, ), neuroscience (Hu2012motifs, ; Battiston2017motif, ), and computational biology (Mangan2003ffl, ; Przulj2007comparison, ). Furthermore, the counts of motifs have also been used to automatically uncover fundamental design principles in complex systems (Milo2002motifs, ; Mangan2003ffl, ; Milo2004superfamilies, ).The scale of graph datasets has led to a number of algorithms for estimating the frequency of motif counts (Ahmed2014graphsample, ; Elenberg2016profiles, ; Bressan2017graphlets, ; Jain2017cliques, ; Wang2018moss, ). For example, just the task of estimating the number of triangles in a graph has garnered a substantial amount of attention (Tsourakakis2009doulion, ; Avron2010counting, ; Lim2015mascot, ; Seshadhri2013wedge, ; Eden2017triangles, ; Stefani2017triest, ). Many of these algorithms are based on sampling procedures amenable to streaming models of graph data (Feigenbaum2005streaming, ; McGregor2014streaming, ). At this point, there is a reasonably mature set of algorithmic and statistical tools available for approximately counting motifs in large graph datasets.
While graphs have become large enough to warrant frequency estimation algorithms, graph datasets have, at the same time, become richer in structure. A particularly important type of rich information is time (Kossinets2008pathways, ; Holme2012temporal, ; Farajtabar2015coevolve, ; Gaumont2016dense, ; Scholtes2017networks, ). Specifically, in this paper, we consider datasets where edges are accompanied by a timestamp, such as the time a transaction was made with a cryptocurrency, the time an email was sent between colleagues, or the time a packet was forwarded from one IP address to another by a router. Accordingly, motifs have been generalized to incorporate temporal information (Zhao2010communicationmotifs, ; Kovanen2011motifs, ; Paranjape2017motifs, ) and have already been used in a variety of applications (Lahiri2007structure, ; Shao2013temporal, ; Meydan2013prediction, ; Kovanen2013motifs, ). However, we do not yet have algorithmic tools for estimating frequencies of temporal motifs in these large temporal graphs. This is especially problematic since including timestamps increases the size of the stored data; for example, a traditional email graph would only record if one person has ever emailed another person, whereas the temporal version of the same network would record every time there is a communication.
To exacerbate the problem, counting temporal motifs turns out to be fundamentally more difficult in a computational complexity sense. In particular, we prove that counting basic temporal star motifs is NPcomplete. This contrasts sharply with stars in traditional static graphs, which are generally considered trivial to count (the number of noninduced edge stars with center node is simply , where is the degree of ). Thus, our result highlights how counting problems in temporal graphs involve fundamentally more challenging computations, thus further motivating the need for approximation algorithms.
Here we develop the first frequency estimation algorithms for counting temporal motifs. We focus on the definition of temporal motifs from Paranjape et al. (Paranjape2017motifs, ), but our methodology is general and could be adapted for other definitions. Our approach is based on sampling that employs as a subroutine any algorithm (satisfying some mild conditions) that exactly counts the number of instances of temporal motifs. Thus, our methodology provides a way to accelerate existing algorithms (Paranjape2017motifs, ; Mackey2018chronological, ), as well as better exact counting algorithms that could be developed in the future.
At a basic level, our sampling framework partitions time into intervals, uses some algorithm to find exact motif counts in a subset of the intervals, and weights these counts to get an estimate of the number of temporal motifs. A key challenge is that the time duration of a temporal motifs can cross interval boundaries, which makes it challenging to obtain an accurate frequency estimator since motifs of larger duration are more likely to be omitted. At its core, our sampling framework uses importance sampling (mcbook, )
in two different ways. First, we use importance sampling as a way to design an unbiased estimator by appropriately scaling the exact counts appearing in some intervals. Second, we use importance sampling as a way to (probabilistically) choose which intervals to sample, which reduces the variance of our unbiased estimator.
In addition to the scalability advantages offered by sampling, our framework has two other important features. First, the sampling requires a smaller amount of memory. We show an example where this enables us to count a complex motif on a large temporal graph when existing exact counting algorithms run out memory. Second, the sampling procedure has builtin opportunity for parallel computation, which provides a path to faster computation with exact counting algorithms that do not have builtin parallelism.
As discussed above, our sampling framework employs an exact counting algorithm as a subroutine. The constraints on the algorithm are that it must provide the exact counts along with the socalled duration of the motif instance (the difference in the earliest and latest timestamp in the edges in the motif instance; for example, the duration in the top left motif instance in Figure 1 is 32  16 = 16). This constraint holds for some existing algorithms (Mackey2018chronological, ) but not for others (Paranjape2017motifs, ). An additional contribution of our work is a new exact counting algorithm for a class of star motifs that is compatible with our sampling framework. As an added bonus, this new exact counting algorithm actually outperforms existing algorithms.
We test our sampling procedure on several temporal graph datasets from a variety of domains, ranging in size from 60,000 to over 600 million temporal edges and find that our sampling framework can improve run time by one to two orders of magnitude while maintaining a relative error tolerance of 5% in the counts. The variance analysis of our error bounds tends to be pessimistic, since we make no assumptions on the distribution of timestamps within our datasets. Thus, we also show empirically that our worstcase bounds are far from what we see in the data.
2. Preliminaries on temporal motifs
We first review some basic notions of temporal motifs. There are a few types of temporal motifs, which we discuss in the context of related work in Section 6. Here we review the definitions used by Paranjape et al. (Paranjape2017motifs, ), which is one of the more flexible definitions that also poses difficult computational challenges.
A temporal edge
is a timestamped directed edge between an ordered pair of nodes. A collection of temporal edges is a
temporal graph (see Figure 1). Formally, a temporal graph on a node set is a collection of tuples , , where each and are elements of and each is a timestamp in . There can be many temporal edges from to (one example is in email data, where one person sends an email to another many times). We assume that the timestamps are unique so that the temporal edges in a graph can be ordered. This assumption makes the presentation of the paper simpler, but our methods can handle temporal graphs with nonunique timestamps.If we ignore time and duplicate edges, the temporal graph induces a standard (static) directed graph. Formally, the static graph of a temporal graph on a node set is a graph , where . Edges in are called static edges.
Next, we formalize temporal motifs (illustrated in Figure 1).
Definition 2.1 (Temporal motif (Paranjape2017motifs, )).
A node, edge temporal motif consists of a multigraph with nodes and edges and an ordering on the edges of .
We often find it convenient to represent by an ordered sequence of edges . Definition 2.1 is a template for a temporal graph pattern, and we want to count how many times the pattern appears in a temporal network. Furthermore, we are interested in how often the motif occurs within some time span . A collection of edges in a temporal graph is a instance of a temporal motif if it matches the same edge pattern of the multigraph , the temporal edges occur in the specified order , and all of the temporal edges occur within a time window (see Figure 1). We now formalize this definition.
Definition 2.2 (Motif instance (Paranjape2017motifs, )).
A timeordered sequence , , of unique temporal edges from a temporal graph is a instance of the temporal motif if

There exists a bijection on the vertices in such that and , ; and

The edges all occur within the time span, i.e., .
With this definition, motif instances are defined by just the existence of edges (a general subgraph) and not the nonexistence of edges (an induced subgraph).
We are interested in counting how many motifs appear within a maximum time span of time units. Our sampling framework will also make use of the actual duration of motif instances, or the difference in the latest and earliest timestamp of a motif instance. We formalize this notion in the following definition.
Definition 2.3 (Motif duration).
Let , , be an instance of a motif as per Definition 2.2 with . Then the duration of the instance, denoted , is .
3. Counting temporal stars is hard
Star motifs as in Figure 1 are one of the fundamental small graph patterns and are used in, e.g., anomaly detection (Akoglu2010oddball, ) and graph summarization (Koutra2015summarizing, ). In static graphs, counting noninduced instances of stars is simple. Given the degree of node , is the center of node stars. Thus, there is a simple polynomialtime algorithm for computing the total number of stars.
In contrast, once we introduce temporal information, it turns out that stars become hard to compute. Specifically, we show in this section that counting temporal stars is NPcomplete, and even determining the existence of a temporal star motif is NPcomplete. This result serves two purposes. First, it highlights that the computational challenges with temporal graph data are fundamentally different from those in traditional static graph analysis. Second, the computational difficulty in such a simple type of temporal motif motivates the need for scalable approximation algorithms, which we develop in the next section. We begin with a formal definition of a temporal star motif.
Definition 3.1 ().
A temporal star is a temporal motif where the multigraph is connected and has node set with edges , , where either or is 0, .
The restriction that either or is 0 means that each edge either originates from node 0 or enters node 0. The ordering of the edges in the multigraph needed by Definition 2.1 is arbitrary—we only need the star structure of the multigraph. We will show that determining the existence of an instance of a temporal star in a temporal graph is NPcomplete and then generalize our result to an even more restricted class of star motifs. We begin with the formal decision problem.
Problem 1 ().
Given a temporal graph , a temporal star , and a time span , the kStarMotif problem asks if there exists at least one instance of in .
To establish NPcompleteness, we reduce kClique to kStarMotif. A kClique problem instance is formalized as follows: given an undirected graph and an integer , the kClique problem asks if there exists at least one clique of size in .
Theorem 3.2 ().
kStarMotif is NPcomplete.
Proof.
Our input is an instance of kClique on a vertex set . Assume that the nodes in are numbered from to (Figure 2). We construct an instance of kStarMotif:
The timestamps come from the set , and we think of the timestamps as partitioned into blocks, with timestamps in each block. If the timestamp of an edge lies in , then we say that the edge belongs to block . Each block then corresponds to a node in , with the first and last timestamp in each block reserved for the backward edges we add to . For each node in the original graph, we add the two backward edges in block to node in , and for each neighbor of , we add a forward edge using the timestamp in the th position of block . Figure 2 is a schematic of the construction. Observe that if there is a clique in , then by construction the star motif occurs in .
Intuitively, the backward edges added to and serve as “bookends”. If the two backward edges corresponding to a node are found to be part of , then each of the other nodes in the motif has to contribute a forward edge with timestamps between the two backward edges of . By construction of , an edge connected to can only have a timestamp in block if is connected to in . This implies that is connected to the other nodes selected in the motif . Applying this argument to each node in the motif , there must be a clique in the original graph . ∎
The result does not depend on having edges in two directions. We call a star motif unidirectional if all of the edges in either originate from or enter the center node (node 0, in our notation).
Theorem 3.3 ().
kStarMotif is NPcomplete even when restricted to unidirectional stars.
Proof.
(Sketch.) Instead of using two backward edges for bookkeeping, we can expand the size of each block to and use the first and last timestamps within the block as the bookends. Thus, the graph in the previous proof is modified by connecting the center node to each node with forward edges using the timestamps reserved for bookkeeping in block . The motif is modified by requiring the same forward edges as before, plus an additional forward edges, with timestamps in and for . By using forward edges for each bookend, we ensure that any occurrence found in must include at least one edge from each bookend of the chosen nodes. This allows us to again argue that forward edges must be between the bookends of block , implying that there is a clique. ∎
These hardness results illustrate the computational difficulties in counting temporal graph patterns, which motivates scalable approximation algorithms for counting such patterns. We next present a general sampling framework for scalable estimation of the number of instances of temporal motifs.
4. Algorithmic sampling framework
Suppose we are given a motif , a time span , and a temporal graph . In this section, we develop a sampling framework to estimate the number of instances of in , which we denote by . Our sampling framework will employ some algorithm that can compute exactly the number of instances of on temporal subgraphs of . The requirements on the algorithm are that, given a temporal graph , a motif , and a time span , the algorithm outputs a sequence of the countduration pairs , where is the number of instances of the motif with duration . We denote this output by . We work from these assumptions in this section, and Section 5 discusses compatible algorithms.
Intervals and the count vector
. We begin with some definitions and technical lemmas that will later be used to develop our estimator. Let be a random integer uniformly drawn from for some input integer that controls the size of the sampling windows. We call a shift, and we will eventually make use of multiple shifts within our sampling framework. We consider the set of intervals of width with shift :(1) 
For an instance of the motif with duration
, it is easy to see that the probability (over a random choice of shift
) that is completely contained within an interval in is(2) 
Next, for an interval , let
be an indicator random variable which equals 1 if
is completely contained in and 0 otherwise. For each interval , we associate a weighted count of the number of instances of motif completely contained in the interval :(3) 
Let be the vector of such counts:
(4) 
(here, denotes the th coordinate of ). Next, let be an indicator random variable that equals 1 if the motif instance is completely contained in an interval in and 0 otherwise. Then . The following lemma says that is an unbiased estimator for the motif count for any value of .
Lemma 4.1 ().
.
Proof.
Since , . ∎
The next lemma bounds the variance of .
Lemma 4.2 ().
.
Proof.
First, we have that
where and range over the instances of the motif . Using the bounds and ,
Putting everything together,
∎
Our sampling framework estimates in order to estimate the number of motif instances . The basic idea of our approach is to use importance sampling to speed up this estimation task, by picking a set of intervals in and computing their weights. Here, computing the weight for an interval uses an exact motif count restricted to the interval . Equivalently, we (i) sample a subset of coordinates of , (ii) compute their values exactly, and (iii) combine them to estimate . We describe this procedure next.
Importance sampling for an estimator. Let denote the th coordinate of , corresponding to interval . Our estimator is a random variable defined as follows. First, we sample interval (independently) with some probability . These values will be based on simple statistics of the intervals; we will specify choices for later but note for now that they do not necessarily sum to 1. Second, let be an indicator random variable corresponding to interval , where equals 1 if is picked and 0 otherwise. Finally, our estimator is
(5) 
Our first result is that is an unbiased estimator for , the number of instances of the motif .
Theorem 4.3 ().
The random variable in Eq. 5 is an unbiased estimator for the number of motif instances, i.e., .
Proof.
First, note that . For any , . Hence, by Lemma 4.1. ∎
Next, we work to bound the variance of our estimator . To this end, it will be useful to define a scaled version of :
(6) 
The following lemma provides a useful equality on the variance of our estimator in terms of and , conditioned on the shift .
Lemma 4.4 ().
.
Proof.
By independence of the ,
Therefore, . ∎
We are now ready to bound the variance of .
Theorem 4.5 ().
.
Proof.
Our analysis thus far has been for a single shift . If we repeat the above computations for randomly chosen shifts and report the mean of the estimates, then the variance is reduced by a factor of . Algorithm 1 outlines the the overall sampling procedure, assuming that the sampling probabilities are given along with the exact counting algorithm . In the algorithm, we use to denote the subgraph restricted to interval and to denote the output of the exact counting algorithm on the interval, which is a sequence of the countduration pairs of motif instances contained in the interval. The algorithm also explicitly states that the parallelism that can be performed over the samples.
Choosing the sampling probabilities. In order to get average squared error , we need to set the parameters as follows:
(7)  
(8) 
The first term in the lefthand side of Eq. 8 combines (i) a natural measure of sparsity of the distribution of motifs with (ii) the extent of correlation between the sampling probabilities and the (weighted) motif counts for intervals . In order to understand this, let denote the dimension of and consider the simple uniform setting of (so one interval is sampled in expectation). In this case, the term becomes
(9) 
Equation 9 is a natural measure of sparsity of the vector . In the extreme case where is a vector with only one nonzero coordinate, the value is , and in the other extreme where is a uniform vector, the value is bounded above by 1. In the sparse case, we need to increase the sampling probabilities—thus sampling more intervals—to compensate for the large variance. In the worst case, this would require looking at all the intervals, i.e., we get no running time savings from sampling (however, we will see in our experiments that the data is far from the worst case in practice). Nonetheless, the ability of the algorithm to pick sampling probabilities gives flexibility to mitigate the dependence on the sparsity of . To illustrate this point, in the extremely favorable case when is proportional to , i.e., (so one interval is sampled in expectation), the first term on the lefthand side of Eq. 8 is less than 1. This analysis suggests that a good choice of sampling probabilities roughly balances the two terms:
A priori, we do not know or . What we can do is choose the sampling probabilities by some easily measured statistic that we think is correlated with , such as the number of temporal edges or number of static edges. For this paper, we simply choose to be proportional to the number of temporal edges in the interval, i.e.,
(10) 
where is a small constant (in practice on the order of 10–100). This leads to substantial speedups, as we will see in the next section. There are certainly more sophisticated approaches one could take to choose the , and we leave this as an avenue for future research.
Streaming from sampling. When memory is at a premium, the sampling framework above can be made memory efficient. By considering the windows in chronological order, the edges of past windows do not need to be stored. By running several estimators in parallel, we can achieve any accuracy we want while only needing to store edges in an interval of at most at a time. As we will see in our experiments, the memory savings allows us to processes larger temporal graphs than we could with an exact algorithm.
5. Computational experiments
In this section, we use our sampling framework from Section 4 and various exact temporal motif counting algorithms to count temporal motifs on realworld datasets. By exploiting sampling and the ability sample in parallel, we obtain substantial speedups with modest computational resources and only a small error in the estimation.
Data. We gathered 10 datasets for our experiments. Paranjape et al. analyzed seven of them (Paranjape2017motifs, ), and we collected three larger datasets to better analyze the performance of our methodology. Table 1 lists summary statistics of the datasets, and we briefly describe them below. Each dataset is a collection of timestamped directed edges. The time resolution of each dataset is 1 second, except for the EquinixChicago dataset, where the time resolution is 1 microsecond.
dataset  # nodes  # static  # temporal  time 

edges  edges  span  
CollegeMsg  1.9K  20.3K  59.8K  194 days 
emailEucore  986  24.9K  332K  2.20 years 
MathOverflow  24.8K  228K  390K  6.44 years 
AskUbuntu  157K  545K  727K  7.16 years 
SuperUser  192K  854K  1.11M  7.60 years 
WikiTalk  1.09M  3.13M  6.10M  6.24 years 
StackOverflow  2.58M  34.9M  47.9M  7.60 years 
Bitcoin  48.1M  86.8M  113M  7.08 years 
EquinixChicago  12.9M  17.0M  345M  62.0 mins 
RedditComments  8.40M  517M  636M  10.1 years 
CollegeMsg (Panzarasa2009patterns, ).
A network of private messages sent on an online social network at the University
of California, Irvine.
emailEucore (Paranjape2017motifs, ).
A collection of internal email records from a European research institution.
MathOverflow, AskUbuntu, SuperUser, and StackOverflow (Paranjape2017motifs, ).
These datasets are derived from user interactions on Stack Exchange question and
answer forums. A temporal edge represents a user replying to a question,
replying to a comment, or commenting on a question.
WikiTalk (Leskovec2010Governance, ; Paranjape2017motifs, ).
A network of Wikipedia users making edits on each others’ “talk pages.”
Bitcoin (Kondor2014inferring, )
. A network representing timestamped transactions on Bitcoin. The addresses were partially aggregated by a deidentification heuristic
(Reid2012analysis, ) implemented by Kondor et al., using all transactions up to February 9, 2016 (Kondor2014inferring, ). Timestamps are the creation time of the block on the blockchain containing the transaction. We will release this dataset with the paper.EquinixChicago (CAIDAdata, ). This dataset was constructed from passive internet traffic traces from CAIDA’s monitor in Chicago on February 17, 2011. Each edge represents a packet sent from one anonymized IP address to another. Data was collected from the “A direction” of the monitor.
RedditComments (Hessel2016science, ). This dataset was constructed from a large collection of comments made by users on https://www.reddit.com, a popular social media platform. A comment from user to user at time induces a temporal edge in our dataset.
5.1. Exact counting algorithms
Our sampling framework is flexible since it can use any algorithm that exactly counts temporal motifs as a subroutine, provided that this algorithm can be transformed to output the countduration pairs , where is the number of instances of the motif with duration . A recently proposed “backtracking” algorithm satisfies this constraint (Mackey2018chronological, ). The fast algorithms for 2node, 3edge star motifs introduced by Paranjape et al. do not satisfy these requirements, since the algorithm uses an inclusionexclusion rule that cannot output the durations. However, we still use this algorithm as a baseline in our experiments. We also create a new exact counting algorithm that is compatible with our sampling framework, which we describe below.
Backtracking algorithm (BT, (Mackey2018chronological, )). The backtracking algorithm examines the edges of the input graph in chronlogical order and matches one edge of the motif at a time. The software was not released publicly, so we have reimplemented it with some optimizations. The algorithm is compatible with our sampling framework. The algorithm is inherently sequential, so our parallel sampling framework is especially useful with this method.
Fast 2node, 3edge algorithm (F23, (Paranjape2017motifs, )). Paranjape et al. introduced a collection of algorithms for counting motifs with at most 3 edges. Here we use their specialized algorithm for motifs with two nodes and three edges, which produces exact counts in time linear in the number of edges. This algorithm is incompatible with our sampling framework, but we still use it for comparative purposes.
dataset  BT  BT+S  BT+PS  EX23  EX23+S  EX23+PS  F23  F23+P  error (%) 

CollegeMsg  0.076  0.072  0.038  0.017  0.016  0.009  0.056  0.054  0.33 
EmailEu  0.339  0.305  0.191  0.073  0.078  0.027  0.217  0.165  1.44 
MathOverflow  0.545  0.361  0.143  0.233  0.148  0.097  0.998  0.878  1.74 
AskUbuntu  1.414  1.305  0.500  0.592  0.311  0.176  2.534  2.371  1.84 
SuperUser  2.590  1.446  0.483  1.097  0.194  0.104  4.595  4.129  1.69 
WikiTalk  15.92  14.88  5.463  4.737  3.645  0.876  20.46  18.23  0.89 
StackOverflow  198.9  160.8  79.58  108.1  69.50  17.81  299.2  230.1  1.95 
Bitcoin  514.0  520.7  102.3  494.4  233.5  88.66  10348  10135  3.59 
EquinixChicago  480.4  180.3  37.64  382.7  56.33  24.64  477.3  383.8  0.00 
RedditComments  7301  7433  2910  1563  3154  367.4  6602  5036  4.83 
A new exact algorithm for 2node 3edge motifs (EX23). We devised a new algorithm for 2node, 3edge motifs that is compatible with our sampling framework. In this case, each pair of nodes in the input graph forms an independent counting problem. For each pair of nodes in the input that are neighbours in the static graph, we gather all temporal edges between the two nodes. Then we fix the first and last edge of the 3edge motif by iterating over all pairs of gathered temporal edges. By maintaining an additional counter of the number of edges between the two fixed edges, we can count the number of temporal motifs that begins on the first fixed edge and ends on the second fixed edge (the procedure is outlined in Algorithm 2). Overall, this procedure takes time, where the sum iterates over all pairs of nodes in the graph, and is the number of temporal edges between nodes and . With additional code complexity from special tree structures, the running time can be improved to and still be compatible with our sampling framework. However, this optimization is not crucial for the main focus of our paper, which is the acceleration of counting algorithms with sampling. Thus, we use the simpler unoptimized algorithm, which we will see actually outperforms the other exact counting algorithms.
5.2. Performance results
We now evaluate the performance of several algorithms: (i) the three baseline exact counting algorithms described in the previous section (BT, F23, EX23); (ii) the F23 baseline with parallelism enabled (F23+P); (iii) our sampling framework on top of backtracking and our new exact counting algorithm (BT+S, EX23+S); and (iv) our parallelized sampling framework on top of backtracking and our new exact counting algorithm (BT+PS, EX23+PS). As explained above, the F23 algorithm is incompatible with our sampling framework; we include the algorithm and its parallelized version as baselines for fast exact counting.
All algorithms were implemented in C++, and all experiments were executed on a 16core 2.20 GHz Intel Xeon CPU with 128 GB of RAM. The algorithms ran on a single thread unless explicitly stated to be parallel. The parallel algorithms used 16 threads. In the case of the sampling algorithm, parameters are set so that the approximations are within 5% relative error of the true value.
Experiments on a 2node, 3edge motif. Table 2 reports the performance of all algorithms on the 2node, 3edge temporal motif in Figure 3 (we chose this motif to allow us to compare against one of the fast algorithms of Paranjape et al. (Paranjape2017motifs, )). The time span was set to 86400 seconds = 1 day in all datasets except EquinixChicago, where was 86400 microseconds (these are the same parameters used in exploratory data analysis in prior work (Mackey2018chronological, )).
We highlight three important findings. First, our new EX23 algorithm with parallel sampling is the fastest algorithm on every dataset. Comparing our algorithm against the previous state of the art, we see speedups up to 120 times faster than the slowest exact algorithm (see the results for Bitcoin). This is in part due to the fact that our EX23 algorithm is actually faster the than the backtracking algorithm (BT) and the fast algorithm of Paranjape et al. (F23). In other words, our proposed exact algorithm already outperforms the current state of the art.
Second, in all cases, parallel sampling provides a substantial speedup over the exact baseline algorithm. Speedups are typically on the order of 2–6x improvements in running time. We used 16 threads but did not optimize our parallel algorithms; there is ample room to improve these results with additional software effort.
Third, the running time of the backtracking algorithm with sampling (BT+S) is often comparable to simple backtracking (BT). In these cases, we hypothesize that the backtracking algorithm has enough overhead and is pruning enough edges to make simple sampling not worthwhile under our parameter settings. However, parallel sampling with the backtracking algorithm (BT+PS) can yield substantial speedups (see, e.g., Bitcoin, EquinixChicago, and RedditComments). This illuminates an important feature of our sampling framework, namely, we get parallelism for free. The backtracking algorithm is inherently sequential, but parallel sampling can achieve substantial speedups. Thus, future research in the design of fast exact counting algorithms can largely leave parallelism to be handled by our sampling framework. Finally, although not reported, the sampling framework requires a smaller amount of memory than the exact algorithms; thus, if no parallelism is available, we can at least gain in terms of memory, if not in speed.
We used the heuristic in Eq. 10 to determine the sampling probabilities in these experiments. To understand why this heuristic worked for these datasets (i.e., the relative errors are small), we measured the correlation of and the coordinates of the vector used in the sampling framework (Figure 4). In datasets such as CollegeMsg and WikiTalk, the correlation between is large, and consequently, the relative error in the estimates is small (Table 2).
Experiments on a 4node, 4edge bifan motif. Next, we show the results of the backtracking algorithm on a socalled “bifan” motif (Figure 3). This motif has four nodes and four temporal edges. The static version of the motif appears is important in networks from a variety of domains, including sociology (Zhang2013potential, ), neuroscience (Benson2016hoo, ), gene regulation (Dobrin2004aggregation, ), and circuit design (Milo2002motifs, ). Since this motif has four nodes and four edges, our new fast algorithm and the fast algorithm of Paranjape et al. cannot be used. Thus, we focus on accelerating the backtracking algorithm with (parallelized) sampling. For these experiments, we set , as running the algorithm with exceeds the alloted memory on the Bitcoin and RedditComments dataset. Table 3 shows the results.
dataset  BT  BT+S  BT+PS  error (%) 

CollegeMsg  0.081  0.076  0.069  1.02 
EmailEu  0.353  0.307  0.120  0.85 
MathOverflow  0.528  0.362  0.041  3.60 
AskUbuntu  1.408  0.909  0.078  4.52 
SuperUser  2.486  1.269  0.164  2.33 
WikiTalk  51.85  35.35  11.99  2.01 
StackOverflow  221.7  93.10  5.208  4.88 
Bitcoin  1175  985.9  269.3  3.09 
EquinixChicago  481.2  45.50  5.666  1.33 
RedditComments  ✗  6739  2262  – 
Again, the parallel sampling procedure provides a substantial speedup over the baseline algorithm. We emphasize that our simple parallelization technique is a property of the sampling procedure and not a property of the exact algorithm. In fact, the exact algorithm is inherently sequential, and the sampling framework enables parallelism in a trivial way with minimal loss in accuracy. Moreover, since the sampling algorithm only examines a portion of the graph at a time, it uses much less memory than the exact counting algorithm. For example, with the RedditComments dataset, the exact algorithm ran out of memory, while the sampling algorithms completed successfully (thus no relative error is reported). This feature is useful in streaming applications, where memory is limited.
6. Additional related work
Our sampling framework relies on importance sampling, which is used for finding motifs in gene sequence analysis (siddharthan2008phylogibbs, ; chan2010importance, ; gupta2007variable, ; liang2008sequential, ); here “motif” refers to short string patterns in DNA. (The term “motif” in the context of network analysis is borrowed from this domain (shen2002network, ).) We have already covered much of the related work in sampling algorithms for pattern counting in static graphs, so we summarize additional research related to various definitions of temporal motifs here. Some of these are for sequences of static snapshot graphs (Zhang2014dynamic, ; Lahiri2007structure, ; Jin2007trend, ), which is a different data model than the one in this paper, where edges have timestamps in a continuum. For the data in our paper, there are motifs based on “adjacent events” that require each new edge in a sequence to be within a certain timespan of each other (Zhao2010communicationmotifs, ; Gurukar2015commit, ; Hulovatyy2015exploring, ). These definitions are slightly more restrictive than the one by Paranjape et al. analyzed here; however, the principles of our techniques could also be adapted to these cases as well. Kovanen et al. use the same notion of event adjacency but also restrict motif instances to cases where the events are consecutive for all nodes involved (i.e., within the span of the motif instance, there can be no other temporal edge adjacent to one of the nodes) (Kovanen2011motifs, ). This definition is even more restrictive in the events that it captures but it does allows for much faster exact counting algorithms, e.g., triangle motifs can be counted in linear time in the size of the data. Thus, speeding up computation with sampling is less appealing for this definition. Finally, there is also a line of research in finding dense subgraphs in datasets similar to our model, which is a specific type of motif (Gaumont2016dense, ; Viard2018enumerating, ; Viard2016computing, ).
7. Discussion
We have developed a sampling framework for estimating the number of instances of temporal motifs in temporal graphs. Overall, our sampling framework is flexible in several ways. First, the framework is built on top of exact counting algorithms. Improvements in these algorithms can be used directly within our framework, provided it meets the conditions needed by our framework. In fact, we demonstrated this to be the case with the specialized algorithm we developed for 2node, 3edge motifs in Section 5, which was faster than existing exact counting methods for that particular motif. Second, our sampling framework provides natural parallelism, which allowed us to achieve orders of magnitude speedups on algorithms that do not have obvious parallel implementations. Finally, the sampling is inherently less memory intensive, which allowed us to estimate motif counts on datasets on which exact algorithms cannot even run (Table 3); thus, our framework makes knowledge discovery feasible in new cases.
An extraordinary amount of research has gone into scalable estimation algorithms for counting patterns in static graphs. Our paper takes this line of research in a new direction by considering richer patterns that arise when temporal information is incorporated into the graph. We anticipate that our work will open new challenges for algorithm designers while simultaneously providing a solution for domain scientists working with largescale temporal networks.
Acknowledgements
This research was supported in part by NSF Award DMS1830274. We thank Eugene Y.Q. Shen for donating the computing resources that made some of the larger datasets possible.
References
 (1) The CAIDA UCSD Anonymized Internet Traces – February 17, 2011. http://www.caida.org/data/passive/passive_dataset.xml.
 (2) N. K. Ahmed, N. Duffield, J. Neville, and R. Kompella. Graph sample and hold: a framework for biggraph analytics. In Proceedings of KDD, 2014.
 (3) L. Akoglu, M. McGlohon, and C. Faloutsos. OddBall: Spotting anomalies in weighted graphs. In Advances in Knowledge Discovery and Data Mining. 2010.

(4)
H. Avron.
Counting triangles in large graphs using randomized matrix trace estimation.
In Workshop on Largescale Data Mining: Theory and Applications, volume 10, 2010.  (5) F. Battiston, V. Nicosia, M. Chavez, and V. Latora. Multilayer motif analysis of brain networks. Chaos, 27(4):047404, 2017.
 (6) A. R. Benson, D. F. Gleich, and J. Leskovec. Higherorder organization of complex networks. Science, 353(6295):163–166, 2016.
 (7) M. Bressan, F. Chierichetti, R. Kumar, S. Leucci, and A. Panconesi. Counting Graphlets: Space vs Time. In Proceedings WSDM. ACM Press, 2017.
 (8) S. Chakrabarti, M. Ester, U. Fayyad, J. Gehrke, J. Han, S. Morishita, G. PiatetskyShapiro, and W. Wang. Data mining curriculum: A proposal. Intensive Working Group of ACM SIGKDD Curriculum Committee, 140, 2006.
 (9) H. P. Chan, N. R. Zhang, and L. H. Chen. Importance sampling of word patterns in dna and protein sequences. Journal of Computational Biology, 2010.
 (10) R. Dobrin, Q. Beg, A.L. Barabási, and Z. Oltvai. Aggregation of topological motifs in the Escherichia coli transcriptional regulatory network. BMC Bioinformatics, 2004.
 (11) T. Eden, A. Levi, D. Ron, and C. Seshadhri. Approximately counting triangles in sublinear time. SIAM Journal on Computing, 46(5):1603–1646, jan 2017.
 (12) E. R. Elenberg, K. Shanmugam, M. Borokhovich, and A. G. Dimakis. Distributed estimation of graph 4profiles. In Proceedings of WWW, 2016.
 (13) 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 coevolution. In Advances in Neural Information Processing Systems, pages 1954–1962, 2015.
 (14) J. Feigenbaum, S. Kannan, A. McGregor, S. Suri, and J. Zhang. On graph problems in a semistreaming model. Theoretical Computer Science, 348(23):207–216, 2005.
 (15) N. Gaumont, C. Magnien, and M. Latapy. Finding remarkably dense sequences of contacts in link streams. Social Network Analysis and Mining, 6(1), sep 2016.
 (16) M. Gupta and J. G. Ibrahim. Variable selection in regression mixture modeling for the discovery of gene regulatory networks. Journal of the American Statistical Association, 2007.
 (17) S. Gurukar, S. Ranu, and B. Ravindran. COMMIT: A Scalable Approach to Mining Communication Motifs from Dynamic Networks. In Proceedings of SIGMOD, 2015.
 (18) J. Han, J. Pei, and M. Kamber. Data mining: concepts and techniques. Elsevier, 2011.
 (19) K. Henderson, B. Gallagher, T. EliassiRad, H. Tong, S. Basu, L. Akoglu, D. Koutra, C. Faloutsos, and L. Li. RolX: structural role extraction & mining in large graphs. In Proceedings of KDD, 2012.
 (20) J. Hessel, C. Tan, and L. Lee. Science, askscience, and badscience: On the coexistence of highly related communities. In ICWSM, 2016.
 (21) P. Holme and J. Saramäki. Temporal networks. Physics Reports, 2012.

(22)
Y. Hu, J. Trousdale, K. Josić, and E. SheaBrown.
Motif statistics and spike correlations in neuronal networks.
BMC Neuroscience, 13(Suppl 1):P43, 2012.  (23) Y. Hulovatyy, H. Chen, and T. Milenković. Exploring the structure and function of temporal networks with dynamic graphlets. Bioinformatics, 2015.
 (24) S. Jain and C. Seshadhri. A fast and provable method for estimating clique counts using Turán’s theorem. In Proceedings of WWW. ACM Press, 2017.
 (25) R. Jin, S. McCallen, and E. Almaas. Trend motif: A graph mining approach for analysis of dynamic complex networks. In Proceedings of ICDM, 2007.
 (26) D. Kondor, I. Csabai, J. Szüle, M. Pósfai, and G. Vattay. Inferring the interplay between network structure and market effects in Bitcoin. New J. of Physics, 2014.
 (27) G. Kossinets, J. Kleinberg, and D. Watts. The structure of information pathways in a social communication network. In Proceeding of KDD, 2008.

(28)
D. Koutra, U. Kang, J. Vreeken, and C. Faloutsos.
Summarizing and understanding large graphs.
Statistical Analysis and Data Mining: The ASA Data Science Journal
, 8(3):183–202, may 2015.  (29) L. Kovanen, M. Karsai, K. Kaski, J. Kertész, and J. Saramäki. Temporal motifs in timedependent networks. J. of Stat. Mech.: Theory and Experiment, 2011.
 (30) L. Kovanen, K. Kaski, J. Kertesz, and J. Saramaki. Temporal motifs reveal homophily, genderspecific patterns, and group talk in call sequences. PNAS, 2013.
 (31) M. Lahiri and T. Y. BergerWolf. Structure prediction in temporal networks using frequent subgraphs. In IEEE Symposium on Computational Intelligence and Data Mining. IEEE, 2007.
 (32) J. Leskovec, D. Huttenlocher, and J. Kleinberg. Signed networks in social media. In Proceedings CHI, 2010.
 (33) J. Leskovec, D. P. Huttenlocher, and J. M. Kleinberg. Governance in social media: A case study of the wikipedia promotion process. In Proceedings of ICWSM, 2010.
 (34) K.c. Liang, X. Wang, and D. Anastassiou. A sequential monte carlo method for motif discovery. IEEE transactions on signal processing, 56(9):4496–4507, 2008.
 (35) Y. Lim and U. Kang. MASCOT: Memoryefficient and Accurate Sampling for Counting Local Triangles in Graph Streams. In Proceedings of KDD, 2015.
 (36) P. Mackey, K. Porterfield, E. Fitzhenry, S. Choudhury, and G. Chin Jr. A chronological edgedriven approach to temporal subgraph isomorphism. arXiv, 2018.
 (37) S. Mangan and U. Alon. Structure and function of the feedforward loop network motif. PNAS, 2003.
 (38) A. McGregor. Graph stream algorithms: A survey. ACM SIGMOD Record, 2014.
 (39) C. Meydan, H. H. Otu, and O. Sezerman. Prediction of peptides binding to MHC class I and II alleles by temporal motif mining. BMC Bioinformatics, 2013.
 (40) R. Milo, S. Itzkovitz, N. Kashtan, R. Levitt, S. ShenOrr, I. Ayzenshtat, M. Sheffer, and U. Alon. Superfamilies of evolved and designed networks. Science, 2004.
 (41) R. Milo, S. ShenOrr, S. Itzkovitz, N. Kashtan, D. Chklovskii, and U. Alon. Network motifs: Simple building blocks of complex networks. Science, 2002.
 (42) C. C. Noble and D. J. Cook. Graphbased anomaly detection. In Proceedings of KDD, 2003.
 (43) A. B. Owen. Monte Carlo theory, methods and examples. 2013.
 (44) P. Panzarasa, T. Opsahl, and K. M. Carley. Patterns and dynamics of users’ behavior and interaction: Network analysis of an online community. J. of the American Society for Information Science and Technology, 60(5):911–932, 2009.
 (45) A. Paranjape, A. R. Benson, and J. Leskovec. Motifs in temporal networks. In Proceedings of WSDM. ACM Press, 2017.
 (46) N. Przulj. Biological network comparison using graphlet degree distribution. Bioinformatics, 23(2):e177–e183, jan 2007.
 (47) N. Przulj, D. G. Corneil, and I. Jurisica. Modeling interactome: scalefree or geometric? Bioinformatics, 20(18):3508–3515, jul 2004.
 (48) F. Reid and M. Harrigan. An analysis of anonymity in the bitcoin system. In Security and Privacy in Social Networks, pages 197–223. 2012.
 (49) K. Rohe and T. Qin. The blessing of transitivity in sparse and stochastic networks. arXiv, 2013.
 (50) R. A. Rossi and N. K. Ahmed. Role discovery in networks. IEEE Transactions on Knowledge and Data Engineering, 27(4):1112–1131, apr 2015.
 (51) I. Scholtes. When is a network a network?: Multiorder graphical model selection in pathways and temporal networks. In Proceedings of KDD, 2017.
 (52) C. Seshadhri, A. Pinar, and T. G. Kolda. Triadic measures on graphs: The power of wedge sampling. In Proceedings of SDM, 2013.
 (53) H. Shao, M. Marwah, and N. Ramakrishnan. A temporal motif mining approach to unsupervised energy disaggregation: Applications to residential and commercial buildings. In Proceedings of AAAI, 2013.
 (54) S. S. ShenOrr, R. Milo, S. Mangan, and U. Alon. Network motifs in the transcriptional regulation network of escherichia coli. Nature Genetics, 2002.
 (55) R. Siddharthan. Phylogibbsmp: module prediction and discriminative motiffinding by gibbs sampling. PLoS computational biology, 2008.
 (56) L. D. Stefani, A. Epasto, M. Riondato, and E. Upfal. TRIÈST: Counting Local and Global Triangles in Fully Dynamic Streams with Fixed Memory Size. ACM Transactions on Knowledge Discovery from Data, 11(4):1–50, jun 2017.
 (57) J. Sun, C. Faloutsos, S. Papadimitriou, and P. S. Yu. GraphScope: Parameterfree Mining of Large Timeevolving Graphs. In Proceedings of KDD, 2007.
 (58) C. E. Tsourakakis, U. Kang, G. L. Miller, and C. Faloutsos. DOULION: counting triangles in massive graphs with a coin. In Proceedings of KDD, 2009.
 (59) C. E. Tsourakakis, J. Pachocki, and M. Mitzenmacher. Scalable motifaware graph clustering. In Proceedings of WWW, 2017.
 (60) J. Ugander, L. Backstrom, and J. Kleinberg. Subgraph frequencies: Mapping the empirical and extremal geography of large graph collections. In Proceedings of WWW, 2013.
 (61) T. Viard, M. Latapy, and C. Magnien. Computing maximal cliques in link streams. Theoretical Computer Science, 609:245–252, jan 2016.
 (62) T. Viard, C. Magnien, and M. Latapy. Enumerating maximal cliques in link streams with durations. Information Processing Letters, 133:44–48, may 2018.
 (63) P. Wang, J. Zhao, X. Zhang, Z. Li, J. Cheng, J. C. Lui, D. Towsley, J. Tao, and X. Guan. MOSS5: A fast method of approximating counts of 5node graphlets in large graphs. IEEE Transactions on Knowledge and Data Engineering, 2018.
 (64) Q.M. Zhang, L. Lü, W.Q. Wang, and T. Z. and. Potential theory for directed networks. PLoS ONE, 2013.
 (65) X. Zhang, S. Shao, H. E. Stanley, and S. Havlin. Dynamic motifs in socioeconomic networks. Europhysics Letters, 108(5):58001, dec 2014.
 (66) Q. Zhao, Y. Tian, Q. He, N. Oliver, R. Jin, and W.C. Lee. Communication motifs: a tool to characterize social communications. In Proceedings of CIKM, 2010.