三代测序 - Oxford Nanopore (ONT) 数据分析 - 数据质控和过滤

原创

三代测序说

修改于 2024-03-05 10:23:42

3.8K8

修改于 2024-03-05 10:23:42

对于二代测序，可以使用fastQC软件对原始下机数据质量进行全面的统计及绘图，多样本时可进一步结合使用multiQC。但是由于二代测序（如 illumina 平台）独特的特性，例如读长长度一致，可能包含测序接头（Adapter) 和PCR 引起的片段重复（Duplication）等。因此，同是fastq格式文件，fastQC并不适合用来处理三代测序数据，包括Oxford Nanopore（ONT）测序平台原始下机的fastq文件。

一、ONT下机数据质量

当前ONT测序质量虽然有很大的改善，但准确性依然不及二代测序，例如illumina或者BGIseq等。2018-2019年主流芯片R9.4 准确率对于2D reads为94%，1Dreads仅为86% ，如下图Fig.2b所示(1)。

现在最新的R10.4.1及配套试剂（2），使用Dorado v0.2.5进行高准确度碱基转换（High accuracy ，HAC）原始数据平均碱基质量可达Q20（99.0%），进行超高准确度碱基转换（Super accuracy ，SUP）原始数据平均碱基质量可达Q23（99.5%），如果进行Duplex测序，原始数据平均碱基质量可达Q30（99.9%）。

当前Nanopore测序下级数据仍然以Q7作为数据过滤的一个标准。PacBio三代测序平台由于有CCS校准模式（也有参数可以设置），一般来说HiFi数据个人觉得不用做质控，但对于现在的ONT数据个人建议做个质控查看一下数据质量，如果有需要可以做一下reads过滤。相信随着ONT的测序准确度越来越高以及Duplex测序方式，今后也可能弱化质控过程。

二、演示ONT数据（fastq）

PRJNA809646：https://www.ebi.ac.uk/ena/browser/view/PRJNA809646

Whole-Genome Sequence analysis of Clinically Isolated Carbapenem Resistant Escherichia coli from Iran (IRANCOLI ONT, SRR21721846，SRR21721848）。

#下载 MinION sequencing: IRANCOLI ONT, SRR21721846
$ wget -c -t 0 ftp://ftp.sra.ebi.ac.uk/vol1/fastq/SRR217/046/SRR21721846/SRR21721846_1.fastq.gz    #数据大小521M

#下载 MinION sequencing: IRANCOLI ONT, SRR21721848
$ wget -c -t 0 ftp://ftp.sra.ebi.ac.uk/vol1/fastq/SRR217/048/SRR21721848/SRR21721848_1.fastq.gz    #数据大小543M

#fastq文件改名
mv SRR21721846_1.fastq.gz SRR21721846.fastq.gz 
mv SRR21721848_1.fastq.gz SRR21721848.fastq.gz

三、数据质控软件

直接上结论：NanoComp 用以对下机序列（reads）数据质控，Chopper用以数据过滤，剪切和去污染序列足以。 1. NanoComp 包含了 NanoStat，NanoPlot 和 NanoQC 的内容，而且多样本之间的比较和整合性比较好，所以 多样本 可以直接使用 NanoComp ，单个样本可以运行 Nanoplot 。 2.对下机原始序列进行过滤，剪切和污染序列的去除可以使用 Chopper 。 3.如果运行 NanoPlot ，则会产生 NanoStat.txt 文件，结果和单独运行 NanoStat 的结果一样。 直接运行Nanoplot就行 。

Wouter De Coster， Rosa Rademakers 实验室的博后，编写了一系列软件对长度长测序原始下机数据进行统计质控和对污染序列进行去除的软件，也是现在常用的ONT质控软件。

NanoComp: comparing multiple runs 比较多组长度长测序和比对数据。
NanoStat: statistic summary report of reads or alignments 长度长测序和比对数据的统计报告。
NanoFilt: filtering and trimming of reads 过滤和修剪长度长序列。
NanoLyse: removing contaminant reads (e.g. lambda control DNA) from fastq 从fastq文件中去除污染序列（e.g. lambda对照/参照DNA序列）。
NanoPlot: Plotting tool for long read sequencing data and alignments. 对长度长测序和比对数据的统计报告进行可视化。

以上所有软件都存放在Nanopack的github里面，可以通过一行安装命令进行安装:

#一键安装ONT系列质控软件
$ pip install nanopack

# NanoComp        NanoFilt        NanoLyse        NanoPlot        NanoStat

Nanopack: https://github.com/wdecoster/nanopack

备注：

1. Nanostat 已经被运行速度更快的 Cramino 替代了：

CRAMINO: A tool for quick quality assessment of cram and bam files, intended for long read sequencing 对长度长序列的cram和bam文件进行质量评估。

2. NanoFlit 和 NanoLyse 被运行速度更快的 chopper 整合替代了：

Chopper: A rust implementation combining NanoLyse and NanoFilt into one faster tool for filtering, trimming, and removing contaminants。

3.nanoQC需要单独安装。

下面我们介绍每一个软件的具体安装和使用方法：

1. NanoPlot

主页网址：：https://github.com/wdecoster/NanoPlot

软件安装

#pip install
$ pip install NanoPlot
$ pip install NanoPlot --upgrade

#conda install
$ conda install -c bioconda nanoplot

以前通过conda安装过，在新买的云服务器上，两种方法都报错。于是使用singularity的方法下载镜像, 如果大家也有安装问题可以尝试使用singularity镜像。

#从官方docker镜像地址下载nanoplot的singularity镜像
$ singularity pull nanoplot.sif docker://nanozoo/nanoplot

#使用
$ singularity exec nanoplot.sif NanoPlot

软件使用

#官方例子
$ NanoPlot --summary sequencing_summary.txt --loglength -o summary-plots-log-transformed  
$ NanoPlot -t 2 --fastq reads1.fastq.gz reads2.fastq.gz --maxlength 40000 --plots dot --legacy hex
$ NanoPlot -t 12 --color yellow --bam alignment1.bam alignment2.bam alignment3.bam --downsample 10000 -o bamplots_downsampled


#实际使用
$ NanoPlot -t 12 --N50 --dpi 300 \
--fastq samplename.fq.gz \
--title samplename \
-o samplename

$ singularity exec nanoplot.sif NanoPlot -t 12 --N50 --dpi 300 \
--fastq samplename.fq.gz \
--title samplename \
-o samplename \

# -N50 在片段长度N50处划线 
# --dpi 设置图片dpi
# -t 线程数
#--title 所有表的标题

演示数据

$ NanoPlot -t 24 --N50 --dpi 300 \
--fastq SRR21721846.fastq.gz \
--title SRR21721846 \
-o SRR21721846_NanoPlot

$ NanoPlot -t 24 --N50 --dpi 300 \
--fastq SRR21721848.fastq.gz  \
--title SRR21721848 \
-o SRR21721848_NanoPlot

结果展示

$ cd SRR21721846_NanoPlot/

#结果文件：
NanoPlot-report.html  主要报告文件
LengthvsQualityScatterPlot_dot.html  
LengthvsQualityScatterPlot_kde.html  
NanoStats.txt      
Yield_By_Length.html                     
Non_weightedHistogramReadlength.html  
WeightedHistogramReadlength.html
Non_weightedLogTransformed_HistogramReadlength.html
Non_weightedLogTransformed_HistogramReadlength.png 
WeightedLogTransformed_HistogramReadlength.html
NanoPlot_20231120_2259.log

每个结果目录中都有NanoPlot-report.html文件，用浏览器打开即可查看结果报告索引，报告中包括所有统计数据和交互图。

摘要统计Summary statistics

LengthvsQualityScatterPlot_dot

点状展示每条read的长度和质量的分布，两侧加柱状图进一步呈现长度和质量的分布情况。可以看出长度分布在2K以下，质量集中在Q13.

LengthvsQualityScatterPlot_kde

核密度估计(Kernel Density Estimation，kde）展示每条read的长度和质量的分布，两侧加柱状图进一步呈现长度和质量的分布情况。可以看出长度分布在2K以下，质量集中在Q13.

Nonweighted Histogram of Read Lengths 无加权长度分布直方图

Nonweighted LogTransformed Histogram of read lengths 无加权及取对数后的长度分布直方图

Non-weighted LogTransformed Histogram Read Length

测序数量较大，且长度分布极不均匀且偏短，只在底部看到一条线，或一个峰。此时就需要将数据进行以10为底对数转换再观察。

Weighted Histogram of Read Lengths 加权后的长度分布

Weighted LogTransformed Histogram of Read Lengths 加权及取对数后的长度分布

Weighted Histogram of Read Lenghts after Log Transformation

X轴为长度，Y轴是碱基数量，更好地看出不同长度上的碱基数量分布。

加权(weighted)的含义: They're weighted by the number of nucleotides per bin. As such longer reads count 'heavier'. The height of the bar, therefore, does not represent the number of reads, but the number of nucleotides in that bin.

Yield by length 长度产出图

X轴为长度，Y轴为产量的频率。一般为越长越少。

2. NanoQC

主页网址：https://github.com/wdecoster/nanoQC

软件安装

#pip 安装
$ pip install nanoQC

#conda安装
$ conda install -c bioconda nanoqc

软件使用

nanoQC [-h] [-v] [-o OUTDIR] fastq

positional arguments:
  fastq                 Reads data in fastq.gz format.

optional arguments:
  -h, --help            show this help message and exit
  -v, --version         Print version and exit.
  -o, --outdir OUTDIR   Specify directory in which output has to be created.
  -l, --minlen int      Minimum length of reads to be included in the plots
                        This also controls the length plotted in the graphs
                        from the beginning and end of reads (length plotted = minlen / 2)

-l 参数制定最短的reads长度限制

演示数据

$ nanoQC SRR21721846.fastq.gz \
-o SRR21721846_NanoQC

$ nanoQC  SRR21721848.fastq.gz \
-o  SRR21721848_NanoQC

结果展示

cd SRR21721846_NanoQC/

#结果文件：
nanoQC.html
NanoQC.log

输出结果包含log文件和html报告文件（主要看html）：

报告中包含read长度、碱基含量和碱基质量，个人觉得这个可以作为头尾碱基过滤的一个参考。

3. NanoStat

主页网址：https://github.com/wdecoster/nanostat

如果运行NanoPlot，则会产生NanoStat.txt文件，结果和单独运行NanoStat的结果一样。

软件安装

#pip 安装
$ pip install nanostat

#conda安装
$ conda install -c bioconda nanostat

软件使用

#官方例子
$ NanoStat --fastq reads.fastq.gz --outdir statreports
$ NanoStat --bam alignment.bam alignment2.bam

演示数据

$ NanoStat --fastq SRR21721846.fastq.gz -o SRR21721846_NanoStat -n SRR21721846.txt
$ NanoStat --fastq SRR21721848.fastq.gz -o SRR21721848_NanoStat -n SRR21721848.txt

# -o 要和 -n 一起连用，才能生成 SRR21721848_NanoStat 文件夹以及文件夹里的 SRR21721848.txt 文档（结果）

统计结果

$ cat SRR21721848.txt 
#SRR21721848 统计结果
General summary:         
Mean read length 平均read长度:                2,390.9
Mean read quality 平均碱基质量 :                  11.9
Median read length 长度中位数 :              1,691.0
Median read quality 碱基质量中位数:                12.8
Number of reads 读长总数:               247,524.0
Read length N50 累计半总长的片段大小(N50) :                 3,819.0
STDEV read length:               2,335.6
Total bases 总碱基数:               591,806,555.0
Number, percentage and megabases of reads above quality cutoffs
>Q5:    247524 (100.0%) 591.8Mb
>Q7:    247524 (100.0%) 591.8Mb
>Q10:   229057 (92.5%) 553.8Mb
>Q12:   157978 (63.8%) 404.2Mb
>Q15:   29039 (11.7%) 87.7Mb
Top 5 highest mean basecall quality scores and their read lengths
1:      20.2 (199)
2:      19.6 (1538)
3:      19.5 (1930)
4:      19.4 (3240)
5:      19.0 (2429)
Top 5 longest reads and their mean basecall quality score
1:      37983 (15.3)
2:      35820 (12.8)
3:      33832 (11.2)
4:      33532 (12.4)
5:      33501 (10.2)

4. Cramino

主页网址：https://github.com/wdecoster/cramino

软件安装

#conda安装
$ conda install -c bioconda cramino

软件使用

$ cramino --ubam samplename.bam -t 12

需要比对后或者未比对的bam文件，这里就不做实际演示了，个人认为跑一个NanoPlot就足以。

结果展示

# Histogram for read lengths:
     0-2000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
  2000-4000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
  4000-6000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
  6000-8000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
 8000-10000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
10000-12000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
12000-14000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
14000-16000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
16000-18000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
18000-20000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
20000-22000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
22000-24000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
24000-26000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
26000-28000 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
28000-30000 ∎∎∎∎∎∎∎∎∎∎∎∎
30000-32000 ∎∎∎∎∎∎∎∎∎
32000-34000 ∎∎∎∎∎∎
34000-36000 ∎∎∎∎
36000-38000 ∎∎
38000-40000 ∎
40000-42000 ∎
42000-44000 ∎
44000-46000 
46000-48000 
48000-50000 
50000-52000 
52000-54000 
54000-56000 
56000-58000 
58000-60000 
     60000+ 


# Histogram for Phred-scaled accuracies:
  Q0-1 
  Q1-2 
  Q2-3 
  Q3-4 
  Q4-5 
  Q5-6 ∎∎∎
  Q6-7 ∎∎∎∎∎∎∎∎∎∎∎∎
  Q7-8 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
  Q8-9 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
 Q9-10 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q10-11 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q11-12 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q12-13 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q13-14 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q14-15 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q15-16 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q16-17 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q17-18 ∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎∎
Q18-19 ∎∎∎∎
Q19-20 ∎
Q20-21 
Q21-22 
Q22-23 
Q23-24 
Q24-25 
Q25-26 
Q26-27 
Q27-28 
Q28-29 
Q29-30 
Q30-31 
Q31-32 
Q32-33 
Q33-34 
Q34-35 
Q35-36 
Q36-37 
Q37-38 
Q38-39 
Q39-40 
  Q40+

三、数据质控多样本整合

NanoComp

主页网址：https://github.com/wdecoster/nanocomp

Compare multiple runs of long read sequencing data and alignments. Creates violin plots or box plots of length, quality and percent identity and creates dynamic, overlaying read length histograms and a cumulative yield plot.

软件安装

#pip安装
$ pip install NanoComp

#This script is written for Python3.

软件使用

#官方示例
$ NanoComp --bam alignment1.bam alignment2.bam alignment3.bam --outdir compare-runs
$ NanoComp --fastq reads1.fastq.gz reads2.fastq.gz reads3.fastq.gz reads4.fastq.gz --names run1 run2 run3 run4

#实际样本
$ nohup NanoComp -t 24 -f pdf  \
--fastq SRR21721846.fastq.gz SRR21721848.fastq.gz \
--names SRR21721846 SRR21721848 \
-o NanoComp &

# -f 图片以pdf的格式输出
# -t 运行线程数

结果展示

每个结果目录中都有 NanoComp-report.html 文件，用浏览器打开即可查看结果报告索引，报告中包括所有统计数据和交互图。内容其实和NanoPlot产生的结果一样，只不过有多样本之间的比较。

摘要统计Summary statistics

总read数，总碱基数，读长N50数据：

Read读长，对数转换的Read读长，碱基质量分布：

长度分布直方图，标准化的长度分布直方图，加权后的长度分布：

取对数后的长度分布，标准化及取对数后的长度分布，加权及取对数后的长度分布：

四、数据过滤软件

chopper

官方网站: https://github.com/wdecoster/chopper

Rust implementation of NanoFilt+NanoLyse, both originally written in Python. This tool, intended for long read sequencing such as PacBio or ONT, filters and trims a fastq file. Filtering is done on average read quality and minimal or maximal read length, and applying a headcrop (start of read) and tailcrop (end of read) while printing the reads passing the filter. 基于平均Read质量，最小最长Read读长，头尾碱基修剪

安装方法
直接下载二进制文件，加入路径。
Conda安装。

#conda安装
$ conda install -c bioconda chopper

使用方法

#使用说明
FLAGS:
    -h, --help       Prints help information
    -V, --version    Prints version information

OPTIONS:
        --headcrop      Trim N nucleotides from the start of a read [default: 0] 序列前端剪切
        --maxlength     Sets a maximum read length [default: 2147483647] 最大长度
    -l, --minlength     Sets a minimum read length [default: 1] 最小长度
    -q, --quality       Sets a minimum Phred average quality score [default: 0] 最小平均碱基质量（Phred评分）
        --tailcrop      Trim N nucleotides from the end of a read [default: 0] 序列末端剪切
        --threads       Number of parallel threads to use [default: 4] cpu线程数
        --contam        Fasta file with reference to check potential contaminants against [default None] 去除污染序列

#使用示例
$ gunzip -c reads.fastq.gz | chopper -q 10 -l 500 | gzip > filtered_reads.fastq.gz

实际使用

根据质控结果和实际需求调整

$ gunzip -c SRR21721846.fastq.gz | chopper -q 15 -l 500 | gzip > SRR21721846_filtered_reads.fastq.gz

参考文献：

Wang, Y., Zhao, Y., Bollas, A., Wang, Y., & Au, K. F. (2021). Nanopore sequencing technology, bioinformatics and applications. Nature Biotechnology
Raw read accuracy
NanoPlot：三代纳米孔测序数据质量评估 - 刘永鑫
基因组组装---Nanopore数据评估（nanoqc和NanoPlot套件工具）

原创声明：本文系作者授权腾讯云开发者社区发表，未经许可，不得转载。

如有侵权，请联系 cloudcommunity@tencent.com 删除。

生物基因

原创声明：本文系作者授权腾讯云开发者社区发表，未经许可，不得转载。

如有侵权，请联系 cloudcommunity@tencent.com 删除。

登录后参与评论

0 条评论

热度