Copula回歸模型與應(yīng)用

王 靜，劉新紅，吳 萍（.中國(guó)中醫(yī)科學(xué)院廣安門醫(yī)院，北京0005；2.中國(guó)人民大學(xué) 統(tǒng)計(jì)學(xué)院，北京00872；.北京石油化工學(xué)院 數(shù)理系，北京0267）

Copula回歸模型與應(yīng)用

王 靜1,2，劉新紅3，吳 萍1
（1.中國(guó)中醫(yī)科學(xué)院廣安門醫(yī)院，北京100053；2.中國(guó)人民大學(xué) 統(tǒng)計(jì)學(xué)院，北京100872；3.北京石油化工學(xué)院 數(shù)理系，北京102617）

文章探討Copula回歸模型應(yīng)用于具有相關(guān)關(guān)系的指標(biāo)時(shí)的優(yōu)勢(shì)，通過(guò)示例分析，揭示了用Copula回歸模型與普通回歸模型對(duì)數(shù)據(jù)擬合的不同，表明了當(dāng)不同指標(biāo)間存在相關(guān)關(guān)系時(shí)，Copula回歸模型能更加客觀準(zhǔn)確地反映數(shù)據(jù)背后的關(guān)系。

Copula函數(shù)；回歸模型；相關(guān)

1 Copula回歸模型

Copula是一種通過(guò)單個(gè)變量的邊緣分布構(gòu)造多個(gè)變量的聯(lián)合分布的一種數(shù)學(xué)方法，1959年，Sklar提出了Copula函數(shù)，將多元隨機(jī)變量的邊緣分布和它們之間的相關(guān)結(jié)構(gòu)分開(kāi)研究，相關(guān)結(jié)構(gòu)不受邊緣分布的限制。隨后，很多學(xué)者發(fā)現(xiàn)了Copula理論在研究相關(guān)性方面的價(jià)值。Copula不僅可以反映線性相關(guān)，也能描述非線性相關(guān)。而最常用的Pearson相關(guān)系數(shù)只能反映變量間的線性相關(guān)程度，無(wú)法捕捉到非線性的關(guān)系。其他常用的一些相關(guān)性系數(shù)如Kendall's τ系數(shù)、Spearman系數(shù)和Gini關(guān)聯(lián)系數(shù)等雖然能在一定程度上反映變量間的非線性相關(guān)性，但都不能全面完整地刻劃變量間的相關(guān)結(jié)構(gòu)。

根據(jù)Sklar定理[1]，若二維隨機(jī)向量(X，Y)的聯(lián)合分布函數(shù)為H(x,y)，邊緣分布函數(shù)分別為F(x)和G(y)，則存在一個(gè)Copula函數(shù)C(u,v)，滿足：

如果F和G是連續(xù)的，則Copula函數(shù)C是唯一確定的。

在實(shí)際應(yīng)用中，常用的二元Copula函數(shù)有橢圓Copula函數(shù)和阿基米德Copula，橢圓Copula主要包括Gaussian Copula與t Copula，而阿基米德Copula除了包括常用的Gumbel Copula，Clayton Copula，F(xiàn)rank Copula，還包括Joe Copula、BB1 Copula等。這里只列出Gaussian Copula、Clayton Copula、Gumbel Copula和Frank Copula的具體形式，其他Copula可參考相關(guān)文獻(xiàn)。

（1）Gauss Copula函數(shù)

其中Φ(·)為標(biāo)準(zhǔn)正態(tài)分布的分布函數(shù)。

（2）Clayton Copula函數(shù)

當(dāng)θ→0時(shí)，隨機(jī)變量獨(dú)立；當(dāng)θ→+∞時(shí)，隨機(jī)變量完全相關(guān)。Clayton Copula的密度函數(shù)也具有非對(duì)稱性，其密度分布呈“L”字型，即上尾低下尾高。

（3）Gumbel Copula函數(shù)

當(dāng)θ=1時(shí)，隨機(jī)變量獨(dú)立；當(dāng)θ→+∞時(shí)，隨機(jī)變量完全相關(guān)。Gumbel Copula的密度函數(shù)具有非對(duì)稱性，其密度分布呈“J”字型，即上尾高下尾低。

（4）Frank Copula函數(shù)

當(dāng)θ→0時(shí)，隨機(jī)變量獨(dú)立；Frank Copula的密度函數(shù)具有對(duì)稱性，其密度分布呈“U”字型。

對(duì)于兩個(gè)存在相關(guān)關(guān)系的隨機(jī)變量X和Y，在各自回歸模型的基礎(chǔ)上，將這種相關(guān)關(guān)系用Copula函數(shù)刻畫，建立兩個(gè)因變量的聯(lián)合分布函數(shù)，即可建立Copula回歸模型：

其中，式（2）表示了隨機(jī)變量X的分布及回歸模型，式（3）表示了隨機(jī)變量Y的分布及回歸模型，式（4）表示了隨機(jī)變量X與Y的聯(lián)合分布。Peter X.-K.Song等[2]在2009年用Gaussian Copula對(duì)人體燒傷面積與存活率這兩個(gè)相關(guān)的結(jié)局指標(biāo)進(jìn)行了聯(lián)合回歸分析，論證了與分別單獨(dú)回歸分析相比，聯(lián)合回歸分析由于考慮了指標(biāo)間的相關(guān)性，并且能得到一個(gè)結(jié)局指標(biāo)關(guān)于另一個(gè)結(jié)局指標(biāo)的條件分布，從而具有更高的估計(jì)效能和推斷效能[4]。

例1 (2018年四川達(dá)州)如圖1，二次函數(shù)y=ax2+bx+c的圖象與x軸交于點(diǎn)A(-1，0)，與y軸的交點(diǎn)B在(0，2)與(0，3)之間(不包括這兩點(diǎn))，對(duì)稱軸為直線x=2．下列結(jié)論：①abc<0；②9a+3b+c>0；③若點(diǎn)點(diǎn)是函數(shù)圖象上的兩點(diǎn)，則y1

2 應(yīng)用示例

2．1 數(shù)據(jù)來(lái)源及變量選擇

數(shù)據(jù)來(lái)源于某藥物療效評(píng)價(jià)試驗(yàn)，樣本量為402例，因變量為西醫(yī)量表評(píng)分差值（變量名為y1）和中醫(yī)證候評(píng)分差值（變量名為y2）。自變量有藥物（drug）、醫(yī)院（cn）、年齡（age）、病程（bch）和基線評(píng)分（c0及zz0），除基線評(píng)分為連續(xù)變量外，其他均為分類變量。

2．2 建立普通回歸模型

對(duì)因變量 y1和 y2進(jìn)行正態(tài)分布檢驗(yàn)，Kolmogorov-Smirnov檢驗(yàn)的P值分別為0.001和0.000，均不服從正態(tài)分布。它們的偏度系數(shù)分別為0.0557和0.2245，具有明顯的右偏特點(diǎn)，故選取逆Gumbel分布（Reverse Gumbel），這是一種極值分布，概率密度函數(shù)為：

其數(shù)學(xué)期望和方差為：μ+0.5772σ和1.6449σ2。

西醫(yī)量表評(píng)分差值（y1）的回歸模型：

全模型：

中醫(yī)證候評(píng)分差值（y2）的回歸模型：

全模型：

經(jīng)過(guò)模型選擇，剔除不顯著自變量（P＞0.05），最后選定模型為：

兩個(gè)回歸模型參數(shù)估計(jì)結(jié)果如表1所示，逆Gumbel分布的參數(shù)σ估計(jì)值分別為5.0542和4.7871。自變量drug在兩個(gè)模型中都不顯著，也就是說(shuō)，兩種藥物在西醫(yī)量表評(píng)分和中醫(yī)證候評(píng)分改善方面，均沒(méi)有顯著差異。殘差服從正態(tài)分布，可見(jiàn)，選取的分布是合適的。

表1 μ1和μ2普通回歸模型參數(shù)估計(jì)結(jié)果

2．3 建立Copula回歸模型

在不考慮自變量的情況下，y1和y2的Pearson相關(guān)系數(shù)為0.6340，Kendall's τ相關(guān)系數(shù)為0.4818，表明這兩個(gè)因變量是相關(guān)的。將這種相關(guān)關(guān)系用Copula函數(shù)刻畫，建立兩個(gè)因變量的聯(lián)合分布函數(shù)，即可建立Copula回歸模型。

其中，H(y1，y2)表示(y1，y2)的聯(lián)合分布函數(shù)，F(xiàn)(y1)和G(y2)分別為y1和y2的邊緣分布函數(shù)，y1，y2均服從逆Gumbel分布，密度函數(shù)見(jiàn)式（5），分布函數(shù)易知，C(u,v)表示Copula函數(shù)。

C(u,v)有多種選擇，使用常用的Gauss Copula、Gumbel Copula、Clayton Copula和Frank Copula函數(shù)。參數(shù)估計(jì)采用極大似然法，通過(guò)R軟件的GAMLSS包和CDVINE包[3,4]實(shí)現(xiàn)。使用4種Copula函數(shù)的回歸模型的AIC值分別為2329.911、2376.840、2310.876和2318.165，可見(jiàn)，Clayton Copula回歸模型擬合效果最好。

在Clayton Copula回歸模型中，Copula函數(shù)中的參數(shù)估計(jì)值為1.3534，Kendall's τ相關(guān)系數(shù)為0.4036。 y1和y2逆Gumbel分布中的參數(shù)σ估計(jì)值分別為5.0005和5.0824，模型參數(shù)的具體估計(jì)結(jié)果見(jiàn)表2所示。

表2 μ1和μ2的Clayton Copula回歸模型參數(shù)估計(jì)結(jié)果

將Copula回歸模型估計(jì)結(jié)果與普通回歸模型參數(shù)估計(jì)結(jié)果比較，可以看出自變量回歸系數(shù)的點(diǎn)估計(jì)值相差不多，但是，中醫(yī)證候評(píng)分差值回歸模型中，自變量藥物（drug）和年齡（age）由不顯著因素變成了顯著因素，也就是說(shuō)，兩種藥物在中醫(yī)證候評(píng)分改善方面有顯著性差異，這也印證了Peter X.-K.Song等在2009年得出的聯(lián)合回歸比單獨(dú)普通回歸具有更高統(tǒng)計(jì)推斷效能的結(jié)論。

3 討論

Copula是分析相關(guān)關(guān)系的有力工具，近十年來(lái)在金融、保險(xiǎn)、生物和醫(yī)藥等領(lǐng)域得到廣泛的應(yīng)用。Copula回歸模型是一種考慮了因變量之間相關(guān)關(guān)系的聯(lián)合回歸模型[5]，當(dāng)因變量之間存在相關(guān)關(guān)系時(shí)，聯(lián)合回歸分析比普通回歸分析具有更高的推斷效能，所以，用普通回歸模型分析時(shí)不顯著的因素，用Copula回歸模型分析時(shí)可能變成顯著因素?？梢?jiàn)，當(dāng)指標(biāo)間存在相關(guān)性時(shí)，Copula回歸模型能夠揭示普通回歸模型所揭示不了的差異，能夠更充分挖掘數(shù)據(jù)背后的關(guān)系，從而更加客觀準(zhǔn)確地揭示隱藏在數(shù)據(jù)背后的規(guī)律。

[1]Nelsen R B.An Introduction to Copulas[M].New York:Springer,2006.

[2]Song P,Li M,Yuan Y.Joint Regression Analysis of Correlated Data Us?ing Gaussian Copulas[J].Biometrics,2009,（65）.

[3]Stasinopoulos M,Rigby B.Generalized Additive Models for Location Scale and Shape(GAMLSS)in R[J].Journal of Statistical Software, 2007,23(7).

[4]Brechmann E,Schepsmeier U.Modeling Dependence With C-and D-Vine Copulas:The R Package CDVINE[J].Journal of Statistical Software,2013,52(3).

[5]Kolev N,Paiva D.Copula-based Regression Models:A Survey[J].Jour?nal of Statistical Planning and Inference,2009,(139).

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

O212

A

1002-6487（2016）24-0079-03

科技部重大新藥創(chuàng)制課題(2013ZX09303301)；中國(guó)中醫(yī)科學(xué)院廣安門醫(yī)院所級(jí)科研基金課題（2011S264）

王 靜（1978—），女，河北唐山人，助理研究員，博士研究生，研究方向：數(shù)理統(tǒng)計(jì)。

劉新紅（1978—），女，河北保定人，博士，講師，研究方向：風(fēng)險(xiǎn)管理與非壽險(xiǎn)定價(jià)。

（通訊作者）吳 萍（1960—），女，江西萍鄉(xiāng)人，碩士，研究員，研究方向：臨床試驗(yàn)管理。