H1N1 甲型流感全球航空传播与早期预警研究

2010 年 第 55 卷 第 12 期：1128 ~ 1133 www.scichina.com csb.scichina.com 英文版见: Chang C Y, Cao C X, Wang Q, et al. The novel H1N1 influenza A global airline transmission and early warning without travel containments. Chinese Sci Bull, 2010, 55, doi: 10.1007/s11434-010-0180-9 论 文 《中国科学》杂志社 SCIENCE CHINA PRESS H1N1甲型流感全球航空传播与早期预警研究 常超一①⑥, 曹春香①*, 王桥③, 陈玉⑤⑥, 曹志东④, 张颢①, 董磊①⑥, 赵坚①⑥, 徐敏①⑥, 高孟绪①⑥, 钟少波①, 何祺胜①⑥, 王劲峰②, 李小文① ① 遥感科学国家重点实验室, 中国科学院遥感应用研究所, 北京 100101; ② 中国科学院地理科学与资源研究所, 北京 100101; ③ 中国环境保护部卫星环境应用中心, 北京 100012; ④ 中国科学院自动化研究所复杂系统与智能科学重点实验室, 北京 100190; ⑤ 中国科学院对地观测与数字地球中心, 北京 100101; ⑥ 中国科学院研究生院, 北京 100049 * 联系人, E-mail: cao413@irsa.ac.cn 2009-10-31 收稿, 2010-01-14 接受 国家重点基础研究发展计划(编号: 2007CB714404)、国家自然科学基金(批准号: 40871173)和国家科技重大专项课题(编号: 2008ZX10004-012) 资助项目 摘要 一种新型流感(H1N1 甲型流感)已经全球蔓延. 早期的研究显示, 在没有采取防控措施 的情况下, 国际航空旅行可能引发一次严重的全球流感大流行. 在本研究中, H1N1 甲型流感早 期在墨西哥和美国的爆发被运用于 SEIR 流行病传播模型, 用以估计该病毒的传染性及疫情初 步参数, 构建了一个基于实时航空数据的新的空间传播模型以评估 H1N1 甲型流感从墨西哥传 播到世界各地的潜在风险. 结果显示, H1N1 病毒的基本再生数 R0 约为 3.4, 有效再生数随防控 措施下降明显, 研究还发现西班牙、加拿大、法国、巴拿马和秘鲁是最有可能发生严重的流感 大流行的国家. 关键词 H1N1 甲型流感 航空传播 早期预警 基本再生数 防控措施 近几个月, 由一种新的甲型流感病毒 H1N1 病毒 所引发的传染病, 开始在世界各地蔓延[1]. 2009 年 4 月 29 日, 世界卫生组织(WHO)发出对 H1N1 病毒全 球传播的警告. 大规模的传播已被证实(卫生组织猪 流感 报告 软件系统测试报告下载sgs报告如何下载关于路面塌陷情况报告535n,sgs报告怎么下载竣工报告下载 , http://www.who.int/csr/don/2009_05_05/en /index.html). 航空旅行被怀疑是一个长距离传播的 关键因素 . 迅速上升的感染者数量以及现代十分便 捷的交通网络引发流行病学家的担心——这场流感 是否会成为另一场类似 1918 年的全球流感大流行? H1N1 甲型流感是一种人际间传播的传染病, 通 过直接的身体接触或飞沫传播[2]. 飞机旅行是一种密 切接触的理想的封闭环境 . 虽然目前没有明确证据 表明该病毒是通过飞机向世界传播的 , 但全球航空 网络是病毒传播的重要场所. 例如, 美国、加拿大、韩 国和中国香港报告的许多案例都是从墨西哥旅行不 久后确诊 , 这意味着航空旅行正在帮助该病毒向世 界各地蔓延. 世界卫生组织 2009 年 4 月 27 日宣布墨 西哥成为疫区以后 , 许多国家颁布措施限制到墨西 哥的国际旅行. 我们需要构建一个合适的数学模型描述该病的 流行病学特征, 包括感染率、平均感染力时间及基本 再生数. 经典的流行病学 SIR 模型, 通常被用来模拟 缺乏有效的控制策略下的流行病学过程[3]. SIR 模型简化了复杂的传染病过程. 目标人群被 划分为 3 类人群: 易感者、感染者和移出者, 并用微 分方程模拟三者之间的动态变化. SIR 模型假设所有 群体完全均匀混合 , 因此易感者具有同样的感染概 率. 但是, 当甲流爆发初期感染者有限的情况下, 易 感者亦相对有限, 这些假设即会失效. 而易感者-暴 露者-感染者-移出者模型(简称 SEIR)将有限的暴露 1129 论 文 者考虑进 SIR 模型[4]. 通过粗略估计接触率和死亡/ 移出率, 该 SEIR 模型即可模拟甲型流感的流行病学 过程. 然而, 病毒的病理学特征因人而异, 因此, 鲜 有研究分析了这一新的新型病毒的空间分布特征 . 而研究感染者的空间分布对制定公共卫生防控措施 相当重要. 传染病防控措施应该是模型的反馈 , 它对疾病 的感染率有影响. 2009 年 5 月 2 日, 世界卫生组织将 甲型流感全球流行的警戒级别提升到第 5 级(共 6 级), 这意味着本地传播在墨西哥以外的其他国家发生 . 如果没有有效的防控措施遏制, 病毒蔓延即会失控. 2009 年 6 月 11 日, 世卫组织再次将全球流行的警戒 级别升至最高的 6 级. 在这 40天里, 许多国家采取了 严格的防控措施. 因此, 当地的传输过程中一定程度 上被抑制 , 而航空旅行则是病毒全球传播的主要载 体. 这期间, 全球各国的防控策略发生了变化. 例如, 中国改变了将防疫工作的重点从密切健康转而以社区 防控为主. 在不考虑防控措施的复杂性的前提下进行 全球疫情模拟是没有意义的. 因此, 研究病毒的这种 大的空间尺度的扩散取决于如何建立一个通用的模型 来估计国际航空旅行对病毒传播的影响. 即使没有明 确的证据显示, H1N1 病毒是从墨西哥起源的, 许多的 案例将甲型流感的爆发和从墨西哥城旅行回国的乘客 联系起来[2,5]. 如果缺乏以墨西哥城为始发地的航空旅 行信息, 疫情动态将难以估算. 本文的研究目的是模拟正在流行的 H1N1 甲型 流感早期阶段的全球传播及对其空间扩散进行早期 预警. 我们从互联网上收集了 2009 年 7 月 3 日以前 的所有数据, 引入一个 logistic 模型到 SEIR 模型, 用来 有效地估算防控措施对疫情的影响 . 公共卫生防疫 措施被认为是控制疫情传播的一个有力举措. 此外, 墨西哥航空旅行被认为是初步估计潜在空间扩散的 首要因素. 在这里, 我们结合各地的航空可达性, 对 H1N1 甲型流感的流行病学特征, 空间分布以及航空 旅行扩散进行探索性的分析、估计和预测. 1 数据整理及模型 1.1 基于网络的数据整理 全球甲型 H1N1 的发病数据来自于互联网 . Henry Niman 博士通过 Rhiza 实验室和谷歌跟踪疫情 (Henry Niman, http://www.flutracker.com). 甲型 H1N1 流感病毒是一种全新病毒, 人类对其没有免疫力, 这 表明 H1N1 病毒“人-人”的传播可能导致数百万人死 亡的大流行. 因此, 及早发现疾病活动并迅速有效地 采取措施可以减少季节性和大流行性流感的影响. 通过在线搜索引擎来检索用以监测健康行为是 一种改善早期预警的方法[6], 这些检索结果由世界各 地的数以百万计的用户每天提交 . 数据来源包括官 方来源, 新闻报道和用户的贡献, 每天更新多次. 这 些数据中不仅有实验室确诊病例 , 而且有通过强大 的搜索引擎处理得到的具有类似症状的敏感案例 . 幸运的是, 数据中包含有关 记录 混凝土 养护记录下载土方回填监理旁站记录免费下载集备记录下载集备记录下载集备记录下载 的详细说明, 并提供 相对准确的所有感染的空间位置 , 这是我们研究的 潜在基础. 航空旅行在传播疾病方面起作用的概率很高 . Brownstein 等人[7]的研究利用 1996~2005 年美国多个 城市的流感死亡的数据(由美国疾病控制和预防中心 (CDC)提供), 连同国内国际航线的旅客估计人数做 统计分析 . 这项研究表明基于计算机模型模拟的结 果, 航空旅行可能在传播疾病方面起作用. 由于航空旅行是全球传播的重要因素 , 我们从 实时航空旅行互联网(http://feeyo.com/)收集大型机场 航线的相关信息. 就我们所知, 各国的航空统计数据 通常在年底由政府计算, 由于季节性差异, 实时航空 统计数据很难估计. 然而, 我们的研究中应用的数据 是在 5 月份获取的 2009 年 4 月实时售票数据, 航空 公司根据售票情况调整航线 , 通过提高上座率以降 低航空公司的成本 . 飞往墨西哥的旅客数量与全世 界甲型 H1N1 确诊病例频率具有相关性[2]. 世界各地 的 135个城市被选入航空监测系统[8], 其中 60 个位于 美国, 其余的分布其他疫区以表示潜在的传播趋势. 基于互联网的信息提供航班号, 起飞、到港时间以及 飞机的类型 (http://www.carnoc.com). 出于经济方面 的考虑, 假定所有的航线都为经济最优化的满载, 因 此每个航班的乘客数量可以通过乘以该类型飞机的 最大载荷计算. 1.2 建立流行病时间动态模型 在克马克-麦肯德里克模型是一个固定人群的经 典的 SIR模型. 它被用来描述流行病中传染病感染者 数量随时间快速上升和下降的过程, 比方说鼠疫(伦 敦 1665~1666 年, 孟买 1906 年)和霍乱(伦敦 1865)[9]. 它假定在流行病爆发期间人口规模是固定的(即人口 总数不变), 没有潜伏期 , 感染周期和疾病流行时间 2010 年 4 月 第 55 卷 第 12 期 1130 一致. 同时, 它还假定该人群的年龄、空间分布, 或 社会构成完全同质. d ( ) ( )( ) , d S S t I tt t N β= − (1) d ( ) ( )( ) ( ), d I S t I tt E t t N β σ= − (2) d ( ) ( ), d R t I t t γ= (3) 其中, t 是时间, ( )S t 表示易感者数量, ( )I t 为感染者 数量, ( )R t 为移出者数量(康复或者免疫或者死亡者 数量), ( )tβ 为感染率, ( )tγ 为移出率. 基本再生数 R0, 等于接触率除以移出率 , 即在 平均感染期内一个感染者所感染的易感者数量 . 基 本再生数表示该传染病的传染力 . 在传统的无控制 措施的情况下, 当 R0<1, 每个感染者将会在传染期 内传染少于一人给下一代, 因此疫情将逐渐消亡. 当 R0>1, 每个感染者将会在传染期内将传染给超过一 人, 因此疫情会继续蔓延. Kermack-McKendrick 模型 现在又重新得到重视 . 更多种不同的 Kermack- McKendrick 模型能更好地反映实际的特定疾病的生 物学特征. 然而 , 不同人群均匀混合的假设在现实中几乎 不存在. 此外, 有效的控制措施以及接种计划会改变 疾病的扩散. 事实上, 控制措施会在某种程度上降低 基本再生数, 因此, 平均每个感染者传染的实际的第 二代感染者的数量被称为有效再生数. 当有效再生数 >1 的时候, 传染的趋势会愈演愈烈, 而当有效再生数 <1, 传播链会不可避免最终消亡. 接种、抗病毒药物 以及隔离也能阻止病毒传播. 因此, 当我们考虑控制 措施的时候, LSEIR 模型可以被用来更好地模拟 H1N1 甲型流感在时间维上的传播和扩散[10]. 在 LSEIR 模型 中, 一个 logistic 模型被用来调控模型中的接触率. 它 用 4 个参数分别来响应第二稳态值、阈值、曲率以及 感染事件, 能够很好地拟合流行病的动态[11]. d ( ) ( )( ) , d S S t I tt t N β= − (4) d ( ) ( )( ) ( ), d E S t I tt E t t N β σ= − (5) d ( ) ( ) ( ), d I E t t I t t σ γ= − (6) d ( ) ( ), d R t I t t γ= (7) 2 1 3 4 ( ) , 1 exp[ ( )] p t p p t p β = + + × − (8) 其中暴露者人群与易感者人群在时刻 t 的有限接触次 数用函数 ( )tβ 来表示. 参数 p1, p2, p3, p4 分别代表基 台值、值域、曲率以及位移. 这些参数具有流行病学 特征: (p1+p2)/γ可模拟基本再生数, 而 p1/γ表示采取控 制措施后的有效再生数. 1.3 建立流行病空间动态模型 我们收集了 2009 年 4 和 5 月全球 144 个重要机 场的航班数据(通过国际航空运输协会的机场编码生 成网页搜索 , 然后记录搜索结果). 判定是否进入机 场搜索名单的依据是根据其在国际贸易中的排名以 及是否是某个区域或者国家机场来筛选的 . 在北美 洲和南美洲, 我们增加了机场的样本数, 以便模拟能 覆盖感染者的全球旅行的范围. 部分感染者甚至可能会通过航空旅行从墨西哥 城前往其他国家. 由于乘机率很难被测算, 我们假设 从墨西哥城乘机出发的感染者占总感染者的比例是α. 而考虑到部分乘客可能会二次中转, 转机率设为μ. 那 么假设航线 AB1 通过 C 型飞机连接了 A 市和 B 市, 那 么 AB 两市之间通过航线 AB1 旅行的乘客数则为 C AB1 ,7 n S P = (9) 其中, Cn 表示 C 型飞机的最大载客量, S 表示该航班 每周的班次. 因此, 通过 n 次航班往来 A, B 两市的乘 客数为 AB AB 1 . n iP P= ∑ (10) 那么, 从墨西哥城传播到 A 城的感染者数为 A 1 A A 1 ( ) , n i n N i N I t P I P P α = ∑ ∑ ∑ (11) 其中 , A 1 n N i N P P∑ ∑ 表示从墨西哥城出发的总的乘客 数, A 1 n iP∑ 则表示前往 A 城的乘客数. 如果 B 城没有 1131 论 文 来自墨西哥城的直达航班而需通过墨西哥城—A城—B 城中转, 那么从墨西哥城传播到 B 城的感染者数为 A AB 1 1 A 1 ( ) . n n i i B n N i AN N N I t P P I P P P α μ= × × ∑ ∑ ∑ ∑ ∑ (12) 所有的两次以内的中转航班都考虑进了我们的 模型 , 这其中包括直接中转的航班(例如 , 从墨西哥 城—B 城—A 城), 也包括中转两次的航班. 在考虑中 转两次的航班的时候 , 我们去除了目的地重复出现 在中转地的航线. 例如从墨西哥城—B 城—A 城—C 城 算作墨西哥城—C 城, 但墨西哥城—B 城—A 城—B 城 不考虑进墨西哥城—B 城以内. 根据目的地, 所有从墨 西哥城出发通过一次中转和二次中转前往其他城市 的感染者被计算出来. 2 讨论 正如研究中使用 LSEIR 模型所能得到的一样 , 我们拟合了它的 6 个参数. 从 2009 年 4 月 22 日到 5 月 11 日墨西哥城的 H1N1 甲型流感大爆发被用来估 计这 6 个基本参数, 而之后到 7 月 3 日的数据被用来 检验模型结果 . 这期间墨西哥实际累计爆发病例与 模拟爆发病例的比对如图 1 所示. 因为 LSEIR 模型中 4 个微分方程存在于成对的 解的集合, 所以估计最优解集合会较为困难. 我们就 采用 Runge-Kutta 公式来解这些方程 . 一般来说 , Runge-Kutta 公式对于普通微分方程来说是一个有效 的数学最优解. 然而, 因为这些参数存在于一个六维 实数空间, 所以我们必须尽量减小误差. 在本研究中, 一个计算最优化解法——遗传算法被用来解这些方 图 1 墨西哥城实际累积爆发病例与模拟爆发病例比对 程. 利用基因遗传和变异类似的生物过程的原理, 遗 传算法被广泛用来解决最优化估计. 在本研究中, 我 们采用一个常规的遗传算法估计解参数 , 置信水平 为 90%. 2009 年 7 月 3 日以后, 世界卫生组织不再更 新甲流数据, 所以之后的参数很难被正确估计. 根据墨西哥城的城市规模 , 疫区人口设定为 1900 万人. 最初的感染者设定为 1 人, 乘机率和中转 率设为 0.1. 通过多次模拟, 我们计算出所有参数解. 根据 LSEIR 模型的定义, 可以确定出 H1N1 的流行病 重要参数. 其中基本再生数 R0为 3.379415, 表示在传 染初期的时候平均大约有 3.4 个易感者被 1 个感染者 传染. 但是在实际防疫措施的控制下, 有效再生数回 落至 1.8 左右. 不仅如此, 从一个传染个体到一个大 爆发的可能性可能在不同的情况下变化非常大 . 尽 管超级传播事件还没有像 SARS 这样被报道过, 但是 根据我们早期发现 , 被不同传染个体所传染的数目 视不同情况变化明显. H1N1 甲型流感的基本再生数 是低于于 SARS 和其他一些传染性疾病的, 但是高于 Fraser 的估计(R0=1.4~1.6, 基于截止到 4 月 30 日的墨 西哥数据 )[3]. 潜伏期和感染期分别是 1.000007 和 2.96765 d, 因此传染期共为 4 d(根据美国疾病控制中 心对个案研究的报告, H1N1 甲型流感的潜伏期通常 是 24~48 h)这个结果也显示, 在有效的防疫措施的控 制下, 有效再生数最终会回落到 0.417797, 意味着墨 西哥城的 H1N1 甲流疫情最终会得到控制. 如果我们 将乘机率和转机率调整为 0.2, 即可得到二次模拟的 结果. 通过全球航空旅行的传播, H1N1 甲型流感将会 传播到其他国家. 在 WHO 迅速的警告后, 很多国家 已经采取了控制措施. 在我们的研究中, 全球航空旅 行对于传染病的分布起了重要作用 . 考虑到感染持 续期不会超过 7 d(在我们研究中近似为 4 d), 并且在 传染后, 飞机旅行通常会花费 1~2 d, 传染性在国际 旅行中会降低. 并且在其他国家, 数目更多的易感者 将会被从墨西哥返回的二代感染者和三代感染者所 感染 . 航空传播减缓了传播速度但是增加了他国旅 客感染的风险. 在我们的研究中, 大多数城市都是这 个国家的首都, 或者是国际性城市. 所有城市与墨西 哥城通过航空旅行的中转城市不超过 2 个, 乘机率和 中转率最初设为 10%. 我们考虑了全球传播的 3 个不 同阶段: (1) 最初阶段: 最初在墨西哥的疫情开始向 全球传播 ; (2) 疫情失控阶段 : 墨西哥局部大爆发 , 2010 年 4 月 第 55 卷 第 12 期 1132 超过 800 人被感染(根据 WHO 报告, 在 4 月 28 日, 墨 西哥城有 800 名感染者同时还有超过 1000 人的疑似 病例需要被确认); (3) 最终阶段: 疫情逐渐消退. 最 初, H1N1 的病症和病理学都是不清楚的, 直到美国 疾病控制中心实验室确诊了美国加利福尼亚州的两 个病例. 因此大爆发并未被注意到, 潜在的传播持续 了好几天 . 有限的在本地居民和游客之间的交叉传 染有可能发生 , 但是并没有引起对于潜在的全球大 流行的国际关注. 在 4月 27 日, 加拿大和西班牙报道 病例后, WHO 将全国警告等级提高到 4 级, 这意味着 这个病毒可能导致“社区级大爆发”. 两天后, 警告等 级被调整到 5 级, 意味着全国性的大流行即将到来. 我们评估了从墨西哥到世界可能性传播. 结果显示, 墨西哥的邻国和发达国家, 拥有更高的 H1N1 入侵风 险(图 2). 这可能是由于发达国家人们热衷于前往墨 西哥旅游的原因. 由于 WHO 没有公布 2009 年 5 月 11 日世界各地爆发的甲型流感病例数, 因此我们用 5 月 22 日公布的数据来比对模拟的结果(图 3). 这个研究的误差可能来自于不同的防御策略的 干扰. 因为采用了防疫政策, H1N1 甲型流感将倾向 于灭绝 , 而且传染的随机性造成了潜在的地方性爆 发的差异. 3 结论 本研究采用实时航空数据对全球传播的甲型 H1N1 流感进行粗略的估计. 大尺度的空间扩散清楚 的描绘出病毒是通过航空旅行传播的 . 研究证实 H1N1 和航空旅行之间存在确实的关系, 基于这些发 现可以做出预测. 此外, 对这种新型 H1N1 流感病毒 的流行特征, 包括基本的生殖数量, 传染性时间进行 了估计并对控制策略进行了讨论 . 伴随世界卫生组 织的警告, 各国采取一些措施来控制 H1N1 病例数, 包括航空公司的限制 , 严格的流行病学监测和相关 的公众意识运动, 这些策略似乎已经控制了 H1N1 病 毒的传播 . 本研究特别分析了从墨西哥通过航空旅 行的病毒传播, 结果表明在国际航线上控制 H1N1 病 毒的传播是避免疫病大流行的首选 . 卫生和抗病毒 药物应该首先用在飞机上 , 因为狭窄的车厢是病毒 图 2 根据航空路线模拟 2009 年 5 月 11 日甲流感世界范围内爆发的分布图 1133 论 文 图 3 2009 年 5 月 22 日甲流感实际爆发案例分布图 传播的理想地点 , 而交叉感染延长了病毒的传染时 间. 由于控制措施对 H1N1 病毒扩散的具有很强影响, 今后的研究将集中于找出特定遏制政策与甲型 H1N1 流感的传播之间的关系. 参考文献 1 陈继明, 孙映雪, 刘朔, 等. 2009 年北美出现的甲型 H1N1 流感病毒的来源和去向分析. 科学通报, 2009, 54: 1657—1660 2 Fraser C, Donnelly C A, Cauchemez S, et al. Pandemic potential of a strain of influenza A (H1N1): Early findings. Science, 2009, 324: 1557—1561 3 Kermack W O, McKendrick A G. A contribution to the mathematical theory of epidemics. Proc Roy Soc Lond A, 1927, 115: 700—721 4 Anderson R M, May R M. Infectious Diseases of Humans: Dynamics and Control. London: Oxford University Press, 1991 5 Daniel C. Swine flu jumps continents. Nature, 2009, 458: 405 6 Ginsberg J, Mohebbi M H, Patel R S, et al. Detecting influenza epidemics using search engine query data. Nature, 2009, 457: 1012—1014 7 Brownstein J S, Wolfe C J, Mandl K D. Empirical evidence for the effect of airline travel on inter-regional influenza spread in the United States. PLoS Med, 2006, 3: e401 8 Declan B. Swine flu goes global. Nature, 2009, 458: 1082—1083 9 Hethcote H W. The mathematics of infectious diseases. SIAM Rev, 2000, 42: 599—653 10 Wang J F, McMichael A J, Meng B. Spatial dynamics of an epidemic of severe acute respiratory syndrome in an urban area. Bull World Health Org, 2006, 84: 965—968 11 曹志冬. 时空数据的动态预测模型与空间分析技术研究——以广州 SARS 流行为例. 博士学位论文, 北京: 中国科学院研究生院, 2008. 72—80