- 在线时间
- 1630 小时
- 最后登录
- 2024-1-29
- 注册时间
- 2017-5-16
- 听众数
- 82
- 收听数
- 1
- 能力
- 120 分
- 体力
- 563408 点
- 威望
- 12 点
- 阅读权限
- 255
- 积分
- 174245
- 相册
- 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传播模型6 r3 }4 C" ]0 [8 n
#SI疾病传播模型的原理: f, n' K& U2 s5 B4 C0 i
在经典的传染病模型中,种群(Population)内N个个体的状态可分为如下几类
1 i" p( O9 i9 U7 s- G8 y+ c
, h- m( i' u1 v7 P8 B易感状态(Susceptible)。一个个体在感染前是处于易感状态的,即该个体有可能被邻居个体感染。+ h) |% j' a, j
易感状态I(Infected)。一个感染上某种病毒的个体就称为是处于感染状态。,即该个体还会以一定概率感染其邻居个体。, G$ E L5 Q4 z4 H& D
移除状态(Remove,Refractory或者Recovered)。也成为免疫状态或恢复状态,当一个个体经历过一个完整的感染周期后,该个体就不再被感染,因此就可以不再考虑改革提。
( i4 w' n3 T9 hSI传播模型是最简单的疾病传播模型,模型中的所有个体都只可能处于两个状态中的一个1 x; w6 B: W$ v& c# M4 n; f- |# L8 q
即易感(S)状态或感染(I)状态。SI模型中的个体一旦被感染后就永远处于感染状态。8 a( Z$ h6 ], V3 v0 b0 q0 ?5 P
在给定时刻t,令S(t)与I(t)分别代表该时刻处于易感和感染状态的个体数目,显然有) g2 U) s; g: V
S(t)+I(t)恒等于N,这里,N是个体总数。随着时间t的增长,易感个体与感染个体的接触9 Q9 ^0 L0 u# C5 _8 z
会导致感染个体数量的增加。加入由于个体之间的接触而导致疾病传播的概率为β,疾病仅在
- e! ?6 F1 W; E' K) B5 _感染个体和易感个体之间进行接触时才会以概率β将疾病传染给易感个体。在时刻t,易感个体的比例为S(t)/N,感染个体的数量为I(t),一次,易感个体的数量将以如下变化率减少
: V. \( f+ N% w* E! Ads/dt = -β*S(t)I(t)/N
1 O2 r# ?7 L9 w) |同时,感染个体的数量会以与易感个体相反的变化率增加,+ X2 |& o( O2 o; I! m: @
ds/dt = βS(t)*I(t)/N
2 Y7 L0 k% ]& ^2 _8 k分别将时刻t处于易感状态和感染状态的个体所占比例记为,
7 p: L" ?' _/ F, j/ D/ js(t)=S(t)/N
: r: q7 E- N" U* ?4 ^i(t)=I(t)/N& ^ O6 X2 k: |" F+ M) w" C9 V
显然有,
3 E$ F& x5 s9 ^; f! D/ is(t)+i(t)恒等于1,此时之前的公式可以记做; d) X/ g' L& C* [: W
ds/dt=-βsi" J* U. l3 v8 d
di/dt=βsi
6 G% `- u* x1 f- F# X即( E9 u4 s7 y* S
di/dt=βi(1-i)3 B; t, _0 X$ B( A1 ^
上式也成为Logistic增长方程式(Logistic growth equation),% g8 O2 P: F" j8 ^
方程的解和图像如图# m5 X8 @+ Q6 W
8 O9 U) n: [ }4 v/ ?* e5 O
代码和相关文件以及环境链接:链接:https://pan.baidu.com/s/1JSfHuTPaglFimeEBLdSDyQ
2 ]- X/ p' j# t j. _- [0 ~2 A& t提取码:z448
j% F) o" ?1 P9 \4 O! S0 r k8 D- X k% c4 r7 s; D" `% k
- D* e+ I r8 Y! s' r" I# i9 O
'''
" W$ j4 c( ^- ]实验环境Python2.7.13,igraph包,cairo包,numpy包
7 a$ n+ g- F+ r. _% S'''8 ~6 t0 T, f: N K( G$ g
# -*- coding:utf8 -*% `. b" t3 W7 g* Q! |) S
from igraph import *# ~3 y" p% N, J
import numpy as numpy6 q: d5 d* x! U- X
from numpy import *" x# ^- \; g) y1 p4 u1 A, p
import random$ F) q3 k; j8 ?/ x3 F: U( ~3 K2 J
- o7 i, U' E% l) O: v+ [6 Gdef len_arr(infected_array,nodes_num):#获取感染数组长度
% z0 R& ]% ~6 F% q1 T6 p3 E len_value=0#初始化长度
/ ]5 b# M- c9 n2 A) E- W0 Q len_value=nodes_num-infected_array.count(-1)#被感染数量是结点总数减去未感染节点数(未感染的结点被标记为-1)
# `3 _- N0 |) P2 `. J return len_value: L3 Z9 H( c" o2 r* _
; g7 [, a- @3 F" n% r" y8 E9 rg=Graph.Read_GML("C:\python27\e1.gml")#将本地保存的网络数据读入变量g(生成图)
5 m; c% W6 @, y% hsummary(g)& K( u( B6 P! g
nodes_num=g.vcount()#统计图中的结点个数
# X2 N b3 C. t( K* Knet_mat=g.get_adjacency(type=GET_ADJACENCY_BOTH)#将网络数据转换为邻接矩阵存储在变量net_mat
. l, K. j% a. I# `9 Ag.vs["color"]=["white"]#给图的顶点序列颜色赋值白色
! H- j/ _0 ?8 C5 F# Q! B5 Y0 ]0 W a: M9 ba=[arange(nodes_num)+1]*3#声明一个N行3列的数组a
* I3 P: \. g/ e* Rnodes_state=matrix(a).T#nodes_state通过转置a矩阵创建,用于存放每个节点的状态信息以及其被感染的时间(这个是理解算法的重中之重!!!)
: c3 _3 d8 V6 G+ T7 J5 o #第一列是节点编号,第二列是节点状态,感染状态用-2表示,第三列是节点感染的时间3 {+ y! L/ ]& n6 e4 {/ j8 ~1 a6 M' ?1 R% X
print(nodes_state)
' c+ i- k, h% M/ o" v6 D- Tinfected_array=[-1]*34#用于存放本轮被感染的结点, 这些结点将参与下一次感染 34代表网络节点数
, b! s/ W: `+ }) L6 b) \4 tprint(infected_array)
4 S6 e, M* o' |9 ^( d: h8 R y* `' Y
infe_rate=1#传播率(感染率) 1代表邻接点100%被感染
! S, i* F6 C9 f4 @8 N' r! u& }set_time=2#传播次数(感染次数) 2次7 d& \1 A0 u# @+ Z. n/ {
source_seed=1#感染源位置
; l% Q% N; Y# n+ Z# J6 onodes_state[0:nodes_num,2]=-1#给所有节点初始化感染时间为-1
3 A0 \: B2 H) p* z8 Z ^" v. N8 Ynodes_state[source_seed-1,1]=-2#设置第一个感染源感染状态 -2代表感染状态
/ O, D7 {! _/ xnodes_state[source_seed-1,2]=1#设置第一个感染源的感染时间为1
7 X$ S! z% N: d; yg.vs[source_seed-1]["color"]="red"#将感染的顶点颜色标红
1 R* b8 a4 y1 i Y$ finfected_array[0]=source_seed#将感染源的位置存入被感染节点列表
4 k$ x! u8 X5 Pplot(g)#绘制
7 n# d" U' x0 U) ^8 C9 j! X
# \, C- S0 Z1 T) a/ w! Istop=False#感染过程结束的标记
- J u. \% z# B) {, o) _temp_time=0#第几次感染
; P% V# ?: o( Gtemp_len=0#本轮的感染源数量初始化
" f# x, f- ~0 m* N
' t) j2 F! Y ~( `1 i4 ]while not stop:
1 j( q+ w9 e1 F6 w. |; m" Z: ] i=0#记录让每个感染源都传播一次2 a. s. O* K7 c% O# I* O( L
if len_arr(infected_array,nodes_num)>0 and len_arr(infected_array,nodes_num)<=nodes_num:#感染可以进行
% R# X% e2 f: h/ N& n" c temp_len=len_arr(infected_array,nodes_num)#获取本轮的感染源数量
( _3 P- `1 K" E! `% w' T while i<temp_len:% x+ W, Y! \2 J$ v% m
temp_time=nodes_state[infected_array-1,2]#获取每一个节点的感染时间
& d/ v- d. s, N$ k: o' x nei_count=0#下一轮可以被感染到的节点数量6 \! p% R2 p1 o) ]! {
#生成下一轮可能被感染的节点的集合nei_arr6 e# q- x1 e7 D6 ?6 }2 n( p' M9 ~
for j in range(nodes_num):#遍历节点, G+ V4 |; L& \0 g" j3 X
if net_mat[infected_array-1,j]==1 and nodes_state[j,1]!=-2:#是邻接节点而且未被感染' q: ^* s3 j1 W C9 a, m2 q! U! L
nei_count=nei_count+1#下一轮可以被感染到的节点数量++
" U. m6 A4 R. ^+ s# y7 z nei_arr=[-1]*nei_count#用于临时存放本轮被感染的结点, 这些结点将参与下一次感染. D+ e, F. a7 N1 E" r! U! L7 N
t=05 M6 z% v8 m+ j/ N
for j in range(nodes_num):, ]5 m" P! |$ x5 ?) o _4 @
if net_mat[infected_array-1,j]==1 and nodes_state[j,1]!=-2:" K8 ^$ r( M1 J0 s) {6 `
nei_arr[t]=j+1
& c8 h* a' U+ m" G. B t=t+1
1 f5 q& _( O9 T* d) \9 J ran_infe_arr=random.sample(range(nei_count),int(nei_count*infe_rate))#随机生成会被感染的节点的数组
" H) j' A+ a B4 p5 s #random.simple(arg1,num) 从arg1集合中随机取num个数据生成一个对象
3 l. V% k+ _; G if len(ran_infe_arr)>0:#存在需要被感染的节点
: P7 d+ I. _+ ? t=0#让ran_infe_arr内每个感染源都被感染- ~) @+ P( F# p9 x7 P
while t<len(ran_infe_arr):#对刚才生成的会被感染的数组内的节点进行感染
! }5 l$ g# J7 S! T* G5 @ nodes_state[nei_arr[ran_infe_arr[t]]-1,1]=-2#标记为感染状态- e- \. a" ?9 m1 T) b, m0 C
nodes_state[nei_arr[ran_infe_arr[t]]-1,2]=temp_time+1#记录感染时间
`, U1 C, c8 M8 t& p3 s infected_array[len_arr(infected_array,nodes_num)]=nei_arr[ran_infe_arr[t]]#将此次感染节点放入总的感染节点数组中
3 \ Y% x# J9 ^& w# s g.vs[nei_arr[ran_infe_arr[t]]-1]["color"]="pink"#将此次感染的节点集的所有节点颜色置为粉色9 P/ z% Q0 X9 w0 _, n& K
plot(g)#绘制& R$ J0 t# H9 x' u
t=t+1
+ d3 ^8 b2 `; s; C3 ?1 g y. U' o i=i+10 y/ k! d3 c) T; G6 S* D
if temp_time>set_time-1:#当执行感染的次数等于设置的次数结束感染1 x. I' T4 Y+ ]: N
stop=True , G+ ^: b6 L% ~" \3 i
Q4 {& Q2 J# L
/ f- v k8 B) H
视频演示bilibili传送门
# B! l1 `( b2 v" _; o效果图
) K9 n; g. |, T4 X
y8 j# Y$ n! A, i; ?( p. H8 o* T! f* Z( e1 y
7 m: U& K3 c; I8 V- c0 E/ d9 A4 _! ?3 |. u9 M% `' L
: M% I: V2 ]6 M5 @: v6 {* N& n- J& {
3 N* M& c9 _ v; \
* Y4 f! i- C$ C; X* F; @
) g7 Q9 _* o5 _8 Q! B {
8 X2 S6 F! ^! _# {: M+ K0 x b2 l
& F2 h; c, y3 l, \" x2 e
. C) g9 t+ o2 }4 T* |
————————————————! m, J- v' _8 a) y* q3 [& x
版权声明:本文为CSDN博主「eck_燃」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
' i- k! B+ i( t) |( ^3 A原文链接:https://blog.csdn.net/wdays83892469/article/details/80878862
, z8 r. d$ p2 U f
/ z/ v r5 Y7 O6 m1 p# q" H% Y$ P4 U+ r; o W+ x6 J
|
zan
|