{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Introduction\n",
"\n",
"In this workbook, we demonstrate the usage of Graphs.jl to execute diffusion simulations on random graphs. The broad aim is to replicate the results from the influential paper, \"A simple model of global cascades on random networks\".\n",
"\n",
"Watts, D. J. (2002). A simple model of global cascades on random networks. Proceedings of the National Academy of Sciences, 99(9), 5766-5771."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Outline\n",
"\n",
"Through this model, Watts demonstrated that a simple threshold rule could be used as a basis to explain global cascades, i.e., events where small initial random shocks spread to a significant proportion of the network. \n",
"\n",
"One rather surprising result from this study is that even when the average degree of the network is low, a global cascade can be triggered by randomly chosen seeds. This throws a spaner in the usual conception that large scale diffusions on a network are triggered only by \"special people\", i.e., influencers. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Model\n",
"\n",
"A node in a graph is faced with a binary decision - to engage or to not engage (e.g., with new products or discussions with friends). They make this decision based on a simple rule - they compute the fraction of their neighbors that have engaged, compare it with their personal threshold, and engage if the fraction of engaged neighbors exceeds the threshold. In other words, the thought process is on the lines, \"If at least 18% of my friends bought the new iPhone, I would want to buy it too\".\n",
"\n",
"The execution of this model is random and asynchronous, i.e., each node checks the status of its neighbors in a random order and decides whether to engage or not depending on the fraction of neighbors engaged at that point"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m registry at `C:\\Users\\Thibaut\\.julia\\registries\\General`\n",
"\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m git-repo `https://github.com/JuliaRegistries/General.git`\n",
"\u001b[32m\u001b[1m Resolving\u001b[22m\u001b[39m package versions...\n",
"\u001b[32m\u001b[1m No Changes\u001b[22m\u001b[39m to `C:\\Users\\Thibaut\\.julia\\environments\\v1.8\\Project.toml`\n",
"\u001b[32m\u001b[1m No Changes\u001b[22m\u001b[39m to `C:\\Users\\Thibaut\\.julia\\environments\\v1.8\\Manifest.toml`\n",
"\u001b[32m\u001b[1m Resolving\u001b[22m\u001b[39m package versions...\n",
"\u001b[32m\u001b[1m No Changes\u001b[22m\u001b[39m to `C:\\Users\\Thibaut\\.julia\\environments\\v1.8\\Project.toml`\n",
"\u001b[32m\u001b[1m No Changes\u001b[22m\u001b[39m to `C:\\Users\\Thibaut\\.julia\\environments\\v1.8\\Manifest.toml`\n",
"\u001b[32m\u001b[1m Resolving\u001b[22m\u001b[39m package versions...\n",
"\u001b[32m\u001b[1m No Changes\u001b[22m\u001b[39m to `C:\\Users\\Thibaut\\.julia\\environments\\v1.8\\Project.toml`\n",
"\u001b[32m\u001b[1m No Changes\u001b[22m\u001b[39m to `C:\\Users\\Thibaut\\.julia\\environments\\v1.8\\Manifest.toml`\n"
]
}
],
"source": [
"using Pkg\n",
"for pkg in [\"Graphs\", \"Plots\", \"GR\"]\n",
" Pkg.add(pkg)\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"using Distributions\n",
"using Graphs\n",
"using Plots; gr()\n",
"using Random: shuffle, seed!\n",
"\n",
"# make sure this tutorial is reproducible\n",
"seed!(20130810);"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"fraction_engaged (generic function with 1 method)"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"function fraction_engaged(node::Int,\n",
" G::Graphs.SimpleGraphs.SimpleGraph,\n",
" node_status::BitVector)\n",
" \"\"\"\n",
" Computes the fraction of neighbors engaged within the neighborhood\n",
" of a given node. It uses the node status to check the engagement status of\n",
" each of the nodes neighbors\n",
" \"\"\"\n",
" num_engaged_neighbors = 0\n",
" for nbr in neighbors(G, node)\n",
" if node_status[nbr] == true\n",
" num_engaged_neighbors += 1\n",
" end\n",
" end\n",
" return num_engaged_neighbors/length(neighbors(G, node))\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"update_node_status! (generic function with 1 method)"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"function update_node_status!(G::Graphs.SimpleGraphs.SimpleGraph,\n",
" node_status::BitVector,\n",
" threshold::Vector{Float64})\n",
" \"\"\"\n",
" This function executes the random asynchronous updates of the entire network\n",
" at each time step. In this conceptualization, each time step comprises mini\n",
" time steps where a randomly shuffled node list updates at each iteration\n",
" \"\"\"\n",
" \n",
" for node in shuffle(vertices(G))\n",
" if node_status[node] == false\n",
" if fraction_engaged(node, G, node_status) > threshold[node]\n",
" node_status[node] = true\n",
" end\n",
" end\n",
" end\n",
"\n",
" return nothing\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"diffusion_simulation (generic function with 1 method)"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"function diffusion_simulation(n::Int,\n",
" z::Int,\n",
" threshold::Vector{Float64},\n",
" T::Int,\n",
" n_realizations::Int)\n",
" \"\"\"\n",
" Executes the diffusion simulation.\n",
" It creates a new Watts-Strogatz graph at each realization and seeds a single node\n",
" (initialization). It then updates the network for the defined tnumber of time\n",
" steps.\n",
"\n",
" The idea is to run the diffusion simulation a very large number\n",
" of times and count the instances where we observe a\n",
" global cascade, i.e., number of nodes engaged after the simulation\n",
" process is a sizeable proportion of the network\n",
"\n",
" Hyper Parameters of the model\n",
" ----------\n",
" 1. Number of nodes in the Watts-Strogatz graph (n)\n",
" 2. Average degree (z)\n",
" 3. Threshold (distribution or a specific value)\n",
" 4. Time steps for simulation to be run\n",
" 5. Number of realizations\n",
"\n",
" Output\n",
" -----------\n",
" A vector of number of engaged nodes at the end of each realization\n",
" of the simulation\n",
"\n",
" Intended usage of results\n",
" -----------\n",
" Plot the cascades on the z-ϕ phase space; replicate results from\n",
" \"A simple model of global cascades on random networks\", Watts (2002)\n",
" \"\"\"\n",
"\n",
" output = Vector{Int}(undef, n_realizations)\n",
" beta = z/n\n",
"\n",
" for r in 1:n_realizations\n",
" G = watts_strogatz(n, z, beta)\n",
" # Select a single random node from the network and seed it\n",
" node_status = falses(nv(G))\n",
" node_status[sample(vertices(G))] = true\n",
"\n",
" # Update the network for predefined number of time steps\n",
" for _ in 1:T\n",
" update_node_status!(G, node_status, threshold)\n",
" end\n",
" output[r] = sum(node_status)\n",
" end\n",
"\n",
" return output\n",
"\n",
"end"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Results\n",
"\n",
"We have two parameters that we wish to explore - the average degree of the network and the threshold for engagement. While we leave the full exploration of the parameter space to the reader, we demonstrate how this exploration could be done below. We take the example of two particular values of (z, $\\phi$), i.e., (5, 0.18) and (6, 0.18) at which a phase transition occurs, i.e., at average network degree of 5 global cascades occur while at degree 6 no cascades occur irrespective of threshold!"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" 2.289443 seconds (2.01 M allocations: 526.882 MiB, 3.20% gc time, 0.33% compilation time)\n"
]
}
],
"source": [
"const N = 10^4\n",
"z, threshold = 5, 0.18\n",
"\n",
"@time data = diffusion_simulation(N, z, fill(threshold, N), 50, 100);"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"histogram(data, xlab=\"Number of engaged nodes\", ylab=\"Frequency\", legend=false)"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" 3.093600 seconds (2.01 M allocations: 526.875 MiB, 2.25% gc time)\n"
]
}
],
"source": [
"z, threshold = 6, 0.18\n",
"\n",
"@time data = diffusion_simulation(N, z, fill(threshold, N), 50, 100);"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"