NetworkX上带有代理的SIRS模型

Question

我想在NetworkX图上模拟和可视化SIRS模型（恢复/删除的节点可能再次变得敏感）。每个节点也作为代理运行，并且可以在每个时间步选择以概率p隔离，因此不能在该时间步中感染

我用EoN.Gillespie_simple_contraction模拟了一个SIRS模型，我试图了解我是否可以修改这些方法以包括代理行为，或者我是否需要编写一个定制方法

我发现可以修改一些EoN方法：Modified SIR model

以下是我正在使用的代码：

import networkx as nx
import matplotlib.pyplot as plt 
import EoN
from collections import defaultdict

# parameters required for the SIRS model
a = 0.1
b = 0.01
y = 0.001
d = 0.001

# Simple contagions
# the below is based on an example of a SEIR disease (there is an exposed state before becoming infectious)
# from https://arxiv.org/pdf/2001.02436.pdf

Gnp = nx.gnp_random_graph(500, 0.005)

H = nx.DiGraph() #For the spontaneous transitions
H.add_edge('I', 'R', rate = b)  # an infected node can be recovered/removed
H.add_edge('I', 'S', rate = y)  # an infected node can become susceptible again
H.add_edge('R', 'S', rate = d)  # a recovered node can become suscepticle again

J = nx.DiGraph() #for the induced transitions
J.add_edge(('I', 'S'),('I', 'I'), rate = a)  # a susceptible node can become infected from a neighbour
IC = defaultdict(lambda: 'S')

# set all statuses except one to susceptible. only one node shall begin infected
for node in range(500):
    IC[node] = 'S'
IC[0] = 'I'

return_statuses = ('S', 'I', 'R')
print('doing Gillespie simulation')

t, S, I, R = EoN.Gillespie_simple_contagion(Gnp, H, J, IC, return_statuses, tmax = 500)

print('done with simulation, now plotting')
plt.plot(t, S, label = 'Susceptible')
plt.plot(t, I, label = 'Infected')
plt.plot(t, R, label = 'Recovered')
plt.xlabel('$t$')
plt.ylabel('Number of nodes')
plt.legend()
plt.show()

Answer 1

为此，我们将引入一个新的节点状态'X'，以表示易受影响的隔离节点。我在'H'中添加了一条从'S'到'X'的边，另一条从'X'到'S'的边。我已经将'X'添加到return_statuses中，并绘制了它

为了简化代码，我删除了将所有内容最初分配为'S'的循环。因为它使用了一个default_dict，所以它的开头都是隐式的'S'

您提到了下一个“时间步”并以概率转移p。这不是代码中会发生的事情。代码是连续工作的。所以它需要一个利率。它并没有具体说明节点保持隔离的时间，而是给出了节点在任何时刻离开的可能性。如果这对你来说没有意义，你应该读一读指数分布

虽然我没有在这里使用它，但您可能想知道，EoN附带了允许您animate the simulation的工具

import networkx as nx
import matplotlib.pyplot as plt 
import EoN
from collections import defaultdict

a = 0.1
b = 0.01
y = 0.001
d = 0.001
to_isolation_rate = 0.05
from_isolation_rate = 1

# Simple contagions
# the below is based on an example of a SEIR disease (there is an exposed state before becoming infectious)
# from https://arxiv.org/pdf/2001.02436.pdf

Gnp = nx.gnp_random_graph(500, 0.005)

H = nx.DiGraph() #For the spontaneous transitions
H.add_edge('I', 'R', rate = b)  # an infected node can be recovered/removed
H.add_edge('I', 'S', rate = y)  # an infected node can become susceptible again
H.add_edge('R', 'S', rate = d)  # a recovered node can become suscepticle again
H.add_edge('S', 'X', rate = to_isolation_rate)
H.add_edge('X', 'S', rate = from_isolation_rate)

J = nx.DiGraph() #for the induced transitions
J.add_edge(('I', 'S'),('I', 'I'), rate = a)  # a susceptible node can become infected from a neighbour
IC = defaultdict(lambda: 'S')

IC[0] = 'I'

return_statuses = ('S', 'I', 'R', 'X')
print('doing Gillespie simulation')

t, S, I, R, X = EoN.Gillespie_simple_contagion(Gnp, H, J, IC, return_statuses, tmax = 500)

print('done with simulation, now plotting')
plt.plot(t, S, label = 'Susceptible')
plt.plot(t, I, label = 'Infected')
plt.plot(t, R, label = 'Recovered')
plt.plot(t, X, label = 'Isolating')
plt.xlabel('$t$')
plt.ylabel('Number of nodes')
plt.legend()
plt.show()

Answer 2

可以为所需功能编写定制方法：

def SIRS_simulation_with_quarantine(graph, params, tmax):
    a, b, y, d, q = params
    t=0
    graph_current_timestep = graph.copy()
    graph_previous_timestep = graph.copy()
    Infected_nodes = [x for x in graph_current_timestep.nodes() if graph_current_timestep.nodes[x]['state'][0]=='I']
    S = []
    I = []
    R = []
    T = []
    
    while(t<tmax):
        
        for node in graph.nodes():
            
            # initialise random probablilty for that node
            p = np.random.rand()
            # determine neighbours for that node
            neighbors = [n for n in graph_previous_timestep.neighbors(node)]
            
            # depending on the state, check for any state changes and update node attributes
            if graph_previous_timestep.nodes[node]['state'][0] == 'S':
                if (q >  np.random.rand()):                                             # if the node didnt quarantine...
                    if len(set.intersection(set(neighbors), set(Infected_nodes)))>0:    # if any neighbors are infected...
                        if (p < a):                                                     
                            attrs = {node: {'state': 'I'}}
                            nx.set_node_attributes(graph_current_timestep, attrs)                   
            elif graph_previous_timestep.nodes[node]['state'][0] == 'I':                
                if (p < b):
                    attrs = {node: {'state': 'R'}}
                    nx.set_node_attributes(graph_current_timestep, attrs)                    
                elif (p < y):
                    attrs = {node: {'state': 'S'}}
                    nx.set_node_attributes(graph_current_timestep, attrs)

            elif graph_previous_timestep.nodes[node]['state'][0] == 'R':                
                if (p < d):
                    attrs = {node: {'state': 'S'}}
                    nx.set_node_attributes(graph_current_timestep, attrs)

        # after the changes to states have occured for all nodes, record number of nodes at each state, increment timestep
        Infected_nodes = [x for x in graph_current_timestep.nodes() if graph_current_timestep.nodes[x]['state'][0]=='I']
        S.append(len([x for x in graph_current_timestep.nodes() if graph_current_timestep.nodes[x]['state'][0]=='S']))
        I.append(len(Infected_nodes))
        R.append(len([x for x in graph_current_timestep.nodes() if graph_current_timestep.nodes[x]['state'][0]=='R']))
        T.append(t)      
        graph_previous_timestep = graph_current_timestep.copy()
        t = t+1
    
    return T, S, I, R