Chapter 6

Efficient Simulation of Network

Models

Next, we will deepen the discussion of the implementation of graph models includ-
ing technically more advanced graph models. This discussion includes alternative
implementations of the Erdős Rényi graph and different implementations of the
configuration model. Moreover, we will discuss the time and space complexity of
these simulation algorithms.

6.1 Time and space complexity of algorithms and data-

structures
In order to analyze and discuss the efficiency of algorithms and data-structures, the
notations of complexity analysis provide an useful tool. Therefore, before turning
to particular algorithms on graphs, we review briefly some essential definitions and
notations for the time and space complexity of algorithms and problems are given.
We also look at the non-standard multi-parametric input size case, which will be
important in order to discuss algorithms on the important case of sparse graphs
and networks. For a more comprehensive introduction, see [8].
The running time of an algorithm is measured by the number of elementary op-
erations required to produce the desired output. Elementary operations can be ad-
dition, multiplication, comparison, logical operators (and/or) and, depending on
the computational model, we might also be allowed to perform nonlinear function
computations, such as (sin(.), cos(.), log()., etc.). Note, that in computer science, the
model computer is typically a Turing machine, named after the British computer
science pioneer Alan Turing. It is a kind of minimalistic automaton which can per-
form all computations that can be performed on a conventional (von Neumann)
computer. To keep the discussion simple, we will here assume that the elementary
operations have about the same running time and process at most a single item of
the input and that this running time is given by some constant.
An algorithm that solves a computational problem, computes a solution to a
problem using a sequence of elementary operations and making use of some mem-
ory. In order to describe the performance of algorithms independent from the par-
ticular hardware and programming language, it is common to categorize their per-

119
120 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

Figure 6.1: Running time for algorithms of different time complexity

formance by means of computational complexity classes. The computational (time

and space) complexity of algorithms is typically provided in asymptotical (Landau)
notation (also called ‘Big O’ notation). Given a function T : N → R, for instance the
running time of an algorithm, we write

• T ∈ O(h), if and only if ∃n 0 ∈ N, c ∈ R+ : ∀n ≥ n 0 : T (n) ≤ c · h(n)

• T ∈ Ω(h), if and only if ∃n 0 ∈ N, c ∈ R+ : ∀n ≥ n 0 : T (n) ≥ c · h(n)

• T ∈ Θ(h), if and only if T ∈ O(h) and T ∈ Ω(h)

Alternatively, we can state

• T ∈ O(h) ⇔ limn→∞ T (n)/h(n) < ∞,

• T ∈ Ω(h) ⇔ limn→∞ T (n)/h(n) > 0

• T ∈ Θ(h) ⇔ 0 < limn→∞ T (n)/h(n) < ∞

In Figure 6.1, we compare running time for different time complexities.

The following exercise serves to highlight and understand some often used re-
sults:

Exercise 6.1 Show that the following results hold:

1. For constant q ∈ N: Θ( f (n))) = Θ(q f (n))

2. log(n!) ∈ O(n log n);

3. For constant q ∈ N its holds that: Θ(log(n q )) = Θ(log(n));

4. For constant q ∈ N its holds that: log(n)q ∈ O(n 2 ); (Hint: use the rule of L’Hospital)
p
n
5. 2 ̸∈ Ω(2n ).

The asymptotic notation is also used to discuss the space used by the algorithm.
Consider an algorithm with running time T (n) which uses memory of size S(n) to
perform the computation for input n. We say that

• its time complexity is in O(h(n)), if and only if T ∈ O(h), and

• its space complexity is in O(h(n)), if and only if S ∈ O(h).

6.2. DATA STRUCTURES AND BASIC ALGORITHMS FOR GRAPHS 121

Importantly, we distinguish between the worst case time/space complexity (=maxium

time/space required over all problem instances of length n) and the average case
time/space complexity (=average over all problem instances given a probability mea-
sure on the space of instances of length n). If it is not indicated otherwise, in these
lecture notes we will refer to the worst case time/space complexity.
When implementing random graph models, the outcome of an algorithm can
be different for each run and it depends on some random numbers. Algorithms that
depend on random numbers (which follow a given distribution) are called proba-
bilistic algorithms. Typically, their output and running time also depends on the
random numbers and will follow some stochastic distribution.
In the discussion of data structures used to store data (e.g. a specification of a
network), the space complexity the data structure requires. Closely related is the
question of the running time of basic update and query operations, given the data is
stored in a certain data structure. We will discuss this for graphs or ‘network data’ in
the following section.
In case of graphs, the number of edges, say m = |E |, and the number of vertices,
say n = |V |, needs both to be considered in order to analyse the performance of data
structures and algorithms. For instance, in the case of thin graphs with m ∈ O(n)
other data structures and algorithms are advantageous than in the case of dense
graphs with m ∈ Θ(n 2 ).
In such multi-parametric cases, we can use the following definition :
h : N × N → R and define T ∈ O(h) if and only if ∃n 0 ∈ N, m 0 ∈ N, c ∈ R+ : ∀n ≥ n 0 , m ≥
m 0 : T (n, m) ≤ c · h(n, m).

6.2 Data structures and basic algorithms for graphs

There are three common data structures used to store graphs:

• Edge list: The graph is stored as a set of vertices V , and a set of edges (pairs of
vertices) E ⊆ V × V .

• Adjacency list: The graph is stored as an adjacency list, that is an array of ver-
tices, and for each vertex a list of nodes that are adjacent to it.

• Adjacency matrix: The graph is stored in a double indexed array a[i ][ j ], i ∈

{1, ..., n}, j ∈ {1, ..., n} and n = |V |, where a[i ][ j ] = 1 if node i is connected to
node j and a[i ][ j ] = 0 otherwise.

The space complexity of these data structures is as follows:

• Edge list: Θ(|V | + |E |). (Assuming a list of nodes is provided besides the list of
edges)

• Adjacency list: Θ(|V | + |E |). (Every node requires one storage cell and every
edge is represented once by a node, in the list of the predecessor.)

• Ajacency matrix: Θ(|V |2 ). (Every pair is represented by means of a boolean

variable, even if there is no edge.)

In practise, the above data structures can be augmented by auxiliary data struc-
tures, such as search trees, used to quickly locate elements in a list. In the context of
databases and data structures such auxiliary data structures are called index. They
122 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

1 1 2 null
from↓ to→ 1 2 3 4
2 3 null 1 1 1 0 0
2 0 0 1 0
3 4 null 3 0 0 0 1
4 0 0 0 0

4 null

Figure 6.2: Adjacency list representation (left), and adjacency matrix representation
(right) of the graph ({1, 2, 3, 4}, {(1, 1), (1, 2), (2, 3), (3, 4)})

serve to make computations more efficient, but they require additional storage. Be-
sides balanced search trees, hash indices are also often used for this purpose. Alter-
native to using in index, we might also keep lists sorted which makes certain oper-
ations more efficient, e.g., checking the existence of an edge, while others might be
less efficient, e.g., insertion of an edge. The implementation of an balanced tree in-
dex, e.g. B+ trees, AVL trees, or red-black trees, is however not trivial and for details
we refer the interested reader to the literature, e.g., [2].
The time complexity of some basic operations are reviewed next:

• Querying the existing of an edge: For an adjacency matrix this can be accom-
plished in O(1). For an edge list the time complexity is in O(|E |). If an tree
index is used or sorting is maintained it reduces to O(log(|E |)). (Binary search
can be used to find the edge). For adjacency lists, the time complexity is in
O(|V |). (The first node needs to be found in the outer list, and then the sec-
ond node in the list connected to this node.). Using a tree index or sorting, this
can be further reduced to O(log n).

• Inserting and deleting edges: For the adjacency matrix the time complexity
of insertion and deletion of an edge is both in O(1). For edge list, the inser-
tion is O(1) and the deletion is O(n). (insertion at the end of the list, deletion
requires to move all subsequent edges in the list, at most |E | − 1, in memory).
If a tree index is used the time complexity of insertion will be in O(log(|E |))
and that of deletion in O(log(|E |)). (Rebalancing the tree is required). Keep-
ing sorted lists will yield a time complexity in O(|E |) for both insertion and
deletion. (Subequent elements will have to be moved in memory). The time
complexity of insertion and deletion of edges in the adjacency list can be ob-
tained by replacing |E | by |V | and are also governed by the efforts required for
list processing, this time a list of length |V |.

According to the above, it matters for the time and space complexity of graph
algorithms, in which format the data is represented. In the complexity analysis of
algorithms of graphs, this needs to be adressed. The adjacency matrix is a good rep-
resentation if the number of nodes stays fixed and there are frequent queries and
updates required. For sparse graphs, the better space efficiency of edge list and ad-
jacency list is however significant. Moreover, these data structures have advantages
for models where nodes are dynamically inserted or deleted.
6.3. ERDŐS-RÉNYI RANDOM GRAPH IMPLEMENTATION 123

In the context of random graph models, the adjacency matrix representation is

often uses, as it allows quick update of edges based. To implement the preferential
attachment model, one might consider also adjacency list representations, in par-
ticular if nodes are inserted one by one. Moreover, the adjacency list representation
allows for a more efficient implementation of shortest path algorithms, as they are
required in when computing higher order graph properties.

*Alternative graph representations There are also many alternative representa-

tions of graphs that are either used in domain specific context or that are mainly
of theoretical interest. We will briefly summarize geometrical representations on
graphs and succint circuit representations.
One example are geometrical representations of graphs such as the unit disc
graph, where an edge exists if the discs of radius r around the vertices have a non-
empty intersection. This graph just requires the storage of the coordinates of the
spheres and a radius. This provides a space complexity of O(|V |). However, not
every graph can be represented as a unit disc graph. There are other geometrical
representations of graphs with important practical applications. In so called visibil-
ity graphs an edge exists, if an object, e.g., a line segment) can be seen from another
object, e.g., another line segment. Visibility graphs are frequently used in computer
vision and robot motion planning. Only in 3-D it is possible to represent all graphs
by means of visibility graphs, see also [8].
To summarize, many of the alternative graph representation show different trade-
off in terms of time or in space complexity. Often they cannot represent all possible
graphs. In particular application domains, however, they can be very useful.

6.3 Erdős-Rényi random graph implementation

The goal is to write a program that samples instances of random graph models. The
generation of the instance can be viewed as a stochastic simulation of a process that
generates the graph. The simulation of the graph should ideally be such, that the
likelihood of generating instances follows a specified probability distribution. Be-
fore implementing the more advanced configuration model, we will start with im-
plementing the random graph model for testing oursoftware and for the purpose of
comparison.

6.3.1 Simulating Erdős-Rényi Graph models

In the literature on random graph models, there are two variants of the Erdős-Rényi
Graph model. With G p (n, p), we will denote the model that was discussed in previ-
ous chapters and denoted with E R, and with G m (n, m) a model that pre-defines the
exact number of edges in the graph. More precisely:

1. G p (n, p): An instance of the Erdős-Rényi Graph G p (n, p) is a graph with n

nodes and it is decided with equal probability p for each pair i and j whether
there is an edge from i to j .

2. G m (n, m): An instance of the Erdős-Rényi Graph G m (n, m) is a graph with m

edges. The edges are assigned randomly and without repetition among the
pairs i to j with i ̸= j .
124 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

The simulation of the G p (n, p) model was already discussed in Chapter 3. One way
to generate a random graph of type G(n, m) is to use Algorithm 1.

Algorithm 1 Generating edges for G m (n, m) random graph model.

1: function G ENERATE R ANDOM E DGES ((m: no. of edges, n: no. of vertices))
2: k ←0 ▷ Generate all potential n(n − 1) Edges
3: for i = 1 : n do
4: for j = 1 : n do
5: if i ̸= j then
6: k ← k +1
7: e[k] ← (i , j )
8: for i = 1 : n(n − 1) − 1 do ▷ Shuffle the edge list by random permutation.
9: z ← uniform random integer between i and n(n − 1)
10: Swap content of e[i ] and e[z]
11: return e[1], . . . , e[m] ▷ Return first m edges

The algorithm above uses the Fisher Yates shuffle to generate a random permu-
tation, in computer science also known as Algorithm 235 [3].

Exercise 6.2 Discuss the time and space complexity of function G ENERATE R ANDOM E DGES
in terms of n = |V | and m = |E |. Can the time and space complexity be improved? If
so, how?

Exercise 6.3 (Simulation of Erdös-Renyi graph with a fixed number of edges) In this
exercise you can use the example implementation of G(n, m) from previous exercise
(or implementing a new graph model if you prefer):

1. Implement the G m (n, m) as a random graph model (you can also use the code
from the previous exercise).

2. Compare the size of the biggest connected component for G p (n, m/(n(n − 1)))
and G m (n, m) for different settings of m: m = n/2, m = n, m = n(n − 2) + 1 and
n = 5 and n = 20. Run the graph model repeatedly (say 10 times) and, explain
your observation.

Exercise 6.4 Discuss the time and space complexity of Stochastic Block Model (SBM),
that you implemented in the last programming homework.

6.4 Watts-Strogatz Model

The Watts-Strogatz Model [?] is a mathematical model used to describe small-world
networks, and exhibit both high clustering and short path lengths between nodes.
Let’s assume that you have to create a network of n nodes where each node might
connect with k nearest neighbors on both sides and the rewiring probability is p.
The model works as follows:

1. Lattice Structure: Start with a ring lattice of n nodes where each node is con-
nected to its k neighbors within a certain range.
6.5. CONFIGURATION MODEL IMPLEMENTATION 125

2. Rewiring: Randomly rewire some of the edges with a probability p, intro-

ducing randomness to the network. This rewiring decreases the average path
length while maintaining a high clustering coefficient.

For p = 0: The network is a regular lattice, with high clustering but large average
path lengths.
For p = 1: The network becomes a random graph, with low clustering but short av-
erage path lengths.
For intermediate values of p: The model produces a small-world network with both
high clustering and short path lengths, which mirrors many real-world networks.
An example of Watts-Strogatz model is shown in Figure ??.

Figure 6.3: Watts-Strogatz Graph

Exercise 6.5 This exercise focuses on a better understanding of the Watts-Strogatz

model:
1. Write down the algorithm to implement Watts-Strogatz model.
2. Implement this algorithm in Python.
3. Discuss the time and space complexity of your algorithm.

Exercise 6.6 Check if your algorithm from Exercise 6.5 lead to disconnected network.
If yes, update your algorithm so that it always returns a connected graph.

6.5 Configuration model implementation

With this section we return to the so-called Configuration Model [15, 1, 2] (CM for
short) introduced in Chapter 4. We first recall the abstract idea behind the model
and make general considerations. We then focus on the specific implementations
that have been proposed to realize this idea.
Let us consider undirected networks first. The goal of the CM model is to assign
each vertex i a desired degree k i and then generate an ensemble of graphs com-
patible with the resulting degree sequence ⃗k, by drawing links at random between
vertices in such a way that the desired degree sequence is realized. All graphs com-
patible with the desired degree sequence should be realized with the same proba-
bility. The degree sequence ⃗ k can be drawn from any desired degree distribution
126 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

P (k), e.g. a scale-free one with desired exponent. However, as discussed in Chap-
ter 4 (Subsection 4.1.3), the degree sequence must be graphical, i.e. realizable by at
least one graph.
For directed networks, the CM is easily generalized by assigning each vertex i
a given in-degree k ii n and a given out-degree k iout , or in other words by specifying
the in-degree sequence ⃗ k i n and the out-degree sequence ⃗ k out . Now the goal is to
generate an ensemble of random directed graphs in such a way that the desired in-
and out-degree sequences are simultaneously realized.
Note that, unlike the preferential attachment model (see Chapter 4), the CM
does not make explicit hypotheses on how networks organize themselves in a given
structure. It rather is a null model, generating a random(ized) ensemble of graphs
once some low-level information is assumed as an input. As we will see in detail in
Chapter 10, the CM can also be used as a benchmark for empirical data: a compar-
ison between the CM and a real-world network allows us to check whether some of
the higher-order properties observed in the real-world network are consistent with
those generated by the CM using the same degree sequence ⃗ k ∗ as the real network
∗
G . If this is the case, then one can conclude that the observed higher-order prop-
erties are a mere outcome of the specified form of the degree distribution (which is
a first-order property, according to our discussion in Chapter 1), being consistent
with a random assignment of links compatible with the degree sequence. If this
is not the case, then the observed deviations from the model indicate interesting
structural patterns that cannot be traced back to the null hypothesis (i.e. they are
not explained by the degree sequence alone).

6.5.1 Link stub reconnection

Now, we will discuss the implementation of the CM model, i.e., also called ‘link stub’
connection process. Here, as many link stubs (half links) as the prescribed degree
k i are initially attached to each vertex i . All these stubs are then randomly matched,
with the aim of realizing a random graph with the desired degree sequence ⃗ k, or, in
the directed case ⃗ k out ,⃗
kin.
In principle, iterating this process generates as many realizations of the network
as desired, and samples the graph ensemble defined by the given degree sequence.
The connection of edges should be done in such a way that every possible network
with the prescribed degree distribution occurs with the same probability. A high
level description of this algorithm is given in Algorithm 2. To accomplish this the
algorithm first attaches stubs according to ⃗ k and then connects them randomly. A
snapshot of this algorithm in the third iteration is given in Figure 6.4.

Homework 6.1 Write a graph model that generates a graph of the configuration
model by starting from the high level algorithm description above (Algorithm 2). Write
this code in the Link_Stub_Reconnection function. (20%)

Homework 6.2 How can the algorithm be made more time efficient by using a Fisher
Yates shuffle as in Algorithm 1? Discuss the time/space complexity of Algorithm 2 and
the time/space complexity of the algorithm using the Fisher Yates shuffle. Use the ‘Big
O’ notation and express the time complexity in terms of the number of nodes and
edges. (15%)
6.5. CONFIGURATION MODEL IMPLEMENTATION 127

Figure 6.4: After three iterations of the simulation algorithm. (see http://
[Link]/generating_networks_desired_degree_distribution)

Algorithm 2 ‘Link stub’ reconnection algorithm for generating a configuration

model with prescribed degree distribution.
1: function G ENERATE C ONFIGURATION M ODEL (k[1], ..., k[n])
2: Input k[1], ..., k[n]: graphical degree sequence
3: Initialize the list of edges to null
4: Attach link stubs to vertices in the graph according to the prescribed degree
5: sequence.
6: while Not all link stubs are connected do
7: Find a ‘link stub’ that is not yet connected
8: Choose randomly another stub that is not yet conncted and connect
9: the two link stubs
10: Add the resulting edge to the list of edges
return list of edges
128 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

Homework 6.3 Run the program and produce repeatedly (say thirty times) graphs
with the following degree sequences: 12 vertices, each:

⃗
k (1) = (2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2)

and
⃗
k (2) = (7, 4, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1)
Report on the number of self loops and repeated edges found in each of the experi-
ments. (10%) Is there a significant difference in the number of self loops for graphs of
the degree sequence ⃗
k (1) as compared to ⃗k (2) ? If so, what could be a (plausible) expla-
nation? (5%)

6.6 Problems with the implementation of the Configu-

ration Model
We are now going to discuss some limitations of the (repeated) ‘link stub’ reconnec-
tion algorithm and suggest ways of how to overcome them.
A big problem with the above implementations of the CM is that it either gives
rise to undesired self-loops and multiple links between two vertices or in case of the
repeated configuration model, it might require a large number of repetitions.
If these extra links are not discarded, the resulting graph ensemble cannot be
consistently compared with real networks that do not admit such kinds of links,
since this comparison might highlight patterns that are merely due to the differ-
ences between the two topological classes. At a mathematical level, one can identify
degree sequences for which these problems can be avoided (e.g. for sparse graphs
in the asymptotic limit n → ∞; see Chapter 4). However, in practice, if one wants to
define the CM using the empirical degree sequence ⃗ k ∗ of a real-world network G ∗ ,
then the above problems cannot be easily avoided, as empirical degree sequences
typically violate the above assumptions.
If double links and self-loops are deliberately removed, then the original degree
sequence can no longer be realized. Computationally, if one rejects the attempted
matching of two link stubs that would result in double links or self-loops, the al-
gorithm will typically ‘get stuck’ in configurations where large-degree vertices have
no more eligible partners to connect to. So, operatively, the link stub reconnection
method is typically unsatisfactory when an empirical degree sequence ⃗ k ∗ is taken
as input.

6.6.1 Repeated Configuration Model

One way to address this issue is the Repeated Configuration Model, which runs the
configuration model multiple times until a network with no self loops and multi-
edges is created. There is a mathematical prove that the probability of creating a
graph without self loops and multi-edges, given it exists, is positive. Therefore the
expected time for the repeated configuration model to produce a feasible network is
finite; it might however be high. In the homework we will test how many repetitions
it takes for different graphs to generate a feasible instance.

Homework 6.4 Repeated Configuration model:

1. Create a new copy of your program. Change your program so that when a graph
6.6. PROBLEMS WITH THE IMPLEMENTATION OF THE CONFIGURATION MODEL129

Figure 6.5: Elementary step of the local rewiring algorithm for a) undirected and b) directed
networks. Two edges, here (A, B ) and (D,C ), are randomly chosen from graph G 1 and the
vertices at their ends are exchanged to obtain the edges (A,C ) and (D, B ) in graph G 2 . Note
that the degree of each vertex is unchanged (in the directed case, the in- and out-degrees are
separately conserved).

with self-loops and multi-edges is discovered, it is discarded, and the program starts
over to try with a new graph. Implement this method in repeated_configuration_model
function.

2. Run your program repeatedly (30 times) to generate graphs with degree distri-
bution ⃗k (1) = (2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2) and ⃗
k (2) = (7, 4, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1). Based
on your experiments, how often does your program have to restart on average? What
could be the difference in the average number for the two degree distributions? (try to
state a brief answer, no proof required). (15%)

6.6.2 The local rewiring algorithm

In an attempt to overcome the above limitation, some important variants of the CM
have been introduced. Maslov, Sneppen and Zaliznyak [15] proposed to start from
a whole network G ∗ (not just its degree sequence ⃗ k ∗ ) and iteratively randomize the
∗
topology of G in a degree-preserving manner, thereby generating an ensemble of
random graphs with the same degree sequence ⃗ k ∗ as the original network. Their
randomization process consists in what they call the local rewiring algorithm, which
deliberately avoids the occurrence of multiple links and self-loops in the random
networks.
The local rewiring algorithm consists of the iteration of the elementary step shown
in Fig. 6.5a and Fig. 6.5b for undirected and directed networks, respectively: two
links are randomly chosen from the initial graph G 1 and the vertices at their ends
are exchanged in such a way that the new graph G 2 has the same degree sequence
of the initial one. If the ‘new’ links already exist in the network, then the step is
aborted and two different links are randomly chosen again. In this way an ensemble
of random networks is generated, having the same degree sequence as the original
network and no multiple links or self-loops.
Note that, since any two vertices cannot be connected more than once, in the
local rewiring algorithm the presence of links between high-degree vertices is sup-
pressed, determining a certain degree of ‘spurious’ disassortativity which is not due
to a ‘basic’ anticorrelation between vertex degrees. This important point highlighted
by Maslov, Sneppen and Zaliznyak led them to show that much of the disassortativ-
ity observed in the Internet (see Section 1.2.2) can be accounted for in this way.
130 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

Recently, it has been proved mathematically that the local rewiring algorithm
is biased, i.e., it does not explore the space of graphs compatible with the degree
constraints uniformly [6, 7]. Roughly speaking, the root of the problem is the fact
that the algorithm explores with higher probability the graph configurations that are
‘closer’ to the orginal network. In order to overcome this problem, one should intro-
duce a suitable ‘acceptance probability’ for each attempted configuration. However,
the calculation of this probability is computationally demanding because it depends
on the current configuration, and should therefore be repeated at each step of the
algorithm.

Exercise 6.7 (Local Clustering Coefficient) Implement a function that computes the
average local clustering coefficient of a graph - please do not use the built-in cluster-
ing coefficient computation of NetworkX (or any other library) - (see Section 1). You
can use the program below, which computes the local clustering coefficients of nodes,
as a basis. Compute the average local clustering coefficient for graphs G 1 and G 2 in
Figure 6.6.
Note: The local clustering coefficient of a node is the fraction of triangles that ac-
tually exist over all possible triangles in its neighborhood. The average local clustering
coefficient of a graph G is the mean of local clustering coefficient of all nodes in G.

import numpy as np

# Function to compute the local clustering coefficient of each

node in an undirected graph
def local_clustering_coefficient(adj_matrix):
n = adj_matrix.shape[0]
clustering_coeff = [Link](n) # Initialize an empty numeric
array to store clustering coefficients

for i in range(n):
neighbors = [Link](adj_matrix[i, :] == 1)[0] # Find the neighbors of node i
k_i = len(neighbors) # Degree of node i

if k_i < 2:
clustering_coeff[i] = 0 # If degree is less than 2,
clustering coefficient is 0
continue

# Create the sub-adjacency matrix for the neighbors

sub_matrix = adj_matrix[neighbors[:, None], neighbors]

# Calculate the actual number of edges (E_i) among the neighbors

E_i = [Link](sub_matrix) / 2 # Divide by 2 because
each edge is counted twice

# Calculate the local clustering coefficient for node i

clustering_coeff[i] = (2 * E_i) / (k_i * (k_i - 1))

return clustering_coeff
6.7. SIMULATING THE BARABÁSI-ALBERT MODEL OF PREFERENTIAL ATTACHMENT131

# Test the function using a sample adjacency matrix

adj_matrix = [Link]([[0, 1, 1, 0],
[1, 0, 1, 1],
[1, 1, 0, 1],
[0, 1, 1, 0]])

# Compute the local clustering coefficients

clustering_coeff = local_clustering_coefficient(adj_matrix)
print("Local clustering coefficients: ", clustering_coeff)
Homework 6.5 (Local rewiring algorithm) 1. Implement the local rewiring algo-
rithm in a procedure that can start from a given graph with correct (or prescribed)
degree distribution.
2. Next, run the local rewiring algorithm on the two initial graphs G 1 and G 2 depicted
in Figure 6.6. Show in a plot how the average local clustering coefficient changes over
time. Describe the differences between the plot for G 1 and G 2 briefly in your own
words. (15%)

Figure 6.6: Graphs for initializing the rewiring algorithm. G 1 and G 2 that both have
the degree sequence (2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2).

6.7 Simulating the Barabási-Albert model of preferen-

tial attachment
Network models can also be inspired by the process that determines network grows.
In the models we simulated so far, the nodes were given and only connections were
added by the simulation. In a network growth model, nodes might appear one-by-
one and attach to some random nodes. A typical example for such a network growth
model is the preferential attachment model and we will look at the simulation and
analysis of this model in the following.
The Barabási-Albert (BA) model is a preferential attachment model. It serves as
a contrast to the Erdős-Rényi model, focusing on the generation of scale-free net-
works in which the degree distribution is not concentrated around a mean value,
132 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

but there appear nodes of all orders of magnitudes including large nodes which
form hubs in terms of network connectivity. While the Erdős-Rényi model produces
graphs G = (V, E ) where the likelihood of an edge between any two nodes is con-
stant, the BA model incorporates preferential attachment to create networks that
better mimic real-world structures.
Let us discuss next the simulation of the BA model. In the BA model, a graph
starts with m 0 initial nodes. New nodes are then sequentially added to the graph,
each establishing m ≤ m 0 edges to existing nodes. The probability Pr(k) of connect-
ing to a node with degree k is given by:

kq
Pr(k) = P q
j (k j )

where the sum is taken over all degrees k j of existing nodes, and q is a constant,
that is in the canonical BA model set to q = 1. This preferential attachment results
in a network whose degree distribution follows a power law, thus characterizing it
as scale-free. High-degree nodes, often called “hubs," are a hallmark of such net-
works [9].
For implementation, an adjacency matrix A = a i j can represent the graph G =
(V, E ). The Python source code and a plot illustrating the BA model are provided at
the end of this section (see Figure 6.7). We can see some interesting characteristics
of the network from the simulated data:

• Clearly, the network is connected, which should be obviously the case due to
the construction

• There is a wide range of degrees in the networks with some nodes with high
degree (so-called hubs) and others with a small degree.

• By further investigation of the data, we will find that the distribution is scale-
free

Note, that the fact that the BA is scale-free and forms a small world network is
quite sensitive to the fact that preferential attachment works with a linear prefer-
ential attachment probability (q = 1). In the homework of this chapter you will be
invited to change the preferential attachment probability to a super- or sub-linear
proportionality q < 1, q > 1 and observe how this affects the scale free distribution
properties.

import numpy as np
import random
import networkx as nx
import [Link] as plt

# Barabasi-Albert model function as before

def barabasi_albert(N, m0, m):
adj_matrix = [Link]((N, N))

for i in range(m0 - 1):

for j in range(i + 1, m0):
adj_matrix[i, j] = 1
adj_matrix[j, i] = 1
6.7. SIMULATING THE BARABÁSI-ALBERT MODEL OF PREFERENTIAL ATTACHMENT133

Figure 6.7: Barabási-Albert Network Plot

for new_node in range(m0, N):

degree = [Link](adj_matrix[:new_node, :new_node], axis=1)
probabilities = degree / [Link](degree)
selected_nodes = [Link](range(new_node), size=m, p=probabilities)

for node in selected_nodes:

adj_matrix[new_node, node] = 1
adj_matrix[node, new_node] = 1

return adj_matrix

# Parameters and generate network

N = 100
m0 = 5
m = 3
adj_matrix = barabasi_albert(N, m0, m)

# Create graph object

134 CHAPTER 6. EFFICIENT SIMULATION OF NETWORK MODELS

Figure 6.8: Degree vs. Frequency of each degree in a BA network having n=100, and
m=3.

g = nx.from_numpy_array(adj_matrix)

# Calculate degrees for each node for sizing

deg = [Link]([d for n, d in [Link]()])

# Plotting the graph

[Link](figsize=(10, 7))
pos = nx.spring_layout(g)
[Link](g, pos, with_labels=False, node_size=deg * 20, node_color=’skyblue’, edge_colo
[Link]("Barabasi-Albert Network")
[Link]()

Below find a visualization of the degree distribution on log-log plot in 6.8. It can
be observed that the distribution follows a linear downward trend for the logarith-
mically scaled data.

#code to plot degree distribution

import [Link] as plt
import networkx as nx
import collections

def plot_degree_dist(G):
degrees = [[Link](n) for n in [Link]()]
6.7. SIMULATING THE BARABÁSI-ALBERT MODEL OF PREFERENTIAL ATTACHMENT135

freq=[Link](degrees) #compute unique degrees and their frequencies

unique_deg=[Link]()
deg_freq=[Link]()
print(unique_deg)
print(deg_freq)
[Link](unique_deg, deg_freq, ’ro’) # code to make a scatterplot
[Link](’log’)
[Link](’log’)
[Link](’Degree’)
[Link](’Frequency’)
[Link]()

plot_degree_dist(nx.barabasi_albert_graph(100, 3))

Exercise 6.8 Think of an alternative method to implement Barabási-Albert model.

Hint: Create a list of nodes to pick uniformly at random.

Homework 6.6 1. Modify the Barabási-Albert model in the source code given above
to sample from a sublinear probability distribution and a superlinear probability dis-
tribution. Implement both in BA_sublinear_PA and BA_superlinear_PA functions, re-
spectively.
2. Add the log-log plots for the degree distributions for a sample of these two network
models.
The sublinear preferential attachment model is given by:
p
k
Pr(k) = P q
j kj

and the superlinear preferential attachment model is given by:

k2
Pr(k) = P 2
j (k j )
.
Describe briefly your findings. Does the plot resemble a plot of a network with a
scale free distribution? (20%)
Bibliography

[1] Russel, S.J. and P. Norvig (1995): Artificial intelligence: a modern approach.
Series in Artificial Intelligence, Prentice Hall.

[2] Skiena, S. S. The Algorithm Design Manual. New York: Springer-Verlag, pp. 177
and 179, 1997.

[3] Durstenfeld, R. (July 1964). "Algorithm 235: Random permutation". Commu-

nications of the ACM. 7 (7): 420. doi:10.1145/364520.364540

[4] Galperin, H., and Wigderson, A. (1983). Succinct representations of graphs. In-
formation and Control, 56(3), 183-198.

[5] Sergei Maslov, Kim Sneppen, Alexei Zaliznyak, Detection of topological pat-
terns in complex networks: correlation profile of the internet, Physica A: Sta-
tistical Mechanics and its Applications, Volume 333, 15 February 2004, Pages
529-540

[6] A.C.C Coolen, A. De Martino, A. Annibale, J. Stat. Phys. B 136, 103567 (2009).

[7] E.S. Roberts, A.C.C. Coolen, Phys. Rev. E 85, 046103 (2012).

[8] Helmut Alt, Michael Godau, Sue Whitesides, Universal3-dimensional visibility

representations for graphs,Computational Geometry, Volume 9, Issue 1, 1998,
Pages111-125

[9] Barabási, Albert-László and Albert, Réka, Emergence of scaling in random net-
works, Science, 286(5439), 509–512, 1999.

[10] Piva, G. G., Ribeiro, F. L., & Mata, A. S. (2021). Networks with growth and prefer-
ential attachment: modelling and applications. Journal of Complex Networks,
9(1), cnab008.

[wattsstrogatzmodel]

136