基于R语言ARIMA模型在慢阻肺急性加重患者发病预测中的应用

成都市第三人民医院信息部(610031)

郭慧敏 杜 军△ 黄路非

基于R语言ARIMA模型在慢阻肺急性加重患者发病预测中的应用

成都市第三人民医院信息部(610031)

郭慧敏 杜 军△黄路非

目的 建立慢阻肺急性加重入院人次的自回归积分滑动平均模型(ARIMA)，科学预测慢阻肺急性加重入院人次，为该病的诊治以及合理利用医疗资源提供理论依据。方法 使用R语言(v.3.2.3)做模型的识别、模型的参数估计与检验，建立ARIMA模型，对某院2013-2015年慢阻肺急性加重出院人次进行模型拟合，用2016年1～3月的预测值与实际值作比较，检验模型的预测能力，并且预测2016年4-6月慢阻肺急性加重入院人次。结果 经过多次检验，确定ARIMA(2,2,1)(1,1,1)12模型预测能力最佳，其残差序列是白噪声。用2016年1～3月数据来检验模型，其MAPE的绝对值均小于10%，说明模型的拟合优度相对较好，预测能力可靠，根据该模型预测2016年4～6月该院慢阻肺急性加重入院人次分别为162、160、159。结论 ARIMA模型能够很好的拟合慢阻肺急性加重的入院人次并进行短期预测，模型显示2016年该院的急性支气管炎的入院人次将有所上升，为医院合理利用医疗资源提供了有力依据。

R语言 慢阻肺急性加重 预测 回归滑动平均混合模型

慢性阻塞性肺疾病(简称慢阻肺)是一种逐渐削弱患者呼吸功能的破坏性慢性肺部疾病，被称为呼吸道疾病中最“不动声色”的隐形杀手。患者在日常“稳定期”时疾病特征不明显，却在“急性加重期”时症状骤然出现或原有症状急剧恶化[1]，对患者的生活质量产生巨大影响， 加速肺功能恶化， 也是患者住院和死亡的主要原因。根据世界卫生组织估计，到2020年，慢阻肺将成为全球第三大死亡原因，预防慢阻肺急性加重的重要性正日益得到关注，最新版慢阻肺全球创议(GOLD指南)首次强调了慢阻肺急性加重风险和症状的评估[2]。本文建立某院慢阻肺急性加重的自回归积分滑动平均模型(autoregressive integrated moving average model，简称ARIMA)，探讨慢阻肺急性加重入院人次的变化规律。

R语言是一种为统计计算和绘图而生的语言和环境，它是一套开源的数据分析解决方案，由一个庞大且活跃的全球性研究型社区维护[3]，针对不同的场景和统计功能R语言的用户贡献了大量优秀的包(package)。本文采用R语言作为统计分析的工具。

资料与方法

1.资料来源

从某院病案首页信息管理系统中检索2013-2015年出院第一诊断为慢阻肺急性加重(ICD-10编码为J44.1)的病人为研究对象[4]，以每月的入院人次构成时间序列。

2.原始数据获取

R语言通过RODBC包连接数据库[5]，并且提取2013-2015年的研究数据，建立时间序列(stats包中的st函数)，以月为单位绘制原始数据的时间序列图(图1)。

图1 某院慢阻肺急性加重患者住院人次时间序列

3.ARIMA模型建立

ARIMA模型主要步骤[6]为：(1)序列平稳化：据图1，采用季节性ARIMA(p,d,q)(P,D,Q)12模型，对原始数据进行差分并检验序列平稳性，可以确定d、D；(2)模型识别：对处理后的序列做自相关和偏相关处理，发现低阶的自相关系数较大，但随着时滞长度的增加而衰减，可认为处理后的序列符合运用ARIMA模型的条件，并且确定p值和q值；(3)模型参数估计与检验：通过Ljung-Box检验等方法检测模型残差，判断模型的适合性；(4)预测：利用该院2016年1～3月的慢阻肺急性加重的实际入院人次，检验建立模型的预测效果，并预测2016年4～6月慢阻肺急性加重的入院人次，产生可信区间。

4.编程与实现

用R语言(v.3.2.3)作为统计分析的工具，用RODBC中的函数提取原始数据，采用stats包中的函数对原始数据进行时间序列处理，用funitRoot包中的函数进行单位根检验，用forecast包中的函数进行预测。

结 果

1.序列平稳化

对原始数据进行一阶差分和二阶差分(timeSeries包中的diff函数)，经单位根检验(fUnitRoot包中的unitrootTest函数)一阶差分为非稳态序列，二阶差分为稳态序列(P<0.05)，同时白噪声检验说明序列不是纯随机性检验，序列具有值得我们提取的相关信息。由此得到d=2。

2.模型识别

对二阶差分处理后的数据进行自相关和偏相关(stats包中的acf函数和pacf函数)处理，得到ACF图(图2)和PACF图(图3)。观察图2得，二阶自相关系数后都未超出±2倍估计标准差，即自相关系数1阶以后截尾，初步确定q=1；观察图3得，三阶偏相关系数后都未超过±2倍估计标准差，即偏相关系数2阶以后截尾，初步确定p=2。

图2 自相关系数图

图3 偏相关系数

3.模型的参数估计与检测

由以上两步我们大致可以确定p、d、q的大致范围，在p、d、q确定的范围，对模型反复调试和检验(stats包中的arima函数)，季节模型的参数采取0、1、2从低阶到高阶逐个尝试的办法，根据模型参数检验结果和参数间的相关系数对模型反复调试和检验，以赤池信息准则(Akaike Information Criterion，AIC准则)和决定系数R2作为依据确定最优模型(表1)。

根据比较发现模型ARIMA(2,2,1)(1,1,1)12的AIC=335.78最小并且R2=85.23最大，则我们可以确定模型ARIMA(2,2,1)(1,1,1)12拟合较好，对残差序列做自相关图(图4)；图4表明一阶残差序列的自相关系数截尾，据此认为残差序列是白噪声。同时做Ljung-Box检验(stats包中的Box.test函数)，其统计量无统计学意义(P>0.05)。

表1 备选ARIMA模型拟合优度统计量

图4 残差序列的自相关系数图

4.预测

利用ARIMA(2,2,1)(1,1,1)12模型对该院慢阻肺急性加重入院人次以及95%可信区间进行检验(forecast包中的forecast.Arima函数)，结果见表2；由表2得，预测值均在95%可信区间，而且其平均绝对百分误差(MAPE)也非常小，预测模型精度高。据此模型外推，可以预测2016年4～6月慢阻肺急性加重入院人次分别为162、160、159，结果见图5。

表2 2016年1～3月慢阻肺急性加重入院人次预测值与实际值

图5 慢阻肺急性加重入院人次拟合曲线

讨 论

ARIMA模型是由Box和Jenkins(1970)提出的一种时间序列的建模方法，其预测精确度较高，且不需要知道影响预测变量的相关因素,可将各种因素包括未知因素的综合效应统一蕴含在时间变量中， ARIMA模型在医疗卫生领域有广阔的应用前景[7]。

本文通过对某院的2013-2015年慢阻肺急性加重入院人次建立ARIMA模型，序列经过平稳化、模型识别、参数估计及检测等步骤发现ARIMA(2,2,1)(1,1,1)12模型可以有效的拟合每月的入院人次，并且做短期预测；通过对模型的检验，其月份的误差率较低，有一定的实用价值。通过预测可知，2016年4～6月该院的慢阻肺急性加重的入院人次为分别为162、160、159，而且95%可信区间均包含每月的预测值，同时我们对比历史数据发现，该病种的入院人次呈逐年上升的趋势，在配置急性支气管炎的医疗资源时我们应该参考该预测值，优化资源分布，合理利用有限的医疗资源。

数据显示，慢阻肺急性加重患者每年都有所增长，笔者认为由于我国经济不断发展，人们的生活环境污染加剧，该医院所处城市空气质量日渐恶化，呼吸系统疾病的发病率逐年上升，人们应积极的锻炼身体，养成良好的生活习惯，戒除香烟或自觉躲避二手烟、雾霾天等污染源，定期做针对性的医疗检查，是积极面对这种高发病率疾病的有效措施。政府应当对环境治理加大力度，提高国民身体素质。

[1]Faustini A,Stafoggia M,Colais P,et al.EpiAir Collaborative Group.Air pollution and multiple acute respiratory outcomes.Eur Respir J,2013,42(2):304-313.

[2]陈亚红,王辰.2015 年更新版GOLD 慢性阻塞性肺疾病诊断、治疗和预防的全球策略简介.《中国医学前沿杂志(电子版)》,2015,7(2):34-39.

[3]Robert I.Kabacoff著,高涛,肖楠,陈钢译.R语言实战.北京:人民邮电出版社,2013:4-5.

[4]董景五.疾病和有关健康问题的国际统计分类,第1版.人民卫生出版社,2008.

[5]孙振球,徐勇勇.医学统计学,第2版．北京:人民卫生出版社,2002:351-371.

[6]郭慧敏,杜军,练正秋.ARIMA模型应用于临床悬浮红细胞用量预测的探讨.中国输血杂志,2014,8(30):829-832.

[7]熊志斌.基于ARIMA与神经网络集成的GDP时间序列预测研究.数理统计与管理,2011,30(2):306-314.

(责任编辑：刘 壮)

△通信作者：杜军，E-mail:yixi_370785@163.com