QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 4921|回复: 0
打印 上一主题 下一主题

数学建模之传染病SIR模型(新冠真实数据)

[复制链接]
字体大小: 正常 放大
杨利霞        

5273

主题

82

听众

17万

积分

  • TA的每日心情
    开心
    2021-8-11 17:59
  • 签到天数: 17 天

    [LV.4]偶尔看看III

    网络挑战赛参赛者

    网络挑战赛参赛者

    自我介绍
    本人女,毕业于内蒙古科技大学,担任文职专业,毕业专业英语。

    群组2018美赛大象算法课程

    群组2018美赛护航培训课程

    群组2019年 数学中国站长建

    群组2019年数据分析师课程

    群组2018年大象老师国赛优

    跳转到指定楼层
    1#
    发表于 2021-6-22 15:35 |只看该作者 |倒序浏览
    |招呼Ta 关注Ta

    ) Q9 V  A7 @7 o0 Y, `% _+ U( H数学建模之传染病SIR模型(新冠真实数据)$ l4 z; v/ V5 A
    传染病模型的基本问题
    ) a2 |/ L2 t9 F描述传染病的传播过程
    ! e. U% _: ]) e0 P分析受感染人数的变化规律
    ; B  Q4 S- J+ E3 x/ j! |( H预报传染病高潮到来的时刻
    * x: q4 O$ f- ]7 b8 a; D; Q. `2 \预防传染病蔓延的手段
    ! n! `& F6 \! x" t: c! W# }0 c# Y5 m按照传播过程的一般规律用机理分析方法建立模型, |# v$ B  X& G  i! z
    注:我们这里是介绍数学医学领域中基本的传染病模型。不从医学角度分析各种传染病的特殊机理,按照传播过程的规律建立微分方程模型.
    6 ^, ]+ g7 S0 D2 I9 G! D1 f% Z( I0 }4 ^4 _% o$ _
    ; ?/ l6 X) \& l: @# a) i
    建立模型* Z1 R% g7 I0 g. O) k) e
    模型一1 }; k8 u. b( ~% R/ u! M
    假设:
    / K0 N- ^  b+ J, Y% C" F- R) {; a' }6 v# w# c; N, J8 [/ a: Z
      m2 `, s) T5 p: B  A' s0 ~+ @
    设已知感染人数为i ( t ) i(t)i(t)(病人数量随时间变化)0 R; P5 z3 x2 e" H: I0 c
    设每个病人(单位时间)每天有效接触(足以使人治病)人数为λ \lambdaλ. {$ B" v. B/ r5 t% u3 L/ I
    模型:
    % m7 R% e4 F1 m2 c% D; _; n单位时间Δ t \Delta{t}Δt内,新 增 的 人 数 ( 现 有 − 原 有 ) = 原 有 的 × λ 新增的人数(现有-原有)=原有的 \times \lambda新增的人数(现有−原有)=原有的×λ,即( X2 K( P. \! I9 i% |

    ( [1 F- _( H9 k) F7 F

      X% b0 m, L) x! Mi ( t + Δ t ) − i ( t ) = λ i ( t ) Δ t i(t+\Delta{t})-i(t)=\lambda i(t)\Delta{t}i(t+Δt)−i(t)=λi(t)Δt  C1 c- N. b: Y  a& m
    一开始的感染人数为i 0 i_0i & X8 j, x% m$ a9 N. s  y* K
    0
    9 Y' ^8 F: g6 u8 E​       
    ' E0 k- z/ _9 O ( E; y6 p% {/ Y' ~
    i ( 0 ) = i 0 i(0)=i_0i(0)=i 5 u2 q) U& c$ X+ f3 t
    0
    * i* A( e5 g  x4 A2 |/ o) d​        " z) _0 S  T& b# t

    8 ^' p* m$ S8 M7 y! F1 U1 q- Q9 ^2 L解微分方程可以得到
    # d  Z" T+ P* B4 r& j+ [i ( t ) = i 0 e λ t i(t)=i_0e^{\lambda t}i(t)=i + z0 Y! t4 `3 q+ w3 N
    0
    " B, j- I3 D# S- r) M​       
    # H% i( M( l& Z( M1 U5 Q2 C2 t e   D8 R+ F: G6 C& Z, Z; N, i" k
    λt: z; c5 E; X" E3 A. a0 B9 N

    ( b- f; b9 A* p0 X# S0 d  K; V- o所以可以可到当λ → ∞ \lambda \rightarrow \infinλ→∞时i ( t ) → ∞ i(t) \rightarrow \infini(t)→∞
    7 _5 S+ w% F; t# y& k当然这是不可能的,因为我们考虑的因素太少了,首先一个是,若有效接触的是病人,则不能使病人数增加,所以必须区分已感染者(病人)和未感染者(健康人)看模型二来解决这个问题1 X0 v$ s9 R% C/ A2 T

    8 f) J6 ?# h& X" i" Y  |+ K% H

    ) ^! E4 p( H- w模型二- I5 T, ]/ l$ m& T6 I
    假设:9 O, C1 Y6 X7 f7 ]2 w  Y

    0 P1 N: r' K% g6 m" z- X" c3 H

    % P: b; l  A# ]4 {+ i! @! ^将人群分为两类:易感染者(Susceptible,健康人)和已感染者(Infective, 病人).
    1 m6 s$ \' d7 k$ h" ]2 n9 K总人数N不变,时刻t健康人和病人所占比例分别为s ( t ) s(t)s(t)和i ( t ) i(t)i(t), 有s ( t ) + i ( t ) = 1 s(t)+i(t)=1s(t)+i(t)=1
    / M; E7 s) m+ W每个病人每天有效接触人数为λ \lambdaλ(日接触率),且使接触的健康人致病." s+ ^/ r  O3 i1 b
    建模:
    ( z/ \* G3 ]8 e' ~1 T每天新增的总人数为原有的人数乘以每个人可以传染的健康的人数,再乘Δ t \Delta tΔt0 ?$ b. B: B& w( ~, h1 q: B+ b
    3 y# Y! }) P# v- O
    & N. D/ B4 _+ ~  }- g! o
    Δ t \Delta tΔt除过去,两遍N约分得到下面,: W1 v- H7 z: I1 @0 f9 ~

    " H. w/ W- ^' r9 I: c$ n
    - U8 r2 }# m( [, j2 _7 |3 B! p
    MATLAB解一下这个微分方程9 d/ E0 O( B' Z$ l7 U
    6 C' A& L' f% ^* m8 o) ^0 e% X
    ; W3 ~: q9 h0 g( W8 Q, k
    y=dsolve('Dy=n*y*(1-y)','t');# m0 V2 p/ u! ?# _1 t

    . [! U) T# c9 h1 `8 c5 @
    3 B6 \( I+ f5 I% ~, n9 W# t2 T, j
    y =# r: X5 ?0 p& i. K/ E  s1 b
    -1/(exp(C1 - n*t) - 1)$ Z) G) C5 H" r( ^
                          0
    . i4 b: q) t9 x                      1  x( t* x- z7 ?2 G; J. j5 g( f
    1
    $ o, ?5 D" @0 f1 g& @5 z; |" }2- `# u2 b( V4 n  R8 b
    3
    ) D! ]5 s8 d3 e# C% q7 l" Q4. F% ~7 h4 N& L" d
    5' A6 Y0 _' H0 V2 z
    6( W3 r% \& \' {* g! f% v
    写规范点就是这个函数
    ! g# N4 V& _0 A
    % l* u: e) A, N" G) ]; J. e( L
    4 [5 j1 V& w) P- T6 |: s5 Z" _
    函数图像大致为( S) v! o+ X+ M) K- |1 @

    4 |9 w$ j8 c" W6 u7 L& s1 P

    # v5 X. ~) ^# k6 R, F, r可以看出t = t m t=t_mt=t 4 X" u  ~! v9 J6 n2 K* i7 n$ u
    m1 P+ \5 L  A6 l  D
    ​       
    , T% Q$ C! h0 o% ]/ _$ m 时这里图像的斜率有个最大值,其也就是传染的最快的时候,即传染病的高潮时刻,当然t m t_mt
    / I# T' T5 Y& v5 Qm
    2 K3 ]3 s7 C) F/ U8 l  V7 n) O​        3 P9 i. V; `- W: W% O; d! R; f3 x
    是可以求出来的
    2 d# k, m/ }) f( A' A/ v" l6 }1 {" Z" ?! e0 o
    0 \, H1 E3 [1 G% E( e
    再看原式,当t → ∞ t\rightarrow \infint→∞时i → 1 i\rightarrow 1i→1% h- [0 ~+ }* k/ H' X' Q9 g. v8 F
    病人的比例为1,当然这也是不可能的,因为我们还没有考虑有没有可能治愈,看模型三2 Y% R5 U5 {: U9 C

    , d! V' d4 w0 b: B5 ^' j6 q% h

    ! w& I' O$ B0 M3 [模型三( Q5 o$ b  ?) T
    假设:$ l& P( i% q6 a+ e) Z
    6 y1 L$ k) r8 H& D
    & b; B$ f4 W8 z8 R8 Q9 D9 O
    传染病无免疫性如伤风、痢疾等——病人治愈成为健康人,健康人可再次被感染。& D$ }0 c% v" ^( t7 W5 W. X0 V
    病人每天治愈的比例为μ \muμ (日治愈率),1 μ \frac{1}{\mu} ( c* _! o- y" u, H1 Z
    μ
    9 Y) r, A  p& t- u3 X$ S' K- v# H' ~3 W1" l5 j4 e: S: W
    ​        9 @  J5 w& x: M8 [0 ~
    为感染期,
    $ R! n' {# @) s0 o4 z3 I模型% m, }$ T' g0 p* `: ]. }
    这是减去了治愈人数之后的新增人数* [. L* c7 V; d
    & }6 S) n( M. h: G- x& u* ]+ Y

    ( d: B/ G; J6 m0 a9 b. V. O7 `
    : _# X; B2 p$ _
    & _2 K. e: Y- f6 L- Y% ?
    σ \sigmaσ 为一个感染期内每个病人的有效接触人数,称为接触数  J3 }+ C: T3 w5 F4 ^

    % @- w, ~/ ~8 U2 L

    ( U! G+ g# H3 B' U1 o可以画出上面的图形分析下
    8 {! k$ [; B. z5 q8 V$ Z: Z
    * r$ m$ s6 _5 l6 W) m3 z; D% j
    . r1 q! z  C. m' ~" ~6 o
    对上面的公式进行分析,可以得到,当i = 1 − 1 σ i=1-\frac{1}{\sigma}i=1−
    # y2 m! s3 @( Z0 ^, x6 t5 jσ
    6 C2 s" l! |2 T* W# G13 i- e! r! X( o" ?0 f/ t9 h
    ​       
    + D) e1 p. I! E. U5 |1 g& F 时,i ii对t的导数为0这也就到了i ii的最大值;当0 < i < 1 − 1 σ 0<i<1-\frac{1}{\sigma}0<i<1−
    8 \( c8 E$ V! qσ
    ( m" d: Z% L# p8 J( C  t11 g3 o0 U/ g- M  y: a; r
    ​       
    ; h2 Y0 M; Q/ U& i# C2 k# o/ X: K 时,d i / d t < 0 di/dt<0di/dt<0,i单调递增,且在d i / d t di/dtdi/dt最大时,i的斜率最大,增速最快;当i > 1 − 1 σ i>1-\frac{1}{\sigma}i>1−
    2 E: Y9 @# L1 w' m/ z# bσ* g9 O$ x9 G/ X0 z1 n4 A5 i1 L( B% @* C
    1& |/ p3 e! N8 r& n7 Z3 S/ F5 X% b
    ​       
    6 B- a. x% a# b2 w% K8 e ,d i / d t < 0 di/dt<0di/dt<0,i是单调递减的。3 E5 s, {  P$ _/ f- h; E

    : {5 h) U8 U  I* }5 R" @
    5 g% V& x3 T1 |' F
    当然我们也可以画出i ii随t的函数图像
    # K, R2 B+ E7 C9 W" T
    8 h! T1 U; E* x2 t+ C

    . O' H7 s+ @6 y: O) S先看红线,若初始条件i 0 > 1 − 1 σ i_0>1-\frac{1}{\sigma}i
    7 @8 Q$ P3 }+ ]) I- m7 Y) F6 s$ p0
    + o4 k! F$ Q- ?4 }/ C​        8 O0 d% L$ w% D
    >1− * g/ f2 E4 k) _8 p' w/ R7 `
    σ
      b& [2 \! ~2 k  i/ ~  T1 e* J* c  p1
    1 b6 h) T% g* g  s2 s​        / U* g5 Y! U' R* P( g" x
    d i / d t < 0 di/dt<0di/dt<0,i就是单调递减的,& L. H1 o2 Z6 n) G( u. z
    若若初始条件i 0 < 1 − 1 σ i_0<1-\frac{1}{\sigma}i ! _0 \; E( @2 C$ s( K+ \% k
    0
    & {3 y3 e) |- F5 s​        & H, j5 |# P. O
    <1− % }( r3 h  g: v8 t* [8 i( {5 z
    σ5 z  |/ J% Z# s9 m/ P! A
    14 w" ]! ?; Z" D; T8 B+ z& @
    ​       
    , g8 u& A* W0 N: s ,i就是递增的,可以看到i对t的导数图像有一个最大值,下面的黑线就有一个增加速率最快的一个值,按S形曲线增长. D4 R' @/ _; M5 W# N- d) k! r* {5 ]

    & y! ^# C* T) P* B4 x
    $ t$ O: T4 G. u( q% O
    σ = < 1 \sigma =<1σ=<1时d i / d t < 0 di/dt<0di/dt<0 i肯定是单调下降的,最终降到0
    7 d$ R! a  Q* M1 [( b% W8 M7 k9 `! S, M' I

    # d3 I/ ], b$ a% |0 E1 ^" n. A+ F% k+ x7 ?* b- p3 Q

    ' }7 n0 A( [& c) ?9 _) t& I; [7 F; v综上:
    : k; ?# F0 F( v6 T. E6 t5 _想让患病者越来越少,σ \sigmaσ必须小于等于1,即感染期内有效接触使健康者感染的人数不超过原有的病人数.
    : p, r0 x1 N6 O* C0 v; Q
    6 v- s  s9 k8 d% k# D  y
    , Q* ~1 I3 v. v+ Q# o! u
    这里我们分析的是感染之后还能感染的情况,但有些病毒感染之后会在体内生成抗体,就不会再被感染了,下面我们分析这种情况。$ ^0 \7 S( N) a/ z- J
    / X/ `) }' `. y0 I5 Z" G# ?) q
    0 k3 q3 j; d. m
    模型四 SIR模型! z- F% i1 s+ g' f/ ~
    SIR模型是常见的一种描述传染病传播的数学模型,其基本假设是将人群分为以下三类:1 c9 [6 x/ H0 H* c, m% q
    ! A% d% a7 ~$ E) q7 W/ E+ e) }- [

    ( \  D1 L2 J( Y5 Z0 ^7 ^. y; j1 易感人群(Susceptible):指未得病者,但缺乏免疫能力,与感病者接触后容易受到感染。
    / m* b  M0 L6 s2 V; j
    ' }4 a  u4 ?  `( y% Y( s

    - V( Z$ I# Q4 R) S4 Q2 z2 感染人群(Infective):指染上传染病的人,他可以传播给易感人群。
    6 T- U5 u1 V) Z3 d+ F; L0 ?5 R, r+ Z% g9 `- v
    9 t1 H' t* _5 {8 X  l' H
    3 移除人群(Removed):被移出系统的人。因病愈(具有免疫力)或死亡的人。这部分人不再参与感染和被感染过程。
    ! Z" `1 H4 P9 R$ T0 l) U" H8 |5 ?0 ]/ f2 ^9 F

      p' s# M; b0 |假设:
    9 V* M5 f: C+ @' @! Q% U5 Z; v6 z4 h5 H& [; P  A; p: h- ]4 G

      {/ t; G# G. `1 q传染病有免疫性如天花、麻疹等——病人治愈后移出感染系统,称移出者(Removed).; y# E# ~, o3 v: |) w  T3 L6 |
    总人数N不变,健康人、病人和移出者的比例分别为s ( t ) , i ( t ) , r ( t ) s(t), i(t), r(t)s(t),i(t),r(t).8 C4 H+ a  g7 h" `' e
    病人的日接触率为λ \lambdaλ , 日治愈率为μ \muμ, 接触数 σ = λ μ \sigma=\frac{\lambda}{\mu}σ= $ t' Y# l" E& ?+ F
    μ+ J8 X  e) ]; n; \
    λ# P5 D/ K) a* f7 }8 j# W
    ​        4 L$ y5 A( p/ i- V

    9 U6 C8 s% @( o  z* [建模:" t" g: n- ^$ A
    s ( t ) + i ( t ) + r ( t ) = 1 s(t)+ i(t)+ r(t)=1s(t)+i(t)+r(t)=15 f8 z+ S9 E, A. G
    这个就是病人减去治愈的人,和上一个模型是一样的, }  q' w( O  a! c* C, q- v

    9 C- H5 e9 R6 D1 C5 R2 i

    ; b7 u5 ^" O6 w- A3 {3 S因为有治愈后是有免疫性的,所以可能被感染的总人数要减少,减去移除者就是
    " L- n, ]/ d# ?. J9 ~' a6 Q- i1 M/ O$ O2 `
    2 P! I3 ~" {1 s" {
    将上式化简为:3 C9 H$ g* o# M/ b( k" k3 V
    9 {6 s# d3 N6 o4 ~9 T$ ?. T9 T
    8 T8 j$ @: K0 U* e3 p5 n7 S9 _! d
    i 0 + s 0 ≈ 1 i_0+s_0\approx 1i : j# D2 i$ O6 J
    0  H$ T0 m5 s0 U; \+ s* D$ I! n$ w
    ​       
    5 l4 A$ e3 S& T6 ^2 S2 S +s
    $ i: D8 n8 k2 u0 I" ?" ~( `0
    9 V0 _7 h2 R1 M% \​       
    " B; S4 d( {! k. @3 m  Y: | ≈1(通常r ( 0 ) = r 0 r(0)=r_0r(0)=r
    2 {0 I6 Q& _/ A9 ]2 n! R+ e( r# [0! B( `, n3 X5 U' w1 Q1 R$ f
    ​       
    ; B8 N! H, M  B$ I  g9 M 很小)* }2 [5 |4 ?: V* f
    # c9 Z1 d& T  [( U1 q) D

    ; d- c8 u! Z+ p+ H7 i) l关于i(t) , s(t) 的非线性微分方程组,没有解析解,只能通过数值计算得到s(t), i(t), r(t)的曲线,下面来看下曲线的数值解的MATLAB程序0 u1 p0 e: f  t' |4 x' A7 a

    , C, U( B) M* p- g! c* {8 H9 q

    ( `+ V2 r& Z  b7 C# T  V5 h这里我们先设λ = 1 , μ = 0.5 , i 0 = 0.01 , s 0 = 0.99 \lambda =1, \mu=0.5, i_0=0.01,s_0=0.99λ=1,μ=0.5,i
    - M" W; {" l  O9 c" S3 c09 @. \! G6 c/ |1 }! h$ t7 l% C7 p
    ​        $ Y9 e# g) ^3 @; P' t! [0 R
    =0.01,s ( r8 v+ i% R* n/ |: [
    0
    + o* E8 ^8 m4 n2 Q​       
    6 o5 A0 g& G9 k- m7 c3 S =0.99
    0 m$ J4 g2 I6 b) B6 N  @. l也就是平均一个病人人传染一个正常人,治愈率为0.5;开始的病人比例为0.01,正常人为0.99,设没有天生带有病毒抗体的人,所以r 0 = 0 r_0=0r
    4 s' W6 a/ l1 K- n+ ?% h0' i* ?* H$ z8 ~  M4 i
    ​        . q1 M7 z% I' m, z
    =0,之后若果病人被治愈,则具有抗体了,有抗体的人为:r = 1 − i − s r=1-i-sr=1−i−s# y: T7 M9 m8 P& h# v$ {
    : m; k" x- @* t- W8 U& B

    ) D' B2 E3 N, W2 its=0:40;
    # T5 r. K: f9 N. |x0=[0.01, 0.99];. ^6 Q/ Z% ^5 [4 N- ]2 j
    [t,x]=ode45('ill',ts,x0);
    2 ?3 F5 d1 \5 l7 Ir=1-x(:,1)-x(:,2);6 M8 n- G2 ~, Z9 u
    plot(t,x(:,1),t,x(:,2),ts,r),grid
    + O4 [* c. p+ K; nlegend('i(t)','s(t)','r(t)')1 I" X# o8 x' @+ I  D) V
    # o4 I4 c  w: N% A* G! _
    5 ]- B4 y4 W* J" `9 W
    function y=ill( t,x). A  a. R% K2 `2 J
    a=1;
    / h: J( e3 X- ]( K1 ib=0.5;6 i) ^. v+ c( [% g. c
    y=[a*x(1)*x(2)-b*x(1);-a*x(1)*x(2)];6 y" E; u0 d- N5 c5 w+ p5 `% b
    1
    $ `! ?8 i; ?3 @2
    4 C! q. l; a6 E/ P) k* _! H, q* m5 q3' f' D+ l9 Z2 R# Y/ A' l3 v/ H
    4% z# i6 F* Y: d* [9 C! J
    57 w7 l' o+ X7 I% S
    6% W$ j) G1 q- h# {$ P. L4 ~
    7
    - n7 q- r! X& V+ Y) Z6 D  `8" k9 y0 F( @( \& L8 o& _$ Z" ~& M
    9
    - `: m) {8 i( m/ C0 Q# U10
    " G& }8 }5 c' ]4 I4 z% r2 ]11# C% W+ I3 v8 f+ j0 e( ^4 H7 v
    0 T0 t3 y3 a9 Q( C+ N
    * r' z$ _; B6 K" U3 t1 g) x
    可以看出:s(t)单调减,r(t)单调增,都趋于稳定, i(t)先增后减趋于0.
    8 v; i6 ^5 `" L2 N4 C' J结果分析5 }, }  k+ N+ E8 k
    先回顾一下参数+ |; G) e' W: v3 ]8 g: s
    接 触 率 λ ; 治 愈 率 μ ; 1 / μ   平 均 传 染 期 ( 病 人 治 愈 所 需 平 均 时 间 ) ; σ = λ / μ   接 触 数 ( 感 染 期 内 每 个 病 人 有 效 接 触 人 数 ) 接触率 \lambda;治愈率 \mu ; 1/ \mu~平均传染期 (病人治愈所需平均时间);\sigma =\lambda/\mu~接触数 (感染期内每个病人有效接触人数)接触率λ;治愈率μ;1/μ 平均传染期(病人治愈所需平均时间);σ=λ/μ 接触数(感染期内每个病人有效接触人数)  g& z  r+ J% g
    可以分析出:
    1 t, P: X( h2 K8 v* E% W9 n# v3 K, p3 E
    . [! v7 O( E. F) Q3 U# T
    随着卫生健康思想水平高,接触率λ \lambdaλ变小: G! [- l% W5 I; L
    随着医疗水平的提高,治愈率μ \muμ增大
    ( T) Q, o2 ~+ E接触数σ = λ / μ \sigma =\lambda/\muσ=λ/μ减小——有助于控制传播.% W5 k& _3 k! K; g
    我们可以试试稍微减少一下λ \lambdaλ,增大μ \muμ,来看下效果# B4 p/ n; L  m) Z
      E6 z" z% Q/ f% \% L, o1 V. @  F
    4 {7 T! B/ j5 e* z! O8 q9 j
    ts=0:40;
    9 x7 ]2 b1 H9 ?" I. h4 ?9 mx0=[0.01, 0.99];  D, T' Y  T8 n5 x
    [t,x]=ode45('ill',ts,x0);- C) p" \; Y5 X% ]% A! l- H
    r=1-x(:,1)-x(:,2);) u8 j2 L, t3 ~6 U
    plot(t,x(:,1),t,x(:,2),ts,r),grid
    " v. x8 U% S! Dlegend('i(t)','s(t)','r(t)')
    , b: e) U' M. i3 M2 T8 R  Y' O1 O
    0 ^* {. g) i5 b9 a* b# Q

    ' i& \! L, F1 _* V4 a1 i$ Jfunction y=ill( t,x), S1 x$ n2 M  }9 M0 @; n, e& v% ^. _
    a=0.8;/ w8 {6 W+ d5 D% {5 r: @  N: s
    b=0.6;
      Q) P/ q- p. z1 J' r6 Xy=[a*x(1)*x(2)-b*x(1);-a*x(1)*x(2)];
    ! K2 \5 O2 J+ _3 d& O4 u+ y1* ~3 ?+ Z2 j3 T) G
    2
    # m1 `: k; M1 w+ `8 E6 q3 P33 k+ Y, n9 _% y/ h8 N2 W
    4
    % ]$ R1 _. h- A) w3 Q6 Q* ]5' o4 `- F; g& f! y+ j$ w
    6/ j" G2 m$ T/ J
    7
    , w8 `" h" I8 c3 j6 X& v8
    ( d5 @) }1 w9 N% b* f9
    7 i# g# w! z$ v" M, w0 Q105 w. N1 o1 m1 [3 W$ @
    11. V) P0 ]6 O2 c1 I/ j
    : @3 R0 ~8 K) `/ ^
    $ d- ]9 F8 C8 \( T
    综上我们可以得出结论:想要减少传染病的传播,我们就要在接触数σ \sigmaσ上下功夫。
    5 h" v  ]6 z3 v6 f$ x8 Q# I9 }' \+ l
    ' U5 [9 A% ?8 C' z
    实战建模
    3 I6 V1 S; Y6 j7 y数据处理
    * w0 _  a3 f  M
    ! M4 L) Y* _; A# i

      u) W3 L, k6 [% y- F8 o5 j首先,我用python爬虫爬取了丁香医生官方数据,一共5534条数据 特征包括感染、死亡、治愈的总数,当日感染、死亡、治愈新增,疑似病例,时间,省份等14个特征
    % i. g9 v5 e: T7 C* f1 s8 s* Y! z& a' g, F' ~9 ?
    & Y! M4 \4 ]1 S

    . t# M9 v4 q4 t# p

    6 ~. ^; S- }% c3 I* q7 z5 u" J  ^( K% d7 _然后用python进行数据提取,提取了较为典型的湖北省的数据作为我的参考依据
    " H+ E8 B+ x& y
    3 M7 o7 J" V% E8 c7 m+ Z- Y

    % i; a2 Y. z+ _! \: U) G* e! T+ [( ]" D5 z

    ! D8 Y1 {# A0 |* h然后用python对数据进行清洗,提取出了患病总数,现存患者总数,死亡总数,治愈总数,时间,省份这几个特征
    + A) N. l( V: Z' f2 m" r9 C$ H: P. |: r# n- Z

    ' B/ z( N* |# x2 D' Q0 m! e对日期格式进行修改,值保留月和日,并与死亡人数的位置交换
    ) ^7 G, G  k5 u6 d
    " k* {" n$ v3 g! y3 u# W
    ; u) v/ t8 z* x8 Z& O8 `
    这里我用python对提取的四个特征分别进行了数据分析(主要包括计算最值,平均值等,),并把1.20日作为第一天,7.02日作为最后一天也就是第165天,做了可视化可视化处理。
    9 A! S, E* U( R5 d- }! v感染人数示意图
    1 L+ [2 V0 B9 S$ C8 M- [$ _
    1 u1 H2 ]( N$ p* Y2 P: l4 i1 ^

    $ ^9 x2 V# E$ i3 t9 W6 s& _- w& f治愈人数示意图) Q0 o; v& D, n5 m  C. C6 M  L
    / A/ d# P& r  O4 H4 F7 u5 U# D

    8 f: f& @" c4 j8 ?* b; n9 x% l# ?# N* D9 L
    2 D* w% B9 T+ K: [# I
    现存患者数量图; z$ }) ?5 u4 ]  m9 J. s

    9 L& H& w2 }  O. `8 d7 E9 W! [
    8 q% p% D" k; C! Z
    死亡人数示意图
    , w1 r; c4 v! [1 U8 Y& A( }  h7 p2 }  R. s/ q& I. M7 S1 q4 S+ S
    + f8 o* Y8 d, p( i  H* Y, A

    - n6 F5 l- s$ v. |% n& m9 b
    7 Q; ]7 ~6 L. l& L+ B; H
    经过上面的图片与describe数据分析,我们发现有一天是异常的,患者多出了平时的十倍左右,经过查阅资料,这天因加强了检测标准,所以增多了很多。为了避免这个数据的影响我们选择将这一天删去(或者用平均数或中位数代替也可)( ~, o0 V. D7 J
    将上面清理过的数据存放到csv文件中
    3 u+ q: _3 I' a2 j+ e
    3 u  ]+ m7 j. c$ Z* I4 F( W

    0 o% E9 x7 @1 k' e- _/ h4 @模型建立% x# a5 Q+ t2 Z$ f% d, J8 s
    模型假设
    1 V; y2 Z2 X% N7 \经过上面数据的分析,我们大体可以进行如下假设:
    5 Y- R& i" U  M1.由于不存在封闭情况,考虑开放体系。- ?! E( G1 \$ S6 Y4 s
    2.目前数据以天为单位发布,因此不考虑连续变化情况,只考虑离散的方程。
    , e3 m. M3 D* }6 {- S3 a3 q( g4 R3.新型冠状病毒的治愈人数和死亡人数相对较 小,因此只考虑 Susceptible(易感)和 Infected(感染) 两类人群。设易感人群总数为N) H- |7 m( w! P8 E
    4.经专家鉴定新冠病毒患者治愈后至少六个月之内不会再被感染,所以设治愈后移出易感人群。
    ) @# @2 t/ y# P$ W/ L5.设每个病人每天有效接触人数为 λ \lambdaλ(日接触率),且使接触的健康人致病.9 \( C  T# Y) x5 y) P
    6.设病人每天治愈的比例为 μ \muμ(日治愈率)
    + t* t1 n" o' k0 J7.时刻t健康人、病人和移出者的数量分别为 s(t), i(t), r(t).
    + E+ A2 X8 c% Y# x  k+ g/ n! g4 \  N# f8 d

    ) O/ d% `& I4 W2 [- P5 [; ]8 ]模型一
    , ^+ A7 q' l* q2 D
    ) P! y$ X& U, u

    * s/ H* ?3 U1 m: L分析可以得到移出者r(t)=治愈人数+死亡人数0 Y- h5 J3 Q; ~; B
    通过python数据处理,我们算出了r(t)的值,并将其可视化
    * ]! E6 i4 e7 ~; q6 i% B2 @
    & v2 i- e0 `6 X4 y/ `
    ) K0 y( S) \" G" |7 e, S3 U

    7 I7 O$ o4 H# t1 _1 A% O
    1 K' v( z' u8 Z# J
    我用MATLAB对其进行了拟合,拟合图像为  Q  M# Q' X0 a3 h8 h

    . A3 U# _2 w+ U. ]" k4 l

    & v. g; j$ ]! v) U* ^- S
    0 a4 A5 I9 m3 }3 _" r4 ~

    9 n+ {# \, Y' i9 P& b' |5 @) \* j# F2 b/ k- H$ R

    ; k$ k/ m; L& i& S# w2 N3 |, M! S3 `分析可以得到患者 i(t)=患病总数-移出者$ j3 X) f8 }. ~. ^, N% S9 O
    可以通过csv文件的currentConfirmedCount 直接获得i(t)数据,当然也可以通过 i(t)=confirmedCountv - r(t)获得,对此我也做了可视化展示
    6 C" N% T7 E0 s6 [* @5 t% o& x' s0 e7 H( A
    9 q. D3 N& r3 }+ ]5 N, W& @
    通过MATLAB程序对其进行拟合,可以得到r(t)的函数图像大致为
    ' B/ h2 }6 R8 V6 I8 p$ v) j2 Y9 E& P( |/ F: P- l2 }
    9 s2 ~( C1 Y, t, }9 q

    * ~3 O) v' P3 l' h# z" n

    0 i- R) V% f4 `2 \% A7 `
    $ t5 l* f- [9 b
    # I4 L8 f! {; e. ?; W2 H' X
    为了方便,利于公式推导,我们先设时刻t健康人、病人和移出者的数量分别为 s(t), i(t), r(t). 所以有
    + a% Y* c# y+ y* w7 m5 C) z; s) k, [& F. y6 b6 R: I# B
    # _* c/ I% w$ r5 U. Q) d- b5 Y" u
    可以推导出每日新增病例的表达式: L7 t6 w# Y& B
    4 u, Y2 _& T# x. W# d9 z

    / o, a1 C: `/ q) k4 z9 A" ^6 g& b2 V" P; E

    / D% i2 X6 p& C( x  I  E1 }
    ) l/ n+ C5 i+ \

    + G% [; R; p+ K9 l* O0 C& f5 p由以上两个公式可以推导出以下两个微分方程; z" F; y- Z/ {( t3 L: ~
    . X: G+ H) }7 m
    4 m7 |& \, m1 J5 H5 I3 N) v

    * ^* X( c  C2 x  D! V0 F
    7 b" |( V/ k/ r8 f
    可以知道初值  f0 U" j4 l* }/ j
    i ( 0 ) = i 0 . s ( 0 ) = s 0 i(0)=i_0.s(0)=s_0i(0)=i
    6 q7 M% E. H$ N6 O) W6 M/ x0
    4 @2 g+ h. Y$ J8 U* ?& ?- c​        0 r4 P5 M1 f2 ~4 t
    .s(0)=s
    6 E- s: r9 P) S. f; A" j; i0
    & R. ~# l9 p' c) W0 X​        9 @. U) A+ K$ j; G9 T; t* a0 b

    ' y0 p" r* S# L2 y! U& K  l因为一开始治愈的和死亡的肯定很少,所以r0可以看为0,于是就有:5 D0 C. c* r/ e+ J; V
    i 0 + s 0 = 1 i_0+s_0=1i - _# {# \1 f' P. ?: ]0 \
    0
    5 G9 z8 a$ D7 s  {8 g​        ; F& E' F2 U: ?2 e& S0 L; r
    +s
    3 l! A3 i* Q3 r& L; g: Q" y$ T0- l3 d8 r9 z# X9 H
    ​        9 q- \% Q6 W; g. e
    =1
    ( e$ |- b0 N3 x% x4 Y  d; s) O通过解以上微分方程我们可以根据经验假设λ \lambdaλ (日接触率)和 μ \muμ(日治愈率)的值分别为1和0.5(也就是每个患者可能使1个正常人患病,患者可能有0.5的概率被治愈);由于一开始患者肯定比正常人少很多,所以我们设i0=0.01,s0=0.99。对其求解可以得到s(t), i(t), r(t),的变化图像
    1 ?  W! ~4 d3 T7 o9 P: s. d; [1 N' f3 q

    ' {6 e1 n- O7 J2 [1 I/ K
    $ b' d8 p2 v7 k& c% W9 J4 w

    + V$ H; b  h1 m( h, hMATLAB程序如下* ~% Y' b( ]) K, a0 }6 x
    ts=0:40;
    0 ~. u" Q; {& D5 L) v6 m- @) Ax0=[0.01, 0.99];$ Y& Y, W/ G* t9 J, C' X
    [t,x]=ode45(‘ill’,ts,x0);. D( f: l$ T. g* |" i) W
    r=1-x(:,1)-x(:,2);
    / [) c: I! i% X! w! P5 iplot(t,x(:,1),t,x(:,2),ts,r,ts,x(:,1)/x(:,2))" s6 v$ [( H/ I; O/ V% }
    legend(‘i(t)’,‘s(t)’,‘r(t)’)
    . Q5 O1 N2 b) x2 M6 W2 \) M. Q* T8 k
    7 f# k, I/ P3 u
    % C) y" i6 j3 j3 V, u
    function y=ill( t,x)
    2 q! n2 H# U8 w3 ha=1;
    7 y6 N0 ?6 v3 I% ]- u/ Gb=0.5;
    " u# x* X2 ?2 U% l4 W6 Ey=[ax(1)x(2)-bx(1);-ax(1)*x(2)];
    7 B. v% w% A$ t  g) J& d, J1 i  l5 {7 `
    2 ^/ n& x( R* _/ A& P6 y1 b1 J7 V) g( E
    结果分析:患病人数肯定有个高潮,但之后高潮就会减弱,并逐步降低为0。随着医疗卫生条件的不断提升,患者的 λ \lambdaλ(日接触率)肯定降低,μ \muμ (日治愈率)肯定上升,所以我们可以把λ \lambdaλ调一点为0.8,μ \muμ调高一点为0.6,可以得到以下趋势图。所以应对传染病很关键的一点是我们要提高医疗卫生条件
    - n: N8 d6 J) G: c5 i0 R8 ~: e, A9 e# [4 l+ X- _, }1 L- L" {
    . S4 y7 B6 v( ^* Z3 a  A6 X9 l
      ]) F: Q* N* m% J& \
    模型二
    + |+ x$ X+ v6 @2 E# k. V8 I( Y9 b. f
    & y( [4 U& e- Y* f, M/ k1 A3 L- z

    % E! `' Z$ V, b( ]- N实际上,λ \lambdaλ (日接触率)和 μ \muμ(日治愈率)都是随着时间变化的,这里我们设s(t), i(t), r(t) 为第t天健康人、病人、移除者(病愈与死亡之和)的数量, s(t)+ i(t)+r(t)=N..) m7 Q" j6 R9 L2 d+ m
    (t), (t) ~第t天感染率, 移除率(治愈率与死亡率之和)
    1 l/ w8 A+ l) b有 d i / d t = λ ( t ) s ( t ) i ( t ) − ( t ) i ( t ) di/dt=\lambda (t)s(t)i(t) - (t)i(t)di/dt=λ(t)s(t)i(t)−(t)i(t)
    5 c% o/ K  d( e& [$ `4 f因为s远大于i, r,s(t)视为常数,所以有* `0 B' x) L9 g7 z" Y

    9 t9 B  D" Z) y, b5 P- ]
    & [8 H+ U+ s8 E  ]+ l
    : e4 S' E* k2 Y# F. t

    4 _- W, i8 P! U6 y% n: O2 \5 x  c取差分近似导数: ^3 p% H- }( w) ]8 p* A. p

    " q4 |5 u* C1 q/ M* m/ y7 l

    . M+ x  m" J1 |0 X7 X
      L! d" L* c! {/ Y3 H2 g

    6 y+ R9 a1 Z" I* J1 M6 M: h% m我们可以先用真实数据对(t)进行展示并进行拟合( a8 U4 R$ ]6 u

    , \- r6 m( Z( ~  L5 g+ \

    1 d; J6 T8 t9 S/ C8 S
    : a# H1 J2 p' H" k4 v# H0 k, H. A

    9 X  r" D: F0 W9 \4 Z. J当然同样的方法对(t)进行拟合1 W9 S7 Z5 m( l; C6 A+ I% F6 m7 Q: s
    , m, Q2 U+ K0 ~# I! e2 P
    ( n5 k  ^' x' t' g5 w( e) b" X
    做不出来了,好难,光这些东西就弄了四天,到了数学建模国赛得多难多累啊,哎,让我这个小白手足无措。毕竟还没有正规的培训,这个模型等期末考完试一定好好做做!!!7 g6 e  y2 e, c& B* q* G" p) W7 w
    冲国奖
    ( `0 I/ A) J4 r( Z4 n# D* @冲国奖1 k8 Q" t' F3 Z. C; B
    冲国奖
    0 R/ q; k9 p4 _8 v————————————————
    - r" M, w! L7 S版权声明:本文为CSDN博主「小白不白嘿嘿嘿」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
    ) ]" G8 P: Q7 |4 f8 S原文链接:https://blog.csdn.net/weixin_45755332/article/details/1070946309 k: ]4 q1 w+ N8 l  c/ i+ n) K

    ' }$ H; c6 X( E' t6 |
    9 R, k% m# u- D- f
    zan
    转播转播0 分享淘帖0 分享分享0 收藏收藏0 支持支持0 反对反对0 微信微信
    您需要登录后才可以回帖 登录 | 注册地址

    qq
    收缩
    • 电话咨询

    • 04714969085
    fastpost

    关于我们| 联系我们| 诚征英才| 对外合作| 产品服务| QQ

    手机版|Archiver| |繁體中文 手机客户端  

    蒙公网安备 15010502000194号

    Powered by Discuz! X2.5   © 2001-2013 数学建模网-数学中国 ( 蒙ICP备14002410号-3 蒙BBS备-0002号 )     论坛法律顾问:王兆丰

    GMT+8, 2025-7-9 22:49 , Processed in 0.478304 second(s), 50 queries .

    回顶部