CN114566215B - 一种双端成对的剪接位点预测方法 - Google Patents
一种双端成对的剪接位点预测方法 Download PDFInfo
- Publication number
- CN114566215B CN114566215B CN202210178009.6A CN202210178009A CN114566215B CN 114566215 B CN114566215 B CN 114566215B CN 202210178009 A CN202210178009 A CN 202210178009A CN 114566215 B CN114566215 B CN 114566215B
- Authority
- CN
- China
- Prior art keywords
- sequence
- splice site
- sample
- representing
- model
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 29
- 238000012549 training Methods 0.000 claims abstract description 24
- 238000013527 convolutional neural network Methods 0.000 claims abstract description 17
- 239000013598 vector Substances 0.000 claims abstract description 9
- 238000011176 pooling Methods 0.000 claims description 18
- 238000012795 verification Methods 0.000 claims description 17
- 238000012163 sequencing technique Methods 0.000 claims description 15
- 108020004999 messenger RNA Proteins 0.000 claims description 13
- 238000012545 processing Methods 0.000 claims description 13
- 238000011144 upstream manufacturing Methods 0.000 claims description 13
- 108020004414 DNA Proteins 0.000 claims description 9
- 238000012360 testing method Methods 0.000 claims description 9
- 229930024421 Adenine Natural products 0.000 claims description 6
- GFFGJBXGBJISGV-UHFFFAOYSA-N Adenine Chemical compound NC1=NC=NC2=C1N=CN2 GFFGJBXGBJISGV-UHFFFAOYSA-N 0.000 claims description 6
- 108700024394 Exon Proteins 0.000 claims description 6
- 229960000643 adenine Drugs 0.000 claims description 6
- OPTASPLRGRRNAP-UHFFFAOYSA-N cytosine Chemical compound NC=1C=CNC(=O)N=1 OPTASPLRGRRNAP-UHFFFAOYSA-N 0.000 claims description 6
- UYTPUPDQBNUYGX-UHFFFAOYSA-N guanine Chemical compound O=C1NC(N)=NC2=C1N=CN2 UYTPUPDQBNUYGX-UHFFFAOYSA-N 0.000 claims description 6
- 239000002773 nucleotide Substances 0.000 claims description 6
- 125000003729 nucleotide group Chemical group 0.000 claims description 6
- 230000035945 sensitivity Effects 0.000 claims description 5
- 108091092195 Intron Proteins 0.000 claims description 4
- 238000010276 construction Methods 0.000 claims description 4
- 108091026890 Coding region Proteins 0.000 claims description 3
- 108091028043 Nucleic acid sequence Proteins 0.000 claims description 3
- 239000003795 chemical substances by application Substances 0.000 claims description 3
- 229940104302 cytosine Drugs 0.000 claims description 3
- -1 i.e. Proteins 0.000 claims description 3
- 239000011159 matrix material Substances 0.000 claims description 3
- 241000894007 species Species 0.000 claims description 3
- 239000000126 substance Substances 0.000 abstract description 4
- 238000000605 extraction Methods 0.000 abstract description 2
- 230000006870 function Effects 0.000 description 11
- 238000004422 calculation algorithm Methods 0.000 description 10
- 238000013135 deep learning Methods 0.000 description 6
- 108090000623 proteins and genes Proteins 0.000 description 6
- 238000010801 machine learning Methods 0.000 description 5
- 238000005457 optimization Methods 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 238000011161 development Methods 0.000 description 3
- 238000003559 RNA-seq method Methods 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 230000008569 process Effects 0.000 description 2
- 102000004169 proteins and genes Human genes 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- ORILYTVJVMAKLC-UHFFFAOYSA-N Adamantane Natural products C1C(C2)CC3CC1CC2C3 ORILYTVJVMAKLC-UHFFFAOYSA-N 0.000 description 1
- 241000244203 Caenorhabditis elegans Species 0.000 description 1
- 108091027974 Mature messenger RNA Proteins 0.000 description 1
- 102100029469 WD repeat and HMG-box DNA-binding protein 1 Human genes 0.000 description 1
- 101710097421 WD repeat and HMG-box DNA-binding protein 1 Proteins 0.000 description 1
- 210000000349 chromosome Anatomy 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 238000003066 decision tree Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000003475 lamination Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 231100000590 oncogenic Toxicity 0.000 description 1
- 230000002246 oncogenic effect Effects 0.000 description 1
- 230000001717 pathogenic effect Effects 0.000 description 1
- 238000010827 pathological analysis Methods 0.000 description 1
- 238000007637 random forest analysis Methods 0.000 description 1
- 238000009877 rendering Methods 0.000 description 1
- 238000012216 screening Methods 0.000 description 1
- 230000003068 static effect Effects 0.000 description 1
- 238000007794 visualization technique Methods 0.000 description 1
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
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
- G16B20/30—Detection of binding sites or motifs
-
- 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
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/02—Neural networks
- G06N3/08—Learning methods
- G06N3/084—Backpropagation, e.g. using gradient descent
-
- 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
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Theoretical Computer Science (AREA)
- General Health & Medical Sciences (AREA)
- Biophysics (AREA)
- Software Systems (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Computation (AREA)
- Data Mining & Analysis (AREA)
- Molecular Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Medical Informatics (AREA)
- Biotechnology (AREA)
- Computing Systems (AREA)
- Biomedical Technology (AREA)
- General Physics & Mathematics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Mathematical Physics (AREA)
- Evolutionary Biology (AREA)
- Computational Linguistics (AREA)
- General Engineering & Computer Science (AREA)
- Bioinformatics & Computational Biology (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Public Health (AREA)
- Epidemiology (AREA)
- Databases & Information Systems (AREA)
- Bioethics (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Genetics & Genomics (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种双端成对的剪接位点预测方法,该方法通过获取双端成对的剪接位点样本序列作为基准数据集和独立数据集;将碱基序列通过基于序列本身、物理化学性质等多种特征提取方式进行编码;组合多种特征作为一个多通道多维的向量表示;训练卷积神经网络模型;最后进行评估。这种预测方法可以结合样本多种特征表示方式,帮助卷积神经网络充分学习样本内在模式,提高了双端成对的剪接位点预测的准确率。
Description
技术领域
本发明涉及基因的剪接位点识别预测技术领域,具体是一种双端成对的剪接位点预测方法。
背景技术
随着测序技术的发展,研究人员拿到越来越多的测序下机数据。然而在现阶段,生物体参考基因组上的剪接位点注释尚不完整,还有许多人们未发现的新剪接位点。剪接位点不仅是外显子和内含子边界的分割位置,同时也对外显子之间的连接起到关键作用。外显子连接后的序列为成熟的mRNA,mRNA将经过翻译修饰后表达为蛋白质。如果在错误的位置发生剪接,可能导致基因错误地表达了致病蛋白,致使机体无法完成正常生命活动,甚至可能致癌。因此,正确识别剪接位点是十分关键的研究,不仅能够清晰地认识正常的机体生命活动,还促进基因注释、病理分析、可变剪接、剪接变异等下游分析的研究。
目前,针对剪接位点的研究中,主要包含将供体剪接位点和受体剪接位点单独分为两个模型的传统机器学习方法和深度学习方法,以及近两年提出的同时包含供体剪接位点和受体剪接位点的序列作为一个样本的传统机器学习方法和深度学习方法。在将供体剪接位点和受体剪接位点单独分为两个预测模型的传统机器学习方法中,研究者通过提取剪接位点上游和下游的部分碱基作为数据集,然后提取特征并采用机器学习算法学习样本序列内在信息,构建模型并成功预测。例如,Pertea等人采用了决策树算法,并通过马尔科夫算法对其增强以捕获剪接位点周围信息,开发出模型GeneSplicer。Zhang等人利用带有贝叶斯核的线性SVM算法来区分真假剪接位点。Pashaei等人提出一种结合AdaBoost和FDDM编码方法的混合算法以针对剪接位点问题进行预测。该方法经过实验验证,性能略差于近来兴起的深度学习方法。其中主要原因在于特征提取步骤,研究者未能够输入有效特征,致使模型无法学习关键特征以区分正负样本。
近年来,随着深度学习技术迅猛发展,并伴随着深度学习不需要研究者手动提取特征的巨大优势,研究者引入深度学习技术对剪接位点进行预测。例如,Du等人基于卷积神经网络构建了DeepSS模型,针对人类和秀丽隐杆线虫数据集预测剪接位点。Zuallaert等人基于CNN构建了SpliceRover模型以预测剪接位点,通过算法解释了作者提出的五个假设。Albaradei等人通过改进序列编码方法结合卷积神经网络构建了模型Splice2Deep,并在五个样本集中获得较好的准确率和泛化性。Dutta等人利用双向LSTM算法和可视化技术对剪接位点进行处理,并提供SpliceVisual独立工具。
但是,上述研究针对剪接位点问题的研究中,把供体剪接位点序列作为样本训练为一个独立的模型,受体剪接位点序列作为样本训练为另一个独立的模型。这在很大程度上分裂了供体剪接位点和受体剪接位点的关系,当用户输入一条样本序列进模型,模型只能判断该序列位点是否为供体剪接位点或者是否为受体剪接位点,当拿到预测的供体位点结果时,却无法知道与之成对的受体位点位置,这无法为研究者提供成对的剪接位置,导致研究人员并不能在合适位置把内含子切掉。所以,近年来,有研究者把包含成对的供体位点和受体位点作为样本序列进行预测。
在以包含供体位点和受体位点的序列作为训练样本的研究中,可以看到Mapleson等人基于RNA-Seq和随机森林方法构建模型有效地识别真假剪接序列,Zhang等人基于卷积神经网络训练DeepSplice模型,通过RNA-seq比对后的数据来发现新的剪接位点。这些方法有效地解决了单一供体(受体)模型缺失对应受体(供体)信息的问题,不过上述文章还存在不少缺陷,例如数据样本没有包括非经典剪接位点数据集、数据编码方式单一等。
因此,虽然在剪接位点预测问题上,研究者已经得到了许多突破和改进,但还存在不少问题。第一个问题、数据样本分裂了供体剪接位点和受体剪接位点的关系;第二个问题、模型不能对非经典剪接位点进行处理;第三个问题、数据编码方式单一,忽略了碱基间物理化学性质。
发明内容
本发明的目的在于针对现有剪接位点预测问题存在的缺陷,而提供一种双端成对的剪接位点预测方法。
实现本发明目的的技术方案是:
一种双端成对的剪接位点预测方法,包括如下步骤:
1)以人类参考基因组序列为来源,依据参考基因组序列文件和参考基因组注释文件收集剪接位点序列数据,所述剪接位点序列数据包括经典剪接位点序列和非经典剪接位点序列,对收集到的剪接位点序列数据进行数据处理,包括对数据的长度、内含子和外显子进行区域识别处理,以及正负样本划分处理后,将剪接位点序列数据分为训练集、验证集和测试集;
所述的剪接位点序列数据,包括真剪接位点序列即正样本、假剪接位点序列即负样本,每个数据集中的正负样本数量相等;
所述的数据处理,为了在一个样本序列中同时包含供体剪接位点和受体剪接位点,采用测序对比后的sam文件作为DNA序列注释的文件;所述sam文件为读段对比到参考基因组的输出结果;所述读段是对mRNA测序后的结果,测序结果为多个短序列;由于mRNA是经过DNA转录而来,即mRNA包含的是DNA上的外显子序列,则对mRNA测序后的多个短序列读段为DNA上的外显子序列,将读段对比回到基因组后,确定外显子区域,从而推断外显子区域与相邻未对比到的区域边界为剪接位点;再通过未比对到区域与两头相邻外显子区域的位置信息推断中间未比对到区域和两头相邻外显子区域的边界为对应的供体剪接位点和受体剪接位点;最后,通过获取供体剪接位点上下游序列和受体剪接位点上下游序列作为一个样本序列;样本序列长度为204,包括供体剪接位点GT两个碱基和上下游各50个碱基,包括受体剪接位点AG两个碱基和上下游50个碱基;
2)对步骤1)获得的训练集、验证集和测试集样本进行特征编码,样本序列由腺嘌呤A、腺嘧啶T、胞嘧啶C、鸟嘌呤G和未知N组成,N表示可能为A、T、C、G中任意一个;采用序列位置信息、顺序信息,以及物理化学性质对样本序列进行编码,将序列字符编码为数值格式,采用Mismatch、Kmer、RevKmer、IDKmer、Subsequence、DAC、DCC、DACC、TAC、TCC、TACC、MAC、GAC、NMBAC、PseDNC、PseKNC、PC-PseDNC-General、PC-PseTNC-General、SC-PseDNC-General、SC-PseTNC-General多种特征表示方式,得到样本的特征编码向量;
所述Mismatch特征表示方式,设α是长度为k的子串,(k,m)-mismatch特征图在α上定义为:
其中A表示有限的字母,包括A、C、G或T碱基;若β∈N(k,m)(α),其中β是与α最多不匹配的k-mer子串的集合,然后输入序列的特征图x是k-mer子串的特征向量之和为:
(k,m)-mismatch核定义为特征空间中对应的特征映射的点积:
3)构建卷积神经网络模型,模型的网络构建表达式为:
Lable of class=ffcn(fconv2(fconv1(Sequence nucleotide signal)))
其中Lable of class表示卷积神经网络模型最终的分类,Sequence nucleotidesignal表示碱基序列所对应的输入特征编码,fconv1表示第一层卷积层,fconv2表示第二个卷积层,ffcn表示将输入特征经过卷积步骤后的中间结果传入全连接层;
对于输入x,每个通道上都有一个过滤器ω(1,c),第一个卷积层的第一个过滤器点积运算结果z1,(i,j,k)表示为:
z1,(i,j,k)=(x*ω1,c)i,j,k+b1,(k,1)
其中i,j和c分别表示该卷积层输出的行、列和通道,k为当前层的过滤器,b1,(k,1)表示第一个卷积操作过滤器k的偏置值;
对于通道i,有z(1,i):
z(1,i)=xi*ω1,c(i)+b1,i
三个通道的卷积层输出结果z1,(i,j,c)为:
其中i,j和c分别表示最终输出的行、列和通道;l,m,n分别代表过滤器的行、列和通道,k为当前层使用的过滤器符号表示;
接着进入池化层,池化层分为平均池化、最小池化和最大池化,主要聚合特征映射的空间信息,减小网络内信息传输的向量大小。采用最大池化,保留突出特征,最终通过softmax函数,输出该样本属于每类的预测得分,公式如下:
fi(z)=exp(zi)/∑jexp(zj)
其中fi(z)表示样本属于第i个分类的总预测分数,zj表示属于第j个分类节点的得分,zi表示属于第i个分类节点的得分;
基于预测分值和训练集中标签的实际值计算损失,通过反向传播不断缩小差距,以使得模型性能得到提高,最终获得剪接位点预测模型;每次训练从训练集中获取128个训练样本训练,直到将训练集所有样本输入并训练模型。模型以交叉熵损失函数更新反向传播,进行30次迭代,每次迭代则输入验证集样本对每次迭代训练后的模型进行性能验证,验证集的使用可以提前避免模型进入过拟合和挑选更优的超参数。交叉熵损失函数对于每个类别我们的预测得到的概率为p和1-p,此时交叉熵损失函数L的表达式为:
其中,Li表示样本i的损失函数,N表示样本总数,yi表示i的label,正类为1,负类为0;pi表示样本i预测为正类的概率;
4)将步骤2)进行特征编码后的测试集输入步骤3)得到的训练好的卷积神经网络模型中,获取模型的预测分值并构建混淆矩阵,最终在准确率(Accuracy,ACC)、特异性(specificity,SP)、敏感性(sensitivity,SN)、F分数(F-score)、MCC马修斯相关系数(Matthews correlation coefficient,MCC)和受试者工作特征曲线下面积(area undercurve,AUC)评估五个物种的供体剪接位点和受体剪接位点的性能,表达式如下:
其中TP、TN、FP和FN分别表示真阳性、真阴性、假阳性和假阴性的样本数目。
本发明提供的一种双端成对的剪接位点预测方法,可以将该方法中完成训练和测试好的卷积神经网络模型植入服务器中,并搭建一个剪接位点服务平台,方便其他研究人员可将需要预测的剪接位点序列,直接在网站中可视化拖拽上传而不需冗余下载代码和模型,为剪接位点问题的研究增加便捷性。本发明采用碱基序列的物理化学性质编码结合卷积神经网络对两端成对的剪接位点预测模型性能,与现有技术不同相比,有如下优点:
1、本发明能够有效提升双端成对的剪接位点预测性能。
2、本发明创新性地将基于序列信息,物理化学性质等多种特征编码作为剪接位点特征表示方式,能够更全面获取序列的相关信息。
3、本发明提供双端成对的剪接位点预测平台,能够极大地方便研究者开展双端成对的剪接位点研究。
附图说明
图1为一种双端成对的剪接位点预测方法的总体框架图;
图2为本发明实施例的样本示意图;
图3为本发明实施例的经典剪接位点序列示意图;
图4为本发明实施例的非经典剪接位点序列示意图。
具体实施方式
下面结合附图和实施例对本发明内容做进一步阐述,但不是对本发明的限定。
实施例:
如图1所示,一种双端成对的剪接位点预测方法,包括如下步骤:
1)以人类参考基因组序列为来源,依据参考基因组序列文件和参考基因组注释文件收集剪接位点序列数据,具体的是收集人类剪接位点数据集首先需要从NCBI数据库下载人类参考基因组序列,接着从GenCode数据库下载参考基因组注释文件,结合参考基因组序列和注释文件获取需要的信息。
所述剪接位点序列数据包括经典剪接位点序列和非经典剪接位点序列如图2和图3所示,对收集到的剪接位点序列数据进行数据处理,包括对数据的长度、内含子和外显子进行区域识别处理,以及正负样本划分处理后,将剪接位点序列数据分为训练集、验证集和测试集;
所述的剪接位点序列数据,包括真剪接位点序列即正样本、假剪接位点序列即负样本,每个数据集中的正负样本数量相等。
所述的数据处理,为了在一个样本序列中同时包含供体剪接位点和受体剪接位点,采用测序对比后的sam文件作为DNA序列注释的文件;所述sam文件为读段对比到参考基因组的输出结果;所述读段是对mRNA测序后的结果,测序结果为多个短序列;由于mRNA是经过DNA转录而来,即mRNA包含的是DNA上的外显子序列,则对mRNA测序后的多个短序列读段为DNA上的外显子序列,将读段对比回到基因组后,确定外显子区域,从而推断外显子区域与相邻未对比到的区域边界为剪接位点;再通过未比对到区域与两头相邻外显子区域的位置信息推断中间未比对到区域和两头相邻外显子区域的边界为对应的供体剪接位点和受体剪接位点;最后,通过获取供体剪接位点上下游序列和受体剪接位点上下游序列作为一个样本序列;样本序列长度为204,包括供体剪接位点GT两个碱基和上下游各50个碱基,包括受体剪接位点AG两个碱基和上下游50个碱基;
具体是由于剪接操作从内含子两端剪接,因此剪接位点通常成对出现,如图4所示。通常,注释信息并不是直接给出两端剪接位点的位置,而是给出外显子的起始位置和终止位置。因此,需要依据外显子的位置和基因的位置计算处理得到内含子序列的起始位置,然后利用bedtools工具包从参考基因组序列中提取序列样本,该部分序列样本为正样本,而负样本将从同一染色体不包含剪接位点的序列中提取。
本实施例中,将额外提供对sam文件的处理。sam文件包含的是测序读段比对回到参考基因组的起始位置,比对质量等信息。将依据上述信息进行处理,获取需要预测的序列样本。数据处理过程包括sam文件关键信息识别与提取,samtools工具与bedtools工具结合从参考基因组上获取序列,数据过滤与筛选等。
2)对步骤1)获得的训练集、验证集和测试集样本进行特征编码,样本序列由腺嘌呤A、腺嘧啶T、胞嘧啶C、鸟嘌呤G和未知N组成,N表示可能为A、T、C、G中任意一个;采用序列位置信息、顺序信息,以及物理化学性质对样本序列进行编码,将序列字符编码为数值格式,采用Mismatch、Kmer、RevKmer、IDKmer、Subsequence、DAC、DCC、DACC、TAC、TCC、TACC、MAC、GAC、NMBAC、PseDNC、PseKNC、PC-PseDNC-General、PC-PseTNC-General、SC-PseDNC-General、SC-PseTNC-General多种特征表示方式,得到样本的特征编码向量;
所述Mismatch特征表示方式,设α是长度为k的子串,(k,m)-mismatch特征图在α上定义为:
其中A表示有限的字母,包括A、C、G或T碱基;若β∈N(k,m)(α),其中β是与α最多不匹配的k-mer子串的集合,然后输入序列的特征图x是k-mer子串的特征向量之和为:
(k,m)-mismatch核定义为特征空间中对应的特征映射的点积:
3)构建卷积神经网络模型,模型的网络构建表达式为:
Lable of class=ffcn(fconv2(fconv1(Sequence nucleotide signal)))
其中Lable of class表示卷积神经网络模型最终的分类,Sequence nucleotidesignal表示碱基序列所对应的输入特征编码,fconv1表示第一层卷积层,fconv2表示第二个卷积层,ffcn表示将输入特征经过卷积步骤后的中间结果传入全连接层;
对于输入x,每个通道上都有一个过滤器ω(1,c),第一个卷积层的第一个过滤器点积运算结果z1,(i,j,k)表示为:
z1,(i,j,k)=(x*ω1,c)i,j,k+b1,(k,1)
其中i,j和c分别表示该卷积层输出的行、列和通道,k为当前层的过滤器,b1,(k,1)表示第一个卷积操作过滤器k的偏置值;
对于通道i,有z(1,i):
z(1,i)=xi*ω1,c(i)+b1,i
三个通道的卷积层输出结果z1,(i,j,c)为:
其中i,j和c分别表示最终输出的行、列和通道;l,m,n分别代表过滤器的行、列和通道,k为当前层使用的过滤器符号表示;
接着进入池化层,池化层分为平均池化、最小池化和最大池化,主要聚合特征映射的空间信息,减小网络内信息传输的向量大小。采用最大池化,保留突出特征,最终通过softmax函数,输出该样本属于每类的预测得分,公式如下:
fi(z)=exp(zi)/∑jexp(zj)
其中fi(z)表示样本属于第i个分类的总预测分数,zj表示属于第j个分类节点的得分,zi表示属于第i个分类节点的得分;
基于预测分值和训练集中标签的实际值计算损失,通过反向传播不断缩小差距,以使得模型性能得到提高,最终获得剪接位点预测模型。每次训练从训练集中获取128个训练样本训练,直到将训练集所有样本输入并训练模型。模型以交叉熵损失函数更新反向传播,进行30次迭代。每次迭代则输入验证集样本对每次迭代训练后的模型进行性能验证,验证集的使用可以提前避免模型进入过拟合和挑选更优的超参数。交叉熵损失函数对于每个类别我们的预测得到的概率为p和1-p,此时交叉熵损失函数L的表达式为:
其中,Li表示样本i的损失函数,N表示样本总数,yi表示i的label,正类为1,负类为0;pi表示样本i预测为正类的概率。
此外,模型优化算法能够加速模型收敛,优化算法的每个步骤中,将逐步更新参数猜测值,以减少训练样本中的预测误差。在每个新的猜测中,还会为验证样本构建预测,并且当验证样本错误开始增加时,优化将终止。通过尽早结束参数搜索,参数将朝着初始猜测方向缩小。修正梯度估计可以有效缓解梯度估计值随机性的方式,进而提升优化效率。Adam算法需要计算梯度平方g的指数加权平均和梯度&的指数加权平均。它的参数更新差值Δθn为:
其中,和/>是修正后的加权平均值,α为设定的超参数,学习率ε通常设为0.001。
4)将步骤2)进行特征编码后的测试集输入步骤3)得到的训练好的卷积神经网络模型中,获取模型的预测分值并构建混淆矩阵,最终在准确率(Accuracy,ACC)、特异性(specificity,SP)、敏感性(sensitivity,SN)、F分数(F-score)、MCC马修斯相关系数(Matthews correlation coefficient,MCC)和受试者工作特征曲线下面积(area undercurve,AUC)评估五个物种的供体剪接位点和受体剪接位点的性能,表达式如下:
其中TP为真阳性,代表预测类别为正且真实类别为正的样本数量,TN为真阴性,代表预测类别为负且真实类别为负的样本数量,FP为假阳性,代表预测类别为正,但真实类别为负的样本数量,FN为假阴性,代表预测类别为负但真实类别为正的样本数量。特别是,MCC马修斯相关系数是应用在机器学习中,用以测量二分类的分类性能的指标。该指标考虑了真阳性、真阴性和假阳性和假阴性,通常认为该指标是一个比较均衡的综合指标,即使是在两类别的样本含量差别很大时,同样适用。MCC本质上是一个描述实际分类与预测分类之间的相关系数,它的取值范围为[-1,1],取值为1时表示对受试对象的完美预测,取值为0时表示预测的结果还不如随机预测的结果,-1是指预测分类和实际分类完全不一致。
模型预测步骤中,本发明将训练好的模型放到服务器中,搭建一个高并发、高可用、高性能的剪接位点服务器平台。本发明构建的数据分析平台的界面是由JSP、CSS、JQuery、bootstrap及其扩展包实现。界面后台操作基于Java服务器开发套件,包括Struts2和hibernate。所有剪接位点数据样本和其注释信息存储在MySQL数据库(https:// www.mysql.com/)或者静态文件中,样本数据集方便用户在数据库中查询想要的序列。预测分析功能需要用户提交序列到后端,然后经过一系列耗时的操作,最后将预测或分析结果反馈给用户。本发明引入基于Perl和CGI(https://metacpan.org/pod/cgi)的分布式框架,并行处理预测分析任务,以减少用户提交任务后等待的时间。
Claims (2)
1.一种双端成对的剪接位点预测方法,其特征在于,包括如下步骤:
1)以人类参考基因组序列为来源,依据参考基因组序列文件和参考基因组注释文件收集剪接位点序列数据,所述剪接位点序列数据包括经典剪接位点序列和非经典剪接位点序列,对收集到的剪接位点序列数据进行数据处理,包括对数据的长度、内含子和外显子进行区域识别处理,以及正负样本划分处理后,将剪接位点序列数据分为训练集、验证集和测试集;
2)对步骤1)获得的训练集、验证集和测试集样本进行特征编码,样本序列由腺嘌呤A、腺嘧啶T、胞嘧啶C、鸟嘌呤G和未知N组成,N表示可能为A、T、C、G中任意一个;采用序列位置信息、顺序信息,以及物理化学性质对样本序列进行编码,将序列字符编码为数值格式,采用Mismatch、Kmer、RevKmer、IDKmer、Subsequence、DAC、DCC、DACC、TAC、TCC、TACC、MAC、GAC、NMBAC、PseDNC、PseKNC、PC-PseDNC-General、PC-PseTNC-General、SC-PseDNC-General、SC-PseTNC-General多种特征表示方式,得到样本的特征编码向量;
所述Mismatch特征表示方式,设α是长度为k的子串,(k,m)-mismatch特征图在α上定义为:
其中A表示有限的字母,包括A、C、G或T碱基;若β∈N(k,m)(α),其中β是与α最多不匹配的k-mer子串的集合,然后输入序列的特征图x是k-mer子串的特征向量之和为:
(k,m)-mismatch核定义为特征空间中对应的特征映射的点积:
3)构建卷积神经网络模型,模型的网络构建表达式为:
Lable of class=ffcn(fconv2(fconv1(Sequence nucleotide signal)))
其中Lable of class表示卷积神经网络模型最终的分类,Sequence nucleotidesignal表示碱基序列所对应的输入特征编码,fconv1表示第一层卷积层,fconv2表示第二个卷积层,ffcn表示将输入特征经过卷积步骤后的中间结果传入全连接层;
对于输入x,每个通道上都有一个过滤器ω(1,c),第一个卷积层的第一个过滤器点积运算结果z1,(i,j,k)表示为:
z1,(i,j,k)=(x*ω1,c)i,j,k+b1,(k,1)
其中i,j和c分别表示该卷积层输出的行、列和通道,k为当前层的过滤器,b1,(k,1)表示第一个卷积操作过滤器k的偏置值;
对于通道i,有z(1,i):
z(1,i)=xi*ω1,c(i)+b1,i
三个通道的卷积层输出结果z1,(i,j,c)为:
其中i,j和c分别表示最终输出的行、列和通道;l,m,n分别代表过滤器的行、列和通道,k为当前层使用的过滤器符号表示;
接着进入池化层,池化层分为平均池化、最小池化和最大池化,采用最大池化,保留突出特征,最终通过softmax函数,输出该样本属于每类的预测得分,公式如下:
fi(z)=exp(zi)/∑jexp(zj)
其中fi(z)表示样本属于第i个分类的总预测分数,zj表示属于第j个分类节点的得分,zi表示属于第i个分类节点的得分;
基于预测分值和训练集中标签的实际值计算损失,通过反向传播不断缩小差距,以使得模型性能得到提高,最终获得剪接位点预测模型;每次训练从训练集中获取128个训练样本训练,直到将训练集所有样本输入并训练模型;模型以交叉熵损失函数更新反向传播,进行30次迭代,每次迭代则输入验证集样本对每次迭代训练后的模型进行性能验证,验证集的使用可以提前避免模型进入过拟合和挑选更优的超参数;交叉熵损失函数对于每个类别我们的预测得到的概率为p 和1-p,此时交叉熵损失函数L的表达式为:
其中,Li表示样本i的损失函数,N表示样本总数,yi表示i的label,正类为1,负类为0;pi表示样本i预测为正类的概率;
4)将步骤2)进行特征编码后的测试集输入步骤3)得到的训练好的卷积神经网络模型中,获取模型的预测分值并构建混淆矩阵,最终在准确率ACC、特异性SP、敏感性SN、F分数F-score、马修斯相关系数MCC和受试者工作特征曲线下面积AUC评估五个物种的供体剪接位点和受体剪接位点的性能,表达式如下:
其中TP、TN、FP和FN分别表示真阳性、真阴性、假阳性和假阴性的样本数目。
2.根据权利要求1所述的一种双端成对的剪接位点预测方法,其特征在于,步骤1)中,所述的剪接位点序列数据,包括真剪接位点序列即正样本、假剪接位点序列即负样本,每个数据集中的正负样本数量相等;
所述的数据处理,为了在一个样本序列中同时包含供体剪接位点和受体剪接位点,采用测序对比后的sam文件作为DNA序列注释的文件;所述sam文件为读段对比到参考基因组的输出结果;所述读段是对mRNA测序后的结果,测序结果为多个短序列;由于mRNA是经过DNA转录而来,即mRNA包含的是DNA上的外显子序列,则对mRNA测序后的多个短序列读段为DNA上的外显子序列,将读段对比回到基因组后,确定外显子区域,从而推断外显子区域与相邻未对比到的区域边界为剪接位点;再通过未比对到区域与两头相邻外显子区域的位置信息推断中间未比对到区域和两头相邻外显子区域的边界为对应的供体剪接位点和受体剪接位点;最后,通过获取供体剪接位点上下游序列和受体剪接位点上下游序列作为一个样本序列;样本序列长度为204,包括供体剪接位点GT两个碱基和上下游各50个碱基,包括受体剪接位点AG两个碱基和上下游50个碱基。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210178009.6A CN114566215B (zh) | 2022-02-25 | 2022-02-25 | 一种双端成对的剪接位点预测方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210178009.6A CN114566215B (zh) | 2022-02-25 | 2022-02-25 | 一种双端成对的剪接位点预测方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114566215A CN114566215A (zh) | 2022-05-31 |
CN114566215B true CN114566215B (zh) | 2024-03-22 |
Family
ID=81716797
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210178009.6A Active CN114566215B (zh) | 2022-02-25 | 2022-02-25 | 一种双端成对的剪接位点预测方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114566215B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115579060B (zh) * | 2022-12-08 | 2023-04-04 | 国家超级计算天津中心 | 基因位点检测方法、装置、设备及介质 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110010201A (zh) * | 2019-04-16 | 2019-07-12 | 山东农业大学 | 一种rna选择性剪接位点识别方法及系统 |
CN113178227A (zh) * | 2021-04-30 | 2021-07-27 | 西安交通大学 | 多组学融合剪接位点的识别方法及系统、设备和存储介质 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
SG11201912745WA (en) * | 2017-10-16 | 2020-01-30 | Illumina Inc | Deep learning-based splice site classification |
-
2022
- 2022-02-25 CN CN202210178009.6A patent/CN114566215B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110010201A (zh) * | 2019-04-16 | 2019-07-12 | 山东农业大学 | 一种rna选择性剪接位点识别方法及系统 |
CN113178227A (zh) * | 2021-04-30 | 2021-07-27 | 西安交通大学 | 多组学融合剪接位点的识别方法及系统、设备和存储介质 |
Non-Patent Citations (1)
Title |
---|
基于卷积神经网络的基因剪接位点预测;李国斌;杜秀全;李新路;吴志泽;;盐城工学院学报(自然科学版);20200630(第02期);第20-24页 * |
Also Published As
Publication number | Publication date |
---|---|
CN114566215A (zh) | 2022-05-31 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CA2424031C (en) | System and process for validating, aligning and reordering genetic sequence maps using ordered restriction map | |
KR20200010488A (ko) | 심층 학습 기반 변이체 분류자 | |
US6625545B1 (en) | Method and apparatus for mRNA assembly | |
WO2015061099A1 (en) | Systems and methods for transcriptome analysis | |
CN114566215B (zh) | 一种双端成对的剪接位点预测方法 | |
CN111180013B (zh) | 检测血液病融合基因的装置 | |
CN116959585B (zh) | 基于深度学习的全基因组预测方法 | |
CN115206423A (zh) | 基于标签指导的蛋白质作用关系预测方法 | |
CN113764031A (zh) | 一种跨组织/物种rna中n6甲基腺苷位点的预测方法 | |
CN116710922A (zh) | 变异文献解读知识库的构建方法、解读方法及电子设备 | |
CN116994645B (zh) | 基于交互式推理网络的piRNA与mRNA靶标对的预测方法 | |
CN115995262B (zh) | 基于随机森林及lasso回归解析玉米遗传机理的方法 | |
Zhou et al. | Generalizable prediction of potential miRNA-disease associations based on heterogeneous graph learning | |
Barrabés et al. | Genomic Databases Homogenization with Machine Learning | |
Alavi et al. | scQuery: a web server for comparative analysis of single-cell RNA-seq data | |
CN117976040A (zh) | 变异致病性注释方法、预测变异效应图谱构建方法及系统 | |
Strauch | Improving diagnosis of genetic disease through computational investigation of splicing | |
Denti | Algorithms for analyzing genetic variability from Next-Generation Sequencing data | |
CN116631512A (zh) | 基于深度分解机的piRNA与疾病关联关系预测方法 | |
JP2023550539A (ja) | オリゴヌクレオチド配列における欠失の検出 | |
Kumar et al. | Sequence-Based Algorithms in Computational Biology | |
JP2023540803A (ja) | 多数のノイズのある配列からからコンセンサス配列を生成する深層学習ベースの技法 | |
CN117476252A (zh) | 一种基于知识图谱的病因病理预测方法 | |
CA3200803A1 (en) | Methods for genomic identification of phenotype risk | |
Thulasi | Integrative Computational Analysis of Complex Human Diseases |
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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |