李望晨，于貞杰，王在翔，張利平

（濰坊醫(yī)學(xué)院“健康山東”重大社會風(fēng)險預(yù)測與治理協(xié)同創(chuàng)新中心，山東濰坊261053）

醫(yī)療統(tǒng)計(jì)周期性預(yù)測問題的比較研究

李望晨，于貞杰，王在翔，張利平

（濰坊醫(yī)學(xué)院“健康山東”重大社會風(fēng)險預(yù)測與治理協(xié)同創(chuàng)新中心，山東濰坊261053）

文章立足醫(yī)療統(tǒng)計(jì)領(lǐng)域周期性預(yù)測問題，以門診人次為例進(jìn)行建模方法比較，為類似研究提供方法參考。采用方式為：方式1：門診人次季度數(shù)據(jù)由X11法時序分解后提取季節(jié)指數(shù)消除季節(jié)波動，用ARIMA法和拋物線擬合外推。方式2：周期差分提取季節(jié)信息，二階差分提取長期趨勢，直接用季節(jié)ARIMA法建模。季節(jié)波動和短期相關(guān)間交互影響不清，試用簡單季節(jié)或乘積季節(jié)模型。得出結(jié)論：方式1對周期時序資料分解后，由組合方法依次直觀反映；方式2直接用季節(jié)模型擬合預(yù)測，乘積季節(jié)模型性能并非有替代性，簡單季節(jié)模型更優(yōu)。算例為載體用幾類方法比較研究，對于醫(yī)院統(tǒng)計(jì)領(lǐng)域周期性預(yù)測問題有借鑒意義。

時間序列;組合模型；X11法；ARIMA；乘積季節(jié)模型

0　引言

醫(yī)療統(tǒng)計(jì)工作領(lǐng)域經(jīng)常見時間序列預(yù)測問題，如傳染病或慢性病發(fā)病數(shù)、發(fā)病率、門診或住院人次、醫(yī)療營業(yè)收入、衛(wèi)生費(fèi)用支出等。這些事物受到社會、經(jīng)濟(jì)、環(huán)境等復(fù)雜因素的影響，對其解讀分析也是衛(wèi)生政策制定者、公共衛(wèi)生干預(yù)者、衛(wèi)生行政管理者關(guān)心的問題。

由于影響因素復(fù)雜性，傳統(tǒng)回歸模型結(jié)構(gòu)模式確定、參數(shù)選擇困難，指標(biāo)難以量化。時序資料本身隨時間延伸表現(xiàn)規(guī)律性，適合于時間序列分析。時間序列方法在數(shù)學(xué)、統(tǒng)計(jì)學(xué)、計(jì)量經(jīng)濟(jì)和信息科學(xué)推動下，理論體系和應(yīng)用實(shí)現(xiàn)較為成熟[1,2]，近年來在醫(yī)療衛(wèi)生領(lǐng)域得到推廣應(yīng)用，如ARIMA法、曲線擬合法、灰色系統(tǒng)法、數(shù)據(jù)挖據(jù)技術(shù)，其中ARIMA法有代表性，對于隨機(jī)時序資料擬合較好，對于趨勢性、周期性特點(diǎn)資料也有很強(qiáng)適應(yīng)性。醫(yī)療統(tǒng)計(jì)領(lǐng)域常見周期性數(shù)據(jù)資料[3](如門診人次、收入、費(fèi)用)，多采集為月度或季度資料，觀察期長、數(shù)據(jù)豐富、周期明顯、長期平緩，這類事物規(guī)律性好于傳染病疫情時序資料。對此有必要引入適當(dāng)方法建立模型，借助算例實(shí)現(xiàn)為類似研究提供方法借鑒。

第一種方式，鑒于月度或季度周期變化時序資料周期性、趨勢性較明顯，可用時間序列分解法提取信息，X-11過程法可用于提取季節(jié)指數(shù)[4]，趨勢因素可用曲線擬合法或ARIMA法。時序分解后從兩種視角建立擬合模型，外推預(yù)測及合并信息后再預(yù)測未來情況。第二種方式，傳染病預(yù)測為當(dāng)前熱點(diǎn)問題，乘積ARIMA法對傳染病月度數(shù)據(jù)表現(xiàn)應(yīng)用代表性。對于季節(jié)性變化資料，首先用步長差分方法提取季節(jié)信息，低階差分提取長期趨勢信息，季節(jié)方法與ARIMA模型結(jié)合建模型，必要時比較簡單季節(jié)模型與乘積季節(jié)模型研究效果。

1　對象資料

以醫(yī)療統(tǒng)計(jì)領(lǐng)域時序資料為立足點(diǎn)建立模型和進(jìn)行實(shí)證研究，方法學(xué)視角下對于門診或住院人次、醫(yī)療收入等周期性指標(biāo)有普適借鑒意義。某二級甲等醫(yī)院2004—2014年按季度統(tǒng)計(jì)門診人次數(shù)據(jù)為例，以2004—2013年建立擬合模型，對2014年進(jìn)行預(yù)測驗(yàn)證。門診人次數(shù)據(jù)隨季度變化呈周期波動，年度間有平緩增長特點(diǎn)，見表1和圖1所示。

表1　某醫(yī)院門診人次數(shù)據(jù)

2　建模方式一

原始序列有明顯趨勢變化特點(diǎn)，對2004—2013年時序分解，提取季節(jié)指數(shù)后剩余序列有趨勢增長及隨機(jī)變化特點(diǎn)，用曲線擬合法或ARIMA法擬合趨勢規(guī)律并外推預(yù)測未來[5]，與季節(jié)指數(shù)相乘后還原為2014年情況。組合方法在SAS軟件環(huán)境下計(jì)算實(shí)現(xiàn)方便。

2.1季節(jié)指數(shù)提取

由于門診人次數(shù)據(jù)有明顯季節(jié)變化特點(diǎn)，可經(jīng)時間序列分解提取長期趨勢、周期性季節(jié)指數(shù)、不規(guī)則變化部分等。X-11法由美國國情調(diào)查局建立，由于需要用到11次移動平均而得名，可據(jù)此對門診人次數(shù)據(jù)提取季節(jié)指數(shù): 92.882%、102.616%、99.544%、104.963%，消除季節(jié)因素后長期趨勢見圖2、剩余隨機(jī)序列見圖3所示。

圖1　原始數(shù)據(jù)

圖2　長期趨勢

圖3　不規(guī)則波動

2.2ARIMA法

ARIMA表達(dá)式Φ(B)▽dxt=Θ(B)εt，其中φi為自回歸系數(shù),θj為移動平均系數(shù)。Φ(B)=(1-φ1B-…-φpBp)為自回歸多項(xiàng)式,Θ(B)=(1-θ1B-…-θqBq)為移動平均多項(xiàng)式, Bkxt=xt-k為k步延遲算子,▽d=(1-B)d為d階差分算子。

根據(jù)數(shù)據(jù)資料確定模型結(jié)構(gòu)、識別參數(shù)、檢驗(yàn)?zāi)Ｐ惋@著性，否則要不斷調(diào)整優(yōu)化模型。趨勢變化序列要經(jīng)低階差分后轉(zhuǎn)化為平穩(wěn)序列。由AIC、SBC最小準(zhǔn)則選擇最優(yōu)模型，由最小二乘法識別參數(shù)。由殘差純隨機(jī)性檢驗(yàn)后，由模型進(jìn)行短期外推預(yù)測。

正試期第2、第4周采用全糞收集方式采集糞樣，每組用蛇皮袋固定于隨機(jī)選取的4只試驗(yàn)羊尾后，每日更換兩次并稱重，將飼料樣與糞樣帶回實(shí)驗(yàn)室進(jìn)行常規(guī)分析計(jì)算養(yǎng)分表觀消化率。正試期第2、第4周每組隨機(jī)選取4只試驗(yàn)羊，晨飼后3 h利用負(fù)壓原理經(jīng)口腔采集瘤胃液各30 ml，四層紗布過濾后，采用上海雷磁pHS-3C型精密pH計(jì)立即測定pH值，之后將濾液分裝在3個10 ml滅菌離心管中，置于-20℃下保存待測NH3-N濃度和MCP含量。

門診人次序列消除季節(jié)指數(shù)，長期趨勢序列用ARIMA法擬合。趨勢序列經(jīng)二階差分消除趨勢。在p，q≤6中自動尋優(yōu)定階，BIC(1,0)=11.757最小，認(rèn)為AR(1)模型最合適，這也符合簡化原則。條件最小二乘法識別參數(shù)，說明參數(shù)顯著(P=0.0001)，(1+0.57965B)(1-B)2xt=εt為表達(dá)式。經(jīng)模型顯著性檢驗(yàn)，殘差為純隨機(jī)序列，延遲階數(shù)6,12,18,24時，P=0.7355,0.7495,0.8347,0.0672＞0.05，模型對序列擬合顯著有效。對2014年各季度趨勢外推預(yù)測，得到70202,71347,72614,73810，由季節(jié)指數(shù)合并計(jì)算2014年門診人次65205,73213,72283,77473，相對誤差絕對值為0.59%,1.60%,1.21%,1.81%。

2.3曲線擬合法

門診人次序列消除季節(jié)因素后有平穩(wěn)長期趨勢特點(diǎn)，可用曲線擬合法建立模型，如二次或三次拋物線、指數(shù)曲線、Logistic曲線等。經(jīng)試選，本例適合三次拋物線反映趨勢，對近期數(shù)據(jù)變化規(guī)律擬合更好，表達(dá)式為xt=52710+457×t-23.344×t2+0.555×t3。決定系數(shù)為R2=0.982，F(xiàn)檢驗(yàn)P＜0.001，模型擬合不錯，對2014年各季度趨勢外推預(yù)測，得到70424,71808,73286,74860，由季節(jié)指數(shù)合并計(jì)算2014年門診人次65411,73687,72952,78576，相對誤差絕對值為0.91%,2.26%,2.15%,3.26%。

門診人次月度數(shù)據(jù)消除季節(jié)波動后，ARIMA法和趨勢外推法均可用于建立模型。ARIMA法為隨機(jī)時序方法，將歷史序列和殘差序列納入線性模型。不易對事物規(guī)律直觀解釋，優(yōu)點(diǎn)是擬合性能和短期預(yù)測精度高。趨勢外推曲線能解釋事物隨時間變化規(guī)律，須有曲線變化特點(diǎn)，序列波動明顯存在時不太適合。經(jīng)計(jì)算比較，ARIMA效果優(yōu)于多項(xiàng)式曲線。趨勢擬合與外推效果比較見圖4—圖5，組合模型對原始序列擬合及外推效果見圖6所示。

圖4　ARIMA模型

圖5　三次拋物線模型

圖6　三種組合模型擬合外推效果比較

3　建模方式二

考慮門診人次數(shù)據(jù)隨季度呈周期波動變化、隨年度長期平緩趨勢增長特點(diǎn)，借鑒當(dāng)前傳染病預(yù)測文獻(xiàn)中常見季節(jié)ARIMA模型[6～7]，用以建立擬合模型及外推預(yù)測。

3.1簡單季節(jié)AR IMA模型

由于門診人次序列不平穩(wěn)，用周期步長差分消除季節(jié)波動，用低階差分消除平緩增長趨勢，將數(shù)據(jù)轉(zhuǎn)化為平穩(wěn)序列。模型結(jié)構(gòu)如下：

門診人次序列隨季節(jié)周期波動、隨年度有長期趨勢，季節(jié)變化用4步差分，經(jīng)一階差分或二階差分消除趨勢特點(diǎn)，過度差分會損失信息，差分序列若平穩(wěn)化則可以建立模型。經(jīng)計(jì)算分析，以2階差分消除趨勢、以4步差分消除季節(jié)波動，經(jīng)白噪聲檢驗(yàn)得知，延遲6階時，P=0.1452＞0.05，該序列可以建模分析。觀察ACF、PACF是否落入兩倍標(biāo)準(zhǔn)差內(nèi)，調(diào)試階數(shù)、確定結(jié)構(gòu)和識別參數(shù)。模型結(jié)構(gòu)以簡單為原則，發(fā)現(xiàn)ACF一階以后、PACF一階及四階時落入兩倍標(biāo)準(zhǔn)差范圍并繼續(xù)拖尾，考慮低階疏系數(shù)模型，模型ARIMA((4)，1，1)表達(dá)式如下：

(1+0.38669B4)(1-B4)(1-B)2xt=(1+0.90475)εt。經(jīng)延遲6,12,18,24階時，P=0.4509,0.8271,0.9703,0.5971，說明殘差為白噪聲序列。對2014年門診人次各季度外推預(yù)測，得到65127,72420,71452,76067；計(jì)算相對誤差絕對值為0.47%,0.50%,0.05%,0.03%，預(yù)測效果非常好。

3.2乘積季節(jié)AR IMA模型

如果季節(jié)波動與長期趨勢間有復(fù)雜交互影響，要用乘積季節(jié)ARIMA模型。該方法在傳染病發(fā)病率、發(fā)病數(shù)預(yù)測中是當(dāng)前熱點(diǎn)，可能是傳染病月度或季度資料特點(diǎn)決定的。對于本例門診人次預(yù)測來說，簡單季節(jié)模型經(jīng)計(jì)算預(yù)測效果不錯，可認(rèn)為季節(jié)波動與短期相關(guān)間并無復(fù)雜關(guān)系，以下仍用乘積季節(jié)模型驗(yàn)證和進(jìn)行比較。

以ARMA(p，q)提取短期相關(guān)性，以ARMA(P，Q)提取季節(jié)相關(guān)性。假設(shè)二者存在交互乘積關(guān)系，構(gòu)造為乘積模型: ARIMA(p，d，q)×(P，D，Q)S:Φ(B)ΦS(B)▽D▽dxt=Θ(B)ΘS(B)εt;

其中Θ(B)=(1-θ1B-…-θqBq),

前面得知長期趨勢用2階差分，季度周期用四步差分。調(diào)試階數(shù)及識別參數(shù)，模型結(jié)構(gòu)ARIMA(1，2，0)×(1，1，0)4，表達(dá)式為(1-0.18884B)(1-0.32985B4)(1-B4)(1-B)2xt。=εt殘差序列延遲6,12,18,24時，P=0.1178,0.5641,0.6726,0.1108＞0.05，認(rèn)為模型顯著。對2014年門診人次各季度外推預(yù)測，得到65061,72547,71978,76881，計(jì)算相對誤差絕對值0.37%, 0.68%,0.78%,1.04%。然后根據(jù)ACF和PACF特點(diǎn)，繼續(xù)選擇可能模型結(jié)構(gòu)如ARIMA(1，2，1)×(1，1，0)4，表達(dá)式為(1+0.40425B)(1-0.25856B)(1-B4)(1-B)2xt=(1+0.9655B)εt。延遲6,12,18,24時，殘差序列P=0.1878,0.6384,0.7670,0.2152＞0.05，認(rèn)為模型顯著。對2014年門診人次各季度外推預(yù)測，得到64988,72311,71860,76605，計(jì)算相對誤差絕對值0.25%, 0.35%,0.62%,0.67%。

乘積季節(jié)ARIMA模型預(yù)測精度稍微改善，又反復(fù)試取了多種模型結(jié)構(gòu)，但相對簡單季節(jié)ARIMA模型并未明顯改善預(yù)測精度，認(rèn)為季節(jié)波動和短期相關(guān)無交互關(guān)系，比較見圖7所示。

圖7　兩種季節(jié)模型擬合外推效果

以上從兩種方式、四種模型對門診人次季度數(shù)據(jù)進(jìn)行預(yù)測。第一種方式：X11過程法與ARIMA、多項(xiàng)式曲線設(shè)計(jì)組合模型。第二種方式：假設(shè)季節(jié)波動和短期相關(guān)性是否有交互關(guān)系，構(gòu)造簡單季節(jié)ARIMA模型、乘積季節(jié)ARIMA模型。2004—2013年季度數(shù)據(jù)用于確定模型結(jié)構(gòu)、識別參數(shù)，對2014年情況預(yù)測情況比較見表2所示。

表2　模型預(yù)測精度比較

4　討論

門診、收入等醫(yī)療統(tǒng)計(jì)工作領(lǐng)域多數(shù)指標(biāo)為例，月度、季度數(shù)據(jù)隨時間變化往往有規(guī)律性。影響因素復(fù)雜而無法確定模型結(jié)構(gòu)或識別參數(shù)，可以嘗試由時序數(shù)據(jù)由時間序列分析方法建模，時序數(shù)據(jù)可分解為長期趨勢、季節(jié)波動和隨機(jī)無序部分，由數(shù)據(jù)資料計(jì)算季節(jié)指數(shù)，消除周期波動后分析長期趨勢，也可直接應(yīng)用季節(jié)ARIMA模型。

一種方式是，鑒于季節(jié)指數(shù)較明顯，以X11法作為時間序列分解法提取信息，原始序列不再有季節(jié)波動特點(diǎn)，考慮這類序列有平滑且類似曲線變化特點(diǎn)，以某類曲線擬合時序規(guī)律，也可由ARIMA法直接建模。經(jīng)外推預(yù)測后將預(yù)測值與季節(jié)指數(shù)合成未來預(yù)測值。另一種方式是，序列季節(jié)波動與長期趨勢不給予分解提取，而是直接建立季節(jié)ARIMA模型，若季節(jié)波動與短期相關(guān)交互不明顯則用簡單季節(jié)模型，若交互影響存在則用乘積季節(jié)模型。

以某醫(yī)院門診人次預(yù)測為例，季度數(shù)據(jù)周期性變化有規(guī)律，歷史數(shù)據(jù)中季節(jié)波動且平穩(wěn)增長，觀測期長、數(shù)據(jù)連貫、規(guī)律性強(qiáng)，明顯有季節(jié)波動和平穩(wěn)延續(xù)特點(diǎn)，后期增幅明顯，長期較平緩且無轉(zhuǎn)折因素影響。采用以上幾種模型演示實(shí)施流程，為醫(yī)療統(tǒng)計(jì)工作提供方法借鑒。兩種方式又分別給出兩種算法，比較模型原理特點(diǎn)、擬合性能和預(yù)測精度。以2004—2013年數(shù)據(jù)擬合模型，預(yù)測2014年情況并與真實(shí)值計(jì)算相對誤差，X11-ARIMA法1%左右，X11-三次多項(xiàng)式2%左右，簡單季節(jié)ARIMA0.5%以下，乘積季節(jié)ARIMA0.5%左右。

第一種方式為組合建模方法，適用于季節(jié)指數(shù)較明確情況，消除季節(jié)波動以后的長期趨勢規(guī)律，ARIMA法好于曲線擬合法。第二種方式直接用季節(jié)ARIMA模型，經(jīng)算例驗(yàn)證表明季節(jié)ARIMA模型擬合性能、預(yù)測精度均好于第一種方式。預(yù)測領(lǐng)域多數(shù)文獻(xiàn)表明，月度或季度數(shù)據(jù)基本上用乘積季節(jié)ARIMA模型，這種思路未必千篇一律，經(jīng)本例驗(yàn)證，門診人次數(shù)據(jù)中季節(jié)波動和短期相關(guān)性無交互影響，簡單季節(jié)ARIMA模型甚至優(yōu)于乘積季節(jié)ARIMA模型，具體資料情況下可以進(jìn)行試選優(yōu)化，以便改善模型對資料擬合及預(yù)測適應(yīng)性。

[1]王燕.應(yīng)用時間序列分析[M].北京:中國人民大學(xué)出版社,2012.

[2]徐國祥.統(tǒng)計(jì)預(yù)測與決策[M].上海:上海財(cái)經(jīng)大學(xué)出版社,2012.

[3]耿娟.ARIMA模型在醫(yī)院門診量預(yù)測中的應(yīng)用[J].中國衛(wèi)生統(tǒng)計(jì), 2014,31(8).

[4]申銅倩,劉文東,胡建利等.X11-ARIMA過程在痢疾疫情預(yù)測中的應(yīng)用研究[J].中國衛(wèi)生統(tǒng)計(jì),2014,31(3).

[5]李望晨.基于增長特征法與ARIMA的人均衛(wèi)生事業(yè)費(fèi)趨勢預(yù)測比較研究[J].中國衛(wèi)生統(tǒng)計(jì),2014,31(3).

[6]陸波,閔思韜,閔紅星等.應(yīng)用ARIMA模型預(yù)測麻疹發(fā)病率的可行性研究[J].中國衛(wèi)生統(tǒng)計(jì),2014,32(1).

[7]張愛紅,周培,申銅倩等.乘積季節(jié)ARIMA模型在食源性疾病預(yù)測中的應(yīng)用[J].中國衛(wèi)生統(tǒng)計(jì),2014,31(3).

（責(zé)任編輯/浩天）

C913.4

A

1002-6487（2016）19-0078-03

教育部人文社科基金資助項(xiàng)目（15YJCZH087；14YJAZH101；14YJA630098）；山東自然科學(xué)基金資助項(xiàng)目（ZR2015HL101）；山東統(tǒng)計(jì)科研項(xiàng)目（KT16230；KT16231）

李望晨（1980—），男，山東濰坊人，碩士，副教授，研究方向：衛(wèi)生管理統(tǒng)計(jì)。于貞杰（1971—），女，山東濰坊人，博士，教授，研究方向：衛(wèi)生事業(yè)管理。王在翔（1963—），男，山東濰坊人，碩士，教授，研究方向：衛(wèi)生管理統(tǒng)計(jì)。（通訊作者）張利平（1980—），男，山東濰坊人，博士研究生，講師，研究方向：職業(yè)衛(wèi)生評價。