论文部分内容阅读
摘要：针对目前关键基因预测不准确和预测算法缺乏等问题，本文提出一种基于控制理论的关键基因预测算法。首先，从TCGA数据库收集结直肠癌数据，使用计算机工具预处理数据，并利用结直肠癌数据和LncMAP数据库数据构建lncRNA-TF-gene调控网络。然后，设计一种新的筛选方法，基于控制理论中的最小驱动节点集思想和可控性动态分类理论，筛选得到关键节点基因集;将突变得分和网络拓扑分析方法得分融合分析，得到潜在关键基因集。最后，对关键节点基因集和潜在关键基因集取交集，得到关键基因集。结合相关文献和CGC数据库对关键基因集进行验证，证实了该预测算法的有效性，为预测结直肠癌关键基因提供了一种新的思路和方法。
关键词：控制理论;关键基因;结直肠癌;最小驱动节点集
中图分类号：TP301 文献标识码：A
文章编号：1009-3044（2021）30-0028-05
开放科学（资源服务）标识码（OSID）：
Key Gene Prediction Algorithm Based on Control Theory
SONG Zi-jian， YUE Xin-lei， LI Jian-wei
（School of Artificial Intelligence and Data Science， Hebei University of Technology， Tianjin 300130， China）
Abstract：Aiming at the problems of inaccurate prediction of key genes and lack of prediction algorithm， this paper proposes a key gene prediction algorithm based on control theory. Firstly， the data of colorectal cancer were collected from TCGA database， and preprocessed by computer tools. The lncRNA-TF-gene regulatory network was constructed using colorectal cancer data and LncMAP database data. Then， a new screening method is designed， based on the idea of minimum driven node set in control theory and controllable dynamic classification theory， the key node gene set is screened; the mutation score and network topology analysis score are fused to get the potential key gene set. Finally， the intersection of key node gene set and potential key gene set is obtained. Combined with the relevant literature and CGC database to verify the key gene set， the effectiveness of the prediction algorithm is confirmed， which provides a new idea and method for predicting the key genes of colorectal cancer.
Key words：control theory; key genes; colorectal cancer; minimum driven node set
1 引言
癌癥是严重危害中国人民健康的公共卫生问题，且近几年来罹患癌症人数处于不断增长状态。而作为第三大恶性肿瘤的结直肠癌由于具有恶性程度高、病程进展迅速、易复发和转移等特点，对人类健康和生命安全构成重大威胁。结直肠癌的形成始于正常上皮粘膜向过度增生上皮的转化，这些过度增殖的肠上皮细胞由于失去其组织和结构，所以具有形成腺瘤的能力。少部分的结直肠癌病例是遗传性的，而绝大部分是散发性的[1]。研究表明[2]，结直肠癌的发生发展与多基因及多因素有关。而目前对于结直肠癌中的分子机制和关键基因尚未完全了解，导致治疗结直肠癌的研究受到了阻碍。因此，如何更好地研究结直肠癌中的分子机制和关键基因成了一个亟待解决的问题，开发新的关键基因预测算法对于癌症的诊断、治疗和预后都具有重要意义。
关键基因在人类生理过程和疾病发生过程中具有不可忽视的调控作用，了解关键基因在疾病中的功能和作用，对疾病的调控方式、复杂通路、发生发展、治疗和预后等方面具有重大意义[3]。由于现有技术方面的限制，仅通过生物实验测定基因功能所耗费的时间成本和经济成本难以估量。针对这一难题，利用当前现有的数据发现结直肠癌关键基因的预测算法亟待开发。目前，已有许多数据库收集已被实验证实的癌症的关键基因，经整理后供广大的研究人员使用，也有依据单一的算法预测癌症的关键基因。已有文献表明[4]，可以通过基因的突变信息来筛选癌症的关键基因，也有文献研究表明[5]，可以基于生物网络筛选癌症的关键基因，这些理论和研究成果为进一步开发关键基因预测算法提供了帮助。
综上所述，本章提出了一种新的关键基因预测方法，以结直肠癌数据为研究对象进行分析。在该方法中，引入了控制理论的有关思想和方法，首次将突变得分和网络拓扑分析方法得分相结合。实验结果表明，本研究提出的预测算法能有效地预测结直肠癌关键基因，该工作将在预测关键基因方面发挥重要作用。 2 lncRNA-TF-gene网络构建
基于控制理论的关键基因预测算法的总体框架如图1所示。
如图1所示，预测算法的主要流程如下：首先，收集TCGA数据库中结直肠癌的数据。针对以往预测算法中忽略lncRNA的问题，构建lncRNA-TF-gene调控网络，可以更为全面地研究lncRNA与基因在结直肠癌间的作用机制和功能。然后，引入控制理论中的可控性动态分类理论和最小驱动节点集思想，通过分析调控网络建立关键节点基因集，并提出基于融合信息的关键基因筛选方法，将突变得分和网络拓扑分析方法得分相结合，筛选得到关键基因集。最后，结合文献和数据库数据验证关键基因集与结直肠癌之间的关系，以检验预测算法的有效性。
2.1 TCGA结直肠癌数据收集
TCGA计划[6]由国家癌症研究所与国家人类基因组研究所共同展开，对多种癌症的差异基因图谱进行全面绘制，并构建了TCGA数据库对数据进行整理，供广大研究人员使用研究。结直肠癌的数据来源TCGA数据库，数据集预处理过程主要包括如下步骤：1）收集TCGA数据库中结直肠癌数据，利用TCGA官方提供gdc-client软件进行数据下载，得到COAD数据和READ数据。其中，COAD指结肠癌，READ指直肠癌，综合两个数据集进行整合得到完整的结直肠癌数据;2）编写R语言代码对数据采用RMA算法进行背景校正和归一化处理;3）利用Ensembl数据库[7]下载gtf压缩文件，对表达谱矩阵Ensembl ID进行基因名的转化。
2.2 lncRNA-TF-gene调控网络构建
转录因子（transcription factor，TF）是一群能与基因上游特定序列专一性结合，从而保证目的基因以特定的强度在特定的时间与空间表达的蛋白质分子。转录调控对疾病的影响是研究的热点之一，科学家发现lncRNA与转录因子存在转录调控[8]。通过lncRNA-TF-gene的调控网络，可以分析结直肠癌中基因间的调控关系。
而目前常用的lncRNA调控关系数据库是LncMAP数据库[9]，它主要探索人类lncRNA的调控功能，并收集公开的各种测序数据以及lncRNA的相互作用关系数据，并利用线性回归方法预测转录因子与基因之间的关联性。
线性回归[10]是回归分析中第一种经过严格研究并在实际应用中使用的类型，广泛应用于基因调控网络的分析中。通过构建损失函数，来求解损失函数最小时的参数w和b。可以表达成如下公式：
[y=wx+b] （1）
利用LncMAP数据库获取lncRNA-TF-gene调控关系数据，依据阈值为P<0.05，score>0.95进行筛选，得到lncRNA-TF-gene调控关系数据集A，并结合python语言遍历文献后得到lncRNA-TF-gene调控关系数据集B。编写R语言代码整理数据集A和数据集B后得到交集调控关系数据集。利用Cytoscape3.7.2软件构建lncRNA-TF-gene调控网络并进行初步分析，以供后续研究。
3 关键节点集构建
3.1 控制理论介绍
控制理论[11]是以数学方法和计算机技术为主要工具，研究各种控制策略及控制系统的理论、方法和技术。它通过对网络拓扑结构及网络环境下先进控制理论与方法的研究，充分利用网络资源，实现从决策到控制的全过程优化。基于控制理论的复杂网络研究为生物信息学提供了新的思路，复杂网络是从真实系统中抽象得到的一类具有复杂拓扑结构的网络模型。现实世界的真实系统大多是非线性系统，但非线性关系很难被刻画，因此借鉴线性系统的动力学，可以研究非线性的复杂网络[12]。根据控制理论，包含 N 个节点的线性系统动力学方程如下：
[dxtdt=Axt+But ]（2）
其中，向量[xt]表示网络内的N个节点在t时刻的状态。[ut]表示外部控制的M个节点在[t]时刻的状态。若存在一个连续的输入[ut]，使得在有限的时间内，能使得线性时不变系统能够从任意的初始状态达到任意的终止状态，则称这个系统是可控的。
3.2 最小驱动节点集思想
Liu等[13]结合图论方法和结构可控性定理，提出了一个基于图模型的最大匹配算法求解最小驱动节点集的复杂网络可控性分析框架。該框架提出通过最小驱动节点集可以实现控制整个网络的目的，从理论上证明了最小驱动节点集就是最大匹配中的未匹配节点集。最小驱动节点集中节点数量越多，意味着复杂网络的可控性越低，反之亦然。此外，最小驱动节点集的节点数量由网络中的度决定，网络的度越大，最小驱动节点集中节点数量越少，复杂网络的可控性就越低。
3.3 可控性动态分类理论
可控性动态分类理论基于Kalman判定条件[14]和网络的拓扑结构判断整个网络是否满足可控性。核心思想是将网络看作矩阵，将矩阵中的所有元素设为0或已确定的常数。若一个网络是可控的，那么网络中边的权重不影响网络的可控性。Kalman判定条件如下：
设A为网络的连接矩阵，B为外部输入矩阵，那么定义矩阵C为：
[C=B，AB，A2B，…，AN?1B] （3） 当矩阵C为满秩时，即
[rankC=rankB，AB，A2B，…，AN?1B=N]时，网络是可控的。
Vin等人基于可控性动态分类理论分析了蛋白质互作网络的可控性[15]，并对Liu等人的研究进行了改进，将蛋白质相互作用网络中的每一个节点逐个删除，检验最小驱动节点集中节点数量的变化，若节点数量增加，则说明删除的节点为关键节点;若节点数量保持不变，则说明删除的节点为普通节点;若节点数量减少，则说明删除的节点为冗余节点。实验结果表明，关键节点代表的蛋白质是患病的主要因素。而普通节点和冗余节点这两种类型的蛋白质的代表性不足。
3.4 关键节点基因集构建流程
基于上述理论，借助Liu[16]的框架代码，结合readr包后导入结直肠癌的表达谱数据，编写R语言代码分析构建lncRNA-TF-gene调控网络，阈值设为P<0.05，分析后构建最小驱动节点集。基于可控性动态分类理论，移除最小驱动节点集的每个节点，依据移除后对调控网络的影响，将最小驱动节点集的节点进行分类，分为关键节点、普通节点和冗余节点。利用write.csv函数导出分类结果和关键节点基因，以便后续分析。
4 基于融合信息筛选关键基因
针对目前筛选关键基因方法单一的问题，本节研究将突变频率和网络拓扑分析方法得分进行融合，利用突变得分和网络拓扑分析方法得分的融合信息对lncRNA-TF-gene调控网络进行分析，得到潜在关键基因，以便后续研究。
4.1 突变得分
突变得分通常指每一个配子发生突变的概率[17]，即用一定数目配子中的突变配子数表示。在无性繁殖的细菌中，突变得分代表用每一个细胞世代中每个细菌发生突变的概率，即用一定数目的细菌在一次分裂过程中发生突变的次数表示。不同生物或同一生物个体的不同基因的自发突变率是不相同的。突变得分的估算方法因生物生殖方式的不同而不同。在有性生殖的生物中，突变率通常为每一个配子发生突变的概率，即用一定数目配子中的突变配子数表示。在无性繁殖的细菌中，突变率代表每一个细胞世代中每个细菌发生突变的概率，即用一定数目的细菌在一次分裂过程中发生突变的次数表示。
通过TCGA数据库下载结直肠癌对应突变数据集，基于突变的变体分类，仅选择具有功能的突变，例如splice_site，in_frame_del和frame_shift_del等，并計算突变频率。基因的变异越频繁，在排名列表中就越高。根据TCGA数据的突变频率计算调控网络中基因的突变得分，利用Python语言排序后整理，得到基因突变得分列表。
4.2 网络拓扑分析方法得分
网络拓扑统计参数反映网络的结构和动力学特性，对网络结构的分析一般包括以下几个方面：将交互关系复杂的系统建模为网络，然后计算网络拓扑参数来研究网络的拓扑结构特性，根据网络演化过程中网络拓扑参数的变化来分析网络演化对网络结构的影响[18]。
本研究通过3个网络拓扑方法对潜在mRNA关键基因进行联合分析，即degree、MNC和MCC网络拓扑分析方法[19]。degree代表节点的度，v为图中节点数目，MNC为最大邻域分量，其公式为：
[MNCv=VMCv] （4）
其中[V ]指网络中节点的总和，[MCv]指子图与子图间的最大连接分量。
MCC网络拓扑分析方法可以很好地分析节点数目较多的复杂网络。给定一个节点v，MCC的定义为：
[MCCv=C∈SvC?1！] （5）
其中，[C]为属于[Sv]的一个团，[Sv]为包含节点[v]的最大团的集合，[C?1]表示小于[C]的所有正整数的阶乘。
利用Cytoscape3.7.2软件的cytoHubba插件中的MNC、MCC和Degree网络拓扑分析方法分析lncRNA-TF-gene调控网络，根据网络拓扑分析方法得分进行排序，利用python语言整理后得到网络拓扑分析方法得分列表。
4.3 潜在关键基因筛选
若一个基因可以同时出现在突变得分列表和网络拓扑分析方法得分列表中，则认为该基因是潜在关键基因。因此，利用python语言的numpy程序包[20]分别统计突变得分列表和网络拓扑分析方法得分列表中排名前50名的基因，并通过R语言对两个列表下排名前50名的基因取交集，整理后得到潜在关键基因集。
4.4 关键基因筛选及验证
结合上文的关键节点基因集和潜在关键基因集，通过R语言找到两个基因集的共同基因，整理后命名为关键基因集，关键基因在结直肠癌的发病过程中起到了重要的调控作用。为了验证实验结果的有效性，利用python遍历文献中的实验数据，结合CGC数据库[21]中的数据对关键基因集中的基因进行验证。利用两方面数据对关键基因集进行验证，以便对预测算法的效果进行评估。
5 实验结果与分析
5.1 lncRNA-TF-gene调控网络
利用TCGA数据库共得到60477条基因表达谱数据，经背景校正、归一化处理、基因名转化后，整理为csv格式，以便后续分析。
利用LncMAP数据库共获取到102111条结直肠癌中lncRNA-TF-gene调控关系数据，整理后命名为lncRNA-TF-gene调控关系数据集A，利用python语言共收集到32011条结直肠癌中lncRNA-TF-gene调控关系数据，整理后命名为lncRNA-TF-gene调控关系数据集B。将调控关系数据集A与调控关系数据集B取交集后得到10521条调控关系数据，整理后命名为交集调控关系数据集，部分数据如表1所示。结果表明，lncRNA-TF-gene调控关系显著相关并且具有显著差异性，可作为调控网络进行研究。利用Cytoscape3.7.2软件构建lncRNA-TF-gene调控网络，并采用NetworkAnalyzer功能对调控网络进行分析。 5.2 關键节点基因集
基于最小驱动节点集思想，对lncRNA-TF-gene调控网络进行筛选分析，共得到1021个驱动节点基因，整理后命名为最小驱动节点集。依据可控性动态分类理论，对最小驱动节点集进行分类，分类后的结果如表2所示，扇形图如图2所示。可控性动态分类筛选出的14%的关键节点基因，可作为研究潜在的关键基因。关键节点基因集的部分结果如表3所示，其中，K为lncRNA-TF-gene调控中节点的度，Kin为节点的入度，Kout为节点的出度，TypeI为冗余和普通节点基因，TypeII为关键基因，0表示否，1表示是。
5.3 潜在关键基因集
利用R语言下载TCGA数据库中结直肠癌的临床数据，保存为RData格式。编写R语言代码计算突变得分，利用python语言的numpy程序包整理后得到突变得分列表，得分前10名的基因如表4所示。利用Cytoscape3.7.2软件的cytoHubba插件的MNC、MCC和Degree网络拓扑指标联合分析lncRNA-TF-gene调控网络，利用Python语言的numpy程序包整理后得到网络拓扑分析方法得分列表，得分前10名的基因如表5所示。
利用Python语言的numpy程序包分别统计突变得分列表和网络拓扑分析方法得分列表中排名前50名的基因，对两个列表中排名前50名的基因集取交集，得到36个潜在关键基因。
5.4 关键基因集及验证
结合上文的关键节点基因集和潜在关键基因集，对两个基因集取交集，共得到17个关键基因（CASP3、JAK1、CDK2、CCL5、FN1、PCNA、MAPK8、SOD1、PPARG、VEGFA、TGFB1、CDKN1A、AKT1、TP53、BCL2、MAPK3和MAPK1）。为了验证实验的有效性，通过CGC网站和大量文献对关键基因进行验证，发现17个关键基因中有13个已经通过CGC网站和相关文献得到证实，与结直肠癌有密切关系。这表明其余4个关键基因有很大可能与结直肠癌也密切相关。这些关键基因对未来结直肠癌机理研究和开发靶向药物均具有重要参考价值。
6 结论
本文提出了一种基于控制理论的结直肠癌关键基因预测算法，本算法基于突变得分和网络拓扑分析方法得分的融合信息进行关键基因筛选。实验结果证实，关键基因的预测准确性较以往的预测算法得到了显著提高。本算法不仅可以应用在结直肠癌的关键基因预测，还可以广泛应用在各类癌症的关键基因预测中，为今后发现癌症的关键基因提供新的一种思路和选择。
参考文献：
[1] Cancho V G，Bazán J L，Dey D K.A new class of regression model for a bounded response with application in the study of the incidence rate of colorectal cancer[J].Statistical Methods in Medical Research，2020，29（7）：2015-2033.
[2] Maduni? K，Zhang T，Mayboroda O A，et al.Colorectal cancer cell lines show striking diversity of their O-glycome reflecting the cellular differentiation phenotype[J].Cellular and Molecular Life Sciences，2021，78（1）：337-350.
[3] Baratti D，Kusamura S，Pietrantonio F，et al.Progress in treatments for colorectal cancer peritoneal metastases during the years 2010-2015.A systematic review[J].Critical Reviews in Oncology/Hematology，2016，100：209-222.
[4] Reimand J，Bader G D.Systematic analysis of somatic mutations in phosphorylation signaling predicts novel cancer drivers[J].Molecular Systems Biology，2013，9（1）：637.
[5] Hou J P，Ma J.DawnRank：discovering personalized driver genes in cancer[J].Genome Med，2014，6（7）：56.
[6] Network C G A R，Weinstein J N，Collisson E A，et al.The cancer genome atlas pan-cancer analysis project[J].Nature Genetics，2013，45（10）：1113-1120.
[7] Flicek P，Amode M R，Barrell D，et al.Ensembl 2012[J].Nucleic Acids Research，2012，40（database issue）：D84-D90.
[8] Ji Y，Gu Y Q，Hong S，et al.Comprehensive analysis of lncRNA-TF crosstalks and identification of prognostic regulatory feedback loops of glioblastoma using lncRNA/TF-mediated ceRNA network[J].Journal of Cellular Biochemistry，2020，121（1）：755-767. [9] Li Y S，Li L L，Wang Z S，et al.LncMAP：Pan-cancer atlas of long noncoding RNA-mediated transcriptional network perturbations[J].Nucleic Acids Research，2018，46（3）：1113-1123.
[10] Calinski T，Steel R G D，Torrie J H.Principles and procedures of statistics：a biometrical approach[J].Biometrics，1981，37（4）：859.
[11] Narendra K S.Neural networks for control theory and practice[J].Proceedings of the IEEE，1996，84（10）：1385-1406.
[12] 王婉宁.基于复杂网络结构可控性的致病基因识别方法研究[D].西安：西安理工大学，2019.
[13] McDowell J J.Matching theory in natural human environments[J].The Behavior Analyst，1988，11（2）：95-109.
[14]WELCH G， BISHOP G. An Introduction to the Kalman Filter[M]. University of North Carolina at Chapel Hill， 1995.
[15] Vinayagam A，Gibson T E，Lee H J，et al.Controllability analysis of the directed human protein interaction network identifies disease genes and drug targets[J].Proceedings of the National Academy of Sciences of the United States of America，2016，113（18）：4976-4981.
[16] Liu Y Y，Slotine J J，Barabási A L.Controllability of complex networks[J].Nature，2011，473（7346）：167-173.
[17] Devarakonda S，Rotolo F，Tsao M S，et al.Tumor mutation burden as a biomarker in resected non-small-cell lung cancer[J].Journal of Clinical Oncology，2018，36（30）：2995-3006.
[18] 周庚.复杂网络节点中心性度量算法的研究及应用[D].兰州：兰州理工大學，2020.
[19] Chin C H，Chen S H，Wu H H，et al.cytoHubba：identifying hub objects and sub-networks from complex interactome[J].BMC Systems Biology，2014，8（Suppl 4）：S11.
[20] van der Walt S，Colbert S C，Varoquaux G.The NumPy array：a structure for efficient numerical computation[J].Computing in Science & Engineering，2011，13（2）：22-30.
[21] Sondka Z，Bamford S，Cole C G，et al.The COSMIC Cancer Gene Census：describing genetic dysfunction across all human cancers[J].Nature Reviews Cancer，2018，18（11）：696-705.
【通联编辑：唐一东】
关键词：控制理论;关键基因;结直肠癌;最小驱动节点集
中图分类号：TP301 文献标识码：A
文章编号：1009-3044（2021）30-0028-05
开放科学（资源服务）标识码（OSID）：
Key Gene Prediction Algorithm Based on Control Theory
SONG Zi-jian， YUE Xin-lei， LI Jian-wei
（School of Artificial Intelligence and Data Science， Hebei University of Technology， Tianjin 300130， China）
Abstract：Aiming at the problems of inaccurate prediction of key genes and lack of prediction algorithm， this paper proposes a key gene prediction algorithm based on control theory. Firstly， the data of colorectal cancer were collected from TCGA database， and preprocessed by computer tools. The lncRNA-TF-gene regulatory network was constructed using colorectal cancer data and LncMAP database data. Then， a new screening method is designed， based on the idea of minimum driven node set in control theory and controllable dynamic classification theory， the key node gene set is screened; the mutation score and network topology analysis score are fused to get the potential key gene set. Finally， the intersection of key node gene set and potential key gene set is obtained. Combined with the relevant literature and CGC database to verify the key gene set， the effectiveness of the prediction algorithm is confirmed， which provides a new idea and method for predicting the key genes of colorectal cancer.
Key words：control theory; key genes; colorectal cancer; minimum driven node set
1 引言
癌癥是严重危害中国人民健康的公共卫生问题，且近几年来罹患癌症人数处于不断增长状态。而作为第三大恶性肿瘤的结直肠癌由于具有恶性程度高、病程进展迅速、易复发和转移等特点，对人类健康和生命安全构成重大威胁。结直肠癌的形成始于正常上皮粘膜向过度增生上皮的转化，这些过度增殖的肠上皮细胞由于失去其组织和结构，所以具有形成腺瘤的能力。少部分的结直肠癌病例是遗传性的，而绝大部分是散发性的[1]。研究表明[2]，结直肠癌的发生发展与多基因及多因素有关。而目前对于结直肠癌中的分子机制和关键基因尚未完全了解，导致治疗结直肠癌的研究受到了阻碍。因此，如何更好地研究结直肠癌中的分子机制和关键基因成了一个亟待解决的问题，开发新的关键基因预测算法对于癌症的诊断、治疗和预后都具有重要意义。
关键基因在人类生理过程和疾病发生过程中具有不可忽视的调控作用，了解关键基因在疾病中的功能和作用，对疾病的调控方式、复杂通路、发生发展、治疗和预后等方面具有重大意义[3]。由于现有技术方面的限制，仅通过生物实验测定基因功能所耗费的时间成本和经济成本难以估量。针对这一难题，利用当前现有的数据发现结直肠癌关键基因的预测算法亟待开发。目前，已有许多数据库收集已被实验证实的癌症的关键基因，经整理后供广大的研究人员使用，也有依据单一的算法预测癌症的关键基因。已有文献表明[4]，可以通过基因的突变信息来筛选癌症的关键基因，也有文献研究表明[5]，可以基于生物网络筛选癌症的关键基因，这些理论和研究成果为进一步开发关键基因预测算法提供了帮助。
综上所述，本章提出了一种新的关键基因预测方法，以结直肠癌数据为研究对象进行分析。在该方法中，引入了控制理论的有关思想和方法，首次将突变得分和网络拓扑分析方法得分相结合。实验结果表明，本研究提出的预测算法能有效地预测结直肠癌关键基因，该工作将在预测关键基因方面发挥重要作用。 2 lncRNA-TF-gene网络构建
基于控制理论的关键基因预测算法的总体框架如图1所示。
如图1所示，预测算法的主要流程如下：首先，收集TCGA数据库中结直肠癌的数据。针对以往预测算法中忽略lncRNA的问题，构建lncRNA-TF-gene调控网络，可以更为全面地研究lncRNA与基因在结直肠癌间的作用机制和功能。然后，引入控制理论中的可控性动态分类理论和最小驱动节点集思想，通过分析调控网络建立关键节点基因集，并提出基于融合信息的关键基因筛选方法，将突变得分和网络拓扑分析方法得分相结合，筛选得到关键基因集。最后，结合文献和数据库数据验证关键基因集与结直肠癌之间的关系，以检验预测算法的有效性。
2.1 TCGA结直肠癌数据收集
TCGA计划[6]由国家癌症研究所与国家人类基因组研究所共同展开，对多种癌症的差异基因图谱进行全面绘制，并构建了TCGA数据库对数据进行整理，供广大研究人员使用研究。结直肠癌的数据来源TCGA数据库，数据集预处理过程主要包括如下步骤：1）收集TCGA数据库中结直肠癌数据，利用TCGA官方提供gdc-client软件进行数据下载，得到COAD数据和READ数据。其中，COAD指结肠癌，READ指直肠癌，综合两个数据集进行整合得到完整的结直肠癌数据;2）编写R语言代码对数据采用RMA算法进行背景校正和归一化处理;3）利用Ensembl数据库[7]下载gtf压缩文件，对表达谱矩阵Ensembl ID进行基因名的转化。
2.2 lncRNA-TF-gene调控网络构建
转录因子（transcription factor，TF）是一群能与基因上游特定序列专一性结合，从而保证目的基因以特定的强度在特定的时间与空间表达的蛋白质分子。转录调控对疾病的影响是研究的热点之一，科学家发现lncRNA与转录因子存在转录调控[8]。通过lncRNA-TF-gene的调控网络，可以分析结直肠癌中基因间的调控关系。
而目前常用的lncRNA调控关系数据库是LncMAP数据库[9]，它主要探索人类lncRNA的调控功能，并收集公开的各种测序数据以及lncRNA的相互作用关系数据，并利用线性回归方法预测转录因子与基因之间的关联性。
线性回归[10]是回归分析中第一种经过严格研究并在实际应用中使用的类型，广泛应用于基因调控网络的分析中。通过构建损失函数，来求解损失函数最小时的参数w和b。可以表达成如下公式：
[y=wx+b] （1）
利用LncMAP数据库获取lncRNA-TF-gene调控关系数据，依据阈值为P<0.05，score>0.95进行筛选，得到lncRNA-TF-gene调控关系数据集A，并结合python语言遍历文献后得到lncRNA-TF-gene调控关系数据集B。编写R语言代码整理数据集A和数据集B后得到交集调控关系数据集。利用Cytoscape3.7.2软件构建lncRNA-TF-gene调控网络并进行初步分析，以供后续研究。
3 关键节点集构建
3.1 控制理论介绍
控制理论[11]是以数学方法和计算机技术为主要工具，研究各种控制策略及控制系统的理论、方法和技术。它通过对网络拓扑结构及网络环境下先进控制理论与方法的研究，充分利用网络资源，实现从决策到控制的全过程优化。基于控制理论的复杂网络研究为生物信息学提供了新的思路，复杂网络是从真实系统中抽象得到的一类具有复杂拓扑结构的网络模型。现实世界的真实系统大多是非线性系统，但非线性关系很难被刻画，因此借鉴线性系统的动力学，可以研究非线性的复杂网络[12]。根据控制理论，包含 N 个节点的线性系统动力学方程如下：
[dxtdt=Axt+But ]（2）
其中，向量[xt]表示网络内的N个节点在t时刻的状态。[ut]表示外部控制的M个节点在[t]时刻的状态。若存在一个连续的输入[ut]，使得在有限的时间内，能使得线性时不变系统能够从任意的初始状态达到任意的终止状态，则称这个系统是可控的。
3.2 最小驱动节点集思想
Liu等[13]结合图论方法和结构可控性定理，提出了一个基于图模型的最大匹配算法求解最小驱动节点集的复杂网络可控性分析框架。該框架提出通过最小驱动节点集可以实现控制整个网络的目的，从理论上证明了最小驱动节点集就是最大匹配中的未匹配节点集。最小驱动节点集中节点数量越多，意味着复杂网络的可控性越低，反之亦然。此外，最小驱动节点集的节点数量由网络中的度决定，网络的度越大，最小驱动节点集中节点数量越少，复杂网络的可控性就越低。
3.3 可控性动态分类理论
可控性动态分类理论基于Kalman判定条件[14]和网络的拓扑结构判断整个网络是否满足可控性。核心思想是将网络看作矩阵，将矩阵中的所有元素设为0或已确定的常数。若一个网络是可控的，那么网络中边的权重不影响网络的可控性。Kalman判定条件如下：
设A为网络的连接矩阵，B为外部输入矩阵，那么定义矩阵C为：
[C=B，AB，A2B，…，AN?1B] （3） 当矩阵C为满秩时，即
[rankC=rankB，AB，A2B，…，AN?1B=N]时，网络是可控的。
Vin等人基于可控性动态分类理论分析了蛋白质互作网络的可控性[15]，并对Liu等人的研究进行了改进，将蛋白质相互作用网络中的每一个节点逐个删除，检验最小驱动节点集中节点数量的变化，若节点数量增加，则说明删除的节点为关键节点;若节点数量保持不变，则说明删除的节点为普通节点;若节点数量减少，则说明删除的节点为冗余节点。实验结果表明，关键节点代表的蛋白质是患病的主要因素。而普通节点和冗余节点这两种类型的蛋白质的代表性不足。
3.4 关键节点基因集构建流程
基于上述理论，借助Liu[16]的框架代码，结合readr包后导入结直肠癌的表达谱数据，编写R语言代码分析构建lncRNA-TF-gene调控网络，阈值设为P<0.05，分析后构建最小驱动节点集。基于可控性动态分类理论，移除最小驱动节点集的每个节点，依据移除后对调控网络的影响，将最小驱动节点集的节点进行分类，分为关键节点、普通节点和冗余节点。利用write.csv函数导出分类结果和关键节点基因，以便后续分析。
4 基于融合信息筛选关键基因
针对目前筛选关键基因方法单一的问题，本节研究将突变频率和网络拓扑分析方法得分进行融合，利用突变得分和网络拓扑分析方法得分的融合信息对lncRNA-TF-gene调控网络进行分析，得到潜在关键基因，以便后续研究。
4.1 突变得分
突变得分通常指每一个配子发生突变的概率[17]，即用一定数目配子中的突变配子数表示。在无性繁殖的细菌中，突变得分代表用每一个细胞世代中每个细菌发生突变的概率，即用一定数目的细菌在一次分裂过程中发生突变的次数表示。不同生物或同一生物个体的不同基因的自发突变率是不相同的。突变得分的估算方法因生物生殖方式的不同而不同。在有性生殖的生物中，突变率通常为每一个配子发生突变的概率，即用一定数目配子中的突变配子数表示。在无性繁殖的细菌中，突变率代表每一个细胞世代中每个细菌发生突变的概率，即用一定数目的细菌在一次分裂过程中发生突变的次数表示。
通过TCGA数据库下载结直肠癌对应突变数据集，基于突变的变体分类，仅选择具有功能的突变，例如splice_site，in_frame_del和frame_shift_del等，并計算突变频率。基因的变异越频繁，在排名列表中就越高。根据TCGA数据的突变频率计算调控网络中基因的突变得分，利用Python语言排序后整理，得到基因突变得分列表。
4.2 网络拓扑分析方法得分
网络拓扑统计参数反映网络的结构和动力学特性，对网络结构的分析一般包括以下几个方面：将交互关系复杂的系统建模为网络，然后计算网络拓扑参数来研究网络的拓扑结构特性，根据网络演化过程中网络拓扑参数的变化来分析网络演化对网络结构的影响[18]。
本研究通过3个网络拓扑方法对潜在mRNA关键基因进行联合分析，即degree、MNC和MCC网络拓扑分析方法[19]。degree代表节点的度，v为图中节点数目，MNC为最大邻域分量，其公式为：
[MNCv=VMCv] （4）
其中[V ]指网络中节点的总和，[MCv]指子图与子图间的最大连接分量。
MCC网络拓扑分析方法可以很好地分析节点数目较多的复杂网络。给定一个节点v，MCC的定义为：
[MCCv=C∈SvC?1！] （5）
其中，[C]为属于[Sv]的一个团，[Sv]为包含节点[v]的最大团的集合，[C?1]表示小于[C]的所有正整数的阶乘。
利用Cytoscape3.7.2软件的cytoHubba插件中的MNC、MCC和Degree网络拓扑分析方法分析lncRNA-TF-gene调控网络，根据网络拓扑分析方法得分进行排序，利用python语言整理后得到网络拓扑分析方法得分列表。
4.3 潜在关键基因筛选
若一个基因可以同时出现在突变得分列表和网络拓扑分析方法得分列表中，则认为该基因是潜在关键基因。因此，利用python语言的numpy程序包[20]分别统计突变得分列表和网络拓扑分析方法得分列表中排名前50名的基因，并通过R语言对两个列表下排名前50名的基因取交集，整理后得到潜在关键基因集。
4.4 关键基因筛选及验证
结合上文的关键节点基因集和潜在关键基因集，通过R语言找到两个基因集的共同基因，整理后命名为关键基因集，关键基因在结直肠癌的发病过程中起到了重要的调控作用。为了验证实验结果的有效性，利用python遍历文献中的实验数据，结合CGC数据库[21]中的数据对关键基因集中的基因进行验证。利用两方面数据对关键基因集进行验证，以便对预测算法的效果进行评估。
5 实验结果与分析
5.1 lncRNA-TF-gene调控网络
利用TCGA数据库共得到60477条基因表达谱数据，经背景校正、归一化处理、基因名转化后，整理为csv格式，以便后续分析。
利用LncMAP数据库共获取到102111条结直肠癌中lncRNA-TF-gene调控关系数据，整理后命名为lncRNA-TF-gene调控关系数据集A，利用python语言共收集到32011条结直肠癌中lncRNA-TF-gene调控关系数据，整理后命名为lncRNA-TF-gene调控关系数据集B。将调控关系数据集A与调控关系数据集B取交集后得到10521条调控关系数据，整理后命名为交集调控关系数据集，部分数据如表1所示。结果表明，lncRNA-TF-gene调控关系显著相关并且具有显著差异性，可作为调控网络进行研究。利用Cytoscape3.7.2软件构建lncRNA-TF-gene调控网络，并采用NetworkAnalyzer功能对调控网络进行分析。 5.2 關键节点基因集
基于最小驱动节点集思想，对lncRNA-TF-gene调控网络进行筛选分析，共得到1021个驱动节点基因，整理后命名为最小驱动节点集。依据可控性动态分类理论，对最小驱动节点集进行分类，分类后的结果如表2所示，扇形图如图2所示。可控性动态分类筛选出的14%的关键节点基因，可作为研究潜在的关键基因。关键节点基因集的部分结果如表3所示，其中，K为lncRNA-TF-gene调控中节点的度，Kin为节点的入度，Kout为节点的出度，TypeI为冗余和普通节点基因，TypeII为关键基因，0表示否，1表示是。
5.3 潜在关键基因集
利用R语言下载TCGA数据库中结直肠癌的临床数据，保存为RData格式。编写R语言代码计算突变得分，利用python语言的numpy程序包整理后得到突变得分列表，得分前10名的基因如表4所示。利用Cytoscape3.7.2软件的cytoHubba插件的MNC、MCC和Degree网络拓扑指标联合分析lncRNA-TF-gene调控网络，利用Python语言的numpy程序包整理后得到网络拓扑分析方法得分列表，得分前10名的基因如表5所示。
利用Python语言的numpy程序包分别统计突变得分列表和网络拓扑分析方法得分列表中排名前50名的基因，对两个列表中排名前50名的基因集取交集，得到36个潜在关键基因。
5.4 关键基因集及验证
结合上文的关键节点基因集和潜在关键基因集，对两个基因集取交集，共得到17个关键基因（CASP3、JAK1、CDK2、CCL5、FN1、PCNA、MAPK8、SOD1、PPARG、VEGFA、TGFB1、CDKN1A、AKT1、TP53、BCL2、MAPK3和MAPK1）。为了验证实验的有效性，通过CGC网站和大量文献对关键基因进行验证，发现17个关键基因中有13个已经通过CGC网站和相关文献得到证实，与结直肠癌有密切关系。这表明其余4个关键基因有很大可能与结直肠癌也密切相关。这些关键基因对未来结直肠癌机理研究和开发靶向药物均具有重要参考价值。
6 结论
本文提出了一种基于控制理论的结直肠癌关键基因预测算法，本算法基于突变得分和网络拓扑分析方法得分的融合信息进行关键基因筛选。实验结果证实，关键基因的预测准确性较以往的预测算法得到了显著提高。本算法不仅可以应用在结直肠癌的关键基因预测，还可以广泛应用在各类癌症的关键基因预测中，为今后发现癌症的关键基因提供新的一种思路和选择。
参考文献：
[1] Cancho V G，Bazán J L，Dey D K.A new class of regression model for a bounded response with application in the study of the incidence rate of colorectal cancer[J].Statistical Methods in Medical Research，2020，29（7）：2015-2033.
[2] Maduni? K，Zhang T，Mayboroda O A，et al.Colorectal cancer cell lines show striking diversity of their O-glycome reflecting the cellular differentiation phenotype[J].Cellular and Molecular Life Sciences，2021，78（1）：337-350.
[3] Baratti D，Kusamura S，Pietrantonio F，et al.Progress in treatments for colorectal cancer peritoneal metastases during the years 2010-2015.A systematic review[J].Critical Reviews in Oncology/Hematology，2016，100：209-222.
[4] Reimand J，Bader G D.Systematic analysis of somatic mutations in phosphorylation signaling predicts novel cancer drivers[J].Molecular Systems Biology，2013，9（1）：637.
[5] Hou J P，Ma J.DawnRank：discovering personalized driver genes in cancer[J].Genome Med，2014，6（7）：56.
[6] Network C G A R，Weinstein J N，Collisson E A，et al.The cancer genome atlas pan-cancer analysis project[J].Nature Genetics，2013，45（10）：1113-1120.
[7] Flicek P，Amode M R，Barrell D，et al.Ensembl 2012[J].Nucleic Acids Research，2012，40（database issue）：D84-D90.
[8] Ji Y，Gu Y Q，Hong S，et al.Comprehensive analysis of lncRNA-TF crosstalks and identification of prognostic regulatory feedback loops of glioblastoma using lncRNA/TF-mediated ceRNA network[J].Journal of Cellular Biochemistry，2020，121（1）：755-767. [9] Li Y S，Li L L，Wang Z S，et al.LncMAP：Pan-cancer atlas of long noncoding RNA-mediated transcriptional network perturbations[J].Nucleic Acids Research，2018，46（3）：1113-1123.
[10] Calinski T，Steel R G D，Torrie J H.Principles and procedures of statistics：a biometrical approach[J].Biometrics，1981，37（4）：859.
[11] Narendra K S.Neural networks for control theory and practice[J].Proceedings of the IEEE，1996，84（10）：1385-1406.
[12] 王婉宁.基于复杂网络结构可控性的致病基因识别方法研究[D].西安：西安理工大学，2019.
[13] McDowell J J.Matching theory in natural human environments[J].The Behavior Analyst，1988，11（2）：95-109.
[14]WELCH G， BISHOP G. An Introduction to the Kalman Filter[M]. University of North Carolina at Chapel Hill， 1995.
[15] Vinayagam A，Gibson T E，Lee H J，et al.Controllability analysis of the directed human protein interaction network identifies disease genes and drug targets[J].Proceedings of the National Academy of Sciences of the United States of America，2016，113（18）：4976-4981.
[16] Liu Y Y，Slotine J J，Barabási A L.Controllability of complex networks[J].Nature，2011，473（7346）：167-173.
[17] Devarakonda S，Rotolo F，Tsao M S，et al.Tumor mutation burden as a biomarker in resected non-small-cell lung cancer[J].Journal of Clinical Oncology，2018，36（30）：2995-3006.
[18] 周庚.复杂网络节点中心性度量算法的研究及应用[D].兰州：兰州理工大學，2020.
[19] Chin C H，Chen S H，Wu H H，et al.cytoHubba：identifying hub objects and sub-networks from complex interactome[J].BMC Systems Biology，2014，8（Suppl 4）：S11.
[20] van der Walt S，Colbert S C，Varoquaux G.The NumPy array：a structure for efficient numerical computation[J].Computing in Science & Engineering，2011，13（2）：22-30.
[21] Sondka Z，Bamford S，Cole C G，et al.The COSMIC Cancer Gene Census：describing genetic dysfunction across all human cancers[J].Nature Reviews Cancer，2018，18（11）：696-705.
【通联编辑：唐一东】