数学建模社区-数学中国
标题:
数学建模中的传染病模型及其编程求解
[打印本页]
作者:
浅夏110
时间:
2020-5-17 09:59
标题:
数学建模中的传染病模型及其编程求解
问题的提出
# y( r( L7 M9 b7 j
医生们发现，在一个民族或地区，当某种传染病流传时，波及到的总人数大体上保持为一个常数。即既非所有人都会得病也非毫无规律，两次流行（同种疾病）的波及人数不会相差太大。如何解释这一现象呢？试用建模方法来加以证明。
2 @" \/ x2 b1 s: ~& A2 o7 Y
6 l3 z2 v& X" X( ]. c/ v- o
指数模型
# l* E9 ~8 f: w# r- }1 b
定义已感染人数为i(t) i(t)i(t)，假设每个病人单位时间有效接触（足以使人致病）的人数为λ \lambdaλ，那么，在时间段Δt \Delta tΔt内，病人的增量可以用如下的公式进行计算
- @+ h" x9 A: H$ _, g. S
i(t+Δt)−i(t)=λi(t)Δt i(t+\Delta t) - i(t) = \lambda i(t)\Delta t
$ x$ P! O% n% Q6 Q2 d8 {9 F+ q
i(t+Δt)−i(t)=λi(t)Δt
6 ]5 {7 |- s( b7 r7 A
9 Y9 M% S; g" f- ^! m& e9 [$ t
将i(t) i(t)i(t)移到等式的右边，我们得到如下的递推公式
; g1 t3 u0 h) I$ T% n3 h( {
i(t+Δt)=i(t)+λi(t)Δt i(t+\Delta t) = i(t) + \lambda i(t)\Delta t
! E. `/ M# G& k2 o0 N
i(t+Δt)=i(t)+λi(t)Δt
9 ~: n4 C- t6 t& C7 o
& Y) @ m$ z+ A# T e) c. x; m
以上递推公式意味着，我们可以通过当前时刻的病人人数和致病参数λ \lambdaλ，计算得到Δt \Delta tΔt时间后的病人人数，将以上思想在Python中进行实现，代码如下。
1 f# b+ j0 {# S( N: i
4 \. ~1 k, T8 |, j& y
import matplotlib.pyplot as plt
& y) ]+ g/ E' T$ P3 c- \
%matplotlib inline
. _4 g$ {7 V3 C4 H
deltaT = 0.01
; t7 o8 X5 A0 ~3 O# b
lamb = 2
5 p1 I5 Z; }( S, g7 c5 ]
i_list = []
# H$ h/ C( s& u3 @, }6 m* d% o
i0 = 0.08; # 初始有8%的人患病
# i7 d8 D: V' l5 e1 c6 b
i_list.append(i0)
, H& v1 E$ N6 r/ [7 p
Tot_Time = 10
+ h8 H4 V* u8 W2 g3 r
TotStep = int(Tot_Time/deltaT)
& w4 r- ~/ f( V
##
) @% x: I F, G' Q- Q
for i in range(TotStep):
& V _4 r+ b: G- U- ?2 ] ~
i_new = i_list[-1] + lamb * i_list[-1] * deltaT
. c$ @, {! ?# L2 N2 a* y: B3 L6 C
i_list.append(i_new)
7 p* Y# [, r/ F+ J8 I: ]7 j
plt.plot(i_list)
$ g$ }3 E/ V0 N! I2 \3 d9 K
9 ^. j. m) W, b/ A& l) V+ a
将以上代码在Jupyter Notebook中运行，得到病人人数的变化趋势见下图，从中我们可以看到病人的增长是指数级的，在短短十天后，已经有3000万人患病！这显然不符合实际情况的，那么问题出在哪里了呢？
! g$ a/ K9 A* \9 a
$ [6 ?3 _# f+ r; K
* v& C* D7 a( j
实际上，若病人解除的是病人，并不能够使病人再次患病，实际上以上的算法导致了重复计数现象的发生。解决办法：必须区分已感染者和未感染者。
* ^ Q( R0 \- @1 U0 B5 r
& u0 ^$ \0 g+ q
SI模型
; ]1 A# }" C% K
现在我们将人群分成两个群体：已感染者（病人，Infected）和未感染者（健康者，Suspect），该模型称为SI模型，模型假设：
, L2 w# m/ j* H2 t7 I/ K7 N
6 N2 ^' l W" w/ l' B K# x
在研究时间内，不考虑死亡率和出生率，即总人数N NN不变，病人和健康人的比例分别为i(t) i(t)i(t)和s(t) s(t)s(t)
/ R q ]( m x" M- s% L
每个病人在单位时间内有效接触并致病的人数为λ \lambdaλ，且只有接触健康人才会致病，称λ \lambdaλ为日接触率
y% ]6 n$ k% \/ H3 y/ K9 ^- W
仿照指数模型里面的建模方法，在时间段Δt \Delta tΔt内，病人的增量可以用如下的公式进行计算
4 q- Z% s- d% g6 t7 ^6 ^
+ T+ Z& x+ S7 ~' y$ Z* d' Y
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
3 S5 U9 w0 O- t
N[i(t+Δt)−i(t)]=[λs(t)]Ni(t)Δt
7 {% E4 {/ J T0 x# e( {
2 _7 E" o! a! \; M
消去N NN，再将i(t) i(t)i(t)移到等式的右边，我们得到如下的递推公式
/ E+ u% I! ^3 m; f0 o% Q
i(t+Δt)=i(t)+λi(t)s(t)Δt i(t+\Delta t) = i(t) + \lambda i(t)s(t)\Delta t
7 Q [ Q4 i+ V9 E7 D7 ]! g) o; B9 M
i(t+Δt)=i(t)+λi(t)s(t)Δt
) m' D% H1 f3 F0 w
* P) l# {9 u: c
同样地，我们可以通过当前时刻的病人人数和致病参数λ \lambdaλ，计算得到Δt \Delta tΔt时间后的病人人数，将以上思想在Python中进行实现，代码如下：
- D. K# c9 S. }' X
6 Y E( c( G' e: h
import matplotlib.pyplot as plt
) I- E1 x: z$ y3 r8 O
%matplotlib inline
6 W" c8 Z1 R8 o0 h
deltaT = 0.01
) }* }' S1 v$ N# p! P& H9 Q9 E
lamb = 2
* o3 {3 k8 W1 x4 ]" F2 z2 N, c- M
i_list = []
" X1 H! n* q) P5 \
s_list = []
0 n' @7 v4 ~0 B
i0 = 0.08; # 初始有8%的人患病
/ i. I% A: c; x7 U3 S3 @ G
i_list.append(i0)
1 M- ]! B" c8 l! E2 T8 r0 ?
s_list.append(1 - i0)
/ `( l/ _2 [& D
Tot_Time =5
$ X% T/ U ~$ I: W' a
TotStep = int(Tot_Time/deltaT)
W9 i$ o( ]0 R# @8 e z- d
##
7 u! U# E6 O8 k$ I, d
for i in range(TotStep):
; Z; \) s A, a
i_new = i_list[-1] + lamb * i_list[-1] * deltaT * s_list[-1]
/ Z- @2 r7 ^ I9 `5 N! Y7 u. ?/ k
i_list.append(i_new)
* h/ J h+ z3 X/ \7 E( ?4 R
s_list.append(1- i_new)
% S3 W" B" U, z9 X/ Y, \- k
Time = [i * deltaT for i in range(TotStep + 1)]
* E7 U( f$ y9 ~: w( c, ^" M1 W
plt.plot(Time,i_list)
$ T! ^* j6 _1 ~! B( a4 S" h* ^
plt.plot(Time,s_list)
) Q0 {: u3 i1 d" X
plt.title("SI",fontsize = 20)
9 j: g6 W8 e' V
plt.xlabel("Time")
$ K+ H2 X0 U% g4 g! y
plt.ylabel('i(t)')
2 C( ^) b2 U/ m/ X9 ^+ @+ F: _
9 M3 R- a7 C! Y1 x
从SI模型我们可以看到，病人比例不再会出现"指数爆炸"的情况，在t→∞ t \rightarrow \inftyt→∞时最大患病比例为1。在SI模型中，病人数量的增长曲线是一个典型的S型曲线，又称为Logistic曲线，该曲线在生物学上经常被用来描述物种的增长模。
3 w! X% l% A* {6 ]
1 ]/ v s) F0 |
然后，SI模型的结论告诉我们，无论λ \lambdaλ多么小，最终人群都会患病，这显然也是不符合实际情况的。
. | K6 K- \ H% g5 y) C" L8 C# G+ }
————————————————
9 i) Y. `+ ~( ?+ {( S0 a
版权声明：本文为CSDN博主「任公子ha」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
- R9 U- t, ]$ t k3 e- F) S4 L
原文链接：https://blog.csdn.net/baidu_26746963/article/details/93918383
% q) C* p7 f6 A, |: P% c
' D1 i: E; B. J! _9 D& s9 s7 o
4 X7 u6 X! N) C: C/ h X* |
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5