張昌爽

（重慶大學(xué)資源與安全學(xué)院，重慶 400044）

工傷事故灰色預(yù)測模型展現(xiàn)工傷事故時間狀態(tài)序列數(shù)據(jù)及其相應(yīng)變化總體波動趨勢。建立灰色馬爾可夫預(yù)測模型，用馬爾可夫預(yù)測法確定未來狀態(tài)轉(zhuǎn)移規(guī)律，馬爾可夫灰色預(yù)測具備顯著優(yōu)勢。馬爾可夫鏈為狀態(tài)可列、時間齊次馬爾可夫過程[1-3]。

1 灰色權(quán)重馬爾可夫鏈預(yù)測步驟

1.1 狀態(tài)劃分

由有序聚類法構(gòu)造工傷人數(shù)分級標(biāo)準(zhǔn)，構(gòu)建馬爾可夫鏈狀態(tài)指標(biāo)值分級空間，確定各時段指標(biāo)所對應(yīng)狀態(tài)[4-6]。變量x1，x2，x3歸類為{xi,…,xj},j≥i，定義均值向量：

式中，D(i,j)為量化變量段內(nèi)部各變量間的差異值，以徑向值{xi，…，xj}(j≥i)進(jìn)行確定；徑向量數(shù)值越小，表示區(qū)間向量差異越小，彼此越接近。

將n個變量序列分為K類，分類誤差函數(shù)如下：

式中，k未知分類數(shù)，計算k=1，2，…，n所有最佳分割。

誤差函數(shù)e[P(n,k)]通過線性規(guī)劃，對應(yīng)最優(yōu)分類標(biāo)準(zhǔn)通過在min{e[P(n,k)]}中得以確定，確定分類數(shù)K變化趨勢，作出誤差函數(shù)e[P(n,k)]與K的關(guān)系圖，函數(shù)曲線拐點對應(yīng)的K值為最優(yōu)分類數(shù)。

1.2 確定狀態(tài)轉(zhuǎn)移概率矩陣

轉(zhuǎn)移概率以{x1，x2，…，xn}定義馬爾可夫鏈指標(biāo)序列，包含m個狀態(tài)，指標(biāo)值序列從狀態(tài)i經(jīng)過1步轉(zhuǎn)移到達(dá)狀態(tài)j的頻數(shù)以fij表示，i，j∈E，轉(zhuǎn)移步長是1個時間單位，2個甚至m個時間單位。

式中，fij為轉(zhuǎn)移頻數(shù)矩陣(fij)組成單元，轉(zhuǎn)移頻數(shù)矩陣的第i行第j列元素fij除以對應(yīng)各行的單元值總和所得，表示為Pij，i，j∈E。

1.3 馬氏性檢驗

用χ2統(tǒng)計量檢驗離散型序列隨機(jī)變量。M表示隨機(jī)變量指標(biāo)值序列存在的演變狀態(tài)，fij為指標(biāo)值序列{x1，x2，…，xn}從狀態(tài)i經(jīng)過1步轉(zhuǎn)移到達(dá)狀態(tài)j的頻數(shù)值，i，j∈E。

式中，p.j為邊際概率，轉(zhuǎn)移頻數(shù)矩陣第j列之和除以矩陣所有行和列總和。

當(dāng)n的基數(shù)足夠大，統(tǒng)計量驗證指標(biāo)：

式中，pij為轉(zhuǎn)移概率。

在顯著性水平a，服從自由度(m-1)2的χ2分布，查表得分位點的客觀數(shù)值，計算后得研究指標(biāo)統(tǒng)計量χ2值。若，則xi符合馬氏性檢驗標(biāo)準(zhǔn)，否則該數(shù)據(jù)序列不適用于馬爾可夫鏈來預(yù)測。

1.4 各階自相關(guān)系數(shù)計算

自相關(guān)系數(shù)rk，k∈e；rk為第k階（滯時為k季度）自相關(guān)系數(shù)；xl為第l 季度工傷人數(shù)；x季度平均受傷人數(shù)均值；n季度平均受傷人數(shù)序列長度。

式中，rk為各階自相關(guān)系數(shù)。

規(guī)范各階自相關(guān)系數(shù)，即：

式中，wk步長馬爾可夫鏈權(quán)重；m預(yù)測最大階數(shù)。

以狀態(tài)特征值法、線性插值法予以推斷，將預(yù)測數(shù)據(jù)加入原矩陣序列之中，重復(fù)上述步驟。

式中，μ為狀態(tài)特征值。

其中，i≤μ＜i+1，β＞0為狀態(tài)調(diào)整因子，狀態(tài)區(qū)間下限為，上限為，預(yù)測指標(biāo)值通過上下限兩類算法，結(jié)果一致。

式中，x為指標(biāo)預(yù)測值；k=0，狀態(tài)值增加與指標(biāo)值一致；否則，k=1。

非零向量πj={π1,π2,…,πE}，使得πjPij=πj，其中Pij為狀態(tài)轉(zhuǎn)移概率矩陣，有πj={π1,π2,…,πE}為馬爾可夫鏈平穩(wěn)分布，處于平穩(wěn)狀態(tài)。當(dāng)存在m＞0，pm中諸元素皆非負(fù)非零，則p為正規(guī)概率矩陣。

2 模型應(yīng)用

以2010年1季度至2019年第二季度共38季度工傷人數(shù)預(yù)測2019年3季度的受傷人數(shù)，將2019年3季度實際預(yù)測數(shù)據(jù)加入序列模型中，預(yù)測2019年4季度受傷人數(shù)，并與社保局實際備案的工傷人數(shù)對比，驗證預(yù)測效果。

2.1 分級準(zhǔn)則確定狀態(tài)

用有序聚類法將序列分5個人數(shù)區(qū)間比較，工傷人數(shù)等級標(biāo)準(zhǔn)以1～5表示，狀態(tài)依次由低到高（見表1）。

表1 工傷人數(shù)狀態(tài)劃分表

2.2 馬氏性檢驗

2010年1 季度至2019年第4季度共40個季度某企業(yè)工傷事故人數(shù)序列作馬氏性檢驗。每季度工傷人數(shù)頻數(shù)轉(zhuǎn)移矩陣(fij)5×5、步長為1的1步轉(zhuǎn)移概率矩陣(pij)5×5如下：

由(fij)5×5及(pij)5×5得邊際概率、統(tǒng)計量χ2值（見表2、表3）：

表2 邊際概率

表3 統(tǒng)計量χ2計算表

計算χ2為54.5556，在顯著性水平a=0.05下，查表得分位點由于因此工傷人數(shù)序列滿足馬氏性檢驗。

2.3 馬爾可夫鏈權(quán)重確定

2010 年1 季度至2019 年第二季度共38 季度工傷數(shù)據(jù)，受傷人數(shù)序列各階自相關(guān)系數(shù)為：r1=0.8633，r2=0.6697，r3=0.4693，r4=0.2874，r5=0.0665。自相關(guān)系數(shù)歸一平均化，作為各滯時馬爾可夫鏈權(quán)重：w1=0.3664，w2=0.2842，w3=0.1992，w4=0.1220，w5=0.0282。

2.4 不同步長轉(zhuǎn)移概率矩陣

根據(jù)2010年1季度至2019年第二季度共38個季度工傷數(shù)據(jù)，馬爾可夫工傷序列統(tǒng)計得不同滯時轉(zhuǎn)移概率矩陣：

2.5 創(chuàng)建預(yù)測表

以2018年第二季度至2019年第二季度共5個季度工傷人數(shù)、狀態(tài)轉(zhuǎn)移概率矩陣對2019年第三季度工傷人數(shù)進(jìn)行預(yù)測分析（見表4）。

表4 2019年第3季度工傷事故人數(shù)預(yù)測表

2019 年第三季度實際工傷人數(shù)為24 人/季度，由表可知，max(pi,i∈E)=0.5776，i=3，即預(yù)測2021年第三季度工傷人數(shù)狀態(tài)為3，工傷人數(shù)x滿足：26≤x＜39，預(yù)測數(shù)值與實際受傷人數(shù)相近，預(yù)測有效。取調(diào)整因子β=1，得狀態(tài)特征值μ=3.131。下限法得2019年第三季度工傷人數(shù)預(yù)測值：

x=+此值與實際相對誤差僅有1.25%，以2010年第1季度至2019年第3季度共39個季度工傷人數(shù)序列預(yù)測2019年第4季度的工傷狀態(tài)，得各階自相關(guān)系數(shù)、各步長馬爾可夫鏈權(quán)重（見表5）。

表5 2010.1-2019.3季度序列各階自相關(guān)系數(shù)及權(quán)重

2010年第一季度至2019年第3季度共39個季度工傷事故受傷人數(shù)，不同步長（滯時）轉(zhuǎn)移概率矩陣：

以2018年3季度至2019年第3季度共5個季度工傷人數(shù)預(yù)測2019年第4季度的工傷人數(shù)，預(yù)測見表6。

表6 2019.4季度工傷事故人數(shù)預(yù)測

2019 年 第4 季 度 實 際 工 傷 人 數(shù) 為22 人/ 季 度，max(pi,i∈E)=0.4248，i=2，即預(yù)測2019年第4季度工傷人數(shù)狀態(tài)為2，受傷人數(shù)x滿足：13≤x＜26，預(yù)測結(jié)果與實際情況完全一致，預(yù)測有效。

取調(diào)整因子β=0.95，狀態(tài)特征值μ=2.0231。下限法得2019年第4季度人工傷預(yù)測值：

此值與實際相對誤差僅有1.36%

2.6 不同步長馬爾可夫鏈特征分析

此鏈狀態(tài)關(guān)聯(lián)互通，對任意i,j∈E,i←→j,(i≠j)，非周期，故為不可約正常返鏈。設(shè)平穩(wěn)分布為(πj,j∈E)：以步長為2的轉(zhuǎn)移概率矩陣p(2)，平穩(wěn)分布、極限分布、各狀態(tài)重現(xiàn)期見表7。

表7 平穩(wěn)分布、極限分布及各狀態(tài)復(fù)現(xiàn)期

狀態(tài)j的復(fù)現(xiàn)期為Tj，其對應(yīng)概率pj=1/Tj=πj。各狀態(tài)的復(fù)現(xiàn)期為T1=19.2678（季度），T2=3.8536（季度），T3=2.7525（季度），T4=5.6721（季度），T5=6.2696（季度）。在中長期預(yù)測中，企業(yè)工傷人數(shù)為狀態(tài)3的概率最大，平均每2.7525季度重現(xiàn)，復(fù)現(xiàn)概率為0.3633；其次復(fù)現(xiàn)的為狀態(tài)2，每3.8536季度復(fù)現(xiàn)1次，概率為0.2595；復(fù)現(xiàn)概率最低的為狀態(tài)1為0.0519，復(fù)現(xiàn)期19.2678季度。

3 結(jié)論

1）加權(quán)馬爾可夫模型狀態(tài)的復(fù)現(xiàn)期T1=19.2678（季度），T2=3.8536（季度），T3=2.7525（季度），T4=5.6721（季度），T5=6.2696（季度）。長期預(yù)測中，企業(yè)工傷人數(shù)狀態(tài)3的概率最大，平均每2.7525季度重現(xiàn)，復(fù)現(xiàn)概率0.3633；其次復(fù)現(xiàn)的為狀態(tài)2，每3.8536季度復(fù)現(xiàn)1次，概率為0.2595；復(fù)現(xiàn)概率最低的為狀態(tài)1，為0.0519，復(fù)現(xiàn)期19.2678季度。

2）實際應(yīng)用價值在于以某城市溫室氣體排放序列來預(yù)測未來溫室氣體年排放狀況，以某國家煤炭資源消耗量序列來預(yù)測未來煤炭年消耗狀況等。