趙立戎

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

0 引言

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

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

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

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

1 廣義線性濾波模型

1.1 貝葉斯GLM的二階近似法

GLM 的標準結構為：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階設計矩陣，h為聯(lián)結函數(shù)，ε是一個隨機誤差向量。可見，Yi具有相同的分布形式，但位置參數(shù)θi和尺度參數(shù) 1/λi可以取不同的值。 E（Y）=μ（θ）=h-1（Xβ）。

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

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

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

考慮參數(shù)向量β隨時間隨機變化，則GLM模型形式變?yōu)椋篩（t）=h-1（X（t）β（t））+v（t），其中 Y（t）是 Y 在 t時刻的觀測數(shù)據(jù)向量，設計矩陣 X（t）、參數(shù)向量 β（t）和誤差向量 v（t）都表示t時刻的相應狀態(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ù)向量的維度隨時間變化。這樣就構成了動態(tài)廣義線性模型（DGLM）的基本形式。在實際應用中，t可以選為事故發(fā)生年，也可以選為日歷年。

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

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

為了估計 E［h-1（M（t|t）β（t）|Y（t））］和 E［h-1（M（t+1|t）β（t+1）|Y（t）］的值及對應的方差，在估計過程中應用貝葉斯GLM二階近似的結論，隨著觀測數(shù)據(jù)向量 y（1）,y（2），…逐漸引入模型，對β（?）的均值及方差的估計值逐步遞推調整，即產生一個濾波，濾波的起始點是對 E［β（1）］和 Var［β（1）］的先驗估計值。在應用廣義線性濾波模型評估非壽險準備金時，初始值的選定依賴于對經驗數(shù)據(jù)等先驗信息的分析和精算師的主觀判斷，隨著迭代步驟的增加，初始值對整個計算結果的敏感性將越來越低。在濾波迭代的過程中，每一步迭代都是從 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通常并不屬于同一分布族，只有在迭代過程中應用貝葉斯二階近似法，用同屬EDF自然共軛分布族的二階近似式代替原后驗似然，貝葉斯GLM濾波迭代才能夠實現(xiàn)遞推。DGLM的濾波迭代的邏輯結構如圖1所示。

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

2 伽瑪廣義線性濾波模型

在GLM中假設誤差項服從伽瑪分布，使用對數(shù)聯(lián)結函數(shù)（伴隨典則聯(lián)結），即 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）

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

Θ（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ù)結構隨時間變化的規(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）］，用當前觀測值與前期預測值之差 Y（t）-（t|t-1）作為狀態(tài)參數(shù)更新的基礎。

近似計算：Γ（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），作為下一輪迭代的基礎。

根據(jù)近似公式：

3 實例分析

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

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

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

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

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

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

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

[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]張連增.未決賠款準備金評估的隨機性模型與方法[M].北京：中國金融出版社，2008.