Randomly Evolving Idiotypic Networks: Structural Properties and Architecture

0.87 MB·English

by  Holger Schmidtchen

Randomly Evolving Idiotypic Networks:
Structural Properties and Architecture
Holger Schmidtchen,1,2 Mario Thu¨ne,1 and Ulrich Behn1,2,∗
1Institut fu¨r Theoretische Physik, Universita¨t Leipzig, POB 100 920, D-04009 Leipzig, Germany
2International Max Planck Research School Mathematics in the Sciences,
Inselstraße 22, D-04103 Leipzig, Germany
(Dated: January 19, 2012)
WeconsideraminimalisticdynamicmodeloftheidiotypicnetworkofB-lymphocytes. Anetwork
noderepresentsapopulationofB-lymphocytesofthesamespecificity(idiotype),whichisencodedby
abitstring. Thelinksofthenetworkconnectnodeswithcomplementaryandnearlycomplementary
2 bitstrings, allowing for a few mismatches. A node is occupied if a lymphocyte clone of the corre-
1 spondingidiotype exists, otherwise it is empty. There is a continuous influxof new B-lymphocytes
0 of random idiotype from the bone marrow. B-lymphocytes are stimulated by cross-linking their
2 receptors with complementary structures. If there are too many complementary structures, steric
hindrance prevents cross-linking. Stimulated cells proliferate and secrete antibodies of the same
n
idiotypeas their receptors, unstimulated lymphocytesdie.
a
Depending on few parameters, the autonomous system evolves randomly towards patterns of
J
highly organized architecture, where the nodes can be classified into groups according to their
7
statisticalproperties. Weobserveanddescribeanalyticallythebuildingprinciplesofthesepatterns,
1
whichallowtocalculatenumberandsizeofthenodegroupsandthenumberoflinksbetweenthem.

Thearchitectureofallpatternsobservedsofarinsimulations canbeexplainedthisway. Atoolfor
]
B real-time pattern identification is proposed.
C
PACSnumbers: 87.18.-h,87.18.Vf,87.23.Kg,87.85.Xd,64.60.aq,05.10.-a,02.70.Rr
.
o
i
b I. INTRODUCTION toireandtheautonomousdynamicsofanimmunesystem
-
q not exposed to foreign antigen. It provides a mechanism
[ B-lymphocytes play a crucial role in the adaptive of immunological memory. Imagine that an antigen Ag

1 immune system. They express Y-shaped receptor isrecognizedbyanantibodyAb1. Thus,thecloneofAb1
v molecules, antibodies, on their surface. Antibodies have expandsandpossiblymeetsanothercloneofcomplemen-
8 very specific binding sites (idiotopes) which determine tary idiotype Ab2. Both mutually stimulate each other
1 their idiotype. All receptors of a given B-cell have the and they persist even after Ag has disappeared. Ab2 is
6 same idiotype. B-cells are stimulated to proliferate if structurallysimilartoAg andcanbeconsideredasinter-
3 their receptors are cross-linked by structures which are nalimageofAg. Furthermore,the networkisthoughtto
1. complementary to the idiotopes, e.g. by foreign antigen. control autoreactive clones. All these issues are beyond
0 StimulatedB-cells thus survivewhereasunstimulated B- the concept of clonal selection.
2 cells die, this process is called clonal selection [1]. The network paradigm got an immediate enthusias-
1
Afterafewgenerations,stimulatedB-cellsdifferentiate tic response and idiotypic interactions were considered
:
v to plasma cells which secrete soluble antibody molecules as the major regulating mechanism of the immune sys-
Xi of the same idiotype, which may bind to complementary tem. However,the rapidprogressofmolecular immunol-
sites on antigen and mark them for further processing. ogy, difficulties in the direct experimental verification,
r
a B-lymphocytes of different idiotype are continuously andthe discoveryofother regulatingmechanisms letthe
produced in the bone marrow in a remarkable diversity. interest of experimental immunologists decay. Yet, for
The variety of the potential idiotype repertoire, created systembiologiststhenetworkparadigmalwaysremained
by somatic reshuffling of gene segments and mutations attractive. Several aspects of the original concept were
[2], was combinatorially estimated [3] to exceed most revisedinduecourse. Mostnotably,VarelaandCoutinho
likely 1010. [5–7] suggested second generation networks with an ar-
Complementary structures may be found on antigen chitecture comprising a strongly connected central part
but could also be situated on other antibodies of com- withautonomousdynamicsandasparselyconnectedpe-
plementary idiotype. B-lymphocytes can stimulate each ripheralpartforlocalizedmemoryandadaptiveimmune
other and, thus, they form a functional network, the id- response. Inasense,theyreconcilebothparadigmsofid-
iotypic network [4]. iotypicnetworksandclonalselection. Areadablehistory
Jerne’s concept of the idiotypic network explains in a ofimmunologicalparadigmscanbefoundin[8]. Reviews
naturalwaythediversityoftheexpressedidiotypereper- withfocusonidiotypicnetworksare[9]withanemphasis
onmodelingapproachesandmorerecently [10]withem-
phasis on new immunological and clinical developments.
∗ [email protected] Todaythemainactivitiesareinclinicalresearch. Idio-
2
typicinteractionsarethebaseofalltherapieswithmono- the reader familiar with the systems behavior. Consid-
clonal antibodies [11]. New experimental techniques ering global and local network characteristics the exis-
make large-scale studies of the idiotypic repertoire fea- tence of groups of nodes which share statistical proper-
sible[12]whicharenecessarytoinferthenetworksarchi- ties [13] is confirmed and more details are revealed. In
tecture. Sec. IV we describe certain regularities in the bitstrings
In this paper we consider a minimalistic model of the of nodes which belong to the same group. The observed
idiotypic network, which was first formulated and inves- patterns can be constructed from pattern modules [17].
tigated in [13]. In this model a node represents lympho- The construction principle is explained first for the sim-
cytesandantibodiesofagivenidiotype. Lymphocytesof plest pattern and generalized afterwards. These findings
complementary idiotype can stimulate each other. The allow to calculate the number of groups,the groupsizes,
corresponding nodes are connected by links. and the linking between groups. A new observable, the
Idiotypesarerepresentedbybitstrings[14]oflengthd, center of mass, is introduced that proves very useful in
b b b with b 0,1 . Ideally, d is chosen such real-time pattern identification. In Sec. V we apply this
d d−1 1 i
that 2d·i·s·the size of t∈he{pote}ntial repertoire. The nodes conceptconsideringspecific patterns observedin simula-
ofthenetworkarelabeledbythesebitstrings. Twonodes tions, among them a dynamic pattern with core groups,
v and w are linked if their bitstrings are complementary peripheral groups, stable holes, and singletons that re-
allowing for up to m mismatches, i.e. their Hamming sembles in some aspects the biological network [5, 6]. In
distance is d (v,w) > d m. The corresponding undi- the appendix we calculate the scaling of the relative size
H
−
rected graph is the base graph G(m). Each node has the of these groups for systems of biological size.
d
same number of neighbors, κ = m d . It represents
k=0 k
the potential idiotypic repertoire with all possible inter-
P (cid:0) (cid:1) II. CONTEXT AND RELATED MODELS
actions.
Not all idiotypes are expressed in the real network.
In our minimalistic model we only account whether an The topic falls into several scientific disciplines. It
idiotype is present or not, correspondingly the node is is natural to place our investigations in the context of
occupied, n(v) = 1, or empty, n(v) = 0. The subgraph network theory. Network theory has applications in a
of G(m) induced by the occupied nodes represents the plethora of different, multidisciplinary fields [18–20] and
d has received great attention in the community of statis-
expressed idiotypic network at a given time.
tical physicists in the last decade.
We describe the temporal evolution in discrete time.
Amajorbodyofresearchdealswithgrowingnetworks,
The influx of new idiotypes from the bone marrow is
where new nodes are attached to the existing nodes ran-
modeled by occupying empty nodes with probability p
domly, or depending on properties of the existing nodes.
[15]. For survival a B-lymphocyte needs stimulation by
In this context deletion of nodes is only considered to
complementary structures. The dose-response curve is
study the resilience of the network against random or
log-bell shaped, cf. [16] and Refs. therein. If there
targeted attacks [21]. For recent reviews see [22–24], cf.
aretoomanycomplementarystructures,crosslinkingbe-
also [19, 25].
comes less likely due to steric hindrance, the stimulus is
Natural networks however do not grow without limit
reduced. In our model an occupied node survives if the
but stay finite and evolve towards a functional archi-
numberofoccupiedneighborsisbetweentwothresholds,
tecture. There are several modes to enable evolution:
t and t . The rules of parallel update are
L U
(i) Adding nodes but keeping the growth balanced by
(i) Occupy empty nodes with probability p deletion[26–31]ormergingofnodes[32–34]. (ii)Keeping
(ii) Count the number of occupied neighbors n(∂v) of the nodes unaffected but add, delete, or reorganize the
node v. If n(∂v) is outside the window [t ,t ], set links [35–40]. (iii) Addition and deletion of both nodes
L U
the node v empty and links [41–47]. Again, this can be done randomly or
depending on the properties of the nodes and its neigh-
(iii) Iterate.
bors.
Driven by the random influx of new idiotypes the net- Generic observables to characterize networks include
work evolves towards a quasi-stationary state of non- thedegreedistribution,centrality,betweenness,cliquish-
trivial, functional architecture in which groups of nodes ness, modularity, clustering coefficients, and diameter.
can be identified according to their statistical proper- Forinstance,growingnetworksusingpreferentialattach-
ties. Crucial for that is besides the random occupation menthaveapowerlawdegreedistributionlikemanyreal
ofemptynodes,thatoccupiednodesareemptiediflinked world networks. However,the characteristic exponent of
with too few or too many occupied nodes. preferential attachment networks is larger than the one
The paper is organized as follows. In the next section found in naturalnetworks. This was a major motivation
we discuss the model in its scientific context in other to study evolving networks.
disciplines and its relation to other models of idiotypic Inmanynaturalnetworksnodeshaveindividualprop-
networks. In Sec. III we provide simulation results. We erties which control their potential linking. Clearly, in
sketcha typicalrandomevolutionof the systemto make ourcasethisistheidiotype. Alsoproteinnetworks,tran-
3
scription networks and generally signaling networks be- dow update rule to simulate the internal dynamics and
long to this class. a zero/one clone population. However, while we con-
Nodes may have an internal state which can change sideradiscreted-dimensionalhypercubicshapespace,in
depending on their neighborhood in the network, or on theirmodelthecomplementaryidiotypesliveondifferent
external influences. This dynamics has a typical time sheets of a 2D continuous shape space [65]. A summary
scale that is shorter than the time scale for evolution of ofresultsobtainedinmodelswithcontinuousshapespace
the network’s architecture. The interplay of these pro- is given by Bersini [66]. Several aspects of modeling in
cesses came into the focus of research only in the last continuousanddiscreteshapespacearediscussedin[67].
few years. For a recent review and a status report see An early network model inspired by spin glass physics
[48, 49]. Our model is a very early example where this was proposed by Parisi [68] to describe immunological
interplayisstudied[13]. Inthepresentpaperwedescribe memory. Theinteractionbetweenidiotypesinthemodel
the building principles of the architecture. ofBarraandcollaborators[69–71]isalsotakenfromspin
Our network model is a Boolean network [50], since glass physics. Their model describes a given number of
eachnode can be only in one of two states, empty or oc- idiotypepopulationseachwithaconstantnumberoflym-
cupied. Thenodesareupdatedinconsecutivetime steps phocytes. Each cell can be in a firing or quiescent state.
depending on its own and its neighborhood occupation. The strength of the ferromagnetic coupling between the
The model is also a cellular automaton, for a compre- idiotypes(alsoencodedbybitstrings)modelstheaffinity,
hensive monograph see [51]. More precisely, since the whichisrelatedtothecomplementarityofthebitstrings.
update depends only on the sum of the neighbor states Barra and Agliari [69, 70] compute the degree distribu-
and the state of the node itself, it is a totalistic cellular tion, the type and number of loops, and consider the
automaton. Cellularautomatanaturallyinvolveunoccu- scaling behavior. They describe primary and secondary
pied nodes. In our model unoccupied nodes, holes, play immuneresponsesandunderstandlowandhighdosetol-
animportantrole. Wedistinguishholesthatcouldbeoc- erance as a network phenomenon. In our model the log-
cupied from stable holes which can not be occupied due bell-shaped response is integrated as a B-cell property
to overstimulation. A simple network picture disregards [6, 16]. Essential for our approach is a random influx of
nodes which are not occupied. new cells and a selection mechanism. Hence, we can de-
scribe the evolution of the network towards a functional
A famous example of a totalistic two-dimensional de-
architecture.
terministic cellular automatonis Conway’sGame of Life
A model which distinguishes between antibody
[52]. Both,the survivalofanoccupiedcellandthe occu-
molecules and lymphocytes with a dynamics including
pation of an empty cell, are governed by window rules.
an influx similar to ours is proposed by Ribeiro et al.
Many interesting patterns, depending on the initial con-
[72]. They compare the dynamics on given architectures
ditions static or dynamic, have been described in detail
of random and scale-free networks.
and classified. The Game of Life was transfered to a va-
riety of lattices, e.g. to three dimensional cubic lattices IMMSIM, invented by Celada and Seiden [73, 74], is
[53,54],totriangular,pentagonalandhexagonaltessela- alsoamodifiedcellularautomaton. Itincorporatesmany
tions,toPenrosetilings[54–56]andtosmallworldgeom- immunologic agents, including antigen presenting cells,
etry [57]. Larger than Life [58–60] increases the radius B-cells, T-cells, antigens, antibodies and cytokines. It
of the neighborhood. Also, a probabilistic version on a describesbothhumoralandcellularresponses. Idiotypes
2d square lattice has been proposed, where stochastic are also characterized by bitstrings. The vast amount
deviations from the deterministic update are permitted of interacting agents increases drastically the number of
[61]. The mean occupation of cells, a global order pa- parameters. The model is intended to be as realistic as
rameter, undergoes a sharp phase transition for increas- possible and to provide experimental and practical im-
ingstrengthofstochasticity. Ontheoccasionofthe40th munologistswith a tool to test hypotheses in silico. Our
anniversary appeared a comprehensive collection of re- approachisin asensecomplementary. We aimatanun-
centresultsontheGameofLifeanditsdescendants[62]. derstanding of the principles governing the autonomous
Our model can be considered as a further version of the evolutiontowardafunctionalarchitecture. Thereforewe
Game of Life on a high dimensional graph, where empty investigate a minimalistic model with a small number of
nodesarerandomlyoccupied,while the survivalofoccu- parameters,whichneverthelessexhibitsessentialfeatures
pied nodes is governed by a deterministic window rule. of the biological network.
Starting from an empty graph we observe an evolution Theconceptofidiotypicnetworksinspiredapplications
toward a complex, highly organized architecture. in computer sciences, e.g. artificial immune systems for
the detection of intrusion, of spam or viruses [75, 76].
Themodelcanalsobecategorizedasastochastic,non-
linear dynamical system.
There exists a variety of models for B-cell networks.
III. SIMULATION RESULTS
References [9, 16] give comprehensive surveys of mod-
elling approaches.
For instance, Stewart and Varela [63, 64] proposed In this section we report on simulations on the base
a model, which also has a random influx and a win- graph G(2), which consists of 212 = 4096 nodes, each of
12
4
which has κ=79 links to other nodes. The window rule system size such catastrophes become less likely. In the
parameters are [t ,t ] = [1,10]. The lower threshold is thermodynamic limit we expect a breaking ofergodicity.
L U
biologically motivated, a node needs at least one occu- In our finite system there can be transitions between
pied neighbor to survive. The upper threshold is chosen differentpatternsonaroutewithoutextinctioncatastro-
to enable a non-trivial, dynamic pattern of complex ar- phe but via the formation of an intermediate, unstable
chitecture. As discussed below in more detail, higher t giant cluster, for more details see [13].
U
would allow a broader variety of static patterns. How- Modeling biological systems we should keep in mind
ever, the concept developed in this paper still applies. that they are finite and have a finite life expectation.
The simulations start with an empty base graph. The The quasi-stationary state could persist for times longer
influx p varies from 0 to 0.11. This covers the range in than the life span of the individual but transitions be-
which we find interesting patterns, above p ' t /κ = tween different states can not be excluded. Moreover,
U
0.127 there is only trivial random behavior. the parameters could vary during the individual’s life.
A. Random Evolution B. Global Characteristics
We start the simulations with an empty base graph Afirstcharacterizationofthedifferentpatternscanbe
occupying nodes with probability p. In the first time obtainedconsidering globalquantities. They include the
step only those nodes survive which have at least one number of occupied nodes on the base graph
occupiedneighbor(havingmorethant occupiedneigh-
U
bors is unlikely in the beginning). The surviving nodes n(G)= n(v), (1)
represent seeds the neighbors of which will survive if oc- v∈G
X
cupied. Hence, we observe a rapid growth of a giant
where n(v) 0,1 is the occupation of node v, the size
clusteruntilmoreandmorenodeshavemorethant oc-
U ∈{ }
of the largest cluster in the set of present clusters
cupiedneighbors. Exceedingthe upper thresholddeletes
C
a node. Thus, the giant cluster decays and many stable
Cmax =max(C ), (2)
holesarecreated,i.e. nodeswiththenumberofoccupied | | C∈C | |
neighbors above the upper threshold t . Fig. 1 shows a
U
time series for the evolution to a stationary state where and the average size of the clusters
the number of stable holes increases up to its stationary
1
value. C = C , (3)
C
The empty base graph is a highly symmetric object. h| |i NC | |
C∈C
X
Due to the random influx the symmetry is broken and
the system falls into a network configuration of lower where N is the current number of clusters, and
C S
h · i
symmetry. The typicalresult of the evolutionis a quasi- denotes the average over the elements of some set .
S
stationarypatternofmutuallydependentoccupiednodes Clusters, i.e. connected parts of the occupied subgraph,
andstableholes. Properlysituatedoccupiednodescreate are very characteristic for patterns. Finally, we mention
stable holes and in return stable holes in the neighbor- the number of stable holes h∗(G), i.e. empty nodes with
hood of an occupied node may prevent its overstimula- n(∂v)>t .
U
tion. Figure 1 shows a generic time series of these global
Inthe genericcasefora givenp a certainpatterntype characteristics for a parameter setting where the system
isfoundmostfrequently. Occasionally,depending onthe evolves to a stationary 8-cluster pattern.
history of the driving process, also other patterns occur. In the stationarystate we can considertemporally av-
Once established they all can live for a long time. (This eragedglobalquantities. Thetemporalaverageisdefined
can be proved in simulations preparing the pattern as as
initial configuration.)
1
In principle, the system is ergodic [77]. This becomes x= x , (4)
t
T T
immediately clear considering the following unlikely but 1 0
− t∈(XT0,T1]
possible event. For any choice of p > 0 the occupation
of nodes by the influx can be such that the application where T should be larger than the relaxation time
0
of the window rule leads to an empty base graph. After in which the system reaches the stationary state, and
this extinction catastrophe the further realization of the (T ,T ] is the averagingperiod.
0 1
driving process, the influx, determines to which pattern Table I gives results for three patterns which occur
the system evolves. For any p > 0 an infinite trajectory for different influx p. The static patterns are named ac-
contains partial trajectories leading to any possible pat- cording to the size of the characteristic clusters. With a
tern. In practice, in our simulations we have never seen deeper understanding of the architecture we shall name
suchanextinctioncatastrophebutonly quasi-stationary them by the number of groups of nodes, see Sec. IV be-
patterns usually living for a long time. For increasing low.
5
3000 2425
2424 2429
2500
Occ. nodes
Largest cluster 3782
des 2000 Avg.S ctalubslete hr osliezes 2553 3783 3778 3449
o
of n 1500 3654 2758
ber 2361
m 3814
Nu 1000 2409 3798 2393
3790
3270
500
3780
2417 4038
2937
0
0 200 400 600 800 1000 1200 1400
Time 2169 2427
FIG. 1. Time series of the number of occupied nodes n(G),
the size of the currently largest cluster |Cmax|, the average FIG. 2. A selection of 2-clusters as they appear for moder-
cluster size h|C|iC, and the number of stable holes h∗(G) on ate influx. We see ten 2-clusters of occupied nodes (black).
the base graph G(122) with [tL,tU]= [1,10] and p=0.01. The 2T-hcleuystaerre(aellmcpotnynceicrtceldest)hrboeucgohmaeshoucbcu(gpriaeyd),.tIhfethheubelewveonutldh
system evolvesto a stationary 8-cluster pattern.
not survive the next update. There is one link (solid) with
onemismatch,theotherlinks(dashed)havetwomismatches.
The nodes are labeled with the decimal expressions of their
TABLE I. Temporal averages and standard deviations of
bitstrings. Figure produced using yEd [78].
global characteristics for three typical patterns. Data from
500,000 iterations.
pied graph. We see one large cluster of about 400 nodes
Pattern 8-cluster 2-cluster dynamic and about 120 isolated nodes. These nodes had at least
p 0.005 0.015 0.025 one occupied neighbor which was removed in latest up-
n(G) 1024.6± 0.8 1023.3± 0.9 533.4± 27.7 date. Further,approximately1800stable holesand1800
h∗(G) 2816.0± 0.0 3055± 20 1828± 114 unstable holes have been observed, which are not shown
in the figure. We can clearly distinguish a central and a
|Cmax| 8.0± 2.8 4.5± 2.1 410± 20
peripheral part, which are supposed to be of functional
h|C|i 7.96± 0.04 2.01± 0.01 4.49± 0.96
C importance in the biological idiotypic network.
NC 128.6± 0.8 510.4± 3.7 123.1± 22.0
C. Local Characteristics
The 8-cluster pattern at p=0.005 has 128 clusters of
size8,whichresultsinatotalaveragepopulationof1024 Besides global quantities, time averagesof local quan-
occupiednodes,i.e. onefourthofallnodes. Theremain- tities characterizingeverysingle node canbe considered.
ingnodesareholes,2816ofwhicharestable(n(∂v)>tU), Elucidatingarethemeanoccupationn(v),thenumberof
and 256 are not stable (n(∂v)6tU). Empty nodes with occupied neighbors n(∂v), and the mean life time τ(v),
n(∂v) < tL are not counted as stable holes, since they which is defined as
could easily become occupied and sustained by the ran-
1
dom influx. Figure 1 shows the temporal evolution from τ(v)= n (v), (5)
t
an empty base graph towards an 8-cluster pattern. b(v)+nT0 t∈(XT0,T1]
Inthe 2-clusterpatternatp=0.015the 510clustersof
size2togetheroccupyaboutonefourthofthebasegraph. whereb(v)isthenumberofbirthsduringtheobservation
Theremainingthreequartersarestableholes. Defectsin time, i.e. the number of new occupations of the node by
theperfect2-clusterpatternallowthatoccasionallysome the influx. Of course, b(v)+n = 0 must be fulfilled,
T0 6
ofthenodepairsbecomeconnectedviaacentralhuband otherwise τ(v) has no meaning.
alargerstar-shapedclusterisformed,seeFig.2. t =10 We can identify groups of nodes sharing statistical
U
was chosen to exclude the occupation of the hub in the propertiesasproposedin[13]. Figure4showsmeanoccu-
perfect pattern. Both the 8- and the 2-cluster patterns pation,meanlifetime,andthenumberofoccupiedneigh-
are quasi-static, the temporal fluctuations are small. bors vs. influx probability. The groups appear as peaks
A more complex, dynamic pattern evolves for larger in the histograms. The number of occupied neighbors
p ' 0.03. The standard deviation of the temporal av- proves most suitable to distinguish the different groups.
erages is one order of magnitude larger than in the two For small and moderate influx a clear group struc-
static patterns. Figure 3 shows a snapshot of the occu- ture is visible. The data for the mean life time and the
6
1.1
1 4000
0.9 static dynamic random
3500
n 0.8
o
ati 0.7 3000
p
cu 0.6 2500
c
O
n 0.5 2000
a
e 0.4
M 1500
0.3
1000
0.2
0.1 500
0 0
0 0.02 0.04 0.06 0.08 0.1
Influx p
100000
4000
static dynamic random
10000 3500
e 3000
m 1000
e ti 2500
FIG. 3. (Color online) Snapshot of the occupied graph of n lif 2000
a 100
thecomplexconfigurationofadynamicpatternforp=0.025. Me 1500
The nodes are colored according to their membership in dif-
ferent groups, see Fig. 15 (top). Nodes within a group have 10 1000
similar statistical characteristics, see text. Figure produced 500
using yEd [78].
1 0
0 0.02 0.04 0.06 0.08 0.1
Influx p
80
mean occupation suggest that the patterns are static for
p/0.03. These patterns have groups of occupied nodes 70 static dynamic random 4000
with a high mean life time, other groups are stable holes s 60 3500
or sparsely occupied nodes. For 0.03 / p / 0.08 the bor 3000
h 50
patterns are dynamic, but still stationary. Also in dy- g
ei 2500
ntcihamamneigcoefspoacctoctneutrpinniesudotuhnseolryde.esWarieshislsemtaiabnllle,sttahhtoeilceosac.ncudTpdiheyednamsmuebiacgnrpalapifthe- ccupied N 3400 12500000
O 20
terns all nodes remain in their groups, for high influx 1000
p ' 0.08 the patterns become short-lived, groups dis- 10 500
solve and reemerge in a different configuration. Their
0 0
distinctionbytemporalaveragesbecomesmoreandmore 0 0.02 0.04 0.06 0.08 0.1
difficult. For very high influx p 0.12 the dynamics is Influx p
≫
entirely random.
FIG. 4. Time averages of the mean occupation (top), the
mean life time (center), and the mean number of occupied
neighbors (bottom) of each node for increasing values of p,
IV. CONCEPT OF PATTERN MODULES ∆p = 5/4096. Each graph shows isolines of a histogram in
top view. The grey scale reflects the number of nodes with
the corresponding property. Regimes of different temporal
We find regularities in the bitstrings encoding nodes
behavior are separated by vertical lines. For each p data is
belonging to the same group. This allows to identify
from time series of 500,000 iterations in the stationary state,
generalbuildingprinciples. Patternsarebuildfrommore starting from an empty base graph G(2).
elementary objects, so called pattern modules. With this 12
concept we can derive all structural properties, the size
of the groups and their linking, of the patterns observed
A. Simplest Case: The 2-Cluster Pattern
so far in the simulations.
We first explain the principles on the example of the
2-cluster pattern, and develop then the detailed concept We introduce the building principles of patterns con-
forthegeneralcase. Attheendofthesection,exploiting sideringthesimplestpattern,the2-clusterpatternwhich
this concept we introduce the center of mass vector, a appears for moderate influx. In simulations starting
tool which allows to identify in simulations the majority from empty base graphs it is rare and only found for
of patterns in real time. 0 < p / 0.03, but if prepared as initial condition it is
7
very stable for a larger range of p up to 0.045. By their 2−Cluster
statistical characteristics, cf. Table II, we can distin-
guish three groups of nodes, frequently occupied nodes 1 1
(S1) with a highmeanlife time, permanently empty sta- 22
ble holes (S ), and rarely occupied potential hubs (S )
3 2 11
which link together up to tU 2-clusters if occupied. 57
2 56 56
TABLE II. Local characteristics of the three groups in the
2-cluster pattern for p = 0.025 compared with the ideal
Hubs
pattern. Data from a timeseries of 500,000 iterations.
22
11
S1 S2 S3 3 1
mean occupation hn(v)iSi 0.993 0.0004 0.000 Stable holes
n 1 0 0
ideal
FIG. 5. The three groups of the 2-cluster pattern and their
occ. neighbors hn(∂v)iSi 1.002 10.95 55.64 linking. Thecirclesizescorrespondtothegroupsizes. Thear-
n(∂v) 1 11 56 rowsandthenumbersnexttothemindicatehowmanynodes
ideal
ofagroupexertinfluenceonthenodesofanothergroup,e.g.
eachnodeofS isstimulatedby11nodesfromS . Thenum-
mean life time hτ(v)iSi 6923 0.016 0.000 ber of links can2 be counted in simulations or tak1en from the
link matrix that is derived in Sec. IVB.
Looking at the node indices i in decimal representa-
v
tion we observed that the sum of the two indices in a S =2d−2 occupiednodesand S =2 2d−2 potential
1 2
2-cluster is constant in a realization. In a different real- | | | | ×
hubs.
ization the index sum can be different. For instance, in
Wefurthercancomputethenumberofoccupiedneigh-
Fig. 2 the index sum within all 2-clusters is 6207. This
bors n(∂v) of a node v of any group. Since all nodes of
indicates regularities at the level of the bitstrings.
S are occupied in the ideal pattern, n(∂v) is given by
1
We found that all occupied nodes are identical in ex-
the number of links between v and nodes in S . A link
1
actlytwobits,sayatpositionkandl. Themembersofa
between two nodes exists if their bitstrings are comple-
2-clusterare complementary inall other bits, in symbols
mentary except for up to two mismatches. If v S ,
1
we write ∈
it has two bits in common with all other nodes in S ,
1
···b ···b ··· is linked with ···b ···b···, namely bk and bl. Thus, all remaining bits must be ex-
k l k l
actly complementary. There is only one node w S ,
1
∈
where the bar denotes the bit inversion. The bitstrings which obeys this constraint. If v S or v S , there is
2 3
∈ ∈
of all stable holes are also equal in the same two bit po- one pre-determinedmismatch or none, respectively. The
sitions k and l. However, they are inverse to b and b remainingmismatchescanbedistributedamongthed 2
k l
−
of the occupied nodes. Potential hubs have exactly one non-determinant bits. Thus a node in S has
i
inverseandoneequalbitinthesepositions. Asthesebits
i−1
play a crucial role, we call them determinant bits. The d 2
n(∂v) = − (6)
regularities are summarized by ideal j
j=0(cid:18) (cid:19)
X
occupied nodes S ···b ···b ···
1 k l occupied neighborsin the ideal pattern. For small influx
···b ···b ···
potential hubs S k l this is in goodagreementwith the simulations,cf. Table
2 ( ···bk···bl··· II. In a similar way for all nodes the number of links to
stable holes S ···b ···b ··· . nodes in different groups can be calculated, the result is
3 k l
visualized in Fig. 5. The derivation for the general case
The example in Fig. 2 has the determinant bits in posi- is given in the next subsection.
tions 7 and 12, b =b =1. Thisregularityencouragedthefollowingconcept. Con-
7 12
This allows to explain all structural properties of the sidering the two determinant bits as coordinates of a
pattern observed in the simulations. We can construct two-dimensional space, they define the corners of a two-
an ideal 2-cluster pattern, a configuration in which all dimensional hypercube. The corner with coordinates
nodes of group S are occupied and the others remain (b ,b ) represents an occupied node, the opposite cor-
1 k l
empty. It is ideal in the sense that there are no defects ner (b ,b ) is a stable hole, and the neighboring corners
k l
but also no hubs. of (b ,b ) are potential hubs. We call this structure a
k l
Since all other bits can take all possible combinations, pattern module, because it is the building block for the
the size of the groups can be calculated, e.g. there are entire regular configuration. In a different picture, we
8
_
B10 patternmodule. OneofthemislabelledS . Allnodesin
B10 1
groupS havethesamedeterminantbitsb b . The
_ 1 1··· dM
other groups are ordered such that the determinant bits
B00
B00 of S differ in j 1 positions from the determinant bits
j
−
of S . It is clear that there are d +1 groups. Groups
1 M
S and S are equivalent, they have the same size
j dM+2−j
and linking properties, see below. In typical patterns,
occupation breaks the symmetry. We usually label the
smallest group of the more occupied half as S , which is
1
_ of course arbitrary.
B11
B11 The size of the groups can be obtained by elementary
_
combinatorics. Since, the i-th groupdeviates in i 1 out
B01 of d determinant bit positions of S , there are− dM
B01 M 1 i−1
choices. This number is multiplied by the number of
_ pattern modules on the base graph given by 2d−dM.(cid:0)Th(cid:1)e
module B module B
group size is
FIG. 6. Two pattern modules with complementary non-
determinant bit chains B and B, respectively, on a two mis- S =2d−dM dM , i=1,...,d +1. (7)
i M
match base graph with a 2-cluster configuration. The two- | | i 1
(cid:18) − (cid:19)
dimensional modules are congruently occupied, each consist-
ing of one occupied node (black, ·00), two potential hubs The factor dM is called the relative group size, since it
(gray, ·01 and ·10) and one stable hole (white, ·11). The is the groupis−iz1e normalized by the size of S . It is inde-
1
positions andvaluesofthedeterminantbitsarechosenwith- (cid:0) (cid:1)
pendent of the base graph dimension d and the number
out loss of generality. The links have no mismatch (bold),
of mismatches m.
one(solid),ortwomismatches(dashed). Thereare2d−2 pat-
As anexample,we canconstruct2-clusterpatternson
tern modules, each pair of which with complementary non-
determinant bits contributea pair of occupied nodes. a base graph G(m) by means of pattern modules with
d
exactly one occupied node (S ). The dimension of the
1
pattern module d then has to equal the number of al-
M
can also understand an ideal configuration as consisting lowed mismatches m. The number of qualitatively dis-
of 2d−2 congruently occupied ‘parallel worlds’. Figure 6 tinguishable groupsis m+1,etc. A 2-clusterpatterncan
illustrates the concept of pattern modules. emerge if the lower threshold is tL 6 1 and the upper
Any choice of the two determining bits is of course threshold obeys 1 6 tU 6 d m. The 2-cluster pattern
−
on 1-mismatch graphs described in [13] is an instance of
possible, all corresponding patterns are equivalent, the
2-cluster pattern is 22 d -fold degenerated, where the suchapattern. However,inthe1-mismatchcaseonehalf
× 2 of all nodes are occupied, the other half are stable holes.
first factor represents the choice of the two determinant
(cid:0) (cid:1)
bits, and the second factor gives the number of possible
positions of these bits in the bitstring of length d. It is
2. Linking
the individual history (the realization of the random in-
flux), which selects the determining bits andthus breaks
the symmetry. Each node on G(m) has κ = m d links, κ is con-
d k=0 k
strained by the allowed number of mismatches m. We
P (cid:0) (cid:1)
considera patternwithd determinantbits. Eachnode
M
in group S is linked to L neighbors in group S . The
B. General Case i ij j
L are the entries of the link matrix L. L defines the
ij
architecture of a pattern built of modules of dimension
1. Groups
d .
M
The dynamics of a node depends on the number of
Many results for 2-cluster patterns on the G(2) base occupied neighbors. The mean occupation is a typical
12
graph can be generalized to more complex architectures commonpropertyof nodes belonging to the same group.
andotherchoicesofdandm. Thisincludesthe8-cluster Thus, the knowledge of the group membership of the
pattern mentioned in Sec. III and other static patterns node’s neighbors is of crucial importance to understand
as well as the dynamic pattern. Their structure is cor- its statistical properties.
rectly described by pattern modules with more than two Within the concept of pattern modules L can be de-
determinant bits. rived combinatorially. Recall that the determinant bit-
In the same way as for the 2-cluster pattern, we de- string of a node in S deviates in l 1 bits from the de-
l
−
finethepatternmoduleasahypercubeofdimensiond , terminant bitstring of a node v in S . In the following
M 1 1
where d is the number of determinant bits. There are we consider a node v(i) chosen such that the i 1 bits
M
−
twogroupswhicharerepresentedbyonlyonenodeinthe inverse to the corresponding bits of v are left-aligned,
1
9
its non-determinant bits are denoted by B. This choice sincethenumberofbitscomplementarytov ischanged.
1
is without loss of generality, because all the arguments Wehavetoinvertasecondbitsothatthenumberofmis-
do not depend on the labeling of the bit position. From matcheswithv remainsconstant. The numberofnodes
1
the nodes in S we choose v(j) such that the j 1 bits in S with ∆ +2k mismatches to v(i) is computed as
j j ij
− | |
inverse to the corresponding bits of v are right-aligned follows.
1
and the non-determinant bits are B. There is no other Incase(ii)weinvertk bitsinthei 1left-alignedbits
−
node in Sj with less mismatches to v(i). ofv(j). Thereare i−k1 possibilities. Tocompensatethis,
We have to distinguish whether or not the partial bit- wehavetoinvertkbitsamongthej 1right-alignedbits
strings of length i 1 and j 1 do overlap. These cases of v(j). Altogether(cid:0), th(cid:1)ere are now ∆−+k mismatches to
− − ij
aTrheerdeisacrreimtihnraeteedcasbeys.the value of ∆ij = dM−i−j+2. vd(iM)−lii+v1ingpoosnsibthileitideMs t−oid+is1trriibguhtte-atlhigenmed. Tbihtus.sTwheehraevaere
(i)∆ =0. Alldeterminantbits ofv(i) andv(j) arecom- ∆ij+k
ij
plementary. (cid:0) (cid:1) i 1 d i+1
M
(ii)∆ >0. v(i) andv(j) share∆ determinantbitswith − −
ij ij k ∆ +k
v . This case is illustrated below. (cid:18) (cid:19)(cid:18) ij (cid:19)
1
nodes in S with ∆ +2k mismatches to v(i).
v b ... b b ... b b ... b j ij
1 1 i−1 i dM−j+1 dM−j+2 dM Case (iii) is similar. We invertk bits in the dM i+1
−
v(i) b1 ... bi−1 bi ... bdM−j+1 bdM−j+2 ... bdM rTioghcto-amlipgennesdatbeittshoisf,vw(je),htahveeretoarienvdeMrt−kki+b1itspoasmsiobnilgititehse.
v(j) b ... b b ... b b ... b d j +1 left-aligned bits of v(j(cid:0)). Toge(cid:1)ther, ∆ +k
1 i−1 i dM−j+1 dM−j+2 dM M − ij
mismatchestov(i)liveonthei 1left-alignedbits. There
i−1 bits ∆ij=dM−i−j+2 bits j−1 bits are i−1 possibilities. This−leads to
| {z } | {z } | {z } |∆ij|+k
(iii) ∆ij <0. v(i) and v(j) share ∆ij determinant bits (cid:0) (cid:1)
| | d i+1 i 1
which are inverse to the corresponding bits of v1, see M− −
diagram below. (cid:18) k (cid:19)(cid:18)|∆ij|+k(cid:19)
v1 b1 ... bdM−j+1 bdM−j+2 ... bi−1 bi ... bdM noTdehseirnesSujltwiinthca|∆seij(|i)+is2kobmtaisinmeadtcshetetsintogv∆(i).=0, and
ij
v(i) b1 ... bdM−j+1 bdM−j+2 ... bi−1 bi ... bdM the result for redistributing only |∆ij| mismatches de-
rived before is reproduced for k =0.
v(j) b1 ... bdM−j+1 bdM−j+2 ... bi−1 bi ... bdM Since the total number of mismatches between linked
nodes can not exceed m, l and k are constrained by l+
−∆ij=i+j−dM−2 bits dM−i+1 bits
2k+ ∆ 6m. Summarizing we obtain for ∆ >0
ij ij
i−1|bits {z }| {z } | |
d d i 1 d i+1
The|number of mism{aztches between v}(i) and v(j) is Lij = − M − M−
l k ∆ +k
∆ij . If ∆ij 6mthereisalinkbetweenv(i)andv(j). In lX,k=0(cid:18) (cid:19)(cid:18) (cid:19)(cid:18) ij (cid:19)
| | | |
this case there are further nodes in Sj which link to v(i), ×11(l+2k+∆ij 6m) , (8)
thenumberofwhichiscalculatedcombinatorially. There
are m ∆ additionally allowed mismatches which can and for ∆ij 60
ij
−| |
appearamongnon-determinantand/ordeterminantbits.
d d i 1 d i+1
Incases(ii)and(iii)therearefurthernodesinSj with Lij = − M − M−
l ∆ +k k
|∆ij|mismatchestov(i). Theseareobtainedbydistribut- l,Xk=0(cid:18) (cid:19)(cid:18)| ij| (cid:19)(cid:18) (cid:19)
ingthe ∆ mismatchesamongthed i+1right-aligned
| ij| M− 11(l+2k+ ∆ij 6m) . (9)
bits in case (ii), or among the i 1 left-aligned bits in × | |
−
case (iii). This leads to dM−i+1 and i−1 nodes in the In some cases the evaluation of the sums in Eqs. (8, 9)
∆ij |∆ij|
respective cases. leads to simple rules, see [79].
(cid:0) (cid:1) (cid:0) (cid:1)
Now we consider additional mismatches. Among the Of course, the total number of links of a node is
d d non-determinant bits we can distribute l mis-
− M dM+1
matches in L =κ. (10)
ij
j=1
d d X
− M The link matrix L=(Lij) has a symmetry given by
l
(cid:18) (cid:19)
L =L . (11)
different ways. ij dM+2−i,dM+2−j
Among determinant bits additional mismatches can If the link matrix entries of the ith row are multiplied
onlyappearinpairs. Thisisrelevantform−|∆ij|>2. If by the size of groupSi, the resulting matrix with entries
we invert one bit in v(j), the resulting node is not in S Λ = S L additionally obeys Λ =Λ .
j ij i ij ij ji
| |
10
m
1 2 3
|D | = 0
ij 4
no self−
|D |
coupling ij = 1 7 6 5
...
|D | = m
ij
core
1 2 3 4 5 6
no self−
coupling
12 11 10 9 8 7
FIG. 8. Groups and their linking for even (above), and odd
d +1
M dM (below) on a two-mismatch base graph. The lines show
possible links between nodes of the groups. The circles indi-
catetheself-coupling within thecore groups. Itis suggestive
FIG. 7. The general structure of the (dM+1) ×(dM+1)
to call arrangement of thegroups caterpillar representation.
link matrix. Matrix entries along the 2m+1 solid diagonal
lines are non-zero, cf. Eqs. (8, 9). All other entries are zero.
Groups in the center (gray square) couple to themselves and
are naturally called core groups. Among the groups without element in L. Such an element exists if dM is an even
self-couplingwecanfurtherdistinguishgroupsthatcoupleto number. Thus,
the core from those which do not. See further discussion in
thetext. m if either m or d odd
#core groups= M . (13)
( m+1 otherwise
Itmaybeofinteresttoknowthenumberoflinksfrom
Groups without self-coupling exist if the number of
a node of S to nodes of S with a given number of mis-
matches µ,iwhich is denotjed by Lµ. It is obtained re- groups dM+1 is larger than the number of core groups.
ij Such groups are related to the square submatrix of L
placinginEqs.(8,9)11(l+2k+|∆ij|6m)byδl+2k+|∆ij|,µ, with all entries equal to zero, cf. Fig. 7.
which gives, e.g. for ∆ >0
ij Thenumberofcoregroupsforagivenmonlydepends
on whether d is even or odd but is independent on d.
d d i 1 d i+1 M
Lµ = − M − M− . (12) We can represent the content of Fig. 7 in an alternative
ij µ 2k ∆ k ∆ +k
k (cid:18) − − ij(cid:19)(cid:18) (cid:19)(cid:18) ij (cid:19) way exploiting that groups with index i have the same
X
properties as those with index d + 2 i. In Fig. 8
M
Obviously, it holds Lij = mµ=0Lµij. these groups are on the same position i−n the two par-
allel strands. The core groups are on the right end of
P
the diagram. We distinguish two cases, with even and
3. Architecture odd number of core groups, respectively, which form the
“head” of a caterpillar. If we increase d by two, this
M
For the link matrix and, thus, for the architecture of does not change, but only the “tail” of the caterpillar
patterns we can state a number of general properties. gains an additional segment. This observation can be
Because ∆ = d i j+2 6 m is the necessary used to determine the scaling of the groupsizes for large
ij M
| | | − − |
andsufficientconditionfortheexistenceoflinksbetween d, cf. the appendix of this paper.
groupsiandj,allnon-zeromatrixentriesaresituatedon Note, that all these general properties in this section
a band along the secondarydiagonal,see the illustration arestructuralinformationaboutthepatterns. Theyonly
of the generalstructure of the link matrix in Fig. 7. The depend on the parameters d, m and on the choice of a
widthofthis bandis 2m+1,i.e. m controlsthe rangeof pattern module dM.
influence ofagroup. GroupS interactswith the 2m+1
i
groups S ,... ,S .
dM+2−i−m dM+2−i+m
We always find groups with self-coupling, because 4. Center of Mass
there is always a quadratic block of non-zero matrix en-
tries in the center of L. Such groups are called core In the simulations time series of 2d nodes are gener-
groups. The size of the core, i.e. the number of core ated, e.g. to calculate the mean occupation. In or-
groups, depends on m and on the existence of a central der to identify patterns in real time it has proven use-