劉 娟 高 潔

(江南大學(xué)理學(xué)院，無錫 214122)(2010年4月16日收到;2010年8月4日收到修改稿)

甲型流感病毒DNA序列的長記憶ARFIMA模型*

劉 娟 高 潔

(江南大學(xué)理學(xué)院，無錫 214122)(2010年4月16日收到;2010年8月4日收到修改稿)

流感病毒分為三類:甲型(A型)，乙型(B型)，丙型(C型).在這三種類型中甲型(A型)流感病毒是最致命的流感病毒，對(duì)人類引起了嚴(yán)重疾病.本文對(duì)甲型流感病毒DNA序列建立了一種新的時(shí)間序列模型，即CGR(Chaos Game Representation)弧度序列.利用CGR坐標(biāo)將甲流病毒DNA序列轉(zhuǎn)換成CGR弧度序列，且引入長記憶ARFIMA模型去擬合此類序列，發(fā)現(xiàn)隨機(jī)找來的10條H1N1序列，10條H3N2序列都具有長相關(guān)性且擬合很好，并且還發(fā)現(xiàn)這兩種序列可以嘗試用不同的 ARFIMA模型去識(shí)別，其中 H1N1可用 ARFIMA(0，d，5)模型去識(shí)別，H3N2可用ARFIMA(1，d，1)模型去識(shí)別.

甲型流感，時(shí)間序列模型，CGR，ARFIMA(p，d，q)模型

PACS:87.10.Vg，02.50.Fz

1.引 言

流感是一種反復(fù)出現(xiàn)的傳染病，在全球引起了高發(fā)病率和高死亡率[]1.流感病毒分為三類:甲型(A型)，乙型(B型)，丙型(C型).在這三種類型中甲型(A型)流感(以下簡(jiǎn)稱甲流)病毒是最致命的流感病毒，給人類帶來了嚴(yán)重的疾病.甲流病毒根據(jù)其表面的血凝素(hemagglutinin，HA)和神經(jīng)氨酸酶(neuraminidase，NA)基因的不同又可分成16個(gè)HA亞型(H1-H16)和9個(gè)NA亞型(N1—N9)，不同的HA和NA形成了甲流病毒的許多亞型，如H1N1，H3N2，H5N1 等等［2—4］.筆者參看了許多文獻(xiàn)幾乎沒有看到用時(shí)間序列模型來挖掘甲型流感病毒的特性的，因而本文采用時(shí)間序列模型來分析甲型流感病毒.

1992年，Peng等［5］提出了 DNA一維游走模型.同年Voss等［6］提出了不同的觀點(diǎn)，他們發(fā)現(xiàn) DNA序列的譜密度顯示的1/fβ噪聲無處不在，意味著當(dāng)0<β<1存在長相關(guān)性，認(rèn)為不僅在非編碼區(qū)序列中在編碼序列中也存在長相關(guān)性.另一方面Buldyrev 等［7，8］設(shè)計(jì)了一個(gè)廣義-Lévy 游走模型去生成一個(gè)模型序列，使用所有可用的DNA序列發(fā)現(xiàn)主要在非編碼序列中呈現(xiàn)長相關(guān)性.基于該模型，Tai等［9］提出了一個(gè)二維修正-Lévy游走模型.為區(qū)分 C和T，A和 G，Lou等建立了二維和三維游走模型［10］，Yu 等則建立了圖譜［11，12］，來研究 DNA 序列的相關(guān)性.2006 年，Lopes和 Nunes［13］引入長記憶ARFIMA(0，d，0)模型去擬合 DNA序列的一維游走序列.2009年，Gao等［14］基于 CGR(chaos game representation)坐標(biāo)提出了一種DNA序列轉(zhuǎn)換成一個(gè)時(shí)間序列(CGR-游走序列)的方法，并引入長記憶ARFIMA(p，d，q)模型來分析.

本文對(duì)甲流病毒DNA序列提供了一種新的時(shí)間序列模型，即CGR弧度序列.利用CGR坐標(biāo)將甲流病毒DNA序列轉(zhuǎn)換成CGR弧度序列，且引入長記憶ARFIMA模型去擬合此類序列，發(fā)現(xiàn)隨機(jī)找來的10條H1N1序列，10條H3N2序列都具有長相關(guān)性且擬合很好，并且還發(fā)現(xiàn)這兩種序列可以嘗試用不同的 ARFIMA模型去識(shí)別，其中 H1N1可用ARFIMA(0，d，5)模型去識(shí)別，H3N2可用 ARFIMA(1，d，1)模型去識(shí)別.

2.ARFIMA模型

如果隨機(jī)過程 {xt}是平穩(wěn)的，且滿足方程Φ(B)Δdxt= Θ(B)εt，其中，－ 0.5

則稱{xt}服從 －0.5

因此，{xt}可看作是分?jǐn)?shù)差分噪聲導(dǎo)出的 ARMA(p，q)過程.當(dāng)2d－1=－1時(shí)，d=0即為短記憶過程;所以當(dāng)2d－1>－1時(shí)，d∈(0，0.5)具有長記憶的特征.

3.基于CGR的時(shí)間序列模型

1990年Jeffrey提出了一種 DNA序列可視化的方法即 CGR方法［15］.CGR是一種迭代映射技術(shù)，它把序列中的每個(gè)單元，如蛋白質(zhì)序列中氨基酸，DNA中的核苷酸，映射到一個(gè)連續(xù)的坐標(biāo)空間中去.

正方形的四個(gè)頂點(diǎn)對(duì)應(yīng)四種核苷酸.在這里，用DNA序列代替隨機(jī)數(shù)，每一個(gè)堿基的坐標(biāo)都可以來確定下一個(gè)堿基的位置.我們?nèi)?A(0，0)，T(1，0)，G(1，1)，C(0，1)，并且取點(diǎn)(0.5，0.5)為起始點(diǎn).

下面給出DNA迭代函數(shù)，也可以認(rèn)為是 CGR算法 的 公 式 化 形 式［15，16］. 對(duì) 于 一 個(gè) 序 列 S =

其中 gi={(0，0)，(1，0)，(1，1)，(0，1)}，gi和 si相對(duì)應(yīng).

對(duì)于一個(gè)DNA序列，定義

其中yn是CGRn的y坐標(biāo)值，xn是 CGRn的 x坐標(biāo)值.則得到一個(gè)數(shù)據(jù)序列 {Rn:n=1，2，…，N}，我們把它作為一個(gè)時(shí)間序列，并稱它為“CGR弧度序列”.

以甲流病毒 H1N1序列 CY056890為例，數(shù)據(jù)來自 NCBI網(wǎng)站，其網(wǎng)址:http://www.ncbi.nlm.nih.gov/.

它的CGR弧度序列“游走圖”如下表1.

表1 CY056890序列所選部分前8個(gè)ACATGGTA游走結(jié)果

4.甲流 H1N1型病毒 CY056890的數(shù)據(jù)分析

圖1(a)是CY056890序列CGR弧度序列圖(位置380—2170)，樣本容量1791.這些數(shù)據(jù)變動(dòng)較大，呈現(xiàn)非平穩(wěn)特征.考慮對(duì)此過程作d階差分.先對(duì)原序列作對(duì)數(shù)變換然后再做一階差分結(jié)果如圖1(b)所示，可見除少數(shù)地方呈現(xiàn)異方差外，基本呈現(xiàn)平穩(wěn)性.

圖2(a)(ACF)和圖2(b)(PACF)為樣本取對(duì)數(shù)再一階差分后的自相關(guān)函數(shù)圖形和偏自相關(guān)函數(shù)圖形.可見ACF衰減迅速，而PACF衰減緩慢，這意味著原序列具有長記憶特征.

圖1 (a)甲流H1N1型病毒CY056890的弧度序列圖;(b)取對(duì)數(shù)再一階差分圖

圖2 (a)取對(duì)數(shù)再一階差分的樣本自相關(guān)圖;(b)取對(duì)數(shù)再一階差分的樣本偏自相關(guān)圖

圖3給出了方差圖［17］是一個(gè)估計(jì)長記憶參數(shù)d的有用工具.對(duì)于一個(gè)長記憶時(shí)間序列{Rn}，它的均值珔Rk的方差滿足作log［Var(珔Rk)］關(guān)于log(k)的散點(diǎn)圖，對(duì)散點(diǎn)圖線性擬合，可估計(jì)得到線性方程的斜率為 －0.6877，令該斜率為2d－1= －0.6877，即可得 d的估計(jì)值0.156.

根據(jù)上述理由我們可選擇CGR弧度序列顯示長記憶特征.目的是利用上述特點(diǎn)為序列建立一個(gè)合適的模型.因此，可以考慮長記憶 ARFIMA(p，d，q)模型(d∈ (0，0.5))，p，q定階時(shí)為考慮實(shí)用性，僅考慮 p，q均小于等于5的 ARFIMA(p，d，q)模型.由 Akaike 信 息 判 別 準(zhǔn) 則［18，19］，可 選 ARFIMA(0，0.156，5)模型來擬合.

為檢驗(yàn)該模型的合理性，選擇了一個(gè)合適的檢驗(yàn)統(tǒng)計(jì)量 LB 檢驗(yàn)統(tǒng)計(jì)量［20，21］

其中rk是滯后k的樣本自相關(guān)函數(shù)，n是樣本容量，M是一個(gè)取定的比n小的正整數(shù).

表2顯示了對(duì)于各滯后階數(shù)，LB統(tǒng)計(jì)量的p值均顯著大于0.1，意味著擬合模型的殘差序列應(yīng)為白噪聲(純隨機(jī))，因而可以認(rèn)為 ARFIMA(0，0.156，5)模型能很合理地?cái)M合 CY056890序列的CGR-游走序列.

圖3 CY056890DNA序列的CGR弧度序列方差圖

表2 殘差的自相關(guān)檢驗(yàn)

表3給出了被選擇的 ARFIMA(0，0.156，5)模型的參數(shù)估計(jì)，5個(gè)參數(shù)的 T檢驗(yàn)統(tǒng)計(jì)量的p值均顯著小于 0.005.這意味著 ARFIMA(0，0.156，5)模型能有效地?cái)M合這個(gè)CGR弧度序列.

表3 條件最小二乘估計(jì)

5.其余9條H1N1序列和10條H3N2序列數(shù)據(jù)分析

表4和表5分別給出了隨機(jī)選的9條H1N1序列和10條 H3N2序列的數(shù)據(jù)信息、被選擇的ARFIMA(p，d，q)模型及參數(shù)估計(jì).從計(jì)算結(jié)果可得d均位于(0，0.5);對(duì)于各滯后階數(shù)，LB統(tǒng)計(jì)量的 p值除極個(gè)別外其余均顯著大于0.1;且每個(gè)被選擇的模型中各參數(shù)的T檢驗(yàn)統(tǒng)計(jì)量的p值均顯著小于0.01.所有這些結(jié)果都顯示ARFIMA(p，d，q)模型能很合理很有效地?cái)M合這些不同的CGR弧度序列且還發(fā)現(xiàn)所選H1N1序列均為ARFIMA(0，d，5)模型，所選 H3N2序列均為 ARFIMA(1，d，1)模型.所以我們可嘗試用 ARFIMA(0，d，5)模型，ARFIMA(1，d，1)模型分別去識(shí)別H1N1序列，H3N2序列.

表4 9條H1N1序列的數(shù)據(jù)信息、被選擇的ARFIMA模型和參數(shù)估計(jì)

表5 10條H3N2序列的數(shù)據(jù)信息、被選擇的ARFIMA模型和參數(shù)估計(jì)

6.結(jié) 論

本文基于CGR坐標(biāo)提出了一種將甲流病毒DNA序列轉(zhuǎn)換成時(shí)間序列(CGR弧度序列)的方法，并引入長記憶模型ARFIMA模型來分析，首先分析了甲流H1N1型病毒CY056890序列，從圖1到圖3可知弧度序列顯示長記憶特征，并選擇了ARFIMA(0，0.156，5)模型去擬合它，從表 2到表 3發(fā)現(xiàn)擬合合理有效.

然后又分析了隨機(jī)找來的19條序列的CGR弧度序列，從表4和表5可知所有 ARFIMA(p，d，q)模型都有效合理.并且從表4中還發(fā)現(xiàn)所選H1N1序列均為 ARFIMA(0，d，5)模型，表5所選 H3N2序列均為 ARFIMA(1，d，1)模型.

由此可見，DNA序列的CGR弧度序列能由長記憶ARFIMA(p，d，q)模型有效合理地?cái)M合，并且還可嘗試用 ARFIMA(0，d，5)模型，ARFIMA(1，d，1)模型分別去識(shí)別H1N1序列、H3N2序列.作為具有完善算法的經(jīng)典時(shí)間序列模型，不僅可以幫助我們得到甲流病毒DNA序列清晰的結(jié)構(gòu)，而且還可幫助我們有效識(shí)別甲流中的兩種亞型.

本文僅對(duì)甲流中的兩種亞型進(jìn)行了研究分析，后面我們將研究分析甲流中的其他亞型以及乙型丙型流感病毒.

［1］ Morens D，F(xiàn)olkers G，F(xiàn)auci A 2004 Nature 430 242

［2］ Chen J M，Sun Y X，Liu S 2009 Chinese Science Bulletin 54 1657(in Chinese)［陳繼明、孫映雪、劉 朔 2009科學(xué)通報(bào)54 1657］

［3］ Webster R G，Bean W J，Gorman O T 1992 Microbiol Rev.56 152

［4］ Shi X M，Shi L，Zhang J F 2010 Chin.Phys.B 19 038701

［5］ Peng C K，Buldyrev S，Goldberg A L，Havlin S，Sciortino F，Simons M，Stanley H E 1992 Nature 356 168

［6］ Voss R F 1992 Phys.Rev.Lett.68 3805

［7］ Buldyrev S V，Goldberger A L，Havlin S，Peng C K，Simon M，Stanley H E 1993 Phys.Rev.E 47 4514

［8］ Buldyrev S V，Goldberger A L，Havlin S，Mantegna R N，Matsa ME，Peng C K，Simon M，Stanley H E 1995 Phys.Rev.E 51 5084

［9］ Tai Y Y，Li P C，Tseng H C 2006 Physica A 369 688

［10］ Luo L F，Lee W J，Jia L J，Ji F M，Tsai L 1998 Phys.Rev.E 58 861

［11］ Yu Z G，Chen G Y 2000 Theor.Phys.33 673

［12］ Yu Z G，Anh V，Gong Z M，Long S C 2002 Chin.Phys.11 1313

［13］ Lopes S R C，Nunes M A 2006 Physica A 361 569

［14］ Gao J，Xu Z Y 2009 Chin.Phys.B 18 370

［15］ Jeffrey H J 1990 Nucleic Acid Res 18 2163

［16］ AlmeidaJonas， carrico Joao A， Maretzek António 2001 Bioinformatics 17 429

［17］ Beran J 1994 Statistics for long-memory Processes(New Work:Chapman Hall)

［18］ Hosking J R M 1984 Water Resour.Res.20 1898

［19］ Crato N，Ray B K 1996 Journal of Forcasting 15 107

［20］ Ljung G M，Box G E P 1978 Biometrika 65 297

［21］ Li W K，Mcleod A I 1986 Biometrika 73 217

Long-memory ARFIMA model for DNA sequences of influenza A virus*

Liu Juan Gao Jie
(School of Science，Jiangnan University，Wuxi 214122，China)(Received 16 April 2010;revised manuscript received 4 August 2010)

Influenza viruses are divided into three types:A，B and C.Among them，type A virus is the most virulent human pathogen and causes the most severe disease.In this paper，we propose a new time series model for influenza A virus DNA sequence，i.e.chaos game representation(CGR)radians series.The CGR coordinates are converted into a time series model，and a long-memory ARFIMA(p，d，q)model is introduced to simulate the time series model.We select randomly 10 H1N1 sequences and 10 H3N2 sequences in analysis.we find in these data a remarkably long-range correlation and fit the model reasonably by ARFIMA(p，d，q)model，and also find that we can use different ARFIMA models to identify the two kinds of sequences，i.e.ARFIMA(0，d，5)model and ARFIMA(1，d，1)model that can identify H1N1 and H3N2 respectively.

influenza A virus，time series model，chaos game representation(CGR)，ARFIMA(p，d，q)model

.E-mail:ezhun6669@sina.com

*江南大學(xué)創(chuàng)新團(tuán)隊(duì)發(fā)展計(jì)劃(批準(zhǔn)號(hào):2008CX002)中央高?；究蒲袠I(yè)務(wù)經(jīng)費(fèi)專項(xiàng)資金(批準(zhǔn)號(hào):JUSRP21117)資助的課題.

.E-mail:ezhun6669@sina.com

*Project supported by the Innovative Research Team of Jiangnan University(Grant No.2008CX002)the Foundamental Research Founds for the Central Universities(Grant No.JUSRP21117).

PACS:87.10.Vg，02.50.Fz