PaSS：用于 PacBio 测序的测序模拟器

张文敏, 1 本佳, 1魏朝春1 , 2

作者信息文章注释版权和许可信息免责声明

本文已被PMC 中的其他文章引用。

抽象的

背景

第三代测序平台，如 PacBio 测序，近年来发展迅速。PacBio 测序比第二代测序（或下一代测序，NGS）技术产生更长的读取，并且它具有独特的测序错误模式。有效的读取模拟器对于评估和促进用于 PacBio 测序数据分析的新生物信息学工具的开发至关重要。

结果

我们开发了一种新的Pa cBio S测序S模拟器 (PaSS)。它可以从目前可用的 PacBio 测序数据中学习序列模式。除了读取长度和错误率的分布之外，我们还包括了一个特定于上下文的测序错误模型。与现有的 PacBio 测序模拟器如 PBSIM、LongISLND 和 NPBSS 相比，PaSS 在许多方面表现更好。组装测试还表明，PaSS 模拟的读取与实验测序数据最相似。

结论

PaSS 是用于 PacBio 测序的有效序列模拟器。它将促进第三代测序数据新分析工具的评估和开发。

电子补充材料

本文的在线版本 (10.1186/s12859-019-2901-7) 包含补充材料，授权用户可以使用。

关键词：三代测序，二代测序，PacBio测序，测序模拟器，测序错误，序列模式

去：

背景

包括 PacBio 或 SMRT（单分子实时）测序和纳米孔测序在内的第三代测序技术正在引发基因组学研究的一场革命，因为它们为研究人员提供了前所未有的测序读取长度 [ 1 ]。Pacific BioSciences 开发的 SMRT 测序是应用最广泛的第三代测序技术之一[ 2 ]。

越来越多的生物信息学工具和算法，例如序列比对程序 BLASR [ 3 ] 和 GraphMap [ 4 ]、基因组组装程序 canu [ 5 ] 和 miniasm [ 6 ] 以及结构变异调用程序 PBHoney [ 7 ] 和 Sniffles [ 8 ]新兴的 SMRT 数据分析。此外，PacBio 测序已迅速开发出多个版本。使用针对特定版本的 PacBio 技术的测序模拟器模拟的读数对这些工具进行基准测试和评估至关重要。PacBio 数据的模拟有助于指导用户为自己的研究项目选择最合适的分析工具或方法 [ 9]。此外，生成计算机数据可以显着降低改进下游分析工具所需的成本和时间[ 10 ]。

PacBio reads 的特点与二代测序 reads 有很大的不同。它能够产生大约 10-15 kb 的读数，比现有的第二代测序方法长得多。长读数可用于跨越重复或复杂区域，例如大的结构变异，因为可以更精确地确定基因组中读数的映射位置。因此，长读长在分析重复区域和大的结构变异方面显示出优越性。相比之下，二代测序难，容易出现组装错误和空缺。然而，与第二代测序技术中的不到 1% 相比，每个碱基的错误率可能在 15% 左右，并且错误以插入缺失为主 [ 11]。然而，高错误率可以通过单分子循环测序或多通道测序来缓解。在测序过程中，可以使用环状模板对目标分子的正向和反向链进行多次测序[ 12 ]。称为聚合酶读数的输出序列可以分成多个称为子读数的读数。通过生成这些子读数的共识，可以提高最终输出的测序读数质量。尽管吞吐量仍然很低，但最新的测序仪 Sequel 可以生成比旧的测序仪 PacBio RSII 多 7 到 10 倍的序列。它可以产生 5-10Gb 碱基，每次运行大约 365 k~ 500 k 读数 [ 13 ]。此外，与 NGS 方法相比，PacBio 测序速度更快且没有 GC 偏差 [14 ]。

目前，有几种工具可以模拟 PacBio 读取，例如 PBSIM [ 15 ]、LongISLND [ 16 ] 和 NPBSS [ 17]。所有这些模拟器都可以估计读取长度分布，但只有 LongISLND 考虑了 PacBio 平台的多通道测序。PBSIM 可以使用基于模型或基于采样的方法来模拟读取。但是 PBSIM 的读取长度分布与当前数据并没有很好的匹配。LongISLND 采用称为extended-kmer 的序列上下文敏感方法来处理均聚物依赖性偏差，并且可以以多种文件格式输出。NPBSS 可以使用实际错误率和质量值 (QV) 之间的关系，但模拟时间较长。对于来自最新测序仪 Sequel 的序列，使用了固定质量值 (QV)，因此 QV 不代表实际错误率，而 PBSIM 和 NPBSS 模拟测序错误的方法是基于 QV。此外，LongISLND 无法处理 Sequel 数据的文件格式。最重要的是，这三个模拟器仅根据比对结果中的比对区域构建了它们的测序错误模型，因此缺少有关测序错误的一些信息，尤其是那些质量低的区域。

为了抓住测序技术的创新和改进现有方法，我们提出了一种新的 PacBio 序列模拟器 PaSS。PaSS 可以从真实的 PacBio 数据生成定制的测序模式模型，并使用定制或经验的测序模型为输入参考基因组生成子读数。最后，将 PaSS 和一些流行的现有模拟器进行了比较。结果和组装测试表明，PaSS 可以高保真地模拟 PacBio 读取。

去：

执行

通常，PaSS 可以使用先前为某些给定参考基因组构建的测序错误模型生成计算机读取。测序错误模型也可以从真实的 PacBio 测序数据中重新估计。所涉及的方法介绍如下。

从真实测序数据估计测序错误模型

为了更好地模拟 PacBio 测序，研究了 PacBio 测序仪的多通道测序。我们注意到目标读取长度和较长模板循环次数减少之间存在折衷。PaSS 可以从真实测序数据中学习读取序列模式。通过数的分布及其相应的读取长度分布记录在模型中以用于序列生成。

为了了解错误在读取中的分布情况，我们将 PacBio 读取与参考序列对齐。在我们尝试使用多种比对工具进行长读长后，我们采用 BLASR [ 3 ] 将测序读长与参考基因组或高质量的从头组装进行比对。分析真实测序数据的比对结果，提取测序误差模型，作为仿真阶段的输入。由于这些区域的高错误率，某些读取的头部和尾部区域可能无法与参考序列对齐。估计整个聚合酶读数的未对齐部分的比率可以获得更完整的模型（参见附加文件 1：图S1）。整个读取的平均测序质量因读取而异。PaSS 通过基于 kmer 的分析来学习不同错误类型（匹配/插入/删除/替换）的比率及其对应的序列上下文模式。表 S1 显示了真实测序数据的不同错误事件中的 64 k-mers (k = 3) 频率。每个事件都记录有其对应的 3 碱基序列作为参考，连续错误被视为一个事件。一些 k-mer 的错误率相对较高，尤其是前两个碱基相同的 k-mer。大肠杆菌和秀丽隐杆线虫数据集中的趋势（图 1和附加文件 1：图 S2）看起来很接近，这是合理的，因为这两个数据集来自同一个序列器 RSII P6-C4。误差大小分布也来自对齐结果。虽然我们观察到读数中相对位置的测序错误偏差，但我们没有在当前版本的 PaSS 中包含这种模式。

图。1

基于 64 k-mers (k = 3) 在大肠杆菌 K12、秀丽隐杆线虫和拟南芥的真实测序数据集中的测序错误率分布

PacBio 多通道测序读数的模拟

图 2说明了模拟过程。首先，根据 pass-number 的分布估计正反向循环的数量，读取长度由该 pass-number 的相应长度分布确定。然后，PaSS 从用户指定的参考基因组序列中随机抽取一个无错误读数。如果所选序列包含 Ns，则这些 Ns 将在 read 中随机替换为 ACGT。收集到的读数被视为序列模板，其子读数在正向和反向链之间交替。最后，引入错误以读取输出。标记为来自同一模板的reads根据聚合酶reads内部的相对位置分为推定未比对部分和比对部分。对于假定的未对齐部分，我们使用预设的高错误率。 1：表 S2），我们选择 0.4 作为默认值。对于对齐的区域，根据模型中记录的特定于上下文的 bin 随机绘制事件类型。当发生错误时，然后从模型中得出错误的长度。从真实的 PacBio 数据中，我们发现插入的核苷酸取决于序列上下文。因此，如果错误是插入，则插入的核苷酸也由上下文决定。如果错误是替换，则根据十二种替换类型的分布引入替换模式。

图 2

PaSS的系统图。测序概况（或测序错误模型）可以从真实测序数据及其与参考基因组的比对生成。可以根据参考基因组和测序图谱（或错误模型）模拟读数。对于每次读取，首先从参考基因组中选择一个序列片段，然后根据配置文件添加测序错误，其中包括读取两端的预设错误率、通过次数、读取长度分布、基于上下文误差模型和误差大小

真正的 PacBio 测序数据集

为了评估 PaSS 的性能，选择了三个用于大肠杆菌、秀丽隐杆线虫和拟南芥的真实 PacBio 测序数据集进行基准测试。附加文件1：表 S3 显示了这些数据集的简要统计数据，可以从附加文件1：表 S3 中列出的网站免费下载。为了对测序模拟器的性能进行全面评估，我们包含了来自两个不同平台 RSII 和 Sequel 的真实测序数据。大肠杆菌和秀丽隐杆线虫的测序数据来自 RSII 测序平台，而拟南芥测序数据来自最新的 Sequel 平台。

模拟方法比较

为了进行公平的比较，我们首先尝试从真实的测序数据中估计所有方法的测序模型，并使用为同一基因组生成的测序模型模拟读数。由于 NPBSS 程序只能模拟单个染色体，因此我们每次模拟一条染色体并将C.elegans和A. thaliana的读数混合。LongISLND 无法从 Sequel 数据生成配置文件，我们没有使用 LongISLND模拟拟南芥的读取。

去：

结果和讨论

实施了一种名为 PaSS 的用于 PacBio 测序的新测序模拟器。我们使用三个序列数据集将 PaSS 与三种现有的流行方法 PBSIM、LongISLND 和 NPBSS 进行了比较（有关更多详细信息，请参见方法）。

模拟结果和比较

比较模拟读数的长度分布和真实测序数据的长度分布，结果如图 3（A）所示。所有模拟器都获得与真实测序数据相似的长度分布。PBSIM 中定义的最大读取长度的默认值已过时，无法重新配置。

图 3

模拟数据与真实测序数据的比较。( a ) 读取长度分布的比较。子图 a、b 和 c 显示了分别为大肠杆菌 K12、秀丽隐杆线虫和拟南芥这三种生物收集的读取长度分布。( b ) 误差大小分布的比较。将来自不同方法的模拟读数与生物体大肠杆菌 K12的错误大小分布的真实测序数据进行比较。三个子图分别显示了插入、删除和替换的概率密度条形图

然后我们评估了错误基础的长度分布。图 3 (B) 显示了大肠杆菌错误碱基的长度分布（附加文件1：图 S3。对于秀丽隐杆线虫和拟南芥）。尽管大多数插入缺失的长度都是一个碱基，但大约有 15-20% 的插入和 7-10% 的缺失包含多个碱基。与插入缺失不同，绝大多数替换的长度是一个碱基。PBSIM 和 NPBSS 读取在这方面与实际测序数据更加不同，因为它们的模型中仅包含单碱基错误。

接下来，我们使用 Kolmogorov-Smirnov 检验（KS 检验）来确定真实测序数据和模拟数据的两个概率分布是否不同。该检验的原假设是两组数据来自同一分布。对读长分布和错误碱基分布进行KS检验。表 S4 和 S5 中的结果p值（参见附加文件1) 拒绝原假设，这表明所有模拟数据和真实测序数据之间的两个分布是可区分的。然而，在大多数情况下，真实测序数据与 PaSS 模拟数据之间的测试统计量 D 在几个模拟器中是最小的。检验统计量 D 是两个分布之间差异的最大值。因此，它表明来自 PaSS 的模拟数据分布与真实数据的分布之间的距离是最接近的。此外，它与附加文件1中所示的一致：图 S4 和 S5。

表 1（附加文件1：表S6-S7）显示了比对结果的统计数据，从中我们可以看到PaSS中插入、删除和替换的比对率和错误率比现有方法更符合真实测序数据。PBSIM、LongISLND和NPBSS模拟reads的99%以上的bases可以与reference对齐，而真实测序reads和PaSS reads的对齐率更一致，三个数据集的对齐率在89%到94%之间. 因为只有对齐的区域被分析并包含在估计的轮廓中，所以这三个模拟器忽略了未对齐的区域。如前所述，Sequel 测序数据中的质量值 (QV) 并不代表实际的错误率。所以，PBSIM 无法从 Sequel 平台的真实测序数据中获得合理的参数。如果我们确实使用 PBSIM 从真实的 Sequel 测序数据中重新估计测序错误模型，错误率可能低于 1%，远低于应有的水平。附加文件1：图S6显示了整个聚合酶读数的平均准确度（1-错误率）分布。测序读数的质量并不统一，PaSS 提供了比其他工具更真实的模拟结果。一般来说，PaSS 可以比其他模拟器更好地模拟 PacBio 数据，尤其是对于新的 Sequel 数据。

表格1

PBSIM、LongISLND、NPBSS、PaSS 模拟读数的统计数据和大肠杆菌 K12基因组的真实测序数据。读数与大肠杆菌K12基因组对齐

方法	对齐速率（读取）	对齐率（基础）	错误率	插入	删除	替代
真实数据	96.71%	91.74%	14.54%	9.42%	3.86%	1.27%
PSIM	99.99%	99.73%	12.27%	7.22%	3.16%	1.89%
长岛	99.90%	99.92%	11.07%	7.09%	2.77%	1.20%
NPBSS	100.00%	99.93%	11.48%	2.67%	6.05%	2.76%
经过	95.84%	92.53%	14.39%	8.97%	3.80%	1.62%

在单独的窗口中打开

我们进一步研究了读取中错误率与碱基相对位置之间的相关性。我们将每个聚合酶读数平均分成十个片段，然后计算每个片段的平均错误率。如附加文件1：图S7所示，真实测序数据的错误率在前一个或两个片段处迅速降低，然后在读取的末端片段处略有增加。我们尝试将reads分成10个均分区间和1个区间，发现10区间模型和1区间模型的模拟结果相似。因此，我们最终采用了一种区间模型。

速度比较

为了比较四个模拟器的速度，我们分别模拟了测序深度为 170 和 50的基因组大肠杆菌 K12和秀丽隐杆线虫的读数。我们在附加文件1中报告了所有模拟器的计算时间：表 S8。PaSS 可以与多线程并行运行。因此，PaSS 的不同运行时间列出了不同的线程数。PBSIM 是最快的工具，而 NPBSS 是最慢的工具。如果 PaSS 使用超过 4 个线程，则 PaSS 比 PBSIM 和 LongILSND 快。

使用装配结果评估模拟

我们对 PBSIM、NPBSS、LongISLND 和 PaSS 模拟的读取以及真实测序数据进行了组装。用 5、10、15、20、25、30、35 和 40 的测序深度模拟每个基因组。另外模拟了 45 测序深度的秀丽隐杆线虫，另外模拟了 45 和 50 测序深度的大肠杆菌。我们使用 canu 进行了从头组装，canu 是一种为嘈杂的长读序列器设计的组装器 [ 18 ]。夸斯特 [ 19] 用于将程序集与参考基因组进行比较，并根据某些特征评估模拟器。比较来自模拟数据集和真实测序数据的组装的重叠群数量、基因组分数、每 100 kb 的插入缺失、每 100kbp 的错配、N50 和插入缺失长度的重叠群。组装结果如图 4所示，附加文件1：图 S8 和 S9 是针对E. coli K12、C. elegans和A. thaliana分别。装配结果评估是模拟器之间的间接比较。与其他模拟器相比，PaSS 的结果显示出与真实测序数据结果更相似的模式。就组装中的重叠群数量和正在组装的基因组部分而言，曲线趋向于稳定在 25 倍的 大肠杆菌和 35 倍的有机体秀丽隐杆线虫. 这表明这些测序深度可以在实际实验中采用。从对齐的重叠群数量覆盖参考的比例来看，PaSS的组装结果比其他模拟器更接近真实测序数据。更重要的是，关于每 100kbp 对齐碱基的平均插入缺失和错配数的指标也表明，PaSS 模拟读取比其他模拟器更接近真实数据。尽管从高测序深度 (>30X) 模拟数据组装的 contigs 可以覆盖所有模拟器的大部分参考基因组，但较低测序深度 (5-30X) 的结果确实显示了真实数据和模拟数据之间的差距，并且来自 PaSS 的模拟读取比其他模拟器更类似于真实数据。

图 4

比较真实读取的组装结果和不同方法模拟的组装结果。六个子图显示了测序深度与结果 ( a ) 重叠群数量、( b ) 基因组分数、( c ) 每 100kbp 插入缺失数、( d ) 每 100kbp 错配数、( e ) N50 和 ( f )大肠杆菌 K12组装结果中 contigs 插入缺失的长度。X轴代表测序深度

去：

结论

在本文中，我们建议使用 PaSS 来模拟 PacBio 测序读取，以跟上最新的测序技术。我们将序列上下文结合到 PaSS 的测序模型中。与现有方法相比，在 PaSS 中也考虑了由于高错误率而无法对齐回参考的序列部分。根据我们的评估，与现有的模拟系统相比，PaSS 可以模拟 PacBio 测序读数更类似于真实的 PacBio 数据。总体而言，PaSS 是一种有效的序列模拟器，可以生成具有已知基本事实的基准数据集，从而有助于评估最新的生物信息学工具。此外，由于没有可用于数据分析的金标准，因此它可以作为研究人员的指导。

但是，PaSS 仍然可以在各个方面进行改进。首先，kmer的长度是有限的。其次，基于对齐结果估计误差模型的方法并不完善。对齐工具的算法和性能会影响估计的误差模型，并可能带来额外的偏差。第三，测序模拟器无法针对不同物种轻松定制或更新。

去：

可用性和要求

项目名称： Pass。

项目主页： http ://cgm.sjtu.edu.cn/PaSS

操作系统： Linux。

编程语言： Perl 和 C。

其他要求： Perl（5.10.1 或以上）、gcc（4.8.0 或以上）。

许可证： GNU GPL。

非学术人员使用的任何限制：无。

去：

附加文件

附加文件 1：（4.7M，docx）

PaSS 的补充材料（包括补充数字和表格）。(DOCX 4780 KB)

去：

致谢

感谢上海交通大学高性能计算中心 (HPCC) 的计算。

资金

这项工作得到了国家自然科学基金（61472246）、国家基础研究计划（2013CB956103）、国家高技术研发计划（863）（2014AA021502）和中国科学院跨院科研基金的资助。上海交通大学（YG2017ZD01）。资助者在研究设计、数据收集和分析、发表决定或手稿准备方面没有任何作用。

数据和材料的可用性

本文所有数据均可在PaSS: Sequencing simultor for PacBio Sequencing获取。

去：

作者的贡献

CW构思并设计了这项研究。WZ 和 BJ 实现了流水线，WZ 测试了流水线。WZ、BJ 和 CW 撰写了手稿。

去：

伦理批准和同意参与

不适用。

去：

利益争夺

作者声明他们没有竞争利益。