用建模仿真看疫情数据
来源：《北京市新冠肺炎疫情数据初步分析》
作者：华如科技 季秀才
导 读
突如其来的新型冠状病毒肺炎疫情已影响到我们每个人的工作生活。疫情变化数据及趋势是大家密切关注的焦点。疫情发展到了什么阶段，病毒传播强度具体有多大，目前的应对措施是否得当，疫情是否已得到了控制，什么时候会向好的方向发展，拐点是否出现，还需多久结束，生活和工作何时能够恢复正常等等，这些问题会持续萦绕在我们每个人的心头。在关注北京疫情数据的同时，我尝试通过数学模型和方法对疫情数据进行分析，来定量解答上述疑问，并初步整理了相关分析过程与结论。本文是对个人数据分析工作的初步总结，文中的分析及有关论述均为本人学术观点，欢迎有兴趣的同仁给予批评、指正，并一起开展更为深入的分析工作。
一
数据及概念说明
(一) 数据来源
本次分析所用基础数据由北京市卫健委在官网上发布的北京市新冠疫情通报整理获得，数据的时间跨度为2020年1月20日至2020年2月14日，共26天（14日以来每日所公布数据基本符合本文模型推算趋势，对此来源数据无需再追加补入），如表1所示，其中空缺数据未公布。
（表1：北京卫健委公布北京市新冠肺炎疫情数据汇总表）
北京市卫健委官网 ：http://wjw.beijing.gov.cn/wjwh/ztzl/xxgzbd/gzbdyqtb
北京市卫健委每天发布的北京市疫情通报中会公布前一天的北京市新冠疫情数据，各天公布的疫情数据项类型并不完全相同。其中，每天都公布的疫情数据项主要包括：
累计确诊病例数量；
每日新增确诊病例数量；
新增确诊病例为湖北及其他省份接触者数量；
每日新增确诊病例为其他确诊病例密切接触者数量；
每日新增确诊病例既有湖北及其他省接触史又是密切接触者数量；
每日新增确诊病例为无湖北接触史者数量；
累计治愈病例数量；
累计死亡病例数量。
每日疑似病例数量、累计密切接触医学观测人数、累计解除医学观察人数等数据是从2月6日起开始公布。除表 1所列数据外，疫情通报中还会公布部分病例的性别、年龄、由湖北省返京时间、开始有呼吸道症状时间、就诊时间，以及确诊病例在北京各区分布数量等信息。本文主要以表 1所列数据为基础数据开展相关分析工作。
（二）基本概念
1. 确诊病例和疑似病例
确诊病例、疑似病例通过流行病学调查和临床表现两方面来确定。国家卫健委在其官网发布的《新型冠状病毒肺炎诊疗方案（试行第五版 修正版）》中，明确给出了新型冠病毒肺炎的确诊病例和疑似病例的认定标准：
上述方案中也规定了新型冠病毒肺炎疑似病例的处置和排除条件：
数据来源：http://www.nhc.gov.cn/xcs/zhengcwj/202002/d4b895337e19445f8d728fcaf1e3e13a.shtml
在当前全社会对疫情极度重视程度下，北京市的相关医院会很快对就诊患者进行流行病学调查和临床检验，因此在分析时我们可以认为北京市公布的新冠肺炎确诊病例和疑似病例为当天就诊病例，并在当天对其实施了隔离治疗。新冠肺炎病毒感染者在被认定为确诊病例和疑似病例前具备对外的传染性，在认定后由于被隔离不具备对外传染性。
另外，需要明确的一点是疑似病例可能并不是新冠肺炎病毒感染者，会在持续的治疗中持续的进行筛查识别，被确认的将成为新的确诊病例；排除的将排除疑似，认定为新冠肺炎病毒非感染者。
2. 密切接触者及医学观察
国家卫健委在其发布的《新型冠状病毒肺炎防控方案（第四版）》的附件3《新型冠状病毒肺炎病例密切接触者管理方案（第四版）》中，给出了新型冠状病毒肺炎密切接触者的判定标准：
同时，该方案也规定了密切接触者应采取的处理措施：
数据来源：
http://www.nhc.gov.cn/xcs/zhengcwj/202002/573340613ab243b3a7f61df260551dd4.shtml
上述规定明确新型冠病毒肺炎确诊病例和疑似病例的密切接触者都要实施医学观察。因此在分析中暂且认为北京市公布的医学观察人员全部为新冠肺炎确诊病例和疑似病例的密切接触者。
二
直观分析
（一）人员类别及比例分析
1. 人员类别分析
北京市公布的新冠肺炎疫情数据相关的人员类别包括：确诊病例（IC）、疑似病例（SC）、死亡病例（DC）和实施医学观察的密切接触者（CC）。
各类人员中，确诊病例比较复杂，公布的疫情信息中进行了进一步的细分，包括四类：
具有湖北及其他省份接触史的确诊病例（IC_HC）；
为其他确诊病例密切接触者的确诊病例（IC_CC）；
既有湖北及其他省接触史又是密切接触者的确诊病例（IC_HC&CC）；
无湖北及其他省份接触史同时也不是密切接触者的确诊病例（IC_N）
各类确诊病例间的关系可以用如下的集合图说明。
（图 1：各类新冠肺炎间关系图）
IC_HC和IC_CC的交集为IC_HC&CC。C_HC和IC_CC是可以通过流行病学调查确定病毒传播链的病例，而IC_N很可能是流调难以确定传播链的病例。可以比较肯定地认为：
1) IC_HC与C_HC&CC的差集，即IC_HC\C_HC&CC(图1中的棕色部分)，为外部输入到北京的新冠肺炎感染者，暂且称之为输入确诊病例；
2) IC_CC与C_HC&CC的差集与IC-N的合集，即（IC_CC\C_HC&CC）∪IC_N为在北京当地被感染的病例，暂且称之为当地确诊病例；
C_HC&CC有可能是输入病例也有可能是当地病例，在相对乐观的估计中可以将其作为输入病例对待，在相对保守的分析中可以将其作为当地病例对待，本文持相对乐观态度，即认为在严密的防控下北京的就地病例相对较少，后续将IC_HC全部作为输入病例进行分析。
另外，此次疫情发生在春节期间，以及北京的人口流动特色等方面也值得剖析：
1) 北京人口在春节前主要向外流，春节前输入北京的人口相对少，春节放假期间输入北京的人口相对较少，因此相对于疫情发生于其他时间而言，此时发生疫情北京的输入病例相对少；另外，也可以认为大部分输入病例于产生前来到北京；
2) 基本上是在除夕前后全国疫情形势紧张起来，而此时开始春节放假，大多数人不上班，市内人员流动也少了，使得人们密切接触的大多为家人，而工作、交通、朋友聚会等密切接触相对较少，客观上减少了人员的接触数量，同时也有利于人员的管控，这可能是密切接触人员总量相对较少的重要原因；
3) 复工后会有大量人员返京，有引进新的输入病例的风险。
2. 人员比例分析
每天各类人员的累计数量在不断变化，截止2月14日，最新的各类确诊病例的数量和比例如表 2所示，当地比例中各类人员的比例如表 3所示。
（表2：2月14日各类新冠肺炎确诊人员数量和比例类型）
（表3 ：2月14日各类新冠肺炎当地确诊人员数量和比例类型）
从表 2和表 3中的数据可以看出以下几点：
1) 到目前为止北京地区输入病例占确诊病例的一半以上；
2) 大部分的北京地区当地病例为密切接触者，说明北京地区的防控效果显著，未造成大范围的不可追踪的传播；
3) 目前的死亡率也比较低，一方面说明确诊病例得到的有效的救治，另一方面也说明在医疗资源相对充足的情况下新冠肺炎的死亡率较低，当然该病的死亡率还要看疫情结束后的统计结果。
（二）基本概念
我们基于基础数据计算了累计输入确诊病例、累计当地确诊病例、累计治愈病例等数据，具体如表 4所示；各数据项的数时间发展曲线图如图 2、图 3和图 4所示。
（表4：累计输入确诊病例、累计当地确诊病例、累计治愈病例列表）
（图 2：累计输入确诊病例数量数据随时间发展曲线图）
（图 3：累计当地确诊病例数量数据随时间发展曲线图）
（图 4累计治愈病例数量数据随时间发展曲线图）
由上述图表我们分析获得以下几点：
1) 2月7日至2月14日连续8天，输入病例数量每天只增加1或2位，大多数时间不增加，说明北京本轮的输入确诊病例已基本全部识别完了；
2) 当地确诊病例的增长在2月12日以后每天以逐渐减小的个位数增长，增长幅度越来越少，很快我们也会将北京当地确诊病例全部识别出来，北京本轮拐点似乎已经到来了；
3) 2月12日以来治愈率增长迅速，按照发展发展趋势相对保守估计在3月8日以前，已经确诊的大部分病例将会痊愈出院。
如果当前在北京的人们保持密切接触情况不变，当连续14天无新增确诊病例，北京的疫情可视为基本结束。
三
模型分析
应用数学模型和方法对流行病传播过程进行分析可以获得更深层次的认识。基于北京市新冠疫情实际，本文对经典的“流行病仓室模型”（网上可以找到很多有关论述）进行了改造，并结合实际疫情数据，对北京当前的新冠疫情的发展过程进行了分析。
（一）基本假设
为了简化问题，本文基于前述分析提出如下假设：
1.新冠肺炎病毒携带者在认定为确诊病例和疑似病例后将被隔离，不具备对外的传染性；
2.新冠肺炎病毒携带者在潜伏期和发病期具备相同的传染性：
3.北京市公布的医学观察人员全部为新冠肺炎确诊病例和疑似病例的密切接触者；
4.只考虑新冠肺炎病毒携带者对其密切接触者的传染；
5.北京市的输入确诊病例全部是在春节前进入北京，北京卫健委在疫情通报中透露，最早的一例有武汉旅行史的确诊病例1月4日返京，1月9日出现发热症状，1月14日在医疗机构就诊，因此假定输入确诊病例进入北京的时间段为1月4日至1月23日；
6.只要是新冠感染者，其就诊被认定为疑似病例后能够很快检测认定为确诊病例，疑似转确诊的时间忽略不计；
7.治愈出院病人不会再被传播感染新冠病毒。
（二）模型描述
1. 人群分类
传统的仓室模型把人群分成了5类：易感人群、潜伏期人群、发病人群、治愈人群和病死人群。
本文对人群进行如下分类：
S 类，易感者人群(Susceptible)，指未得病者，但缺乏免疫能力，与感染者接触后容易受到感染；
H类，潜藏者人群（Hidden），指未就诊被识别为确诊病人或疑似病人的新冠病毒感染者人群染上新冠病毒的人群，可以传染S 类成员，其成员即可能处于潜伏期、也可能处于发病期；进一步细分为输入潜藏者人群—IH和本地潜藏者人群—LH；
I 类，确诊病例人群 (Infectious)，被隔离，不会传染S 类成员；进一步细分为输入潜伏者II和本地潜伏者人群LI；
R 类，康复者 (Recovered)，指因病愈出院而具有免疫力的人群，不会再被感染；
D类，病逝者（Dead）；
E类，外部感染源（External）。
2. 建模
上述6类人群的成员数量时间变化过程应用以图 5所示的图形描述。
（图 5各类人群转移过程图）
图 5中各个参数的意义如下：
（1） β，新冠病毒感染者对其接触人群的日接触感染率；
（2） γ，新冠病毒感染者的日确诊率；
（3） κ，确诊病例的日治愈率；
（4） μ，确诊病例的日死亡率；
（5） λ，日输入病例率。
图 5描述的过程应用微分方程模型表达如下:
可以简写成如下表达式：
从公开的基础数据的基本趋势上我们基本能够分析疫情的走势。但是对于疫情的深入认识，包括传播的强度、传染的周期等都潜藏于模型的系数当中，其中：
1/γ可以表征传染的周期，即一个病毒感染者在感染病毒到确诊被隔离的平均时长；
R0=βS0/γ，称为传染病的基本再生数，是指一个病毒感染者进入这个易感人群到他被确诊隔离平均会传染的人数，表明了这个传染病在这个群体中的传播强度，R0<1,传染病逐渐消亡，R0>1传染病将不断扩大；传染病控制就是要将R0将控制在1以下，而且越小越好。
（三）模型求解
1. 初值设定
（1） 输入病例总量：201例
以2月14日累计的有湖北及外省接触史确诊病例为输入病例总量，共计201例。这些病例全部于2月20日前进入北京。
（2） 易感人群总量：1304人
以确诊病例在确诊病例和疑似病例总和的占比为比例，从2月14日累计医学观察人员总量中计算，1843*375/（175+375）=1304人。
2. 求解结果
我们应用公开的疫情数据对模型参数进行可反演计算，近似的估计获得的模型参数取值如下：
接触感染率：β=4.13E-5;
确诊率：γ=0.11
基本再生数：R0=βS0/γ=4.14E-5×1304/0.11=0.4908
治愈率：κ=0.026,治疗周期38.5天
（四）分析
通过对模型参数的解读可以获得北京2月14日前的疫情基本特征：
1. 1/γ=9.09，平均传染期大约为9天，即一个新冠病毒感染者从感染到被隔离平均时间为9天左右。到2月14日，北京的疫情经过一个月左右，即3个左右的传染波次。
2. 近期，有报道称钟南山院士的团队研究表明新冠病毒潜伏期的中位数为3.0天。如果按照本文传染周期的计算结果，9.09-3.0=6.09，说明北京的新冠病例从发病到确诊平均时长为6天。
3. R0=0.4908，一个病毒感染者平均传染0.49个人。而近期一些报道称，新冠病毒在不受控的情况下的基本再生数至少在2以上，这说明北京的疫情防控措施将病毒的基本再生数控制在了较小于1，防疫控制措施还是非常有效的。
4. 2月14截止，北京的输入病例大约200人，按照 1个人平均传染大约0.5个人计算， 3个传染周期，确诊病例总量大约为375个，这与2月14日的累计病例数正好吻合。这一点也印证了本文模型的有效性。
5. 1/κ=38.5,确诊病例的平均治疗时间为38.5天。这说明北京已有病例的治疗周期是比较长的。2月1日和2月2日的疫情通报中公布了当天的出院病例的平均住院治疗时间分别为10.25天和9.9天，与模型计算结果差距很大。但是，在近期制约率不断提高的情况下，根据2月11日到14日的疫情数据平均日治愈率大约为0.04，据此计算的平均1/0.05=25天，与公布的数据差距也比较大。分析原因，公布的住院时间是前期轻症患者的住院时间。另外，疫情前期治愈率很低，很长时间是0或0.01左右，平均住院时间很长，本模型计算疫情以来的平均治疗时间为38.5是有合理性的。
四
一些建议
如果北京新冠病例从发病到确诊平均时为6天，说明很多感染者发病后没有尽快的就诊。可能是担心牵扯身边人被隔离，或担心到医院就诊被传染，或受其他恐惧心理等影响，很多出现呼吸道症状的患者先是进行了自我隔离下的感冒治疗，发现症状无法控制后，自我确认感染病毒的可能性很高的条件下才去就诊，所以使得发病到确诊时间间隔较长。但是，这样做客观上极大地增加了密切接触者的感染风险。如果能第一时间就诊并确诊，使得平均传染期接近潜伏期，基本再生数将降到4.14E-5×1304×3.0=0.162，能够非常有效地减小病毒的传播强度。
因此，为了自己和与家人、同事等密切接触人群考虑，也为了疫情能够尽早结束考虑，建议大家出现新冠相关症状后能够尽快就诊排除或确认。如果不及时就诊很可能将轻症拖成重症，重症的治疗代价、病死率以及后遗症等完全不是轻症可比的（https://mp.weixin.qq.com/s/CCuMLwTs7bnm_f4ZHSU1CQ）。其次，建议政策上采取更多措施确保人们到医院就诊被感染的风险极小，并加强宣传消除大家的疑虑，并鼓励有症状的人员能够尽早就诊。
另外，虽然北京上一轮的病毒感染过程已经接近尾声，但是复工带来的新一轮的返京潮确实带来新的风险。复工之前，人们的密切接触人员只是家人，但是复工之后，同事也成为了密切接触人员，一个人的密切接触人员会成倍的增长，即S0成倍扩大，可能使得基础再生数R0=βS0成倍上涨，甚至疫情会出新一轮的反复，因此北京最近针对返京人员的隔离要求是非常有必要的：
一方面，建议要返京的人员能够在相对短的时间内集中返京，使得大家尽量同一时间隔离，能够在较短时间通过把感染者识别筛选出来，尽量缩短疫情持续时间；
另一方面，近期返京人员途中尽量做好防护，返京后切实做好个人隔离，同时其他人也要尽量控制自己的对外接触，保持警惕，再忍耐一段时间，使得疫情能够一直保持目前的下降势头，胜利会更早地来到。
至2月18日的最新数据表明，已经连续三天未有新增输入确诊病例，并且新增本地病例也是在个位数变化，并且18日治愈率出现了一个跃升，进一步说明疫情正在加速减弱。
（图 6：截止至2月18日的北京公布疫情数据）
五
结束语
以上是对北京近期疫情数据浅显分析结果的分享，受时间、经验、专业能力等方面限制，加之建模与分析计算相对简略，难免存在很多不足，甚至是错误；如果考虑时间滞后等因素，应用差分动力方程或概率模型等其他模型可能会得到更为精确的结果。欢迎大家批评指正。借用我们华如科技内部通知中提到的一句话，“没有一个冬天不可逾越，也没有一个春天不会到来”，只要我们科学认知、遵从法纪、严密防控，春风必将如期而至！