長沙市金屬和類金屬職業病危害因素檢測超標率時間序列分析*

何 俊 胡平成 金若剛 羅 磊 李繼猛

?

長沙市金屬和類金屬職業病危害因素檢測超標率時間序列分析*

何 俊1胡平成2金若剛1羅 磊1李繼猛1

目的 了解長沙市金屬和類金屬職業病危害因素檢測超標率變化趨勢，建立季節性ARIMA模型進行短期預測，為控制工作場所金屬和類金屬危害提供措施建議。方法 利用2008-2015年檢測數據為訓練集進行建模，利用2016年上半年數據為驗證集檢驗，并進一步預測2016-2017年的超標率。結果 金屬和類金屬超標率長期趨勢為逐年下降，且具有一定的季節效應，其中二、三季度較高，一、四季度較低。建立ARIMA(0,1,1)×(0,1,1)4模型，經檢驗殘差為白噪聲序列，預測的平均絕對誤差(MAE)為4.54%。結論 長沙市金屬與類金屬危害的控制近年來取得了一定的成效，超標率呈逐年下降趨勢。2017年第二季度是金屬與類金屬危害控制的重點。季節性ARIMA模型在金屬與類金屬檢測超標率預測方面有一定的可行性，今后可以和其他模型聯合降低預測誤差。

金屬與類金屬 職業病危害因素 超標率 ARIMA模型

長沙市工業企業以機械制造、金屬加工、汽車制造、電氣電子、化工等產業為主導，存在錳、鉛、錫、鋅、銅、鎘、鉻、磷、砷等多種金屬和類金屬職業病危害因素。勞動者在生產過程中可能因意外事故、燃燒或焊接吸入高濃度金屬煙霧導致急性中毒，低劑量長時間接觸金屬和類金屬也可引起慢性中毒[1]。預防和控制金屬和類金屬職業病危害因素對保護勞動者身體健康至關重要。

自回歸滑動平均混合模型簡稱ARIMA模型，目前已在公共衛生領域得到廣泛應用，如流感、結核、乙肝、猩紅熱、甲肝、呼吸道傳染病、手足口病、梅毒的發病預測[2-9]，住院人數、門診量、門診收入變化趨勢的季節調整等方面[10-12]，但在職業衛生領域的應用還鮮有報道。本文通過分析2008-2016年上半年長沙市工業企業金屬和類金屬職業病危害因素的檢測結果，利用季度超標率建立季節性ARIMA模型，對未來金屬和類金屬危害的變化趨勢進行預測，以制定降低工作場所有害物質濃度和保護勞動者健康的策略與措施。

資料與方法

1.資料來源

本研究資料來自2008-2016年上半年長沙市疾病預防控制中心依法對長沙市用人單位每年一次的定期檢測，這里對用人單位的定義與《中華人民共和國職業病防治法》的闡述保持一致，即存在粉塵、放射性物質和其他有毒、有害因素的企業、事業單位和個體經濟組織。長沙市疾病預防控制中心定期檢測的用人單位主要來自職業衛生服務市場化前的主管企業，分布在全市各行業和地區，可基本代表長沙市存在職業病危害因素用人單位的危害狀況。

2.分析方法

ARIMA模型是Box-Jenkins方法中成熟的時間序列分析預測方法，模型主要分析時間序列的隨機性、平穩性和季節性，根據Akaike′s information Criterion(AIC)和Schwarz′s Bayesian Criterion(BIC)選擇合適的參數進行預測，并通過Box-Ljung檢驗、平均絕對誤差(MAE)、平均相對誤差(MAPE)對模型進行評價。ARIMA建模主要包括三個內容：自回歸、滑動平均、差分求和，季節性模型一般表示為ARIMA(p,d,q)×(P,D,Q)s，其中p、P表示模型的普通、季節自回歸階數，d、D分別表示普通、季節平穩化的差分階數，q、Q表示普通、季節移動平均階數[13]。模型公式為：

φ(B)=1-φ1B-φ2B2-…-φpBp

θ(B)=1+θ1B+θ2B2+…+θqBq

ΦP(Bs)=1-Φ1Bs-Φ2B2s-…-ΦPBPs

ΘQ(Bs)=1+Θ1Bs+Θ2B2s+…+ΘQBQs

▽d=(1-B)d

本研究采用excel 2016進行數據錄入匯總，利用R(3.3.1)進行ARIMA建模、檢驗、預測和繪圖。

結 果

1.金屬和類金屬超標率序列特征

本研究數據為2008-2016年上半年金屬和類金屬檢測數據，選擇2008-2015年數據為訓練集，2016年上半年數據為驗證集進行建模。對訓練集數據用stl函數進行分解，圖1展現了原始序列、季節性、長期趨勢、不規則變化四個部分，序列長期趨勢明顯，基本呈下降趨勢，且具有一定的季節效應，高峰位于每年的第二、三季度，而第一、四季度較低。

圖1 超標率時間序列分解圖

2.季節性ARIMA模型的建立

(1)序列平穩化

金屬和類金屬超標率時間序列有明顯的長期趨勢，提示可能為非平穩序列，由于原始序列部分季度數據為0，為降低序列波動，采用加1后自然對數變換的方式進行統計變換。對新序列進行普通差分并繪制ACF(自相關系數)圖和PACF(偏相關系數)圖，其中ACF圖顯示序列存在一定的周期性，第1、5、9、13階相關系數顯著增大，對序列再進行季節性差分，經ADF單位根檢驗，統計量為-6.23，P=0.01<0.05，說明兩次差分后的序列無單位根，為平穩序列，見圖2。

(2)模型的識別與建立

圖2兩次差分后的序列仍存在一定的周期性，其中ACF圖第1、5、13階相關系數增大，可視為季節性1階后截尾或拖尾，PACF圖第1、4、5、8、9、12、13、16階顯著增大，可視為季節性1階后截尾、3階后截尾或拖尾，則P、Q有(0,1)、(1,0)、(1,1)、(3,0)、(3,1)五種可能。

圖2 兩次差分的ACF圖和PACF圖

周期內ACF圖可看做1階后截尾或拖尾，PACF圖可看做拖尾，則p、q有(0,1)、(1,1)兩種可能。利用P、Q、p、q可能的參數組合分別建模，由于是差分序列，各模型均未包含常數項，計算得AIC和BIC值見表1。其中ARIMA(0,1,1)×(0,1,1)4的AIC和BIC值均最小，t檢驗表明模型參數不為零，ARIMA(0,1,1)×(0,1,1)4為最優模型，表達式為：(1-B4)(1-B)xt=(1-0.68B4)(1-0.95B)wt，各參數見表2。

表1 模型的AIC、BIC值

表2 ARUNA(0,1,1)×(0,1,1)4模型參數

(3)模型的診斷與驗證

對ARIMA(0,1,1)×(0,1,1)4模型回歸殘差進行Ljung-Box檢驗，滯后3階到16階的自相關系數相應Ljung-Box檢驗統計量的概率值都大于0.05，可以認為模型序列不存在自相關，為白噪聲序列，見圖3。

圖3 Ljung-Box檢驗結果

為驗證ARIMA(0,1,1)×(0,1,1)4模型的效果，利用模型對長沙市2016年第一、二季度金屬類職業病危害因素檢測超標率預測，預測結果分別為1.48%、4.89%，實際超標率為0、12.5%，所選模型的平均絕對誤差(MAE)為4.54%。

(4)預測

將2016年上半年實際數據加入ARIMA(0,1,1)×(0,1,1)4重新建模，經Box-Ljung檢驗，3-16階P值均大于0.05，殘差為白噪聲序列。利用新模型對2016-2017年長沙市金屬類職業病危害因素檢測超標率進行預測，6個季度超標率分別為2.63%、1.95%、0、5.89%、0.93%、0.26%，超標率峰值在第二季度，見圖4。

圖4 長沙市金屬類職業病危害因素檢測超標率預測

討 論

金屬和類金屬是長沙市用人單位存在的主要職業病危害因素之一，控制其工作場所濃度和降低超標率是急慢性中毒類職業病防治的重中之重[15]。本研究采用季節性ARIMA模型對金屬和類金屬超標率歷史數據進行了分析，對短期內超標率進行了預測，旨在發現近年來長沙市金屬和類金屬職業病危害變化的特點和規律，從而制定有力的控制措施。

時間序列的長期趨勢顯示檢測超標率呈逐年下降趨勢，表明長沙市職業病防治工作近年來取得了一定的成效。時間序列同時也存在一定的季節性，二、三季度較一、四季度有更高的超標率，一方面這可能與某些毒物在氣溫較高時容易揮發和擴散有關，另一方面長沙市疾病預防控中心檢測行業的季度分布并不均勻，機械加工和電氣電子等金屬和類金屬工藝較多的用人單位在二、三季度檢測較多，從而導致金屬類職業病危害超標率更高[15]。

本研究利用長沙市金屬與類金屬職業病危害因素的季度超標率建立ARIMA(0,1,1)×(0,1,1)4模型，對未來超標率進行了短期預測。模型2016年上半年超標率的驗證表明預測值與實際值的發展趨勢一致，通過診斷和檢驗說明該模型對長沙市金屬與類金屬超標率的預測具有一定的可行性，但2016年第二季度還存在較大預測誤差，在以后的研究中可試驗ARIMA與其他模型聯合提高預測能力。

長沙市金屬與類金屬類職業病危害因素有較明顯的時間趨勢和季節趨勢，根據序列的規律對未來超標率進行預測，可預先了解金屬與類金屬類職業病危害情況的發展狀況，以提前發現高風險企業，采取適當的控制措施降低金屬與類金屬對勞動者身體健康的損害。預測顯示2017年第二季度金屬與類金屬的危害較大，提示我們在工作中尤其要關注此時期內檢測企業的金屬與類金屬危害，發現問題需及時向安全生產監管等政府部門進行反饋，密切注意工作場所金屬與類金屬濃度的變化，從而防患于未然，降低職業性急性金屬與類金屬中毒的發病率。模型以歷史超標率數據建立，在今后的預測中還需加入每年的新數據進行修正，提高金屬與類金屬超標率預測的準確性[9]。

[1]孫貴范.職業衛生與職業醫學.第7版.北京:人民衛生出版社,2012.

[2]Song X,Xiao J,Deng J,et al.Time series analysis of influenza incidence in Chinese provinces from 2004 to 2011.Medicine,2016,95(26):e3929.

[3]Li XX,Wang LX,Zhang H,et al.Seasonal variations in notification of active tuberculosis cases in China,2005-2012.Plos One,2013,8(7):e68102.

[4]陳正利,陳偉,許汴利.應用ARIMA模型對河南省1991-2011年乙型肝炎發病趨勢分析.中國衛生統計,2013,30(3):401-402.

[5]朱玉.單純ARIMA模型和ARIMA-GRNN組合模型在猩紅熱月發病率中的預測效果比較.合肥:安徽醫科大學,2011.

[6]朱奕奕,馮瑋,趙琦,等.ARIMA乘積季節模型在上海市甲肝發病預測中的應用.復旦學報(醫學版),2012,39(5):460-464.

[7]時照華.探索ARIMA模型在呼吸道傳染病疫情預測中的應用.合肥:安徽醫科大學,2013.

[8]董晶.長沙市手足口病疫情分布特征及預測模型的建立與評估.長沙:中南大學,2014.

[9]陳偉,陳正利,李少芳,等.ARIMA模型在河南省梅毒月發病率預測中的應用.中國衛生統計,2013,30(4):604-606.

[10]薛允蓮,張晉昕,劉貴浩,等.醫院住院人數序列的春節效應調整.中國衛生統計,2010,27(5):473-476.

[11]王玖,韓春蕾,欒奕昭.組合預測在醫院門診量預測中的應用.中國衛生統計,2012,29(6):881-883+886.

[12]郭志武,蒲繼紅,滕國召.基于ARIMA模型的春節因素調整方法研究.中國衛生統計,2009,26(6):573-576+579.

[13]孫振球,徐勇勇.醫學統計學.第4版.北京:人民衛生出版社,2014.

[14]Shumway RG,Stoffer DS.Time series analysis and its applications with R examples.Third Edition.Springer,2010:83-161.

[15]何俊.2008-2013年長沙市工業企業職業病危害因素縱向研究.長沙:中南大學,2014.

(責任編輯：郭海強)

2016年湖南省衛生計生委科研項目(B2016231)

1.長沙市疾病預防控制中心職業病防治所(410004)

2.中南大學湘雅公共衛生學院