基于Copula函數(shù)的中藥有效成分群譜效分析*

李智慧 陸 濤，3 楊中林 黃 園 陶禹希 言方榮△

基于Copula函數(shù)的中藥有效成分群譜效分析*

李智慧1，2陸 濤1，2，3楊中林1，4黃 園1，2陶禹希1，2言方榮1，2△

目的 通過構(gòu)建適合的二元Copula函數(shù)模型，對中藥有效成分群指紋圖譜與藥物活性進(jìn)行分析，探索簡單易行的中藥質(zhì)量控制方法。方法 本文首先分別模擬分析二元Copula函數(shù)和二元正態(tài)分布對中藥有效成分群譜效關(guān)系的擬合情況，然后選取合適的Copula函數(shù)模型分析18批次懷牛膝藥材指紋圖譜與成骨細(xì)胞增殖活性的相關(guān)關(guān)系，最后利用選擇的Copula函數(shù)模型對中藥懷牛膝進(jìn)行質(zhì)量判別。結(jié)果 五種常見的Copula函數(shù)中，二元單參數(shù)Gumbel Copula函數(shù)對所研究的中藥懷牛膝有效成分群譜效關(guān)系有較好的擬合效果，其相關(guān)系數(shù)τ=0.5547，ρs=0.7422，λup=0.6384，λlo=0。結(jié)果表明此函數(shù)上尾部呈現(xiàn)出較強的相關(guān)性而下尾部漸近獨立，即質(zhì)量好的中藥其細(xì)胞增值活性也較強，質(zhì)量差的對細(xì)胞藥效影響則不顯著。結(jié)論 Copula函數(shù)能夠用于分析中藥有效成分群譜效關(guān)系，尤其對尾部相關(guān)性的描述為質(zhì)量辨別提供了新方法，為中藥譜效關(guān)系研究提出新思路。

Copula 中藥 有效成分群 譜效分析 質(zhì)量控制

1.中國藥科大學(xué)天然藥物活性物質(zhì)與功能國家重點實驗室(210009)

2.中國藥科大學(xué)理學(xué)院數(shù)學(xué)教研室

3.中國藥科大學(xué)分子設(shè)計與藥物發(fā)現(xiàn)實驗室

4.中國藥科大學(xué)中藥學(xué)院

△通信作者::言方榮，Email:f.r.yan@163.com

隨著中藥現(xiàn)代化的發(fā)展，探索符合中醫(yī)藥特點的中藥質(zhì)量分析與評價技術(shù)，是推進(jìn)中藥現(xiàn)代化和國際化發(fā)展的關(guān)鍵任務(wù)之一。由于中藥是多個基本單元組成的系統(tǒng)，其臨床療效不是各組成部分療效簡單的線性加和〔1〕，檢測中藥一個或幾個化學(xué)成分含量難以準(zhǔn)確評價中藥質(zhì)量優(yōu)劣，以往的中藥質(zhì)量檢測方法實現(xiàn)不了根據(jù)“量”控制“質(zhì)”的目的。尋找與藥效相關(guān)的成分群，建立反映中藥內(nèi)在質(zhì)量的藥效指紋圖譜，通過中藥化學(xué)指紋圖譜來評價中藥產(chǎn)品質(zhì)量已成為中藥質(zhì)量分析領(lǐng)域的前沿研究方向〔2〕。中藥化學(xué)指紋圖譜是能夠標(biāo)識藥物內(nèi)多種化學(xué)成分特性的多維多息譜圖，它是對中藥物質(zhì)基礎(chǔ)的一種整體表達(dá)，能較好地體現(xiàn)中藥成分的復(fù)雜性和相關(guān)性〔3〕。因此可借助其找到與藥效相關(guān)的藥材的特征，進(jìn)而得到“譜-效”關(guān)系。

Copula函數(shù)是一類將聯(lián)合分布函數(shù)與他們各自的邊緣分布函數(shù)連接在一起的函數(shù)，也稱為連接函數(shù)。它可以捕捉變量間非線性、非對稱以及尾部相關(guān)關(guān)系，因此用Copula函數(shù)描述中藥“譜-效關(guān)系”復(fù)雜的非線性特點是適宜的。Copula函數(shù)模型不僅可以用于研究一般情況下變量之間的相關(guān)關(guān)系，還可以用于研究極值相關(guān)關(guān)系〔4〕，所以在中藥質(zhì)量優(yōu)劣的譜效分析中更為實用和有效。本文借助二元Copula函數(shù)建立中藥懷牛膝有效成分群指紋圖譜與成骨細(xì)胞增殖活性之間的譜效關(guān)系分析模型并進(jìn)行懷牛膝質(zhì)量優(yōu)劣判別。

Copula函數(shù)理論

Copula函數(shù)是把隨機向量x1，x2，…，xN的聯(lián)合分布函數(shù)F(x1，x2，…，xN)與各自的邊緣分布函數(shù)FX1(x1)，F(xiàn)X2(x2)，…，F(xiàn)XN(xN)相連接的連接函數(shù)，即函數(shù)C(u1，u2，…，uN)使F(x1，x2，…，xN)=C［FX1(x1)，F(xiàn)X2(x2)，…，F(xiàn)XN(xN)］。

1.二元Copula函數(shù)

定義1(Nelsen，2006) 二維 Copula函數(shù)是指滿足以下性質(zhì)的函數(shù)C(u，v):

(1)C(u，v)的定義域為［0，1］×［0，1］;

(2)至少存在一個u0∈［0，1］和一個v0∈［0，1］，使得C(u0，v)=0=C(u，v0);對任意 0≤u1≤u2≤1 和0≤v1≤v2≤1 有C(u2，v2)－C(u2，v1)－C(u1，v2)+C(u1，v1)≥0;并且對任意的u，v∈［0，1］，滿足C(u，1)=u，C(1，v)=v。

根據(jù)Sklar定理若二元聯(lián)合分布函數(shù)H(x，y)存在邊緣分布F(x)和G(y)，則存在一個 Copula函數(shù)C(u，v)，滿足H(x，y)=C［F(x)，G(y)］。連續(xù)邊緣函數(shù)確定唯一的Copula函數(shù)，當(dāng)為離散邊緣分布函數(shù)時有Copula函數(shù)使之成立。這不僅提供了在不研究邊緣分布的情況下分析變量之間相關(guān)結(jié)構(gòu)的途徑，同時也為求取聯(lián)合分布函數(shù)提供了方法。

2.基于Copula函數(shù)的相關(guān)性度量

變量間的相關(guān)關(guān)系通常用相關(guān)系數(shù)來度量。目前，有很多種相關(guān)系數(shù)，如線性相關(guān)系數(shù)、秩相關(guān)系數(shù)和尾部相關(guān)系數(shù)等。線性相關(guān)系數(shù)能反映變量間的線性關(guān)系;秩相關(guān)系數(shù)善于描述變量一致性;尾部相關(guān)系數(shù)適用于極值理論下相關(guān)性測度。Copula函數(shù)是一種更靈活穩(wěn)健的相關(guān)性分析工具，有其自身的特點:如果對變量進(jìn)行嚴(yán)格單調(diào)增變換，由Copula函數(shù)導(dǎo)出的相關(guān)性度量值不會改變，反映的是嚴(yán)格單調(diào)增變換下的相關(guān)性，比線性相關(guān)系數(shù)的適用范圍寬泛〔5〕。

定義2(Kendall秩相關(guān)系數(shù)) 令(x1，y1)和(x2，y2)為獨立同分布的隨機向量，定義

為Kendall秩相關(guān)系數(shù)，記為τ。

尾部相關(guān)系數(shù)是一個廣泛用于極值理論的測度，用來表示當(dāng)一個觀測變量的實現(xiàn)值為極值時另一個變量也為極值的概率。

定義3(尾部相關(guān)系數(shù)) 令X、Y為兩個連續(xù)的隨機變量，具有邊緣分布F(x)和G(y)，Copula函數(shù)C(u，v)，分別定義

為上尾相關(guān)系數(shù)和下尾相關(guān)系數(shù)。其中u=F(x)，v=G(y)〔6〕。

3.常見的二元Copula函數(shù)

(1)正態(tài)Copula函數(shù)

其中，ρ為變量間的線性相關(guān)系數(shù);Φ－1為標(biāo)準(zhǔn)正態(tài)分布函數(shù)的逆函數(shù)。

(2)t-Copula函數(shù)

其中，ρ為變量間的線性相關(guān)系數(shù);k為自由度;tk－1為自由度為k的一元t分布的分布函數(shù)的逆函數(shù)。

3.阿基米德Copula函數(shù)

其中，φ(u)稱為阿基米德 Copula函數(shù)C(u1，u2，…，uN)的生成元，φ－1(u)是 φ(u)的反函數(shù)，在區(qū)間［0，∞)上連續(xù)并且非增。

阿基米德族Copula函數(shù)眾多，常用的二元單參數(shù)Copula函數(shù)有Gumbel Copula、Clayton Copula和Frank Copula(見表1)。

4.常見Copula函數(shù)與相關(guān)性分析

Gaussian、t和Frank Copula函數(shù)具有對稱的尾部，其中t-Copula對隨機變量之間的尾部相關(guān)的變化較為敏感，能更好地捕捉到隨機變量之間的對稱的尾部相關(guān)關(guān)系;Gaussian和Frank Copula在分布的尾部兩個變量是漸進(jìn)獨立的。Gumbel Copula的密度函數(shù)具有非對稱性，其密度函數(shù)呈“J”字形(圖1)，即上尾高下尾低，對變量在分布上尾部的變化十分敏感，能夠快速捕捉到上尾相關(guān)的變化;而對變量在分布下尾部的變化不敏感，難以捕捉到下尾相關(guān)的變化。Clayton Copula與Gumbel Copula函數(shù)密度函數(shù)圖相反，其密度分布呈“L”字形，即上尾低下尾高，能反映出變量在下尾部的相關(guān)關(guān)系而無法描述在上尾部的變化情況。

表1 常用的二元單參數(shù)阿基米德Copula函數(shù)

模擬分析

目前，中藥質(zhì)量控制中整體系統(tǒng)的觀點已被接受，以單一或幾個成分來判斷藥材質(zhì)量具有一定的局限性，不能體現(xiàn)中藥多成分、多靶點的特征〔7〕;因此本文擬以18批中藥懷牛膝有效成分群指紋圖譜實驗數(shù)據(jù)和成骨細(xì)胞增值率為研究對象，借助二元Copula函數(shù)得到“譜-效”相關(guān)性，并獲得藥材質(zhì)量判斷的新方法。以往文獻(xiàn)常假定中藥譜效關(guān)系成多元正態(tài)分布，然而這樣的模型假設(shè)不一定適合復(fù)雜的中藥系統(tǒng)。為檢驗二元正態(tài)函數(shù)和Copula函數(shù)對懷牛膝的譜效關(guān)系的擬合情況，本研究中各產(chǎn)生1000組二元正態(tài)隨機數(shù)和5種滿足常見Copula函數(shù)形式的隨機數(shù)模擬分析。將有效成分群指紋圖譜值記作X，藥效值記作Y，用U和V表示X和Y的邊緣分布函數(shù)。模擬結(jié)果(表2)表明:Gumbel Copula函數(shù)關(guān)于U和V的數(shù)學(xué)期望計算結(jié)果E(U)Gum=0.4952，E(V)Gum=0.4857比二元正態(tài)分布E(U)Gau=0.7902，E(V)Gau=0.7269更接近實驗結(jié)果E(U)exp=0.4881，E(V)exp=0.5306，并且比其他Copula函數(shù)結(jié)果更優(yōu)。

表2 二元正態(tài)分布和常見Copula函數(shù)的模擬情況

實例分析

1.數(shù)據(jù)處理

數(shù)據(jù)來源于文獻(xiàn)〔7〕并采用其數(shù)據(jù)處理方法。將有效成分群指紋圖譜正、負(fù)系數(shù)相關(guān)峰峰面積比值記作AR(area rate)，成骨細(xì)胞增殖活性記作AT。

2.Copula函數(shù)模型選擇

Copula函數(shù)在實際應(yīng)用中的關(guān)鍵是函數(shù)形式的選擇，不同Copula函數(shù)模型可能導(dǎo)致不同的分析結(jié)果，因此選擇合適的Copula函數(shù)模型十分重要。本文中采用比較相關(guān)系數(shù)和解析法〔8〕結(jié)合的方式選擇合適的Copula函數(shù)。解析法借助于經(jīng)驗Copula函數(shù)，經(jīng)驗分布具有較好的統(tǒng)計性質(zhì)，為使用方法提供了保證。此外，經(jīng)驗分布還可以減少假設(shè)所帶來的誤差〔9〕。

對于任意Copula函數(shù)集合Ck，最優(yōu) Copula函數(shù)的選擇準(zhǔn)則是考慮它們與經(jīng)驗Copula函數(shù)CN(u，v)之間的平方歐式距離，則有

即Copula函數(shù)選擇的解析法〔10〕。

本研究在Matlab R2010a中算出AR和AT之間的線性相關(guān)系數(shù)ρ，Kendall秩相關(guān)系數(shù)τ，Spearman秩相關(guān)系數(shù)ρs，并估算出5個常見Copula函數(shù)的相關(guān)系數(shù)值。根據(jù)Copula函數(shù)選擇的解析法，計算出平方歐式距離并參考相關(guān)系數(shù)值，平方歐氏距離越小相關(guān)系數(shù)與真實值越接近，則相應(yīng)的Copula函數(shù)模型較為理想，能較好的反應(yīng)譜效相關(guān)關(guān)系;反之，則用于描述譜效關(guān)系效果不佳。

結(jié)果分析

根據(jù)樣本值估計AR和AT各自分布函數(shù)，并進(jìn)行假設(shè)檢驗，結(jié)果表明AR～N(0.8280，0.17772)，AT～N(236.7704，113.82032)。相關(guān)系數(shù)和平方歐式距離計算結(jié)果(見表2)表明:Gumbel Copula函數(shù)平方歐式距離=0.0351在五種Copula函數(shù)計算結(jié)果中較小并且Kendall秩相關(guān)系數(shù)估算值 τGum=0.5547與真實值τ=0.5686最為接近。觀察AR和AT的分布函數(shù)相關(guān)關(guān)系圖(圖2)，綜合分析以上結(jié)果，選取二元單參數(shù)Gumbel Copula函數(shù)作為AR和AT的聯(lián)合分布的擬合函數(shù)，參數(shù)α=2.2455，尾部相關(guān)系數(shù)λup=0.6384，λlo=0。這表明這18批次中藥懷牛膝正負(fù)相關(guān)峰面積比值較大者對其成骨細(xì)胞增值率影響較大，而比值較小的即使發(fā)生較大的改變對藥效也幾乎沒有什么影響。所以在進(jìn)行懷牛膝質(zhì)量控制的時候要特別注意對AR值較大的藥材的選擇，因為它們的變化能引起藥效較大變化。

圖1 Gumbel Copula密度函數(shù)(α=2．2455)

圖2 AR和AT分布函數(shù)關(guān)系圖

表3 經(jīng)驗Copula函數(shù)與常見Copula函數(shù)的平方歐式距離

討 論

本文引入尾部相關(guān)性度量指標(biāo)λup和λlo，得到了譜效關(guān)系在尾部的變化趨勢，并將其運用于中藥懷牛膝質(zhì)量優(yōu)劣的判別中。常用的相關(guān)系數(shù)實際上是線性變換下不變的一種相關(guān)性指標(biāo)，涉及到非線性函數(shù)的相關(guān)性，會導(dǎo)出錯誤的結(jié)論，而由Copula函數(shù)導(dǎo)出的相關(guān)性度量可以更準(zhǔn)確的描述出變量間的非線性相關(guān)關(guān)系，因此應(yīng)用范圍更廣。

表4 (u，v)和(1－u，1－v)計算及排序結(jié)果

表4 (u，v)和(1－u，1－v)計算及排序結(jié)果

C(u，v) 序號 C(1－u，1－v) 序號 AR排序 AT排序0.0157 13 0.9571 13 13 13 0.0323 7 0.8901 7 7 7 0.0747 12 0.7359 12 10 12 0.0787 8 0.7321 8 12 8 0.1504 10 0.4062 4 8 6 0.2739 11 0.3894 15 11 5 0.3282 16 0.3889 5 16 15 0.3363 6 0.3737 14 4 4 0.3444 4 0.3679 6 14 14 0.3549 5 0.2880 11 15 9 0.3842 14 0.2667 16 5 11 0.3865 15 0.2432 10 6 16 0.5610 9 0.2371 9 9 3 0.6422 3 0.1815 3 3 10 0.7153 18 0.1358 18 18 18 0.8008 17 0.0668 17 17 17 0.8241 2 0.0544 1 1 2 0.8519 1 0.0350 2 2 1

在本研究中Copula函數(shù)模型表現(xiàn)良好，但仍存在以下兩個問題有待進(jìn)一步研究:其一，本文中選擇了5種常見的Copula函數(shù)對中藥有效成分群譜效關(guān)系的尾部相關(guān)性進(jìn)行了實證研究，從中選擇出最優(yōu)的Copula作為擬合函數(shù)。但是由于Copula函數(shù)族的龐大和選擇方法的多樣性，并沒有論證所選擇的Copula是符合數(shù)據(jù)特征的最優(yōu)形式，采用何種方法選取最優(yōu)Copula函數(shù)模型有待進(jìn)一步探究。其二，中藥的多樣性導(dǎo)致不同藥材之間的差異顯著。一個Copula函數(shù)往往只能適用于某一個中藥，對其他則不一定合適，因此本文中選擇出的函數(shù)不能適用于所有中藥。

1．秦華珍，劉磊，王曉倩，等．中藥劑量與量效關(guān)系的思考．四川中醫(yī)．2011，(6):48-49．

2．李云飛，程翼宇，范驍輝．中藥多維譜效關(guān)系研究思路探討．中國天然藥物，2010，(3):167-170 ．

3．齊方，蓉蓉，薛付忠．中藥藥性特征標(biāo)記的PLS統(tǒng)計模式識別模型．中國衛(wèi)生統(tǒng)計，2011，(6):628-637．

4．韋艷華，張世英．Copula理論及其在金融分析上的應(yīng)用．北京:清華大學(xué)出版社，2008．

5．朱新玲．相關(guān)系數(shù)與Copula函數(shù)相關(guān)性比較研究．武漢科技大學(xué)學(xué)報，2009，32(6):664-668．

6．Nelson RB．An Introduction to Copulas．New York:Springer，1998:214-216．

7．周培培，言方榮，張春鳳，等．基于成分群動態(tài)變化探索藥材質(zhì)量優(yōu)劣判斷方法初步研究．中醫(yī)藥學(xué)報，2012，40(1):63-68．

8．于波，陳希鎮(zhèn)，杜江．Copula函數(shù)的選擇:方法與應(yīng)用．?dāng)?shù)理統(tǒng)計與管理，2008，27(6):1027-1033．

9．于波．Copula函數(shù)模型的選擇．統(tǒng)計與決策，2009(14):153-154．

10．閆寶偉，郭生練，肖義，等．基于兩變量聯(lián)合分布的干旱特征分析．干旱區(qū)研究，2007，24(4):538-541．

The Effective Components Analysis in Traditional Chinese Medicine Based on Copula

Li Zhihui，Lu Tao，Yang Zhonglin，et al．Department of Mathematics，China Pharmaceutical University(210009)，Nanjing

ObjectiveWe analyze the activity of the effective component group by constructing suitable binary Copula model，and then get a simple and convenient method to control the traditional Chinese medicine quality.Methods First，simulations on Copula models and the binary normal functions are taken respectively to fit the dose-response relationship.Second，Achyranthes bidentata BI.from 18 batches are chosen for study.The correlation analysis between the HPLC fingerprints of samples and proliferation activity of osteoblasts are carried out with suitable bivariate Copula function model.Results Gumbel Copula function is the most suitable model forAchyranthes bidentata BI.spectrum-response relationship.It has the correlation coefficients τ=0.5547，ρs=0.7422，λup=0.6384andλlo=0.It is a strong link in the upper tail，but asymptotic independence in the lower tail.That is to say that the good quality traditional Chinese medicine has stronger cell proliferation activity，however，the poor quality ones has insignificant effect.Conclusion Copula function can be used to analyze the relationship between spectrum and activity of traditional Chinese medicine active ingredients group，especially the tail dependence and drug activity.It provides a new method for TCM quality discrimination，and puts forward new ideas for TCM spectrum activity relationship research.

Copula;Traditional Chinese medicine(TCM);Effective components group;Spectrum activity relationship;Quality control

中央高校專項業(yè)務(wù)經(jīng)費(JKQ2011032，JKPZ2013015);國家自然科學(xué)基金重點項目(NSFC 81130068)

(責(zé)任編輯:劉 壯)