| 传染病的微分方程模型 | 您所在的位置：网站首页 › 疾病的流行曲线 › 传染病的微分方程模型 |
传染病的微分方程模型
诸如霍乱、天花、疟疾等曾经肆虐全球的传染性疾病已经得到有效的控制.一些新的、不断变异着的传染病毒却悄悄向人类袭来20世纪80年代，十分险恶的艾滋病毒开始肆虐全球，至今仍在蔓延.2003年春来历不明的SARS病毒突袭人间，给人们的生命财产带来极大的危害.截至2019 年8月21日，刚果（金）已累计报告埃博拉（Ebola）确诊和可能病例2934例，其中死亡病例1965例. 从3000万到零，中国创造了无疟疾的未来建立传染病的数学模型，描述传染病的传播过程，分析受感染人数的变化规律，探索制止传染病蔓延的手段一直是各国有关专家和官员关注的课题.不同类型传染病的传播过程有其各自不同的特点，弄清这些特点需要相当多的病理知识，这里不可能从医学的角度一一分析各种传染病的传播，而只是按照一般的传播机理建立几种模型.模型1——自由增长模型模型假设设时刻t的病人人数x(t)是连续、可微函数.每天每个病人有效接触，即足以使人致病的接触，人数为常数\lambda>0.模型构建考察t到t+\Delta t病人人数的增加，有x(t+\Delta t)-x(t)=\lambda x(t)\Delta t. \\ 再设t=0时有x_0个病人.得到微分方程\frac{dx}{dt}=\lambda x,~~x(0)=x_0~~~~~~~~(1). \\ 模型求解方程(1)的解为x(t)=x_0e^{\lambda t}~~~~~(2). \\ 说明：随着t的增加，病人人数x(t)无限增长.评价不符合实际情况.原因：在病人有效接触的人群中，有健康人也有病人，而其中只有健康人才可以被传染为病人，所以在改进的模型中必须区别这两种人.模型2——SI模型模型假设在疾病传播期内所考察地区的总人数不变，既不考虑生死，也不考虑迁移.人群分为易感染者（Susceptible）和已感染者（Infective）两类（取两个词的首字母，称之为SI模型），以下简称健康者和病人.每个病人每天有效接触的平均人数是常数，称为日接触率.当病人与健康者有效接触时，使健康者受感染变为病人.基本符号N\sim总人数.s(t)\sim时刻t易感染者在总人数中所占的比例.i(t)\sim时刻t已感染者在总人数中所占的比例.\lambda\sim日接触率.模型构建根据假设，每个病人每天可使\lambda s(t)个健康者变为病人，因为病人数为Ni(t)，所以每天共有(Ni(t))\cdot(\lambda s(t))=\lambda Ns(t)i(t) \\ 个健康者被感染. 故\lambda Nsi就是病人数Ni的增加率，有N\frac{di}{dt}=\lambda Nsi~~~~~~~(3). \\ 又因为s(t)+i(t)=1~~~~~~~~~~(4). \\ 再记初始时刻（t=0）病人的比例为i_0.综上得\frac{di}{dt}=\lambda i(1-i),~~i(0)=i_0~~~~~~~(5). \\ 方程(5)是logistic模型. 模型求解分离变量法：\frac{di}{i(1-i)}=\lambda dt.它的解为i(t)=\frac{1}{1+(\frac{1}{i_0}-1)e^{-\lambda t}}~~~~(6). \\ i(t)\sim t和\frac{di}{dt}\sim i的图形（取\lambda=1,i_0=\frac{1}{4}）： SI模型i~t曲线 SI模型di/dt~i曲线当i=1/2时，\frac{di}{dt}达到最大值(\frac{di}{dt})_m，这个时刻为t_m=\frac{1}{\lambda}\ln(\frac{1}{i_0}-1)~~~~(7). \\ 时刻t_m病人增加得最快，可以认为是医院的门诊量最大的一天，预示着传染病高潮的到来，是医疗卫生部门关注的时刻.因为t_m 与\lambda成反比，而日接触率\lambda表示该地区的卫生水平，\lambda越小卫生水平越高. 所以改善保健设施、提高卫生水平可以推迟传染病高潮的到来.当t\rightarrow\infty时，i\rightarrow1，即所有人终将被传染，全变为病人，这不符合实际情况.原因：模型中没有考虑到病人可以治愈，人群中的健康者只能变成病人，病人不会再变成健康者.模型3——SIS模型有些传染病如伤风、痢疾等愈后免疫力很低，可以假定无免疫性.于是病人被治愈后变成健康者，健康者还可以被感染再变成病人，所以这个模型称SIS模型.模型假设在疾病传播期内所考察地区的总人数不变，既不考虑生死，也不考虑迁移.人群分为易惑染者和已感染者两类，以下简称健康者和病人.每个病人每天有效接触的平均人数是常数，称为日接触率. 当病人与健康者有效接触时，使健康者受感染变为病人.每天被治愈的病人数占病人总数的比例为常数\mu，称为日治愈率.病人治愈后成为仍可被感染的健康者.\frac{1}{\mu}是这种传染病的平均传染期. 模型构建考虑到假设4，SI模型的（3）式应修正为N\frac{di}{dt}=\lambda Nsi-\mu Ni~~~~~(8). \\ （4）式s(t)+i(t)=1不变.于是（5）式应改为\frac{di}{dt}=\lambda i(1-i)-\mu i,~~i(0)=i_0~~~~~~~~(9). \\ 模型求解伯努利方程：\frac{di}{dt}+(\mu-\lambda) i=-\lambda i^2.令y=\frac{1}{i}，则\frac{dy}{dt}-(\mu-\lambda)y=\lambda,y(0)=\frac{1}{i_0}. \\ 常系数变易法，得y=(\frac{1}{i_0}+\frac{\lambda}{\mu-\lambda})e^{(\mu-\lambda)t}-\frac{\lambda}{\mu-\lambda}. \\ 最后i=\frac{1}{y}.图形分析定义\sigma=\frac{\lambda}{\mu}~~~~(10). \\ 注意到\lambda和\frac{1}{\mu}的含义，可知\sigma是整个传染期内每个病人有效接触的平均人数，称为接触数.利用\sigma，方程（9）可以改写作\frac{di}{dt}=-\lambda i[i-(1-\frac{1}{\sigma})]~~~~~~(11). \\ 分离变量法：\frac{di}{i[i-(1-\frac{1}{\sigma})]}=-\lambda dt，两边求不定积分，得i=\frac{1}{a-(a-\frac{1}{i_0})e^{-\lambda t/a}},a=\frac{1}{1-\frac{1}{\sigma}}. \\ 由方程（11）容易先画出\frac{di}{dt}\sim i的图形.当\sigma>1时，比如取\sigma=2,\lambda=1.当\sigma\leq1时，比如取\sigma=\frac{1}{2},\lambda=1. SIS模型di/dt~i曲线（sigma;amp;amp;amp;amp;gt;1） SIS模型di/dt~i曲线（sigma;amp;amp;amp;amp;lt;=1）再画出i\sim t的图形.当\sigma>1时，比如取\sigma=2,\lambda=1,i_0=\frac{1}{4}.当\sigma\leq1时，比如取\sigma=\frac{1}{2},\lambda=1,i_0=\frac{1}{4}. SIS模型i~t曲线（sigma;amp;amp;amp;amp;gt;1） SIS模型i~t曲线（sigma;amp;amp;amp;amp;lt;=1）接触数\sigma=1是一个阈值.当\sigma>1时，i(t)的增减性取决于i_0 的大小（仅画了增加的情形），但其极限值i(\infty)=1-\frac{1}{\sigma}随\sigma 的增加而增加.当\sigma\leq1时，病人比例i(t)越来越小，最终趋于0.这是由于传染期内经有效接触从而使健康者变成的病人数不超过原来病人数的缘故.模型4——SIR模型大多数传染病如天花、流感、肝炎、麻疹等治愈后均有很强的免疫力，所以病愈的人既非健康者（易感染者），也非病人（已感染者）.他们已经退出传染系统，这种情况比较复杂.模型假设总人数N不变.人群分为健康者、病人和病愈免疫的移出者（Removed）三类，称SIR模型.时刻t三类人在总人数N中占的比例分别记作s(t)，i(t)和r(t).病人的日接触率为\lambda，日治愈率为\mu（与SI模型相同），传染期接触数为\sigma=\frac{\lambda}{\mu}.模型构成由假设1，有s(t)+i(t)+r(t)=1~~~~~(12). \\ 根据条件2，方程（8）仍成立，即N\frac{di}{dt}=\lambda Nsi-\mu Ni. \\ 对于病愈免疫的移出者而言，有N\frac{dr}{dt}=\mu Ni~~~~~~~~~~~~(13). \\ 记初始时刻的健康者和病人的比例分别是s_0>0和i_0>0.由(8)，(12)，(13)式，SIR模型的方程为\left\{ \begin{array}{ll} \frac{di}{dt}=\lambda si-\mu i, & \hbox{$i(0)=i_0$,} \\ \frac{ds}{dt}=-\lambda si, & \hbox{$s(0)=s_0$.}~~~~~~~~(14) \end{array} \right. \\ 模型求解从以上方程组(14)无法求出s(t)和i(t)的解析解（用公式显示表示），先作数值计算.在方程组(14)中设\lambda=1,\mu=0.3,i(0)=0.02,s(0)=0.98.用MATLAB软件编程，新建函数：function y=ill(t,x) a=1;b=0.3; y=[ax(1)x(2)-bx(1),-ax(1)*x(2)]'; end 绘图：ts=0:50; x0=[0.02,0.98]; [t,x]=ode45('i11',ts,x0);[t,x] plot(t, x(:, 1),t, x(:, 2)),grid,pause plot(x(:,2),x(:,1)),grid   输出的简明计算结果列入： i(t),s(t)的数值计算结果 i(t)（蓝色），s(t)（黄色）图形病人比例i(t)由初值增长至约t=7时达到最大值，然后减少，t\rightarrow\infty,i\rightarrow0. \\ 健康者比例s(t)则单调减少，t\rightarrow\infty,s\rightarrow0.0398. \\  i~s图形（横坐标为s）i\sim s的图形称为相轨线.初值i(0)=0.02,s(0)=0.98相当于图中最右下角的点(0.98,0.02).随着t的增加，(s,i)沿轨线自右向左运动.设i=i(s)，其中s=s(t). 由复合函数求导公式\frac{d i}{dt}=\frac{di}{ds}\cdot\frac{ds}{dt}，并将方程（14）代入，得\frac{di}{ds}=\frac{\mu}{\lambda}s-1,i(s_0)=i_0. \\ 解之得i=\frac{\mu}{\lambda}\ln s-s+c, \\ 其中常数c=i_0+s_0-\frac{\mu}{\lambda}\ln s_0. 因为\lambda=1,\mu=0.3,i(0)=0.02,s(0)=0.98，故c=1.0061，且i=0.3\ln s-s+1.0061. \\ 相轨线分析 SIR模型的相轨线s\sim i平面称相平面，相轨线在相平面上的定义域D=\{(s,i)|s\geq0,i\geq0,s+i\leq 1\}~~~~~（15）. \\ 在方程(14)中消去dt，用\sigma=\frac{\lambda}{\mu}，得\frac{di}{ds}=\frac{1}{\sigma s}-1,~~i|_{s=s_0}=i_0~~~（16）. \\ 解之得i=(s_0+i_0)-s+\frac{1}{\sigma}\ln\frac{s}{s_0}~~~（17）. \\ 在定义域D内，(17)式表示的曲线即为相轨线，其中箭头表示了随着时间t的增加s(t)和i(t) 的变化趋向.记s(t)，i(t)和r(t)在t\rightarrow\infty时的极限值分别记作s_\infty，i_\infty和r_\infty.不论初始条件i_0,s_0如何，病人终将消失，即i_\infty=0~~~~(18). \\ 从图形上看，不论相轨线从P_1或从P_2点出发，当t充分大时，它都会与s轴相交.2. 最终未被感染的健康者的比例是s_\infty. 在(17)式中令i=0 得到，s_\infty是方程s_0+i_0-s_\infty+\frac{1}{\sigma}\ln\frac{s_\infty}{s_0}=0~~~~(19). \\ 在(0,\frac{1}{\sigma})内的根. 在图形上，s_\infty是相轨线与s轴在(0,\frac{1}{\sigma}) 内交点的横坐标.3. 若s_0>\frac{1}{\sigma}，则i(t)先增加；当s=\frac{1}{\sigma}时，i(t)达到最大值i_m=s_0+i_0-\frac{1}{\sigma}(1+\ln \sigma s_0)~~~~(20). \\ 然后当s>\frac{1}{\sigma}时，i(t)减小且趋于0.s(t)则一直单调减小至s_\infty.如图中由P_1(s_0,i_0)出发的轨线.4. 若s_0\leq\frac{1}{\sigma}，则i(t) 单调减小至零，s(t)单调减小至s_\infty.如图中由P_2(s_0,i_0)出发的轨线.如果仅当病人比例i(t)有一段增长的时期才认为传染病在蔓延，那么\frac{1}{\sigma}是一个阈值.当s_0>\frac{1}{\sigma}（即\sigma>\frac{1}{s_0}）时传染病就会蔓延.而减小传染期接触数\sigma，即提高阈值\frac{1}{\sigma}，使得s_0\leq\frac{1}{\sigma}（即\sigma\leq \frac{1}{s_0}），传染病就不会蔓延（健康者比例的初始值s_0是一定的，通常可认为s_0 接近1）.即使s_0>\frac{1}{\sigma}，从(19)，(20)式可以看出，\sigma减小时，s_\infty 增加，i_m降低，也控制了蔓延的程度.在\sigma=\frac{\lambda}{\mu}中，人们的卫生水平越高，日接触率\lambda越小；医疗水平越高，日治愈率\mu越大，最终\sigma越小，所以提高卫生水平和医疗水平有助于控制传染病的蔓延.另一方面\sigma s=\lambda s\cdot \frac{1}{\mu}是传染期内一个病人传染的健康者的平均数，称为交换数，其含义是一个病人被\sigma s个健康者交换.所以当s_0\leq \frac{1}{\sigma}，即\sigma s_0\leq 1时，必有\sigma s\leq 1. 既然交换数不超过1, 病人比例i(t)绝不会增加，传染病不会蔓延.群体免疫和预防根据对SIR模型的分析，当s_0\leq\frac{1}{\sigma}时传染病不会蔓延.所以为制止蔓延，除了提高卫生和医疗水平，使阈值\frac{1}{\sigma} 变大以外，另一个途径是降低s_0，这可以通过比如预防接种使群体免疫的办法做到.忽略病人比例的初始值i_0，设移出者比例的初始值r_0，则有s_0=1-r_0，于是传染病不会蔓延的条件s_0\leq \frac{1}{\sigma}可以表为r_0\geq1-\frac{1}{\sigma}~~~~(21). \\ 这就是说只要通过群体免疫使初始时刻的移出者比例（即免疫者比例）r_0满足(21)式，就可以制止传染病的蔓延.以上内容参考自：姜启源，谢金星，叶俊，《数学模型》.分享以上正文课件的LaTeX原代码： 链接：https://pan.baidu.com/s/1tzH95_eZYLEn1lqyuqTQ4Q 提取码：eq74 下载后解压，用LaTeX相关编辑器双击lec5correct.tex打开，运行后，可得同名pdf文档，内容如下： 链接：https://pan.baidu.com/s/1Ie4htVDofKfHoZT8TnR9hQ 提取码：bua1 |
【本文地址】