CN112735524A - 一种基于神经网络的真实纳米孔测序信号滤波方法及装置 - Google Patents

一种基于神经网络的真实纳米孔测序信号滤波方法及装置 Download PDF

Info

Publication number
CN112735524A
CN112735524A CN202011583559.3A CN202011583559A CN112735524A CN 112735524 A CN112735524 A CN 112735524A CN 202011583559 A CN202011583559 A CN 202011583559A CN 112735524 A CN112735524 A CN 112735524A
Authority
CN
China
Prior art keywords
sequencing
signal sequence
trained
sequence
sequencing signal
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Withdrawn
Application number
CN202011583559.3A
Other languages
English (en)
Inventor
陈为刚
张鹏
韩昌彩
赵毅强
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hefei Institute Of Innovation And Development Tianjin University
Original Assignee
Hefei Institute Of Innovation And Development Tianjin University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hefei Institute Of Innovation And Development Tianjin University filed Critical Hefei Institute Of Innovation And Development Tianjin University
Priority to CN202011583559.3A priority Critical patent/CN112735524A/zh
Publication of CN112735524A publication Critical patent/CN112735524A/zh
Withdrawn legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/04—Architecture, e.g. interconnection topology
    • G06N3/045—Combinations of networks
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/02—Preprocessing
    • G06F2218/04—Denoising
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/22—Source localisation; Inverse modelling

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Biophysics (AREA)
  • General Physics & Mathematics (AREA)
  • Chemical & Material Sciences (AREA)
  • Evolutionary Computation (AREA)
  • Computational Linguistics (AREA)
  • Molecular Biology (AREA)
  • Computing Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Biomedical Technology (AREA)
  • Mathematical Physics (AREA)
  • Software Systems (AREA)
  • Artificial Intelligence (AREA)
  • Data Mining & Analysis (AREA)
  • Analytical Chemistry (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Biotechnology (AREA)
  • Evolutionary Biology (AREA)
  • Medical Informatics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
  • Apparatus Associated With Microorganisms And Enzymes (AREA)

Abstract

本发明公开了一种基于神经网络的真实纳米孔测序信号滤波方法及装置，所述方法包括以下步骤：在纳米孔测序K‑mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；构建基于双向门控循环单元神经网络的真实测序信号处理模型；构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理；本发明的优点在于：能够准确地滤掉待测真实测序信号序列中与实际测序信号序列无关的高频分量并保留有用的高频分量。

Description

一种基于神经网络的真实纳米孔测序信号滤波方法及装置
技术领域
本发明涉及纳米孔测序信号处理领域，更具体涉及一种基于神经网络的真实纳米孔测序信号滤波方法及装置。
背景技术
新一代测序技术的提出使研究人员能够以高通量的方式对DNA和RNA进行测序，这促进了在基因组学、转录组学和表观基因组学等领域的众多突破。目前市场上最受欢迎的新一代测序技术主要包括Illumina、PacBio和Nanopore等测序平台。与其他测序技术不同的是，纳米孔测序技术主要是通过嵌入在分隔两个电解液室的薄膜上的纳米孔来直接检测碱基序列。通过在薄膜上施加电压差，将单链核苷酸序列通过纳米孔穿过薄膜，孔中的碱基会影响孔的电阻，因此可以由随时间变化的电流信号检测到通过孔的碱基序列。与目前常用的短读段测序技术(例如Illumina MiSeq测序平台)相比，纳米孔测序技术具有多个优势。首先，纳米孔测序技术可以从单个核苷酸分子中实时产生测序读段，结合快速的文库制备工作，有效地缩短了从测序生物样本预处理到测序数据分析之间的时间。此外，纳米孔测序技术还可以直接用于RNA测序，而无需提前进行逆转录或扩增工作。任意长度的DNA或RNA分子都可以使用纳米孔测序技术来进行测序，长读段是非常有价值的，它们简化了基因组组装和结构变异检测等研究工作。与此同时，纳米孔测序平台(如MinION测序仪)要比当前的短读段测序平台更加轻便，这使得它能够在传统实验室环境之外进行测序工作。
尽管纳米孔测序技术的优势众多，但复杂的测序环境会导致低信噪比的测序信号序列，这对于进一步的测序数据分析提出了挑战。在纳米孔测序过程中，待测序碱基序列移动通过纳米孔时输出的电流信号测量值存储为16位整数值。对于当前的MinION孔化学，单个DNA链以平均450bp/s的速度穿过纳米孔，而测序电流信号的采样频率为4kHz，这意味着平均每个k-mer有9个离散的测量信号值。尽管由于运动蛋白易位速度的波动，这些值可能是不同的。为了将输出的原始测序电流信号序列转换成碱基序列，需要复杂的碱基识别软件来实现。在纳米孔测序过程中，有以下几个因素会导致低信噪比的原始电流信号序列：通过纳米孔的DNA或RNA的四种碱基结构较为相似，导致不同碱基产生的原始电信号差异较小；原始电流信号主要受到同时占据纳米孔的5或6个碱基的影响，所以一个电信号测量值对应45或46个可能的k-mers；牵引DNA或RNA序列移动的运动蛋白的随机性导致碱基序列的移速不均匀；均聚物通过纳米孔时的电信号没有发生变化，从而导致多个相同碱基的长链被错误解释。上述几个因素使得纳米孔测序技术产生的测序读段精度较低，从而限制了纳米孔测序技术在实际中的应用。
随着纳米孔测序技术的快速发展，针对纳米孔测序数据的分析方法和工具也迅速涌现。例如，用于纳米孔测序技术的基因序列比对工具Graphmap，Minimap2和MashMap2，用于组装纳米孔测序长读段的基因组组装工具Canu和Racon，用于纳米孔测序信号的可视化工具BulkVis和SquiggleKit。可以预见，在不久的将来研究人员将开发出更多的用于纳米孔测序的数据分析方法和工具。研究人员通常使用实验测序数据或仿真测序数据对这些新算法和新工具进行性能测试，与实验测序数据相比，仿真测序数据可以大大节约研发成本，降低数据分析难度，提高研发效率，因此对纳米孔测序数据进行更准确的仿真有助于为纳米孔测序技术开发出更多性能更优的方法和工具。
纳米孔测序数据的仿真算法可以分为纳米孔测序读段仿真和纳米孔测序信号仿真两种。ReadSim、SiLiCO和NanoSim三个纳米孔测序读段仿真软件利用输入的核苷酸序列和配置文件生成仿真测序读段，其中配置文件包含一组预先设置的参数，例如插入率、删节率、替代率、读段长度和质量分数等参数。这三个纳米孔测序读段仿真软件的不同之处在于ReadSim使用固定的配置文件，SiLiCO使用用户提供的配置文件，NanoSim使用用户提供的实际测序数据来学习将在仿真阶段使用的配置文件。尽管纳米孔测序读段仿真软件可以生成高质量的仿真读段，但是纳米孔测序电信号才是纳米孔测序技术的本质，这些纳米孔测序读段仿真软件无法生成纳米孔测序仿真信号，从而限制了这些测序读段仿真软件的应用。
DeepSimulator是第一个完整地模拟纳米孔测序技术整个流程的仿真软件，它可以同时生成纳米孔测序仿真信号和仿真读段。纳米孔测序技术整个流程主要包括三个阶段：首先通过样本制备过程产生用于测序实验的核苷酸样品；然后使用纳米孔测序设备(如MinION)测量核苷酸序列通过纳米孔时输出的电流信号序列，采集到的测序电流信号序列通常存储在fast5文件中；最后使用碱基识别软件将测序电流信号序列转换为碱基序列。相应地，DeepSimulator的主要工作框架包括序列生成模块，信号生成模块和碱基识别模块。首先，序列生成模块随机选择输入参考基因组序列上的起始位置，生成满足实际测序读段长度分布的较短核苷酸序列。然后，信号生成模块根据已知的纳米孔测序6-mer孔模型生成序列生成模块输出的核苷酸序列所对应的纳米孔测序仿真信号序列。最后，碱基识别模块使用碱基识别软件(Albacore、Guppy)将仿真测序信号序列转换为仿真测序读段。在DeepSimulator主要工作框架中的三个模块中，信号生成模块是DeepSimulator软件的核心模块，其首先根据已知的纳米孔测序6-mer孔模型和实际测序信号序列的重复次数分布生成输入核苷酸序列所对应的真实测序信号序列，然后使用低通滤波器滤掉嵌入在真实测序信号序列中与真实测序信号序列无关的高频分量，最后在滤波信号序列上添加高斯噪声输出最终的仿真测序信号序列。
为了生成高质量的纳米孔仿真测序信号序列，DeepSimulator信号生成模块中的一个关键步骤是使用低通滤波器对真实测序信号序列进行滤波处理。真实测序信号序列是由一系列方波组成，其频谱是无限正弦波的组合。为了更加准确地对纳米孔测序信号序列进行仿真，必须对嵌入在真实测序信号序列方波中的高频分量进行滤波处理。DeepSimulator通过将windowed-sinc函数和真实测序信号序列进行卷积来实现对真实测序信号序列的低通滤波处理。考虑到单链核苷酸序列通过纳米孔的速度约为450bp/s，因此低通滤波器的截止频率应大于450Hz。当低通滤波器的截止频率设为950Hz时，可以生成与实际测序信号序列最相似的仿真测序信号序列。低通滤波器是一种可以通过频率低于所选截止频率的信号并衰减频率高于截止频率的信号的滤波器。由于低通滤波器会衰减真实测序信号序列中所有高于截止频率的高频分量，无法保留其中与实际测序信号序列相关的高频分量，这对于某些输入核苷酸序列来说可能会造成仿真测序信号序列与真实测序信号序列之间存在着较大的差异，这给关心仿真测序信号序列输出的用户带来不便。
发明内容
本发明所要解决的技术问题在于传统低通滤波器无法准确地滤掉真实测序信号序列中与实际测序信号序列无关的高频分量的问题。
本发明通过以下技术手段实现解决上述技术问题的：一种基于神经网络的真实纳米孔测序信号滤波方法，所述方法包括以下步骤：
步骤一：在纳米孔测序K-mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；
步骤二：构建基于双向门控循环单元神经网络的真实测序信号处理模型，该信号处理模型包括三层双向门控循环单元神经网络和一层全连接层，其中信号处理模型的输入为输入到第一层双向门控循环单元神经网络的经中值归一化的真实测序信号序列，信号处理模型的输出为全连接层输出的与输入的待测真实测序信号序列等长的滤波信号序列；
步骤三：获取纳米孔测序平台输出的待训练实际测序信号序列并根据待训练实际测序信号序列计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列，将待训练实际测序信号序列和与其对应的待训练真实测序信号序列作为信号处理模型参数训练所需的监督训练数据，构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；
步骤四：将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理。
本发明基于已知的纳米孔测序K-mer孔模型将输入核苷酸序列转换为预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；然后建立真实测序信号序列处理模型，使用监督训练数据完成信号处理模型参数的训练；最后使用完成训练的信号处理模型对待测真实测序信号序列进行滤波处理，利用神经网络学习到实际测序信号序列的时频特性，进而能够准确地滤掉待测真实测序信号序列中与实际测序信号序列无关的高频分量并保留有用的高频分量。
进一步地，所述步骤一包括：
步骤101：在纳米孔测序K-mer孔模型中输入一条具有T个碱基的核苷酸序列X＝x1,x2,...,xT，xT表示第T个碱基，从输入核苷酸序列的第一个碱基开始，每次移动一个碱基取出一个K-mer，直到移动到倒数第K个碱基位置结束，T个碱基共得到T-K+1个K-mer；
步骤102：纳米孔测序K-mer孔模型包含每个K-mer所对应的预期电流信号值，对照纳米孔测序K-mer孔模型依次查找T-K+1个K-mer所对应的预期电流信号值，生成输入核苷酸序列X对应的预期测序信号序列Y＝y1,y2,...,yT-5，其中yi表示从X中的位置i开始的K-mer所对应的预期电流信号值；
步骤103：根据实际测序信号序列中的信号值重复次数分布，将预期测序信号序列中的每个信号值重复多次，生成与实际测序信号序列之间具有相似长度分布的待测真实测序信号序列。
进一步地，所述步骤二包括：
步骤201：将中值归一化后的待测真实测序信号序列I＝i1,i2,...,iT作为信号处理模型中第一层双向门控循环单元神经网络的输入向量，其中双向门控循环单元神经网络的基本单元由一个前向传播的门控循环单元和一个后向传播的门控循环单元组成，首先计算t时刻的前向输出向量
Figure BDA0002864906830000051
计算公式为：
Figure BDA0002864906830000052
Figure BDA0002864906830000053
Figure BDA0002864906830000054
Figure BDA0002864906830000055
其中，it为t时刻的输入信号值，
Figure BDA0002864906830000056
为t-1时刻的前向输出向量，σ为sigmoid函数且σ(z)＝1/(1+e-z)，
Figure BDA0002864906830000057
表示第一中间变量，
Figure BDA0002864906830000058
表示第二中间变量，
Figure BDA0002864906830000059
表示第三中间变量，
Figure BDA00028649068300000530
表示两个等维向量的对应元素相乘，
Figure BDA00028649068300000510
为第一权重矩阵，
Figure BDA00028649068300000511
为第二权重矩阵，Wf为第三权重矩阵，
Figure BDA00028649068300000512
为第一偏置向量，
Figure BDA00028649068300000513
为第二偏置向量，然后计算t时刻的后向输出向量
Figure BDA00028649068300000514
计算公式为：
Figure BDA00028649068300000515
Figure BDA00028649068300000516
Figure BDA00028649068300000517
Figure BDA00028649068300000518
其中，
Figure BDA00028649068300000519
为t+1时刻的后向输出向量，
Figure BDA00028649068300000520
表示第四中间变量，
Figure BDA00028649068300000521
表示第五中间变量，
Figure BDA00028649068300000522
表示第六中间变量，
Figure BDA00028649068300000523
为第四权重矩阵，
Figure BDA00028649068300000524
为第五权重矩阵，Wr为第六权重矩阵，
Figure BDA00028649068300000525
为第三偏置向量，
Figure BDA00028649068300000526
为第四偏置向量，最后t时刻的输出向量ht是由t时刻的前向输出向量
Figure BDA00028649068300000527
和t时刻的后向输出向量
Figure BDA00028649068300000528
连接得到，即
Figure BDA00028649068300000529
其中||表示向量的连接符号；
步骤202：将第一层双向门控循环单元神经网络的输出向量作为第二层双向门控循环单元神经网络的输入向量，第二层双向门控循环单元神经网络的输出向量作为第三层双向门控循环单元神经网络的输入向量，其中第二层双向门控循环单元神经网络、第三层双向门控循环单元神经网络的计算过程与第一层双向门控循环单元神经网络的计算过程相同，不同层之间的权重矩阵与偏置向量参数不同；
步骤203：将三层双向门控循环单元神经网络中最后一层的输出向量作为全连接层的输入向量进行处理，全连接层的输出向量为与输入的待测真实测序信号序列等长的滤波信号序列O＝o1,o2,...,oT。
进一步地，所述步骤三包括：
步骤301：信号处理模型的监督训练数据包含待训练实际测序信号序列和与其对应的待训练真实测序信号序列，首先读取纳米孔测序平台输出的fast5文件得到待训练实际测序信号序列，然后使用连续小波动态时间规整算法计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列；
步骤302：构建信号处理模型的损失函数为
Figure BDA0002864906830000061
其中，cosh表示双曲函数，其中，oi表示信号处理模型的第i个输出信号，ri表示实际测序信号序列中第i个信号；
步骤303：将每层双向门控循环单元神经网络中权重矩阵以及偏置向量初始化为区间
Figure BDA0002864906830000062
内的均匀分布，其中ni表示每层双向门控循环单元神经网络的输入向量维度，no表示每层双向门控循环单元神经网络的输出向量维度；
步骤304：使用Adam自适应优化器迭代计算得到最小化损失函数的模型参数，完成信号处理模型的参数训练过程。
更进一步地，所述步骤301包括：
步骤3011：读取纳米孔测序平台输出的fast5测序文件得到待训练实际测序信号序列，其中fast5文件为HDF5文件格式；
步骤3012：对待训练实际测序信号序列进行碱基识别得到测序读段，采用基因序列比对算法将测序读段比对到参考基因组序列，根据比对结果得到与待训练实际测序信号序列对应的参考基因组序列片段，然后使用参考基因组序列片段和纳米孔测序K-mer孔模计算出待训练预期测序信号序列；
步骤3013：使用连续小波动态时间规整算法完成待训练实际测序信号序列与待训练预期测序信号序列之间的点到点映射，根据映射结果，求得待训练实际测序信号序列中每个信号值所对应的预期信号值，进而计算出待训练实际测序信号序列所对应的待训练真实测序信号序列。
本发明还提供一种基于神经网络的真实纳米孔测序信号滤波装置，所述装置包括以下步骤：
待测真实测序信号序列生成模块，用于在纳米孔测序K-mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；
信号处理模型构建模块，用于构建基于双向门控循环单元神经网络的真实测序信号处理模型，该信号处理模型包括三层双向门控循环单元神经网络和一层全连接层，其中信号处理模型的输入为输入到第一层双向门控循环单元神经网络的经中值归一化的真实测序信号序列，信号处理模型的输出为全连接层输出的与输入的待测真实测序信号序列等长的滤波信号序列；
信号处理模型训练模块，用于获取纳米孔测序平台输出的待训练实际测序信号序列并根据待训练实际测序信号序列计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列，将待训练实际测序信号序列和与其对应的待训练真实测序信号序列作为信号处理模型参数训练所需的监督训练数据，构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；
滤波处理模块，用于将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理。
进一步地，所述待测真实测序信号序列生成模块还用于：
步骤101：在纳米孔测序K-mer孔模型中输入一条具有T个碱基的核苷酸序列X＝x1,x2,...,xT，xT表示第T个碱基，从输入核苷酸序列的第一个碱基开始，每次移动一个碱基取出一个K-mer，直到移动到倒数第K个碱基位置结束，T个碱基共得到T-K+1个K-mer；
步骤102：纳米孔测序K-mer孔模型包含每个K-mer所对应的预期电流信号值，对照纳米孔测序K-mer孔模型依次查找T-K+1个K-mer所对应的预期电流信号值，生成输入核苷酸序列X对应的预期测序信号序列Y＝y1,y2,...,yT-5，其中yi表示从X中的位置i开始的K-mer所对应的预期电流信号值；
步骤103：根据实际测序信号序列中的信号值重复次数分布，将预期测序信号序列中的每个信号值重复多次，生成与实际测序信号序列之间具有相似长度分布的待测真实测序信号序列。
进一步地，所述信号处理模型构建模块还用于：
步骤201：将中值归一化后的待测真实测序信号序列I＝i1,i2,...,iT作为信号处理模型中第一层双向门控循环单元神经网络的输入向量，其中双向门控循环单元神经网络的基本单元由一个前向传播的门控循环单元和一个后向传播的门控循环单元组成，首先计算t时刻的前向输出向量
Figure BDA0002864906830000081
计算公式为：
Figure BDA0002864906830000082
Figure BDA0002864906830000083
Figure BDA0002864906830000084
Figure BDA0002864906830000085
其中，it为t时刻的输入信号值，
Figure BDA0002864906830000086
为t-1时刻的前向输出向量，σ为sigmoid函数且σ(z)＝1/(1+e-z)，
Figure BDA0002864906830000087
表示第一中间变量，
Figure BDA0002864906830000088
表示第二中间变量，
Figure BDA0002864906830000089
表示第三中间变量，
Figure BDA00028649068300000822
表示两个等维向量的对应元素相乘，
Figure BDA00028649068300000810
为第一权重矩阵，
Figure BDA00028649068300000811
为第二权重矩阵，Wf为第三权重矩阵，
Figure BDA00028649068300000812
为第一偏置向量，
Figure BDA00028649068300000813
为第二偏置向量，然后计算t时刻的后向输出向量
Figure BDA00028649068300000814
计算公式为：
Figure BDA00028649068300000815
Figure BDA00028649068300000816
Figure BDA00028649068300000817
Figure BDA00028649068300000818
其中，
Figure BDA00028649068300000819
为t+1时刻的后向输出向量，
Figure BDA00028649068300000820
表示第四中间变量，
Figure BDA00028649068300000821
表示第五中间变量，
Figure BDA0002864906830000091
表示第六中间变量，
Figure BDA0002864906830000092
为第四权重矩阵，
Figure BDA0002864906830000093
为第五权重矩阵，Wr为第六权重矩阵，
Figure BDA0002864906830000094
为第三偏置向量，
Figure BDA0002864906830000095
为第四偏置向量，最后t时刻的输出向量ht是由t时刻的前向输出向量
Figure BDA0002864906830000096
和t时刻的后向输出向量
Figure BDA0002864906830000097
连接得到，即
Figure BDA0002864906830000098
其中||表示向量的连接符号；
步骤202：将第一层双向门控循环单元神经网络的输出向量作为第二层双向门控循环单元神经网络的输入向量，第二层双向门控循环单元神经网络的输出向量作为第三层双向门控循环单元神经网络的输入向量，其中第二层双向门控循环单元神经网络、第三层双向门控循环单元神经网络的计算过程与第一层双向门控循环单元神经网络的计算过程相同，不同层之间的权重矩阵与偏置向量参数不同；
步骤203：将三层双向门控循环单元神经网络中最后一层的输出向量作为全连接层的输入向量进行处理，全连接层的输出向量为与输入的待测真实测序信号序列等长的滤波信号序列O＝o1,o2,...,oT。
进一步地，所述信号处理模型训练模块还用于：
步骤301：信号处理模型的监督训练数据包含待训练实际测序信号序列和与其对应的待训练真实测序信号序列，首先读取纳米孔测序平台输出的fast5文件得到待训练实际测序信号序列，然后使用连续小波动态时间规整算法计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列；
步骤302：构建信号处理模型的损失函数为
Figure BDA0002864906830000099
其中，cosh表示双曲函数，其中，oi表示信号处理模型的第i个输出信号，ri表示实际测序信号序列中第i个信号；
步骤303：将每层双向门控循环单元神经网络中权重矩阵以及偏置向量初始化为区间
Figure BDA00028649068300000910
内的均匀分布，其中ni表示每层双向门控循环单元神经网络的输入向量维度，no表示每层双向门控循环单元神经网络的输出向量维度；
步骤304：使用Adam自适应优化器迭代计算得到最小化损失函数的模型参数，完成信号处理模型的参数训练过程。
更进一步地，所述步骤301包括：
步骤3011：读取纳米孔测序平台输出的fast5测序文件得到待训练实际测序信号序列，其中fast5文件为HDF5文件格式；
步骤3012：对待训练实际测序信号序列进行碱基识别得到测序读段，采用基因序列比对算法将测序读段比对到参考基因组序列，根据比对结果得到与待训练实际测序信号序列对应的参考基因组序列片段，然后使用参考基因组序列片段和纳米孔测序K-mer孔模计算出待训练预期测序信号序列；
步骤3013：使用连续小波动态时间规整算法完成待训练实际测序信号序列与待训练预期测序信号序列之间的点到点映射，根据映射结果，求得待训练实际测序信号序列中每个信号值所对应的预期信号值，进而计算出待训练实际测序信号序列所对应的待训练真实测序信号序列。
本发明的优点在于：本发明基于已知的纳米孔测序K-mer孔模型将输入核苷酸序列转换为预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；然后建立真实测序信号序列处理模型，使用监督训练数据完成信号处理模型参数的训练；最后使用完成训练的信号处理模型对待测真实测序信号序列进行滤波处理，利用神经网络学习到实际测序信号序列的时频特性，进而能够准确地滤掉待测真实测序信号序列中与实际测序信号序列无关的高频分量并保留有用的高频分量。
附图说明
图1为本发明实施例所公开的一种基于神经网络的真实纳米孔测序信号滤波方法的示意图；
图2为本发明实施例所公开的一种基于神经网络的真实纳米孔测序信号滤波方法的实际测序信号序列中信号值重复次数分布示意图；
图3为本发明实施例所公开的一种基于神经网络的真实纳米孔测序信号滤波方法中门控循环单元神经网络结构示意图；
图4为本发明实施例所公开的一种基于神经网络的真实纳米孔测序信号滤波方法中三层双向门控循环单元神经网络结构示意图；
图5为本发明实施例所公开的一种基于神经网络的真实纳米孔测序信号滤波方法中信号处理模型参数训练过程中损失函数下降示意图；
图6为本发明实施例所公开的一种基于神经网络的真实纳米孔测序信号滤波方法中信号处理模型的监督训练数据示意图；
图7为现有技术的低通滤波器和本发明信号处理模型输出滤波信号波形对比示意图，图7中(A)为低通滤波器输出信号示意图，图7中(B)为本发明信号处理模型输出信号示意图；
图8为真实测序信号、DeepSimulator仿真测序信号、本发明仿真测序信号波形对比示意图，图8中(A)为真实测序信号波形示意图，图8中(B)为DeepSimulator仿真测序信号波形示意图，图8中(C)为本发明仿真测序信号波形示意图；
图9为DeepSimulator仿真测序信号、本发明仿真测序信号与真实测序信号的DTW距离对比示意图；
图10为真实测序信号、DeepSimulator仿真测序信号和本发明仿真测序信号连续小波变换时频图，图10中(A)为真实测序信号时频图，图10中(B)为DeepSimulator仿真测序信号时频图，图10中(C)为本发明仿真测序信号时频图；
图11为DeepSimulator仿真测序信号、本发明仿真测序信号与真实测序信号连续小波变换时频图之间的PCC比较示意图；
图12为真实测序读段、DeepSimulator仿真测序读段、本发明仿真测序读段的错误特性对比示意图；
图13为真实测序读段、DeepSimulator仿真测序读段、本发明仿真测序读段和肺炎克雷伯菌定制模型仿真测序读段的错误特性对比示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
如图1所示，一种基于神经网络的真实纳米孔测序信号滤波方法，主要以牛津纳米孔技术公司的MinION测序平台(R9.4测序芯片)测序得到的4000条Lambda病毒测序样本、4000条大肠杆菌测序样本、4467条不动杆菌测序样本、15178条肺炎克雷伯菌测序样本、16742条粘质沙雷氏菌测序样本和11047条金黄色葡萄球菌测序样本为例，将Lambda病毒测序样本和大肠杆菌测序样本作为训练数据集也即本发明所述的待训练真实测序信号序列，剩余测序数据作为测试数据集也即本发明所述的待测真实测序信号序列，所述方法包括以下步骤：
步骤S1：在纳米孔测序K-mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；具体过程为：
步骤101：在纳米孔测序K-mer孔模型中输入一条具有T个碱基的核苷酸序列X＝x1,x2,...,xT，xT表示第T个碱基，从输入核苷酸序列的第一个碱基开始，每次移动一个碱基取出一个K-mer，直到移动到倒数第K个碱基位置结束，T个碱基共得到T-K+1个K-mer；需要说明的是，给定一个整数值K，从输入核苷酸序列的的第一个位置开始，取一连续K个碱基的短串，称之为K-mer，本发明实施例中，K取6，因此，本发明采用纳米孔测序6-mer孔模型，每次移动一个碱基取出一个6-mer。
步骤102：纳米孔测序K-mer孔模型包含每个K-mer所对应的预期电流信号值，对照纳米孔测序K-mer孔模型依次查找T-K+1个K-mer所对应的预期电流信号值，生成输入核苷酸序列X对应的预期测序信号序列Y＝y1,y2,...,yT-5，其中yi表示从X中的位置i开始的K-mer所对应的预期电流信号值；
步骤103：如图2所示，根据实际测序信号序列中的信号值重复次数分布，将预期测序信号序列中的每个信号值重复多次，生成与实际测序信号序列之间具有相似长度分布的待测真实测序信号序列。
步骤S2：构建基于双向门控循环单元神经网络的真实测序信号处理模型，该信号处理模型包括三层双向门控循环单元神经网络和一层全连接层，其中信号处理模型的输入为输入到第一层双向门控循环单元神经网络的经中值归一化的真实测序信号序列，信号处理模型的输出为全连接层输出的与输入的待测真实测序信号序列等长的滤波信号序列；具体过程为：
步骤201：将中值归一化后的待测真实测序信号序列I＝i1,i2,...,iT作为信号处理模型中第一层双向门控循环单元神经网络的输入向量，其中双向门控循环单元神经网络的基本单元由一个前向传播的门控循环单元和一个后向传播的门控循环单元组成，首先计算t时刻的前向输出向量
Figure BDA0002864906830000121
如图3所示，计算公式为：
Figure BDA0002864906830000122
Figure BDA0002864906830000123
Figure BDA0002864906830000124
Figure BDA0002864906830000131
其中，it为t时刻的输入信号值，
Figure BDA0002864906830000132
为t-1时刻的前向输出向量，σ为sigmoid函数且σ(z)＝1/(1+e-z)，
Figure BDA0002864906830000133
表示第一中间变量，
Figure BDA0002864906830000134
表示第二中间变量，
Figure BDA0002864906830000135
表示第三中间变量，
Figure BDA00028649068300001326
表示两个等维向量的对应元素相乘，
Figure BDA0002864906830000136
为第一权重矩阵，
Figure BDA0002864906830000137
为第二权重矩阵，Wf为第三权重矩阵，
Figure BDA0002864906830000138
为第一偏置向量，
Figure BDA0002864906830000139
为第二偏置向量，然后计算t时刻的后向输出向量
Figure BDA00028649068300001310
计算公式为：
Figure BDA00028649068300001311
Figure BDA00028649068300001312
Figure BDA00028649068300001313
Figure BDA00028649068300001314
其中，
Figure BDA00028649068300001315
为t+1时刻的后向输出向量，
Figure BDA00028649068300001316
表示第四中间变量，
Figure BDA00028649068300001317
表示第五中间变量，
Figure BDA00028649068300001318
表示第六中间变量，
Figure BDA00028649068300001319
为第四权重矩阵，
Figure BDA00028649068300001320
为第五权重矩阵，Wr为第六权重矩阵，
Figure BDA00028649068300001321
为第三偏置向量，
Figure BDA00028649068300001322
为第四偏置向量，最后t时刻的输出向量ht是由t时刻的前向输出向量
Figure BDA00028649068300001323
和t时刻的后向输出向量
Figure BDA00028649068300001324
连接得到，即
Figure BDA00028649068300001325
其中||表示向量的连接符号；
步骤202：如图4所示，将第一层双向门控循环单元神经网络的输出向量作为第二层双向门控循环单元神经网络的输入向量，第二层双向门控循环单元神经网络的输出向量作为第三层双向门控循环单元神经网络的输入向量，其中第二层双向门控循环单元神经网络、第三层双向门控循环单元神经网络的计算过程与第一层双向门控循环单元神经网络的计算过程相同，不同层之间的权重矩阵与偏置向量参数不同；
步骤203：将三层双向门控循环单元神经网络中最后一层的输出向量作为全连接层的输入向量进行处理，全连接层的输出向量为与输入的待测真实测序信号序列等长的滤波信号序列O＝o1,o2,...,oT。
步骤S3：获取纳米孔测序平台输出的待训练实际测序信号序列并根据待训练实际测序信号序列计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列，将待训练实际测序信号序列和与其对应的待训练真实测序信号序列作为信号处理模型参数训练所需的监督训练数据，构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；具体过程为：
步骤301：信号处理模型的监督训练数据包含待训练实际测序信号序列和与其对应的待训练真实测序信号序列，首先读取纳米孔测序平台输出的fast5文件得到待训练实际测序信号序列，然后使用连续小波动态时间规整算法计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列；
步骤302：构建信号处理模型的损失函数为
Figure BDA0002864906830000141
其中，cosh表示双曲函数，其中，oi表示信号处理模型的第i个输出信号，ri表示实际测序信号序列中第i个信号；
步骤303：将每层双向门控循环单元神经网络中权重矩阵以及偏置向量初始化为区间
Figure BDA0002864906830000142
内的均匀分布，其中ni表示每层双向门控循环单元神经网络的输入向量维度，no表示每层双向门控循环单元神经网络的输出向量维度；
步骤304：使用学习率为0.0001的Adam自适应优化器迭代计算得到最小化损失函数的模型参数，完成信号处理模型的参数训练过程，训练期间将batch_size设置为64，迭代次数设置为1000，参数训练中损失函数下降过程参见图5。
其中，所述步骤301包括：
步骤3011：读取纳米孔测序平台输出的fast5测序文件得到待训练实际测序信号序列，其中fast5文件为HDF5文件格式；
步骤3012：对待训练实际测序信号序列进行碱基识别得到测序读段，采用基因序列比对算法将测序读段比对到参考基因组序列，根据比对结果得到与待训练实际测序信号序列对应的参考基因组序列片段，然后使用参考基因组序列片段和纳米孔测序K-mer孔模计算出待训练预期测序信号序列；
步骤3013：使用连续小波动态时间规整算法完成待训练实际测序信号序列与待训练预期测序信号序列之间的点到点映射，参见图6，根据映射结果，求得待训练实际测序信号序列中每个信号值所对应的预期信号值，进而计算出待训练实际测序信号序列所对应的待训练真实测序信号序列。
步骤S4：将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理。
为了验证本发明提出的基于双向门控循环单元神经网络的信号处理模型可以高效且准确地滤掉标准测序信号序列中无用的高频分量，分别使用DeepSimulator软件信号生成模块的低通滤波器和本发明构建的信号处理模型对标准测序信号序列进行了滤波处理，输出的滤波信号波形对比结果参见图7。根据滤波信号波形对比结果可以看到，与传统的低通滤波器相比本发明提出的基于双向门控循环单元神经网络的信号处理模型只对标准测序信号序列中6-mer变化处的信号进行滤波处理，这更符合实际纳米孔测序的原理。本发明所提出的信号处理模型利用神经网络可以准确地学习到标准测序信号序列与真实测序信号之间的内在关系。
通过在输出滤波信号上添加高斯噪声可以生成最终的仿真测序信号，参见图8，DeepSimulator软件与本发明都可以生成与真实测序信号波形相似的仿真测序信号。为了定量分析两种仿真测序信号与真实测序信号之间的相似度，使用动态时间规整算法进行分析，这是比较两个信号之间差异的标准方法。在测试数据集中随机选取1000条测序数据进行了测试，计算结果显示DeepSimulator软件生成的仿真测序信号与真实测序信号的平均归一化DTW距离为0.132，而本发明生成的仿真测序信号与真实测序信号的平均归一化DTW距离为0.121，比DeepSimulator软件降低了7.5％，这表明本发明可以生成更接近于真实测序信号的仿真测序信号。图9显示了1000条测序样本的实验测试结果，图中每一个点代表一条测序数据，对角线以上的点表示本发明仿真效果更好，而对角线以下的点表示DeepSimulator软件更好，根据统计结果，大约89.1％的测试样本表明本发明方法更好。
除了比较仿真测序信号与真实测序信号之间的DTW距离，同时采用连续小波变换算法对三类测序信号进行了时频分析，得到的连续小波变换时频图如图10所示。总体而言，两种仿真测序信号的时频图与真实测序信号的时频图都非常相似，但本发明生成的仿真测序信号时频图更接近于真实测序信号时频图，特别是高频部分。为了比较两种仿真测序信号时频图与真实测序信号时频图之间的相似度，分别计算了两种仿真测序信号时频图与真实测序信号时频图之间的皮尔逊相关系数(Pearson correlation coefficient，PCC)，并且进一步计算了信号时频图低频部分和高频部分处的PCC。使用测试数据集中4000条测序数据计算的平均值如图11所示。整体而言，本发明生成的仿真测序信号时频图与真实测序信号时频图之间的PCC较DeepSimulator提升了9.08％。对于时频图的高频部分，本发明较DeepSimulator提升了18.88％，这进一步说明了本发明提出的基于双向门控循环单元神经网络的信号处理模型可以更准确地滤掉标准测序信号序列中无用的高频分量，保留有用的高频分量，从而生成时域和频域上更接近于真实测序信号的仿真测序信号。
评估仿真测序信号质量的一个重要指标是比较由仿真测序信号生成的仿真读段是否与实际实验生成的测序读段具有相似的错误特性。使用牛津纳米孔技术公司最新的碱基识别软件Guppy分别对两种仿真测序信号和真实测序信号进行碱基识别，然后使用edlib软件计算三种测序读段的错误特性，包括插入率、删节率以及替代率。使用测试数据集提供的四种生物测序样品获得的统计结果如图12所示。由于实际测序读段的准确性与复杂的测序环境相关，因此仿真测序读段的准确性与实际测序读段的准确性有很大不同。尽管通过不同的测序信号仿真方法生成的两种仿真测序读段的准确性相似，但是在仿真读段中三类错误所占的比例并不一致。为了详细比较两种仿真测序读段和四种真实测序读段的错误特性，进一步计算了六种测序读段中三类错误所占的比例。除了肺炎克雷伯菌测序样本外，本发明所产生的仿真测序读段与实际测序读段具有相似的错误分布。然而，对于所有的生物测序样本而言，由DeepSimulator生成的仿真测序读段与实际测序读段之间的差异较大。针对仿真测序读段的分析结果表明，本发明提出的纳米孔测序信号仿真方法大大提高了仿真测序信号的质量。
对于肺炎克雷伯菌测序样本，由DeepSimulator和本发明提出的仿真方法生成的两种仿真测序读段都与实际测序读段之间有很大的差异。为了进一步探索本发明提出的测序信号仿真方法的可扩展性，我们从肺炎克雷伯菌测序样本中随机选择了4000个测序读段来训练提出的基于双向门控循环单元神经网络的信号处理模型，其余的测序读段用作测试数据集以验证肺炎克雷伯菌定制模型的性能。如图13所示，由肺炎克雷伯菌定制模型产生的仿真测序读段具有与真实测序读段相似的插入率和替代率，而仿真测序读段的删节率小于真实测序读段的删节率。与其他两类仿真测序读段相比，肺炎克雷伯菌定制模型产生的仿真测序读段的错误特性更接近于真实测序读段，这表明所提出的信号仿真方法可以有效地学习到真实测序信号的时频特性，并且比DeepSimulator具有更好的可扩展性。
通过以上技术方案，本发明基于已知的纳米孔测序K-mer孔模型将输入核苷酸序列转换为预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；然后建立真实测序信号序列处理模型，使用监督训练数据完成信号处理模型参数的训练；最后使用完成训练的信号处理模型对待测真实测序信号序列进行滤波处理，利用神经网络学习到实际测序信号序列的时频特性，进而能够准确地滤掉待测真实测序信号序列中与实际测序信号序列无关的高频分量并保留有用的高频分量。
实施例2
与本发明实施例1相对应的，本发明实施例2还提供一种基于神经网络的真实纳米孔测序信号滤波装置，所述装置包括以下步骤：
待测真实测序信号序列生成模块，用于在纳米孔测序K-mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；
信号处理模型构建模块，用于构建基于双向门控循环单元神经网络的真实测序信号处理模型，该信号处理模型包括三层双向门控循环单元神经网络和一层全连接层，其中信号处理模型的输入为输入到第一层双向门控循环单元神经网络的经中值归一化的真实测序信号序列，信号处理模型的输出为全连接层输出的与输入的待测真实测序信号序列等长的滤波信号序列；
信号处理模型训练模块，用于获取纳米孔测序平台输出的待训练实际测序信号序列并根据待训练实际测序信号序列计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列，将待训练实际测序信号序列和与其对应的待训练真实测序信号序列作为信号处理模型参数训练所需的监督训练数据，构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；
滤波处理模块，用于将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理。
具体的，所述待测真实测序信号序列生成模块还用于：
步骤101：在纳米孔测序K-mer孔模型中输入一条具有T个碱基的核苷酸序列X＝x1,x2,...,xT，xT表示第T个碱基，从输入核苷酸序列的第一个碱基开始，每次移动一个碱基取出一个K-mer，直到移动到倒数第K个碱基位置结束，T个碱基共得到T-K+1个K-mer；
步骤102：纳米孔测序K-mer孔模型包含每个K-mer所对应的预期电流信号值，对照纳米孔测序K-mer孔模型依次查找T-K+1个K-mer所对应的预期电流信号值，生成输入核苷酸序列X对应的预期测序信号序列Y＝y1,y2,...,yT-5，其中yi表示从X中的位置i开始的K-mer所对应的预期电流信号值；
步骤103：根据实际测序信号序列中的信号值重复次数分布，将预期测序信号序列中的每个信号值重复多次，生成与实际测序信号序列之间具有相似长度分布的待测真实测序信号序列。
具体的，所述信号处理模型构建模块还用于：
步骤201：将中值归一化后的待测真实测序信号序列I＝i1,i2,...,iT作为信号处理模型中第一层双向门控循环单元神经网络的输入向量，其中双向门控循环单元神经网络的基本单元由一个前向传播的门控循环单元和一个后向传播的门控循环单元组成，首先计算t时刻的前向输出向量
Figure BDA0002864906830000181
计算公式为：
Figure BDA0002864906830000182
Figure BDA0002864906830000183
Figure BDA0002864906830000184
Figure BDA0002864906830000185
其中，it为t时刻的输入信号值，
Figure BDA0002864906830000186
为t-1时刻的前向输出向量，σ为sigmoid函数且σ(z)＝1/(1+e-z)，
Figure BDA0002864906830000187
表示第一中间变量，
Figure BDA0002864906830000188
表示第二中间变量，
Figure BDA0002864906830000189
表示第三中间变量，
Figure BDA00028649068300001822
表示两个等维向量的对应元素相乘，
Figure BDA00028649068300001810
为第一权重矩阵，
Figure BDA00028649068300001811
为第二权重矩阵，Wf为第三权重矩阵，
Figure BDA00028649068300001812
为第一偏置向量，
Figure BDA00028649068300001813
为第二偏置向量，然后计算t时刻的后向输出向量
Figure BDA00028649068300001814
计算公式为：
Figure BDA00028649068300001815
Figure BDA00028649068300001816
Figure BDA00028649068300001817
Figure BDA00028649068300001818
其中，
Figure BDA00028649068300001819
为t+1时刻的后向输出向量，
Figure BDA00028649068300001820
表示第四中间变量，
Figure BDA00028649068300001821
表示第五中间变量，
Figure BDA0002864906830000191
表示第六中间变量，
Figure BDA0002864906830000192
为第四权重矩阵，
Figure BDA0002864906830000193
为第五权重矩阵，Wr为第六权重矩阵，
Figure BDA0002864906830000194
为第三偏置向量，
Figure BDA0002864906830000195
为第四偏置向量，最后t时刻的输出向量ht是由t时刻的前向输出向量
Figure BDA0002864906830000196
和t时刻的后向输出向量
Figure BDA0002864906830000197
连接得到，即
Figure BDA0002864906830000198
其中||表示向量的连接符号；
步骤202：将第一层双向门控循环单元神经网络的输出向量作为第二层双向门控循环单元神经网络的输入向量，第二层双向门控循环单元神经网络的输出向量作为第三层双向门控循环单元神经网络的输入向量，其中第二层双向门控循环单元神经网络、第三层双向门控循环单元神经网络的计算过程与第一层双向门控循环单元神经网络的计算过程相同，不同层之间的权重矩阵与偏置向量参数不同；
步骤203：将三层双向门控循环单元神经网络中最后一层的输出向量作为全连接层的输入向量进行处理，全连接层的输出向量为与输入的待测真实测序信号序列等长的滤波信号序列O＝o1,o2,...,oT。
具体的，所述信号处理模型训练模块还用于：
步骤301：信号处理模型的监督训练数据包含待训练实际测序信号序列和与其对应的待训练真实测序信号序列，首先读取纳米孔测序平台输出的fast5文件得到待训练实际测序信号序列，然后使用连续小波动态时间规整算法计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列；
步骤302：构建信号处理模型的损失函数为
Figure BDA0002864906830000199
其中，cosh表示双曲函数，其中，oi表示信号处理模型的第i个输出信号，ri表示实际测序信号序列中第i个信号；
步骤303：将每层双向门控循环单元神经网络中权重矩阵以及偏置向量初始化为区间
Figure BDA00028649068300001910
内的均匀分布，其中ni表示每层双向门控循环单元神经网络的输入向量维度，no表示每层双向门控循环单元神经网络的输出向量维度；
步骤304：使用Adam自适应优化器迭代计算得到最小化损失函数的模型参数，完成信号处理模型的参数训练过程。
更具体的，所述步骤301包括：
步骤3011：读取纳米孔测序平台输出的fast5测序文件得到待训练实际测序信号序列，其中fast5文件为HDF5文件格式；
步骤3012：对待训练实际测序信号序列进行碱基识别得到测序读段，采用基因序列比对算法将测序读段比对到参考基因组序列，根据比对结果得到与待训练实际测序信号序列对应的参考基因组序列片段，然后使用参考基因组序列片段和纳米孔测序K-mer孔模计算出待训练预期测序信号序列；
步骤3013：使用连续小波动态时间规整算法完成待训练实际测序信号序列与待训练预期测序信号序列之间的点到点映射，根据映射结果，求得待训练实际测序信号序列中每个信号值所对应的预期信号值，进而计算出待训练实际测序信号序列所对应的待训练真实测序信号序列。
以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (10)

1.一种基于神经网络的真实纳米孔测序信号滤波方法，其特征在于，所述方法包括以下步骤：
步骤一：在纳米孔测序K-mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；
步骤二：构建基于双向门控循环单元神经网络的真实测序信号处理模型，该信号处理模型包括三层双向门控循环单元神经网络和一层全连接层，其中信号处理模型的输入为输入到第一层双向门控循环单元神经网络的经中值归一化的真实测序信号序列，信号处理模型的输出为全连接层输出的与输入的待测真实测序信号序列等长的滤波信号序列；
步骤三：获取纳米孔测序平台输出的待训练实际测序信号序列并根据待训练实际测序信号序列计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列，将待训练实际测序信号序列和与其对应的待训练真实测序信号序列作为信号处理模型参数训练所需的监督训练数据，构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；
步骤四：将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理。
2.根据权利要求1所述的一种基于神经网络的真实纳米孔测序信号滤波方法，其特征在于，所述步骤一包括：
步骤101：在纳米孔测序K-mer孔模型中输入一条具有T个碱基的核苷酸序列X＝x1,x2,...,xT，xT表示第T个碱基，从输入核苷酸序列的第一个碱基开始，每次移动一个碱基取出一个K-mer，直到移动到倒数第K个碱基位置结束，T个碱基共得到T-K+1个K-mer；
步骤102：纳米孔测序K-mer孔模型包含每个K-mer所对应的预期电流信号值，对照纳米孔测序K-mer孔模型依次查找T-K+1个K-mer所对应的预期电流信号值，生成输入核苷酸序列X对应的预期测序信号序列Y＝y1,y2,...,yT-5，其中yi表示从X中的位置i开始的K-mer所对应的预期电流信号值；
步骤103：根据实际测序信号序列中的信号值重复次数分布，将预期测序信号序列中的每个信号值重复多次，生成与实际测序信号序列之间具有相似长度分布的待测真实测序信号序列。
3.根据权利要求1所述的一种基于神经网络的真实纳米孔测序信号滤波方法，其特征在于，所述步骤二包括：
步骤201：将中值归一化后的待测真实测序信号序列I＝i1,i2,...,iT作为信号处理模型中第一层双向门控循环单元神经网络的输入向量，其中双向门控循环单元神经网络的基本单元由一个前向传播的门控循环单元和一个后向传播的门控循环单元组成，首先计算t时刻的前向输出向量
Figure FDA0002864906820000021
计算公式为：
Figure FDA0002864906820000022
Figure FDA0002864906820000023
Figure FDA0002864906820000024
Figure FDA0002864906820000025
其中，it为t时刻的输入信号值，
Figure FDA0002864906820000026
为t-1时刻的前向输出向量，σ为sigmoid函数且σ(z)＝1/(1+e-z)，
Figure FDA0002864906820000031
表示第一中间变量，rt f表示第二中间变量，
Figure FDA0002864906820000032
表示第三中间变量，
Figure FDA00028649068200000322
表示两个等维向量的对应元素相乘，
Figure FDA0002864906820000033
为第一权重矩阵，Wr f为第二权重矩阵，Wf为第三权重矩阵，
Figure FDA0002864906820000034
为第一偏置向量，
Figure FDA0002864906820000035
为第二偏置向量，然后计算t时刻的后向输出向量
Figure FDA0002864906820000036
计算公式为：
Figure FDA0002864906820000037
Figure FDA0002864906820000038
Figure FDA0002864906820000039
Figure FDA00028649068200000310
其中，
Figure FDA00028649068200000311
为t+1时刻的后向输出向量，
Figure FDA00028649068200000312
表示第四中间变量，
Figure FDA00028649068200000313
表示第五中间变量，
Figure FDA00028649068200000314
表示第六中间变量，
Figure FDA00028649068200000315
为第四权重矩阵，
Figure FDA00028649068200000316
为第五权重矩阵，Wr为第六权重矩阵，
Figure FDA00028649068200000317
为第三偏置向量，
Figure FDA00028649068200000318
为第四偏置向量，最后t时刻的输出向量ht是由t时刻的前向输出向量
Figure FDA00028649068200000319
和t时刻的后向输出向量
Figure FDA00028649068200000320
连接得到，即
Figure FDA00028649068200000321
其中||表示向量的连接符号；
步骤202：将第一层双向门控循环单元神经网络的输出向量作为第二层双向门控循环单元神经网络的输入向量，第二层双向门控循环单元神经网络的输出向量作为第三层双向门控循环单元神经网络的输入向量，其中第二层双向门控循环单元神经网络、第三层双向门控循环单元神经网络的计算过程与第一层双向门控循环单元神经网络的计算过程相同，不同层之间的权重矩阵与偏置向量参数不同；
步骤203：将三层双向门控循环单元神经网络中最后一层的输出向量作为全连接层的输入向量进行处理，全连接层的输出向量为与输入的待测真实测序信号序列等长的滤波信号序列O＝o1,o2,...,oT。
4.根据权利要求1所述的一种基于神经网络的真实纳米孔测序信号滤波方法，其特征在于，所述步骤三包括：
步骤301：信号处理模型的监督训练数据包含待训练实际测序信号序列和与其对应的待训练真实测序信号序列，首先读取纳米孔测序平台输出的fast5文件得到待训练实际测序信号序列，然后使用连续小波动态时间规整算法计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列；
步骤302：构建信号处理模型的损失函数为
Figure FDA0002864906820000041
其中，cosh表示双曲函数，其中，oi表示信号处理模型的第i个输出信号，ri表示实际测序信号序列中第i个信号；
步骤303：将每层双向门控循环单元神经网络中权重矩阵以及偏置向量初始化为区间
Figure FDA0002864906820000042
内的均匀分布，其中ni表示每层双向门控循环单元神经网络的输入向量维度，no表示每层双向门控循环单元神经网络的输出向量维度；
步骤304：使用Adam自适应优化器迭代计算得到最小化损失函数的模型参数，完成信号处理模型的参数训练过程。
5.根据权利要求4所述的一种基于神经网络的真实纳米孔测序信号滤波方法，其特征在于，所述步骤301包括：
步骤3011：读取纳米孔测序平台输出的fast5测序文件得到待训练实际测序信号序列，其中fast5文件为HDF5文件格式；
步骤3012：对待训练实际测序信号序列进行碱基识别得到测序读段，采用基因序列比对算法将测序读段比对到参考基因组序列，根据比对结果得到与待训练实际测序信号序列对应的参考基因组序列片段，然后使用参考基因组序列片段和纳米孔测序K-mer孔模计算出待训练预期测序信号序列；
步骤3013：使用连续小波动态时间规整算法完成待训练实际测序信号序列与待训练预期测序信号序列之间的点到点映射，根据映射结果，求得待训练实际测序信号序列中每个信号值所对应的预期信号值，进而计算出待训练实际测序信号序列所对应的待训练真实测序信号序列。
6.一种基于神经网络的真实纳米孔测序信号滤波装置，其特征在于，所述装置包括以下步骤：
待测真实测序信号序列生成模块，用于在纳米孔测序K-mer孔模型中输入核苷酸序列将其转换为与其对应的预期测序信号序列，将预期测序信号序列中的每个信号值重复多次生成待测真实测序信号序列；
信号处理模型构建模块，用于构建基于双向门控循环单元神经网络的真实测序信号处理模型，该信号处理模型包括三层双向门控循环单元神经网络和一层全连接层，其中信号处理模型的输入为输入到第一层双向门控循环单元神经网络的经中值归一化的真实测序信号序列，信号处理模型的输出为全连接层输出的与输入的待测真实测序信号序列等长的滤波信号序列；
信号处理模型训练模块，用于获取纳米孔测序平台输出的待训练实际测序信号序列并根据待训练实际测序信号序列计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列，将待训练实际测序信号序列和与其对应的待训练真实测序信号序列作为信号处理模型参数训练所需的监督训练数据，构建信号处理模型的损失函数，进行信号处理模型参数的初始化，通过自适应优化器最小化损失函数实现模型参数的训练；
滤波处理模块，用于将待测真实测序信号序列输入到完成参数训练的信号处理模型中实现对待测真实测序信号序列的滤波处理。
7.根据权利要求6所述的一种基于神经网络的真实纳米孔测序信号滤波装置，其特征在于，所述待测真实测序信号序列生成模块还用于：
步骤101：在纳米孔测序K-mer孔模型中输入一条具有T个碱基的核苷酸序列X＝x1,x2,...,xT，xT表示第T个碱基，从输入核苷酸序列的第一个碱基开始，每次移动一个碱基取出一个K-mer，直到移动到倒数第K个碱基位置结束，T个碱基共得到T-K+1个K-mer；
步骤102：纳米孔测序K-mer孔模型包含每个K-mer所对应的预期电流信号值，对照纳米孔测序K-mer孔模型依次查找T-K+1个K-mer所对应的预期电流信号值，生成输入核苷酸序列X对应的预期测序信号序列Y＝y1,y2,...,yT-5，其中yi表示从X中的位置i开始的K-mer所对应的预期电流信号值；
步骤103：根据实际测序信号序列中的信号值重复次数分布，将预期测序信号序列中的每个信号值重复多次，生成与实际测序信号序列之间具有相似长度分布的待测真实测序信号序列。
8.根据权利要求6所述的一种基于神经网络的真实纳米孔测序信号滤波装置，其特征在于，所述信号处理模型构建模块还用于：
步骤201：将中值归一化后的待测真实测序信号序列I＝i1,i2,...,iT作为信号处理模型中第一层双向门控循环单元神经网络的输入向量，其中双向门控循环单元神经网络的基本单元由一个前向传播的门控循环单元和一个后向传播的门控循环单元组成，首先计算t时刻的前向输出向量
Figure FDA0002864906820000071
计算公式为：
Figure FDA0002864906820000072
Figure FDA0002864906820000073
Figure FDA0002864906820000074
Figure FDA0002864906820000075
其中，it为t时刻的输入信号值，
Figure FDA0002864906820000076
为t-1时刻的前向输出向量，σ为sigmoid函数且σ(z)＝1/(1+e-z)，
Figure FDA0002864906820000077
表示第一中间变量，rt f表示第二中间变量，
Figure FDA0002864906820000078
表示第三中间变量，
Figure FDA00028649068200000727
表示两个等维向量的对应元素相乘，Wz f为第一权重矩阵，Wr f为第二权重矩阵，Wf为第三权重矩阵，
Figure FDA0002864906820000079
为第一偏置向量，
Figure FDA00028649068200000710
为第二偏置向量，然后计算t时刻的后向输出向量
Figure FDA00028649068200000711
计算公式为：
Figure FDA00028649068200000712
Figure FDA00028649068200000713
Figure FDA00028649068200000714
Figure FDA00028649068200000715
其中，
Figure FDA00028649068200000716
为t+1时刻的后向输出向量，
Figure FDA00028649068200000717
表示第四中间变量，
Figure FDA00028649068200000718
表示第五中间变量，
Figure FDA00028649068200000719
表示第六中间变量，
Figure FDA00028649068200000720
为第四权重矩阵，
Figure FDA00028649068200000721
为第五权重矩阵，Wr为第六权重矩阵，
Figure FDA00028649068200000722
为第三偏置向量，
Figure FDA00028649068200000723
为第四偏置向量，最后t时刻的输出向量ht是由t时刻的前向输出向量
Figure FDA00028649068200000724
和t时刻的后向输出向量
Figure FDA00028649068200000725
连接得到，即
Figure FDA00028649068200000726
其中||表示向量的连接符号；
步骤202：将第一层双向门控循环单元神经网络的输出向量作为第二层双向门控循环单元神经网络的输入向量，第二层双向门控循环单元神经网络的输出向量作为第三层双向门控循环单元神经网络的输入向量，其中第二层双向门控循环单元神经网络、第三层双向门控循环单元神经网络的计算过程与第一层双向门控循环单元神经网络的计算过程相同，不同层之间的权重矩阵与偏置向量参数不同；
步骤203：将三层双向门控循环单元神经网络中最后一层的输出向量作为全连接层的输入向量进行处理，全连接层的输出向量为与输入的待测真实测序信号序列等长的滤波信号序列O＝o1,o2,...,oT。
9.根据权利要求6所述的一种基于神经网络的真实纳米孔测序信号滤波装置，其特征在于，所述信号处理模型训练模块还用于：
步骤301：信号处理模型的监督训练数据包含待训练实际测序信号序列和与其对应的待训练真实测序信号序列，首先读取纳米孔测序平台输出的fast5文件得到待训练实际测序信号序列，然后使用连续小波动态时间规整算法计算每个待训练实际测序信号序列所对应的待训练真实测序信号序列；
步骤302：构建信号处理模型的损失函数为
Figure FDA0002864906820000081
其中，cosh表示双曲函数，其中，oi表示信号处理模型的第i个输出信号，ri表示实际测序信号序列中第i个信号；
步骤303：将每层双向门控循环单元神经网络中权重矩阵以及偏置向量初始化为区间
Figure FDA0002864906820000082
内的均匀分布，其中ni表示每层双向门控循环单元神经网络的输入向量维度，no表示每层双向门控循环单元神经网络的输出向量维度；
步骤304：使用Adam自适应优化器迭代计算得到最小化损失函数的模型参数，完成信号处理模型的参数训练过程。
10.根据权利要求9所述的一种基于神经网络的真实纳米孔测序信号滤波装置，其特征在于，所述步骤301包括：
步骤3011：读取纳米孔测序平台输出的fast5测序文件得到待训练实际测序信号序列，其中fast5文件为HDF5文件格式；
步骤3012：对待训练实际测序信号序列进行碱基识别得到测序读段，采用基因序列比对算法将测序读段比对到参考基因组序列，根据比对结果得到与待训练实际测序信号序列对应的参考基因组序列片段，然后使用参考基因组序列片段和纳米孔测序K-mer孔模计算出待训练预期测序信号序列；
步骤3013：使用连续小波动态时间规整算法完成待训练实际测序信号序列与待训练预期测序信号序列之间的点到点映射，根据映射结果，求得待训练实际测序信号序列中每个信号值所对应的预期信号值，进而计算出待训练实际测序信号序列所对应的待训练真实测序信号序列。
CN202011583559.3A 2020-12-28 2020-12-28 一种基于神经网络的真实纳米孔测序信号滤波方法及装置 Withdrawn CN112735524A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011583559.3A CN112735524A (zh) 2020-12-28 2020-12-28 一种基于神经网络的真实纳米孔测序信号滤波方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011583559.3A CN112735524A (zh) 2020-12-28 2020-12-28 一种基于神经网络的真实纳米孔测序信号滤波方法及装置

Publications (1)

Publication Number Publication Date
CN112735524A true CN112735524A (zh) 2021-04-30

Family

ID=75606867

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011583559.3A Withdrawn CN112735524A (zh) 2020-12-28 2020-12-28 一种基于神经网络的真实纳米孔测序信号滤波方法及装置

Country Status (1)

Country Link
CN (1) CN112735524A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113851184A (zh) * 2021-09-29 2021-12-28 湖南工商大学 一种基于人工智能的粪大肠杆菌群数的预测方法及装置
CN117831630A (zh) * 2024-03-05 2024-04-05 北京普译生物科技有限公司 为碱基识别模型构建训练数据集的方法、装置及电子设备

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2001151834A (ja) * 1999-06-07 2001-06-05 Toshiba Corp パターン形成材料、多孔質構造体の製造方法、パターン形成方法、電気化学セル、中空糸フィルター、多孔質カーボン構造体の製造方法、キャパシタの製造方法、および燃料電池の触媒層の製造方法
CN101183899A (zh) * 2007-12-19 2008-05-21 天津大学 基于bp网络用于光纤管道泄漏监测装置的管道安全识别方法
WO2018018038A1 (en) * 2016-07-22 2018-01-25 The Regents Of The University Of California System and method for small molecule accurate recognition technology ("smart")
CN110766070A (zh) * 2019-10-22 2020-02-07 北京威信通信息技术股份有限公司 一种基于循环自编码器的稀少信号识别方法及装置
CN110870020A (zh) * 2017-10-16 2020-03-06 因美纳有限公司 利用卷积神经网络(cnns)进行异常剪接检测
CN111060902A (zh) * 2019-12-30 2020-04-24 电子科技大学 一种基于深度学习的mimo雷达波形设计方法
CN111967524A (zh) * 2020-08-20 2020-11-20 中国石油大学(华东) 基于高斯滤波反馈和空洞卷积的多尺度融合特征增强算法
CN112070774A (zh) * 2020-09-16 2020-12-11 西南石油大学 一种用于页岩数字岩心图像分割的神经网络优化方法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2001151834A (ja) * 1999-06-07 2001-06-05 Toshiba Corp パターン形成材料、多孔質構造体の製造方法、パターン形成方法、電気化学セル、中空糸フィルター、多孔質カーボン構造体の製造方法、キャパシタの製造方法、および燃料電池の触媒層の製造方法
CN101183899A (zh) * 2007-12-19 2008-05-21 天津大学 基于bp网络用于光纤管道泄漏监测装置的管道安全识别方法
WO2018018038A1 (en) * 2016-07-22 2018-01-25 The Regents Of The University Of California System and method for small molecule accurate recognition technology ("smart")
CN110870020A (zh) * 2017-10-16 2020-03-06 因美纳有限公司 利用卷积神经网络(cnns)进行异常剪接检测
CN110766070A (zh) * 2019-10-22 2020-02-07 北京威信通信息技术股份有限公司 一种基于循环自编码器的稀少信号识别方法及装置
CN111060902A (zh) * 2019-12-30 2020-04-24 电子科技大学 一种基于深度学习的mimo雷达波形设计方法
CN111967524A (zh) * 2020-08-20 2020-11-20 中国石油大学(华东) 基于高斯滤波反馈和空洞卷积的多尺度融合特征增强算法
CN112070774A (zh) * 2020-09-16 2020-12-11 西南石油大学 一种用于页岩数字岩心图像分割的神经网络优化方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
CHEN WG 等: "Simulation of nanopore sequencing signals based on bigru", 《SENSORS》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113851184A (zh) * 2021-09-29 2021-12-28 湖南工商大学 一种基于人工智能的粪大肠杆菌群数的预测方法及装置
CN117831630A (zh) * 2024-03-05 2024-04-05 北京普译生物科技有限公司 为碱基识别模型构建训练数据集的方法、装置及电子设备
CN117831630B (zh) * 2024-03-05 2024-05-17 北京普译生物科技有限公司 为碱基识别模型构建训练数据集的方法、装置及电子设备

Similar Documents

Publication Publication Date Title
Strelioff et al. Inferring Markov chains: Bayesian estimation, model comparison, entropy rate, and out-of-class modeling
CN112735524A (zh) 一种基于神经网络的真实纳米孔测序信号滤波方法及装置
RU2004101179A (ru) Анализ ямр-данных многократных измерений на основе максимальной энтропии
CN111564179A (zh) 一种基于三元组神经网络的物种生物学分类方法及系统
WO2021213135A1 (zh) 音频处理方法、装置、电子设备和存储介质
CN114822698B (zh) 一种基于知识推理的生物学大样本数据集分析方法及系统
CN114360662A (zh) 一种基于两路多分支cnn的单步逆合成方法及系统
CN110677297A (zh) 一种基于自回归滑动平均模型和极限学习机的组合网络流量预测方法
CN111058840A (zh) 一种基于高阶神经网络的有机碳含量(toc)评价方法
CN114420217A (zh) 一种新型量子化学分子性能预测的方法和系统
CN117434429B (zh) 芯片的稳定性测试方法及相关装置
CN102697491B (zh) 心电图特征波形识别方法和系统
Jorgensen et al. Operator theory of electrical resistance networks
CN110347579B (zh) 基于神经元输出行为模式的深度学习测试用例的选择方法
CN115598714B (zh) 基于时空耦合神经网络的探地雷达电磁波阻抗反演方法
WO2023124779A1 (zh) 基于三代测序数据检测点突变的分析方法和装置
CN112562788A (zh) 一种环状rna-rna结合蛋白关系预测模型构建方法
CN111859241A (zh) 一种基于声传递函数学习的非监督声源定向方法
CN113283584B (zh) 一种基于孪生网络的知识追踪方法及系统
CN116631499A (zh) 一种基于条件离散扩散模型的核酸适体生成方法
CN116434887A (zh) 一种材料配方设计及性能预测方法
CN115982566A (zh) 一种水电机组多通道故障诊断方法
CN113658643B (zh) 一种基于注意力机制对lncRNA和mRNA的预测方法
Ke et al. Discogen: Learning to discover gene regulatory networks
CN113963757A (zh) 基于气体关系和图神经网络的充油电气设备故障诊断方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
CB02 Change of applicant information
CB02 Change of applicant information

Address after: Building 42, Sino German Cooperation and Innovation Demonstration Park, Intersection of Qingtan Road and Zipeng Road, Economic and Technological Development Zone, Hefei City, Anhui Province, 230000

Applicant after: Hefei Institute of innovation and development, Tianjin University

Address before: 14 / F and 16 / F, comprehensive business building, export processing zone, Hefei Economic Development Zone, Anhui Province (south of Yungu road and east of Taozhi Road, Hefei Economic Development Zone)

Applicant before: Hefei Institute of innovation and development, Tianjin University

WW01 Invention patent application withdrawn after publication
WW01 Invention patent application withdrawn after publication

Application publication date: 20210430