劉 娟 高 潔
(江南大學理學院,無錫 214122)(2010年4月16日收到;2010年8月4日收到修改稿)
甲型流感病毒DNA序列的長記憶ARFIMA模型*
劉 娟 高 潔
(江南大學理學院,無錫 214122)(2010年4月16日收到;2010年8月4日收到修改稿)
流感病毒分為三類:甲型(A型),乙型(B型),丙型(C型).在這三種類型中甲型(A型)流感病毒是最致命的流感病毒,對人類引起了嚴重疾病.本文對甲型流感病毒DNA序列建立了一種新的時間序列模型,即CGR(Chaos Game Representation)弧度序列.利用CGR坐標將甲流病毒DNA序列轉換成CGR弧度序列,且引入長記憶ARFIMA模型去擬合此類序列,發現隨機找來的10條H1N1序列,10條H3N2序列都具有長相關性且擬合很好,并且還發現這兩種序列可以嘗試用不同的 ARFIMA模型去識別,其中 H1N1可用 ARFIMA(0,d,5)模型去識別,H3N2可用ARFIMA(1,d,1)模型去識別.
甲型流感,時間序列模型,CGR,ARFIMA(p,d,q)模型
PACS:87.10.Vg,02.50.Fz
流感是一種反復出現的傳染病,在全球引起了高發病率和高死亡率[]1.流感病毒分為三類:甲型(A型),乙型(B型),丙型(C型).在這三種類型中甲型(A型)流感(以下簡稱甲流)病毒是最致命的流感病毒,給人類帶來了嚴重的疾病.甲流病毒根據其表面的血凝素(hemagglutinin,HA)和神經氨酸酶(neuraminidase,NA)基因的不同又可分成16個HA亞型(H1-H16)和9個NA亞型(N1—N9),不同的HA和NA形成了甲流病毒的許多亞型,如H1N1,H3N2,H5N1 等等[2—4].筆者參看了許多文獻幾乎沒有看到用時間序列模型來挖掘甲型流感病毒的特性的,因而本文采用時間序列模型來分析甲型流感病毒.
1992年,Peng等[5]提出了 DNA一維游走模型.同年Voss等[6]提出了不同的觀點,他們發現 DNA序列的譜密度顯示的1/fβ噪聲無處不在,意味著當0<β<1存在長相關性,認為不僅在非編碼區序列中在編碼序列中也存在長相關性.另一方面Buldyrev 等[7,8]設計了一個廣義-Lévy 游走模型去生成一個模型序列,使用所有可用的DNA序列發現主要在非編碼序列中呈現長相關性.基于該模型,Tai等[9]提出了一個二維修正-Lévy游走模型.為區分 C和T,A和 G,Lou等建立了二維和三維游走模型[10],Yu 等則建立了圖譜[11,12],來研究 DNA 序列的相關性.2006 年,Lopes和 Nunes[13]引入長記憶ARFIMA(0,d,0)模型去擬合 DNA序列的一維游走序列.2009年,Gao等[14]基于 CGR(chaos game representation)坐標提出了一種DNA序列轉換成一個時間序列(CGR-游走序列)的方法,并引入長記憶ARFIMA(p,d,q)模型來分析.
本文對甲流病毒DNA序列提供了一種新的時間序列模型,即CGR弧度序列.利用CGR坐標將甲流病毒DNA序列轉換成CGR弧度序列,且引入長記憶ARFIMA模型去擬合此類序列,發現隨機找來的10條H1N1序列,10條H3N2序列都具有長相關性且擬合很好,并且還發現這兩種序列可以嘗試用不同的 ARFIMA模型去識別,其中 H1N1可用ARFIMA(0,d,5)模型去識別,H3N2可用 ARFIMA(1,d,1)模型去識別.
如果隨機過程 {xt}是平穩的,且滿足方程Φ(B)Δdxt= Θ(B)εt,其中,- 0.5 則稱{xt}服從 -0.5 因此,{xt}可看作是分數差分噪聲導出的 ARMA(p,q)過程.當2d-1=-1時,d=0即為短記憶過程;所以當2d-1>-1時,d∈(0,0.5)具有長記憶的特征. 1990年Jeffrey提出了一種 DNA序列可視化的方法即 CGR方法[15].CGR是一種迭代映射技術,它把序列中的每個單元,如蛋白質序列中氨基酸,DNA中的核苷酸,映射到一個連續的坐標空間中去. 正方形的四個頂點對應四種核苷酸.在這里,用DNA序列代替隨機數,每一個堿基的坐標都可以來確定下一個堿基的位置.我們取 A(0,0),T(1,0),G(1,1),C(0,1),并且取點(0.5,0.5)為起始點. 下面給出DNA迭代函數,也可以認為是 CGR算法 的 公 式 化 形 式[15,16]. 對 于 一 個 序 列 S = 其中 gi={(0,0),(1,0),(1,1),(0,1)},gi和 si相對應. 對于一個DNA序列,定義 其中yn是CGRn的y坐標值,xn是 CGRn的 x坐標值.則得到一個數據序列 {Rn:n=1,2,…,N},我們把它作為一個時間序列,并稱它為“CGR弧度序列”. 以甲流病毒 H1N1序列 CY056890為例,數據來自 NCBI網站,其網址:http://www.ncbi.nlm.nih.gov/. 它的CGR弧度序列“游走圖”如下表1. 表1 CY056890序列所選部分前8個ACATGGTA游走結果 圖1(a)是CY056890序列CGR弧度序列圖(位置380—2170),樣本容量1791.這些數據變動較大,呈現非平穩特征.考慮對此過程作d階差分.先對原序列作對數變換然后再做一階差分結果如圖1(b)所示,可見除少數地方呈現異方差外,基本呈現平穩性. 圖2(a)(ACF)和圖2(b)(PACF)為樣本取對數再一階差分后的自相關函數圖形和偏自相關函數圖形.可見ACF衰減迅速,而PACF衰減緩慢,這意味著原序列具有長記憶特征. 圖1 (a)甲流H1N1型病毒CY056890的弧度序列圖;(b)取對數再一階差分圖 圖2 (a)取對數再一階差分的樣本自相關圖;(b)取對數再一階差分的樣本偏自相關圖 圖3給出了方差圖[17]是一個估計長記憶參數d的有用工具.對于一個長記憶時間序列{Rn},它的均值珔Rk的方差滿足作log[Var(珔Rk)]關于log(k)的散點圖,對散點圖線性擬合,可估計得到線性方程的斜率為 -0.6877,令該斜率為2d-1= -0.6877,即可得 d的估計值0.156. 根據上述理由我們可選擇CGR弧度序列顯示長記憶特征.目的是利用上述特點為序列建立一個合適的模型.因此,可以考慮長記憶 ARFIMA(p,d,q)模型(d∈ (0,0.5)),p,q定階時為考慮實用性,僅考慮 p,q均小于等于5的 ARFIMA(p,d,q)模型.由 Akaike 信 息 判 別 準 則[18,19],可 選 ARFIMA(0,0.156,5)模型來擬合. 為檢驗該模型的合理性,選擇了一個合適的檢驗統計量 LB 檢驗統計量[20,21] 其中rk是滯后k的樣本自相關函數,n是樣本容量,M是一個取定的比n小的正整數. 表2顯示了對于各滯后階數,LB統計量的p值均顯著大于0.1,意味著擬合模型的殘差序列應為白噪聲(純隨機),因而可以認為 ARFIMA(0,0.156,5)模型能很合理地擬合 CY056890序列的CGR-游走序列. 圖3 CY056890DNA序列的CGR弧度序列方差圖 表2 殘差的自相關檢驗 表3給出了被選擇的 ARFIMA(0,0.156,5)模型的參數估計,5個參數的 T檢驗統計量的p值均顯著小于 0.005.這意味著 ARFIMA(0,0.156,5)模型能有效地擬合這個CGR弧度序列. 表3 條件最小二乘估計 表4和表5分別給出了隨機選的9條H1N1序列和10條 H3N2序列的數據信息、被選擇的ARFIMA(p,d,q)模型及參數估計.從計算結果可得d均位于(0,0.5);對于各滯后階數,LB統計量的 p值除極個別外其余均顯著大于0.1;且每個被選擇的模型中各參數的T檢驗統計量的p值均顯著小于0.01.所有這些結果都顯示ARFIMA(p,d,q)模型能很合理很有效地擬合這些不同的CGR弧度序列且還發現所選H1N1序列均為ARFIMA(0,d,5)模型,所選 H3N2序列均為 ARFIMA(1,d,1)模型.所以我們可嘗試用 ARFIMA(0,d,5)模型,ARFIMA(1,d,1)模型分別去識別H1N1序列,H3N2序列. 表4 9條H1N1序列的數據信息、被選擇的ARFIMA模型和參數估計 表5 10條H3N2序列的數據信息、被選擇的ARFIMA模型和參數估計 本文基于CGR坐標提出了一種將甲流病毒DNA序列轉換成時間序列(CGR弧度序列)的方法,并引入長記憶模型ARFIMA模型來分析,首先分析了甲流H1N1型病毒CY056890序列,從圖1到圖3可知弧度序列顯示長記憶特征,并選擇了ARFIMA(0,0.156,5)模型去擬合它,從表 2到表 3發現擬合合理有效. 然后又分析了隨機找來的19條序列的CGR弧度序列,從表4和表5可知所有 ARFIMA(p,d,q)模型都有效合理.并且從表4中還發現所選H1N1序列均為 ARFIMA(0,d,5)模型,表5所選 H3N2序列均為 ARFIMA(1,d,1)模型. 由此可見,DNA序列的CGR弧度序列能由長記憶ARFIMA(p,d,q)模型有效合理地擬合,并且還可嘗試用 ARFIMA(0,d,5)模型,ARFIMA(1,d,1)模型分別去識別H1N1序列、H3N2序列.作為具有完善算法的經典時間序列模型,不僅可以幫助我們得到甲流病毒DNA序列清晰的結構,而且還可幫助我們有效識別甲流中的兩種亞型. 本文僅對甲流中的兩種亞型進行了研究分析,后面我們將研究分析甲流中的其他亞型以及乙型丙型流感病毒. [1] Morens D,Folkers G,Fauci A 2004 Nature 430 242 [2] Chen J M,Sun Y X,Liu S 2009 Chinese Science Bulletin 54 1657(in Chinese)[陳繼明、孫映雪、劉 朔 2009科學通報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 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 *江南大學創新團隊發展計劃(批準號:2008CX002)中央高校基本科研業務經費專項資金(批準號: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
3.基于CGR的時間序列模型



4.甲流 H1N1型病毒 CY056890的數據分析







5.其余9條H1N1序列和10條H3N2序列數據分析


6.結 論
(School of Science,Jiangnan University,Wuxi 214122,China)(Received 16 April 2010;revised manuscript received 4 August 2010)