基于深度学习构建时序基因调控网络

1、基于深度学习构建时序基因调控网络摘要：基因调控是生命体中重要的运行机制，运用深度学习来预测P53相关基因的调控关系 对了解基因调控具有重要意义&提出了一种基于卷积神经网络与全连接网络相结合的模 型&在ArrayExpress获得了电离辐射Affymetrix数据集(E-MEXP?549)上，可通过微阵列 表达水平预测基因间调控关系，为保留的验证集提供92. 07%分类准确率，并且该模型的 kappa系数达到0. 84,AUC验证平均精度达到94. 25%，从而构建了带有时延性的p53相 关的基因调控网络&研究结果表明，该模型在筛选出已经论证调控关系的基因对上具有较 好的验证关系，构建出可视化的

2、基因调控网络&在筛选出未知调控关系的基因对上，具有较 好的预测关系与研究价值&关键词：基因调控；高通量数据；深度学习；时间延迟性基因调控网络是控制基因表达的分子相互作用网络&基因间交互错杂的调控关系构成了基因调控网 络(%),这些网络在细胞生命活动的每个阶段都承担着重要的功能，是细胞信号通路的调控者，也是细胞中大 量必需营养物质的控制层&网络中的基因因子之间相互作用，分为促进与抑制两种机制&其中，p53基因 (AB118156)(2系属肿瘤抑制基因家族，是与人类肿瘤基因密切相关的遗传基因，可调节大量目的基因的表 达，进而影响细胞组织的凋亡，分化&随着高通量技术的发展，构建基因调控网络有了更好的

3、数据支撑&高 通量DNA芯片技术进步和几种物种的全基因组序列与基因表达量的可用性，为研究基因调控网络提供了 一种新的策略&例如，Chen等(3)基于动态贝叶斯网络将表观遗传数据(即组蛋白修饰谱)整合到GRN推断 中，提高了 GRN网络的准确性& Khosravi等(4)发现GRN的动态分析可以解释调控因子在不同条件(如癌 症亚型)之间的变化，提出一种最新的逆向工程方法，从每个前列腺状态重建独立的基因调控网络& Khosra? vi等(5)通过观察PcG蛋白通过形成多聚蛋白复合物抑制基因表达，阐明人类PcG蛋白复合物的组成，识别 该复合物结合的靶序列和能够抵消PcG复合物作用的DNA序列&这些方

4、法整合了已知的结合位点信息以 及表达水平、共表达谱，以预测调控相互作用和组装整个网络&近年来，基于神经网络的深度学习模型在各 种分类问题中取得了较好成绩&例如，Li等(6)通过综合分癌症的相关通路，提出了一种新的基于网络的方 法来综合分析DNA甲基化和基因表达数据，以扩展原始的途径& Basavegowda等(7)基于深度学习网络架 构，针对给定的癌症微阵列数据进行分析，对流行的癌症数据进行了高效的分类&这些只使用原始数据的方 法比以前的实验方法获得了更好的性能&与传统的实验方法相比，深度学习方法允许其计算模型得到原始 数据，并自动发现分类所需的复杂表示形式，具有成本低，速度快的优势&本文提出

5、了一种应用于微阵列数 据来预测基因对之间的调控关系的预测方法，该方法采用改进的Fuzzy Nearest Prototype(NFNP)算法动 态筛选数据集，并基于卷积神经网络算法进行了基因间关系的预测，与时间延迟表达谱结合，构建出相应时 延性调控网络&1数据与方法1.1 数据来源采用P53相关基因转录的数据集，构建具有时延性的基因调控网络该数据集含有功能性P53的人 白血病细胞系(M0LT4),通过辐射器每隔#个小时对细胞进行照射，并提取RNA和蛋白质来收集的。实 践过程是对同一个实验独立且细胞制备物相同，同时进行3次，使用Affymetnx U133A微阵列保存，确保 总体的转录反应1.2

6、 方法1.2.1数据预处理 通过BioConduct包中的Limma函数包对微阵列数据集进行预处理，Limma是用于 分析基因表达微阵列数据的软件包，适用于对线性模型分析设计的实验和基因差异表达的评估首先对基 因表达量矩阵进行分组，构建DGEList对象；使用edgeR包中CPM值方法(Count?Per?Million)过滤数据中 的低表达量的基因；利用edgeR中的calcNormFactors()函数对数据标准化；通过estimateDisp()估算离散 值(最终使用Limma包中voom方法进行差异性处理其中删除了信号质量差且在所有时间点上变化较小 的探针从总共22 284个探针中获得

7、了约8 737个探针(.737个基因表达数据)%1.2.2数据筛选采用了成对的Fisher线性判别式来筛选所有类别中具有明显表达差异的基因该方 法通过找到组间平方和与组内平方和的最大比值，寻找组间的良好分离在处理数据时，舍弃一部分数据， 在剩余的数据集里创建模型，进行交叉验证完成该有的功能提出成对Fisher线性判别式(PFLD)：每一次从L基因样本中随机删除部分数据(5%)(对所有的类进行成对两两比较，根据FLD(Fishefs Lmear Discriminant)计算出差异得分d- (+)；整个过程重复P次，计算出最终期望差异d+= .-dp+)/p(4)定义函数F (d(4)定义函数F

8、 (d) =(5)取 F (dE ) 0. 9%这里,P设为100,确保成对Fisher预测比独立预测的特征更加显著为了选出对研究有关的基因范围，设定为了选出对研究有关的基因范围，设定FDR (False Discovery Rate) b 0. 05 和范围0. 5,0. 9# 的 32阈值，同时使用了 maSigPro包处理数据% maSigPro方法也是一种R软件包，专门用于分析 微阵列，该软件包已应用于相同的预处理微阵列数 据集%采用相同的参数设置，使用maSigPro进行 基因选择，筛选出相应基因为了保证选择数据的 鲁棒性，对两种方法交叉论证，发现不同R2中，会 出现基因重叠，如表1

9、%基因筛选数R20.50.60.70.80.9Fisher5232911425013maSigpro6613061394314重合基因数4552631314012重合占比0.690.860.940.930.86表1 Fisher方法与maSigpro方法比较特别是，当R2阈值较高时，maSigPro方法提供的基因与Fisher方法选择的基因重叠较多(* 85%)% 因此，定义前15%为最相关响应探针被认为是一个稳健的选择将经过预处理的探针进一步在每个阵列中 居中并转化为Z分数，筛选对电离辐射最相关的探针选择最相关的响应探针中的前15% (约1312个探 针)作为非线性模型的输入数据1.2.3

10、聚类分析 使用Neural Gas算法12#、NFNP算法13#对1 312个数据进行处理%使用了 Neural Gas 算法和stress function函数14#对数据进行降维，实现了高维输入基因空间到低维子空间的转化%采用 NFNP算法进行最后的数据筛选，算法流程如下所示：开始输入z，将基因向量分类;初始化+ = 1(计算每个原型到Z的距离；计算出C+到Z的距离;将L%到C+的距离段进行排序;初始化+ = 1 ;z为所有的类分配成员;计算距离比E),、1/ M C, ：!愆)=A(1/ Il c, | 1 ,#1/+ ;返回U+(X) ;结束通过上述算法，自动识别最优神经元，最终分配出

11、40组不均匀的具有独特性质共表达基因模块&每个 基因模块代表一组共同表达的基因，这些基因可以通过特定的实验条件或共同的反式调节输入来刺激&对 这些基因模块进行更具体的数据分析，可能会发现它们之间的相互关系以及其在转录调控过程中的复杂机 制&每组内的各个基因相互配对，对其表达水平数据进行皮尔逊相关系数的计算，并筛选相关系数0. 8 的基因对，记录下40组基因模块中表达水平高度相关的基因对作为正样本数据集(5 019对)。同时，筛选记 录出相关系数6 0. 4的基因对，作为负样本数据集(3 930对),进行对照&1.2.4时间延迟基因表达谱在生物体的时间延迟的基因调节是一种常见的现象，多个时间延迟

12、的基因之 间的调控关系存在着一定的联系&这里加入了时间延迟(td)的概念，可以发现多个时间点的基因之间的依 赖关系&对每一组基因构建一个()$ O) X (n X T)的时间延迟基因表达谱TdE(time-delayed expression profiles)矩阵15,其中T列表示每一个基因在t,t ( 1,),t (T 1)时刻的基因表达水平，这样每一行就 是nXT维的向量&当t在T,)-1(范围取值的时候，共产生)T个这样的向量，称作)T个样本& 设定受控基因在t + 1时间点下的表达状态为样本的类别标签&标签的设定为_j2,e, * 0 基因表达水平上调C #$1? , 0基因表达水平

13、下调e 杯=lg(fel) lg (iel)(2)其中，e ,为g,在时间点,下的表达水平，fel为t + 1时刻的最终表达水平，iel为t + 1时刻的初始表达水平& 这样对于每个基因都得到了一个带类别标签的时间延迟基因表达谱D + # (TIE ,C,) &对于每个基因构建好时间延迟基因表达谱后，对预测后的基因对进行时间延迟标签的比对&定义初始 时间延迟td=T + 1,基因对类别标签差异个数为),当基因对中两个基因在t时刻的基因表达水平的调整状态不同时，时间延迟td+m ,直到对比基因对上的每个时刻都对比完；当两个基因在t时刻的基因表达水 平的调整状态一一对应时，记时间延迟td=0&由于

14、T取值不同，所以会出现不同时间点出现不同表达水 平，定义n为T时间点上表达水平类别标签个数,m当)6 n/2(B(为取整函数)时,时间延迟为负&否 则，时间延迟为正&时间延迟td代表了基因对在调控过程中多个时间点上基因之间的依赖关系的延迟&2结果与分析2.1构建基因调控网络选取一组数据集(28个已论证基因)，对其基因对之间的调控关系进行预测&该组基因中，有15个存在 有丝分裂细胞周期过程中，17个存在染色体组织内，分别参与细胞周期的各个过程，具有相互调控的关系& 为了验证结果的生物学意义,考察基因的表达模式和调控关系是否符合细胞周期的特征,对于从模型内提取 预测的调控关系，均与已有的细胞周期基

15、因表达和调控的知识匹配,并定义了更明确的基因间调控的时空关 系。这里对该数据集的基础数据进行了处理，对于每一组基因构建了一个(m-T) X (nXT)的时间延迟基因 表达谱TdE(time-delayed expression profiles)矩阵，并且对于受控基因在7 + 1时刻表达状态进行了类别标 记。对于每个基因构建好时间延迟基因表达谱后，对于预测后的基因对进行时间延迟标签的比对。将基因 时间延迟表达谱与基因模型验证调控关系进行结合，构建基因时延性调控网络。2.2基因调控网络可视化与分析经过分析验证，构建的模型发现了 TP53,JUN,CCNA2这三个基因之间的调控关系&堕是存在时间延

16、迟 的。JUN作为转录因子对CCNA2蛋白的生成具有激活作用，JUN对TP53基因活性具有抑制的作用，而 TP53对CCNA2是具有转录抑制的作用，很明显当JUN表达水平升高，对TP53表达水平会降低，进而影 响TP53对CCNA2的转录抑制的作用，所以CCNA2表达水平会上升。另外一方面,CDCA8作为有丝分裂 关键调节剂，能够诱导细胞有丝分裂，作为着丝粒蛋白的CENPF参与细胞有丝分裂，受CDCA8诱导表达。 当CENPF蛋白增多时，促使了 AURKA激酶的生成。CDCA8调节剂的上升，不仅会间接的促使AURKA 激酶的生成17,同样也会直接影响CENPF蛋白的增多。在一定情况下，进一步说

17、明证实了该模型具有较 高的研究价值。基于不同时间T发现了表达明显的3组基因簇，如图1所示。图1 基因延迟调控网络图(a) CENPE基因调控簇(b) Tp53基因调控簇(c) AURKA基因调控簇图1中所示的每个节点表示一个基因，节点内字母编号为基因的名称；有向边表示基因之间的调控关 系。表示基因间的激活状态，即正调； 表示基因间的抑制状态，即负调。在状态线段边上的数值分别代表T = 1,2,3,4时，不同时刻基因间的调控状态与时间延迟。其中“ + ”号为正调控，“一”号 为负调控，标记的数值为调控反应发生所耗费的时间，即时间延迟。例如, AURKB在调控BIRC5的过程 中，其最终反应调控为

18、正调控。但是，在调控过程中具有不同的调控效果及时间延迟。在时间T = 1,T = 2, T = 3三个时间段内均为正调控，且时间延迟分别为4,3,5个单位时间。在时间T = 4时间段内为负调控, 时间延迟为5个单位时间。3讨论本文开发了一种基于神经网络的模型，该模型用于从微阵列基因表达数据中推断遗传调控机制。模型 的优势在于通过对高通量数据处理和对数据中隐藏特征的提取，可以采用简单的方式深入研究基因表达的 非线性特性。模型加入时间延迟的属性，对于基因调控过程中可能会出现延迟调控现象进行分析，更好的诠 释了基因调控之间的依赖关系，对于p53基因构建新的模型提供了技术支持&18。本研究对于数据集进

19、行了 相关的处理与操作，并对构建的模型进行了相关的模型评估。3.1 模型构建选取了 5 019对相关基因和3 930对不相关基因的表达水平作为数据集，输入的数据是(7,2,1)的向量 (每个基因7次表达水平数据，1对数据为2个基因)。随机打散后选取了其中1/3对基因作为验证集，剩下 作为训练集。由于单元数据维度较小(矩阵维度为7X2X1),使用一个相对较小的点乘网络，建立一个由可 分离式的卷积神经网络(CNN)19和简单的全连接网络(dens)&2。组成的混合/组合网络，如图2所示。Input Convolutional Pooling Fully Connected Output /Acti

20、vation图Input Convolutional Pooling Fully Connected Output /Activation图2 卷积神经网络与全连接网络网络由两层卷积层和两层网络之间的Max Pool层组成#除最后一层网络的激活函数为“sigmoid”外, 其他各层网络的激活函数都为“relu”，模型使用“二值交叉熵(binary cross entropy)”作为损失函数，使用 “RMSprop”作为优化器#32模型优化为了提升模型的收敛速度、提升模型的精度，去除数据的单位限制，转化为无量纲的纯数值，便于不同单 位或量级的指标能够进行比较和加权#采用“sklearn”框架中的

21、“StandardScaler”进行数据归一化处理#神 经网络搭建的过程中涉及到很多超参数，例如,DNN的隐藏层、每层神经元的数量(layersize)；卷积网络的 卷积核(kernelsize)、过滤器大小(filter_size)、Dropout比率。通常，这些超参数通过经验或者在模型训练初期随 机填入，但这样效率不高。现采用超参数搜索的方式确定超参数，由计算机模拟这个过程#使用sklearn框架 中的RandomzedSearchCV方法确定超参数。利用tensorflow2提供的callback进行过拟合处理,记录了训练过 程中各参数的变化，并设置了停止循环的阈值，达到一定准确度自动停

22、止，保证模型准确性#3.3模型评估指标数据TP0.525914FN0.037734FP0.041530TN0.394821表 3 TP、FN、FP、TN 指标为了从侧面客观的评价模型的性能，分别计算了模型的=、 FN、FP、N各项指标数据#其中P、N远大于FN表 3 TP、FN、FP、TN 指标与此同时，为了更好验证模型的统计性能，使用kappa系数作为 展示模型的性能指标，该系数作为衡量分类精度的指标已经取得很 好的验证，能够表示模型整体一致性与分类一致性#最终模型验证 结果为=0. 9207 ： =0. 50.6,2 =0. 83870计算结果证实该模 型具有良好的分类效果与性能#为避免模型的过拟合，增加模型的随机性，将数据集随机打乱，重复三次，最终进行模型验证，提高模型 的鲁棒性#模型数据结果理想值设为0或其中预测值设定范围0%+，设定7偏差7=理想值一预测值, 认为偏差* 0. 5为准确，其