传染病模型(微分方程)

t微分方程建模（传染病模型）的求解。

1、模型1：SI 模型。

假设：（1）t 时刻人群分为易感者（占总人数比例的()s t ）和已感染者（占总人数比例的()y t ） （2）每个病人每天有效接触的平均人数是常数λ，λ称为日接触率，当健康者与病人接触时，健康者受感染成为病人。

分析：根据假设，每个患者每天可以使()s t λ个健康者变为病人，因为病人数为()Ny t ，所以每天共有()()Ns t y t λ个健康者变为病人。

即：dyNNsy dtλ=，且()()1s t y t +=，设初始时刻病人比例为b ，则： (1)(0)dyy y dty bλ⎧=-⎪⎨⎪=⎩，用MATLAB 解此微分方程： >> syms a b>> f=dsolve('Dy=a*y*(1-y)','y(0)=b','t') f =1/(1-exp(-a*t)*(-1+b)/b) %11()1111(1)t ty t b e e b bλλ--==--+- 当0.09,0.1b λ==时，分别在坐标系oty 中作出()y t 的图像，坐标系oyy '中作出(1)y y y λ'=-的图像，>> a=0.1;>> b=0.09;>> h=dsolve('Dy=a*y*(1-y)','y(0)=b','t') h =1/(1-exp(-a*t)*(-1+b)/b) >> f=subs(h) f =1/(1+91/9*exp(-1/10*t)) ()y t 的图像>> ezplot(f,[0,60]) >> grid on>> figure (2)>> fplot('0.1*y*(1-y)',[0,1])>> grid on(1)y y y λ'=-的图像模型分析：（1）当12y =时，dydt达到最大值，则此时病人增速最快。

（2）当t →∞时，()1y t →，即所有的人被传染，全部变为病人，这显然是不符合实际的，其原因是没有考虑到病人可以治愈，人群中的健康者只能变为病人，而病人不会变为健康者。

2、模型2：SIS 模型。

假设：（1）t 时刻人群分为易感者（占总人数比例的()s t ）和已感染者（占总人数比例的()y t ） （2）每个病人每天有效接触的平均人数是常数λ，λ称为日接触率，当健康者与病人接触时，健康者受感染成为病人。

（3）病人每天被治愈的占病人总数的比例为μ，称为日治愈率，显然1μ为这种传染病的平均传染期。

则dyNNsy Ny dtλμ=-。

则建立微分方程模型为： (1)(0)dyy y ydty b λμ⎧=--⎪⎨⎪=⎩用MATLAB 解此微分方程：>> h2=dsolve('Dy=a*y*(1-y)-c*y','y(0)=b','t') h2 =(a-c)/(a-exp(-(a-c)*t)*(-a+c+b*a)/b/(a-c)*a+exp(-(a-c)*t)*(-a+c+b*a)/b/(a-c)*c )>> pretty(h2)/ exp(-(a - c) t) (-a + c + b a) a (a - c)/|a - -------------------------------- \ b (a - c)exp(-(a - c) t) (-a + c + b a) c\ + --------------------------------| b (a - c) /化简：()().().()()()a c t a c t a ce a c ba a e a c ba c a b a c b a c ------++-++-+--2()()()().().()a c t a c t b a c ab a c e a c ba e a c ba c -----=---+++-++ 2()()()().()a c tb ac ab a c c a e a c ba ---=-+---++1()2()().()()a c t ab a c c a e a c ba b a c ---⎡⎤-+---++=⎢⎥-⎣⎦1()1()a c t a a e a c b a c ---⎡⎤=+-⎢⎥--⎣⎦即：1()1()()t y t e b λμλλλμλμ---⎡⎤=+-⎢⎥--⎣⎦。

当（1）λμ≠时，1()1()()t y t e b λμλλλμλμ---⎡⎤=+-⎢⎥--⎣⎦；（2）λμ=时，>> clear>> h2=dsolve('Dy=a*y*(1-y)-a*y','y(0)=b','t') h2 =1/(a*t+1/b)即：11()y t t b λ-⎡⎤=+⎢⎥⎣⎦。

定义λσμ=：一个传染期内每个病人有效接触的平均人数。

则：11,(1)()0(1)y σσσ⎧->⎪∞=⎨⎪≤⎩，用MATLAB 作图像：令0.01λ=，0.05μ=，0.7b =（0.21σ=<） >> clear>> a=0.01;b=0.7;c=0.05;>> h2=dsolve('Dy=a*y*(1-y)-c*y','y(0)=b','t'); >> h22=subs(h2)t3/20/(3/10-3/35 exp(-3/20 t))t3/20/(3/10+1/5 exp(-3/20 t))h22 =-1/25/(1/100-47/700*exp(1/25*t)) >> ezplot(h22,[0,120]) >> grid on0.21σ=<的图像令0.3λ=，0.15μ=，0.7b =或0.3b =分别作图（21σ=>） >> a=0.3;b=0.7;c=0.15;>> h2=dsolve('Dy=a*y*(1-y)-c*y','y(0)=b','t'); >> h23=subs(h2) h23 =3/20/(3/10-3/35*exp(-3/20*t)) >> subplot(2,1,1)>> ezplot(h23,[0,25])>> grid on>> b=0.3; >> h24=subs(h2); >> subplot(2,1,2) >> ezplot(h24,[0,25])grid on21σ=>的图像（上面0.7b =，下面0.3b =）模型分析：（1）1λσμ=≤时，病人比例越来越少，最终趋于零，这是因为传染期内经有效接触从而使健康者变为病人数不超过原来病人数的缘故。

（2）1λσμ=>时，病人比例()y t 增减性是由b 来决定，其极限值1()1y σ∞=-随着σ的增加而增加。

3、模型3：SIR 模型。

假设：（1）人群分为健康者，其比例()s t 、病人()i t 、病愈免疫的移出者()r t 。

（2）病人的日接触率为λ，日治愈率为μ，传染期接触数为λμ。

则()()()1s t i t r t ++=， 对于病愈者而言，dry dtμ=， 设初始时刻的健康者和病人的比例为0s 和0y ，则建立微分方程模型为：00(0),(0)dysy y dt ds sy dt y y s s λμλ⎧=-⎪⎪⎪=-⎨⎪==⎪⎪⎩由于此微分方程组的解析解无法求出，则转为相平面sy 上讨论解的性质。

相轨线的定义域(,)s y D ∈应为：{}(,)0,0,1D s y s y s y =≥≥+≤，由方程组消去dt 并将λσμ=得： 0011s s dy ds s y y σ=⎧=-⎪⎨⎪=⎩用matlb 求解：>> dsolve('Dy=1/cma/s-1','y(s0)=y0','s') ans =1/cma*log(s)-s-1/cma*log(s0)+s0+y0>> pretty(ans)log(s) log(s0)------ - s - ------- + s0 + y0 cma cma 即0001()()lnsy s s y s s σ=+-+（相轨线） 定义域内，1σ=时，00(,)y s 分别取(0.3,0.65)，(0.4,0.35)，(0.5,0.45)，(0.7,0.25)s1-s在同一直角坐标系中作出其图像： >> cma=1;y0=0.3;s0=0.65; >> clear>> f=dsolve('Dy=1/cma/s-1','y(s0)=y0','s'); >> cma=1;y0=0.3;s0=0.65; >> f1=subs(f);>> ezplot(f1,[0,1]) >> hold on>> y0=0.4;s0=0.35; >> f2=subs(f);>> ezplot(f2,[0,1]) >> hold on>> y0=0.5;s0=0.45; >> f3=subs(f);>> ezplot(f3,[0,1]) >> hold onSIR 模型的相轨线>> y0=0.7;s0=0.25; >> f4=subs(f);>> ezplot(f4,[0,1]) >> hold on>> ezplot('1-s',[0,1]) >> grid on模型分析：（1）不论初始条件0s ，0y 如何，病人比例越来越少，最终消失。

（2）最终未被感染的健康者的比例是s ∞，在0001()()lnsy s s y s s σ=+-+中。

令()0y s ∞=时，0001()ln0s s y s s σ∞∞+-+=的单根即为s ∞：最终未被感染的健康者的比例。