问小白 wenxiaobai
资讯
历史
科技
环境与自然
成长
游戏
财经
文学与艺术
美食
健康
家居
文化
情感
汽车
三农
军事
旅行
运动
教育
生活
星座命理

PacBio Iso-Seq三代测序数据分析流程详解

创作时间:
作者:
@小白创作中心

PacBio Iso-Seq三代测序数据分析流程详解

引用
1
来源
1.
https://cloud.tencent.com/developer/article/2382137

随着测序技术的不断发展，长读长测序技术在RNA测序领域展现出巨大潜力。PacBio Iso-Seq技术通过其独特的长读长优势，能够直接获得包含5’UTR、3’UTR、polyA尾的mRNA全长序列，为转录组研究提供了全新的视角。本文将详细介绍PacBio Iso-Seq的实验流程和数据分析方法，帮助读者深入了解这一前沿技术。

全长转录组测序概述

全长转录组（Full-length transcriptome）测序和分析是基于PacBio和Oxford Nanopore三代测序平台，利用其长读长的特性，建库测序时无需对RNA进行打断，可以直接获得包含5’UTR、3’UTR、polyA尾的mRNA全长序列及完整结构信息。这种技术能够准确分析有参考基因组物种的可变剪接及融合基因等结构信息，克服无参考基因组物种转录本拼接组装较短、信息不完整的难题。

最早PacBio关于全长转录组的产品命名为Iso-Seq，全称叫做Isoform-Sequencing，是对自家开发的转录本测序技术的规范化命名。现在使用最新的试剂盒SMRTbell prep kit 3.0进行测序文库的构建。2023年10月初推出的Kinnex全长RNA建库试剂盒，能将5个转录本串联成一条测序read，充分利用其读长长的优势，提高测序通量，配合Revio系统一起使用，使得对全长转录本进行定量变得更为实际。

Iso-Seq方法可对整个cDNA分子（长达10kb或更长）进行测序，无需进行生物信息学转录本组装，因此可以对批量（bulk）和单细胞转录本组中的新基因和异构体进行表征，并进一步：

  • 鉴定可变剪接(AS)事件，包括可变起始位点、终止位点、内含子保留和外显子跳跃事件。
  • 通过开放阅读框(ORF)预测新型同源异构体的功能影响。
  • 检测差异表达的同源异构体和同源异构体的转换事件。
  • 发现肿瘤样本中的基因融合事件。
  • 识别等位基因同源异构体。

PacBio Iso-Seq的实验流程

通过PacBio SMRTbell prep kit 3.0构建Iso-Seq测序文库，适用于PacBio Sequel II和Revio仪器型号。

  1. 通过带有Oligo-dT的引物富集含有polyA尾的mRNA。
  2. 使用Iso-Seq RT enzyme对mRNA进行反转录。
  3. 加入模板转换oligo (Template Switch Oligo, TSO)
  4. 通过PCR扩增富集合成的cDNA，此步骤可加入barcode，用于混样。
  5. 对全长cDNA进行损伤修复、末端修复、末端加A尾。
  6. 连接SMRT哑铃型测序接头，最后结合测序引物，绑定DNA聚合酶，形成完整的SMRT-bell测序文库。

Iso-Seq全长转录组分析流程

Iso-seq基础概念

  • ROI：reads of insert
    ROI，全称reads of insert，可以理解为插入片段。对于全长转录本而言，其ROI reads中包含5’primer和3’primer；而且会出现polyA为结构；（polyA针对mRNA和部分lncRNA）。

  • Artifacts，文库构建过程中可能产生的非正常转录本
    可以理解为，共有两种来源：

  • Artificial Concatemer
    这种序列是由于文库制备阶段，adapter序列错误地将两条转录本的序列链接构成了一个环状分子，这个和adapter浓度有关，通常这种reads产生的比例很少，小于0.5%，在后续的分析中，这部分reads需要去除。

  • PCR Chimera
    在PCR反应中，由于不完全延伸的产物作为了下次扩增反应的引物，导致出现嵌合体序列，直观上看，就是PCR产物来源于两条或者多条reads。PCR产生的嵌合体序列，在PCR反应体系中，这种序列是不可避免的，大约有3%的比例，在后续的分析过程中，可以借助软件去除这部分reads。

  • FL Reads
    FL，Full-length reads，全长转录本。从raw data到ROI，在从ROI去除artifacts reads之后，我们就得到了用于后续分析的clean reads。clean reads就已经是转录本的序列了，我们首先看一下clean reads当中，哪些是全长转录本；哪些不是全长转录本。

Iso-Seq软件

Iso-Seq是PacBio官方开发的一款用PacBio subreads或HiFi数据进行全长转录组分析的软件，最终输出高质量转录本的全长序列。截止到2023年6月7号，最新版本为4.0.0。

Github主页：https://github.com/PacificBiosciences/IsoSeq

软件安装：isoseq和lima

# 使用conda安装isoseq，v4.0.0.
$ conda install -c bioconda isoseq
# lima,用于拆解barcode.
$ conda install -c bioconda lima

Iso-Seq分析流程

整个Iso-Seq的分析流程如下：

  1. 原始下机数据
  • subreads.bam
  • 通过CCS软件获得HiFi Reads,
  • hifi.reads.bam。

  1. cDNA两端引物的去除，拆解barcode，调整转录本5' - 3' (polyA)方向。

  2. 3' polyA尾和嵌合体（concatemer）序列的去除。

  3. 转录本聚类。

  4. Consensus的转录本序列以.fasta格式输出。

单个样本官方示例数据演示

（1）示例数据下载

# Download toy S-read dataset
# This is a toy dataset consisting of ~80k segmented reads (S-reads) from a Kinnex full-length RNA library
$ wget https://downloads.pacbcloud.com/public/dataset/IsoSeq_sandbox/human_80k_Sreads.segmented.bam
# Download the Iso-Seq v2 cDNA primers (from Iso-Seq express 2.0 kit)
$ wget https://downloads.pacbcloud.com/public/dataset/Kinnex-full-length-RNA/REF-primers/IsoSeq_v2_primers_12.fasta

（2）demultiplex，使用lima拆解barcode

$ lima --version
lima 2.9.0
$ lima human_80k_Sreads.segmented.bam IsoSeq_v2_primers_12.fasta human_80k.bam --isoseq --peek-guess

对于iso-seq数据，使用--isoseq加--peek-guess参数来降低假阳性率。可以使用--biosample-csv input.csv添加样本名称,bio sample name。

（3）refine，使用isoseq refine去除poly(A)和嵌合体（concatemer）序列

  • 输入文件为：
  • <primer--pair>.fl.bam
  • 输出文件主要为：
  • <movie>.flnc.bam
$ isoseq refine human_80K.IsoSeqX_bc10_5p--IsoSeqX_3p.bam IsoSeq_v2_primers_12.fasta human_80K.flnc.bam --require-polya
$ ls human_80K.flnc.*
human_80K.flnc.bam
human_80K.flnc.bam.pbi
human_80K.flnc.consensusreadset.xml
human_80K.flnc.filter_summary.report.json
human_80K.flnc.report.csv

--require-polya：只有有polyA尾序列才会被认作full length，并且去除polyA序列。
-j，--num-threads：CPU线程数。
--min-polya-length：最小polyA尾长度，默认值为20。

3' polyA尾和嵌合体（concatemer）序列的去除后得到 Full-Length Non-Concatemer (FLNC) reads。

（4）cluster，转录本聚类

  • 输入文件：
  • <movie>.flnc.bam
  • flnc.fofn。
  • 输出文件：
  • <prefix>.bam
$ isoseq cluster2 human_80K.flnc.bam human_80K.transcripts.bam
$ ls human_80K.transcripts.*
human_80K.transcripts.bam
human_80K.transcripts.bam.pbi
human_80K.transcripts.cluster_report.csv

输出转录本的同源异构体（isoforms）至少有两个或两个以上的FLNC（full-length non-concatemer）序列支持。如果想包含singletons，可以加入--singletons。

关于满足什么条件isoseq cluster算法会将两条序列聚类为同一个转录本序列：

  • 5'端差异小于100bp以内。
  • 3'端差异小于30bp以内。
  • 小于10bp以内的gaps，gaps的数量没有上限。

混样数据官方示例演示

（1）示例数据下载

# This is a 12-plex regular Iso-Seq (non-Kinex) run on Sequel II system consisting of ~3 million HiFi reads.
# Download HiFi reads from a non-Kinnex (regular Iso-Seq) BAM file
$ wget https://downloads.pacbcloud.com/public/dataset/Kinnex-full-length-RNA/DATA-SQ2-UHRR-Monomer/1-CCS/m64307e_230628_025302.hifi_reads.bam
$ wget https://downloads.pacbcloud.com/public/dataset/Kinnex-full-length-RNA/DATA-SQ2-UHRR-Monomer/1-CCS/m64307e_230628_025302.hifi_reads.bam.pbi
# Download the Iso-Seq v2 cDNA primers (from Iso-Seq express 2.0 kit)
$ wget https://downloads.pacbcloud.com/public/dataset/Kinnex-full-length-RNA/REF-primers/IsoSeq_v2_primers_12.fasta

（2）demultiplex，使用lima拆解barcode

$ lima --version
lima 2.9.0
# Demux and primer removal
$ lima --isoseq --peek-guess m64307e_230628_025302.hifi_reads.bam IsoSeq_v2_primers_12.fasta UHRR.bam

每个barcode对输出一个.bam文件，一共12个.bam文件对应12个样本。

（3）合并拆分样本的文件名

# Combine inputs
$ ls UHRR.IsoSeqX*bam > all.fofn
$ cat all.fofn
UHRR.IsoSeqX_bc01_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc02_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc03_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc04_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc05_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc06_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc07_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc08_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc09_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc10_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc11_5p--IsoSeqX_3p.bam
UHRR.IsoSeqX_bc12_5p--IsoSeqX_3p.bam

fofn: "files-of-file-names"的缩写。

（4）refine，使用isoseq refine去除poly(A)和嵌合体（concatemer）序列

# Remove poly(A) tails and concatemer
$ isoseq refine all.fofn IsoSeq_v2_primers_12.fasta UHRR.flnc.bam --require-polya
$ ls UHRR.flnc.*
UHRR.flnc.bam
UHRR.flnc.bam.pbi
UHRR.flnc.consensusreadset.xml
UHRR.flnc.filter_summary.report.json
UHRR_80K.flnc.report.csv

（5）cluster，转录本聚类

$ isoseq cluster2 UHRR.flnc.bam UHRR.transcripts.bam

（6）Polish （可选择）

$ isoseq cluster flnc.fofn clustered.bam --verbose --use-qvs

这里使用isoseq cluster，而不是isoseq cluster2, cluster相比于cluster2比较耗时。Polish会略微提升数据质量，不是必须步骤。运行完成以后获得以下文件：

  • <prefix>.bam
  • <prefix>.hq.fasta.gz with predicted accuracy ≥ 0.99
  • <prefix>.lq.fasta.gz with predicted accuracy < 0.99
  • <prefix>.bam.pbi
  • <prefix>.transcriptset.xml

综上所展示，iso-seq分析流程及每一步输出文件总览，如图16。

参考文献

  1. Iso-seq 必备基础-blog.csdn
  2. pacbio 三代全长转录组数据分析流程
  3. PacBio Iso-Seq Workshop Online
热门推荐
© 2023 北京元石科技有限公司 ◎ 京公网安备 11010802042949号