肖前軍，周金玉，鄧總綱

(湖南理工職業(yè)技術(shù)學(xué)院，湖南湘潭411104)

蛋白質(zhì)序列混沌游戲表示模擬效果的優(yōu)化

肖前軍，周金玉，鄧總綱

(湖南理工職業(yè)技術(shù)學(xué)院，湖南湘潭411104)

蛋白質(zhì)序列的可視化表示——混沌游戲表示呈現(xiàn)出明顯的分形特征.根據(jù)分形的產(chǎn)生機(jī)理，可以用遞歸迭代函數(shù)系統(tǒng)很好地模擬蛋白質(zhì)序列的混沌游戲表示.在實驗中發(fā)現(xiàn)，其模擬效果可以進(jìn)一步優(yōu)化，可能有助于更準(zhǔn)確地進(jìn)行物種親緣關(guān)系的分析.

遞歸迭代函數(shù)系統(tǒng)；混沌游戲表示；模擬；蛋白質(zhì)序列；優(yōu)化

0 引言

蛋白質(zhì)序列是由20種不同的氨基酸組成的一維字符序列，但要直接從中得出更多隱含的生物特性是非常困難的.為此人們設(shè)計出很多方法來進(jìn)行序列分析，例如，統(tǒng)計分析、人工神經(jīng)網(wǎng)絡(luò)技術(shù)和動態(tài)規(guī)劃等方法已經(jīng)成功地應(yīng)用到蛋白質(zhì)序列分析的某些領(lǐng)域，但這些方法不夠直觀且其中很多方法都需要進(jìn)行繁瑣的序列比對.因此，人們探索將蛋白質(zhì)序列轉(zhuǎn)化為數(shù)字信號、曲線或者可視化模型，再利用數(shù)學(xué)方法（如分形理論）進(jìn)行相關(guān)分析，以期從中發(fā)現(xiàn)更多的生物特性.Jeffrey[1]介紹了一種表征蛋白質(zhì)序列的可視化方法——混沌游戲表示，通過在蛋白質(zhì)序列的混沌游戲表示圖中發(fā)現(xiàn)明顯的分形特征[2-5]，而實際的分形是隨機(jī)迭代產(chǎn)生的，因此可以運用遞歸迭代函數(shù)系統(tǒng)[2，5-10]來模擬蛋白質(zhì)序列的混沌游戲表示.實驗發(fā)現(xiàn)模擬效果很好[2，7]，而且其結(jié)果可以進(jìn)行物種進(jìn)化等相關(guān)分析[6]，實驗過程中還發(fā)現(xiàn)其模擬結(jié)果可以進(jìn)一步優(yōu)化.

1 理論模型

1.1 蛋白質(zhì)序列混沌游戲表示（CGR）[1-5]

蛋白質(zhì)序列由20種不同的氨基酸組成[11]，即甘氨酸（G）、丙氨酸（A）、纈氨酸（V）、異亮氨酸（I）、亮氨酸（L）、苯丙氨酸（F）、脯氨酸（P）、甲硫氨酸（M）、色氨酸（W）、半胱氨酸（C）、絲氨酸（S）、蘇氨酸（T）、天冬酰胺（N）、谷氨酸（E）、酪氨酸（Y）、組氨酸（H）、天冬氨酸（D）、谷酰胺（Q）、賴氨酸（K）和精氨酸（R）.根據(jù)detailed-HP模型[5，7，12-13]，可以把以上20種氨基酸分成4類：無極性氨基酸（A、I、L、M、F、P、W和V），負(fù)極性氨基酸（D和E），正極性氨基酸（R、H和K）和不帶電氨基酸（N、C、Q、G、S、T和Y）.對一條長度為l的蛋白質(zhì)序列s=s1s2…sl，其中si是20種氨基酸中的某一種，i=1，2，…，l，定義：這樣就把蛋白質(zhì)的氨基酸序列轉(zhuǎn)換成了一條字符序列X（s）=a1a2…al，其中ai∈{0，1，2，3}.根據(jù)Jeffrey[1]的定義，一條蛋白質(zhì)序列s的混沌游戲表示規(guī)則如下：先在平面上作一正方形[0，1]×[0，1]，以其4個頂點分別代表4個字符0、1、2和3，然后進(jìn)行打點，第一點打在正方形中心與序列X（s）第一個字符對應(yīng)頂點的中點，第i（i≥2）點打在第i-1點與序列X（s）第i個字符對應(yīng)頂點的中點，直到在正方形中作出序列X（s）最后一個字符對應(yīng)的點，這樣得到的圖像稱為蛋白質(zhì)序列s的混沌游戲表示CGR.

1.2 混沌游戲表示測度[5，8，10，12]

根據(jù)蛋白質(zhì)序列混沌游戲表示的作法可知，序列X（s）與其CGR之間具有一一對應(yīng)關(guān)系，即不同的蛋白質(zhì)序列對應(yīng)著不同的混沌游戲表示.為了進(jìn)一步討論的需要，引入蛋白質(zhì)序列的混沌游戲表示的測度這一概念，對給定的蛋白質(zhì)序列，以如下方式來定義其CGR測度：μ（B）＝＃（B）/N，其中＃（B）表示該蛋白質(zhì)序列CGR子集B中點的數(shù)目，N表示蛋白質(zhì)序列的長度.為了便于討論，將正方形[0，1]×[0，1]分割成2k×2k個小正方形網(wǎng)格，對每一種分割，取B為該分割中每一網(wǎng)格與蛋白質(zhì)序列CGR的交集，則＃（B）就是該網(wǎng)格上CGR點的數(shù)目，這樣我們得到2k×2k的測度矩陣.

1.3 遞歸迭代函數(shù)系統(tǒng)（RIFS）[2，5-10]

在完備距離空間X上，一個由壓縮映射集S＝{S1，…，SN}和與其對應(yīng)的概率集P＝（pij）（其中構(gòu)成的系統(tǒng)（S，P），任取起始點A0，按如下的迭代過程產(chǎn)生隨機(jī)點列：An+1=Sσn（An）（n＝0，1，2，…），其中σn在集合{1，2，…，N}中取值且P（σn+1＝i）=pσni，i=1，…，N，則稱（S，P）為遞歸迭代函數(shù)系統(tǒng).

在本文中，我們?nèi)=[0，1]×[0,1]，N=4，壓縮映射為：1，2，3，4，即：

遞歸迭代函數(shù)系統(tǒng)的一個主要結(jié)論是：存在唯一不變測度且其支撐集為遞歸迭代函數(shù)系統(tǒng)的吸引子.由于需要用遞歸迭代函數(shù)系統(tǒng)的吸引子來逼近蛋白質(zhì)序列的混沌游戲表示，所以為了得到給定的蛋白質(zhì)序列對應(yīng)的RIFS的吸引子，我們根據(jù)已知蛋白質(zhì)序列混沌游戲表示的測度，用矩方法來估計遞歸迭代函數(shù)系統(tǒng)中的未知參數(shù)pij.如果μ′與A分別是RIFS的不變測度與吸引子，則μ′的矩為：

記Gmn為蛋白質(zhì)序列的混沌游戲表示測度μ的各階矩，通過解優(yōu)化問題Gmn）2可得參數(shù)pij的估計值，根據(jù)其估計值可以產(chǎn)生上述RIFS吸引子A，記χB為吸引子A的子集B的特征函數(shù)，根據(jù)RIFS遍歷理論可知其不變測度為：μ′（B）＝這樣就可以計算出與上述蛋白質(zhì)序列混沌游戲表示的測度矩陣同樣大小的不變測度矩陣，然后根據(jù)它們的差異程度來判別RIFS的吸引子逼近蛋白質(zhì)序列CGR的好壞.

1.4 CGR模擬效果判斷[2，5-6，10]

1.4.1 相對標(biāo)準(zhǔn)誤差

為了說明不變測度對原始測度的模擬效果，引入相對標(biāo)準(zhǔn)誤差：e=e1/e2，其中e1=為原始測度矩陣的元素，yij為模擬測度矩陣的元素，若e＜1，則說明擬合得比較好，如果e＜0.1，則擬合的效果更好[2].

1.4.2 測度累積行走

測度累積行走定義如下：其中fi是將測度矩陣的k×k網(wǎng)格所有行擴(kuò)展為一行后第i個網(wǎng)格的密度值，是網(wǎng)格上所有密度值的平均值.

2 結(jié)果與分析

數(shù)據(jù)庫Genebank（ftp://ncbi.nlm.nih.gov/genbank/genomes）中有大量的蛋白質(zhì)序列，從中隨機(jī)下載了20種細(xì)菌（見表1）的蛋白質(zhì)序列.以細(xì)菌Thevo（Thermoplasma volcanium，火山熱原體，廣古生菌（Euryarchaeota））為例，根據(jù)1.1節(jié)中的方法，圖1是Thevo蛋白質(zhì)序列的混沌游戲表示，該圖呈現(xiàn)出明顯的分形特征.根據(jù)1.2節(jié)中的方法，圖2是k=7時Thevo蛋白質(zhì)序列混沌游戲表示的測度圖.基于1.3節(jié)中的理論，通過計算得到參數(shù)pij的估計為:

根據(jù)概率矩陣P4×4得到了對應(yīng)的k=7時遞歸迭代函數(shù)系統(tǒng)的吸引子圖3及不變測度圖4.對于Thevo而言，根據(jù)1.4.1節(jié)的方法，計算得到e=0.130 8，根據(jù)1.4.2節(jié)的方法，圖5給出了k=7時Thevo的CGR測度行走表示及其模擬的比較圖，表明用遞歸迭代函數(shù)系統(tǒng)能很好地模擬Thevo的蛋白質(zhì)序列CGR.為了對模擬結(jié)果進(jìn)一步優(yōu)化，我們考慮兩個問題：第一，迭代次數(shù)是否影響迭代結(jié)果？第二，k的不同取值是否影響迭代效果？在文獻(xiàn)［2，5］中，都是選用k=7來進(jìn)行研究.實驗證實迭代次數(shù)并不影響迭代結(jié)果，而k的不同取值可以影響迭代效果.以Thevo為例，圖6～8分別是k=8時Thero的CGR測度，RIFS吸引子及RIFS模擬CGR的不變測度，由圖1和圖3、圖2和圖4的對比結(jié)果與圖1和圖7、圖6和圖8的對比結(jié)果可以直觀地看出，k=8時的模擬效果明顯優(yōu)于k=7時的模擬效果，圖5和圖9分別表示k=7和8時，Thevo的CGR測度行走表示及其模擬的比較圖，也很好地證實了上述結(jié)論.此時，參數(shù)pij的估計為：

而且從表1可以看出，對Thevo而言，k=8時有最好的模擬效果.之后隨機(jī)選取了20種不同細(xì)菌進(jìn)行實驗發(fā)現(xiàn)，當(dāng)k取不同的值時，其模擬效果是不同的，可以選取不同的k值使其模擬效果達(dá)到最優(yōu)（見表1）.其缺點就是過程相對麻煩，需要進(jìn)行反復(fù)的實驗去尋找適合的k值.

圖1 Thevo的CGR

圖2 k=7時Thevo的CGR測度

圖3 k=7時Thevo的RIFS吸引子

圖4 k=7時Thevo的RIFS模擬CGR不變測度

表1 不同k值時RIFS對20種細(xì)菌的CGR模擬結(jié)果

細(xì)菌名稱（簡寫）蛋白質(zhì)序列長度相對標(biāo)準(zhǔn)誤差k=7k=8k=9k=10k=11k=12Cloab1 123 6040.189 50.180 50.092 30.099 20.099 60.099 8 lo640 5460.052 70.037 70.039 10.041 10.043 50.044 8Borbu283 5960.200 80.167 80.152 70.101 40.108 10.115 3Brume950 1280.038 20.036 00.027 30.025 90.026 90.027 8Camje508 8370.080 20.083 80.025 70.026 00.026 70.027 2Caucr1 209 2280.025 90.024 20.015 80.016 20.017 20.018 1Chlmu321 0650.118 30.113 80.092 40.098 10.098 70.113 2ChlpnAYerpek1 273 9040.068 90.058 00.055 10.057 50.058 60.060 5Wigbr201 6721.103 51.116 90.867 70.932 00.961 20.964 0Thevo450 1070.130 80.050 30.053 60.055 20.057 20.063 5Mycpe400 8000.325 30.221 70.230 20.235 20.241 20.251 1Urepa228 2390.405 80.330 30.248 10.262 80.268 90.271 3Bachd1 188 4480.048 70.046 60.034 60.036 70.037 20.038 2Bacsu1 218 6160.068 60.041 80.042 30.043 30.045 10.048 1Bif 364 0780.104 40.098 20.061 90.064 60.074 50.086 7

3 結(jié)語

用分形方法中的遞歸迭代函數(shù)系統(tǒng)對基于detailed-HP模型的蛋白質(zhì)序列可視化方法——混沌游戲表示作了逼近，分析發(fā)現(xiàn)用遞歸迭代函數(shù)系統(tǒng)來逼近蛋白質(zhì)序列的混沌游戲表示時有好的效果，且其模擬結(jié)果可以進(jìn)一步優(yōu)化，為更準(zhǔn)確地進(jìn)行物種親緣關(guān)系的分析提供一種新的方法.

[1] Jeffrey H J.Chaos game representation of gene structure[J].Nucleic Acids Res，1990（18）：2 163-2 170.

[2] 肖前軍.遞歸迭代函數(shù)系統(tǒng)對基于detailed-HP模型的蛋白質(zhì)序列的混沌游戲表示的模擬[J].湘潭師范學(xué)院學(xué)報（自然科學(xué)版），2008，30（1）：14-17.

[3] 肯尼思·法爾科內(nèi).分形幾何中的技巧[M].曾文曲，王向陽，陸夷，等，譯.沈陽：東北大學(xué)出版社，1999.

[4] 肯尼思·法爾科內(nèi).分形幾何——數(shù)學(xué)基礎(chǔ)及其應(yīng)用[M].曾文曲，王向陽，陸夷，等，譯.沈陽：東北大學(xué)出版社，2001.

[5] 肖前軍.功能蛋白質(zhì)序列的混沌游戲表示及其遞歸迭代函數(shù)系統(tǒng)模擬[D].湘潭：湘潭大學(xué)數(shù)學(xué)與計算科學(xué)學(xué)院，2008.

[6] Lau K S，Anh V V，Yu Z G.Recognition of an organism from fragments of its complete genomes[J].Physical Review E，2002（66）：0319101-9.

[7] 喻祖國，Vo Anh，劉家成.迭代函數(shù)系統(tǒng)模型在生物序列分析中的應(yīng)用[J].湘潭大學(xué)（自然科學(xué)學(xué)報），2003（9）：131-139.

[8] Yu Z G,Anh V V,Wanliss J A，et al.Chaos game representation of the dst index and prediction of geomagnetic storm events[J].Chaos,Solitons and Fractals,2007（31）：736-746.

[9] Vrscay E R.Fractal geometry and analysis,chapter iterated function systems：theory，applicationsand the inverse problem[M].Dorsrecht：Kluwer，1991:405-468.

[10] Yu Z G,Shi L，Xiao Q J，et al.Simulation for chaos game representation of genomes by recurrent iterated function systems[J].J Biomedical Science and Engineering,2008（1）:44-51.

[11] 孫嘯，陸祖宏，謝建明，等.生物信息學(xué)基礎(chǔ)[M].北京：清華大學(xué)出版社,2005.

[12] Yu Z G,Anh V V,Lau K S.Chaos game representation of protein sequences based on the detailed HP model and their multifractal and correlation analyses[J].Journal of Theoretical Biology,2004（226）：341-348.

[13] Yu Z G,Anh V V,Lau K S.Fractal analysis of measure representation of large proteins based on the detailed HP model[J].Physica A:Statistical and Theoretical Physics,2004（337）:171-184.

Optimization of Simulation for Chaos Game Representation of Protein Sequence

XIAO Qian-jun，ZHOU Jin-yu，DENG Zong-gang

（Hunan Vocational College of Science and Technology,Xiangtan 411104,Hunan,China）

A fractal pattern is apparent in visual representation of a protein sequence：Chaos Game Representation（CGR）.According to the mechanism of fractal forming，Chaos Game Representation（CGR）of protein sequences can be simulated by Recurrent Iterated Function System（RIFS）.In experiments，a close approximation of CGR is obtained by the method.The results of simulation can be further optimized and the phylogenetic tree can be correctly worked out.

iterated function system；chaos game representation；simulation；protein sequence；optimization

O 29；O 811.4；Q 516

A

1001-4217（2010）01-0035-07

2009-06-24

肖前軍（1978-），男，湖南湘潭人，講師.研究方向：分形幾何與生物信息學(xué).E-mail：xqjxt@126.com

2008年湖南省教育廳項目（08C882）