Construction of edge-ordered multidirected graphlets for comparing dynamics of spatial temporal neural networks
Joshua M. Roldan, Sebastian Pardo G., Vivek Kurien George, Gabriel A. Silva
CConstruction of edge-ordered multidirected graphlets forcomparing dynamics of spatial temporal neural networks
Joshua M. Roldan , Sebastian Pardo G. , Vivek Kurien George , and GabrielA. Silva ∗ Department of Bioengineering, University of California San Diego Department of Neurosciences, University of California San Diego Center for Engineered Natural Intelligence, University of California San Diego
Abstract.
The integration and transmission of information in the brain are dependenton the interplay between structural and dynamical properties. Implicit in any pur-suit aimed at understanding neural dynamics from appropriate sets of mathematicallybounded conditions is the notion of an underlying fundamental structure-function con-straint imposed by the geometry of the structural networks and the resultant latenciesinvolved with transfer of information. We recently described the construction and the-oretical analysis of a framework that models how local structure-function rules giverise to emergent global dynamics on a neural network. An important part of this re-search program is the requirement for a set of mathematical methods that allow us tocatalog, theoretically analyze, and numerically study the rich dynamical patterns thatresult. One direction we are exploring is an extension of the theory of graphlets. Inthis paper we introduce an extension of graphlets and associated metric that maps thetopological transition of a network from one moment in time to another at the sametime that causal relationships are preserved.
The integration and transmission of information in the brain are dependent on the interplay betweenstructural and dynamical properties across many scales of organization. This interplay happens inmolecular and diffusion interactions occurring in neurons and other neural cells, and at physio-logical scales in individual single cells, networks of cells, and eventually in networks of brain ∗ Corresponding author: G.S. Email: [email protected] a r X i v : . [ q - b i o . N C ] J un egions. A critical consideration towards a systems engineering view of the brain is understandinghow mathematically imposed constraints resulting from physical considerations determine neuraldynamics and what the brain is able to do.Implicit in any theoretical or computational pursuit aimed at understanding neural dynamicsfrom appropriate sets of mathematically bounded conditions, whether acknowledged or not, is thenotion of an underlying fundamental structure-function constraint. This fundamental constraintis imposed by the geometry of the structural networks that make up the brain at different scales,and the resultant latencies involved with the flow and transfer of information within and betweenfunctional scales. It is a constraint produced by the very way the brain is wired up, and howits constituent parts necessarily interact, e.g. neurons at one scale and brain regions at a higherscale. The networks that make up the brain, across all the various scales of organization, arephysical constructions over which signals and information must travel. These signals are subjectto processing times and signaling speeds (conduction velocities) that must travel finite distances toexert their effects - to transfer the information they are contributing to the next stage in the system.Nothing is infinitely fast. Furthermore, the latencies created by the interplay between structuralgeometries and signaling speeds is generally at a temporal scale similar to the functional processesbeing considered. So it matters from the perspective of understanding how structure determinesfunction in the brain, and how function modulates structure, for example, in learning or plasticity.Towards such a systems engineering understanding of the brain, on-going efforts in our groupare focused on deriving and studying the mathematical relationships and consequences that thisfundamental structure-function constraint produces. Central to the thesis of our work is the no-tion that the mathematical relationships we discover and prove about information integration andcomputation in biological neural networks should be independent, as much as possible, from un-necessary physiological details. Our objective is to identify and understand the most basic andsimple set of conditions that necessarily result, and to write them down in mathematical formsamenable to deep theoretical analyses. We are attempting to discover the fundamental algorithmsassociated with neural dynamics and neurobiological information representations. By ’unneces-sary biological and physiological details’ we mean that our theoretical constructions should retainonly features deemed essential to the algorithms themselves, while remaining as much as possibleindependent of details responsible for their implementation in the ’wetware’ environment of thebrain. The algorithms and mathematical relationships that underly them should be independent ofthe neurobiological specifics of any particular experimental model, for example.We recently described the construction and theoretical analysis of a framework derived fromthe canonical neurophysiological principles of spatial and temporal summation. This frameworkmodels the competing interactions of signals incident on a target downstream node (e.g. a neuron)along directed edges coming from other upstream nodes that connect into it in a network [4, 18].We considered how temporal latencies produce offsets in the timing of the summation of incomingdiscrete events due to the geometry (physical structure) of the network, and how this results in theactivation of the target node. The framework we constructed models how the timing of differentsignals compete to ‘activate’ nodes they connect into. This could be a network of neurons or anetwork of brain regions, for example. At the core of the model is the notion of a refractory pe-riod or refractory state for each node. This reflects a period of internal processing, or period of2nability to react to other inputs at the individual node level. It is important to note that we did notassume anything about the internal model that produces this refractory state, which could includean internal processing time during which the node is making a decision about how to react. In ageometric network temporal latencies are due to the relationship between signaling speeds (con-duction velocities) and the geometry of the edges on the network (i.e. edge path lengths). We haveshown that the interplay between temporal latencies of propagating discrete signaling events on thenetwork relative to the internal dynamics of the individual nodes - when they become refractoryand for how long - can have profound effects on the dynamics of the network [4], and the definitionof mathematical bounds on conditions required for a formal definition of efficient signaling [18].We have also shown that Basket cell neurons optimize their morphology in order to preserve whatwe call the refraction ratio, a balance between the temporal dynamics of individual nodes relativeto the dynamics of the entire network [14]. Our model allows us to compute and study the localdynamics that govern and give rise to emergent global dynamics on a network. This frameworkand its theoretical analysis are a concrete example of brain invoked algorithms that directly resultfrom the fundamental structure-function constraint, and are independent of any neurobiologicalor biophysical implementation details. We refer the reader to the relevant references for a fulldiscussion [18].An important part of this research program is the requirement for a set of mathematical methodsthat allow us to catalog, theoretically analyze, and numerically study the rich dynamic patternsthat result from network simulations using our framework. One direction we are exploring isan extension of the theory of graphlets, and the technical focus of the work in this paper. Ourmotivation and use of these results is their application to neural dynamics broadly speaking andthe analysis of numerical simulations from our framework applied to the study of the brain and tothe development of a new machine learning architecture we are constructing. But the work itselfif of even broader context and utility, with applications to network theory and the analyses of awide range of dynamic spatial-temporal networks, and to mathematics in the sense that we areextending the theory of graphlets. Specifically, we introduce an extension of graphlets that mapsthe topological transition of a network from one moment in time to another at the same time thatcausal relationships are preserved, including in situations where multiple signals contribute to theinitiation of a downstream activation event in a target node.The remainder of this paper is organized as follows: Section 2 discusses some background andother related work. It provides a short summary of previous work dealing with network motifsand graphlets. Additional information is provided in the Appendix. In Section 3 we give a briefintroduction to different types of networks and how graphlets are derived. Section 4 constructssynaptic signal graphs and introduces graphlets for multidigraphs and edge-ordered multidigraphs.In Section 5 we describe an extension of the graphlet-orbit transition metric to edge-ordered mul-tidigraphs. Section 6 discusses some concluding remarks. Milo et. al. were the first to observe that certain subgraphs, called motifs, occurred in networkswith a much higher frequency than in similar randomized networks [11]. Randomized networks3reserved the same number of edges coming in and out as the corresponding node from the realnetwork, but were assigned random connectivity. Motifs were considered the building blocks oflarger networks. Network motifs characterized classes of networks whose properties are defined byspecific types of elementary structures. Examples of networks motifs include feed- forward loops and
Bi-fans found in transcription networks and neuron synaptic connection networks , respectively(also see [11]).In 2004, Prˇzulj, introduced the concept of graphlets , which are small connected induced sub-graphs of a larger network. Graphlets are similar to network motifs in the sense that both analyzenetworks using subgraphs. However, graphlets are induced subgraphs, whereas network motifsallow partial subgraphs. Graphlets can be understood as a generalization of node degrees. In-stead of just looking at how many other nodes a given node touches, one asks how many triangles,chains, or stars structures does a node participate in. This notion better allows a better answer tothe question ’how many times does each subgraph show up in the graph?’, referred to as the the
Subgraph Counting Problem . Compared to network motifs, graphlets have an advantage in thatthey do not depend on generating corresponding randomized networks, and subgraphs are induced.Different randomizations can produce different network motifs, so consistency is an issue. Requir-ing graphlets be induced subgraphs rather than partial subgraphs yields increased network measureprecisions.Graphlets have found merit in their application to Protein-Protein Interaction (PPI) networks(see [10]). In 2007 Prˇzulj introduced a similarity metric called the graphlet degree distributionagreement (see [13]) which showed to be more efficient in comparing two static networks thannetwork motifs or other existing measures. Graphlets were first introduced for undirected net-works and then extended to directed networks (see [17]); the latter captured more informationand increased the precision by improving the degree of granularity used to describe a network.While we wanted to use graphlets to catalog and study causal signal patterns that resulted from ourmodel, the signal dynamics we need to analyze are intrinsically temporal. Thus, what we need isan inherently temporal or dynamic graphlet construction and associated metric in order to comparedifferent patterns.Milenkovi´c et al. made a first attempt to transform graphlets into a temporal framework inhopes of capturing more information than the original static version of the theory. In [8], Holmeand Saram¨aki defined a temporal network as a set of nodes and a set of events (temporal edges) thatare associated with two time parameters: a start time and a duration time. Any temporal networkcan be modeled as a sequence of snapshots . Each of these snapshots is a static graph whichaggregates the temporal information observed within a fixed time interval. The new temporalversion of graphlets was called dynamic graphlets by Milenkovi´c and colleagues (see [9]). Theseare equivalence classes of isomorphic time-connected temporal subgraphs. More precisely, twoinduced temporal subgraphs are equivalent if they are structurally isomorphic and have the samerelative order of events. Dynamic graphlets were shown to be superior relative to static and static-temporal approaches in their ability to perform network classification and node classification. Theadditional contribution of dynamic graphlets is the acquired ability to capture relationships betweensnapshots.But although dynamic graphlets are able to produce greater accuracy for temporal network4lassification than earlier methods, there persist a number of fundamental limitations. Dynamicgraphlets do not inform how the local topology changes over time, e.g., where poorly connectedsubgraphs become near-cliques. This notion of measuring how the local topology of a networkchanges over time was the motivation behind work by Ribiero, Silva and Aparicio who proposed acomplementary graphlet-based comparison metric of temporal networks by defining graphlet-orbittransitions . Graphlet-orbit transitions account for the possible change of states of a node appearingin different orbits (see [3]); e.g., a node in a periphery of a star moves to the center.For our purposes though, the framework needs to be extended further. Dynamic graphlets onlyallow one event or temporal edge configuring into paths. This makes it incapable of dealing withnetworks that have multiple (signaling) events or multiple edges. But this is critical to our modelbecause multiple signals from the same upstream node activated repeatedly can contribute to thesame running summation that triggers a single activation event in the downstream node it connectsinto. This implies that we need to represent multiple edges between connected nodes pairs, whichdynamic graphlets cannot do. To address this, we extended graphlet-orbit transitions to incorporatedirected edges and multi-directed edges.The relationship between dynamic graphlets and graphlet-orbit transitions is complementary.Dynamic graphlets are essentially causal paths obtained from jumping from one snapshot to thenext. Thus, dynamic graphlets inherently preserve causal order. On the other hand, graphlet-orbittransitions capture information within a snapshot and its transition to the next one, but it blurs themore granular causal information between snapshots, including the order in which the edges thatencode signals arrive due to the aggregation of the dynamics. Here, we take these two graphlet-based temporal metrics and combine them into one framework, so that the topological transitionfrom one snapshot to another can be characterized at the same time that causal relationships arepreserved, including in situations where multiple signals contribute to the initiation of a singledownstream activation event.
In this section we summarize the principle ideas and notation of graphlets and put them in contextwith standard graph theory. The mathematical notation was chosen to be as consistent as possiblewith the existing graphlet literature,and the theoretical construction of our neuro-derived dynamicsignaling framework.
Graphs consist of vertices (or nodes - we use these terms interchangeably here) and edges, wherevertices represent the components of a system and edges the relationships between these compo-nents. Formally, a graph G ( V, E ) is defined by a set of vertices V ( G ) , and a set of edges E ( G ) such that E ( G ) ⊆ V ( G ) × V ( G ) .One can always consider subsets V (cid:48) ( G ) ⊆ V ( G ) and E (cid:48) ( G ) ⊆ E ( G ) . Figure ?? illustrates two subgraphs of G . 5igure 1: A graph G with four vertices.Figure 2: Two subgraphs of G from Figure 1 of size 3.The size of a graph G is defined as the total number of elements in the vertex set | V ( G ) | . Sinceany subgraph is a graph by itself, we also say that the size of a subgraph H is n if | V ( H ) | = n . Wedenote this by H ( n ) . For example, the graph G in Figure 1 has size and both subgraphs are ofsize n = 3 , Figure 2.A special kind of subgraph are induced subgraphs. Induced subgraphs are obtained by ver-tex deletion (as opposed to edge deleted subgraphs constructed by deleting edges but preservingvertices). Given a graph G , let X be a set of vertices deleted from G . The resulting subgraph G − X is an induced subgraph. Note that we are actually interested in the subgraph constructed(i.e. induced) by the set of vertices that are not deleted, not the set of vertices that are. Explicitly, G [ Y ] = G − X is the subgraph induced by the vertices in the subgraph [ Y ] . The edge set of G [ Y ] are all edges in G that have both endpoints in G [ Y ] . For example, given G in Figure 3, H is aninduced subgraph of G constructed by vertex deletion of the two red vertices, but H (cid:48) is not inducedbecause it additionally contains an edge deletion not induced by vertex deletion. They are bothhowever, subgraphs of G .Figure 3: Give the graph G , H is an induced subgraph of G constructed by vertex deletion of thetwo red vertices. However, H (cid:48) is not induced because it additionally contains an edge deletion notinduced by vertex deletion. 6n mathematics, functions are used to relate two objects having the same structure in orderto transfer information. In the graph context, these functions have the feature of preserving theconnections between nodes, also called adjacencies. Given two graphs G and G (cid:48) , a graph homo-morphism f : G −→ G (cid:48) is a function sending the vertex set V ( G ) to the vertex set V ( G (cid:48) ) in sucha way that adjacencies are preserved. That is, if ( u, v ) ∈ E ( G ) , then ( f ( u ) , f ( v )) ∈ E ( G (cid:48) ) . Whenthe function f is also a bijective function, we say that f is a graph isomorphism , wherewith G and G (cid:48) are said to be isomorphic graphs. Intuitively, two graphs G and G (cid:48) are isomorphic if we can gofrom G to G (cid:48) (or vice versa) just by relabeling the vertices without changing their topology.A special case of isomorphism arises when the function f has the same domain and range,that is, when f is an isomorphism from G to itself. These particular functions are called automor-phisms . Automorphisms play a key role in classifying the local topological properties of a vertex.The set of all automorphism under composition defines a group called the automorphism group of G which we denoted by Auto ( G ) (see [5]). The automorphism group Auto ( G ) induces a partitionof G into equivalence classes where two vertices u and v are equivalent if there exists an automor-phism f such that f ( u ) = v . Vertices of the same equivalence class are said to be in the same automorphism orbit or orbit for short. Thus, two vertices in a graph are topologically identical ifand only if they belong to the same automorphism orbit. A similar concept is edge automorphismorbits. One can utilize small connected non-isomorphic induced subgraphs called graphlets andtheir corresponding orbits to characterize the vertices of a larger graph, and hence the larger graphitself. Our work in this paper focuses on vertex graphlets and leaves edge-graphlets for futurework.From a theoretical perspective, the most informative way of applying graphlets to any graph G is to enumerate the induced subgraphs so we can describe the local topology. For example, Figure1 in [10] shows all the induced undirected graphlets up to five nodes, which we include in theAppendix as Figure 13. Observe that within each of these graphs G s , some nodes are colored thesame; this represents that they belong to the same orbit. For instance, from Figure 13 one can seethat G has two orbits. This follows from the fact that G has only two automorphisms, namely, Auto ( G ) = { the identity function, the permutation α of the end points } .Hence, if a, c denote the two endpoints and b the only center point, the orbits are O = { a, c } and O = { b } . Here, the endpoints are colored with black and the center node is colored with white.Likewise, G has only one orbit because there is always an automorphism which sends one vertexto any other. Thus, the vertices are all in the same orbit O signifying they all have equivalenttopological positions.In [13], Prˇzulj introduces a graphlet-based similarity metric for undirected graphs called GraphletDegree Distribution Agreement (GDA) which was first tested in [10] to study protein-protein In-teraction networks. In that work the authors showed that the local network structure was related tospecific biological functions. The GDA is constructed as follows: given a node v of a network G ,we account for the collection of all different graphlets G s in which vertex v belongs to. For eachstep, we indicate the orbit or topological position within the graphlet G s the vertex v is at. Then,we can construct a vector called the Graphlet Degree Vector (GDV) which has dimension equal tothe sum of all the orbits considered, i.e, the number of orbits in graphlet G plus the number oforbits in graphlet G , etc. The entries of the vector reflect the number of times the chosen vertex v i th orbit within G s . The graphlet degree vector of the vertex v , GDV ( v ) , containslocal topological information about this vertex. As an example, consider the graph G in Figure 4.The vertex v that is colored blue and the different orbits contained in the graphlets G , G , G aredrawn.Figure 4: The 4 different orbits within the graphlets G , G , G . Notice that each orbit is coloredby a different color.Thus, the graphlet degree vector associated with the blue vertex is: GDV ( v ) = Repeating the process for for every vertex in the graph, one can construct a matrix whoserows are all the graphlet degree vectors. Furthermore, this matrix can be transformed into anothermatrix by computing the s th -graphlet degree distribution followed by the arithmetic mean of allagreements over all the orbits considered. This matrix is called the graphlet degree distribution . Itgives the number of nodes participating in graphlet G s k times, and has the property of storing notonly the frequencies, but the distribution of the graphlets. For undirected graphs, the tuples ( u, v ) and ( v, u ) represent the same edge, meaning that the orderof the pair ( u, v ) is irrelevant. But when the order matters, in the sense that there is a directionalityto the transfer of information between the vertex pair, then ( u, v ) and ( v, u ) represent differentedges, each one of these representing a different direction. Graphs with this edge property arecalled directed graphs or digraphs . Directed edges represent a way of capturing asymmetric in-formation flow in a network. For example, one can model metabolic reactions in Eukaryotes as anenzyme-enzyme network in which two enzymes-coding genes are connected by a directed line ifthe first enzyme catalyses a reaction whose product is a substrate for a reaction catalysed by thesecond enzyme; a biological example where graphlets have been used (see [17]). Digraphs are anappropriate model for biological neural networks. For example, vertices can represent neurons anddirected edges represent axons and the directional propagation of action potentials between two8eurons, or vertices can be brain regions connected by edges that represent white matter tracts, oreven as a model of a neuron itself from the soma to axonal arborizations with vertices representingpositions along the axon and edges representing axonal segments [14].Figure 5: Directed and bidirected graphs.Two digraphs G and G (cid:48) are isomorphic if there exists an isomorphism f on their underlyinggraphs which preserves the direction of the edges. The automorphism group of a digraph G issimilar to an undirected graph, except now the group is smaller due to the directionality restrictions.This means the orbit sizes, the number of vertices with identical topology, is smaller. Inversely,the number of orbits is greater since there are more different types of vertex topological positions.Nevertheless, the notion of graphlets and metrics in the directed case are still attainable. Given theset of all orbits for directed graphlets of size up to four nodes, Prˇzulj, et. al. extended the GraphletDegree Distribution Agreement (GDA) for directed graphs, showing its efficiency and superiorityover existing network measures (see [17]). Networks can also used to model temporal dynamics. These networks are interchangeably referredto as temporal networks or dynamical networks . A temporal network comprises a set of nodes V and a set E of events (temporal edges) that are associated with a start time and duration (see [9]).Temporal networks who’s dynamic are dependent on their spatial geometry in addition to theirconnectivity form a class of spatial-temporal networks that we have theoretically studied as mod-els of biological neural networks [4, 18]. One way to study global dynamics is to aggregate all thenodes and edges from the temporal information into a single static graph. Another strategy is torepresent the evolving dynamics by a sequence of consecutive network snapshots , each of whichis a static graph gathering all the temporal data observed during a time interval T W . Static net-work metrics, such as GDA, have been used with temporal networks to measure how the topologychanges over time by comparing consecutive snapshots. This is a static-temporal approach. Nev-ertheless, a static-temporal approach, and hence GDA, can toss out temporal information which iscentral for understanding the dynamic’s evolution. In both cases though, the major shortcoming ofthese approaches is that they lose information about how the dynamics evolved over timeIn [3], Aparicio, Ribiero, and Silva introduced a new metric for temporal networks called the orbit-transition agreement (OTA) which considers how local topologies change with time. Theyshowed that OTA was able to show improved clustering accuracies of sets of temporal networkswhen compared to static network motifs and static graphlets.9 Synaptic Signal Graphs
The novel contribution of this paper is an extension of the theory of graphlets in order to allow usto analyze the causal evolution of dynamic patterns formed by our competitive refractory dynamicsmodel. In the Introduction we briefly gave a high level overview of the framework. We refer thereader again to [18] for details on the full development of the model and its theoretical analyses,and [4] for background and results on the effects of network geometry on neural dynamics. In thissection we extend this prior work in order to accommodate the graphlet analysis.
Consider a discrete signal σ sent from some start node s ( σ ) . This signal traverses down its ge-ometric edge until it reaches its target node t ( σ ) . Each signal’s ”lifetime” can be described by afunction D ( ρ, ψ, α ) . These three parameters will facilitate deciphering the signal’s causal effectwithin each of its states, namely if any effect exists. Figure 6 shows the four possible states of asignal. Terminal stateSynaptic state Traversal stateWaiting state
Figure 6: Possible dicrete signal states along an edge connecting two vertices (nodes) in a network.During the ’Waiting state’ a signal has not yet been released by the initiating vertex. The ’Traversalstate’ reflects a signal in transit along the edge. The ’Active state’ reflects the arrival of the signalat its downstream target vertex and a subsequent effect. In the ’Terminal state’ the signal has hada past effect.We consider the first state as a waiting period or waiting state where the signal has not yetbeen sent. The signal necessarily has no effect on the dynamics. The second state begins when thesignal’s start node activates and fires, thus sending the signal towards its target node. We denotethe time when σ is sent from s ( σ ) by ρ ( σ ) . This state is characterized by having σ coming downthe edge. We call this σ ’s traversal state . Observe that the traversal state is bounded by the timewhen σ reaches its target node, which we denote by ψ ( σ ) . Also, during the traversal state, σ doesnot cause any node to activate nor contributes any weight to the target’s running summation. Inother words, traversing has no effect on any nodes, and therefore neither on the dynamics. It isonly when the signal reaches its target node at time ψ ( σ ) does it have an effect and contribute aweight to the running summation. If a signal reaches its target node but the node is refractory, thenthis signal is dropped, which is why its important to also mark σ as entering its synaptic state -10nalogous to the arrival of an action potential at a synaptic terminal in a biological neuron. For therest of the paper, we assume that we are only dealing with signals which do enter their synapticstate, since dropped signals do not have causal effects.When a target node’s running summation surpasses its activation threshold, at time α ( σ ) , thenthe signal weight of σ becomes irrelevant. This is because the target node becomes refractory andresets its running summation to zero, as required by the model. At time α ( σ ) , the signal enters its terminal state . Below we define the possible states of a signal given an observation time T O in thecontinuous time interval [ a, b ] . S T O ( σ ) = waiting if T O < ρ ( σ ) traversal if ρ ( σ ) ≤ T O < ψ ( σ ) synaptic if ψ ( σ ) ≤ T O < α ( σ ) terminal if α ( σ ) ≤ T O We want to investigate the interactions between signals and their collective effect for activatingnetwork nodes. We define the signal set Φ [ t a , t b ] = { σ | t a ≤ ρ ( σ ) ≤ t b } . Φ [ t a , t b ] is the set of all signals sent from any network node within the time interval [ t a , t b ] ⊆ R .By definition, it is clear that if [ t a , t b ] ⊆ [ t a , t c ], then Φ [ t a , t b ] ⊆ Φ [ t a , t c ]. Because there can onlybe a monotonic increase in the set of signals sent within the network as time passes, we can writethe relation Φ[ t a , t c ] = Φ[ t a , t b ] (cid:116) Φ( t b , t c ] (1)In the event that Φ [ t b , t c ] is an empty set, the signal set remains constant from Φ [ t a , t b ] to Φ ( t b , t c ].For our purposes, we notice that the signal set undergoes a change or a transition when a nodeactivates, because the next event involves this node sending signals to its downstream nodes. Also,we assume that there are finitely many time points before a node activates, hence finitely manytransitions in the signal set. We can therefore enumerate these node activation times in increasingorder as { t , t , . . . , t N − , t N } . Let A be an order preserving map A : N −→ R such that A ( j ) = t j for j ≤ N and A ( j ) = t N for j > N . We call the sequence { A ( j ) } j ∈ N : = { A j } j ∈ N the activationseries and a particular A j an activation-step . As a result, we can construct the sequence of signalsets Φ ⊆ : Φ[0 , A ] ⊆ Φ[0 , A ] ⊆ Φ[0 , A ] ⊆ · · · ⊆ Φ[0 , A N − ] ⊆ Φ[0 , A N ] ⊆ . . . (2)This construction preserves all the pertinent causal information in the dynamics, while simulta-neously supporting what will be the requirements of a graphlet-based temporal analysis. Recallthat the definition of a snapshot-based representation introduced in [9] and [3] requires snapshotsto be equally spaced apart by ∆ t . However, here the set of all activation-steps do not necessar-ily need to be separated by a fixed time interval. Thus, we generalize our window sizes to allowthem to change rather than keeping them fixed. Specifically, we chose to take such generalizedsnapshot-based representations of signal dynamics at each activation-step.11efore continuing, we make an important observation. Suppose we had two identical sequencesof signal sets Φ ⊆ and Φ ⊆ . It is not necessary that the actual signal dynamics are the same. Todemonstrate this point, take two signal dynamics D and D , along with their corresponding se-quences such that Φ ⊆ ( D ) = Φ ⊆ ( D ) . We can think of these two signal dynamics as films thatcan be played at three different speeds: slow, normal, or fast. Whether we watch the same film atslow, normal or fast, we still see the same succession of scenes or actions. Thus, the similarity orinvariance between the films is the sequence of scenes and its difference is their respective timelapse . Hence, if two signal dynamics D and D have the same sequence of signal sets, we saythe signal dynamics are equivalent up to time lapse. This is a byproduct of interrogating causalrelationships embedded in possibly variable dynamics. In this section we describe a snapshot-based representation of the competitive refractory dynamicmodel. We do this by constructing synaptic signal graphs G Ψ to play the role of snapshots. In thisway, we will be able to represent signal dynamics as a sequence of synaptic signal graphs.Given a time interval [ t c , t d ] ⊆ [ t a , t b ] , we define the synaptic signal set Ψ[ t c , t d ] = { σ ∈ Φ[ t a , t b ] | t c ≤ ψ ( σ ) ≤ t d ≤ α ( σ ) } .By construction, it follows that Ψ[ t c , t d ] ⊆ Φ[ t a , t b ] . Using the activation series { A j } j ∈ N we canconstruct the sequence of synaptic signal sets: Ψ[0 , A ] , Ψ[0 , A ] , Ψ[0 , A ] , . . . , Ψ[0 , A N − ] , Ψ[0 , A N ] , . . . Let us take a look at a small subsequence: . . . ,
Ψ[0 , A j − ] , Ψ[0 , A j ] , Ψ[0 , A j +1 ] , . . . The synaptic signal set
Ψ[0 , A j ] can be partitioned in two ways, both resulting in two disjoint sets.The first partition classifies signals via their relation to the node(s) activated at time A j . Let theset Ψ act [0 , A j ] consist of all the signals whose target node activates at A j , while the second set Ψ cact [0 , A j ] has all the other synaptic signals whose target node does not activate. Thus, we havethe disjoint union Ψ[0 , A j ] = Ψ act [0 , A j ] ˙ (cid:116) Ψ cact [0 , A j ] . (3)A second partition of Ψ[0 , A j ] is due to the signals in Ψ cact [0 , A j − ] which remain synaptic inthe next activation-step, whereby they are a subset in Ψ[0 , A j ] . Hence, Ψ[0 , A j ] = Ψ cact [0 , A j − ] ˙ (cid:116) Ψ new ( A j − , A j ] , (4)where Ψ new ( A j − , A j ] is the set of new synaptic signals. That is, Ψ new ( A j − , A j ] = { σ | A j − < ψ ( σ ) ≤ A j ≤ α ( σ ) } .12s a consequence, for the subsequence we have: Ψ[0 , A j +1 ] = Ψ cact [0 , A j ] ˙ (cid:116) Ψ new ( A j , A j +1 ]= (cid:16) Ψ[0 , A j ] \ Ψ act [0 , A j ] (cid:17) ˙ (cid:116) Ψ new ( A j , A j +1 ]= (cid:16)(cid:0) Ψ cact [0 , A j − ] ˙ (cid:116) Ψ new ( A j − , A j ] (cid:1) \ Ψ act [0 , A j ] (cid:17) ˙ (cid:116) Ψ new ( A j , A j +1 ] (5)Intuitively, this means that there is subset from Ψ cact [0 , A j − ] which remains in Ψ[0 , A j +1 ] .From any synaptic signal set Ψ[0 , A j ] , we can create a snapshot which is a graph that encodesthe respective set of synaptic signals. We will refer to these graphs as synaptic signal graphs andwrite them as G Ψ [ A j ] . The edge set of any G Ψ [ A j ] corresponds to the synaptic signals in Ψ[0 , A j ] and the vertex set is comprised of all the network nodes.The above synaptic signal set decompositions have a correspondence to their synaptic signalgraph counterparts. It is possible for there to be subgraphs of G Ψ [ A j − ] which are embedded into G Ψ [ A j +2 ] . This is because synaptic signals may remain synaptic through multiple network nodeactivations without their downstream nodes activating. To achieve this, we let H Ψ [ A j ] correspondto the synaptic signal subgraph of the synaptic signal subset Ψ cact [0 , A j ] ⊂ Ψ[0 , A j ] .In the next section we how how synaptic signal graphs are represented by multidigraphs , whichallow us to accurately express a full description of the dynamics by writing down their edge-ordering relating to the order in which the signals became synaptic. A graphlet-based analysis thenfollows by constructing the appropriate graphlet types. In some cases, a node u contributing to the running summation of a node v could could consistof multiple synaptic signals contributing towards a single activation event. This happens when theupstream node u activates multiple times between two consecutive activations of node v such thatat least two signals of u become synaptic before the second activation of v . In order to preserve theproperty of synaptic signals multiplicity, we use multidigraphs. Thus, we will write each synapticsignal graph as a multidigraph.To facilitate notation, for the remainder of the w papere will denote a synaptic signal multidi-graph at A j as M Ψ [ A j ] . Again, for each M Ψ [ A j ] the vertex set V ( M Ψ [ A j ]) equates to the networknodes whereas the edge set corresponds to the set of all synaptic signals in Ψ[0 , A j ] . In this context,there will be directed edges ( u, v ) of degree d in E ( M Ψ [ A j ]) if and only if there exists d synapticsignals σ i from u to v with ≤ i ≤ d and A j < ψ ( σ i ) ≤ A j +1 < α ( σ i ) for i ∈ { , , . . . , d } .We will write the degree of directed edges as d u,v where u denotes the start vertex and v the targetvertex of these edges.Given a graph representation by M Ψ [ A j ] for each snap shot, we can now analyze it by a graphletdecomposition. However, existing graphlets do not capture information about multiple directededges. Rather than make an enumeration of the graphlets (for each of any occurring multipleedges), we will take a different approach. We will construct a vector space as a means to describeand extend directed graphlets to directed multigraphlets in the same framework. We first start byviewing directed graphlets as linear combinations of directed edges, and then generate directed13ulti-graphlets or multidigraphlets by assigning coefficients to these linear combinations whosevalues are equal to the degrees of the multiple directed edges. The advantage of having this al-gebraic construction of multidigraphlets from directed graphlets is two-fold. First, it is a naturalmapping from the automorphism orbits of directed graphlets to those of multidigraphlets. Thisnot only gives a relationship between the orbits, but also topological relationships. The secondadvantage is it provides a fine-tuned enumerating process of multidigraphlets. We begin by describing how to algebraically represent and construct directed graphlets via theiredges. This construction is based on the concept of elementary paths as in [7]. We adapt thefollowing definition to our framework for a multidigraph M : Definition 4.1. [7, Definition 2.2] Given a positive integer r , an elementary r -path in a multidi-graph M is a non-empty sequence a , a , · · · a r − of edges in M such that t ( a i ) = s ( a i +1 ) for i ∈ { , · · · , r − } . We denote this r -path by p = a a · · · a r − , whose start vertex s ( p ) = a andtarget vertex t ( p ) = a r − . Notice that for r = 0 , the -elementary paths have the form p = v , where v ∈ V ( M ) is somevertex of the multidigraph M and when r = 1 , the -elementary paths are comprised by all thedirected edges in E ( M ) . Let K be an integral domain and let us denote by Λ the set of all linearcombinations of -elementary paths and by Λ the set of all linear combinations of -elementarypaths in M over K respectively. Graphically, an element in Λ translates to a finite set of isolatedvertices, whereas directed and multidirected subgraphs in Λ are a linear combination of thosevertices. Thus we define Λ = Λ ⊕ Λ , whose elements are all linear combinations generated byboth -elementary paths and -elementary paths. As we will see, all possible states in which amultidigraph can be found can be represented by a linear combination in Λ .According to the directed graphlets enumeration in [1], graphlet G corresponds to a singledirected edge from node v to node v . Thus, its algebraic representation is v , v ) . So, if westart with the vertex set V = { v , v } , the graphs G , G with edge sets E = { ( v , v ) } and E = { ( v v ) } , respectively, both belong to the graphlet type G . In the first case v belongsto orbit O , and v belongs to orbit O . The second case is the reverse. For graphlet type G (abidirectional edge) the representation is v , v ) + 1( v , v ) , where v , v belong to the same orbit O . Suppose now that we have a graph G whose edge set consists of d edges ( v , v ) . We representthis multiplicity with the coefficient d v ,v ∈ N , such that G is expressed by d v ,v ( v , v ) . Theorbits of this graph are closely linked to the orbits of G , but could be different. Likewise, if wehave a graph G whose edge set had d edges ( v , v ) and d edges ( v , v ) , the orbits will have thesame vertices as elements. However, if d v v (cid:54) = d v v , then the vertices will not be in the same orbitbecause the symmetry breaks as a result of the different degrees of the directed multi-edges. Thus,the graphlet type is expressed as d v v ( v , v ) + d v v ( v , v ) with v and v belonging to distinctorbits. Note that the expression is unique up to reordering. These algebraic expressions can becombined to express any multidigraph, and in particular any multidigraphlet type. This exampledemonstrates the ability of graphlets to be combined together to produce new graphlets.14igure 7: Resulting multidigraphs after applying the node equivalence relation between vertices ofa DAG. Definition 4.2.
We say that two multidigraphs M and M are isomorphic if there exists a bijection f : V ( M ) −→ V ( M ) and a bijection g : E ( M ) −→ E ( M ) such that for every directed multi-edge ( u, v ) ∈ E ( M ) , we have that g ( u, v ) = ( f ( u ) , f ( v )) ∈ M and d u,v = d f ( u ) f ( v ) . In other words, any multidigraph isomorphism is an isomorphisms between the underlyingdirected graphs which preserves the edge multiplicities. In a follow up paper we will introducethe use of this construction as part of an algorithm that will allow us to build up dynamic patternsproduced by the competitive refractory model.
In the same way that the automorphism orbits for directed graphs decreases their number of ele-ments, the vertex orbits of multidigraphs have fewer elements in each orbit equivalence class, butmore classes in total. This is because most of the automorphism groups are trivial; that is, theyonly contain the identity map. Consequently, we lose a lot of the symmetry associated with undi-rected and directed graphlets. Nevertheless, despite the loss of symmetries, orbits can still be usedto indicate which topological position a vertex or edge belongs to within a multidigraphlet.Consider the snapshots taken at A (4) and A (5) in Figure . Since the other snapshots are merelydirected graphs, their graphlets are also directed. For this reason, we work through M Ψ ( A ) and M Ψ ( A ) to explain the new case of multidigraphlets.Figure illustrates the two and three node graphlets of M Ψ ( A ) and M Ψ ( A ) , respectively. Thealgebraic expressions for the two and three nodes graphlets of M Ψ ( A ) can be written algebraicallyas mdG = { ( P, B ); B + R ; 2( R, G ); P + G } mdG = { R + ( P, B ); B + 2( R, G ); P + 2( R, G ); G + ( P, B ) } a) Graphlets of M Ψ ( A ) . (b) Graphlets of M Ψ ( A ) . Figure 8: The two node and three node decomposition of M Ψ ( A ) and M Ψ ( A ) .And the graphlets of M Ψ ( A ) are mdG = { P, B ) + (
B, P ); B + R ; R + G ; P + G } mdG = { R + [2( P, B ) + (
B, P )]; R + B + G ; R + P + G ; G + [2( P, B ) + (
B, P )] } . We introduce a three parameter characterization to catalog submodule spanning multidigraphlets.Let V = | V (cid:48) ( M Ψ ) | be the number of vertices a multidigraphlet has, E = | E (cid:48) ( M Ψ ) | be the number ofnetwork edges selected, and F the number of synaptic signals between the network nodes incidentto the selected edges. Given values for V , E , and F , we can output a set of multidigraphlets M ( V , E , F ) . For example, suppose V = 3 , E = 4 , and F arbitrary. The possible underlyingdirected graphlets (up to isomorphisms) connecting these three vertices are shown in Figure . Thealgebraic expressions for these four multidigraphlets representatives then are given by1. d RB ( RB ) + d BR ( BR ) + d RG ( RG ) + d GR ( GR ) d BR ( BR ) + d GR ( GR ) + d BG ( BG ) + d GB ( GB ) d RB ( RB ) + d RG ( RG ) + d BG ( BG ) + d GB ( GB ) d RB ( RB ) + d BR ( BR ) + d RG ( RG ) + d GB ( GB ) However, it is important to note that the sum of the coefficients d i,j of each expression is alwaysequal to F : 16igure 9: The four possible four underlying directed edges between three nodes which to span theclass of multidigraphlets M (3 , , F ) . The sum of all the edge coefficients must add up to F . The competitive refractory dynamic model treats weight contributions triggered by arriving signalsfrom upstream nodes such that the maximum value of an edge weight to the running summationresults at the time of arrival of the signal, with progressively decaying contributions at subsequenttime steps ( c.f. equations ?? and ?? above). The physiological analog of this process is the decay ofpost-synaptic potentials as a function of the space and time decay constants due to the membranebiophysics. This reflects a critical algorithmic computational component of the neurobiology.Because of this, it is important to be keep track of the order in which the signals becomesynaptic. This motivates an edge-ordering on the multidigraphs M Ψ . As an example, figure compares two cases with the same signals, but shows how the relative order affects how a nodeactivates.The sequence of signal sets is our starting point. Referring back to sequence ( ), for each signalset Φ[0 , A j ] , we apply a total order relation ≤ based off the times they became synaptic, that is,for any σ k , σ l ∈ Φ[0 , A j ] , then σ k ≤ σ l if ψ ( σ k ) ≤ ψ ( σ l ) . We will denote each totally orderedsignal set as Φ ≤ [0 , A j ] for each j ∈ N . Now, given Φ ≤ [0 , A j ] , we will denote by Ψ ≤ [0 , A j ] theinduced totally ordered synaptic signal set. This new condition allows us to re-describe the twosets as monotonic increasing sequences P (Φ ≤ [0 , A j ]) and P (Ψ ≤ [0 , A j ]) . In fact, for the sequence P (Φ ≤ [0 , A j ]) = σ , σ , . . . , σ N , the corresponding synaptic signal sequence P (Ψ ≤ [0 , A j ]) is a sub-sequence σ k , σ k , . . . , σ k S , where N ≥ S are the number of signals in Φ ≤ [0 , A j ] , and Ψ ≤ [0 , A j ] ,respectively [16]. An important observation is that the last synaptic signal in Ψ ≤ [0 , A j ] refers tothe signal(s) which caused its target node(s) to activate.Let us consider now a multidigraph M . We can give two notions of an edge-ordering of M from two perspectives. An absolute edge-ordering is defined as a bijection p abs : E ( M ) −→{ k , k , . . . , k | E ( M ) | } such that k i ≤ k j whenever i < j . A relative edge-ordering is a bijection p rel : E ( M ) −→ { , , . . . , | E ( M ) |} [6]. Observe that one can go from an absolute edge-orderingto a relative one, just by sending k i to i . Now, given a sequence of edges P = e , e , . . . , e l we A total order on a set S is a binary relation ≤ such that the following properties hold: Reflexivity, Antysimmetry,Transitivity and Comparability. P is monotonically increasing sequence of length l if it is a sequence in M such that p ( e i ) ≤ p ( e j ) for all i < j . We denote an edge-ordered multidigraph as an ordered pair ( M, p ) ,where M is a multidigraph and p is an absolute or a relative edge-ordering.To construct an edge-ordered synaptic signal multidigraph ( M Ψ [ A j ] , p ) , we first notice thatfor an absolute edge-ordering, p − abs ( k i ) = ( s ( σ k i ) , t ( σ k i )) ∈ E ( M Ψ [ A j ]) ; that is, the edges getassigned the absolute ordering of the signal. In contrast, for a relative edge-ordering, p − rel ( i ) =( s ( σ k i ) , t ( σ k i )) ∈ E ( M Ψ [ A j ]) - the edge-ordering preserves the relative order in which the signalsbecame synaptic. This is a crucial difference compared to dynamic graphlets. Dynamic graphletshave a relative edge-ordering which represents consecutive causal events. Dynamic graphlets areessentially causal paths where the sequence of edges must have vertices incident to consecutiveedges. In our construction, we loosen this requirement. In fact, no causal relations exist within anedge-ordered multidigraph. This is because every time a signal activates its target node, the repre-sentative edge is removed in the next snapshot since the signal is no longer synaptic by definition.As a consequence, there exists no causal sequence within any ( M Ψ [ A j ] , p ). The causal sequencesare seen between different snapshots ( M Ψ [ A j ] , p ) and ( M Ψ [ A j + c ] , p ) where c > . Definition 4.3.
An edge-ordered multidigraph isomorphism is a multidigraph isomorphism f : V ( M ) −→ V ( M ) such that for the edge-orderings p abs of M and p abs of M , we have p abs ( xy ) = p abs ( f ( x ) f ( y )) . In other words, any edge-ordered multidigraph isomorphism preserves both the topology andthe absolute edge-orderings. This notion lead us to define edge-ordered multidigraphlets as equiv-alence classes of isomorphic edge-ordered multidigraphs where equivalence is with respect to therelative edge-orderings. For instance, consider the edge-ordered multidigraphlets in Figure .Note that the graphlets’ edges carry a relative ordering, not an absolute ordering.Similar to how we created a free K module Λ space for multidigraphlets based on the threeparameters N , E , and F , we will characterize a collection of edge-ordered multidigraphlets in asimilar way. We can consider the relationship between multidigraphlets and edge-ordered mul-tidigraphlets similar to combinations and permutations. For multidigraphlets, the order does notmatter, but for edge-ordered multidigraphlets it does. In this last section we extend a graphlet-based similarity metric in order to compare synaptic signalgraphs and the transitions between the them. We focus on the orbit transition agreement since itprovides a comparison of different temporal networks [3]. We specifically address the question’how do the signal dynamics of networks evolve?’.Related work is exploring what the topologies are for each of the synaptic signal graphs corre-sponding to signal propagations between node activations. We will also analyze the frequency ofgraphlets in order to understand if certain classes of graphlets appear more frequently than others.19 a) Edge-ordered multidigraphs. We suppress the multiple edges to conserve space.(b) Graphlets of ( M Ψ [ A ] , p ) . (c) Graphlets of ( M Ψ [ A ] , p ) . Figure 11: The two and three node edge-ordered multidigraphlet decomposition of M Ψ [ A ] and M Ψ [ A ] . The evolution and time-dependence information of a temporal network is efficiently captured withthe orbit transition agreement (OTA) metric, more so than other metrics based on taking frequen-cies and distributions of graphlets within a network. This is because other metrics do not con-sider the different transitions of the possible states of orbits within a temporal network. Similaritybetween two graphs G and G is given by the average similarity of their graphlet-transition fre-quency for each graphlet-transition.We now extend the orbit transition agreement in [3] for synaptic signal graphs. After classifyingall the graphlets in G Ψ , we can make a matrix. The matrix is a square matrix of size O n × O n ,where O n is the total number of orbits generated by all the synaptic signal graphlets of size n . Theentries tr i,j of this matrix correspond to the number of times an orbit changes from one state toanother.We then enumerate all graphlet-orbit transitions. The similarity metric is based on an arithmeticmean of orbit-transition differences, which is normalized to reduce biases associated with differentsized networks. Normalizing the rows of the orbit transition matrices gives20igure 12: The graphlet-orbit transition matrix corresponding to the possible changes of state of aedge-ordered digraph of size three from one snapshot to another. ntr i,j = tr i,j (cid:80) | O | k =1 tr i,k (6)For any two synaptic signal graphs G and G , we compute their similarity by the average oftheir graphlet-transition frequency for each graphlet-transition ntr i,j . This metric is given by theorbit-transition agreement (OTA): OT A ( G , G ) = 1 | O | × | O | (cid:88) i =1 | O | (cid:88) j =1 (1 − | ntr G i,j − ntr G i,j | ) (7) We introduced an approach for analyzing the similarity between two dynamic patterns computedby our competitive refractory dynamic model operating on an underlying structural network with afixed connectivity topology. The technical approach we took was to create directed graphlets withmultiple edges to encode multiple signals between connected nodes, and then apply edge-orderingin order to account for variable (synaptic) weight contributions to the running summation of arriv-ing signals at a given target node. These constructions extended graphlets to multidigraphlets andedge-ordered multidigraphlets. A challenge presented by the new graphlet types was the enumer-ation process. We overcame this by describing the space of graphlets with a vector space whosecoordinates represented a multidigraphlet class.We studied the graph transitions going from one node activation time to the next. The crucialobservation was the fact that there is a subgraph that is preserved between transitions. This createsa constraint on the graphlet-orbit transition matrices. Another important consideration was toconduct pairwise comparisons by only comparing graphlets observed in one or both synaptic signalgraphs, which saves computational resources. 21uture work will focus on the algorithms that compute these graphlet-based analyses. In ad-dition, we are examining how edge-graphlets perform in contrast to vertex-graphlets, and we areextending graphlets into a persistent framework, thereby extending the work that has been donewith persistent homology and topological data analysis (TDA).22 eferences [1] A
PAR ´ ICIO
D., R
IBIERO
P., S
ILVA
F., (2017)
Extending the Applicability of Graphlets toDirected Networks , IEEE/ACM Transactions on Computational Biology and Bioinformatics,Vol. 14, No. 6.[2] A
PAR ´ ICIO
D., R
IBIERO
P., S
ILVA
F.,
Network comparison using directed graphlets ,arXiv:1511.01964 [cs.SI].[3] A
PAR ´ ICIO
D., R
IBIERO
P., S
ILVA
F.,
Temporal network comparison using graphlet-orbittransitions , arXiv:1707.04572 [cs.SI].[4] B
UIBAS
M., S
ILVA
G.A. (2011),
A framework for simulating and estimating the state andfunctional topology of complex dynamic geometric networks , Neural Computation, 23 183-214.[5] C
HARTRAND
G., L
ESNIAK
L.,
Graphs & Digraphs , 4th Ed.; Chapman & Hall/CDR; UnitedStates, 2005.[6] D E S ILVA
J., M
OLLA
T., P
FENDER
F., R
ETTER
T., T
AIT
M.,
Increasing Paths in Edge-Ordered Graphs: the Hypercube and Random Graph. (2015)[7] G
RIGOR ’ YAN
A., M
URANOV
Y., V
ERSHININ
V., Y AU S.(2018);
Path homology theory ofmultidigraphs and quivers ;[8] H
OLME
P., S
ARAM ¨ AKI , J., (2012)
Temporal networks , Phys. Rep., 519, 97–125.[9] H
ULOVATYY
Y., C
HEN
H., M
ILENKOVI ´ C T.(2015)
Exploring the structure and function oftemporal networks with dynamic graphlets.
Bioinformatics 31(12):i171-i180.[10] M
ILENKOVI ´ C T., P R ˇ ZULJ
N. (2008),
Uncovering biological network function via graphletdegree signatures , Cancer Informatics 6,257-273.[11] M
ILO
R., S
HEN -O RR S., I
TZKOVITZ
S., K
ASHTAN
N., C
HKLOVSKII
D., A
LON
U. (2002)
Network motif: simple building blocks of complex networks . Science 298(5594):824-827.[12] N
EWAZ
K., M
ILENKOVI ´ C T., ”Graphlets in Network Science and Computational Biology”
Analyzing Network Data in Biology and Medicine , edited by Nataˇsa Prˇzulj. (2019), p.193-240.[13] P R ˇ ZULJ
N. (2007)
Biological network comparison using graphlet degree distribution.
Bioin-formatics 23:177–183.[14] P
UPPO
F., G
EORGE
V., S
ILVA
G.A. (2018)
An optimized structure-function design principleunderlies efficient signaling dynamics in neurons
Nature Scientific Reports 8:10460.[15] R
IBIERO
P., S
ILVA
F.,
G-Tries: A Data Structure for Storing and Finding Subgraphs , DataMining and Knowledge Discovery, vol. 28, no. 2, pp. 337-377, 2014.2316] R
UDIN
W.,
Principles of Mathematical Analysis . (1964)[17] S
ARAJLI ´ C A, M
ALOD -D OGNIN
N, N
EBIL Y AVERO ˘ GLU ¨O., P R ˇ ZULJ
N.,
Graphlet-basedcharacterization of directed networks , Scientific Reports 6, 35098 (2016).[18] S
ILVA
G.A.,(2019)
The effect of signaling latencies and node refractory states on the dy-namics of networks , Neural Computation, 31(12).24
Appendix
A.1 Undirected graphlets
Undirected graphs were first used as a suitable model for representing any binary relation. Dueto its simplicity they were used among different research fields showing a high efficiency for de-scribing the complex structure of what they were modeling. After having being fully used forcomparing and classifying networks depending on the frequency and distribution of a small subsetof subgraphs, it was discover that similar networks posses a particular type of signature: these areinduces subgraph which are called graphlets. Below we show all the possible graphs with theirrespective orbits up to 5 vertices. Each orbit represents a different topological property.Figure 13: Note: The induced undirected graphs up to five nodes. Reprinted from [10, Figure 1].
A.2 Directed graphlets
In real-world networks, a directed edge between two nodes represents an interaction of a kind.These interactions can vary from a macroscale to a microscale, reaching social, economical, phys-ical and biological fields: for instance, in biology transcription networks describes all of the reg-ulatory transcription interactions in a cell, whereas in economics one can construct a world tradenetwork to study trade flows.Using directed graphs for representing a system allowed us to gather more information abouthidden relationships within the system itself. This also led to increased complexity of possibletopological states, which were also captured by graphlets. Below we show all the possible graphletswith their respective orbits up to 5 vertices. 25igure 14: Note: The induced directed graphs up to four nodes. Reprinted from [17, Figure 1].
A.3 Orbit transition matrix, undirected case
Each orbits within a graphlet carries a particular topological information. These topological fea-tures are also induced by the automorphism set, that is, for any automorphism f of the graphlet,any node in a j -th orbit goes necessarily to another node which also belongs to the j -th orbit.This implies that any graphlet automorphism induces a permutation of the different orbits sets intothemselves.When comparing two temporal networks, information is lost with the use of techniques basedon the frequency and distribution of graphlets. Neither frequency nor distribution are capable ofconsidering the evolution of the systems. One way to overcome this is taking into account the pos-sible changes of orbit’s states. Below we show the graphlet-orbit transition matrix correspondingto the possible orbit transitions of 3-node in the undirected case.Figure 15: Each entry tr i,j in the matrix represents the number of times orbit i changes to orbit jj