乔 丹
(中国人寿宁夏分公司,宁夏 银川 750001)
结核病在18、19世纪被称为“白色瘟疫”,因为缺少治疗方法,在20世纪早期造成了大量死亡。根据世卫组织(2016)的报告,结核病的主要病因是结核分枝杆菌,它可以通过空气在人群中传播。Chou和Friedman(2016)介绍说,这种细菌攻击肺部,疾病会伴随大量咳嗽和痰,呼吸变得困难,脸颊发红,身体其余部分呈灰白色。此外,健康人在接触结核病人咳嗽、打喷嚏或吐痰时产生的一些具有传染性的细菌,就会有一定的传染概率。Calmette博士和Guerin博士对1例新生儿成功接种卡介苗疫苗后,常规治疗情况有所改善。该疫苗至今已有80年的历史,已成为世界上最常用的疫苗之一;在国家儿童免疫规划的国家,婴儿免疫率超过80%。然而,卡介苗虽然对儿童的脑膜炎和播散性结核病有保护作用,但对预防原发感染和潜伏性肺部感染的再激活没有帮助。后来,由于1946年抗生素“链霉素”的发现,结核病首次得到有效治疗和治愈,病情持续改善。将这种抗生素与结核病联系起来的医生被授予1952年诺贝尔生理学或医学奖。
与过去几十年相比,人口流动在最近几十年随着交通和人民意识及财富的加速发展而迅速增加。以中国为例。国际劳工组织(ILO,2010)的统计数据显示,中国近年来经历了最广泛的内部迁移。截至2009年底,约有2.298亿农民工离开家乡到城市购买更高水平的家庭生活。此外,旅游资源丰富,每年也造成了巨大的人口迁移。
从常识上看,初步的结论是移民在一定程度上影响了结核病的传播。因此,此篇论文打算在SIR模型的基础上构建合适的模型来探索人口流动的影响程度,尝试在详细分析和重复试验后得出结论。
Kermack和McKendrick(1927)这两位数学家提出了一个简单的模型,考虑恒定的人口规模,只包括三类人,即易感个体(S),感染个体(I)和恢复个体(R)[1]。
第一步,提出必要假设并绘制原理图。
S(t),I(t),R(t)代表时间t时每组个体的人数,人口总数N(t)是三类人群之和,接触率“β”和治愈率“γ”是固定的,假设所有个体相同。
第二步,得出微分方程。
由于“β”是两个独立个体进行有效接触的比率,“βN”是单位时间内人群中接触的个体数,“”是易感人群占总人口的比例。
因此,单位时间内,感V染个体可传播个体公式为:
很容易发现,会有γI(t)的人通过治疗或自然治愈疾病。
如图1所示,微分方程如下:
此时引入一个重要参数R0,表示基本SIR模型中由计算出的疾病的基本增值数。Emilia和Richard(2010)写道,这是导致原感染者进入完全易感人群的继发性感染的平均数量。表明该病R0>1时存活,只有R0<1时才会灭绝。根据文献9的结论,时,传染蔓延,只有提高的值,提高公共医疗水平使治愈率“γ”提升,最终使得,使传染不再蔓延。
因此,将欧拉方法用于初始SIR模型,方程如下:
对于一个现实的解决方案,种群大小一定是大于零的。我们需要S(t+h)≥0,所以-βS(t)I(t)h+S(t)≥0,这代表着S(t)*(1-βI(t)h)≥0。因为S(t)≥0,1-βI(t)h≥0。因为I(t)≤0是基础条件,是假设的必要条件。
根据这两个方程,我们假设β和γ的两个实数,并考虑人群中感染人数较少的情况,用EXCEL模拟单位时间内的变化过程。在此变化过程中,当R0>1时,,I(t)有一个最大值且最终趋向于0。
第四步,考虑其他影响因素,调整基础模型。
因SIR模型主要用于描述传染病发展的一般规律,可操作性强,但同时忽略了很多细节,可能不符合现实情况。为了使模型更加真实,模型中考虑了出生/自然死亡率“μ”、疫苗有效接种率“ρ”和疾病致死率“α”。假设所有的新生儿都属于易感人群,现将原理图修改如下。
如图2所示,微分方程如下:
第五步,引入矩阵,用于计算城市间的人口流动。
在所有城市肺结核传染情况都可以用SIR模型表示的背景下,我们考虑了人口迁移的影响。为简单起见,假设每个2个城市之间的人口迁移数量是固定的,这些城市的β和γ相同,定义一个新的参数“m”来表示城市之间的迁移数量。迁移比例不应超过10%,以满足实际情况。下面演示了三种可能的流动方向,箭头表示人口流动情况。
第六步,将人口流动加入模型并进行离散化。
假设共有n个城市,Ni表示i市的总人数,Si表示i市的易感人口,Ii表示i市的感染人口,Ri表示i市的恢复人口。微分方程如下:
离散化后,模型为:
使用如上方程在EXCEL中进行模拟实验,假设实验中参数为实数,我们定义然后,使用这些属性在EXCEL中拟扩散。在实验中,我们定义传染率“β”为0.1,恢复率“γ”为0.45,出生/自然死亡率“μ”为0.1,疾病致死率“α”为0.00015,疫苗有效接种率“ρ”为0.1,初始人口基数“N”为1000,初始感染者“I0”为1,步长“h”为0.01。
为了检验“m”对肺结核传播的影响,我们构建了以上3种不同的城市流动情况。由于在所有情况下,感染人群最终都将成为一个常数,因此,此实验通过记录那个最大值在不同m下的到达时间来评估两者之间的关系。实验结果表明,在开始时,“m”的微小变化会对结果产生很大的影响,所以我们将“m”的实验值设为1、2、4、6、8、10、15、20、25、30、40、50、60、70、80、90和100。假设0时刻只有1个城市有1名感染者,其他城市的感染情况取决于迁移人口。以到达感染峰值的时间作为纵轴,“m”值作为横轴绘制散点图,观察线性关系。绘制图如下:
如图4所示,在最简单的双城案例中,最大值的到达时间在城市1是固定的,但m与城市2中的T(Max(I(t)))呈对数关系。如图5所示,在四个环线相连的城市中,城市2和城市4的趋势线相同,因为它们都是离城市1最近的城市。如图6所示,在最后一种情况下,只有两条清晰的趋势线和方程被证明,因为城市345只与中心城市2进行活动,所以共享了趋势线。
从散点图的分布,我们推测“m”与首次到达传染峰值的时间呈对数关系,通过EXCEl生成趋势线、方程和判定系数来验证[2]。判定系数R square表示趋势线在0到1范围内的拟合程度,判定系数越高,趋势线的可信程度越高。考虑到所有的判定系数都高于0.95,可以说明“m”与首次到达传染峰值的时间“T(Max(I(t)))”之间存在一种对数关系。m越大,达到最大值的时间就越短,对数方程为T(Max(I(t)))=-a*ln(m)+b。
本文首先引入基础的SIR模型,通过增加参数,离散化微分方程等方式进行建模,用软件模拟不同人口流动模式的人群变化过程。实验结果表明,人口流动数量与感染首次到达峰值的时间存在对数关系。从数学角度证明抑制人员流动确实是控制肺结核传播的有效途径之一。但是,因实验在一些不切实际的假设基础上进行,整体实验数据较少,进一步的研究也许可以集中研究这种特殊对数关系产生的原因,并不断改进模型中不现实的方面。