趙立戎

（中國人民大學(xué) 統(tǒng)計學(xué)院，北京100872）

0 引言

對于一個廣義線性模型（GLM），按照貝葉斯統(tǒng)計推斷的理論，假設(shè)GLM的參數(shù)向量服從某個形式已知的先驗分布，其后驗分布就是該參數(shù)向量關(guān)于觀測數(shù)據(jù)的條件分布，通過計算參數(shù)向量的后驗分布解決參數(shù)估計問題。在某些特殊情況如一元離散指數(shù)族下，能夠找到自然共軛的先驗分布，因而容易計算出后驗分布的具體形式，此結(jié)論在多元離散指數(shù)族中并不成立。因此在一般情況下，貝葉斯GLM的后驗分布形式很難直接計算，因而產(chǎn)生了很多近似計算的方法。例如West et al.在保留自然參數(shù)服從自然共軛先驗分布的約束條件下給出了幾種參數(shù)變換的形式。Fahrmeir和Kaufman研究了估計后驗分布的模。近年來還有不少利用MCMC方法處理貝葉斯GLM問題的文獻。在精算領(lǐng)域的應(yīng)用如Scollnik及de Alba等人的研究。MCMC方法便于操作，但得不出后驗分布的具體形式。在實際應(yīng)用中，近似的解析法雖然精度高但計算過于復(fù)雜，而隨機模擬方法雖然便于實現(xiàn)卻又得不到分布形式。

Greg Taylor[4]利用泰勒展開式，在多元離散指數(shù)族下構(gòu)造出后驗分布的二階近似解析表達式，使其滿足與貝葉斯GLM的先驗分布同分布族的性質(zhì)，并能夠進一步產(chǎn)生迭代公式，生成一個與Kalman濾波算法類似的GLM濾波算法，可推廣至動態(tài)廣義線性模型（DGLM）的參數(shù)估計。這種利用二階近似的GLM濾波算法不同于復(fù)雜的傳統(tǒng)解析處理方法，通過使用典則聯(lián)結(jié)或伴隨典則聯(lián)結(jié)函數(shù)能夠簡化計算過程。

在非壽險準(zhǔn)備金評估中，通常使用的確定性模型或隨機性模型都假設(shè)一個不隨時間變化的參數(shù)向量，它們都屬于靜態(tài)模型。當(dāng)保險公司內(nèi)部或外部環(huán)境隨時間發(fā)生變化時，如果考慮模型參數(shù)及參數(shù)結(jié)構(gòu)隨時間發(fā)生變化，就可以使用動態(tài)模型，如精算中常用的Kalman濾波模型，張連增[6]研究了Kalman濾波法在非壽險未決賠款準(zhǔn)備金評估中的應(yīng)用。Kalman濾波模型的應(yīng)用局限性在于假設(shè)觀測數(shù)據(jù)服從正態(tài)分布，不能處理其他常用分布（如泊松分布、伽瑪分布等）類型的數(shù)據(jù)。

本文將從動態(tài)模型的角度出發(fā)，探討非壽險未決準(zhǔn)備金評估問題：應(yīng)用廣義線性濾波模型，基于離散指數(shù)族分布，建立關(guān)于賠付數(shù)據(jù)的動態(tài)廣義線性模型；通過對流量三角形數(shù)據(jù)逐期（事故年或日歷年）迭代，實現(xiàn)模型的參數(shù)估計和準(zhǔn)備金評估。

1 廣義線性濾波模型

1.1 貝葉斯GLM的二階近似法

GLM 的標(biāo)準(zhǔn)結(jié)構(gòu)為：Y=h-1（Xβ）+ε，其中 Y=（Y1，…，Ym）T是觀測數(shù)據(jù)向量，各個Yi（i=1，…m)之間相互獨立且服從離散指數(shù)族（EDF）中的分布，EDF分布的對數(shù)似然函數(shù)形如：L（yi：θi，λi，b，k）=λi［yiθi-b（θi）］+k（λi，yi），β=（β1，… ，βq)T（q≤m)是參數(shù)向量，X是m×q階設(shè)計矩陣，h為聯(lián)結(jié)函數(shù)，ε是一個隨機誤差向量。可見，Yi具有相同的分布形式，但位置參數(shù)θi和尺度參數(shù) 1/λi可以取不同的值。 E（Y）=μ（θ）=h-1（Xβ）。

使用典則聯(lián)結(jié)函數(shù) h=（b′）-1時，模型形式簡化為 θ=Xβ，位置參數(shù)就等于線性響應(yīng)部分。不失一般性，以下討論均在典則聯(lián)結(jié)函數(shù)下進行。忽略不含自然參數(shù)θ的項，Y的對數(shù)似然函數(shù)簡化為 L（y；β，Λ）=yTΛXβ-1TΛb（Xβ），其中 Λ=diag(λ1，…，λm)，1表示所有元素均為1的m階向量。假設(shè)β的先驗對數(shù)似然形式為 π（β）=Mβ-（Mβ），則 β 的后驗對數(shù)似 然 形 式 為 L （β|y，Λ）=（w0TM+yTΛX）β-［n0Tb （Mβ）+1TΛb（Xβ）］，其中M是一個q階非奇異方陣。對后驗似然中的每一項分別按照泰勒展開式在β=β0處展開，最多取到二次項，忽略與β無關(guān)的零次項，再通過一系列正交變換，可以構(gòu)造形如PMβ的近似后驗似然，其中P表示一個線性變換，w1和n1分別對應(yīng)w0和n0。為了方便計算，可選定β0的取值為滿足 h-1（Mβ0）=Eβ［h-1（Mβ）］。 通過二階近似變換后，參數(shù)向量β的先驗似然與后驗似然具有相似的形式，且存在以下替代關(guān)系：w0→w1，n0→n1，M→PM。 因而，在 β 的先驗似然和后驗似然之間形成了一種遞推關(guān)系。值得注意的是，對于非典則聯(lián)結(jié)的其他聯(lián)結(jié)函數(shù)，變換公式要復(fù)雜的多，遞推關(guān)系也未必成立。常用的分布及典則聯(lián)結(jié)函數(shù)有正態(tài)分布和單位聯(lián)結(jié)函數(shù)，伽瑪分布和倒數(shù)聯(lián)結(jié)函數(shù)等。

對于勢方差函數(shù)， 有 h （μ）=μ1-p/（1-p）， 不滿足對任意p∈（-∞，+∞）成立，因而可能使得 b（Xβ）在某些值上沒有定義?？紤]到典則聯(lián)結(jié)函數(shù)存在的這種問題，Taylor[4]定義伴隨典則聯(lián)結(jié)函數(shù)（companion canonical link）h*=b?（b′）-1，并且在伴隨典則聯(lián)結(jié)函數(shù)下得到了與使用典則聯(lián)結(jié)函數(shù)時一致的似然函數(shù)遞推關(guān)系，僅其中涉及參數(shù)的形式發(fā)生一定的變化。常用的分布及伴隨典則聯(lián)結(jié)函數(shù)有伽瑪分布和對數(shù)聯(lián)結(jié)函數(shù)等。

1.2 動態(tài)廣義線性模型

考慮參數(shù)向量β隨時間隨機變化，則GLM模型形式變?yōu)椋篩（t）=h-1（X（t）β（t））+v（t），其中 Y（t）是 Y 在 t時刻的觀測數(shù)據(jù)向量，設(shè)計矩陣 X（t）、參數(shù)向量 β（t）和誤差向量 v（t）都表示t時刻的相應(yīng)狀態(tài)。類似狀態(tài)空間模型的定義，上式為觀測方程，還需定義參數(shù)演化的狀態(tài)方程：β（t+1）=Φ(t+1)β（t）+r（t+1），其中 Φ（t+1）為已知矩陣，r（t+1）是均值為 0 的隨機干擾項，且與 y（1），…，y（t）及 β（0），…，β（t）相互獨立，記Var［r（t+1）］=R（t+1）。允許參數(shù)向量的維度隨時間變化。這樣就構(gòu)成了動態(tài)廣義線性模型（DGLM）的基本形式。在實際應(yīng)用中，t可以選為事故發(fā)生年，也可以選為日歷年。

利用貝葉斯二階近似法得到的先驗似然與后驗似然之間的遞推關(guān)系，可以在保留分布形式的同時通過不斷變換參數(shù)實現(xiàn)多次貝葉斯GLM的模型估計。這就形成了一個與Kalman濾波類似的貝葉斯GLM序列，可以稱之為貝葉斯GLM濾波或廣義線性濾波，利用這一濾波算法實現(xiàn)對DGLM各時刻的參數(shù)估計。

類似狀態(tài)空間模型的定義，用記號t|t-j表示在t時刻利用t-j時刻之前（包t-j含時刻）的全部信息所進行的估計。使用典則聯(lián)結(jié)函數(shù)時，似然函數(shù)的形式轉(zhuǎn)變?yōu)椋?/p>

為了估計 E［h-1（M（t|t）β（t）|Y（t））］和 E［h-1（M（t+1|t）β（t+1）|Y（t）］的值及對應(yīng)的方差，在估計過程中應(yīng)用貝葉斯GLM二階近似的結(jié)論，隨著觀測數(shù)據(jù)向量 y（1）,y（2），…逐漸引入模型，對β（?）的均值及方差的估計值逐步遞推調(diào)整，即產(chǎn)生一個濾波，濾波的起始點是對 E［β（1）］和 Var［β（1）］的先驗估計值。在應(yīng)用廣義線性濾波模型評估非壽險準(zhǔn)備金時，初始值的選定依賴于對經(jīng)驗數(shù)據(jù)等先驗信息的分析和精算師的主觀判斷，隨著迭代步驟的增加，初始值對整個計算結(jié)果的敏感性將越來越低。在濾波迭代的過程中，每一步迭代都是從 E［h-1（M（t|t-1）β（t））|Y（t-1）］和 Var［h-1（M（t|t-1）β（t）|Y（t-1））］通過一系列運算公式遞推到 E［h-1（M（t+1|t）β（t+1））|Y（t）］和 Var［h-1（M（t+1|t）β（t+1）|Y(t））］，t=1，2…，即從利用t-1時刻之前的信息對t時刻的估計值前進到利用t時刻之前的信息對t+1時刻的估計值，直至迭代到流量三角形的最近一期觀測值為止。對于一般的貝葉斯GLM，先驗似然πt|t-1屬于EDF自然共軛分布族，其后驗似然πt+1|t通常并不屬于同一分布族，只有在迭代過程中應(yīng)用貝葉斯二階近似法，用同屬EDF自然共軛分布族的二階近似式代替原后驗似然，貝葉斯GLM濾波迭代才能夠?qū)崿F(xiàn)遞推。DGLM的濾波迭代的邏輯結(jié)構(gòu)如圖1所示。

Taylor[4]將Kalman濾波算法從動態(tài)一般線性模型推廣到DGLM，保留了模型的線性系統(tǒng)成分。應(yīng)注意的是，當(dāng)使用正態(tài)分布和單位聯(lián)結(jié)函數(shù)時，廣義線性濾波就等于Kalman濾波，可見Kalman濾波是廣義線性濾波的一種特殊情況。廣義線性濾波模型的一般形式及其實現(xiàn)過程仍是比較復(fù)雜的，但對于幾種特殊分布（如正態(tài)分布、泊松分布、伽瑪分布）以及使用典則聯(lián)結(jié)或伴隨典則聯(lián)結(jié)函數(shù)的情況，濾波迭代公式能夠大為簡化，因此在解決非壽險準(zhǔn)備金評估問題時具有實際應(yīng)用價值。廣義線性濾波模型對各期參數(shù)的估計是基于迭代的過程，因此各期參數(shù)估計值比一般GLM的參數(shù)估計值較平滑，能夠降低經(jīng)驗數(shù)據(jù)中的異常值對準(zhǔn)備金估計結(jié)果的影響。在早期賠付數(shù)據(jù)缺失或觀測數(shù)據(jù)出現(xiàn)異常值的情況下，應(yīng)用廣義線性濾波模型能夠改善未決賠款準(zhǔn)備金的估計結(jié)果。

2 伽瑪廣義線性濾波模型

在GLM中假設(shè)誤差項服從伽瑪分布，使用對數(shù)聯(lián)結(jié)函數(shù)（伴隨典則聯(lián)結(jié)），即 b（θ）=-log（-θ），h（u）=logu。 建立伽瑪廣義線性濾波模型如下：

觀測方程：Y（t+1）=exp{X（t+1）β（t+1）}+v（t+1)

狀態(tài)方程：β（t+1）=Φ（t+1）β（t）+r（t+1）

為了簡便，將濾波迭代主要關(guān)注的兩個量記為：

Θ（t|t-j）=E［exp{-M（t|t-j）β（t|t-j）}］

Γ（t|t-j）=Var［exp{-M（t|t-j）β（t|t-j）}］

其中矩陣 M（t|t-j）是使參數(shù)方差陣 Var［β（（t|t-j）］對角化的正交矩陣，即使得 M（t|t-j）TQ（t|t-j）M（t|t-j）=Var［β（（t|t-j）］成立，Q（t|t-j）是對角陣。

首先需要選定濾波迭代的初始值，即為E［β（0|0）］和Var［β（0|0））賦先驗值。其他需要選定的輸入值還包括伽瑪分布的參數(shù) Λ（t），可以利用觀測數(shù)據(jù) Y（t）的變異系數(shù)（CV）來計算（λ=1/CV2）。矩陣Φ（t）根據(jù)模型參數(shù)結(jié)構(gòu)隨時間變化的規(guī)律確定。隨機擾動變量的方差R(t)根據(jù)精算師的判斷人為選取，通常選為一個常矩陣，即R（t）=R。

伽瑪廣義線性濾波迭代計算的具體步驟如下：

第一步，計算 E［β（t|t-1）］=Φ（t）E［β（t-1|t-1）］

Var［β（t|t-1）］=Φ（t）Var［β（t-1|t-1）］ΦT（t）+R（t）

近似計算：

Θ（t|t-1）=［I+0.5×Q（t|t-1）］×exp{-E［M（t|t-1）β（t|t-1）］}

Γ（t|t-1）=Q（t|t-1）DIAGexp{-2E ［M（t|t-1）β（t|t-1）］}，

定義 N （t|t-1）=［Γ（t|t-1）］-1［I+0.5×Q（t|t-1)］DIAGexp{-M（t|t-1）E［β（t|t-1）］}

W（t|t-1）=N（t|t-1）Θ（t|t-1）

則 β （t|t-1）=-M（t|t-1）Tlog（N（t|t-1）-1W（t|t-1））。

定義 B（t|t-j）=-DIAGexp{-M（t|t-j）β（t|t-1）}，（j=0，1，…）

G（t）=-DIAGexp{-X（t）β（t|t-1）}

且矩陣 P（t）及 D（t）分別為正交陣和對角陣，滿足

P（t）TD（t）P（t）=N（t|t-1）B（t|t-1）+M（t|t-1）X（t）TDIAG［G（t）Λ（t）Y（t）］X（t）M（t|t-1）T

則 M（t|t）=P（t）M（t|t-1），且

J（t）=P（t）TB（t|t）P（t）［P（t）TD（t）P（t）］-1M（t|t-1）X（t）TΛ（t）

第二步，計算增益矩陣 K（t）=B（t|t）-1P（t）J（t）G（t），代表為新觀測值所賦的信度因子；

第四步，計算 Θ（t|t）=DIAGexp［-M（t|t）β（t|t-1）］［1-K（t）Y（t）-（t|t-1）］，用當(dāng)前觀測值與前期預(yù)測值之差 Y（t）-（t|t-1）作為狀態(tài)參數(shù)更新的基礎(chǔ)。

近似計算：Γ（t|t）=-D（t）-1B（t|t）2

另外，計算 E［β（t|t）］=M（t|t）T{-logΘ （t|t） +0.5 ×Γ （t|t）［DIAGΘ-1（t|t）］21}，Var［β（t|t）］=M （t|t）T［DIAGΘ-1（t|t）］Γ （t|t）［DIAGΘ-1（t|t）］M（t|t），作為下一輪迭代的基礎(chǔ)。

根據(jù)近似公式：

3 實例分析

為了更好地驗證廣義線性濾波模型的估計效果，下面引用 Taylor和 Ashe（1983）[2]的數(shù)據(jù)（表 1）對濾波迭代模型進行實證分析。

根據(jù)流量三角形數(shù)據(jù)的特點，用i表示事故年，j表示進展年，Cij表示事故年i進展j的增量賠款額，Ei表示事故年i的風(fēng)險暴露。以Yij=Cij/Ei作為觀測值。建立伽瑪廣義線性濾波模型，利用Hoerl曲線[1]構(gòu)造GLM的線性響應(yīng)部分。假設(shè)Yij?Gamma（μij，λj），選用對數(shù)聯(lián)結(jié)函數(shù)，μij=exp{βi0+βi1lnj+βi2j}。由于本例中的參數(shù) βi=（βi0，βi1，βi2）T僅與事故年 i有關(guān)， 按照事故年 i（即逐行）進行濾波迭代，則Y（t）=（Yt1，Yt2，…，Yt，11-t）T，β（t）=（βt0，βt1，βt2）T，且

表1 增量賠款額流量三角形

本例中 β（t） t=1，…，10 參數(shù)結(jié)構(gòu)未發(fā)生變化，則 Φ（t）為單位陣。 迭代初始值 E［β（0|0）］和 Var［β（0|0）］選用伽瑪GLM的參數(shù)估計結(jié)果。參數(shù)向量β的濾波估計值E［β（t|t）］如表2所示。由圖2可見由于濾波模型的估計參數(shù)之間存在遞推關(guān)系，各事故年的參數(shù)估計值比較平穩(wěn)，反映了增益矩陣K（t）的平滑作用，降低了觀測數(shù)據(jù)中異常值的影響。

表2 濾波參數(shù)估計值

表3 準(zhǔn)備金估計結(jié)果及穩(wěn)定性比較

協(xié)方差陣Var［β（t|t）］的濾波估計結(jié)果選取了其中 4期的結(jié)果作為示例：

[1]De Jong,P.B.Zehneirth.Claims Reserving State Space Models and the Kalman Filter[J].Journal of the Institute of Actuaries，1983，（110）.

[2]Taylor,G.Ashe,F.Second Moments of Estimates of Outstanding Claims[J].Journal of Economics，1983，（23）.

[3]Taylor,G.Loss Reserving:An Actuarial Perspective[M].Boston:Kluwer Academic，2000.

[4]Taylor,G.Second-order Bayesian Revision of a Generalized Linear Model[J].Scandinavian Actuarial Journal，2008，（4）.

[5]Taylor,G.Gráinne McGuire,Adaptive Reserving Using Bayesian Revision for the Exponential Dispersion Family[J].Variance，2009，（3）.

[6]張連增.未決賠款準(zhǔn)備金評估的隨機性模型與方法[M].北京：中國金融出版社，2008.