文章

Reconstruction of full-length circular RNAs enables isoform-level quantification

作者單位

中國科學院北京生命科學研究院計算基因組學實驗室趙方慶團隊（第一作者：博士研究生鄭毅&助理研究員冀培豐）

發(fā)表期刊

Genome Medicine

發(fā)表時間

2019年01月19日

引言

“君子生非異也，善假于物也?！边@句話很適合該文的方法CIRI-full，它巧妙地應(yīng)用了i)環(huán)狀RNA的本性——環(huán)狀且較短（大部分約300~500bp）; ii)環(huán)化建庫以及更長雙末端測序（PE250, PE300）——環(huán)狀RNA序列全覆蓋。

名詞解釋

FRJ：forward splice junction，前向剪接位點，比如 mRNA 是pre-mRNA 通過剪接體復合物通過前向剪切去除內(nèi)含子將外顯子連在一起的產(chǎn)物。

BRJ：back splice junction，反向剪接位點，比如環(huán)狀RNA通過剪接體復合物將供體3’末端外顯子與受體外顯子5’末端連接在一起的產(chǎn)物。

RO：reverse overlap，反向交疊，如圖灰色陰影所示，表示雙末端測序（paired-end sequencing） Read1 與 Read2 的5’末端或3’末端存在反向一致序列。

Cirexon：circRNA’s exon，構(gòu)成環(huán)狀RNA的外顯子

Isoform：來自同一基因組區(qū)間的環(huán)狀RNA轉(zhuǎn)錄本

FSG：forward splice graph，前向剪接圖

文章亮點

?  提出了識別環(huán)狀RNA全長序列的方法CIRI-full

?  相比基于BSJ的方法，RO方法對低表達豐度的環(huán)狀RNA 更敏感

?  利用FSG的方法精確識別和量化環(huán)狀RNA的轉(zhuǎn)錄本（isoform）

?  6個物種（包括人，恒河猴，小鼠，大鼠，兔，雞）大腦全長環(huán)狀RNA

?  相比BSJ，CIRI-full獲得的轉(zhuǎn)錄本（isoform）能過濾假陽性的差異表達環(huán)狀RNA

數(shù)據(jù)資源

測序數(shù)據(jù)：PRJNA475651

軟件：https:///projects/ciri

環(huán)狀RNA研究的攔路虎——無法準確捕獲全長序列

環(huán)狀RNA 功能的重要性毋庸置疑，然而對其功能的探究仍舉步維艱，主要原因除了人們?nèi)跃窒抻谘芯烤€性RNA分子的思維模式，另一個主要原因就是無法大規(guī)模獲得環(huán)狀RNA的全長序列并準確定量。我們都知道，從序列同源以及物種進化地角度能夠有效地探究一個遺傳分子可能存在地功能；另外，一個遺傳分子表達的準確定量是反應(yīng)其在生物體內(nèi)發(fā)揮功能的有效途徑。然而，大規(guī)模從 RNA-seq 中識別全長 circRNA 的方法仍待開發(fā)。這篇文章就提出了一個新的方法CIRI-full從RNA-seq中探究環(huán)狀RNA轉(zhuǎn)錄本的全長序列以及表達定量。

實驗基礎(chǔ)

當插入片段長度（即測序文庫片段大?。┐笥诃h(huán)狀RNA序列長度，那么就可以觀察到雙末端 reads 的RO事件。

建庫：環(huán)化建庫，最后插入片段大小為300-800bp

測序平臺：Illumina HiSeq 2500 platform

Read長度：250 x 2, 300 x 2

方法概述

該方法總共分為四部分，包括：

?  RO read 的檢測與檢驗（圖A-D）

檢測雙末端 reads（paired-end reads）Read1與Read2的5’末端是否存在反向交疊RO（圖A-B），如果存在則作為候選RO-merged read；隨后分割比對RO-merged read 到參考基因組，最長的作為位置標簽——錨定位點anchor，區(qū)分異常以及無法比對的區(qū)域（圖C）；校正剪切信號GT/AG確定環(huán)狀RNA邊界（圖D）。

?  BSJ 以及 cirexon 檢測

利用 CIRI2 檢測反向剪切事件BSJ，并通過CIRI-AS識別BSJ中的可變剪接事件（single-splice event）。如果構(gòu)成BSJ的兩個reads只落在BSJ內(nèi)部的外顯子即cirexons，那么環(huán)狀RNA的全長序列就能由cirexons線性重構(gòu)。

?  組裝RO與BSJ的reads（圖E）

i） RO-merged read存在3’-RO或兩個末端落在同一個cirexon上

ii） 每個BSJ的所有雙末端reads都落在cirexons

以上兩種情況可以直接確定環(huán)狀RNA全長，而剩余的reads，CIRI-full結(jié)合兩者優(yōu)點——BSJ確定邊界，RO-reads確定內(nèi)部cirexons——進一步組裝全長。

?  轉(zhuǎn)錄本isoform重構(gòu)與定量（圖F）

CIRI-full用BSJ cirexons以及RO-merged reads都構(gòu)建了一個前向剪接圖FSG，利用深度優(yōu)先搜索算法（adapted depth-first search）將FSG分解成通道（path）——候選轉(zhuǎn)錄本（isoform），然后用蒙特卡洛仿真（Monte Carlo simulation）方法模擬每個通道（轉(zhuǎn)錄本）的表達豐度，最后用近乎窮舉的方式尋找每一條通道相對表達豐度（設(shè)置了每個環(huán)狀RNA最多10個轉(zhuǎn)錄本以提高效率）。

方法詳述

1. RO read 的檢測與檢驗

5’-RO識別策略 對于每對 Read1 與 Read2，Read1 與 Read2 其中一個的 5’ 末端的前 10 bp 被分成步長為 1 bp、窗口為 8 bp 的三份（如示例圖，假設(shè)為 Read1 的 5’ 端）。

這些子序列被用來作為種子搜索 Read2 的 5’ 末端，一旦所有種子匹配堿基數(shù) >= 7bp，那么從 Read2 的 5’ 最后一個堿基到 Read1 匹配到的那個位點的序列被提取出來，然后再比對到 Read1 上。如果滿足以下條件，則這對 reads 被認為 5’-RO：reads pair上的長度 >= 13 bp，堿基一致性 >= 95%。隨后，Read1 與 Read2 根據(jù)交疊部分合并為一條長 read 作為候選 RO-merged read以待進一步驗證。3’-RO與此類似。

過濾假陽性 候選RO-merged read 通過BWA-MEM比對確定在參考基因組上的位置，然后將比對長度進行排序，其中最長且比對得分 > 15 的序列被作為anchor（錨定位點）來確定 reads 的位置。如果anchor 的兩頭在基因組上的范圍 < 100kbp，那么將計算ro-merged="" read="" 比對上基因組的長度；長度如若超過了="" read 的一半，將用于再次做局部比對（local="" realignment）來確定ro-merged="" read在基因組上的精確位置——尋找bsj（如果沒有將用anchor），并用動態(tài)規(guī)劃算法計算無法比對以及異常比對的片段。之后將過濾假陽性ro-merged="" read，原則包括：read不包含bsj以及read兩端等長子序列比對位置沒有落在anchor="">

校正GT/AG 由于BWA-MEM無法獲得GT/AG剪切信號以及精確的剪切邊界，另外RO-merged reads可能包含套索結(jié)構(gòu)（lariat structure）。因此，需要對每條候選RO- merged read檢測和校正GT/AG位點，如果不存在將丟棄。另外，對每條read剪接位點（junction site）上下游5-bp比對質(zhì)量進行檢測，如果存在gap或錯配，也將丟棄不用。

2. BSJ 以及 cirexon 檢測

BSJ 用CIRI2檢測，而BSJ中的可變剪切事件（single-splice events）用CIRI-AS來推測。BSJ的檢測方法在所有識別軟件中都比較類似，簡單描述就是截取read的5’與3’末端的子序列比對到參考基因組上，如果比對方向相反將作為候選circRNA。詳細請訪問CIRI 以及CIRI2原文。每個BSJ中的cirexons通過CIRI-AS捕獲的剪切事件推測。

3. 組裝RO與BSJ的reads

基于所識別到的RO-merged reads 以及BSJ邊界信息，利用以下原則構(gòu)建全長環(huán)狀RNA。

?  如果RO-merged read存在3’-RO或兩個末端落在同一個cirexon上，那么該readd被認定為一個全長環(huán)狀RNA；否則等待進一步組裝。

?  如果每個BSJ的所有雙末端reads都落在cirexons上，那么該BSJ的全長環(huán)狀RNA由cirexons線性重構(gòu)；否則BSJ以及cirexons等待進一步組裝。

?  RO與BSJ實際上存在互補信息——BSJ獲得環(huán)狀RNA精確的邊界，而RO識別內(nèi)部結(jié)構(gòu)。

CIRI-full結(jié)合兩者優(yōu)點，將等待進一步組裝的RO-merged reads以及BSJ reads根據(jù)BSJ進行排序聚類，如果兩種read被發(fā)現(xiàn)落在同一個BSJ內(nèi)，那么這些reads將被用來重構(gòu)全長環(huán)狀RNA。另外，RO-merged reads被用來確定BSJ未識別的額外cirexons，如果所有reads都落在cirexons上，那么全長將通過 cirexon 線性組裝；如若不然，BSJ 所識別的 cirexons 將標記為部分重構(gòu)的 circRNA。

4. 轉(zhuǎn)錄本isoform重構(gòu)與定量

對每一個全長環(huán)狀RNA，CIRI-full用BSJ cirexons以及RO-merged reads都構(gòu)建了一個前向剪接圖FSG，其中節(jié)點表示cirexons，邊表示兩個外顯子之間存在前向連接。

理論上，F(xiàn)SG包含了環(huán)狀RNA所有存在的轉(zhuǎn)錄本，然而值得注意的是，由于環(huán)狀RNA本身成環(huán)特性，所以FSG實際上是一個閉環(huán)。利用深度優(yōu)先搜索算法（adapted depth-first search）將FSG分解成通道（path），迭代地從每個節(jié)點開始并在斷點或起始節(jié)點（cirexon）結(jié)束；其中短通道被合并成更長的通道，而冗余的通道被過濾，為了避免過多可能的假陽性，CIRI-full設(shè)置了固定的最大通道數(shù)（默認10）。

詳細步驟：

1> FSG的邊包括4大類：

?  BSJ

?  Phasing FSJ(forward splice junction)，表示只有一個環(huán)狀轉(zhuǎn)錄本

?  Co-occurred FSJs，表示剪切事件的數(shù)目與RO reads的數(shù)目相同

?  剩余的FSJ

其中，包含phasing FSJ和Co-occurred FSJs的通道優(yōu)先級最高，被稱為phased isoforms。而對第四類剩余的FSJs，根據(jù)節(jié)點的測序深度排序，只有top10被保留。所有保留的通道將作為候選轉(zhuǎn)錄本。

2> 為了確定每條通道的相對豐度，用蒙特卡洛仿真（Monte Carlo simulation）方法根據(jù) RNA-seq 插入片段的長度分布（通過 paired-end reads 在參考基因組上的距離推測）模擬每條通道的 BSJ-reads 分布。

3> 為了量化每條通道的相對豐度，CIRI-full 使用了窮舉法（approximate exhaustive search algorithm）：先為每條通道設(shè)置一個相對豐度（正整數(shù)），并且所有通道豐度和等于BSJ-reads總數(shù)?；诔跏枷鄬ωS度以及模擬BSJ-reads，每條通道節(jié)點與邊的累積豐度被計算；根據(jù)mapped BSJ-reads計算節(jié)點與邊的真實豐度；隨后，計算模擬值與真實值的距離。隨后迭代地跟新模擬值，知道距離趨于收斂。

通過以上步驟，最終獲得每個環(huán)狀轉(zhuǎn)錄本的相對豐度。