- 在线时间
- 1630 小时
- 最后登录
- 2024-1-29
- 注册时间
- 2017-5-16
- 听众数
- 82
- 收听数
- 1
- 能力
- 120 分
- 体力
- 563405 点
- 威望
- 12 点
- 阅读权限
- 255
- 积分
- 174244
- 相册
- 1
- 日志
- 0
- 记录
- 0
- 帖子
- 5313
- 主题
- 5273
- 精华
- 3
- 分享
- 0
- 好友
- 163
TA的每日心情 | 开心 2021-8-11 17:59 |
|---|
签到天数: 17 天 [LV.4]偶尔看看III 网络挑战赛参赛者 网络挑战赛参赛者 - 自我介绍
- 本人女,毕业于内蒙古科技大学,担任文职专业,毕业专业英语。
 群组: 2018美赛大象算法课程 群组: 2018美赛护航培训课程 群组: 2019年 数学中国站长建 群组: 2019年数据分析师课程 群组: 2018年大象老师国赛优 |
Python实现简单的SI传播模型4 U* Y$ f% D- q: }2 m" b
#SI疾病传播模型的原理
* U: b& t3 u: A8 l7 H* ]在经典的传染病模型中,种群(Population)内N个个体的状态可分为如下几类
; ?6 {% |: W+ M& m: d* p0 o0 t
4 U0 Y. n0 G5 s6 s9 Z* i易感状态(Susceptible)。一个个体在感染前是处于易感状态的,即该个体有可能被邻居个体感染。
. b/ P2 n6 e' T9 v* @6 c易感状态I(Infected)。一个感染上某种病毒的个体就称为是处于感染状态。,即该个体还会以一定概率感染其邻居个体。
6 z7 L% d! L; F9 J2 `8 ` V. s移除状态(Remove,Refractory或者Recovered)。也成为免疫状态或恢复状态,当一个个体经历过一个完整的感染周期后,该个体就不再被感染,因此就可以不再考虑改革提。
! }3 l8 e4 x4 _SI传播模型是最简单的疾病传播模型,模型中的所有个体都只可能处于两个状态中的一个
$ r' R2 _/ t3 k6 `/ S- ?' i即易感(S)状态或感染(I)状态。SI模型中的个体一旦被感染后就永远处于感染状态。
) R+ E/ j7 V. L9 @# I在给定时刻t,令S(t)与I(t)分别代表该时刻处于易感和感染状态的个体数目,显然有
- \) P: ?2 W8 C6 \' F' AS(t)+I(t)恒等于N,这里,N是个体总数。随着时间t的增长,易感个体与感染个体的接触
: M4 \) X# o9 q0 s2 i会导致感染个体数量的增加。加入由于个体之间的接触而导致疾病传播的概率为β,疾病仅在2 l" {6 h( J( s; a
感染个体和易感个体之间进行接触时才会以概率β将疾病传染给易感个体。在时刻t,易感个体的比例为S(t)/N,感染个体的数量为I(t),一次,易感个体的数量将以如下变化率减少
2 W" c) ]/ z! c* o. I0 t& \6 yds/dt = -β*S(t)I(t)/N
" ]% D; o6 k+ {7 |2 }. ?8 V2 ?同时,感染个体的数量会以与易感个体相反的变化率增加,
+ H6 {( ?4 ]/ n( wds/dt = βS(t)*I(t)/N# D) s: w! d" Q" l+ l, l; e3 P
分别将时刻t处于易感状态和感染状态的个体所占比例记为,1 P" `- \6 ]8 M. v! i; C$ }' a
s(t)=S(t)/N
2 M# N9 l6 \$ N E- Si(t)=I(t)/N* w5 }# B( [. x0 K4 Y
显然有,% p" c7 D$ j; ]
s(t)+i(t)恒等于1,此时之前的公式可以记做1 W4 {1 \& _/ x4 H7 F& J
ds/dt=-βsi0 W% I r; k/ Z$ M
di/dt=βsi
, y: d% L! I" p& H: R" ]2 m; r即
1 _' l9 F. M- {+ I7 Qdi/dt=βi(1-i)
, J8 r. U- t# g6 x* r' T) y上式也成为Logistic增长方程式(Logistic growth equation),0 e! y) K$ ?2 l8 s) c( P* l& K
方程的解和图像如图 `) q) a6 }& C( J7 ]6 v, a% Z; c
$ |/ {8 [ {$ O( V: S) s
代码和相关文件以及环境链接:链接:https://pan.baidu.com/s/1JSfHuTPaglFimeEBLdSDyQ/ D9 P. I1 a; C& ^
提取码:z448
8 P' r0 J7 H, ?" N
/ R) d6 ]1 q( z4 S/ z+ v1 F4 ^. w0 ^
'''
, M" l; G5 L2 \# d; ^实验环境Python2.7.13,igraph包,cairo包,numpy包# ~+ ?# P9 T% z9 Q. A, U
'''; D: p3 x1 F) z0 f; y! x
# -*- coding:utf8 -*! }! i, u" J! }
from igraph import *
' l' m4 M, [' W9 {0 w: G, X& Fimport numpy as numpy7 g# G8 N+ f$ w H1 A! G$ ^
from numpy import *
# b% S0 [5 g! I V) {import random
( q- k8 D" B B5 S7 S" n$ q9 @' O4 j5 R; s
def len_arr(infected_array,nodes_num):#获取感染数组长度! R) t4 F: G) b& \; B
len_value=0#初始化长度* Z, m3 s1 a. ^: }
len_value=nodes_num-infected_array.count(-1)#被感染数量是结点总数减去未感染节点数(未感染的结点被标记为-1)9 @1 n8 U( e9 _3 q
return len_value8 O8 Y% I+ Y' x" ^/ o3 _
3 a, J8 s$ i. o. y% Xg=Graph.Read_GML("C:\python27\e1.gml")#将本地保存的网络数据读入变量g(生成图)3 A, y7 f* E7 F8 \4 d# w, p1 t( Q n6 `" U
summary(g)
; K5 p2 W' C3 dnodes_num=g.vcount()#统计图中的结点个数
6 _- F# e. f$ _2 F& Vnet_mat=g.get_adjacency(type=GET_ADJACENCY_BOTH)#将网络数据转换为邻接矩阵存储在变量net_mat
( q4 N9 ?1 Y5 w* y6 l& K" dg.vs["color"]=["white"]#给图的顶点序列颜色赋值白色
! Z/ m) {& C J5 I# t9 i0 ^; w$ Ra=[arange(nodes_num)+1]*3#声明一个N行3列的数组a
$ @, F/ [! `* f' wnodes_state=matrix(a).T#nodes_state通过转置a矩阵创建,用于存放每个节点的状态信息以及其被感染的时间(这个是理解算法的重中之重!!!)
9 l8 X9 j/ o& U* d #第一列是节点编号,第二列是节点状态,感染状态用-2表示,第三列是节点感染的时间1 S* N! h) K' a
print(nodes_state)9 t% z5 p* B. s- l2 A- [& |
infected_array=[-1]*34#用于存放本轮被感染的结点, 这些结点将参与下一次感染 34代表网络节点数' T8 \) [/ F5 |
print(infected_array)
2 u2 H9 E' p2 p2 b( s2 Z4 g
/ u: W" z; v) T' }+ ~- kinfe_rate=1#传播率(感染率) 1代表邻接点100%被感染6 u/ ?0 j9 N2 Q" }$ \9 [1 \- a, P
set_time=2#传播次数(感染次数) 2次
" A$ e9 G0 k) G) qsource_seed=1#感染源位置
- B% X7 z8 h3 y" unodes_state[0:nodes_num,2]=-1#给所有节点初始化感染时间为-1. j* B; f2 D# D' v! {' {% [# h
nodes_state[source_seed-1,1]=-2#设置第一个感染源感染状态 -2代表感染状态5 B0 T1 }& _+ t5 S# p
nodes_state[source_seed-1,2]=1#设置第一个感染源的感染时间为1: b1 A, Y2 j0 f
g.vs[source_seed-1]["color"]="red"#将感染的顶点颜色标红
7 u1 Z1 x4 n' O: g# p: rinfected_array[0]=source_seed#将感染源的位置存入被感染节点列表
" P$ C# ^$ c6 J$ k rplot(g)#绘制- o. |8 O1 J4 E p0 r3 @6 }
3 Z& s6 x9 p$ Tstop=False#感染过程结束的标记5 \& M- q$ j) Q- ~# {1 Y' T8 _4 z
temp_time=0#第几次感染
/ A. m) ^- `: T3 ttemp_len=0#本轮的感染源数量初始化/ \$ E2 M9 a$ N0 Y
, i: u. X1 ~7 T# B4 e4 U: uwhile not stop:
1 V6 I6 }2 g; d' A$ O3 O i=0#记录让每个感染源都传播一次& b2 x) m/ V5 ^$ f: }& m) n
if len_arr(infected_array,nodes_num)>0 and len_arr(infected_array,nodes_num)<=nodes_num:#感染可以进行
% I& d% E0 w% D6 p& C1 ]2 v temp_len=len_arr(infected_array,nodes_num)#获取本轮的感染源数量
# V( M! p( U- s, X5 E while i<temp_len:# D/ S, r {) b! y* @
temp_time=nodes_state[infected_array-1,2]#获取每一个节点的感染时间
$ ?- @8 E9 v' B7 q- B! s( h nei_count=0#下一轮可以被感染到的节点数量! [5 w- o! C; r. t, Z1 N0 v; y6 F1 E
#生成下一轮可能被感染的节点的集合nei_arr
5 s5 H1 `$ u k# w for j in range(nodes_num):#遍历节点
: e/ d; z% L+ |0 l. Q if net_mat[infected_array-1,j]==1 and nodes_state[j,1]!=-2:#是邻接节点而且未被感染
4 [4 c8 }4 u# f5 h- e% g nei_count=nei_count+1#下一轮可以被感染到的节点数量++
% {, i$ e9 Z& o nei_arr=[-1]*nei_count#用于临时存放本轮被感染的结点, 这些结点将参与下一次感染) I0 j5 R3 `, \* G4 ?* h
t=0
3 j, _) P/ i9 G3 ?: u for j in range(nodes_num):! z- t( N2 P, _- X1 d4 D
if net_mat[infected_array-1,j]==1 and nodes_state[j,1]!=-2:* D+ `5 S" m+ ~' w7 c r
nei_arr[t]=j+1
! Z- C! a+ R) E t=t+1
2 V, w7 I' T, v+ J ran_infe_arr=random.sample(range(nei_count),int(nei_count*infe_rate))#随机生成会被感染的节点的数组
6 m2 H6 \0 k' h' E* S! o #random.simple(arg1,num) 从arg1集合中随机取num个数据生成一个对象. h6 [. ~$ }8 W/ L7 d
if len(ran_infe_arr)>0:#存在需要被感染的节点5 N) c, I) L% l. T8 |
t=0#让ran_infe_arr内每个感染源都被感染
5 R7 z/ x8 k: I! [ while t<len(ran_infe_arr):#对刚才生成的会被感染的数组内的节点进行感染6 j1 G* \% w8 e) m" S
nodes_state[nei_arr[ran_infe_arr[t]]-1,1]=-2#标记为感染状态
. |; |; E* L. `/ s nodes_state[nei_arr[ran_infe_arr[t]]-1,2]=temp_time+1#记录感染时间
$ l3 ]' d& d$ b infected_array[len_arr(infected_array,nodes_num)]=nei_arr[ran_infe_arr[t]]#将此次感染节点放入总的感染节点数组中3 v0 N$ j' D. O7 X! w3 Y; U
g.vs[nei_arr[ran_infe_arr[t]]-1]["color"]="pink"#将此次感染的节点集的所有节点颜色置为粉色
( n. l" \* a- j8 E( O- I' | plot(g)#绘制
7 b( B* Q( a6 {9 o- x) f6 ^, Q t=t+1
. M- ?' c, x4 l! ~ i=i+1
( q# z$ ?, O! N if temp_time>set_time-1:#当执行感染的次数等于设置的次数结束感染, U( C2 {% l5 x
stop=True
+ Z/ R. P6 \' O# J3 m9 E+ V# q4 G3 z
7 J+ c7 L. k$ Q/ G视频演示bilibili传送门; J- j* n$ \- m7 ]* U& @; b
效果图
4 R9 S5 K7 S' B8 y0 M1 E
* h/ N' g/ `! m7 s; Y# t+ l5 ]) T4 F J9 d5 y6 _' `- P
- M4 i$ W4 [5 h" U7 o) p7 [ R. e( Z- W( o7 ]
+ B! D9 q% _. B% s8 q
, ^ p% W6 n% J) f, V! [+ |4 u3 \
' C/ U. x) J( _8 e' i. B3 S& B. z
6 }9 }# A2 P3 A9 L
( x6 U: d7 `0 c7 e4 a
, T! P5 V l- C% E
$ x$ i0 d& `$ t, T
( D! H; u0 N1 L& l& g# d
————————————————
; O$ ]* W- m$ C# C" d1 b( e( \8 ^版权声明:本文为CSDN博主「eck_燃」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。- @8 p* ~& f! h( H" r
原文链接:https://blog.csdn.net/wdays83892469/article/details/808788622 h% a% X/ V# q2 a5 `
% r+ g$ _# e# W; I( A
S3 n8 }+ |5 B1 \+ I, k
|
zan
|