数学建模社区-数学中国

标题: Python实现简单的SI传播模型 [打印本页]

作者: 杨利霞 时间: 2020-4-18 16:16
标题: Python实现简单的SI传播模型
Python实现简单的SI传播模型
#SI疾病传播模型的原理
在经典的传染病模型中，种群(Population)内N个个体的状态可分为如下几类

易感状态(Susceptible)。一个个体在感染前是处于易感状态的，即该个体有可能被邻居个体感染。
易感状态I(Infected)。一个感染上某种病毒的个体就称为是处于感染状态。，即该个体还会以一定概率感染其邻居个体。
移除状态(Remove,Refractory或者Recovered)。也成为免疫状态或恢复状态，当一个个体经历过一个完整的感染周期后，该个体就不再被感染，因此就可以不再考虑改革提。
SI传播模型是最简单的疾病传播模型，模型中的所有个体都只可能处于两个状态中的一个
即易感(S)状态或感染(I)状态。SI模型中的个体一旦被感染后就永远处于感染状态。
在给定时刻t，令S(t)与I(t)分别代表该时刻处于易感和感染状态的个体数目，显然有
S(t)+I(t)恒等于N，这里，N是个体总数。随着时间t的增长，易感个体与感染个体的接触
会导致感染个体数量的增加。加入由于个体之间的接触而导致疾病传播的概率为β，疾病仅在
感染个体和易感个体之间进行接触时才会以概率β将疾病传染给易感个体。在时刻t，易感个体的比例为S(t)/N，感染个体的数量为I(t)，一次，易感个体的数量将以如下变化率减少
ds/dt = -β*S(t)I(t)/N
同时，感染个体的数量会以与易感个体相反的变化率增加，
ds/dt = βS(t)*I(t)/N
分别将时刻t处于易感状态和感染状态的个体所占比例记为，
s(t)=S(t)/N
i(t)=I(t)/N
显然有，
s(t)+i(t)恒等于1，此时之前的公式可以记做
ds/dt=-βsi
di/dt=βsi
即
di/dt=βi(1-i)
上式也成为Logistic增长方程式(Logistic growth equation),
方程的解和图像如图

代码和相关文件以及环境链接：链接：https://pan.baidu.com/s/1JSfHuTPaglFimeEBLdSDyQ
提取码：z448

'''
实验环境Python2.7.13,igraph包，cairo包，numpy包
'''
# -*- coding:utf8 -*
from igraph import *
import numpy as numpy
from  numpy import *
import random

def len_arr(infected_array,nodes_num):#获取感染数组长度
len_value=0#初始化长度
len_value=nodes_num-infected_array.count(-1)#被感染数量是结点总数减去未感染节点数(未感染的结点被标记为-1)
return len_value

g=Graph.Read_GML("C:\python27\e1.gml")#将本地保存的网络数据读入变量g(生成图)
summary(g)
nodes_num=g.vcount()#统计图中的结点个数
net_mat=g.get_adjacency(type=GET_ADJACENCY_BOTH)#将网络数据转换为邻接矩阵存储在变量net_mat
g.vs["color"]=["white"]#给图的顶点序列颜色赋值白色
a=[arange(nodes_num)+1]*3#声明一个N行3列的数组a
nodes_state=matrix(a).T#nodes_state通过转置a矩阵创建，用于存放每个节点的状态信息以及其被感染的时间(这个是理解算法的重中之重！！！)
                                                #第一列是节点编号，第二列是节点状态，感染状态用-2表示，第三列是节点感染的时间
print(nodes_state)
infected_array=[-1]*34#用于存放本轮被感染的结点, 这些结点将参与下一次感染 34代表网络节点数
print(infected_array)

infe_rate=1#传播率(感染率) 1代表邻接点100%被感染
set_time=2#传播次数(感染次数) 2次
source_seed=1#感染源位置
nodes_state[0:nodes_num,2]=-1#给所有节点初始化感染时间为-1
nodes_state[source_seed-1,1]=-2#设置第一个感染源感染状态 -2代表感染状态
nodes_state[source_seed-1,2]=1#设置第一个感染源的感染时间为1
g.vs[source_seed-1]["color"]="red"#将感染的顶点颜色标红
infected_array[0]=source_seed#将感染源的位置存入被感染节点列表
plot(g)#绘制

stop=False#感染过程结束的标记
temp_time=0#第几次感染
temp_len=0#本轮的感染源数量初始化

while not stop:
i=0#记录让每个感染源都传播一次
if len_arr(infected_array,nodes_num)>0 and len_arr(infected_array,nodes_num)<=nodes_num:#感染可以进行
      temp_len=len_arr(infected_array,nodes_num)#获取本轮的感染源数量
      while i<temp_len:
         temp_time=nodes_state[infected_array-1,2]#获取每一个节点的感染时间4 p* Y- l; J( Z! S
         nei_count=0#下一轮可以被感染到的节点数量6 T: d  N' @" U/ n2 L
         #生成下一轮可能被感染的节点的集合nei_arr
8 v' l3 P  ~+ w1 ?; z2 \          for j in range(nodes_num):#遍历节点
, w( R  b- ?' S1 J             if net_mat[infected_array-1,j]==1 and nodes_state[j,1]!=-2:#是邻接节点而且未被感染
% c& ]! V" H- p- u/ ?; m. {                   nei_count=nei_count+1#下一轮可以被感染到的节点数量++) R) m4 z6 |3 I. @2 J# s* X/ m* f
         nei_arr=[-1]*nei_count#用于临时存放本轮被感染的结点, 这些结点将参与下一次感染3 Q/ e! ~2 t* Y6 A$ @& t
         t=0
         for j in range(nodes_num):$ [- S, g% h+ Z" x5 a5 z. _
            if net_mat[infected_array-1,j]==1 and nodes_state[j,1]!=-2:
9 W# t: N" t& @5 Q                   nei_arr[t]=j+15 L, j+ J' t/ j2 K2 f
                  t=t+1
; E- o8 G% e  N2 w' D  }          ran_infe_arr=random.sample(range(nei_count),int(nei_count*infe_rate))#随机生成会被感染的节点的数组  W- ]  g- V+ x, N/ s% X
                                    #random.simple(arg1,num) 从arg1集合中随机取num个数据生成一个对象2 B6 t, r. u) O
         if len(ran_infe_arr)>0:#存在需要被感染的节点9 _. A7 b3 \; X! u* K
            t=0#让ran_infe_arr内每个感染源都被感染
( n( J5 n6 A8 J3 o! f             while t<len(ran_infe_arr):#对刚才生成的会被感染的数组内的节点进行感染9 R8 @1 Z: [$ Y( B& e+ i* O5 t
                  nodes_state[nei_arr[ran_infe_arr[t]]-1,1]=-2#标记为感染状态9 R5 _# D( r7 ?! {  P' l
                  nodes_state[nei_arr[ran_infe_arr[t]]-1,2]=temp_time+1#记录感染时间3 L$ d% ^1 i/ N
                  infected_array[len_arr(infected_array,nodes_num)]=nei_arr[ran_infe_arr[t]]#将此次感染节点放入总的感染节点数组中7 P, B3 r! n+ J  b9 j: x# ]: V
                  g.vs[nei_arr[ran_infe_arr[t]]-1]["color"]="pink"#将此次感染的节点集的所有节点颜色置为粉色* K" X- n! z' d' F
                  plot(g)#绘制
: b: |. W! i" N3 @% h: P  U$ _                   t=t+1
& B/ Z0 l- T5 u          i=i+13 `2 R6 t3 P5 v4 A. x" g
if temp_time>set_time-1:#当执行感染的次数等于设置的次数结束感染4 _2 z# h2 B/ z0 D
      stop=True ' X5 L) v7 t. ~0 g* P9 v, D3 `+ V' M

) |# b  w) r! U/ J/ K# ?7 u4 |3 n% u! _; u; ]5 d6 g. S! N3 _$ s
视频演示bilibili传送门
# B/ Y3 y5 F/ M6 N效果图- G" I+ _% W* ^% m9 v" t5 Z! @

! P* {" Y* H2 I/ A
; }  \( Q/ v4 b: A2 _. n4 ~; F

2020-4-18 16:15 上传
下载附件 (174.22 KB)

. ?2 [  c- S$ Z6 t' _6 z8 D
  j  E3 |) G; q, {- }

2020-4-18 16:15 上传
下载附件 (22.44 KB)

! W! t! {4 |* Y5 E- c7 I0 K9 t9 ?9 ~

% m' M0 e7 g* k3 N

2020-4-18 16:16 上传
下载附件 (42.85 KB)

; z0 t! y' q( C6 V5 f1 p

9 ^( u$ X$ K, C' C# g

2020-4-18 16:16 上传
下载附件 (17.44 KB)

! e+ G% X- ~8 |' U6 `" P. s, L" S$ g  V& n, e4 h

2020-4-18 16:16 上传
下载附件 (19.78 KB)

8 X4 g4 l2 Y+ s* }/ }$ N; m7 F' F

2020-4-18 16:16 上传
下载附件 (21.16 KB)

————————————————. R! l0 g' z* Z( c4 w/ ]1 b
版权声明：本文为CSDN博主「eck_燃」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。) |) p- p4 e, i- {5 S- U
原文链接：https://blog.csdn.net/wdays83892469/article/details/80878862
9 [* [+ P: i- B6 g1 V. X, ~
5 {0 W2 m8 v8 d& S5 u- y/ A2 I, {! Q/ b. W# |0 T% J% ?8 F* t

作者: 尔雅 时间: 2020-4-18 17:50
发表回复棒!- f  U5 l5 B6 ^3 B9 z

欢迎光临数学建模社区-数学中国 (http://www.madio.net/) Powered by Discuz! X2.5