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

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

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

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

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

Copula 中藥 有效成分群 譜效分析 質量控制

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

2.中國藥科大學理學院數學教研室

3.中國藥科大學分子設計與藥物發現實驗室

4.中國藥科大學中藥學院

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

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

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

Copula函數理論

Copula函數是把隨機向量x1，x2，…，xN的聯合分布函數F(x1，x2，…，xN)與各自的邊緣分布函數FX1(x1)，FX2(x2)，…，FXN(xN)相連接的連接函數，即函數C(u1，u2，…，uN)使F(x1，x2，…，xN)=C［FX1(x1)，FX2(x2)，…，FXN(xN)］。

1.二元Copula函數

定義1(Nelsen，2006) 二維 Copula函數是指滿足以下性質的函數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。

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

2.基于Copula函數的相關性度量

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

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

為Kendall秩相關系數，記為τ。

尾部相關系數是一個廣泛用于極值理論的測度，用來表示當一個觀測變量的實現值為極值時另一個變量也為極值的概率。

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

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

3.常見的二元Copula函數

(1)正態Copula函數

其中，ρ為變量間的線性相關系數;Φ－1為標準正態分布函數的逆函數。

(2)t-Copula函數

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

3.阿基米德Copula函數

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

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

4.常見Copula函數與相關性分析

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

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

模擬分析

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

表2 二元正態分布和常見Copula函數的模擬情況

實例分析

1.數據處理

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

2.Copula函數模型選擇

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

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

即Copula函數選擇的解析法〔10〕。

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

結果分析

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

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

圖2 AR和AT分布函數關系圖

表3 經驗Copula函數與常見Copula函數的平方歐式距離

討 論

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

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

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

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

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

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

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

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

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

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

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

8．于波，陳希鎮，杜江．Copula函數的選擇:方法與應用．數理統計與管理，2008，27(6):1027-1033．

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

10．閆寶偉，郭生練，肖義，等．基于兩變量聯合分布的干旱特征分析．干旱區研究，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

中央高校專項業務經費(JKQ2011032，JKPZ2013015);國家自然科學基金重點項目(NSFC 81130068)

(責任編輯:劉 壯)