LncRNA的靶RNA 預測研究

鄭希強，宋子健，李建偉

（河北工業大學人工智能與數據科學學院計算醫學研究所，天津 300401）

隨著生命科學研究的深入和高通量測序技術的飛速發展，研究發現大多數人類基因組均已轉錄[1]。據報道[2]，其中僅有2%的RNA 分子可翻譯為蛋白質，絕大部分RNA 是不能編碼蛋白的非編碼RNA。其中，長度超過200 個核苷酸的非編碼RNA被定義為長鏈非編碼RNAs（long non-coding RNAs，lncRNAs）。近年來，越來越多的研究表明，lncRNAs 參與調控多種重要的生物過程，如調控表觀遺傳染色質修飾、轉錄調控和轉錄后調控等[3]。2012 年，Gong ZJ 等[4]發現HOTAIR 可通過與PRC2復合物相互作用，誘導H3K27 發生甲基化，進而抑制該基因的表達。此外，lncRNAs 與癌癥[5]、心血管疾病[6]和糖尿病[7]等多種復雜疾病的發生發展均存在著密切聯系。因此，深入了解lncRNA的調控功能，對人類復雜疾病的診斷、治療和預后提供新的治療靶點具有重要意義。

1 LncRNA的靶RNA 預測研究

迄今為止，海量的lncRNAs 被不斷發現，但對其調控功能的研究仍不夠充分。通過傳統的生物實驗驗證lncRNA調控功能，結果雖然準確、可靠，但存在實驗周期長、費用高、效率低等問題。單純地依靠傳統的生物實驗已遠遠不能滿足當前需求，借助高效的、基于生物信息學的計算方法預測lncRNAs的調控功能已成為計算醫學研究和對復雜疾病機制研究的熱門方向。研究發現[8]，lncRNAs 常通過堿基互補配對的方式調控靶RNA的表達，進而發揮其調控作用。因此，根據lncRNA 與靶RNA 之間的相互作用關系，通過計算預測lncRNA的靶RNA，利用靶RNA的生物學功能推斷lncRNAs的調控功能已成為一種常見的研究方法。此外，lncRNAs 會折疊成復雜的二級結構來維持它們的穩定性和保守性。lncRNA的二級結構可以幫助預測lncRNAs的靶RNA[9]。lncRNA 二級結構的生物學意義及其堿基序列組成同樣重要，部分lncRNAs的功能活性通常由它們特定的二級結構決定。

基于以上分析，開發lncRNA 靶RNA 預測模型在全基因組范圍預測lncRNA-RNA 相互作用，快速搜索lncRNA 與RNA 之間最穩定的堿基互補配對狀態，成為lncRNA-RNA 相互作用預測的關鍵[10]。隨著生物信息學的不斷發展，涌現出了數量眾多的RNA-RNA 相互作用預測模型，這些模型均可用于lncRNA 靶RNA 預測。常見的RNA-RNA 相互作用預測模型見表1。

表1 常見預測lncRNA 靶RNA 模型的關鍵特性

2 常見的RNA-RNA 相互作用預測模型

2.1 LncTar LncTar 預測模型[11]通過標準化自由能（ndG）輔助科研人員預測lncRNAs的靶RNA。Li J等[11]從RNA 堿基互補配對的角度分析了lncRNAs與靶RNA 之間的相互作用關系，并將多重PCR 技術中引物二聚體的思想運用到篩選lncRNA 靶RNA中。LncTar 預測模型利用“滑動”算法和最近鄰能量模型[15，16]尋找并計算結合區域的自由能，最終通過計算ndG 與自定義閾值實現lncRNA-靶RNA 相互作用關系預測。

LncTar 預測模型包含三個步驟。首先，創建二維矩陣，用于記錄兩個RNAs 分子的互補堿基對。在計算兩個RNAs 形成的結合區域時，采用最近鄰能量模型。其中自由能（ΔG）由給定的焓值（ΔH）、熵值（ΔS）和特定溫度下的熔化溫度T 來計算。具體計算公式：ΔG=ΔH-TΔS。其次，遍歷兩個輸入RNAs 分子形成的全部結合區域，分別計算它們的總自由能（dG），并記錄在哈希表中。遍歷完成后，在哈希表中尋找最小的總自由能（dG）并獲得其結合區域。最后，LncTar 提出一個判定標準——標準化自由能（ndG），用它來反映自由能最小的結合區域中lncRNA-靶RNA 結合的相對穩定性。具體計算公式如下：

其中結合區域長度（bindingregion）是兩個RNAs序列中較短的RNA 序列長度。如果計算出兩個RNAs 形成結合區域的ndG 小于等于人為自定義閾值，就可以判定兩個RNAs 分子有相互作用關系，否則給出相反判定。

LncTar 預測模型的優勢在于對RNA 序列長度沒有任何限制，可適用任意長度的RNA 分子。并且LncTar 提供了標準化自由能，能自動判定兩個RNAs 分子是否有相互作用關系，極大地減輕用戶的判定負擔。此外，LncTar的算法時間復雜度僅為O（n2）。但LncTar 預測模型也有一定的局限性。如因為追求計算速度，LncTar 沒考慮到RNA 二級結構，而lncRNA的二級結構對大部分lncRNA-靶RNA相互作用的判定有著一定的影響。此外，如果兩個RNAs 轉錄本序列長度較長的情況下，預測出的結合區域長度過長，將導致ndG 值過小，具有靶向關系的兩個RNAs 會被誤判。而ndG的閾值如果設置過小，雖然降低了假陽性率，但會導致一部分正樣本預測結果不準。因此，LncTar 預測模型還有很大的提升空間。

2.2 RIsearch2 RIsearch2 模型[12]適用于大規模預測RNAs 間的相互作用。Ferhat A 等[12]在允許G-U 搖擺配對的情況下搜索與種子反向互補匹配的序列片段，并使用后綴數組索引與動態規劃擴展種子的方法快速定位查詢序列和目標序列之間的互補結合區域。RIsearch2 模型分為兩個階段。第一階段：將查詢序列和目標序列各自轉換為后綴數組索引結構，并存儲目標序列后綴數組的反向互補序列。轉換為后綴數組索引結構的算法是基于GUUGle 算法[17]。大部分算法將Watson-Crick 對（A-T，C-G）視為匹配，并且允許存在G-U 搖擺配對。而Ferhat A 等提出一種無需反轉查詢序列來找尋結合區域的方法，將GA 和U-C 匹配視為有效對，然后使用二分查找來確定目標序列后綴數組中與查詢序列匹配的區域。這種方法可以直接匹配目標序列的后綴數組，不需要反轉查詢序列。第二階段：為了減少種子相互重疊造成的冗余，首先將冗余的匹配種子進行過濾，然后通過填充種子上下游區域的動態規劃矩陣擴展剩余的種子區域，接著利用RIsearch[18]中簡化的最近鄰能量模型計算種子匹配區域的自由能，并與用戶定義的分子間自由能閾值進行比較，篩選出滿足條件的匹配區域。最后，通過回溯動態規劃矩陣找到實際的相互作用區域。

RIsearch2 預測模型將RIsearch 模型與GUUGle模型進行深度集成，并且采用多個線程并行處理（OpenMP API）的方式對兩個后綴數組進行匹配。這種高效的計算方式不僅克服了傳統算法一次只能處理一個目標序列的障礙，還提高了計算速度。但RNA-RNA 相互作用預測中忽略RNA 二級結構信息會降低預測的準確性。

2.3 RIblast RIblast 預測模型[13]一種基于種子擴展方法[19]的超快速RNA-RNA 相互作用預測模型。首先，RIblast 使用后綴數組尋找種子匹配區域；其次，通過互補能量模型擴展種子區域的兩端；最后，RIblast 模型使用兩種能量之和作為目標序列和查詢序列相互作用結合區域的輸出結果。其中一種為可訪問能量，這種能量是防止結合區域內形成分子內堿基對所需的能量，通過使用分區函數算法進行計算[20]。另一種能量為雜交能量，它源自兩個RNAs 片段之間分子間堿基配對產生的自由能，可以基于最近鄰能量模型計算。

RIblast 預測模型分為兩個步驟：數據庫構建和RNAs 結合區域搜索。在數據庫構建中，首先是計算目標RNA 數據集中每個片段的可訪問能量，然后將目標RNA 序列反轉，并在未反轉與反轉的序列之間插入分隔符，接著構造串聯序列的后綴數組索引，為了加快RNA 相互作用區域的搜索，可對短字符串預先計算其雜交能量，最后將可訪問能量串聯序列的后綴數組索引和短字符串計算互補能量的結果存儲在數據庫中。在RNAs 結合區域搜索中，首先是計算查詢序列的可訪問能量，并為查詢序列構建后綴數組索引，然后基于查詢和數據庫的兩個后綴數組，查找雜交能量小于閾值能量的種子區域。接著，通過雜交能量和兩個可訪問能量的總和來計算找到的種子區域的相互作用能。繼次，RIblast 在種子區域兩邊進行擴展，沒有間隙。當擴展中不再形成RNA 互補雙鏈體時，RIblast 終止擴展。最后，刪除與其他種子區域完全重疊的種子。另外，也排除超過閾值能量的相互作用結果。

對于RNA-RNA 相互作用預測，可訪問能量的計算時間與序列長度成線性關系，雜交能量的計算時間是種子匹配項的平方，這是提高運行速度的障礙。而RIblast 先進行了大量的預處理工作，然后在查詢序列和目標序列之間找到短匹配區域以加快速度（稱為種子），最后將檢測到的種子兩端擴展，這種方法解決了時間復雜度高的問題。RIblast的不足之處在于最后只給出能量信息，部分用戶無法通過這些能量信息判斷兩個RNAs 分子是否存在相互作用關系。沒有給予用戶一些判定標準，增加了人工判定的負擔。

2.4 ASSA Antonov I 等[14]發現總自由能的值取決于RNA的轉錄本長度和GC 含量兩個因素，這使研究人員很難比較不同特性轉錄本所計算出的RNARNA 相互作用的總自由能。為了解決這個問題，ASSA 預測模型[14]根據兩個相互作用RNAs的轉錄本長度和GC 含量設計出了總自由能的統計顯著性（P值）。ASSA 預測模型包括三個步驟。第一步，使用LASTAL[21]程序包中的局部序列比對工具，搜索查詢序列和目標序列中的局部互補區域。在這項工作中，ASSA 將這些局部互補區域稱為“反義位點”。第二步，首先提取LASTAL 預測的反義位點以及側翼序列組成“假定雙鏈體”，然后根據熱力學模型計算其雜交能量。“假定雙鏈體”中可能存在局部的RNA 二級結構，研究人員通過觀察總結出兩種“假定雙鏈體”類型。第一種類型稱之為“真雙鏈體”，這種類型的RNA 雙鏈體與任一側附近的（側翼）區域不具有強互補性，即單鏈RNA 無法形成二級結構。因此，這種類型的RNA 雙鏈體只形成分子間堿基對。第二種類型稱之為“兩個發夾”。除了RNAs 分子間堿基互補配對以外，兩個RNAs 自身還具有完美的（100%）分子內互補性。為了將“真雙鏈體”與“兩個發夾”兩種類型區分開，ASSA 預測模型使用基于熱力學的工具RNAup[22]，因為它能夠通過計算分子間自由能相對準確地分開兩種類型。ASSA 將RNAup 應用于特定序列組塊（“假定雙鏈體”）而不是全長轉錄本，這樣可以有效并快速地計算RNA 分子間的雜交能量。最后，全部特定序列組塊相加當作兩個RNAs 結合產生的總自由能。第三步，估算獲得的總自由能的統計顯著性。通過將觀察的RNA 序列[總自由能的分布取決于RNA 序列的特征（長度和GC 含量）]與相應隨機序列（ASSA 生成的負樣本）產生的總自由能進行比較，估算總自由能的統計顯著性P值。具體計算公式如下：

其中x 是觀察的RNA 序列通過ASSA 預測的總自由能，Num（ΔGrand≤x|T）是隨機序列產生的總自由能小于或等于x的隨機序列個數，Num（ΔGrand）是隨機序列的總數。統計顯著性P值越小，說明觀察的兩個RNAs 有相互作用的可能性越大。

在ASSA 中的主要創新是數學統計模型，它可以快速計算并估計觀察的RNA-RNA 相互作用產生總自由能的統計顯著性。在RNA-RNA 相互作用預測模型中，ASSA 是唯一估計P值的預測模型。由于僅使用局部RNA 二級結構，對于全序列來說，反義位點上RNA 二級結構可能超出它側翼區域的范圍，忽略該二級結構的自由能會對最后結果的準確性產生一定的影響。

3 總結

不同的lncRNA 靶RNA 預測模型各有其自身的特點和優勢。研究者應根據研究目的和需求，進行合理選擇。但目前lncRNA 靶RNA 預測模型還存在一些不足之處需要深入研究，主要包括以下幾個方面：首先，由于lncRNAs的研究剛剛起步，相關的lncRNA 生物特性認識仍然不足。部分lncRNAs的生物特征如何表征還需進一步探索，如RNA的三級結構，其空間性特征可能在預測lncRNA 靶RNA 中發揮重要的作用。其次，實驗驗證的數據雖然可靠，但數據量較少，通過計算統計的方法來區分兩個RNAs 相互作用的規律還需要進一步研究，以提高預測模型的靈敏度。最后，還應考慮其他因素，包括RNA的染色體定位以及特定RNA 結合蛋白的存在。發現和引入RNA的更多特征是lncRNA 靶基因預測模型發展的一個重要方向。

隨著更多相關數據庫地建立和RNA-RNA 預測方法的改進，相信會有更多、更準確的lncRNARNA 相互作用預測模型被開發并應用于lncRNA調控功能的研究。總之，在全基因組水平上預測lncRNA-RNA 相互作用的關系來尋找lncRNA 靶RNA 仍然有巨大的進步空間。