数学建模社区-数学中国
标题:
数学建模中的传染病模型及其编程求解
[打印本页]
作者:
浅夏110
时间:
2020-5-17 09:59
标题:
数学建模中的传染病模型及其编程求解
问题的提出
' ]8 s, c# o. y5 K
医生们发现,在一个民族或地区,当某种传染病流传时,波及到的总人数大体上保持为一个常数。即既非所有人都会得病也非毫无规律,两次流行(同种疾病)的波及人数不会相差太大。如何解释这一现象呢?试用建模方法来加以证明。
5 F4 c; ]* ^9 c/ e' f1 T; j F
U- R. m9 @( n+ p1 l; I4 d
指数模型
, L' p7 C+ m& c. e: q4 G# N
定义已感染人数为i(t) i(t)i(t),假设每个病人单位时间有效接触(足以使人致病)的人数为λ \lambdaλ,那么,在时间段Δt \Delta tΔt内,病人的增量可以用如下的公式进行计算
6 b6 ~9 |; P2 X
i(t+Δt)−i(t)=λi(t)Δt i(t+\Delta t) - i(t) = \lambda i(t)\Delta t
5 \5 h& _% [+ r! {
i(t+Δt)−i(t)=λi(t)Δt
# B: g* k8 S' n2 T9 h
8 C$ v: R9 _) j, q
将i(t) i(t)i(t)移到等式的右边,我们得到如下的递推公式
6 u. {/ k, x! H. @& n
i(t+Δt)=i(t)+λi(t)Δt i(t+\Delta t) = i(t) + \lambda i(t)\Delta t
& H. `, l' z& `2 k8 Y5 W2 L+ m
i(t+Δt)=i(t)+λi(t)Δt
% a: {$ a6 j O- y' \' T
9 r1 J- n& u1 Z, C7 t$ V5 ?( R
以上递推公式意味着,我们可以通过当前时刻的病人人数和致病参数λ \lambdaλ,计算得到Δt \Delta tΔt时间后的病人人数,将以上思想在Python中进行实现,代码如下。
3 A- |" Y1 I5 z4 Q" x* N
6 Z1 V! ~3 u1 {
import matplotlib.pyplot as plt
9 J1 I# Z4 N% c
%matplotlib inline
' ]/ Z( R' y/ P* ^
deltaT = 0.01
" |' G5 M4 M0 E; q8 _# ]! Q
lamb = 2
4 C0 i* ?( f( l& d) {. R
i_list = []
; U$ B; F& B3 u" m* X
i0 = 0.08; # 初始有8%的人患病
- z' x0 ?7 l! z# b. X
i_list.append(i0)
1 L2 t" U ]7 N Z
Tot_Time = 10
+ b8 [4 n! v6 s4 Y5 [
TotStep = int(Tot_Time/deltaT)
2 U5 f( P* h% r+ I$ L+ o( t
##
. w" x2 C% }0 A& Y' j4 x
for i in range(TotStep):
! z' s4 b! \! e# ^ T' g
i_new = i_list[-1] + lamb * i_list[-1] * deltaT
6 m9 P; M3 |/ U0 a
i_list.append(i_new)
- W- I) m' v! X
plt.plot(i_list)
+ A4 U* y! ]& G; q% Q! e
5 r, v8 ]) c' N. g
将以上代码在Jupyter Notebook中运行,得到病人人数的变化趋势见下图,从中我们可以看到病人的增长是指数级的,在短短十天后,已经有3000万人患病!这显然不符合实际情况的,那么问题出在哪里了呢?
( j. Y9 }; z: b, f% K
- \" X" e& P R) I* [* G5 D
2 D; { r. l6 |3 N5 A
实际上,若病人解除的是病人,并不能够使病人再次患病,实际上以上的算法导致了重复计数现象的发生。解决办法:必须区分已感染者和未感染者。
' G, t6 d: W9 O1 \; Y
& }3 D1 u6 m/ q% h
SI模型
: ~4 u2 l. |" N) d8 X0 W- |
现在我们将人群分成两个群体:已感染者(病人,Infected)和未感染者(健康者,Suspect),该模型称为SI模型,模型假设:
4 K1 Q% Y2 `" R b' {3 f/ L
- N, D7 ~% }+ N8 e4 S
在研究时间内,不考虑死亡率和出生率,即总人数N NN不变,病人和健康人的比例分别为i(t) i(t)i(t)和s(t) s(t)s(t)
& i7 @8 t' m) G6 R! Q+ ?5 B( v6 E
每个病人在单位时间内有效接触并致病的人数为λ \lambdaλ,且只有接触健康人才会致病,称λ \lambdaλ为日接触率
, V: h+ |8 H( W# ^9 f2 n5 v6 H
仿照指数模型里面的建模方法,在时间段Δt \Delta tΔt内,病人的增量可以用如下的公式进行计算
2 a3 P! z; {, E) J% ?
& @( P& j# L K5 Q9 E$ v' I8 i6 s$ x
N[i(t+Δt)−i(t)]=[λs(t)]Ni(t)Δt N[i(t+\Delta t)-i(t)]=[\lambda s(t)] N i(t) \Delta t
, D: F) _4 A/ y" R6 V8 q
N[i(t+Δt)−i(t)]=[λs(t)]Ni(t)Δt
# ]$ H/ d% v3 Z3 I
8 X' \& U- Y4 A4 F; b1 l
消去N NN,再将i(t) i(t)i(t)移到等式的右边,我们得到如下的递推公式
* L) F A' N; u5 W
i(t+Δt)=i(t)+λi(t)s(t)Δt i(t+\Delta t) = i(t) + \lambda i(t)s(t)\Delta t
" `# p5 J S& F; n0 u
i(t+Δt)=i(t)+λi(t)s(t)Δt
* ^% I+ R# \& R8 L' M
6 \) ]% ]' S6 i3 W. ^8 U$ l
同样地,我们可以通过当前时刻的病人人数和致病参数λ \lambdaλ,计算得到Δt \Delta tΔt时间后的病人人数,将以上思想在Python中进行实现,代码如下:
8 X, [" J6 W# j2 E5 x
8 r, k" k- \* ?
import matplotlib.pyplot as plt
' p" [- Y: M1 L0 S# r& S7 o
%matplotlib inline
5 M3 B+ b7 z9 D1 o9 g' j/ f: T
deltaT = 0.01
0 a$ ~7 ?: D2 D! |3 ~
lamb = 2
0 j4 ~2 N5 W% N0 G
i_list = []
0 K4 n6 ]/ x$ p9 {) K n
s_list = []
0 X: v+ Q0 T) H4 ^
i0 = 0.08; # 初始有8%的人患病
3 u# P8 k5 K/ e. r' C3 m
i_list.append(i0)
# r. i6 ~) z: u5 f2 @9 h, D
s_list.append(1 - i0)
. x5 U1 H6 f; s- d4 z
Tot_Time =5
8 w8 |0 U0 |- W0 U4 o# C
TotStep = int(Tot_Time/deltaT)
- W$ Y A5 ~* A+ X# k
##
( J- U( x( Q, N$ j/ h! T
for i in range(TotStep):
; x0 E& O/ O4 B/ U# b
i_new = i_list[-1] + lamb * i_list[-1] * deltaT * s_list[-1]
5 ^) ^' S6 u) u& A8 z& h, ]+ ~
i_list.append(i_new)
0 z7 I5 m7 { Q" q5 E# ? J$ T
s_list.append(1- i_new)
8 Z9 W5 p7 F6 F, @1 {- a2 X3 t) ]( V
Time = [i * deltaT for i in range(TotStep + 1)]
: o+ E; ?6 X* E" u) F
plt.plot(Time,i_list)
8 i! {5 G# i: s# e* i/ t
plt.plot(Time,s_list)
- ~) z# F3 o. {+ h
plt.title("SI",fontsize = 20)
1 k- D4 o( y4 S3 C
plt.xlabel("Time")
& z. O/ @- `, B0 b2 O* h
plt.ylabel('i(t)')
( Z% u. ?: w" t' ?5 I6 K" D
: F0 M+ ]7 ~! l& q4 x
从SI模型我们可以看到,病人比例不再会出现"指数爆炸"的情况,在t→∞ t \rightarrow \inftyt→∞时最大患病比例为1。在SI模型中,病人数量的增长曲线是一个典型的S型曲线,又称为Logistic曲线,该曲线在生物学上经常被用来描述物种的增长模。
( a5 P. ?3 @( F5 Q* z
6 {+ }0 Q1 j7 B3 \
然后,SI模型的结论告诉我们,无论λ \lambdaλ多么小,最终人群都会患病,这显然也是不符合实际情况的。
0 R2 K6 b* I& l4 {3 Y; c% M
————————————————
& q4 e H s% h! @. D8 M
版权声明:本文为CSDN博主「任公子ha」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
6 \8 h: }5 {; l! i
原文链接:https://blog.csdn.net/baidu_26746963/article/details/93918383
, g7 E" E" |) I/ m( S% H" A3 O
+ O5 J4 q/ H7 `; p& V( t
[8 r+ M: T8 v2 D: K
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5