修改后的SIR模型 [英] Modified SIR model

查看:88
本文介绍了修改后的SIR模型的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我正在制作带有增加的疫苗接种参数V的修改后的SIR模型.最初,图中的所有节点都是易受感染的,并且有一些最初感染的人.最初感染邻居的人首先被接种了概率w(这意味着他们不能被感染),然后又被感染了概率b.接种疫苗的总人数由V1控制,而V1是总人口的一小部分.

I am making a modified SIR model with an added vaccination parameter V. InitIALLY all the nodes in the graph are susceptible and there are a few initially infected people. The initial infected people neighbours are first vaccinated with prob w (which means they cant be infected) and then they are infected with prob b. The total number of vaccinated people is controlled by Vl which is a fraction of the total population.

这是我的代码-

import networkx as nx
import random
import scipy
from collections import defaultdict
from collections import Counter
import matplotlib.pyplot as plt
import numpy as np
from statistics import mean
from random import choice
from random import sample

def test_transmission(u, v, p):

    return random.random()<p


def discrete_SIR(G,
                initial_infecteds=None, initial_recovereds = None,beta=0.8,
                rho = None, w=0.5,Vl=1,tmin = 0, tmax = 100000,
                return_full_data = False):


    if G.has_node(initial_infecteds):
        initial_infecteds=[initial_infecteds]           

    N=  G.order()
    t = [tmin]
    S = [N-len(initial_infecteds)]
    I = [len(initial_infecteds)]
    R = [0]
    V = [0]

    susceptible = defaultdict(lambda: True)  
    #above line is equivalent to u.susceptible=True for all nodes.

    for u in initial_infecteds:
        susceptible[u] = False
    if initial_recovereds is not None:
        for u in initial_recovereds:
            susceptible[u] = False

    infecteds = set(initial_infecteds)
    print('len of infected initially',len(infecteds))

    while infecteds and t[-1]<tmax :
        print('len of infected on each iter',len(infecteds))
        new_infecteds = set()
        vaccinated = set()
        #infector = {}  #used for returning full data.  a waste of time otherwise
        for u in infecteds:
            print('u-->' +str(u))
            for v in G.neighbors(u):
                print('v --> '+ str(v))
            ##vaccination
                if len(vaccinated)+V[-1]< (Vl*N)  : #check if vaccination over or not
                #V.append(V[-1]+len(vaccinated))< (Vl*N)
                    #print(len(vaccinated),Vl*N)
                    print("HI")
                    print(V[-1])

                    if susceptible[v] and test_transmission(u, v, w):

                        vaccinated.add(v)
                        susceptible[v] = False
                        print('transmitting vaccination')

                    elif susceptible[v] and test_transmission(u,v,beta):
                        new_infecteds.add(v)
                        susceptible[v]=False
                        print('transmitting infection')
                else:

                    print("BYE")
                    if susceptible[v] and test_transmission(u, v,beta): 
                        new_infecteds.add(v)
                        susceptible[v] = False
                        #infector[v] = [u]



        infecteds = new_infecteds

        R.append(R[-1]+I[-1])
        V.append(len(vaccinated)+V[-1])
        I.append(len(infecteds))
        S.append(N-V[-1]-I[-1]-R[-1])
        #S.append(S[-1]-V[-1]-I[-1])
        t.append(t[-1]+1)

        print('\n')
        print('time is',str(t) +' --> ')
        print('infected is',I)
        print('sum is',R[-1]+V[-1]+I[-1]+S[-1])
        print('R V I S',str(R[-1])+','+str(V[-1])+','+str(I[-1])+','+str(S[-1]))
        print('time t[-1]',t[-1])


    if not return_full_data:
        return scipy.array(t), scipy.array(S),scipy.array(V), scipy.array(I), \
               scipy.array(R)





m=100
G=nx.grid_2d_graph(m,m,periodic=True)

def avg_deg(self,num_nodes):
        return self.number_of_edges() * 2 / num_nodes


def avg_degree(num_nodes,target_deg):

    G=nx.Graph()

    G.add_nodes_from(range(num_nodes))
    while avg_deg(G,num_nodes) < target_deg:
        n1, n2 = sample(G.nodes(), 2)
        G.add_edge(n1, n2, weight=1)

    nx.draw(G)
    plt.show()
    return G    


initial_infections = [(u,v) for (u,v) in G if u==int(m/2) and v==int(m/2)]  

t, S, V, I, R = discrete_SIR(G,initial_infecteds= initial_infections,beta=0.8,w=0.05,Vl=0.1)  

plt.figure()
plt.plot(t,I,ls='-',color='red')
plt.plot(t,R,ls='--',color='orange')
plt.plot(t,S,ls='--',color='green')
plt.plot(t,V,ls='-',color='black')
plt.show()

我的代码存在问题,即S + V + I + R的总数应等于N,并且接种人数最多也要达到5.该数量应该更高.

The problem with my code is the total number of S + V+ I+R should be equal to N and the number of vaccinated people is also at max coming to be 5. It should be higher than that.

推荐答案

您的问题是,V仅衡量该时间步中要接种的个体总数,而不是累计接种量.

Your problem is that V is only measuring the total number of individuals becoming vaccinated in that time step, rather than the cumulative number of vaccinations.

所以 S + I + V + R 不是常数.如果希望 V 为累计疫苗接种数,请执行 V.append(V [-1] + len(vaccinated))

So S+I+V+R is not constant. If you want V to be the cumulative number of vaccinations, then do V.append(V[-1]+len(vaccinated))

如果len(已接种疫苗),我不确定此测试<(Vl * N)也在做您想的.它正在检查在给定的时间范围内接种疫苗的人数是否少于总人口的几分之一.我怀疑您要使用累计疫苗接种数.

I'm not sure that this test if len(vaccinated)< (Vl*N) is doing what you think as well. It is checking whether the number of people vaccinated in a given time step is less than some fraction of the entire population. I suspect you want to use the cumulative number of vaccinations.

这篇关于修改后的SIR模型的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆