Finding Patient Zero: Learning Contagion Source with Graph Neural Networks
Authors: Chintan Shah
Nima Dehmamy
Nicola Perra
Matteo Chinazzi
Albert-Laszlo Barabasi
Alessandro Vespignani
Rose Yu
Venue:  (2021)
Keywords: contagion dynamics
theory of graph neural networks
epidemic modeling

Abstract: Locating the source of an epidemic, or patient zero (P0), can provide critical insights into  the infection's transmission course and allow efficient resource allocation. 
Existing methods use graph-theoretic centrality measures and expensive message-passing algorithms, requiring knowledge of the underlying dynamics and its parameters.
In this paper, we revisit this problem using graph neural networks (GNNs) to learn P0. 
We observe that GNNs can identify P0 close to the theoretical bound on accuracy, without explicit input of dynamics or its parameters.
In addition, GNN is over 100 times faster than classic methods for inference on arbitrary graph topologies.
Our theoretical bound also shows that the epidemic is like a ticking clock, emphasizing the importance of early contact-tracing. 
We find a maximum time after which accurate recovery of the source becomes impossible, regardless of the algorithm used.

Under review as a conference paper at ICLR 2021
FINDING PATIENT ZERO: LEARNING CONTAGION SOURCE WITH GRAPH NEURAL NETWORKS
Anonymous authors Paper under double-blind review
ABSTRACT
Locating the source of an epidemic, or patient zero (PO), can provide critical insights into the infection's 's transmission course and allow efficient resource allo- cation. Existing methods use graph- theorett centrality measures and expensive message-passing algorithms, requiring knowledge of the underlying dynamics and its parameters. In this paper, revisit this problem using graph neural networks (GNNs) to learn PO. We establish a theoretical limit for the identification of PO a class of epidemic models. We evaluate our method against different epidemic models on both synthetic and a real-world contact network considering a disease with history and characteristics of COVID-19. We observe that GNNs can identify PO close to the theoretical bound on accuracy, without explicit input of dynamics its parameters. In addition, GNN is over 100 times faster than classic methods for inference on arbitrary graph topologies. Our theoretical bound also shows that the epidemic like ticking clock, emphasizing the importance of early contact-tracing. We find a maximum time after which accurate recovery of the source becomes difficult, regardless of the algorithm used.
INTRODUCTION
The ability quickly identify the origin of an outbreak, or "finding patient zero", is critically important in the effort to contain an emerging epidemic. The identification of early transmission chains and the reconstruction of the possible paths of diffusion of the virus can be the difference between stopping an outbreak in its infancy and letting an epidemic unfold and affect a large share of a population. Hence, solving this problem would be instrumental in informing and guiding contact tracing efforts carried out by public health authorities, allowing for optimal resource allocation that can maximize the probability of an early containment of the outbreak. Disease spreading modeled as contagion process a network SSroock & Varadhan (007);; Pastor-Satorra, et al. (2015) of human-to- to-hhmon interactions where infected individuals are going to transmit the virus by infecting (with certain probability) their direct contacts. In general, contagion processes capture wide range of phenomena, from rumor propagation on social media to virus spreading over cyber-physical networks Centola & Macy (2007): ): [Baronchelli] (2018); Wang et al. (2013); Mishra & Keshri (2013).
Therefore. learning the source of contagion process would also have broader impact on various domains, from detecting sources of fake news defending malware attacks.
Learning the index case, or patient zero (PO), is a difficult problem. In this paper, we model disease spreading as contagion process (chains of transmissions) over graph. The evolution of an outbreak is noisy and highly dependent on the graph structure and disease dynamics. In addition, in worldw epidemics, there is often a delay from the start of the outbreak to when epidemic surveillance and contact tracing starts. Hence, we might only observe the state of the graph at some intermediate times without access the complete chains of transmission. Furthermore, due to its stochastic nature, the same source node might lead different epidemic spreading trajectories. Finally, learning PO from noisy observations of graph snapshots is computationally intractable and the complexity grows exponentially with the size of the graph (Shah & Zaman 2011..
Most work in learning the dynamics of a contagion process Rodriguez et al. (201)) Mei & Eisner (2017); Li al. (2018a) have focused on inferring the forward dynamics of the diffusion. In epidemiology, for example, PPstor-Satorra & Vespignani| (200T) have studied learning the temporal dynamics of diseases spreading mobility networks. The problem of learning the reverse dynamics

---

Under review as a conference paper at ICLR 2021
and identifying diffusion sources has been largely overlooked due to the aforementioned challenges. Two of the most notable exceptions in the area are "rumor centrality" "Shah & Zaman (2011) for contagion processes trees and Dynamic Message-passingg (DMP) graphs|Lokhov Lokhov et al. (2014)
but both require as input the parameters of the spreading dynamics simulations.
Our goal is to provide fresh perspectives the problem of finding patient zero using graph neural networks (GNNs) Gilmer et al. (2017) First, we conduct a rigorous analysis of learning PO based on the graph structure and the disease dynamics, allowing us find conditions for identifying PO accurately. We test our theoretical results on a set of epidemic simulations on synthetic graphs commonly used in the literature Erdös al. (1959); Albert & Barabási (2002). We also evaluate our method on realistic co-location network for the greater Boston area, finding performance similar to the synthetic data. While collecting labeled data to train GNN to find PO may not possible, training GNN using simulations on real contact-tracing data can provide a fast method for inferring PO and help with planning and resource allocation. To the best of our knowledge, our work is the first to tackle the patient zero problem with deep learning and to test the approach network. In summary, we make the following contributions:
realistic contact
We find upper bounds on the accuracy of finding patient zero in graphs with cycles, indepen- dent of the inference algorithm used.
We show that beyond certain time scale the inference becomes difficult, highlighting the importance of swift and early contact-tracing.
We demonstrate the superiority of GNNs over state-of-the-an message passing algorithms in terms of speed and accuracy. Most importantly our method is model agnostic and does not require the epidemic parameters be known.
We validate our theoretical findings using extensive experiments for different epidemic dynamics and graph structures, including a real-world co-location graph of the COVID-19 outbreak.
RELATED WORK
Learning contagion dynamics Learning forward dynamics of contagion processes on a graph is a well studied problem area. For instance, RRodiiuezz et al. (201T): Du et al. (2013) proposed scalable algorithms to estimate the parameters of the underlying diffusion network, a problem known as network inference. Deep learning has led to novel neural network models that can learn forward dynamics of various processes including neural Hawkes processes Mei & Eisner (2017) and Markov decision processes-based reinforcement learning Li et al. (2018a). Learning forward contagion dynamics have also been intensively studied epidemiology Pastor-Satorras & Vespignani] (200T); Vynnycky & White (2010), social science Matsubara et al. (2012), and cyber-security Prakash et al. (2012). In contrast, research in learning the reverse dynamics of contagion processes is rather scarce. Influence maximization Kempe et al. (2003), for instance, finds small set of individuals that can effectively spread information in a graph, but only maximizes the number of affected nodes in the infinite time limit. Our problem is more difficult we care not just about the number of infected nodes, but which nodes were infected.
Finding patient zero In order to find patient zero, we aim to learn the reverse dynamics of contagion processes. [Shah& & Zaman (201T) were among the first formalize the problem on trees in the context of modeling rumor spreading in a network. Prakash et al. (2012); Vosoughi et al. (2117)
studied similar problems for detecting viruses in computer networks. More recent advances proposed
a dynamic message passing algorithm Lokhov et al. (2014) and belief propagation Altarelli et al.
(2014) to estimate the epidemic outbreak source. Fairly recently, Fanti & Viswanath (2017) reduced the deanonymization of Bitcoin to the source identification problem in an epidemic and analyzes the dynamies properties. On the theoretical side, (Shah & Zaman (201T): Wang et al. (2014) analyzed the quality of the maximum likelihood estimator and rumor centrality, but only for the simple SI model on trees. Antulov-Fannuu model using exact analytical methods and Monte Carlo estimators. Khim & Loh (2016); [ubbeck et al. (2017) proved that it is possible to construct a confidence set for the predicted diffusion source nodes with size independent of the number of infected nodes over regular tree. Our work provides fresh perspectives on the patient zero problem on general graphs based on the recent development of graph neural networks
et al. 2015)) found detectability limits for patient zero in the SIR
2

---

Under review as a conference paper at ICLR 2021
Graph neural networks Graph neural networks have received considerable attention (see several
references in Bronstein et al. (2017); Zhang et al. (2018); Wu et al. (2019); Goyal Ferrara (2018)).
While most research is focused static graphs, few have explored dynamic g Li et al. (2018b);
You et al. 2018); Kipf et al. (2018) Pareja et al. 2019); Trivedi et al. (2019). For example, Kipf et al. (2018) propose deep graph model to learn both the graph attribute and structure dynamics. They use recurrent decoder to forecast the node attributes for multiple time steps ahead. Trivedi et al. (2019) take continuous-tine
modeling approach where they take the node embedding as the input and model the occurrence of an edge point process. Xu et al. (2020) propose a temporal graph attention layer to learn the representations of temporal graphs. However, most research is designed for link prediction tasks and none of these existing studies have studied the problem of learning the
source of the dynamics on a graph.
CONTAGION PROCESS AND PATIENT ZERO
Finding patient zero means tracing the contagion dynamics back to its initial state and identifying the first nodes that started spreading. Here, we describe the disease dynamics on a network using busceptible-Infeeteed Recoverred (SIR) and busceptible-Exposee -ineee--eeeeeeeeeeeeeeeeed & McKendrick (1277) compartmental models that assume that infected individuals develop immunity once they recover from the infections.
(SEIR) (EEIR))Kermack
3.1 CONTAGION PROCESSES ON NETWORKS
In the SIR model, the population is split into three compartments: susceptible (S) who are susceptible to infection by the disease; infected (!) who have caught the disease and are infec- tious; removed (R) who are removed from con- sideration after experiencing the full infectious period.
Continuous time model For a contagion pro- cess a graph G with N nodes, each vertex represents an individual who is in contact only with its neighbors. We can represent the graph using the adjacency matrix A E N where Alú, = 1 if two individuals are connected, otherwise. Let S,, I,, R; be the average probabilities of node i being in each of the states, with S, + + R₁ 1. The SIR dynamics on graph is given by Newman (2018): ds.
Figure I: Visualization of the patient zero prob- lem: uncover the original source (red node, left) given a future state of contagion process (right).
dR, dt = 7
dS,
dl. dR, 0.
dt - Σ SigT
dt
dt
dt
where 3 is the infection rate per contact and y the recovery/death rate. The basic reproductive rate of a discase Ro 3// is defined the number of secondary infections created by index case in fully susceptible population Keeling & Rohani (201I). Here Y₁ is the largest eigenvalue of A. The disease will spread and result in an epidemic if Ro
Discrete time model We can also use an equivalent discrete time SIR model. Let 2 E {S, R} be the state of node time t. For susceptible node its probability to become infected removed at time f + 1 is
p - = -1 II (1- 3A),(())
P(at+i-Ross R =)
The SIR model doesn't account for the incubation period, where an individual is infected but not infectious. This is remedied by introducing "exposed" (E) state, leading the SEIR model. For susceptible node i, the probability to enter the exposed state, and becoming infectious at time + 1 is
P(a = Ela = 5) = 1 II BAg
P(x = M|a = E) = a, (3)
3

---

Under review as a conference paper at ICLR 2021
An infected node eventually enters the removed state with probability y, which is the same as SIR 2 and (3) yield for very small / as II, (1 - B, ze 3Ai I (proof in supp. B
Finding patient zero Finding PO can be formulated as maximum likelihood estimation problem for SIR and SEIR models. Specifically, we observe snapshot of the state of the nodes time step t as x" := (21, ... x),,, with each node's 's state xद E {S, E, I. RJ. The problem of finding PO is
to search for set of nodes 2 {i|x = T, f1,... N}) which led the observed state x².
common approach is to find 2 such that the likelihood of observing x" is maximized:
z argmax 2,12|5h 2 P)
(4)
where Pxx|| 2) the probability of observing x with being the PO set. We assume the number of POs is no larger than k. Estimating the original state of the dynamic system given the future states requires computing the conditional likelihood P(x*) [) exactly, which is intractable due to the combinatories of possible transmission routes.
3.2 FUNDAMENTAL LIMIT OF FINDING PATIENT ZERO
The technical difficulty of finding PO in SIR and SEIR stems from: (1) presence of cycles in graphs (higher-order transmission) (2) the removed state introducing additional uncertainty about temporal order of infections (3) uncertainty of the exact time step of the observed states. For SI dynamics (i.e. infection is permanent) on trees, existing theoretical results (Shh Zaman (20IT); Khim & Loh (2016) have established upper bounds on the detection probability based on estimator called "rumor centrality". For graphs with cycles, finding PO becomes more elusive. We derive the fundamental limit considering the case where time 0 one node, PO, is infected and all of the other nodes susceptible.
Ambiguity of patient zero on cyclic graphs For graph with cycles, if cycle embedded within the infected subgraph, it will reduce the accuracy of predicting PO because multiple spreading scenarios can lead to the same infection pattern in the cycle. For instance, take a 3-regular tree where the infection has started from the root and spread some level. If any two branches connected to make a cycle, both branches become equally likely to have spread the disease downstream. Therefore, having cycles in the infected subgraph reduces the accuracy of finding PO. If O(N) nodes in graph with cycles are infected, cycles will likely interfere with finding PO. Based on this observation, the following theorem estimates the time horizon beyond which finding PO becomes difficult.
Theorem (Time Horizon). Assume SIR dynamics ① connected graph of N nodes, starting with a single patient zero. Denoting the adjacency matrix by A and its largest eigenvalue by X the average infection probability, both over nodes and choice of patient zero, (I()) = (t)/NN) becomes (1) after twax time steps given by
log N
3A1 Ro
max
(5)
7(Ro 1)'
7
Proof:
The proof (in Supp B.5. follows from the exponential behavior of |I at small t Newman
2018
I(t) = Σ exp itBA 1(00 exp [(31 - y)t| (WOD) I(o)) v.
(6)
and the leading eigenvector being positive E, (1 ||0001111 = (Perron-Frooeeiuu - theorem).
The maximum detection accuracy of PO on connected graph would decrease with increasing number of cycles in the infected subgraph, Gy. In special cases we can estimate the penalty to the accuracy due to cycles. For example, in connected Erdös-Rényi (ER) random graphs Erdös et al. (959), where each edge has independent an probability p, the number of cycles can be estimated. To get a conservative estimate, we focus triangles, as they the most prevalent cycles in ER, get an upper bound for the accuracy in the following theorem. This ignores cases where the presence of triangles causes downstream error the arising from other types of cycles.
4

---

Under review as a conference paper at ICLR 2021
Theorem 2 (Detection Accuracy). In contagion process a connected random graph G, with edge probability p and with infected subgraph G1. the prediction accuracy for PO is bounded from above
1/4 3 + (-
1 2
P max
The proof (supp. | follows from estimating number of triangles in subgraph Gy of a dense ER graph and noting each triangle can drop the accuracy of PO to 1/3.
Figure 2 shows an example of how this upper bound behaves for different values of Rg. The graph is uniformly connected ER of N 100, p 2 log N/N and with y 0.4. In conclusion, graphs with cycles, we expect finding PO after time tmax ~ Ollog N) to become difficult. This suggests that find PO contact-tracing must be done promptly and in early stages.
FINDING PATIENT ZERO WITH GRAPH NEURAL NETWORKS
We propose using GNNs for finding PO and show that we can improve significantly upon state-of- the-art methods, e.g. DMP. Moreover, using GNNs gives us the distinct advantage that they are model-agnen aggnostic and do not require access to the epidemic dynamics parameters or the time t of the graph snapshot. Our goal is not to propose novel graph neural network architecture, but to understand the trade-off between different probabilistic inference methods in the context of contagion dynamics. Before we discuss our GNN solution, we briefly review Dynamic Message Passing.
Dynamic Message Passing DMPhookov et al. (2014) estimates the probability of every node being the PO in the SIR model using sessage-passing equations and approximates the joint likelihood with a mean-field field time approach by assuming the following factorization:
2(x x II P(a
II P(a) 23) 1:0n-=1
II Ne
P 2)
(8)
The algorithmic complexity of the DMP equations over a graph with N nodes and / T diffusion steps is O(TN2(k) where (k) is the average degree of the graph. Furthermore, DMP requires providing the SIR epidemic parameters and the time of the graph snapshot before performing inference. For comparison, on a connected random graph, (k) > log N, yielding > O(TN2 log N) time complexity for DMP. trained GNN with layers (2-8 in our case) has complexity O (N2)) in the inference step and does not require inputting the model parameters. This makes it harder to scale DMP for large or dense graphs. DMP is proven to be exact on trees, e.g. Kanoria et al. (011)) and has been used more general graphs with reasonable success. In practice, we find that GNN inference can be 50-100 times faster than DMP, as DMP has to exhaustively infer the marginals of S, I, R states for each node for each time t 2 O. DMP is trying to accurately infer the distribution over the nodes and regr-seegryy infer PO, whereas GNN does not infer node states at each time t and tries to directly infer PO. On the other hand, DMP is unsupervised inference algorithm and does not "learn", while GNN is provided with labels in a standard supervised setting.
Relation between Contagion Dynamics and GNNs Our use of GNNs for finding PO is motivated by the fact that the contagion dynamics ( are a special case of Reaction-Diffusion (RD) processes on graphs Colizza et al. (2007) which structurally equivalent to GNNs, shown next.
Proposition 1. Reaction criftiinndiffusiin dynamics graphs is structurally equivalent the message- passing neural network ansatz.
at Denoting x²(1)= = P(a = (1) of node being in states such as μ E [, or (S, E, time I, Markovian reaction-diffusie dynamics can be written as
(t + 1)= = (
Σ FF(AA, (P(()))))))))) Mo))r)
"t( = Σ WH WHeD + ba
(9)
where 4G (f with A(-) being the step function and ((.) nonlinear function. RD on
graphs is structurally equivalent to Message- passiin Neural Networks (APNN)Giilmed et al. 20117 as RD involves bessage-passing step and node-wise interaction among features (Supp. B.3
5

---

Under review as a conference paper at ICLR 2021
same MPNN. We choose the simpler architecture of GCN as in (T0) instead of general MPNN. Finding PO requires learning the backward dynamics of RD, which seems to require the inverse of the propagation rule (PR). Yet, each node only get infected by its neighbors, backward dynamics requires message passing over the same adjacency matrix and should again have the structure of RD.
the
4.1 LEARNING WITH GRAPH NEURAL NETWORKS
We employ a ofate-of-the-art GNN design, suggested by Dwivedi et al. (2020). We make several modifications to the model architecture to fit our problem. Given one-hot encoded node states of {0, 1)MM as the GNN input, where M is the number of states and where the states are either {S, E, or (S,1,R),
first apply a linear transformation h (0) Uz! with U RCXM
Denote the output of layer by h ( where is the node index. We use graph convolutional network (GCN) Kipf& & Welling (0166) each layer g(h) = og (f(A) h W + b). where W E RCCC b R and f is called the propagation rule GCN. We use f(A) = propagation rule, where Dij = òij Ex Aas is the degree matrix. To include features of the central node, instead of adding self-loops, we use residual connections between GCN layers and notice a significant increase in model performance. The action of these higher GNN layers is given by
for the
((+1) = h1 + (BN(g(h()))) h
y: P ReLu(Q - )
(10)
where is the number of layers and the output layer is parameterized by Q E DxD and P P RIXD to generate E IR representing the probability that node PO. BN(.) denotes Batch Normalization and (-) is leaky-relu nonlinear activation function.
5 EXPERIMENTS
We perform extensive studies on performance of our GCNs in finding PO in SIR and SEIR dynamics over synthetic graphs with various graph topologies. In addition, we generate synthetic epidemic outbreaks on top of a real world co-location network using SEIR compartmental model that is calibrated to simulate a contagion process with characteristics similar to a COVID-19 outbreak.
Experimental Setup We compare the performance of DMP Lokhov al. (014)) and different variants of GCNs, following the architecture we described in sec. 4.1
DMP: Dynamic Message Passing algorithm
Lokhov et al. (2014), we sample graph snapshot O at time & with each node having state xf E {S, 1, RJ, and select the node that has the highest likelihood generating O, that is I).
PO argmax, P(O)z
GCN-S: symmetric GCN Kipf & Welling (2016) f(A) AD/// GCN-R: random walk f(A) D-1A, GCN-N4 M: mixture of propagation rules f(A) A1
CAT: Graph Attention Network Veličković et al. (2017)
We train our models using DGL Wang al. (2019) with PyTorch backend. The task is to predict the probability for each node being PO given a single graph snapshot. We report performance averaged over four random seeds. Details training and hardware given in Supp. A.2 We additionally report inference run times.
Evaluation Metrics We use top-1 accuracy to understand the effectiveness of our method. However, due to the ambiguity of detecting patient zero, elaborated in theoretical analysis, top-1 accuracy may not be the only evaluation measure to be relied upon. Therefore we also calculate the normalized 1
D ΣVED, Ta where D, is the set of test samples time N is the size of the graph and r. is the index of the ground truth PO in the reverse-sorte probability distribution. Normalized rank is retrieval metric that tells us how high the correct patient zero was in the learned output distribution. demonstrates the quality of the output distribution in learning the stochastic dynamics and helps us understand how high PO ranked even if it was ranked the first.
rank defined by RG

---

Under review as a conference paper at ICLR 2021
Performance Comparison on Synthetic Grephs
Inference times
DMP
GCAAS
GEN-R
GCNM
GAT
Accunccy P4
02
Dataset
DMP
GCN
GAT
BA-TTre
14.40 hr
3.89s
3,18s
64
BA-Dense
77.04 hr
4.91s
8.19s
ER-Dense
71.77 hr
4.93s
9.66s
Geometric
70.35 hr
5.34s
10.87s
60
BA-Troe
BA-Dense
DR-Dense
Geometric
Figure 3: Mean Prediction accuracy/speed comparison for different methods for the test set over T = 30 steps and Ro 2.5. The time to perform inference over the test set for different models have been listed above. Note that the time taken by GCN represents the mean time taken by GCN variants. We observe that GNNs beat DMP by large margin both in terms of speed and accuracy. The similar performance of the three GNN models may be due to approaching the theoretical limits.
5.1 EXPERIMENTS WITH SYNTHETIC NETWORKS
We first experiment on three synthetic graph models: ER random graph, Garabási-Alberr (BA) graph Albert & Barabási (202) and Random Geometric Graph (Geometric or RGG) Dall & Chris tensen (2002), and use NDLib (Rossetti et al. 2017) simulate SIR and SEIR epidemic dynamics on the graph (Supp. A.I. We pick a random PO seed node at time t = 0 and run S(E)IR fixed number of steps T. The epidemic parameters (ax, 3, () are chosen such that we can vary Ro to study model performance. We set eigenvalue of the graph. For SEIR, we set a = 0.5. We generate 20, 000 simulations and use 80 - 10 10 train- validation--ess split. For each sample we select & E {1, (1....T' T uniformly at random and try to predict PO at time t = 0 given the graph adjacency matrix A and node features xf.
0.4 and / = Roy/A1 where A is the largest
Empirical Accuracy Theory
1.2
We first compare the top-1 prediction accuracy in SIR and SEIR for different models averaged over < f < T. Fig compares the prediction accuracy and average inference time for different models. We can see that GNN-based models outperform the baseline DMP both in accuracy and efficiency. We also want to note that the training time for GNNs is under 7 hours, sigmificantly less than the fastest DMP run of 14.40 hr. It is also worth emphasizing that DMP requires explicit input of B,7 7 and it while GNNs are model agnostic.
2.cccc)
1.0
Schhea)
Rg 2.5thea)
R.- 3.5
0.6
0.4
0.0
1.5
2.0
2.5
35
40
Time/tear
Figure 2: Theoretical upper bound on accuracy (dashed line) vs ex- perimental (solid line) on ER for varying Rg. While accuracy drops below the theoretical limit at t tmax-
To validate our theory, we plot the theoretical accuracy upper- bound and the empirical accuracy obtained from GNN in Fig. 2 We note that the time scale and drop in accuracy are consistent with our theoretical results on mmax (28) and the upper bound on accuracy P Pmaxx 77). Combined with the fact that all our GNN models have comparable accuracies, this suggests that our GCN-based models may be approaching the fundamental limits we described in 3.2
Fig. 4 shows the trend of accuracy decay over the time steps t for different graph structures and Ro values. As expected, the accuracy is highest on a tree and when is small. In graphs with cycles (BA-Dense, ER-Dense, and Geometric) we also observe nontrivial drop in accuracy which depends both on l and Ro. For SIR we observe a drop in accuracy as a function of Ro and time, consistent with our theoretical upper bound. The decay is slower for SEIR as the latent stage adds delay to the spread of the epidemic. The normalized rank of PO remains high even over longer time horizons, indicating that PO could be narrowed down to small subset of the population with impressive accuracy.
5.2 BOSTON CO-LOCATION NETWORK AND COVID-19 EPIDEMIC TRAJECTORY
Our real-world dataset consists of co-location graph and simulations of an epidemic with the natural progression of COVID-19. The co-location graph is constructed using the Cuebiq data]] for two weeks from 23 March, 2020 to 5 April, 2020 (N 384, 590 nodes). To reduce computational E 30, 376 edges while maintaining the
costs, we sample subgraph with N = 2, 689 nodes and
httpp::/ hwwp///wwwwaa/ cuebi com wttp::/tt./ooo.aaaatooaaa/aaot
for good/ derived from Klein al. 2020a
7

---

Under review as a conference paper at ICLR 2021
[a] BA-Tree
(b) BA Dense
(c) Geometric
(d) Dense
0.88
0.4
Accuracy
0.8
SEI
0.6
0.2
SIR
0.8
0.6
1.00
SEIR
Rank
0.97
0.94
2.5
0.91
17 29 Time Steps
21
13 17 25 Time Steps
5
17 Time Steps
29
Time Steps
Figure 4: Performance of GCN-S for SIR and SEIR epidemic dynamics as top-1 accuracy over the test set. The top-1 recovery accuracy of PO time (first and second row) and normalized rank of PO (third and fourth row) for different graph topologies with varying Rg values. Note that in BA-Tree (a), which is a tree, the accuracy remains fairly high in both SIR and SEIR, consistent with existing literature, and confirming that cycles significantly reduces accuracy of PO.
degree distribution and connectivity patterns the original graph. For the epidemic simulations,
we run a modified SEIR model with asymptomatic infectious states the co-location graph with Ro resembling COVID-19 Chinazzi et al. (2020) and accordingly set Ro 2.5. Each simulation
contains patient zero, selected uniformly at random. The simulation run for 50 days. We create dataset with 10, 000 samples and an 80 10 10 train- validation test split (supp. C
The top-k accuracy performance over differ- ent days when the graph snapshot was ob- served are shown Fig. [a We can see that the top-1 accuracy falls steadily over time, the top-(10, 20) accuracy remains fairly high for the first two weeks suggesting that we can retrieve PO in the most likely 20 nodes out of total 2, 689 candidates.
1.
Raa
0.6
Normeltzes
0.96
0.92
Bo2
0.88
43
15
(a) Top-k accuracy
(b) Normalized rank
Interestingly, while the top-1 accuracy de- creases significantly, degrading by 50% after 14 days, using normalized rank, the model can narrow down the of patient zeros ac- curately even later in the epidemic, as shown Fig [bb, For the normalized rank, PO can be recovered fairly accurately in the first two weeks the epidemic. These results highlight an important trade-off between accurately determining patient zero and retrieving the general infected region.
Figure 5: Performance on Boston co-location net- work with simulations following the natural history of COVID-19, Shown here are the top-k accuracy and normalized rank.
CONCLUSION
We study contagion dynamics on a graph using graph neural networks (GNNs) to learn the reverse dynamics of contagion processes and predict patient zero. We evaluate our method against different epidemic models on both synthetic and real-world contact network with disease with the natural history and characteristics COVID-19. outbreak without explicit input of dynamics parameters. Most notably, GNN accuracy approaches our predicted theoretical upper bound, indicating that further architecture refinements may not improve performance significantly. In addition, GNN is over 100x faster for inference than classic methods for arbitrary graph topologies. Extensions this work may include learning using sequences of graph snapshots, well allowing set of patient zeros.
observe that GNNs can efficiently infer the source an
8

---

Under review as a conference paper at ICLR 2021
REFERENCES
Réka Alhert and Albert-László Barabási. Statistical mechanics of complex networks. Reviews of modern physics, 74(1):47, 2002.
Fabrizio Altarelli, Alfredo Braunstein, Luca Dall' Asta, Alejandro Lage-Castellanos, and Riccardo Zecchina. Bayesian inference of epidemics on networks via belief propagation. Physical review letters, (1((1)):18700) 2014.
Nino Antulov-Fantulin Alen Lančié, Tomislav Smue, Hrvoje Štefančić, and Mile Šikić. Identification of patient zero in static and temporal networks: Robustness and limitations. Physical review letters, 114821)248701, 2015.
Andrea Baronchelli. The emergence of consensus: primer. Royal Society open science, 5(2):172889, 2018.
Lukas Biewald. Experiment tracking with weights and biases, 2020. URL https::/t///tt http::://wwwwaa.CC. com/ Software available from wandb.com.
Michael M Bronstein, Joan Bruna, Yann LeCun, Arthur Szlam, and Pierre Vandergheynst. Geometric deep learning: going beyond euclidean data. IEEE Signal Processing Magazine, 34(4):18-42, 2017.
Sébastien Bubeck, Lue Devroye, and Gábor Lugosi. Finding adam in random growing trees. Random Structures & Algorithms, 50(2):158-11 2017.
Damon Centola and Michael Macy. Complex contagions and the weakness of long ties, American journal of Sociology, 113(3):02--334 2007.
Matteo Chinazzi, Jessica T. Davis, Marco Ajelli, Corrado Gioannini, Maria Litvinova, Stefano Merler, Ana Pastore Piontti, Kunpeng Mu, Luca Rossi, Kaiyuan Sun, Cécile Viboud, Xinyue Xiong, Hongjie Yu, M. Elizabeth Halloran, Ira M. Longini, and Alessandro Vespignani. The effect of travel restrictions on the spread of the 2019 novel coronavirus (covid-19) outbreak. doi: 10. (0.112//cieecceaaaa977 abb77777 URL
Science, 368(6489) 2020, ISSN 0036-8075.
https://science www.sscieeaaa.....e http://cceeeeeeeccccnnnaaaa..o// wwt.e//ccna...oor/oooottttt////e
Vittoria Colizza, Romualdo Pastor-Satoruse and Alessandro Vespignani. fraction-diffusiio pro- cesses and metapopulation models heterogeneous networks. Nature Physics, 3(4):276-282, 2007.
Jesper Dall and Michael Christensen. Random geometric graphs. Physical review E, 66(1):01612,, 2002.
Nan Du, Le Song, Manuel Gomez Rodriguez, and Hongyuan Zha. Scalable influence estimation in continuous-ting time diffusion networks. In Advances in neural information processing systems, pp. 3147-3155, 2013.
Vijay Prakash Dwivedi, Chaitanya K. Joshi, Thomas Laurent, Yoshua Bengio, and Xavier Bresson. Benchmarking graph neural networks, 2020.
Paul Erdös, Alfréd Rényi, et al. On random graphs. Publicationes mathematicae, 6(26))290-297. 1959.
Giulia Fanti and Pramod Viswanath. Deanonymization in the bitcoin p2p network. In Advances in Neural Information Processing Systems, pp. 1364-1311 1373, 2017.
Justin Gilmer, Samuel S Schoenholz, Patrick F Riley, Oriol Vinyals, and George E Dahl. Neural message passing for quantum chemistry. In Proceedings of the 34th International Conference on Machine Learning- Volume 70, pp. 1263-1272. JMLR. org, 2017.
Palash Goyal and Emilio Ferrara. Graph embedding techniques, applications, and performance: A survey. Knowledge-Bann -ased Systems, 151:78-94, 2018.
Yashodhan Kanoria, Andrea Montanari, et al. Majority dynamics on trees and the dynamic cavity method. The Annals of Applied Probability, 21(5):1694 1744118888 2011.

---

Under review as a conference paper at ICLR 2021
Matt J Keeling and Pejman Rohani. Modeling infectious diseases in humans and animals. Princeton University Press, 2011.
David Kempe, Jon Kleinberg, and Éva Tardos. Maximizing the spread of influence through social network. Proceedings of the ninth ACM SIGKDD international conference on Knowledge discovery and data mining, pp. 137-146, 2003.
William Ogilvy Kermack and Anderson G McKendrick. A contribution to the mathematical theory of epidemics. Proceedings of the royal society of london. Series A, Containing papers of mathematical and physical character, 115(772):700-77)
1927.
Justin Khim and Po-Ling Loh. Confidence sets for the source of difffusion in regular trees. IEEE Transactions on Network Science and Engineering, 4(1):27-40, 2016.
Thomas Kipf, Ethan Fetaya, Kuan-Chich Wang, Max Welling, and Richard Zemel. Neural relational inference for interacting systems, arXiv preprint arXiv: /802.:/802.0445 046887, 2018.
Thomas N Kipf and Max Welling. Semi-supervise classification with graph convolutional networks.
arXiv preprint 16Xiv: (6099002907) 2016.
Brennan Klein, LaRocky, S McCabey, L. Torresy, Filippo Privitera, Brennan Lake, Moritz UG Kraemer, John S Brownstein, David Lazer, Tina Eliassi-Rar el. Assessing changes in commuting and individual mobility in major metropolitan areas in the united states during the covid-19 outbreak, 2020a.
Brennan Klein et al. Reshaping a nation: Mobility, commuting, and contact patterns during the covid-19 outbreak, 2020b.
Shuang Li, Shuai Xiao, Shixiang Zhu, Nan Du, Yao Xie. and Le Song. Learning temporal point processes via reinforcement learning. Advances in neural information processing systems, pp. 10781-10791, 201 8a.
Yaguang Li, Rose Yu, Cyrus Shahabi, and Yan Liu. Diffusion convolutional recurrent neural network: Data-driven traffic forecasting. International Conference on Learning Representations (ICLR), 2018b.
Andrey Y Lokhov, Mare Mézard, Hiroki Ohta, and Lenka Zdeborová. Inferring the origin of an epidemic with dynamic message-passing algorithm. Physical Review E, 90(1):01280) 2014.
Naoki Masuda and Petter Holme. Temporal network epidemiology. Springer, 2017.
Yasuko Matsubara, Yasushi Sakurai, Aditya Prakash, Lei Li, and Christos Faloutsos. Rise and fall patterns of information diffusion: model and implications. In Proceedings of the 18th ACM SIGKDD international conference on Knowledge discovery and data mining, pp. 6-14, 2012.
Hongyuan Mei and Jason M Eisner. The neural hawkes process: A neurally self-modulating multivariate point process. In Advances in Neural Information Processing Systems, pp. 6754-676666 6764, 2017.
Bimal Kumar Mishra and Neha Keshri. Mathematical model on the transmission of worms in wireless sensor network. Applied Mathematical Modelling, 4116):103--111,2220 2013.
NetworkX developer team. Networkx, 2014. URL https::/ www..g.itthuioooo 10/
Mark Newman. Networks: introduction. Oxford university press, 2018.
Aldo Pareja, Giacomo Domeniconi, Jie Chen, Tengfei Ma, Toyotaro Suzumura, Hiroki Kanezashi, Tim Kaler, and Charles E Leisersen. Evolvegen: Evolving graph convolutional networks for dynamic graphs arXiv preprint arXiv:1902.10/99 2019.
Romualdo Pastor-Satorras and Alessandro Vespignani. Epidemic spreading in scale-free networks.
Physical review letters, 86(14):3200, 2001.
Romualdo PattrrsSatorras, Claudio Castellano, Piet Van Mieghem, and Alessandro Vespignani.
Epidemic processes in complex networks. Reviews of modern physics, 87(3):925, 2015.
10

---

Under review as a conference paper at ICLR 2021
B Aditya Prakash, Jilles Vreeken, and Christos Faloutsos. Spotting culprits in epidemics: How many and which ones? In 2012 IEEE 12th International Conference on Data Mining, pp. 11-20. IEEE, 2012.
Manuel Gomez Rodriguez, David Balduzzi, and Bernhard Schölkopf. Uncovering the temporal dynamics of diffusion networks. 2011.
Giulio Rossetti, Letizia Milli, Salvatore Rinzivillo, Alina Sirbu, Dino Pedreschi, and Fosca Gi- annotti. Ndlib: python library model and analyze diffusion processes over complex net- works. International Journal of Data Science and Analytics, 5(1):61-79, Dec 2017. ISSN 2364-4168. doi:
URL http:/ /xx doi
141060-017-0086-1
Devavrat Shah and Tauhid Zaman. Rumors in network: Who's the culprit? IEEE Transactions information theory, 57(9):5133-518,, 2011.
Daniel W Stroock and SR Srinivasa Varadhan. Multidimensional diffusion processes. Springer, 2007.
Rakshit Trivedi, Mehrdad Farajtabar, Prasenjeet Biswal, and Hongyuan Zha. Dyrep: Learning representations over dynamic graphs. International Conference on Learning Representations, 2019.
Petar Veličković, Guillem Cucurull, Arantxa Casanova, Adriana Romero, Pietro Lio, and Yoshua Bengio. Graph attention networks. arXiv preprint arXiv: 17xiv:11100 10903, 2017.
Soroush Vosoughi, Mostafa 'Neo' Mohsenvand, and Deb Roy. Rumor gauge: Predicting the veracity of rumors on twitter. ACM transactions on knowledge discovery from data (TKDD), 11(4))1-36, 2017.
Emilia Vynnycky and Richard White. An introduction to infectious disease modelling. OUP oxford, 2010.
Minjie Wang, Lingfan Yu, Da Zheng, Quan Gan, Yu Gai, Zihao Ye, Mufei Li, Jinjing Zhou, Qi Huang, Chao Ma, Ziyue Huang, Qipeng Guo, Hao Zhang. Haibin Lin, Junbo Zhao, Jinyang Li, Alexander J Smola, and Zheng Zhang. Deep graph library: Towards efficient and scalable deep learning on graphs. ICLR Workshop on Representation Learning on Graphs and Manifolds, 2019. URL
https://axxiv.. httpaaxxvv..ggabs///o/ ottp/bs////a////.o// 01315
Yini Wang, Sheng Wen, Yang Xiang, and Wanlei Zhou. Modeling the propagation of worms in networks: survey, IEEE Communications Surveys & Tutorials, 16(2):942-960, 2013.
Zhaoxu Wang, Wenxiang Dong, Wenyi Zhang, and Chee Wei Tan. Rumor source detection with multiple observations: Fundamental limits and algorithms. ACM SIGMETRICS Performance Evaluation Review, 42(1):1-13, 2014.
Zonghan Wu, Shirui Pan, Fengwen Chen, Guodong Long, Chengqi Zhang, and Philip S Yu. A comprehensive survey on graph neural networks. arXiv preprint arXiv:1901.005000 2019.
Da Xu, Chuanwei Ruan, Evren Korpeoglu, Sushant Kumar, and Kannan Achan. Inductive representa- tion learning temporal graphs. arXiv preprini 07X7v:2002.07962, 2020.
Jiaxuan You, Bowen Liu, Rex Ying, Vijay Pande, and Jure Leskovec. Graph convolutional policy network for goal-directed molecular graph generation. arXiv preprint 1E16.080000277. 2018.
Ziwei Zhang, Peng Cui, and Wenwu Zhu. Deep learning on graphs: A survey, arXiv preprint 1eXi::1812.04202 2018.
11

---

Under review as a conference paper at ICLR 2021
A APPENDIX
A.1 DATASET DETAILS
We use three graph models: ER random graph, Barabási-Albert (BA) graph Albert Barabási (2002) and Random Geometric Graph (Geometric or RGG) Dall & Christensen (2002). The density of all three ( E/ (2) is adjustable, but BA produce exact trees. Fixing the number of nodes to N = 1, 000, we first obtain one random instance tree BA, and dense BA, ER and Geometric graphs with F 10, 000 using the NetworkX library (NetworkX developer team 2014) and then use NDLib (Rossetti et af.) 2017) to simulate SIR and SEIR epidemic dynamics on the graph (supp. A. For each sample graph, we pick PO seed node at random to be the patient zero at time f = and then we run SEE)RR a fixed number of steps 7. The epidemic parameters (a, 3,7) are chosen such that can vary Rg study model performance. We set 7 0.4 and 3 Roy/A, A1 where a is the largest eigenvalue of the graph. For SEIR, we set a = 0.5. We generate 20, 000 simulations and use 80 10 10 nain-validation-ter split. For each sample we select {1.- ... T} uniformly at random and try to predict PO at time = 0 given the graph adjacency matrix A and node features 24.
Table T describes the details of the synthetic datasets.
Table 1: Description of the sampled graph statistics
Dataset
of Nodes
of Edges
Density
Diameter
BA-Tree
1,000
999
0.99
19
BA-Dense 1,000
9,900
9.90
4
Geometric
1,000
9,282
9.28
21
ER-Dense
1,000
9,930
9.93
A.2 TRAINING AND HARDWARE
We train the model with an ADAM optimizer for 150 epochs with an initial learning rate of 0.003 and decay the learning rate by 0.5 when the validation loss plateaus with a patience of 10 epochs. We perform hyperparameter tuning over a validation set with a random search strategy. We sweep the hyperparameter space and track experiments using Weights and Biases Biewald (020)) choosing the model with the lowest validation error. We run our experiments on Nvidia 2080T1 GPUs and report performance averaged over 4 random seeds.
A.3 HYPER-PARAMEN PARAMETER DETAILS
Table Description of parameerrsss used. All of our models have been trained with 4 random seeds. The initial learning rate is mentioned in the table below and additionally we decay the learning rate by 0.5 with a patience of 10 epochs when the validation error plateaus. Note that GAT had attention heads and has been trained with 5 layers due to limitation on GPU memory.
Hyperparameters
GCN-S
GCN-R
GCN-M
GAT
Number of Epochs
150
150
150
150
Batch Size
128
128
128
32
GNN Hidden Dim
128
128
128
128
Dropout
0.265
0.265
0.265
0.265
Number of GNN Layers
10
10
10
Initial Learning Rate
0.0033
0.0033
0.0033
0.004
A.4 NOTES ON DMP IMPLEMENTATION
We include DMP
DMP Lokhov et al. (0014) as a baseline against our proposed GNN based method. As DMP does not have code that is publicly available, we implemented DMP using Python for fair
12

---

Under review as a conference paper at ICLR 2021
comparison with GNNs. Accordingly, our implementation of DMP uses DGL Wang et al.) (2019) which enables us to vectorize belief propagation (BP) and marginalization and now it runs in parallel for all nodes.
Given a graph G(V, B), we observe O' as the state of the graph with nodes E V. DMP employs MLE estimation to determine the node iro that may have led to the observed snapshot O. For single sample in our dataset D, we use algorithm [ In order to implement DMP efficiently, we implemented it message-passing on graph using DGL. We sequentially initialize node and edge features for all node and then as we obtaín N = [V| set of graphs with node i acting as PO in G,, argmax, P( P(O The advantage of our implementation then is that
DMP then allows us to obtain we can process all N graphs in parallel as if it were one large graph with N nodes and E² edges thanks to DGL'ss support for batching graphs. passing framework allows to additionally process all the nodes and edges for single time step t parallel. The nature of BP algorithms do not allow us to do away with the for-loop over time and that remains the only sequential aspect of our implementation. Finally, each sample in our test set sequentially. It should be noted that we can further vectorize over a batch of samples in our test set. However, the memory required for DMP is Obb? E²) with b being the size of the batch and so memory requirements quickly blow up. Accordingly, we leave this aspect of implementation for future work.
salient feature of using DGL. that the message
use algorithm Io process
Algorithm 1: Dynamic Message Passing given graph G, snapshot O and time
for i E V do
set node to be PO
initialize node features and edge features with eq (12, 13) in DMP;
for (t 0; < t; = t 1) do
for e E E' do
perform message passing with eq (15, 16, 17) in DMP
for Merjel V do
L marginalize and update node states with eq (18, 19, 20) in DMP. Calculate P(Oli) with eq 21 in DMP.
return ì Panmnar;P(OPp
A.5 EFFECT OF VARYING NUMBER OF GCN-S LAYERS ON T ACCURACY
Fig. shows the t accuracy of PO of the GCN-S model for varying number of layers. We do not observe a significant effect coming from the number of layers. This may be due to the accuracy limitations with mmax and eycles affecting all the models equally, and superseding other effects such as the diameter of the graph. Another possible reason may be that the 20,000 samples on a graph of 1,000 nodes has many repetitions of the same PO, resulting both shallow and deep models memorizing patterns.
(a) BA-Tree
(b) BA-Dense
(c) Geometric
(d) ER-Dense
1.0
0.8
L=2
SIR
1ccura
0.6
L-5
0.4
L=10
0.2
5
13
17
21
25
29
1
9
13
17
21
25
29
5
9
13
17
23
25
29
1
5
9
13
17
21
25
29
Time Steps
Time Steps
Time Steps
Time Steps
Figure 6: Performance of GCN-S for SIR epidemic dynamics as top-1 accuracy over the test set with varying number of layers.
13

---

Under review as a conference paper at ICLR 2021
B THEORETICAL ANALYSIS
B.1 EARLY STAGE EVOLUTION OF SIR AND SEIR
The SIR equation on graph are
dS:
dR; 713,
dS;
dI.
dR;
3 Σ Asj
0.
(11)
dt
dt
dt
dt
dt
In very early stages, when & < 1/7 and E, I R & N, we have S, and have exponential for / because
dS;
dI,
dR:
dI:
-74
dt
dt
dt
dl
dl,
3 Σ
A - 71 = Σ (BA.) - -d6;))
dt
T1(1) Σ (exp(t (BA ,1))), (0)
(12)
Expanding this using the eigen-decoopossitinn = A yields eq. 29
B.2 TRANSITION PROBABILITIES
More generally, when the graph is weighted, the probability of susceptible node getting infected depends on Ai and the probability of node being in the infected state. For brevity, define (t) = P( = p)), with μ E {S, I..,, R}. The infection probability in SIR
2 can written as
P(x /| =) = II BA,p)) AApp 3 [A (
(13)
B.3 REACTION DIFFUSION FORMULATION
For brevity, define (1) = P(a = p). a network diffusion process the assumption that node can only be directly affected by state of node j if there is connection between them, i.e. if Aij 3 0. This restriction means that the general reaction- diffusson process on a graph has the form
Fa(A:)) Σ (aaa)((A))haal)))
(14)
11(1 + 1) =(A:(())) = ( ((Fa(A:P(()))))
(15)
With
halps) = oa (E WH.D + bi
(A))) = ((Ai)))). (A)
(16)
where 0(-) is the step function and a.(-) a nonlinear function. In regular diffusion graph, have two states S, and diffusion is changing the S / state. The probability Pej P(x /|z S) of node getting infected at + 1, given node j was in the infected state at time 1, can be expressed in the form of determined by the adjacency matrix Ag because node j can only infect its neighbors. The infection probability is given by pt (t + 1) = BAisp ( (1) and p for diffusion
p/. Hence,
fz =
91 (A) 3A,
(² = Σ 848.1%.
(17)
In regular diffusion there is no condition on the target node i and even if it is in the I state the dynamics is the same. In the model, however, the infection only spreads to if it is the state. Thus, we have to multiply the dynamics by p = P(x S) which yields
(1 1) (1p
(18)
14

---

Under review as a conference paper at ICLR 2021
This can still be written as by adding the extra functions
fz(z) =
92(A) I,
h2(p,)* - = 0505
(19)
and having
Pr(1 + =FF=))=F(A))())
(20)
where F, = fala ha) in 14 More complex epidemic spreading models such SIR and SEIR can also be written in a similar fashion. In SIR and SEIR the rest of the dynamic equations linear and do not involve the the graph adjacency A at all, meaning ga (a(A) I in the rest of the equations.
B.4 DISCRETE TIME AGENT- BASED SIR A REACTION DIFFUSION SYSTEM
The agent- agent-based models and (3,, which correct for double-coumting of infection from multiple neighbours, are sometimes written as
= -1-
(21)
where & is the total number of neighbors j of which are infected, meaning xf I. We will first
show that this is a special case of the form given in our paper. First, note that in 2 the terms can also be written
(1 -
II (1-
(22)
In the probabilistic model, we have to replace the strict condition of j being in the I state with its probability, so 62g.1 + P(r) = 1) = (t).
P(x = 1| = -1 II
(1-9A;1)
(23)
JEO.
and for small yield
P(att==1|p =) 3 Auph B² Σ AuppAAut 0(33)
(24)
which yields the simplified equation p:(t + 1)' = (t) (t). Note that if the infection
rate per time step large |3 Σ Ayp, can exceed rendering (24 inconsistent with p being
probabilities. Both (23) and (24) both can be written in the form of RD 15 and 9 We utilize the h1.91 and h2, found for diffusion (17) and SI (19)
FA(A:)) = Σ log (1 3A.hh1(p))))
= hha(p)
(25)
and defining the probability as
(1 + =F (1
exp [F2]]] = (1)
II (1 (1-[((-AAsp(())
3p Σ (t)
(26)
B.5 PROOFS
Proposition 2. Reaction Reaction-diffusion dynamics on graphs is structurally equivalent of the message- passing neural network ansatz.
15

---

Under review as a conference paper at ICLR 2021
Proof: Analyzing the full stochastic model requires closely tracking the individual events and varies in each run. Hence, we will work with mean -ield diffusion dynamics using transition probabilities, instead. Denoting Py (t) = P(z μ)) of node being in states such μ {S,I, j at time t, Markovian reaction-diffusion dynamics can be written as
F((As (o,,)))))))))))))))))
11)
ha(p))
Σ Wh a.pyy ba
(27)
where A = ((As))((A); with (()) being the step function and of-) nonlinear function. To see this, note that RD processes on graphs involve a message-paseee passing (MP) step (e.g. an infection signal coming from neighbors of a node), and reaction step where messages of different states / passed to node interact with each other on node RD dynamies such as SIR and SEIR models are also Markovian and the probability y (t) only depends on the probabilities at f 1. These are also the conditions satisfied by MPNN. In (27 AA are a set of propagation rules for the messages, which only nonzero where A is nonzero, same as the aggregation rule in MPNN. To have interactions between states f occurring inside each fixed node i, h(ps) can mix the states but not change the node index leading to the form of h(p)) in (27 which is the general ansatz for neural network with weight sharing for nodes, same as in MPNN, and graph neural networks in general.
B.5.1 PROOF OF THEOREM
Theorem (Time Horizon). Assume SIR dynamics ①T) on a connected graph of N nodes, starting with single patient zero. Denoting the adjacency matrix by A and its largest eigenvalue by X1, the average infection probability, both over nodes and choice of patient zero, (I(t)) = E ((t)/N)po becomes (1) after mmxx time steps given by
log N
321 Ro
Imax
(28)
y(Ro 1)
Proof:
Consider the spectral expansion A = X ((aa with A >... ...> Av- In early
in the disease spreading, all nodes are susceptible, meaning S. s 1, Rs 2 0, and 1, ~ S. Thus, combining the three SIR equations, keeping only 1, the infection spreads as Newman (2018)
T1(1) R
Σ exp [(BA -1),, (00) exp [BB1 -]]1 (1:11) I(0)) v v
(29)
Here, I is the identity matrix, x1 is the largest eigenvalue of A and (()) is the corresponding eigenvector. Averaging over uniform choice of patient zeros, for the average infection probability we have
(T(I) N exp (BA1 - 7)4 (wEO) -TO)
(1) Σ
exp [(32) - 7)4
(30)
N
where we used the inequality between L1 and La norms to get (4) I(0)) po | 1. Connectedness means A is irreducible and by the Setting the lower bound of (30) equal
|
>
theorem its leading
eignvector is positive, hence I w || 2 || and solving for t we get
B.5.2 PROOF OF THEOREM Prri
Proof: If PO is in triangle, we may miss it 2/3 of the times. Thus, the probability of detecting PO is bounded by P < Pers x 2/3, where P.r: is the probability that PO is in a triangle. Since edges in G are uncorrelated, each having probability p. G, is also a connected random graph with the same edge probability p. Hence, in G1 all nodes have degree k xe p|G Pori is one minus the probability
that none of the A neighbors of PO are connected, iec. Piri = (1 - (1) which proves the proposition.
16

---

Under review as a conference paper at ICLR 2021
Iidnggam noldegreen- bink (Drptral Groph with 184, 540 nodelt
Flerogram log Onggeggeee) 10 Nin (Orgiin Graph with nodess
120006
Frequency
Frequeeccy
-30
20000
log(degree)
log(degree)
(a) Degree distribution of original network
(b) Degree distribution of subsampled network
Figure 7: Degree distribution of the original co-location network with 384, 590 nodes and the subsampled network with 2, 689 nodes. We subsample the larger network to find subgraph in order reduce computational costs of our experiments. We observe that the distribution of our subsampled network is similar to the original graph.
C COVID-19 DATA AND SIMULATIONS
Geolocation data Mobility data provided by Cuebiq, a location intelligence and measure- ment platform. Through its Data for Good program (https:/ /ww:://wwwww.w cuwbiq. com/a boutt fata-for-good/ Cuebiq provides access to aggregated and privacy-enhanced mobility data for academic research and humanitarian initiatives. These first-party data are collected from users who have opted in provide access to their GPS location data anonymously, through GDPR-compliant framework. Additionally, Cuebiq provides an estimate of home and work census areas for each user. In order to preserve privacy, noise added to these "personal areas", by upleveling these areas the Census block group level. This allows for demographic analysis while obfuscating the true home location of anonymous users and preventing misuse of data.
Colocation network The method for constructing the co-location graphs is as follows. First, we split each day into five minute time windows, resulting in 288 time bins per day. For every location event, we use its timestamp to assign it to a time bin, then assign the longituee-latiudde coordinate of the observation to an 8-character string known as geohash. geohash defines an approximate grid covering the earth, the area of which varies with latitude. The largest dimensions of an 8-character geohash are 38m x 19m, at the equator. If a user does not have an observation for a given time bin, carry the last observation forward until there is another observation. We finally define two users to be co-located - and therefore to have a timestamped edge in the graph - if they are observed co-location graph is constructed by
in the same geohash in the same time bin. Accordingly. observing the greater Boston area over two weeks from 23 March, 2020 to April, 2020 and results in graph with N 384, 590 nodes. To reduce computational costs, sample a subgraph with N = 2, 689 nodes and [F| = 30, 376 edges with similar degree distribution and connectivity patterns as the original graph and can be observed in Fig 7
Epidemic simulations real data. We run a SEIR model on the real co-location network. In doing so, we select parameters and modify the structure of the model to resemble the natural history of COVID-19 Chinazzi et al. (2020). At each time step nodes, according their health status, can be in one of five compartments: S, E, Ia. or R. Thus, split infectious nodes in two categories. Those that are symptomatic (/) and those that are asymptomatic (Ia). The first category infects susceptible node, with probability À per contact. The second category instead with probability raA. We set ra = 0.5 and consider that probability of becoming asymptomatic as that is the sum of incubation (a 11 and infectious period (y- is set to be 6.5 days. Specifically, we fix a = 2.5 and -- 4 days. In a single, homogencously mixed, population the basic reproductive number of such epidemic model Ro (1 Pa raPa)8/7 where / the per capita spreading rate Keeling & Rohani (201T). Here however, the epidemic model unfolds on top of the real co-location network. Hence, infected nodes are able to transmit the disease only via contacts (with susceptible individuals) established during the observation period. As mentioned above, the value of Rg is defined by the interplay between the disease's 's parameters as well as the structural properties of the approximate = A where (k) is the average number of connections in the network. We obtain = 0.073 after solving for Ro and plugging in (k) 30376/ /2689 11.29. The simulations start
0.5. The generation time,
et al. (2015); Masuda & Holme (2017). For simplicity
17

---

Under review as a conference paper at ICLR 2021
with initial infectious seed selected uniformly at random among all nodes. We then read and store the anneeegggegte network memory. The infection dynamics, which are catalysed by the contacts between infectious and nodes, take place such network. The spontaneous transitions instead (i.c. transition from S to E and the recovery process), take place independently of the connectivity patterns. After the infection and recovery dynamics, we print out the status, with respect to the disease, of each node. Finally, we create a dataset with 10, 000 samples and an 80 - 10 10 train- vuin-vaiiationn-ess split.
18


Reviews:
[{"review_id": "yTBDeoAyfBg", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "An S(E)IR epidemics propagates on a graph, and the goal is to detect its source (P0) only from the observation of the state (S,E,I,R) of every node of the graph at some time $T > 0$. This version of the source detection problem has been studied first by Shah and Zeman (2011) for SI epidemics, as listed in Section 2. The current paper claims to (i) establish new fundamental limits on this problem, showing in particular that after some time the source detection becomes difficult, and (ii) to demonstrate the ability of graph convolutional networks to solve the problem and validate the results on real data.\n\n(i) The theoretical part consists in two short theorems, but their proofs are problematic. \n\n- Theorem 1 builds on an approximate analysis in Newman (2018), which requires a number of assumptions and approximations (for example, a mixing assumption that allows to replace expectations by ensemble averages in ordinary differential equations describing the evolution of $S$, $I$ and $R$; the limitation of time $t \\rightarrow 0$). These assumptions and approximation might not be valid on all graphs (e.g., a line graph has the largest eigenvalue around 2, but it is difficult to imagine that $O(1)$ nodes get infected in $\\Theta(log(N))$ time). The proof cannot rely on implicit assumptions and approximations whose error is not rigorously estimated. Maybe the theorem is true for some class of graphs including Erdős-Rényi graphs, but even in that case a more involved proof is needed, since the current techniques only work for early times $t$ and not for large times as in (5). Some computations could be clarified, for instance if $ \\langle \\psi^{(1)} \\cdot I(0) \\rangle $ denotes the average (over the uniform prior of patient zero over the $N$ nodes) of the scalar product between $\\psi^{(1)}$ and the one-hot vector $I(0)$, why is it equal to $\\lVert \\psi^{(1)} \\rVert_1 $ instead of\n$ \\langle \\psi^{(1)} \\cdot I(0) \\rangle = \\frac{1}{N} \\sum_{i=1}^{N} \\left(  \\psi^{(1)} \\cdot I(0) \\right)_i = \\frac{1}{N} \\lVert \\psi^{(1)} \\rVert_1$ ?\n\n- The proof of Theorem 2 is flawed. It starts with the statement that \"if P0 is in a triangle, we may miss it $2/3$ of the times.\" Why would all nodes of a triangle be equally likely to be the source? With the same rationale, why could not we argue that if P0 is in an edge, we may miss it $1/2$ of the times, which would give an upper bound of $1/2$ (as long as $|G_I|>2$) and would contradict the simulations results? The next statement that `\"in $G_I$ all nodes have degree $k \\approx p|G_I|$\" appears incorrect: suppose that the graph is an E-R graph $G(N,p)$, that $R_0$ is very large and that the infection spreads in a snowball way. Then most of the early infected nodes in $G_I$ have node degree $pN$, and not $p|G_I|$. \n\n- Now, the fact that P0's detection becomes harder over time is an important message, especially these days, but it is not a surprising result that when the infected set is a constant fraction of the population, then it is hard to detect P0. This difficulty was already reported in the initial paper by Shah and Zaman (2011).\n\n- In Section 3, the authors contend that compared to the SI model, the removed state introduces additional uncertainty about the temporal order of infections. Why? Since the state of each node is known, having 3 classes (S,I,R) instead of two (S,I) gives more information, which should ease the task of detecting P0. \n\n(ii) The comparison of the GNNs used by the authors over state-of-the-art message passing algorithms is made only with the DMP method of Lokhov at al (2014), but not with the (in general more accurate) belief propagation method of Altarelli et al (2014), also cited in Section 2. There is no comparison with the rumor centrality method developed by Shah and Zaman (2011) either, in terms of accuracy and speed. Also, it would be interesting to see the comparison at times other than $T=30$. \n\n- In terms of speed, training and inference should clearly be separated. It may be misleading to report in the abstract that GNNs are 100 times faster than state of the art methods: that applies only for inference. It would be more accurate to report that GNNs are 100 times faster for inference and twice faster for training compared to the DMP method (as well described in Section 5.1).  \n\n- It may be unrealistic to assume that all 4 states (S,E,I,R) can be detected for each individual; the exposed state in particular might be very hard to detect. Otherwise the simulations on real data seem well-done. It is unfortunate that with the non-interpretable theoretical results, the simulation results do not give much insight either. For example, it would be interesting to compare the accuracy results to the size of the infected set on the simulations. Maybe it would be more meaningful to normalize the rank by $|G_I|$ instead of $N$. \n\n- In Figure 2, how are the theoretic curves computed? Equation (7) in Theorem 2 depends on $|G_I|$, which is not directly linked to $T$ nor to the epidemic parameters. Is $|G_I|$ computed based on the simulation results? If so, why are then the confidence intervals not given for these curves? \n\n- In Figures 4 and 5, how is the size of the set $G_I$ evolving over time? \n\n- The paper should be proofread, it contains quite a few typos or vague statements, for instance: Theorem 3 in the appendix is actually Theorem 1 in the main paper; Figure 2 caption does not read well (no verb in the sentence : While accuracy drops below...);  Figure 4 caption: cycles significantly reduces accuracy of P0 -> cycles significantly reduce the accuracy of the detection of P0; Bottom of p4: where each edge has independent an probability $p$ -> where each edge has an independent probability $p$; etc", "extra_scores": {}, "extra_text": {}}, {"review_id": "nJYcF9CZdjP", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "\nThis paper proposes to use a graph neural network to infer the source of an epidemic in a network. Given a snapshot of the epidemic, the goal is to determine patient zero without full information of the mechanics of the epidemics but rather learning from historical data.\n\nOverall evaluation:\n\nMy evaluation is borderline. Although the application is timely and of interest, there seems to be little methodological novelty (simply an application of an existing GNN architecture) and, for being an empirical paper, the numerical experiments are not entirely convincing (more details below). The theoretical results add to the contributions.\n\n\nPros:\n\n1 - Timely application of GNNs to a network problem of interest.\n\nNeutral:\n\n1 - The theoretical statements add to the contributions, although quite simple in essence. For example, Theorem 2 is computing the probability that a given node belongs to a triangle in an ER graph. Proposition  1 seems to be more of an observation than a formal statement. Indeed, Reaction-Diffusion dynamics and MPNN are both non-linear processes on networks. It is unclear if Proposition 1 is saying anything more fundamental than this.\n\n\nCons:\n\n1 - No methodological innovation or domain knowledge embedded in the design of the architecture. The authors talk about “several modifications” but these seem to be more accessories (skip connections and batch normalizations) than incorporation of expertise in the architecture or the loss.\n\n2 - For a heavily empirical paper, the experiments are not comprehensive enough. Generalizability across graph types, real-data experiments (instead of simulated data on a real graph), learning from an epidemic mechanism and testing on another one, and learning with just a few observations (how realistic is to observe 20,000 epidemic spreads on the same graph?) are some of the areas that one would expect an empirical paper to cover and are missing in the current version of the manuscript.", "extra_scores": {}, "extra_text": {}}, {"review_id": "9PFUt76mOH0", "reviewer": "AnonReviewer5", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper studies the problem of source detection in an epidemics when one observes the underlying graph and a snapshot of the population at a given time i.e. who is infected or not infected. For a SIR (or SEIR) model, the authors propose to use GNN for this task. The learning procedure is then the following: given a fixed graph G, the authors create a dataset of snapshots by running a SIR on G.\n\nI am not convinced this problem should be solved with a machine learning approach. In most practical cases, we only have access to one snapshot for a given graph and learning is impossible. The authors here solve this issue by simulating many SIR processes but techniques like the one described by Shah and Zaman without any learning seem much more appropriate.\n\nThe authors should compare their results to the results obtained by Shah and Zaman.\n\n[No rebuttal given by the authors] Score unchanged.", "extra_scores": {}, "extra_text": {}}, {"review_id": "8me4farh9X", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Summary: Backtracking source of an epidemic (Patient Zero (P0)) is one of the important research topics of the current era that helps efficient resource allocation. Many of the existing works in this domain use graph-theoretic measures or message passing algorithms to tackle this problem. In contrast, this paper uses recently emerged Graph Neural Networks (GNNs) to learn and efficiently locate P0. It models disease spreading as a contagion process over a graph. While considering cyclic graphs, the paper shows upper bound on the accuracy of finding P0 and presents a bound of time horizon after which inference become difficult. Experimental results on different real-world and synthetic networks show interesting results. \n\nComments: The framing of P0 in terms of learning problem and the proposal of GNNs for solving it is a valuable contribution. Also, the theoretical bounds of time horizon and accuracy is quite interesting. Furthermore, Proposition 1 shows an interesting observation that Reaction Diffusion Dynamics is structurally like Message Passing Neural Networks. Such mapping may open new ways to study information diffusion over networks in a learning perspective. Overall, this paper is well-written and structured and friendly to reviewers. I would suggest the following comments for further improving the quality of the paper\n\n1.\tAlthough the experimental details are given in section 4.1 and Appendix A.1. However, it would be useful to add datasets tagging details that how nodes in each simulation were tagged for training GNNs.\n2.\tSince all the experimental results are based on the randomly chosen P0 and S(E)IR simulations parameters. I encourage the authors adding an ablation study. It would be more helpful to see results on different combination of parameters. Also, the statistical properties such as degree, the number of triangles it lies etc., of P0 would be helpful to see. \n3.\tTo ensure that GNNs work as intended, it would be helpful to see learning curves of GNNs on each type of random graph\n4.\t I encourage the authors to release their code and pre-trained models in a format that is easy to reuse for other researchers\n5.\tThe use of \\citet and \\citep has been mixed throughout the draft. The authors are encouraged to revise the citations accordingly. \n6.\tThe introduction section is lacking in terms of problem motivation. It’s mentioned that the problem is hard however, there is no discussion on the computational complexity of existing graph-theoretic methods and the challenges that need to be addressed. \n7.\tEq. 3 \\alpha is undefined. \n", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper introduces a GNN approach to solve the problem of source detection in an epidemics. While the paper contains some interesting new ideas, the reviewers raised some important concerns about the paper and so the paper should not be accepted in the current form. In particular,

- the paper does not motivate the ML approach to the problem
- the experiments are limited for an empirical paper
- the method used in the paper is not very novel
- the proofs presented in the paper are not formal enough

Decision: Reject