马 盼 王馨梓,2 张莉 唐小新 冉洪雨 王式功
1.成都信息工程大学环境气象与健康研究院，气象环境与健康联合实验室，成都，610225
2.天津市津南区气象局，天津，300350
3.深圳市国家气候观象台，深圳，518000
流感是第一个实行全球监测的传染性疾病，在全球每年导致300 万—500 万人发病、29 万—65 万呼吸道疾病相关死亡病例（张静等，2004；
World Health Organization，2017）。全人群对流感普遍易感，其中孕妇、婴幼儿、老年人和慢性病患者（哮喘、心血管疾病、肾病、糖尿病等）感染后危害更为严重（Nimbalkar，et al，2016）。鉴于流感病毒抗原易变、传染性强，对流感样病例（Influenza Like Illness，ILI）的实时监测是防控季节性流感、应对流感大流行的重要举措（Azziz Baumgartner，et al，2012）。ILI 定义为体温≥38℃，伴有咳嗽或咽痛、全身疼痛等症状的急性呼吸道感染病例；
ILI 已被证实可有效反映流感活动的基本规律（中国国家流感中心，2017；
翟红楠等，2009a，2009b；
王伟，2016）。
研究揭示，流感具有气候依赖性和病毒特异性，其暴发和流行与气象条件关系密切（Ye，et al，2019；
Alonso，et al，2007）。温带地区的流感高峰通常在冬季出现一次，北半球在11 月—翌年3 月，南半球在4—9 月（Finkelman，et al，2017；
Lofgren，et al，2007）；
而热带、亚热带地区的流行模式更加复杂，全年均可能发生、甚至一年内出现多次暴发疫情，流行的局地特征较强（Cheng X W，et al，2013；
Cheng Y H，et al，2016；
Liu T，et al，2018；
Liu X X，et al，2017；
Wang，et al，2017）。此外，在热带和亚热带地区的雨季，流感风险也较高（Soebiyanto，et al，2010；
Moura，et al，2009）。不仅病毒的活性显著受气温、湿度等因素影响，飞沫在空中传播也与多种气象条件有关。低温在多种气候背景下均被认为是影响流感发病的主要因素（Huang，et al，2017；
Soebiyanto，et al，2015）。其他气候因子如降雨量、气温日较差和风速，在不同气候带则具有差异性影响（Gomez-Barroso，et al，2017；
Shaman，et al，2011）。而截至目前，导致流感流行模式多变的机制仍不够清楚（Chong，et al，2015；
Dai，et al，2018）。
深圳地处北回归线以南，属热带向亚热带过渡的海洋性季风气候区，夏季炎热潮湿、伴随热带气旋或季风带来的阶段性强降雨；
冬季短，时常受北方干冷空气影响。作为华南的核心城市之一，深圳具有高度城市化（100%）、人口基数大（2019 年常住人口1343.88 万）、流动人口多、年龄结构特殊等突出特征（深圳市2019 年国民经济和社会发展统计公报，http://www.sz.gov.cn/）。尽管深圳的地理位置毗邻香港，两地气候特征相近，但流感的发病规律差异却较大（Tang，et al，2018）。作为中国的经济特区和重要港口城市，深圳应对城市健康风险的压力大、责任重，对健康气象服务的要求也更高，尤其是相关传染疾病的防控，亟待深入研究。
本研究收集深圳市流感监测哨点医院长序列（2003—2019 年）逐周ILI 数据，首先研究其季节特征与时间趋势；
进而探索ILI 发病风险与多种气象因素的非线性关联与滞后影响，揭示影响流感发病的主导气候因子及其阈值；
最后采用时间序列和多元回归2 种方法进行建模预测和对比验证。以期为相关部门开展健康气象服务和流感防控工作提供有效的科学依据。
2.1 数据来源
本研究所用流感样病例（ILI）资料来自深圳市疾病预防控制中心，包括深圳市流感监测站点2003—2019 年（逐周，共884 周）ILI 就诊人次（合计821305 例），以及每周ILI 人次占就诊总人数的百分比（ILI%），可理解为流感就诊率。除上述原始数据外，考虑到深圳市人口基数波动可能导致门诊人数变化，本研究以每1 万门诊数中的ILI 人数作为标准化的ILI，纳入时间序列模型研究与气象因子的关系。
依据国家流感监测方案，深圳市第一人民医院和妇幼保健院作为流感监测哨点，其内科、儿科、内/儿急诊科和发热科常年进行ILI 监测，专人收集、监测ILI 就诊数及门诊病例总数，并每周上报市疾控中心。两家医院均有20 a 以上的流感监测经验，监测流程完备，故所用数据的准确性、代表性和连续性较好。
相应时段的气象数据来自深圳市气象局，包括日平均/最高/最低气温、相对湿度、水汽压、日照时数、24 h 降水量、日平均风速、日平均/最高/最低气压等要素。气温日较差定义为每日最高气温与最低气温之差，由计算获得。为与疾病资料匹配，对各气象要素同样求其周平均值。由于环境空气质量监测开始的时间较晚（2013 年至今），环境数据暂不包含在本研究中。
2.2 统计方法
2.2.1 分布滞后非线性模型
采用R 软件3.6.1 版本中的dlnm 和mgcv 程序包进行时间序列建模分析，显著性检验水平（α）为0.05。分布滞后非线性模型（Distributed-Lag Nonlinear Models，DLNM）可以较为准确地刻画气象环境要素与某种健康结局的暴露-反应关系，量化滞后和累积效应（Gasparrini，et al，2010，2014）。通过相关分析大致确定影响ILI 的气象因子（P＜0.05），再利用DLNM 分别构建各主要因子的交叉基函数，基于多次敏感性试验，在避免共线性和过拟合的基础上建立ILI 就诊量-气候因子的关联。考虑到对个体而言每周的就诊为小概率事件，采用的广义泊松回归构造模型框架如下
式中，Yt是第t周的ILI 就诊人数，E（Yt）表示模型估计的第t周人数期望值；
cb（*）为气象因子（X）的二维交叉基矩阵，包含气温、相对湿度、风速和气温日较差等变量，lag 表示滞后时间；
s（*）是样条函数，rain 表示降水量的混杂影响，df 为自由度，设为4；
time 表示时间趋势，df 设为51；
holiday 表示节假日的混杂效应，作为哑变量引入。利用赤池信息准则（AIC）调整气象要素的非线性自由度和滞后自由度；
σ表示模型残差。
ILI-气象因子关联强弱用相对危险度（Relative risk，RR）表征，即要素单位变化量导致的ILI 风险与参考条件下（RR=1）的比值。参考值选取各要素累积危险度（CumulativeRR）最小的阈值。
2.2.2 Prophet 时间序列模型
本研究尝试使用Prophet 模型预测ILI 风险，并与传统的线性回归方法做对比。它由Taylor 等（Taylor，et al，2018）在2017 年提出，相较于传统时序模型更加灵活，并有良好的预测性能。其模型构成为
式中，y为被预测变量，t为时间；
g（t）是趋势函数，用来分析序列中的非周期性趋势部分；
s（t）是周期项，表示序列的周期性变化部分，由傅里叶级数拟合；
h（t）是节假日项，代表节假日等少数特定天数的影响；
ɛ为误差项，一般服从高斯分布。
对预测模型进行优度检验的指标包括平均绝对误差（MAE）、均方根误差（RMSE）、平均绝对百分误差（MAPE）、预测准确率（P）4 个参数，其中y和y′分别表示实测值和预测值，n为样本数量，i=1，2，···，n为计数变量，计算公式如下
3.1 数据的基本统计特征
研究数据的描述性统计见表1，可见监测点每周ILI 平均就诊数为925.94，变化范围为134—2823。研究时段深圳的周平均气温、相对湿度、气温日较差、风速的平均值分别是23.35℃、73.39%、6.05℃和2.15 m/s，总体上温、湿度高，风速小。
表1 深圳市2003—2019 年每周流感样病例及气象变量的统计特征Table 1 Statistics of weekly ILI cases and meteorological variables in Shenzhen from 2003 to 2019
由图1 可见，深圳ILI 发病序列具有明显的年际周期，叠加阶段性趋势；
与发病的年周期相对应，大气温、湿度变量的年周期也十分明显，且气候的暖湿特征凸显。气温日较差的变化无序，但整体与气温呈反位相分布，即暖季气温的日内变化更小；
风速总体偏小，并出现缓慢下降趋势。
此外，多数年份ILI 发病人数呈单峰型，个别年份出现多个就诊高峰（如2009 年）。2003 年之后流感发病渐次增加，直至2009 年全球性甲型H1N1 疫情暴发；
此后至2014 年底发病较为平稳，未出现大规模疫情；
而在2015—2016 年全球再次暴发高致病性禽流感（含H5N1、H5N2、H7N9 等亚型）；
之后发病呈显著下降趋势（图1）。考虑到深圳的人口基数仍在增加，说明相关部门对流感的防控见效，市民的预防意识也逐渐提高。
图1 2003—2019 年深圳市每周ILI 就诊人次（a）与气象要素（b—f）的时间序列（a.红色虚线为拟合的ILI 发病年周期，b.平均气温，c.气温日较差，d.相对湿度，e.绝对湿度，f.风速）Fig.1 Time series of weekly ILI cases（a）and meteorological factors（b—f）in Shenzhen from 2003 to 2019（a.the red dashed line indicates the fitted annual cycles of ILI cases，b.temperature，c.daily temperature range，d.relative humidity，e.absolute humidity，f.wind speed）
考察ILI 与ILI%的年内分布，二者均在6 月最高、7 月次之，而11 月—翌年3 月整体较低（图2）。由于当地气候类型较特殊，参照深圳市气象局标准，以5 d 滑动平均气温稳定＞10℃入春、稳定＞22℃入夏为依据，划分以后发现深圳的四季极不均匀（表2），夏季长达196 d，对应ILI 和ILI%也为全年最高；
冬季仅24 d，秋、冬季ILI 发病率偏低。逯建华等（2009）对深圳ILI 疫情的分析显示，3 月以乙型流感为主（72.7%），5、6 月则由甲型流感主导（51.9%）。另有研究发现中国甲型流感的活动与纬度关系密切，33°N 以北的省份呈冬季单高峰模式；
27°N 以南则在4—6 月出现发病高峰；
二者之间的区域则可能同时出现上述双周期（Dai，et al，2018；
Zhang，et al，2020）。而乙型流感呈现较统一的冬、春季流行模式（Yu，et al，2013）。
图2 深圳市2003—2019 年流感样病例数及就诊率的逐月分布Fig.2 Monthly distribution of weekly ILI and ILI% in Shenzhen from 2003 to 2019
表2 深圳的四季划分及其流感样病例特征Table 2 The division of seasons and corresponding ILI characteristics in Shenzhen
对ILI 与气候因子做双变量相关分析（表3），发现发病人数与气温、湿度、降水量呈显著正相关，夏季和过渡季节的相关一致。ILI 与气温日较差和风速总体为显著的负相关，尤以春、秋季最明显。冬季的多个相关系数未通过显著性检验，这与其持续时间短且发病率低均有关（表2、3）。
表3 周ILI 人次与气候因子分季节的斯皮尔曼相关Table 3 Spearman correlations between ILI outpatients and meteorological variables in different seasons
3.2 ILI-气象条件的非线性关联
基于3.1 节的相关分析结果，在时间序列建模中重点考虑气温、相对湿度、气温日较差、风速与流感风险的关联，图3 为相关因子频次分布特征。由于水汽压、比湿、降水等因子与相对湿度共线程度较高，不宜同时纳入模型。
图3 深圳市主要气象要素的频次分布特征（a.平均气温，b.相对湿度，c.风速，d.气温日较差）Fig.3 Frequency distributions of main meteorological factors in Shenzhen（a.air temperature，b.relative humidity，c.wind speed，d.daily temperature range）
图4 为DLNM 拟合的气象因子在滞后0—3 周与ILI 风险的总体关联。发现气温在22℃上下的累积风险最低，相对危险度随气温升高呈指数型增大；
低温的影响相对较弱，但仍然显著。图5 进一步揭示气温的影响在滞后0—3 周的变化，可见高温的效应即时、影响大，而低温的影响滞后1—3 周，这与呼吸系统疾病普遍对气温的响应规律一致（马盼等，2018）。深圳的一项早期研究也支持上述发现，即满足周最高气温30℃、最低气温25℃以上时，极易出现流感发病高峰（翟红楠等，2009b）。文中结果与香港甲型流感发病的最低气温（17.95—24.10℃）（郭貔等，2011）相近。此外，温带地区的流感一般被认为与低温关系密切，相关结果在德国与西班牙（Soebiyanto，et al，2015）、苏格兰（Price，et al，2019）、美国（Barreca，et al，2012）、中国（周艳丽等，2021；
王伟，2016）的典型城市均有报道。
图5 气温对ILI 风险的影响随着滞后时长（a—d.滞后0—3 周）的变化Fig.5 Effects of weekly mean temperature on the risk of ILI at different lag periods（a—d.lag 0—3 week）
相对湿度在70%—75%附近累积风险最高，偏低和过高的相对湿度（＞85%）均可抑制深圳ILI 风险（图4）。Zhang 等（2020）发现上海的低湿和高湿均有助于流感传播，但具体关联性与病毒亚型有关；
也有学者推测低纬度地区湿度与流感发病的正向关联可能是间接的，如降水导致的室内拥挤效应（Soebiyanto，et al，2014）。研究（Lowen，et al，2008；
Alonso，et al，2007）还发现气溶胶的传播效率与湿度呈负相关，因此热带地区更可能通过直接接触传播。此外，温带地区的流感传播与低湿关系更密切，因为干冷天气会增加宿主的易感性，还可提高空气中流感病毒的稳定性（Cheng，et al，2016），而沿海地区常年湿度较高，这类效应不明显。
图4 气象要素（a.平均气温，b.相对湿度，c.风速，d.气温日较差）对ILI 的累积滞后（3 周）危险度及其95%置信区间Fig.4 Effects of meteorological variables（a.air temperature，b.relative humidity，c.wind speed，d.daily temperature range）and their cumulative RR（95%CI）on the ILI at lag of 0—3 weeks
气温与湿度的交互影响通过构造二元响应模型来定性呈现（图6a）。发现ILI 风险最高点出现在气温极高、湿度偏高的情况下，气温为主导因子、湿度协同影响；
ILI 风险最低点对应气温极低而湿度极高的情况。结合图3 及深圳的气候背景，高风险点正好对应高发频率段，因此需要进一步提高预防意识。
不同风速在0—3 周内的效应总体不显著（特别小风速除外，图4c，图7a）。相关结果一致性也较低，山东周平均风速增大不利于流感发生（王伟，2016），而上海的流感风险随风速增加显著增长（Zhang，et al，2020）。一方面，由于风速增大，空气流动加速，不利于病毒颗粒聚集下沉；
但对于适应了暖湿气候的南方人而言，天气突变、风速增大也可能降低人体抵抗力、更易感染呼吸道病毒。
气温日较差具有显著的非线性影响，7—8℃和4℃以下的气温日较差对应ILI 低风险点，过大的日温差会显著增加深圳ILI 风险（图4d）。研究（Zhang，et al，2020）发现气温日较差增大主要影响甲型流感发病，而乙型流感更易在气温日较差较小时传播，一定程度上可以解释本研究的多个阈值点。再者，气温与其日较差的交互效应突出高温的主导性影响（图6b），气温日较差为协变量。最后，气温日较差主要在就诊当周有显著影响，1 周以后作用减弱，影响相对短促（图7b）。
图6 气温分别与（a）相对湿度和（b）气温日较差对ILI 风险交叉效应的三维曲面Fig.6 Three-dimensional curves of cross-effects（a）between temperature and RH，and（b）between temperature and diurnal temperature range
图7 风速（a）和气温日较差（b）在滞后0（a1—b1）和1（a2—b2）周时对ILI 相对危险度的影响Fig.7 Effects of wind speed（a）and daily temperature range（b）on ILI RR at lag of 0（a1—b1）and 1（a2—b2）week
各气象条件的典型值与ILI 发病的滞后性关联见表4。气温、相对湿度、风速和气温日较差的参考值（即RR=1）分别为22℃、34%、2 m/s 和2℃。选取气温5%和95%分位数代表低温、高温效应，发现29.9℃时相对风险高达1.237（95%CI：1.203—1.272）。湿度选取70%为高风险代表值，95%分位数（91%）作为对比，发现适中的湿度更利于ILI 发病。风速分别取1 和3 m/s 代表其极小和极大值，发现风速偏大（3 m/s）会增加ILI 风险，并可滞后1 周。气温日较差适中（5℃）或较大（10℃）均有显著的即时影响，相对危险度量级相当。
表4 各气象要素的典型值对ILI 的单周滞后相对危险度及其95%置信区间Table 4 RR and 95%CI of individual meteorological variables on ILI risk at various lag times
采用Prophet 模型和多元逐步回归分别构建ILI 风险预测模型，用2003—2018 年的资料建模，以2019 年数据做检验，分别进行回报和试预报优度检验（表5）。两种模型的回报准确率均高于86%，暖季（夏季）的回归模型准确率最高；
2019 年Prophet模型的预测准确率接近70%，但低于全年及分季的多元回归模型（约80%），主要源于后者不仅纳入外部的气象因子，还考虑了前一日的ILI 发病人数（表6）。
表5 ILI 预报模型优度检验Table 5 Testing indices for the predictive model of ILI risk
表6 深圳市ILI 风险的逐步回归模型Table 6 The stepwise-regression equations of ILI risk in Shenzhen
图8 显示两种模型对2003—2018 年ILI 发病的年周期和阶段性趋势模拟均较准确，其中回归模型对发病峰值的反映更准确。2019 年Prophet 模型预报曲线在上半年吻合度较好，下半年预测值较实际值整体偏高，可能与多种人为或自然混杂因素的影响有关。据调研，在《中国流感疫苗预防接种技术指南（2019—2020）》发布前后，深圳开始实行每年为中小学在校学生和60 岁以上老人免费接种流感疫苗，这对2019 年发病数的明显下降可能有重要贡献。
图8 基于（a）多元逐步回归和（b）Prophet 时间序列模型的ILI 风险预测效果Fig.8 Predicted results of ILI risk based on（a）the stepwise-regression equations and（b）Prophet time-series model
基于深圳市连续17 a 流感样病例（ILI）的监测资料，通过引入国际上通用的生物-医学时间序列模型，系统性分析了主要气候因子对ILI 发病的影响及其滞后性，并采用2 类方法对ILI 风险进行预测建模研究，主要结论如下：
（1）深圳ILI 就诊人数对流感大流行的响应显著，相应时段会出现就诊人数峰值；
近几年ILI 发病人数呈下降趋势。
（2）夏季ILI 就诊高峰多年稳定维持，6—7 月风险最高；
少数年份在年末出现流感次高峰。
（3）平均气温在22℃上下时，ILI 发病风险最低，高温的影响较强、即时出现，而低温的滞后效应显著；
湿度在70%—75%附近风险最高。
（4）气温日较差和湿度均会与气温产生交叉协同影响，但气温起主导作用；
高温、高湿条件下ILI 风险达到极大值。
（5）深圳市ILI 发病风险具有一定可预测性，Prophet 时间序列模型与多元回归的预测准确率为68% —80%，将前一日发病人数纳入回归模型有助于提高预测准确率。
除了气象因子的多重关联性外，人为因素的影响不可忽视，每年接种疫苗仍是预防流感传播最有效的措施。通过近几年多项公共卫生政策的实施，深圳的流感疫苗接种率在逐步提升（刘刚等，2021；
彭质斌等，2018）。由于深圳气候暖湿、空调使用十分普遍，造成室内空气流通性差，商场、学校、工厂等人员密集处容易发生暴发性疫情；
再者，室内外温差较大也会降低人体抵抗力，导致呼吸道病毒易乘机而入。此外，未来气候变化对相关呼吸道传染病的潜在影响值得关注，将成为下一步研究的重点。
需指明的是，流感样病例并不代表确诊的流感病例，其监测是为了实时监测流感活动水平和流行趋势，帮助发现新型流感病毒。而流感样病例标本的采集和实验室检测均有严格规定，本研究不包括确诊病例和流感病毒毒株的信息，日后还需区分病毒类型（如甲、乙型流感）进行更深入的分析。由于季节性流感的流行模式是由流感病毒、气候条件、社会经济水平和人类活动的复杂相互作用驱动（Tamerius，et al，2019；
Lofgren，et al，2007），提高其预测准确率的难度较高，需加强相关学科的合作研究。