基于卷积神经网络的序列特异性预测研究--云南大学范航恺硕士论文

摘要

基因序列特异性的预测无论是在基因分析领域还是基因调控领域都扮演着重要作用，DNA和RA结合蛋白的特异性模式序列对致病基因的发现也具有指导性作用。

目前生物数据量超级多，传统方法（生物实验、统计分析）耗时耗财耗人，采用深度学习方法可以有效避免这些问题。

本文主要关注点是模体识别问题中模体的序列特异性，针对这个问题，使用卷积神经网络建立了相应的序列预测模型。该模型主要涉及几个阶段：首先将字符串序列转化为数值编码矩阵，再用预测模式的数值矩阵对序列编码矩阵进行卷积运算，从而自动从序列中提取预测模式的特征，卷积输出的矩阵中每一个元素需要与给定的闲值相减，只保留结果大于0的元素传入下一层：再将经过闽值矫正过的输出矩阵输入池化（pooling）层，池化层的作用是对输入数据进行压缩，只保留主要特征：经过池化压缩的矩阵输入一个全连接的神经网络层，这一层通过计算加权和得到一个实数输出分数，将该分数与目标分数相比较，计算出误差，反向传播给前面各层，通过误差调整各层参数，最终拟合得到一个模型，就是序列特异性预测模型。

研究背景

目前生物信息学主要的研究方向有：基因序列分析、致病基因变异的识别、序列特异性的预测及对基因调控方法的建模。 目前最前沿、最热门的方向主要是基因的表达和调控分析。

在基因表达和调控中，基因的中的某些特异性片段对一些特定的结合蛋白具有较强的亲和性，这些生物序列片段也被称为模体（motif），它们在蛋白质合成过程中会起到关键性的作用。因此分析和预测这些特异性片段（即模体），提取其特征模式，对基因分析能够起到指导性的作用。

研究现状

目前国内外对序列特异性预测主要是通过基于各种识别算法来实现的软件评估方法。

模体的序列特异性预测为例，成熟的模体识别算法有Wordup、EM和Gibbs采样等，

其中基于EM（Expectation Maximization）算法Ｗ实现了软件MEME和Imporobizer，基Gibbs采样算法口脚实现了AlignACE、ANN Apec、GLAM、MotifSampler等。

RNA结合位点是一种典型的模体，目前探测这种模体较为常用的是一种称为RAompete的方法。这种方法通过在短RNA序列中探测固定长度的片段模式与相应蛋白质的相互作用，来推断其模式是否具有特异性（即是否是模体）。 借助计算机工具和算法来做序列特异性的发现，是将生物学问题转化为一个字符串处理问题，即在一个大字符串中寻找特定模式的问题。

Wordup算法：与其他位置的基因序列相比，模体序列的出现具有较高的频率。（只适合长度较短的模体，复杂度太高）

MM算法：再EM思想下做了改进。先建立序列的二元有限混合模型，再通过数理统计的方法对参数进行估计。

Gibbs采样：基于假设：每个模体在每一条序列中只出现一次。算法通过对序列进行随机采样，用采样位置的模式不断更新和修正模体在每条序列中出现的位置直到满足收敛条件。

研究意义与所作工作：

降低成本，提高模式发现和预测准确性。

基于CNN做序列特异性预测，处理基因数据并作出预测分析。

改进网络结构

设计针对不同输入序列格式的预处理器，普适性。

第二章 生物学应用背景

模体的概念

在分子生物学领域，转录因子transcription factor指的是能够结合在某特定基因上的、具有特异性核巧酸序列的蛋白质，这些蛋白质的作用是调控基因的转录。 能够与转录因子结合的生物序列称为模体motif，也称为转录因子结合位点(transcription factorbinding sies,TFBSs)它是生物序列特异性的直接物质基础。这些生物序列片段具有独特的序列特异化（sequence specificity）,而这些特异性决定了其具有一些重要的特征。它们一般由数个或数十个碱基组成，位于基因的调控区域中，在转录因子进行基因表达时能与信使RNA(mRNA)结合。不同的模体对不同的蛋白质有其特有的亲和性，即只能与特定的蛋白结合，

在图２－１中，列举了几种模体与特定蛋白质的结合情况，如蛋白质ＣＴＦ只与序列ＧＣＣＡＡＴ发生结合，ＧＣＣＡＡＴ位于整段基因序列的调控区域，它就是一种具有序列特异性的模体。

般把生物进化过程中这些发生了突变的模体称为模体实例（ｍｏｔｉｆｉｎｓｔａｎｃｅ），因为变异了，所以这些模体示例的功能和调控特性也发生了变化。

一般实验获得的是模体实例，无法得到原始模体的实例，所以就引出了模体的表示方法。

模体的表示方法

1.一致模型表示法：将每一个碱基位置上出现次数最多的碱基作为一致序列在该位置上的推断值，生成的一致序列就认为是原始的模体。（丢失信息严重）

PS：

因为模体上发生突变的可能性较小，在生物学中也将其称为保守序列（conservative sequence）。在对模体进行基因特性分析的过程中，需要对模体序列的保守程度进行量化评价，通常使用一致序列与每个模体实例的海明距离（hamming distance）之和来刻画。这个距离之和越小，说明模体的保守程度越高。

当我们需要知道模体在每个碱基位置的保守度和碱基的全局分布情况时，一致模型就不适用了。针对这一问题，我们引出了位置特异性频率－矩阵模型表示法。

2.位置-权重矩阵模型表示法：统计每个碱基位置的每种碱基出现次数，再除Ｗ模体实例的个数，得到每个碱基在模体序列相应位置处出现的概率，将概率值看作对应位置处的权重矩阵。

3.序列logo表示法：碱基所含有的信息量越多，其对应的碱基位置处用来表示碱基的字母就越高。引入序列logo表示法，可以很好地反映各个碱基位点的保守度及该位点对整个模体的重要性。

在logo模型中，输入的是模体的所有模体实例，每个模体实例在相同的碱基位置处碱基不一定相同，反映在logo模型上就是各个模体实例在该位置处所有碱基的堆叠。每个模体实例在某一固定位置上的碱基都对该位置处的信息量有贡献，在该位置处的每个碱基高度都与其在该位置处出现的频率有关，且频率越高排列位置越靠前，在作模体模式预测时该位置处的碱基选取的一般都是最靠前的那一个

模体序列特异性的评价方法

评价序列特异性的质量，确定哪一个是候选模体。基于标准进行打分。

1.Z-scores

2.E-scores

第三章 深度学习方法

相关原理、技术

第四章 序列特异性预测模型的研究与实现

模型框架

感觉这个图画的很乱，个人觉得不需要参考这个图

预处理：输入单链基因序列，将头部留出motif_len-1的位置，推测出另一条链，获得双链数据。

卷积：d个卷积核（motif探测器矩阵）

RELU

池化：maxpooling

全连接：32个神经元，计算模体序列特异性评估分数（这里有点不明白，为啥是32个啊？？？，经过这32个神经元出来32个值，就会有32个分数，在后面分类ud时候，最大是32类？？会影响吗？？？）

softmax:将模体序列特异性评估分数映射成概率值

one-hot编码：将概率值转换为one-hot向量，对序列进行分类。

反向传播：修正参数

数据集选取

RNAcompete数据集（训练和预测RNA结合蛋白序列模式）

包含3部分：

1、sequemces.tsv 含有213130条长度为29nt-38nt的RNA序列文件213130行

2、targets.tsv 含213130条序列对应每个模体的分数，213130行，291列

3、motif 通过RNAcompete方法找到的模体集合文件，291行

共213130条序列，所以假设表格为213130行，291列，每一行代表这每一条序列，有244个分数代表每一行中的序列片段的得分，分数最高的即代表它对应的motif模式存在于该RNA中。

因为数据是通过生物实验得到，所以没有负样本。

通过双碱基随机打乱算法生成负例样本。

（经验：数据集样本>10000，加入的负样本均为10000例）

原始数据（213130个样本）平均划分为AB两个集合，A训练集，B测试集，再给每个集合加10000个负样本，得到最终的数据集A：116565条 B:116565条

校准参数选取

学习率，隐藏层层数......

1、motif探测器长度

2、motif探测器数量

3、学习率[0.0005，0.05]

4、学习势（快速梯度），动态放缩学习率加快参数训练速度。

用的是Nestero势，系数的取值范围为[0.95，0.99]

5、batch_size：64

6、初始化权重：权重值的初始化服从正态分布，对于motif探测器，保证标准值在[1e-7,1e-3]，对于全连接网络，保证标准值在[1e-5,1e-2]

基于卷积神经网络的模型设计

input:一条基因序列s=（s1,...sn），其中si={A,C,G,T,N}，但其实不应该是5阶的吗的，但是很多论文还是选择用4阶的？？？？why

output:score(s)

第五章 预测结果的验证与分析

评价指标

1、皮尔逊积矩相关系数（Pearson product-moment correlation coefficient），度量两个变量X与Y之间的线性相关性，取值[-1,1]

2、斯皮尔曼等级相关系数（Spearman's correlation），衡量两个变量依赖性的非参数指标。一般来说利用单调方程评价两个统计变量的相关性。如果数据中没有重复值且两个变呈完全单调相关时，斯皮尔曼等级相关系数为１或－１。斯皮尔曼相关系数被定义为等级变量之间的皮尔逊相关系数

3、ROC曲线。一般用于二元分类，计算真阳性在结果总数中所占的比例TPR=TP/(TP+FN),伪阳性结果在总数中所占比例为伪阳性率FPR=FP/(FP+TN)，ROC坐标中，伪阳性率FPR定义为X轴，真阳性率TPR定义为Y轴，ROC曲线所围成的面积可以用一个数值表示，即AUC（area under the curve of ROC）,AUC的数值介于0-1，值越大，分类效果越好。

模体预测结果

RNAcompete从213130条RNA序列中找出了291个模体模式，可看作把213130条RNA序列分为291个类，每个类中的元素就是含有相应模体片段的RNA序列，最终我们需要把每一个类（一种模体）中所有元素所含共同的模体片段提取出来。

利用Berkeley开发的基于序列logo表示法的web应用来实现（weblogo.berkeley.edu/logo.cgi），输入同一类RNA序列，输出可视化序列logo。

......

ROC结果：两种方法都比较好。

展望

层数浅

数据集单一