天天看點
傳回

三代全長轉錄組測序群組裝 - BPSO_mynotes

三代全長轉錄組測序群組裝

這些三代全長轉錄組測序的相關問題，可以幫到愛學習的我哦

1.什麼是三代全長轉錄組測序

三代全長轉錄組測序，即利用PacBio三代測序平台對某一物種的mRNA進行測序研究。它以平均超長讀長10-15kb的優勢、結合多片段文庫篩選技術，實作了無需拼接的轉錄本分析，克服了傳統二代轉錄組Unigene拼接較短、轉錄本結構不完整的缺陷，也由于其可直接獲得單個RNA分子從5’端到3’端的高品質全部轉錄組資訊而得名。

2.為什麼要做全長轉錄組測序？

轉錄本非常多樣和複雜，絕大多數基因不符合“一基因一轉錄本”的模式，這些基因往往存在多種剪切形式。通過二代測序，我們可以很準确地進行基因的表達及定量的研究，但是受限于讀長的限制，不能得到全長轉錄本的資訊。

基于二代測序平台的轉錄組産品，首先是把RNA打成小的短片斷進行測序，然後再通過生物資訊的方法進行拼接，将拼接後的序列傳遞給客戶。但是基于二代測序平台的轉錄組，由于讀長的限制（PE150），在轉錄本組裝的過程中會存在較多的嵌合體，并且不能準确地得到完整轉錄本的資訊，進而會大大降低表達量、可變剪接、基因融合等分析的準确性。圖1. 二代和三代轉錄組測序原理及讀長對比

目前基于PacBio的單分子實時測序技術，目前平均讀長已經達到10Kb以上，最長可達80Kb，轉錄組測序不再需要組裝，就可以直接得到全長轉錄本的資訊。

3.二代與三代轉錄組相比兩者分别有哪些優劣勢？三代轉錄組具體優勢可否說明？

表1. 二代和三代轉錄組測序優劣勢對比

從上述對比表格可看出，兩種轉錄組測序技術互有優劣勢，是以在給各位老師在設計課題時，建議老師二+三代轉錄組測序技術同時使用，保證結構準确性、序列完整性及序清單達量準确性，達到資料的最優利用效果以及成本效益最高。

三代轉錄組具體優勢說明如下：

a.超長讀長（平均讀長10-15K，最長讀長80K），可一次将真核生物的全長轉錄本資訊讀取完整；

b.無需進行片段打斷和拼接，避免出現組裝錯誤；

c.基于全長轉錄組測序得到的完整準确的轉錄本資訊，結合二代資料，友善識别特異性表達且做更加精确的基因和轉錄本表達定量。

d.針對有參考基因組的物種，全長轉錄組資訊可以糾正基因組的錯誤組裝、更準确地發現新的轉錄本和基因、分析基因融合事件等。

e.無需鍊特異性建庫，全長轉錄組測序可直接擷取正義鍊、反義鍊及部分LncRNA資訊。

圖2. 三代全長轉錄組測序優勢概覽

4.哪些物種适合做三代全長轉錄組測序？獲得這些全長mRNA資訊有何用處？

無參考基因組物種和有參考基因組物種均适用。

a.對于沒有參考基因組的物種

由于基因組測序成本高，缺乏基因組參考資訊在很大程度上限制了對物種的深入研究。通過三代全長轉錄組測序來建構物種Unigene庫，無需進行序列組裝，就可以獲得該物種轉錄組水準的參考序列（轉錄組水準的參考基因組），為後續研究提供很好的遺傳資訊基礎。

獲得這些全長轉錄本資訊，可以更準确地進行CDS和SSR分析。如果有同一批樣本的二代資料，不但可以提高三代測序資料的使用率，同時可以對這些全長的轉錄本進行更精準的定量分析。

b.對于有不完善參考基因組的物種

參考基因組組裝不完善是普遍存在的問題，特别是多倍體這類物種，給科研工作帶來了極大阻礙。參考基因組組裝不完善，用二代測序會導緻reads比對率低，基因表達定量不準确的問題。用全長轉錄組測序技術可直接獲得轉錄本全長序列，再結合二代測序，會使定量更準确，資料使用率更高，同時基于全長轉錄組資料，可以優化基因結構，進而輔助基因組組裝和注釋。

c.對于具有較好參考基因組的物種

利用三代全長轉錄組測序獲得的資訊是生物體内直接存在的，比基于參考基因組預測到的轉錄組資訊更準确，同時也可準确鑒定基因的可變剪接、融合基因、基因家族和非編碼RNA等資訊。

如果有同一樣本的二代資料，不但可以提高三代資料的使用率，同時還可以深入研究某基因可變剪接形成的不同轉錄本的表達差異。可以确定不同發育階段或不同處理情況下，該基因中高表達轉錄本以及低表達轉錄本。不同樣品的融合基因和等位基因差異，也同樣可以分析。

需要注意的是全長轉錄組測序隻能得到轉錄本全長序列，不可進行基因表達定量。

5. 全長轉錄組測序那麼貴，如何更大程度上降低測序成本？

由于轉錄組資訊呈動态變化且存在組織差異，單一組織得到的全長轉錄本對該物種其他部位組織可能不是很全面或不太适用，是以用一個物種不同部位組織混樣進行高深度測序（針對不同要求及目的，推薦8G、10G和12G等），會得到比較理想的參考轉錄本庫資訊，也是降低測序成本的理想方法。

6. 三代全長轉錄組測序如何選擇測序樣本？

總體原則是根據研究目的進行選擇，舉例說明如下：

a.單個三代轉錄組項目：

①　如果想要獲得該物種相對全面的轉錄本資訊，建議對該物種的不同部位混合取樣；

②　如果隻想研究某個特定的組織部位，建議在不同發育時期對特定組織部位進行取樣；

b. 二+三代轉錄組混合政策項目：

三代轉錄組與二代轉錄組測序取樣部位或時期相對應的同一批樣品，等量RNA混樣測序；

c. 多個三代轉錄組樣品項目：

如果想要研究某物種脅迫處理（其他生物或非生物脅迫都适用）前後變化，建議取對照和處理組（至少兩個樣品）進行對比分析；

①　全長轉錄組混樣測序為了保證資料來源的均一性，一定要等量RNA混合測序，而非等量樣品混合抽取RNA再測序。

②　随着三代轉錄組測序成本逐漸下降，多個三代轉錄組樣品測序的正常時代也即将到來。

7. 全長轉錄本資料量和文庫類型如何确定？

推薦資料量大小需依據物種的複雜程度、基因大小及研究目的來确定。根據已有的項目經驗、資料庫資訊及文章中報道，我們詳細推薦如下：

表2. 推薦性全長轉錄本測序資料量和文庫類型

注：對于全長轉錄組測序，資料量并不是固定的，針對同一物種同一研究目的，測序資料量越多，檢測到的全長轉錄本也會越全面。

8. 全長轉錄組測序為什麼要建3-4個分段文庫？不同文庫資料産出比例如何？

建構分段文庫，是由PacBio平台測序原理所決定。在三代轉錄組測序過程中，建構好的全長文庫需要loading到測序小孔——零模波導孔（ZMWs）中，由于mRNA長度不同，在loading的過程中會出現一定的loading bias，即測序小孔會優先被長度較短的片段占據，每個測序小孔隻能容納一個文庫分子，而大部分長片段則沒有測到。是以為盡量降低loading bias的影響，需要根據測序物種mRNA的長度進行分段，使一個文庫中的序列長度控制在一個較窄的範圍内。故建構分級文庫越多，也會得到更全面的全長轉錄本。

全長轉錄組測序一般推薦至少建構三種文庫類型，1-2Kb、2-3Kb和≥3Kb文庫，資料産出比例為3:2:2。（例如：測8G的資料量，三個文庫分别測3G、3G和2G，也可以根據不同物種調整不同文庫的資料量）；建構1-2Kb、2-3Kb、3-6Kb、≥6Kb四個文庫（例如：測12G的資料量，四個文庫分别測4G、4G，2G和2G。資料量分布一般是2:2:1:1或3:2:2:1。

注：根據甜菜三代全長轉錄組文獻中報道還進一步驗證了一個常識，多數原本轉錄本3\'UTR+5\'UTR長度＞1Kb，是以一般不建議建構＜1Kb文庫，但研究目的是為了獲得較為全面的轉錄本時才會建議建構＜1Kb或0.5-1Kb文庫。

9. 三代全長轉錄組建庫測序的流程是什麼？

圖3. 三代轉錄組建庫測序流程簡圖

簡述以上流程：

a.全長cDNA合成：使用Clontech SMARTer PCR cDNA Synthesis Kit合成全長的cDNA；

b.片段選擇及PCR擴增：采用BluePippinTM儀器直接進行片段篩選并進行擴增；

c.SMRTbell文庫制備：将不同插入片段cDNA加上SMRTbell接頭，并完成文庫建構；

d.測序：文庫進行質控後上機三代平台PacBio測序。

10. 三代全長轉錄組測序的生物資訊分析流程是什麼？具體有哪些分析内容？

圖4. 三代全長轉錄組生物資訊分析流程簡圖

表3. 三代全長轉錄組資訊分析内容

有參考基因組物種

無參考基因組物種

(1) 原始資料處理及過濾；

(2) 測序資料品質評估；

(3)全長轉錄本判定；

(4)轉錄本聚類校正；

(5)與參考基因組序列比對；

(6)全長轉錄本比對注釋；

(7)基因結構優化；

(8)可變剪接鑒定；

(9)新基因預測及CDS預測

(10)LncRNA預測

(11)基因融合鑒定

(1)原始資料處理及過濾;

(2)測序資料品質評估；

(3)全長轉錄本判定；

(4)轉錄本聚類與校正；

(5)全長轉錄本比對注釋；

(6)預測編碼蛋白框（CDS）；

(7) SSR預測；

(8) LncRNA預測；

11. 三代全長轉錄組測序獲得的Clean Reads中，全長序列所占比例是多少？

全長序列所占比例與測序量和建庫長度以及表達量有關。沒有準确的标準，一般全長比例會占到50%左右（與目前文獻報道及官網測試資料水準一緻）。

pacbio 三代全長轉錄組資料分析流程 Iso-Seq 3

Iso-Seq 建庫

Iso-Seq的建庫方案有如下三類：

  1. 整個庫都是一個樣品的全長轉錄組，不需要加barcode區分樣品
  2. 不同樣品的全長轉錄組，加上不同barcode ，可以放在一起進行建庫測序
  3. 一些靶向獲得的部分基因也可以進行全長轉錄組的測序
三代全長轉錄組測序群組裝 - BPSO_mynotes

Iso-Seq 3進行資料分析

Iso-Seq3 進行全長轉錄組的分析，運作流程如下圖所示：

三代全長轉錄組測序群組裝 - BPSO_mynotes

1 ccs(Circular Consensus Calling)

ccs 擷取一緻性的序列，要求每一條測序的reads至少有一端含有引物序列。

ccs test.subreads.bam ccs.bam  --noPolish --minPasses 1 

2 測序引物和barcode的去除

這一步是采用lima 來完成的，這個軟體也是官方最新開發出來的程式，速度和準确度都較以往的算法有了很大的提升。還能夠基于barcode序列區分不同的樣品，

lima ccs.bam barcoded_primers.fasta  demux.ccs.bam --isoseq --no-pbi

3 聚類(cluster)

聚類采用IsoSeq3 軟體來完成，這一步會首先将ployA尾巴去除掉，并将連環結構去除掉，再對相似的序列進行聚類，最好形成全長的reads。

isoseq3 cluster demux.primer_5p--primer_3p.bam  unpolished.bam --verbose

4 抛光(polish)

這一步主要是将聚類的轉錄本，合并成一個完整的一緻性序列。

isoseq3 polish -j 16 unpolished.bam  test.subreads.bam   polished.bam

三代測序技術概述：

  1. PacBio和Oxford Nanopore測序的原理
  2. 三代測序的特點和應用
  3. 三代測序在轉錄組研究的優勢和案例分享
三代全長轉錄組測序群組裝 - BPSO_mynotes

三代測序基本分析流程

三代全長轉錄組測序群組裝 - BPSO_mynotes
  1. 原始測序序列去除接頭和錯誤序列
  2. 提取環形一緻序列讀長(CCS reads)
  3. CCS reads分類(包括全長和非全長CCS reads)
  4. CCS reads聚類(根據CCS reads序列的相似性)獲得最終的轉錄本集合
  5. 最終轉錄本比對回基因組
  6. 轉錄本定量和可變剪接分析

PacBio測序平台基于其獨特的單分子實時測序技術（Single Molecule Real Time，SMRT），通過其超長讀長，均一的覆寫度，高度的一緻性準及确性提供無與倫比的遺傳資訊深度解析。

用 GMAP/GSNAP軟體進行RNA-seq的alignment

首先需要參考基因組：雖然軟體本身提供了一個hg19的參考基因組，并且已經索引好了Human genome, version hg19 (5.5 GB)(http://research-pub.gene.com/gmap/genomes/hg19.tar.gz) ，但是下載下傳很慢，而且不是對所有版本的GSNAP都适用。是以我這裡對我自己的參考基因組進行索引。

gmap_build -D ./ -d  my_hg19.fa

然後取ensemble下載下傳hg19的gtf檔案。

然後還需要把自己下載下傳的gtf檔案也建構索引，需要兩個步驟

cat my_hg19.gtf |  ~/software/gmap-2011-10-16/util/gtf_splicesites > my_hg19.splicesites

cat  my_hg19.splicesites  |   iit_store -o my_hg19.gtf.index

GMAP最早用于講EST/cDNA序列比對到參考基因組上，可以用于基因組結構注釋。後來高通量測序時代，又開發了GSNAP支援高通量資料比對，這篇文章主要介紹GMAP，畢竟高通量轉錄組資料比對大家更喜歡用STAR, HISTA2等軟體。

軟體安裝

下面是我源碼安裝的代碼

wget http://research-pub.gene.com/gmap/src/gmap-gsnap-2018-07-04.tar.gz
tar xf gmap-gsnap-2018-07-04.tar.gz
cd gmap-2018-07-04/
./configure --prefix=$HOME/opt/biosoft/gmap
make -j 20

軟體使用

如下步驟假設你有一個物種的基因組序列和對應的CDS序列，分别命名為"reference.fa"和"cds.fa"

第一步:建構GMAP/GSNAP索引資料庫

GMAP/GSNAP對FASTA檔案中每個記錄下的序列的長度有一定限制, 每一條不能超過4G, 能應付的了大部分物種了。

建構索引分為兩種情況考慮，第一種是一個fasta檔案包含所有的序列

~/opt/biosoft/gmap/bin/gmap_build -d reference reference.fa

第二種則是每個染色體的序列都單獨存放在一個檔案夾裡，比如說你下載下傳人類參考基因組序列解壓後發現有N多個fasta檔案, 然後你就想用其中幾條染色體建構索引

~/opt/biosoft/gmap/bin/gmap_build -d reference Chr1.fa Chr2.fa Chr3.fa ...

注: 這裡的

-d

表示資料K庫的名字，預設把索引存放在gmap安裝路徑下的share裡，可以用

-D

更改.此外還有一個參數

-k

用于設定K-mer的長度, 預設是15, 理論上隻有大于4GB基因組才會有兩條一摸一樣的15bp序列(當然是完全随機情況下)。

第二步:正式比對

建立完索引之後就可以将已有的CDS或者EST序列和參考基因組序列進行比較。

~/opt/biosoft/gmap/bin/gmap -t 10 -d reference -f gff3_gene cds.fa > cds_gene.gff3

其中

-t

設定線程數,

-d

表示參考基因組資料庫的名字, 都是正常參數。我比較感興趣的參數是如何将序列輸出成GFF格式. GMAP允許多種格式的輸出，比如說

-S

隻看聯配的總體情況，而

-A

會顯示每個比對上序列的聯配情況, 還可以輸出蛋白序列(-P)或者是genomic序列(-E). 但是做結構注釋要的gff檔案，參數就是

-f gff3_gene

,

-f gff3_match_cdna

,

-f gff3_match_est

。

參考文獻

要想對一個軟體有更好的認識，最好還是看看他們文章是怎麼說的。

  • GMAP: a genomic mapping and alignment program for mRNA and EST sequences

    Bioinformatics 2005 21:1859-1875 Abstract Full Text, Thomas D. Wu and Colin K. Watanabe

  • Fast and SNP-tolerant detection of complex variants and splicing in short reads

    Bioinformatics 2010 26:873-881 AbstractFull Text, Thomas D. Wu and Serban Nacu

PLEK:轉錄本蛋白編碼潛能預測工具

在之前的文章中，我們介紹過CPC和CNCI這兩款軟體，可以用于預測lncRNA序列。其中CPC基于序列比對的方式，對于注釋資訊相對全面的物種分類效果較好，但是運作速度相對較慢，CNCI基于序列的三聯體堿基組成來區分編碼和非編碼轉錄本，對于注釋資訊缺乏的物種，效果也不錯，但是當序列中存在插入缺失時，其分類效果就變得很差。

在高通量測序産生的資料中，會存在一定的測序錯誤，雖然比例很低，但是基于這樣的序列組裝得到轉錄本然後去預測lncRNA, 對于CNCI這個軟體而言，就會造成相當大的影響。

為了克服上述問題，需要一款運作速度又快，又可以一定程度上降低測序錯誤影響的lncRNA預測軟體，PLEK軟體就是基于這樣的出發點進行開發的。PLEK軟體通過序列的kmer構成來區分編碼和非編碼轉錄本，不需要通過比對來完成，是以運作速度較快，同時其性能受到測序錯誤的影響的機率較低，比較穩定。

可以看到PLEK的運作速度是最快的。該軟體的源代碼托管在sourceforge上，網址如下

https://sourceforge.net/projects/plek/files/

安裝方式如下

wget https://sourceforge.net/projects/plek/files/PLEK.1.2.tar.gz
tar xzvf PLEK.1.2.tar.gz
cd PLEK.1.2
python PLEK_setup.py

基本用法如下

python PLEK.py \
-fasta transcript.fa \
-out output \
-thread 10

隻需要輸入轉錄本對應的fasta格式的檔案就可以了，輸出檔案

output

内容示意如下

三代全長轉錄組測序群組裝 - BPSO_mynotes

第一列代表該轉錄本為

coding

還是

non-coding

, 第二列為打分值，打分值大于0為coding, 小于零為non-coding, 第三列為fasta檔案中的序列辨別符。

預設情況下會調用内置的svm模型，如果你有該物種已知的mRNA和lncRNA轉錄本序列，也可以建構自己的模型，代碼如下

python PLEKModelling.py \
-mRNA mRNAs.fa \
-lncRNA lncRNAs.fa \
-prefix 20190129

運作成功後，會生成字尾為

.model

和

.range

的兩個檔案。在預測時可以通過參數指定svm模型，用法如下

python PLEK.py \
-fasta transcript.fa \
-out output \
-model 20190129.model
-range 20190129.range \
-thread 10

來源：

https://www.cnblogs.com/wangprince2017/p/10852391.html

https://www.omicsclass.com/article/344

http://blog.sina.com.cn/s/blog_182e021ed0102xbqz.html

 https://www.jianshu.com/p/3f331861c364

https://cloud.tencent.com/developer/article/1556485

上一篇: 又有一大老虎官員落馬，太痛心了！
下一篇: RabbitMQ（二）——Work Queues