利用轉(zhuǎn)錄組測(cè)序數(shù)據(jù)分析可變剪接的方法

馬素平

摘 要 可變剪接是調(diào)節(jié)基因表達(dá)和產(chǎn)生蛋白組多樣性的重要因素，同時(shí)參與調(diào)控細(xì)胞分裂、分化及凋亡等重要生物學(xué)過程，異常的可變剪接多與人類疾病有關(guān)。隨著新一代測(cè)序技術(shù)和生物信息學(xué)的快速發(fā)展，以及先進(jìn)計(jì)算方法的提出，使得我們對(duì)可變剪接有了深入的認(rèn)識(shí)。并且基于剪接機(jī)制對(duì)于病的靶向藥物設(shè)計(jì)，已得到了有效的臨床治療效果。本文主要闡述了近年來基于二代測(cè)序技術(shù)開發(fā)的幾種識(shí)別可變剪接的計(jì)算方法，并對(duì)未來的發(fā)展方向進(jìn)行展望。

關(guān)鍵詞 可變剪接;二代測(cè)序技術(shù);生物信息學(xué);分析工具

可變剪接，又稱選擇性剪接（Alternative Splicing，AS），是真核生物基因表達(dá)的普遍調(diào)節(jié)機(jī)制，是指一個(gè)前體mRNA經(jīng)過不同的剪接形式產(chǎn)生多種不同剪接異構(gòu)體的過程。在1978年，Walter Gilbert提出了內(nèi)含子和外顯子命名[1]，不同外顯子組合產(chǎn)生特異的異構(gòu)體。二代測(cè)序技術(shù)的迅速發(fā)展極大地推動(dòng)了人類對(duì)可變剪接的認(rèn)識(shí)?，F(xiàn)有數(shù)據(jù)表明，人類大約有92%-94%的基因都會(huì)經(jīng)歷某種程度的可變剪接行為，并且在20000多種人類蛋白編碼基因中，約37%的基因會(huì)編碼產(chǎn)生不同的蛋白亞型，這表明可變剪接增加了蛋白質(zhì)組的多樣性和復(fù)雜性[2]。

AS對(duì)基因的功能起著重要調(diào)控作用，同一基因的不同亞型可能參與不同的生物學(xué)過程。例如p53抑癌基因（TP53）在DNA受損細(xì)胞的調(diào)控中起著核心作用，然而其Δ133β亞型則可以抑制全長p53β亞型5和6從而誘導(dǎo)腫瘤細(xì)胞的凋亡[3]。另外AS幾乎參與了所有生物學(xué)過程，包括調(diào)節(jié)細(xì)胞的分裂和凋亡、神經(jīng)系統(tǒng)的發(fā)育以及細(xì)胞對(duì)抗多種環(huán)境因素做出的免疫應(yīng)激反應(yīng)等[4]。

另一方面，AS的異常調(diào)節(jié)還與多種遺傳性疾病和惡性腫瘤相關(guān)，包括神經(jīng)退行性疾病、心血管疾病和代謝狀況等。據(jù)報(bào)道，與SNP相關(guān)的遺傳性疾病多達(dá)一半是由于剪接受損引起的[5]。AS的異常調(diào)節(jié)對(duì)癌癥的發(fā)生發(fā)展有重要的作用，為疾病的發(fā)展提供了可能的新穎治療靶標(biāo)和生物標(biāo)志物的來源，而AS位點(diǎn)的預(yù)測(cè)可以為藥物設(shè)計(jì)提供很好的分子基礎(chǔ)。

本文主要闡述了最新的基于二代測(cè)序數(shù)據(jù)計(jì)算可變剪接方法的基本原理及應(yīng)用，并對(duì)未來的研究方向做了初步的展望。

1可變剪接的主要類型

目前大量數(shù)據(jù)研究結(jié)果表明，可變剪接主要包括五種形式（圖1），分別為外顯子跳躍（Skipped Exon，SE），可變5′剪接位點(diǎn)（Alternative 5′ Splice Site，A5SS），可變3′剪接位點(diǎn)（Alternative 3′ Splice Site，A3SS），互斥外顯子（Mutually Exclusive Exons，MXE）和內(nèi)含子保留（Retained Intron，RI）。此外還有兩種不常見的形式：可變的第一個(gè)外顯子（Alternative First Exon，AFE）和可變的最后一個(gè)外顯子（Alternative Last Exon，ALE）。其中SE是人類基因組中最主要的可變剪接模式（約35%），其次是A3SS（16%）與AA5SS（15%），RI發(fā)生最少，僅占到全部事件的1%[6]。

2分析可變剪接的軟件

可變剪接現(xiàn)象從70年代發(fā)現(xiàn)后，其基本的科學(xué)問題聚焦于剪接位點(diǎn)的發(fā)現(xiàn)、差異分析、調(diào)控元件和網(wǎng)絡(luò)的發(fā)現(xiàn)和構(gòu)建等。RNA-seq技術(shù)的發(fā)明，使系統(tǒng)、定量的可變差異剪接分析成為可能。大量測(cè)序數(shù)據(jù)的可變剪接差異分析需要優(yōu)秀的統(tǒng)計(jì)模型和計(jì)算工具，因此一直是生物信息學(xué)領(lǐng)域的重要研究課題之一。

2.1 轉(zhuǎn)錄組讀段回帖

序列比對(duì)軟件STAR（Spliced Transcripts Alignments to a Reference，STAR）[7]是專門針對(duì)RNA-seq數(shù)據(jù)開發(fā)的，該算法使用未壓縮后綴陣列中的連續(xù)最大可比種子進(jìn)行搜索，接著進(jìn)行種子聚類（Seed clustering）和縫合過程（Stitching procedure）。除了典型剪接的無偏從頭檢測(cè)外，STAR還能發(fā)現(xiàn)非典型拼接和嵌合（融合）轉(zhuǎn)錄本，并且能夠比對(duì)全長RNA序列。STAR是目前轉(zhuǎn)錄組測(cè)序分析最常用的軟件之一，由C++編寫，運(yùn)行速度快，同時(shí)也是免費(fèi)開源軟件。

TopHat[8]是基于Bowtie[9]短序列比對(duì)軟件開發(fā)的，它可以實(shí)現(xiàn)對(duì)短序列的快速比對(duì)，根據(jù)比對(duì)的結(jié)果可以檢測(cè)到外顯子間的剪接位點(diǎn)。該算法基本思想是對(duì)Bowtie比對(duì)結(jié)果中大量重疊的片段尋找標(biāo)志性的剪接位點(diǎn)，并列出它們所有可能的組合情況。此外對(duì)于Bowtie沒有比對(duì)到基因組上的片段，TopHat會(huì)重新建立索引進(jìn)行比對(duì)，并與上個(gè)過程得到的所有可能組合進(jìn)行對(duì)比，最終將跨越剪接位點(diǎn)的片段比對(duì)到基因組。

2.2 比對(duì)后的可變剪接分析

rMATS（ replicate multivariate analysis of transcript splicing）是目前最常用的分析RNA-seq數(shù)據(jù)可變剪接事件的軟件[10]。該算法依賴于有重復(fù)的數(shù)據(jù)集，利用分層模型同時(shí)考慮重復(fù)樣本的不確定性和變異性。以SE為例，rMATS首先采用Exon Inclusion level 來定義樣本中可變剪接事件的表達(dá)量，正常的isoform稱為Exon Inclusion Isoform， 發(fā)生了外顯子跳躍的轉(zhuǎn)錄本稱之為Exon Skipping Isoform，比對(duì)到Inclusion? isoform上的reads用L表示，比對(duì)到Skipping Isoform上的reads用S表示（圖2），則該外顯子跳躍的定量如下：

Exon Inclusion level為Inclusion isoform所占的比例，計(jì)算時(shí)用長度（LI和Ls）校正了原始的reads數(shù)。給定閾值c，? 判斷兩個(gè)樣本中相對(duì)應(yīng)的inclusion level是否發(fā)生了變化：c，c用戶可以自定義，取值范圍為0～1，表示兩個(gè)樣本中Inclusion level的差值。如c=0.1，則表示兩個(gè)樣本中該可變剪接事件的Inclusion level相差10%。其他類型的可變剪接事件也可以劃分成上述兩種isoforms，進(jìn)行差異分析比較的是兩組樣本中inclusion level的差異。

Exon Skipping Isoform為Upstream exon 和Downstream exon直接連接形成;Exon Inclusion Isoform為Upstream exon， Alternative exon和Downstream exon連接形成。

該模型以Likelihood-ratio test計(jì)算p值，大大提升了計(jì)算速度。rMATS支持多線程運(yùn)行且支持兩種輸入格式：Fastq或者Bam。根據(jù)計(jì)算時(shí)用到的reads差別，最后會(huì)得到兩組結(jié)果，一種是只用到跨Junction的reads;另一種是比對(duì)到剪接位點(diǎn)上的所有reads。rMATS是目前在RNA-seq數(shù)據(jù)領(lǐng)域應(yīng)用最多的分析可變剪接的工具。

DARTS[11]是2019年3月由邢毅教授團(tuán)隊(duì)在開發(fā)rMATS工具后，新提出的一個(gè)基于深度學(xué)習(xí)檢測(cè)可變剪接的計(jì)算框架，這也是首次將深度學(xué)習(xí)與貝葉斯假設(shè)檢驗(yàn)相結(jié)合，用于RNA的可變剪接分析。DARTS由兩部分構(gòu)成（圖3）：深度神經(jīng)網(wǎng)絡(luò)（DNN）和貝葉斯假設(shè)檢驗(yàn)（BHT）。DNN模型是基于順式序列特征和RNA結(jié)合蛋白（RBP）的mRNA表達(dá)水平預(yù)測(cè)差異剪接的結(jié)果;BHT統(tǒng)計(jì)模型通過將經(jīng)驗(yàn)證據(jù)整合到特定的RNA-seq數(shù)據(jù)集中，并結(jié)合先驗(yàn)的差異可變剪接概率來進(jìn)行預(yù)測(cè)。在訓(xùn)練過程中，DARTS BHT使用無先驗(yàn)概率信息（DARTS BHT（flat），僅將RNA-seq數(shù)據(jù)用于推斷）分析大規(guī)模RNA-seq數(shù)據(jù)，以生成高置信度差異及未發(fā)生剪接的訓(xùn)練標(biāo)簽之間的事件，以此作為DARTS DNN的訓(xùn)練輸入。在應(yīng)用過程中，已訓(xùn)練完成的DARTS DNN用于預(yù)測(cè)用戶特定數(shù)據(jù)集中的差異可變剪接，然后與DARTS BHT（DARTS BHT（info））觀察到的RNA-seq片段的計(jì)數(shù)結(jié)合到一起，以達(dá)到深度學(xué)習(xí)增強(qiáng)的剪接分析。

DARTS BHT（flat）進(jìn)行常規(guī)分析大規(guī)模RNA-seq數(shù)據(jù)中的可變剪接事件，創(chuàng)建帶標(biāo)簽的訓(xùn)練數(shù)據(jù)，用于訓(xùn)練 DNN模型;新的特定RNA-seq經(jīng)DNN模型預(yù)測(cè)作為貝葉斯模型的先驗(yàn)（DARTS BHT（info））;用戶的RNA-seq數(shù)據(jù)則是用于更新先驗(yàn)概率形成后驗(yàn)概率。順式序列特征（Cis-sequence）和反式RBP的mRNA水平（Trans-RBP）：DARTS DNN預(yù)測(cè)差異可變剪接的兩個(gè)因素。先驗(yàn)信息（Prior）：DARTS DNN預(yù)測(cè)的結(jié)果?？赡苄裕↙ikelihood）：特定RNA-seq數(shù)據(jù)集觀察到的讀數(shù)。差異剪接的后驗(yàn)概率（Posterior）由Cis-sequence和Trans-RBP生成的先驗(yàn)概率和Likelihood組成。

DARTS的深度神經(jīng)網(wǎng)絡(luò)部分整合了剪接位置附近的序列信息、進(jìn)化信息、可變剪接產(chǎn)生的RNA二級(jí)結(jié)構(gòu)信息和RBP的表達(dá)信息等，提高了該模型的準(zhǔn)確性。相較于之前的例如MISO[12]、SpliceTrap[13]、SUPPA[14]等軟件，DARTS對(duì)測(cè)序深度的要求降低?？偟膩碚f，DARTS綜合了深度學(xué)習(xí)和貝葉斯假設(shè)檢驗(yàn)統(tǒng)計(jì)模型的優(yōu)點(diǎn)，為低測(cè)序深度的數(shù)據(jù)提供了更好的分析可變剪接的方法，拓展了傳統(tǒng)RNA-seq可變剪接分析的敏感度和準(zhǔn)確度。

3總結(jié)與展望

可變剪接（AS）是基因在轉(zhuǎn)錄后水平的一項(xiàng)復(fù)雜調(diào)控機(jī)制。在最近的幾十年中，我們對(duì)AS的機(jī)理和動(dòng)力學(xué)的了解呈指數(shù)級(jí)增長，人們也越來越認(rèn)識(shí)到AS在人類健康和疾病中發(fā)揮的關(guān)鍵作用。隨著RNA-seq和各種計(jì)算工具的不斷涌現(xiàn)，使得我們可以在全基因組水平檢測(cè)AS現(xiàn)象。同時(shí)隨著近年來單細(xì)胞測(cè)序技術(shù)的發(fā)展，可能會(huì)在單個(gè)細(xì)胞中出現(xiàn)新的剪接現(xiàn)象，這需要我們?nèi)ラ_發(fā)專門針對(duì)單細(xì)胞測(cè)序的更優(yōu)計(jì)算方法。另外機(jī)器學(xué)習(xí)方法將會(huì)加快我們對(duì)可變剪接的認(rèn)識(shí)[15]。最后，可變剪接變體在功能上有助于疾病的耐藥性，因此靶向可變剪接可作為某些疾病的新型療法。目前在人類某些疾病（如SMA）中已取得了很好的臨床治療效果[16]，但在癌癥方面還未有重大突破，相信隨著研究的不斷深入，在不久的將來也能為癌癥提供有效的可選藥物。

參考文獻(xiàn)

[1] Gilbert，Walter. Why genes in pieces？[J]. Nature，1978， 271（5645）：501.

[2] Min-Sik， K. A draft map of the human proteome.[J].Nature， 2014，509（7502）：575-581.

[3] Aoubala M，Murray-Zmijewski F，Khoury M P，et al. P53 directly transactivates Δ133p53α， regulating cell fate outcome in response to DNA damage[J]. Cell death and differentiation，2011，18（2）：248-258.

[4] Ramanouskaya T V，Grinev V V . The determinants of alternative RNA splicing in human cells[J]. Molecular Genetics & Genomics，2017，（3）：263.

[5] Garcia-Blanco M A，Baraniak A P，Lasda E L . Alternative splicing in disease and therapy[J]. Nature Biotechnology，2004，22（5）：535-546.

[6] Wang E T，Sandberg R，Luo S，et al. Alternative isoform regulation in human tissue transcriptomes[J]. Nature，2008，456（7221）：470-476.

[7] Alexander D ， Davis C A ， Felix S ， et al. STAR： ultrafast universal RNA-seq aligner[J]. Bioinformatics，，2013（1）：1.

[8] Trapnell C，Pachter L，Salzberg S L . TopHat： discovering splice junctions with RNA-Seq[J]. Bioinformatics，2009，25（9）：1105-1111.

[9] Langmead B，Trapnell C，Pop M，et al. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome[J]. Genome Biology，2009，10（3）：1-10.

[10] Shihao Shen，Juw Won Park，Zhixiang Lu，et al. rMATS： Robust and flexible detection of differential alternative splicing from replicate RNA-Seq data[J]. proc natl acad sci u s a，2014，111（51）：5593-601.

[11] 佚名. Deep-learning augmented RNA-seq analysis of transcript splicing[J]. Nature Methods，2019，（1）：157.

[12] Yang-Ming S I，Yong-Qiang X，Lu C . Differential splicing event analysis of liver tumor-educated blood platelets RNA-seq data with Hisat2 and MISO[J]. journal of inner mongolia university of science and technology，2016，（3）：73.

[13] Zhang M Q . SpliceTrap： a method to quantify alternative splicing under single cellular conditions[J]. Bioinformatics，2011，27（21）：3010-3016.

[14] more samples. SUPPA： a super-fast pipeline for alternative splicing analysis from RNA-Seq[J]. Biorxiv，2014，（1）：257.

[15] Jaganathan K，Kyriazopoulou Panagiotopoulou S，Mcrae J F，et al. Predicting Splicing from Primary Sequence with Deep Learning[J]. Cell，2019，（3）：176.

[16] Cleveland，D. W. Antisense oligonucleotide therapy for neurodegenerative disease[J].The Journal of clinical investigation，2006， 116（8）：2290.