时间序列模型对细菌耐药率变化趋势预测价值

钱鑫，熊世娟，杨廷秀，胡方芳，谢娟

贵州省人民医院，贵州 贵阳 550002

随着抗菌药物的广泛暴露，病原微生物的不断进化，细菌耐药已成为全球公共健康领域的重大挑战。近年来，随着碳青霉烯类耐药肺炎克雷伯菌（carbapenem-resistant klebsiella pneumoniae，CRKP）检出率呈现明显上升趋势，已严重威胁人类健康，给临床治疗带来巨大困难与挑战。在大数据时代，如何利用所收集的数据为决策提供服务已经逐渐成为这个时代的主题。美国学者 Box 和英国统计学专家 Jenkins 于1976 年提出的自回归移动平均模型（autoregressive moving average model，ARMA）［1］，目前在临床流行病学方面应用较为广泛。本研究旨在探讨ARMA模型对CRKP耐药率趋势的预测价值，为今后医院感染的预防和控制，以及进一步加强抗菌药物管理提供相应的数据支持和理论参考。

1 资料与方法

1.1 研究资料

本次研究的资料来自贵州省人民医院在2013年1 月至2018 年12 月送检的病原学检测结果，去除了同一个患者分离的重复菌株。标本均由各病区护理人员按标准采集，采集流程按照《全国临床检验操作规程》［2］执行。采集的标本均按时送至微生物室，并行标本质量的人工复检，合格的标本在初步培养的基础上，采用法国梅里埃公司 VITEK 2 Compact 全自动微生物鉴定系统进行细菌鉴定及药敏试验，所分离的CRKP 菌株根据美国临床与实验室标准协会（clinical and laboratory standards institute，CLSI）标准对药敏结果进行判定。

1.2 模型介绍

ARMA 模型是时间序列预测分析方法之一［3］，一般常用经差分的ARIMA（p，d，q）模型。AR 是自回归模型，p 为自回归系数；MA 为移动平均模型，q为移动平均系数，d 为成为平稳序列所需差分的阶数。

1.3 模型建立

建模过程主要分3 个步骤进行：（1）确定模型的参数，时间序列模型分析需要保证时间序列的平稳性。如果时间序列表现得不够平稳，则通过差分法可以使其均数与方差趋于平稳化，之后观察所作自相关函数图（autocorrelation function，ACF）和偏自相关函数图（partial correlation function，PACF）的截尾、拖尾情况，初步进行定阶，按低阶到高阶的顺序，提供相关的粗放性模型进行筛选。（2）检验模型的参数，在模型拟合的过程中，如果存在多个模型均通过Box-Ljung 检验，则提示其残差序列均为白噪声序列［4］，根据赤池信息准则（akaike information criterion，AIC）［5］、贝叶斯信息准则（bayesian information criterion，BIC）［6］及许瓦兹贝叶斯准则（schwarz bayesianc riterion，SBC）［7］来选择最优模型。若模型参数无统计学意义，那么分析流程则中止，一般需要返回模型识别阶段重新定阶，再进行模型参数的估计和诊断。（3）评估模型预测价值：比较拟选模型的预测值与实际值，对拟选模型间预测的准确性进行评价。

1.4 数据分析

以月为单位，将耐药率数据录入Excel 2016中，分别导入IBM SPSS Statistics 19.0 和SAS Statistical Discovery 13，将需要分析的原始数据进行平稳化处理后，采用时间序列模型进行拟合分析，以最终确定的模型来预测2019 年的细菌耐药率，并采用绝对误差和相对误差对预测值进行评价，时间序列的平稳性运用EViews10.0，通过扩展迪基-福勒检验（augmented dickey-fuller test，ADF）［8］对单位根进行检验。

2 结果

2.1 全球耐药趋势及医院CRKP 耐药形势

据笔者调研，美国非营利组织“动态疾病、经济和政策研究中心”（center for disease dynamics,economics &policy，CDDEP）对科学文献以及国家和地区监测系统的分析数据提示，中国肺炎克雷伯菌的碳青霉烯类耐药形势已较为严峻，见图1～3。结合院细菌耐药监测情况，近6 年院CRKP 耐药率呈现持续上升趋势，2018 年院美罗培南耐药率为22.62%，稍低于全国年度平均水平28.60%，但该菌2018 年在院的碳青霉烯耐药率较之2017 年出现了大幅度提高，引起了临床的高度关注，特别是碳青霉烯类使用频率较高的重点科室，如重症监护室、血液科等。MDRO千日检出率近6 年来总体变化趋势不大，虽然MDRPA 和MRSA 的千日检出率亦相对稳定，2018 年分别为0.15‰和0.31‰，但CRE 千日检出率出现大幅度抬头趋势，临床面临着较大的压力与挑战。

图1 肺炎克雷伯菌对碳青霉烯类药物的耐药性

图2 肺炎克雷伯菌的细菌耐药率

2.2 模型处理与检验

将耐药率原始数据作时序图，可观察到时间序列的波动比较大，并不能满足时间序列分析的平稳化要求，见图4。该时间序列经ADF 检验，P=0.482，P＞0.05，提示序列存在单位根。对该序列进行一阶差分后序列图显示基本达到平稳，经ADF检验P＜0.05，提示序列无单位根，符合序列平稳化要求，见图5。通过谱密度查看周期性变化，显示序列的季节性波动不明显，因此该时间序列可选择ARIMA（p，1，q）模型，初步筛选ARIMA（1，1，1）、IMA（1，1）、ARI（1，1）三个时间序列模型。初选模型均通过了Ljung-Box Q 检验，提示上述模型观测值是相互独立的，残差是白噪声序列，各备选模型参数信息，见表1。通过方差、AIC、SBC、R 方及均方根误差（root mean square error，RMSE）和平均绝对误差（mean absolute error，MAE）评估模型拟合优度，显示ARIMA（1，1，1）模型权重较高，其Ljung-Box Q 统计量Q=21.348，P=0.166，模型参数估计结果P＞0.05，其残差ACF 和PACF 图，见图6。

图3 MDRO千日检出率变化

图4 原始数据时序图

图5 一阶差分时序图

图6 模型ARIMA（1，1，1）的残差ACF和PACF图

2.3 模型的预测

尽管ARIMA（1，1，1）模型参数权重更高，但与IMA（1，1）相比残差无明显差异，见图7。选取2018 年1 至12 月的真实值分别与两个模型的预测值进行比较，计算模型预测值的相对误差和平均相对误差。结果显示，IMA（1,1）模型的平均相对误差相对较小，IMA（1，1）的拟合效果更接近于真实值，与可视化观察结果一致。基于上述的分析流程和推断，IMA（1，1）模型的拟合度最优。

图7 ARIMA（1，1，1）和IMA（1，1）模型的拟合曲线

3 讨论

近年来，多重耐药菌导致的感染是造成医疗保健相关发病率、死亡率和医疗成本支出的主要来源［9-11］。特别是肺炎克雷伯菌作为肠杆科兼性厌氧菌，是医院感染的重要病原菌［12］。中国的细菌耐药监测网CHINET 显示CRKP 检出率逐年上升［13］。细菌耐药问题已经被G20 峰会和联合国大会正式纳入政治议题。2017 年3 月国卫办医发（2017）10号文出台，进一步强化了碳青霉烯类抗菌药物的管理。因此，制定有效的战略可以对这类病原菌实施监测和预测，同时可以有助于采取有效地干预措施。本研究数学模型的识别方式，将临床收集到以时间为序的细菌耐药数据进行未来值探索与预测，可为我们的公共卫生部门提供重要的决策信息。目前，时间序列模型已广泛应用于患病率的预测，但基于细菌耐药率的预测研究却不多。通过本研究可以看到，贵州省人民医院CRKP 检出率近年大幅增高趋势，预测提示存在持续走高的风险。结合卫健委发布的《碳青霉烯类抗菌药物临床应用评价细则》以及医院实际应用情况分析，需要密切关注以下问题：（1）手卫生不到位使耐药菌株病区内或病区间传播，必然导致碳青霉烯类用量增加，出现耐药菌筛选压力；（2）部分临床医生对经验治疗把握不到位，导致疗程过长，或者对定植菌未仔细评估；（3）不少医务人员存在重症就是耐药菌感染的错误认识，以及碳青霉烯类使用频次不足等问题；（4）高级职称的临床医生对抗菌药物培训参与度不够，持续学习不足，对碳青霉烯类药物应用指征了解不充分；（5）医院抗菌药物TDM 监测手段不完善，对感染性疾病的疗效评估缺失个体化评估手段，一定程度上影响了临床用药的决策。

表1 备选模型的估计参数

表2 2018年1至12月耐药率观测值与预测值的比较

本研究通过原始数据的参数估计和检验，获得了模型IMA（1，1），应用该模型预测的结果显示，拟合值的动态趋势与真实值的升降规律较为相似，提示时间序列模型可较好的分析和预测细菌耐药率的变化趋势。不过，该模型虽然具有不受数据类型束缚、短期预测良好的特点，但考虑到细菌耐药规律的复杂性，模型的长期应用需要保证有足够的数据，并及时补充新的数据，同时对模型类型和参数不断完善修正，才可能确保预测值的精准度，便于“早发现、早干预”，为进一步加强抗菌药物的管理，有效控制细菌耐药作出信息预警。