Demo for SI model of disease on network with rewiring

Consider a network of N individuals, initially randomly connected with average node degree avgdegree. The paper [Gross et al, 2006] studies a disease spreading on the network governed by the following rules: Each node is either infected (I) or susceptible (S). At each time step nodes change their state and their connections according to the following rules:

Gross et al (2006) also formulated a system of three differential equations for the macroscopic quantities i (fraction of infected), ss (density of S-S links) and ii (density of I-I links). They found bistability between an endemic state (a large fraction of nodes is infected) and the disease-free state (the diesase dies out). Note that the disease-free state is absorbing such that, strictly speaking, the endemic state is always only a long-lived transient, even if for the macroscopic model the disease-free state is unstable.

Contents

Provided functions

Data structures

The network is stored in the form of two arrays.

The demo below performs a parameter sweep from $p=0.01$ down to $p=0$, showing an abrupt transition from the endemic state to the disease-free state.

Definition of basic quantities of for the network

The rates for node- and link-changing probabilities are kept in the structure rates. Our initial condition is a random netowrk where a randomly chosen 80% of the nodes are infected.

clear
addpath('../utilities/'); % some auxiliary functions used in the script
%
globaldefinitions; % S=0, I=1
N=1000;            % number of nodes
avgdegree=10;      % average degree
rates=struct('SfromI',0.002,'IfromS',0.004,'rewire',0.04);
network=GenerateNetwork(N,avgdegree);
nodes=zeros(N,1)+S;
nI_ini=round(N*0.8);
nodes(randperm(N,nI_ini))=I;

Setup of range for simulations

p=0.01:-0.001:0;  % parameter range to be scanned
nsteps=1000;      % number of time steps to be simulated for each parameter
fracI=NaN(nsteps,length(p)); % pre-allocated array for fraction of infected

Parameter sweep

clf
for k=1:length(p)
    rates.IfromS=p(k);
    for i=1:nsteps
        [network,nodes,ode]=SIstep(network,nodes,rates,'assert',true);
        fracI(i,k)=sum(nodes==I)/N;
    end
    fprintf('k=%d of %d,p=%g, I/N=%g\n',k,length(p),p(k),mean(fracI(:,k)));
    plot(reshape(fracI(:,1:k),1,k*nsteps));
    xlabel('steps');
    ylabel('fraction of infected nodes');
    title(sprintf(['network of N=%d nodes with step-wise ',...
        'decreasing p from p=%g to p=%g'],N,p(1),p(end)));
    drawnow
end
k=1 of 11,p=0.01, I/N=0.971416
k=2 of 11,p=0.009, I/N=0.971945
k=3 of 11,p=0.008, I/N=0.966145
k=4 of 11,p=0.007, I/N=0.960685
k=5 of 11,p=0.006, I/N=0.954369
k=6 of 11,p=0.005, I/N=0.94582
k=7 of 11,p=0.004, I/N=0.92678
k=8 of 11,p=0.003, I/N=0.900278
k=9 of 11,p=0.002, I/N=0.748936
k=10 of 11,p=0.001, I/N=0.371379
k=11 of 11,p=0, I/N=0.064215

Time profile of another sweep with longer wait for the transients

(nsteps=4000).

time profile of fraction of infected

Corresponding pseudo bifurcation diagram from the parameter sweep

(nsteps=400). Black dots correspond to the mean of the observed fraction of infected. The colour shading represents the histogram (darker= value visited more frequently).

single-parameter bifurcation diagram