回溯搜索丁型肝炎病毒核酶的折疊路徑

鄒彥娟，張小龍，廖 昆

回溯搜索丁型肝炎病毒核酶的折疊路徑

鄒彥娟，張小龍，廖 昆

（萍鄉(xiāng)學(xué)院 機(jī)械電子工程學(xué)院，江西 萍鄉(xiāng) 337000）

RNA生物功能與其折疊的動(dòng)力學(xué)過(guò)程密切相關(guān)。文章采用枚舉法給出丁型肝炎病毒（Hepatitis Delta Virus, HDV）核酶序列的折疊空間，作為回溯搜索的解空間。結(jié)合自由能曲面計(jì)算了空間內(nèi)構(gòu)象之間的轉(zhuǎn)換速率。采用回溯搜索算法，從天然態(tài)開(kāi)始確定了HDV核酶序列的折疊路徑，并利用主方程方法分析了路徑上具有直接轉(zhuǎn)換關(guān)系的結(jié)構(gòu)之間占據(jù)幾率隨時(shí)間的變化。結(jié)果顯示，HDV天然態(tài)的形成過(guò)程呈兩相性，并且主要有3條折疊路徑。另外在折疊過(guò)程中存在一個(gè)滯留時(shí)間較長(zhǎng)的中間態(tài)假結(jié)結(jié)構(gòu)。

折疊路徑；回溯搜索；主方程方法；枚舉法；丁型肝炎病毒核酶

RNA生物功能的實(shí)現(xiàn)與其折疊過(guò)程中形成的結(jié)構(gòu)有關(guān)。RNA在折疊過(guò)程中可能會(huì)沿著不同的路徑形成復(fù)雜的二級(jí)或三級(jí)結(jié)構(gòu)。目前了解RNA結(jié)構(gòu)信息實(shí)驗(yàn)的方法有X射線晶體衍射方法術(shù)[1~3]、核磁共振法[4]和冷凍電鏡技術(shù)[5~6]等。這些實(shí)驗(yàn)方法測(cè)量結(jié)果準(zhǔn)確，可以觀察到RNA的晶體結(jié)構(gòu)及其在折疊過(guò)程中所經(jīng)歷的中間態(tài)，但耗時(shí)長(zhǎng)，成本高。隨著計(jì)算機(jī)的發(fā)展和應(yīng)用，通過(guò)理論計(jì)算模擬，如動(dòng)態(tài)規(guī)劃算法[7]、分子動(dòng)力學(xué)模擬算法[8]、蒙特卡洛算法[9]等，預(yù)測(cè)RNA的動(dòng)力學(xué)過(guò)程成為研究RNA的一個(gè)捷徑。但RNA的空間結(jié)構(gòu)數(shù)目會(huì)隨著鏈長(zhǎng)呈指數(shù)增長(zhǎng)[10]，利用計(jì)算機(jī)很難窮盡所有的構(gòu)象取樣。為此，筆者通過(guò)C語(yǔ)言編程，在個(gè)人計(jì)算機(jī)上運(yùn)行回溯搜索RNA酶的折疊路徑，從而可以分析折疊過(guò)程中的主要結(jié)構(gòu)及其對(duì)折疊過(guò)程造成的影響。

回溯搜索算法是程序設(shè)計(jì)中一種重要的基礎(chǔ)算法，是搜索算法中的一種控制策略?；厮菟阉魉惴ǖ幕驹硎窃谒阉鬟^(guò)程中從某一條搜索路徑一直往前走，直到不能前進(jìn)，然后退回來(lái)選擇另一條路徑再走。利用回溯搜索算法解決問(wèn)題，首先要定義問(wèn)題的解空間，然后利用便于搜索的方法組織解空間，并從問(wèn)題的解空間中搜索得到可行解?；赗NA結(jié)構(gòu)堿基配對(duì)關(guān)系，可以將回溯算法應(yīng)用于尋找RNA折疊過(guò)程任意兩態(tài)之間的折疊路徑。

丁型肝炎病毒（Hepatitis Delta Virus, HDV）核酶是乙型肝炎病毒（Hepatitis B Virus, HBV）的單鏈RNA病毒，在細(xì)胞內(nèi)通過(guò)雙滾環(huán)（double rolling-circle）機(jī)制完成復(fù)制，新產(chǎn)生的RNA再由基因或反基因酶剪切成單體術(shù)[11~12]。Perrotta和Been在實(shí)驗(yàn)中采用化學(xué)電泳法對(duì)樣品分析發(fā)現(xiàn)，正義和反義HDV序列在包含85個(gè)核苷酸時(shí)自剪切活性最佳術(shù)[13~14]。Bevilacqua實(shí)驗(yàn)室對(duì)HDV核酶的研究發(fā)現(xiàn)，HDV核酶的折疊路徑受到假結(jié)（Pseudoknot）的影響，其折疊動(dòng)力學(xué)過(guò)程呈兩相性[15]（折疊路徑分成了兩條）。理論研究方面，Chen等人發(fā)現(xiàn)，HDV核酶在折疊過(guò)程中由于假結(jié)的影響存在快慢兩條路徑[16]。

本文根據(jù)RNA堿基配對(duì)規(guī)則，采用枚舉法定義HDV核酶的構(gòu)象空間，作為回溯搜索的解空間；根據(jù)RNA折疊特點(diǎn)結(jié)合自由能曲面計(jì)算了構(gòu)象之間的轉(zhuǎn)換速率，作為回溯搜索條件；利用主方程方法術(shù)[17~18]計(jì)算各個(gè)構(gòu)象的占據(jù)幾率分布，分析HDV序列的折疊特點(diǎn)；最后回溯搜索計(jì)算了HDV序列由展開(kāi)鏈到天然態(tài)結(jié)構(gòu)的折疊路徑，找出了折疊過(guò)程中出現(xiàn)的重要的中間態(tài)結(jié)構(gòu)，并分析了其在折疊過(guò)程中發(fā)揮的作用。研究結(jié)果顯示，HDV序列折疊過(guò)程存在3條主要的路徑：一條快速路徑，兩條慢速路徑。

1 研究方法

1.1 枚舉生成構(gòu)象回溯空間

RNA分子在細(xì)胞內(nèi)通過(guò)形成連續(xù)的堿基配對(duì)可以折疊形成穩(wěn)定的螺旋結(jié)構(gòu)。筆者以螺旋為基本單元構(gòu)建RNA結(jié)構(gòu)，并限制每個(gè)螺旋至少由三個(gè)連續(xù)的堿基配對(duì)構(gòu)成。在程序中，首先，對(duì)RNA序列進(jìn)行順序搜索，生成所有可能的螺旋組成螺旋數(shù)據(jù)庫(kù)。然后，確定不同螺旋之間的相互關(guān)系（圖1）。螺旋之間存在三種關(guān)系：（a）完全兼容：兩個(gè)螺旋可以共存，無(wú)相互交疊的堿基；（b）完全不兼容：兩個(gè)螺旋不能共存，堿基完全交疊；（c）部分兼容：兩個(gè)螺旋可以部分共存，堿基部分交疊。最后，選取完全兼容或部分兼容的螺旋組裝RNA結(jié)構(gòu)，生成構(gòu)象空間。

圖1 螺旋之間的相互關(guān)系

1.2 計(jì)算構(gòu)象之間的轉(zhuǎn)換速率

在以螺旋為基本單元的模型中，RNA結(jié)構(gòu)之間的躍遷，通過(guò)螺旋的形成、斷開(kāi)和交換來(lái)進(jìn)行。

（1）形成螺旋的速率。RNA形成一個(gè)新的螺旋可從任意一個(gè)“堿基堆積”開(kāi)始，不同的“堿基堆積”對(duì)應(yīng)多條不同的生長(zhǎng)路徑，如圖2（a）所示。對(duì)于任意一個(gè)“堿基堆積”其相鄰的堿基形成速率要遠(yuǎn)遠(yuǎn)快于其他堿基的形成，且其自由能曲面呈逐漸下降的趨勢(shì)，如圖2（b）所示。對(duì)于含有三個(gè)連續(xù)堿基堆積的螺旋，其穩(wěn)定性足以抵抗環(huán)境的熱擾動(dòng)，所以計(jì)算某條路徑的折疊速率可用形成前三個(gè)堿基堆積的速率進(jìn)行估算術(shù)[19~20]，計(jì)算公式是：

（3）螺旋交換速率。如果A與B兩個(gè)螺旋相互交疊，則通過(guò)螺旋A斷開(kāi)、螺旋B再形成的自由能勢(shì)壘的高度，要高于螺旋A斷開(kāi)的同時(shí)螺旋B形成造成的自由能勢(shì)壘的高度?？筛鶕?jù)隧道路徑理論計(jì)算兩個(gè)螺旋之間的交換速率（如圖2（c）所示）：

式中，和分別表示堿基在螺旋A中形成（在螺旋B中打開(kāi)）或堿基在螺旋B中形成（在螺旋A中打開(kāi)）的速率，為玻爾茲曼常數(shù)，為螺旋A與螺旋B所在RNA結(jié)構(gòu)的自由能差值。

1.3 計(jì)算構(gòu)象的占據(jù)幾率

1.4 回溯搜索折疊路徑

2 計(jì)算結(jié)果

2.1 構(gòu)象占據(jù)幾率的變化

本文用于計(jì)算的HDV序列包含95個(gè)核苷酸（CCUGAUGGCCGGCAUGGUCCCAGCCUCCUCGCUGGCGCCGGCUGGGCAACAUUCCGAGGGGACCGUCCCCUCGGUAAUGGCGAAUGGGACGCACA），生成的構(gòu)象空間包含2971個(gè)RNA結(jié)構(gòu)。通過(guò)計(jì)算構(gòu)象空間中各個(gè)態(tài)占據(jù)幾率隨時(shí)間的變化，發(fā)現(xiàn)在折疊過(guò)程中起重要作用的中間態(tài)結(jié)構(gòu)“1627”，在反應(yīng)開(kāi)始后的10-3.5s內(nèi)，大約有50%的HDV序列快速的形成結(jié)構(gòu)“1627”，在10-3.5s到10 s內(nèi)HDV序列形成結(jié)構(gòu)“1627”的速率變得緩慢，10 s到1000 s的時(shí)間內(nèi)，結(jié)構(gòu)“1627”緩慢地轉(zhuǎn)換成其他結(jié)構(gòu)，如圖3（a）所示。我們同時(shí)發(fā)現(xiàn)天然態(tài)結(jié)構(gòu)的形成明顯的分為兩個(gè)階段，在折疊開(kāi)始的18s內(nèi)大約有40%的HDV序列迅速的形成天然態(tài)結(jié)構(gòu)，剩余的大約60%的RNA序列在隨后的35分鐘內(nèi)逐漸形成天然態(tài)結(jié)構(gòu)，如圖3（b）所示。造成天然態(tài)結(jié)構(gòu)的形成分成兩個(gè)階段的原因可能是“1627”在10 s前快速形成，在10s后慢速形成，為了解釋這一原因本文對(duì)HDV的折疊路徑進(jìn)行搜索。

（a）HDV序列折疊動(dòng)力學(xué)過(guò)程，“1”是展開(kāi)鏈，“1627”是一中間態(tài)，“Native”是天然態(tài)結(jié)構(gòu)；（b）天然態(tài)結(jié)構(gòu)占據(jù)幾率隨時(shí)間變化分成兩個(gè)階段

圖3 HDV序列的中間態(tài)和天然態(tài)的占據(jù)幾率隨時(shí)間的變化

2.2 回溯搜索HDV序列的折疊路徑

通過(guò)計(jì)算HDV序列構(gòu)象空間內(nèi)各RNA結(jié)構(gòu)占據(jù)幾率隨時(shí)間的變化，可知HDV序列在折疊過(guò)程中存在滯留時(shí)間較長(zhǎng)的中間態(tài)“1627”。為了解釋“1627”在折疊過(guò)程中發(fā)揮的作用，本文采用回溯搜索算法，計(jì)算了HDV序列的折疊路徑，同時(shí)計(jì)算了在折疊路徑上有直接轉(zhuǎn)換關(guān)系的主要結(jié)構(gòu)之間的轉(zhuǎn)換流量，如圖4（a）、（b）、（c）所示。

HDV天然態(tài)結(jié)構(gòu)由構(gòu)象“1175”與構(gòu)象“1636”轉(zhuǎn)換形成，分別約占87%，5%，見(jiàn)圖4（a）。構(gòu)象“1175”與天然態(tài)結(jié)構(gòu)之間的轉(zhuǎn)換呈現(xiàn)出了二相性，見(jiàn)圖4（a），約40%的HDV序列在折疊開(kāi)始的18s內(nèi)快速轉(zhuǎn)換成天然態(tài)結(jié)構(gòu)，而約50%的HDV序列在35 min內(nèi)緩慢地轉(zhuǎn)換成為天然態(tài)結(jié)構(gòu)，這與天然態(tài)結(jié)構(gòu)的形成分成兩個(gè)階段相一致。

HDV序列從展開(kāi)鏈到天然態(tài)的主要折疊路徑，沿箭頭方向的數(shù)字為轉(zhuǎn)換速率（單位：s-1），括號(hào)中的數(shù)字為構(gòu)象自由能（單位：kcal/mol）

結(jié)論

本文利用回溯搜索算法，通過(guò)搜索HDV核酶由展開(kāi)鏈形成天然態(tài)結(jié)構(gòu)的主要路徑，很好地解釋了HDV折疊過(guò)程中分成兩個(gè)階段的原因。同時(shí)計(jì)算結(jié)果顯示，在HDV序列折疊過(guò)程存在三條主要的折疊路徑。從中間態(tài)結(jié)構(gòu)“1182”開(kāi)始，轉(zhuǎn)換過(guò)程分成了快慢兩個(gè)階段，對(duì)應(yīng)構(gòu)象“1253”的路徑為快速路徑，對(duì)應(yīng)構(gòu)象“1459”的路徑為慢速路徑。導(dǎo)致構(gòu)象“1459”轉(zhuǎn)換形成構(gòu)象“1182”較慢的原因是，在折疊過(guò)程中大約有50%的HDV序列在構(gòu)象“1627”滯留了約100s的時(shí)間，之后這些構(gòu)象才在約35 min的時(shí)間內(nèi)緩慢轉(zhuǎn)換形成構(gòu)象“1459”。本研究進(jìn)一步的工作，可將回溯算法與其他算法結(jié)合，不斷改進(jìn)回溯算法，更好的預(yù)測(cè)RNA序列的折疊過(guò)程。

[1] TERESHKO V, WALLACE ST, USMAN N, et al. X-ray crystallographic observation of “ in-line ” and “ adjacent ” conformations in a bulged self-cleaving RNA / DNA hybrid[J]. RNA, 2001, 7(3): 405~420.

[2] ROBERT R. CRICHTON, RICARDO O. LOURO. Practical approaches to biological inorganic chemistry (Second Edition)[M]. Holland: Elsevier, 2020: 375~416.

[3] JUN S, HIRATA A, KANAI T, et al. The X-ray crystal structure of the euryarchaeal RNA polymerase in an open-clamp configuration[J]. Nat Commun, 2014,5(5132): 1~11.

[4] FüRTIG B, RICHTER C, W?HNERT J, et al. NMR spectroscopy of RNA[J]. ChemBioChem, 2003, 4(10): 936~962.

[5] BARRAUD P, GATO A, HEISS M, et al.Time-resolved NMR monitoring of tRNA maturation[J]. Nature Communications, 2019, 10(1): 3373.

[6] GARMANN RF, GOPAL A, ATHAVALE SS, et al. Visualizing the global secondary structure of a viral RNA genome with cryo-electron microscopy[J]. RNA, 2015, 21(5): 877~886.

[7] AKUTSU T. Dynamic programming algorithms for RNA secondary structure prediction with pseudoknots[J]. Discret Appl Math, 2000, 104(1~3): 45~62.

[8] HASHEM Y, AUFFINGER P. A short guide for molecular dynamics simulations of RNA systems[J]. Methods, 2009, 47(3): 187~197.

[9] CLOTE P, BAYEGAN AH. Mathematical Biology RNA folding kinetics using Monte Carlo and Gillespie algorithms[J]. J Math Biol, 2017, 76(5): 1195~1227.

[10] 鄒彥娟, 石方. 在直角坐標(biāo)系中回溯枚舉RNA結(jié)構(gòu)[J]. 萍鄉(xiāng)學(xué)院學(xué)報(bào), 2018, 35（6）: 57~60.

[11] MACNAUGHTON TB, SHI ST, MODAHL LE, et al. Rolling circle replication of Hepatitis Delta Virus RNA is carried out by two different cellular RNA polymerases[J]. Virology, 2002, 76(8): 3920~3927.

[12] BROWN TS, CHADALAVADA DM, BEVILACQUA PC. Design of a highly reactive HDV ribozyme sequence uncovers facilitation of RNA folding by alternative pairings and physiological ionic strength[J]. J Mol Biol, 2004, 341(3): 695~712.

[13] PERROTTA AT, BEEN MD. A pseudoknot-like structure required for efficient self-cleavage of hepatitis delta virus RNA[J]. Nature, 1991, 350(6317): 434~436.

[14] PERROTTA AT, MICHAEL DB. The self-cleaving domain from the genomic RNA of hepatitis delta virus: sequence requirements and the effects of denaturant[J]. Nucl Acids Res, 1990, 18(23): 6821~6827.

[15] CHADALAVADA DM, SENCHAK SE, BEVILACQUA PC. The folding pathway of the genomic hepatitis delta virus ribozyme is dominated by slow folding of the pseudoknots[J]. J Mol Biol, 2002, 317(4): 559-575.

[16] CHEN JW, GONG S, WANG YJ, et al. Kinetic partitioning mechanism of HDV ribozyme folding[J]. J Chem Phys, 2014, 140(2): 025102.

[17] ZHAO PN, ZHANG WB, CHEN SJ. Cotranscriptional folding kinetics of ribonucleic acid secondary structures[J]. J Chem Phys, 2011, 135(24): 245101.

[18] ZHANG WB, CHEN SJ. Master equation approach to finding the rate-limiting steps in biopolymer folding[J]. J Chem Phys, 2003, 118(7): 3413~3420.

[19] GONG S, WANG YJ, ZHANG WB. The regulation mechanism of yitJ and metF riboswitches[J]. J Chem Phys, 2015, 143(4): 045103.

[20] GONG S, WANG YJ, ZHANG WB. Kinetic regulation mechanism of pbuE riboswitch[J]. J Chem Phys, 2015, 142(1): 015103.

[21] ZHAO PN, ZHANG WB, CHEN SJ. Predicting secondary structural folding kinetics for nucleic acids[J]. Biophys J, 2010, 98(8): 1617~1625.

Backtracking Search for the Folding Pathway of the Hepatitis Delta Virus Ribozyme

ZOU Yan-juan, ZHANG Xiao-long, LIAO Kun

(School of Mechanical and Electronic Engineering, Pingxiang University, Pingxiang Jiangxi 337000, China)

RNA biological functions are closely related to its folding kinetics. In the paper, the folding space of the Hepatitis Delta Virus (HDV) ribozyme sequence is given by enumeration method as a solution space for backtracking search. The conversion rates between conformations within the space were calculated in combination with free energy surfaces. The backtracking search algorithm was used to determine the folding path of the HDV ribozyme sequence from the native state, and the change in the occupancy probability over time between structures with direct transition relationships in the path was analyzed by the master equation method. The results show that the natural state formation process of HDV is biphasic and that there are three main folding paths. In addition, one pseudoknot trap state is formed with a long residence time in the folding process.

folding path; backtracking search; master equation method; enumeration; HDV ribozyme

2020-03-07

萍鄉(xiāng)學(xué)院青年科研基金（2018D0227、2018D0228）

鄒彥娟（1988—），女，河南濮陽(yáng)人，講師，碩士，研究方向：凝聚態(tài)物理學(xué)的研究。

O469；Q61

A

2095-9249（2020）03-0079-06

〔責(zé)任編校：吳侃民〕