基于熵和Fold Change算法的抗前列腺癌药物重定位分析

Prostate Cancer Drugs Repositioning Based on Entropy and Fold Change Algorithm

摘要: 本文针对抗前列腺癌药物重定位问题，首先将R软件中的信息熵和差异表达倍数Fold Change算法结合起来，分析前列腺癌的基因表达数据，并提取含有大量致病基因信息的前列腺癌特征基因，通过FC算法得到2788条上调基因和2222条下调基因；然后对这些特征基因各筛选出200条具有显著信息的特征基因作为检索标签进行cmap (https://portals.broadinstitute.org/cmap/)分析；最后通过cmap负相关分数降序排列，分析确定几种可能对治疗前列腺癌有效的药物。与传统的提取特征基因方法相比，信息熵和Fold Change算法相结合能减少计算步骤和计算量，成本低、耗时少、准确度高。该方法有效可行。

Abstract: Drug repositioning research not only greatly reduces the cycle of new drug development, but also reduces the economic cost. Prostate cancer is a cancer with high mortality in the world. In recent years, the research on single gene is not enough to analyze the complex pathogenesis of cancer. In this paper, information entropy and Fold Change (FC) algorithm in R software were applied to an-alyze the gene expression data of prostate cancer and extract the characteristic genes of prostate cancer with a large number of pathogenic gene information. Then several drugs for prostate cancer were identified by cmap analysis. The gene expression data used in this paper were from TCGA, and 2788 up-regulated genes and 2222 down-regulated genes were obtained through FC algorithm. And characteristics of these genes screened 200 select genes with significant information as retrieval label cmap (https://portals.broadinstitute.org/cmap/) analysis; finally, the negative correlation score of cmap was ranked in descending order, and several possible effective drugs for the treatment of prostate cancer were analyzed and determined. Compared with the traditional method of extracting characteristic genes, the combination of information entropy and Fold Change algorithm can reduce the computational steps and amount of computation, with low cost, low time consumption and high accuracy. This method is effective and feasible.

1. 引言

前列腺癌是发生在男性前列腺组织中的恶性肿瘤，是前列腺腺泡细胞异常无序生长的结果。前列腺癌发病率的高低与地理和种族的差异性有关。在欧美一些发达国家和地区，它是男性最常见的恶性肿瘤，死亡率排在各种癌症的第二位；在亚洲，虽然发病率低于其他西方国家，但是近几年也呈迅速上升趋势。临床上前期主要采用雄激素剥夺疗法(ADT)治疗前列腺癌，然而几乎所有患者最终都会发展为致命性的去势抵抗型前列腺癌(CRPC)。虽然FDA (美国食品药品管理局)批准的第二代抗雄激素药物如Enzalutamide (恩杂鲁胺)和Abiraterone (阿比特龙)等对缓解疾病进展具有一定的功效，但患者很快就会出现临床耐药。因此，临床上迫切需要治疗前列腺癌的特效药。

鉴于国内现有的医疗水平，针对前列腺癌仅能通过常规手术治疗、内分泌及化学药物疗法来提高患者的生活质量，但提高患者的生存期依旧是一个难题。目前，分子靶向治疗已成为肿瘤治疗的研究热点，为前列腺癌的治疗也提供了新的思路和方向。利用基因表达谱等组学技术发现抗前列腺癌的药物靶标可作为一个重要手段。但新药开发是一个耗时费力的高风险过程，充分发掘已有药物的新用途，对药物进行重定位，备受生物医药产业和学者们的青睐 [1] [2] [3] 。

药物重定位又称老药新用，指对曾经用于临床的药物新适应症的发现、确认和应用。包括对处于临床研究阶段或已批准上市的药物进行重定位、重定用途、重评价和重新定位治疗方向等 [4] 。推动一个新药物上市通常需要13~15年，其成本平均需要20~30亿美元，且处于上升趋势。如果对已有药物进行研究，一旦它们拥有不同的医疗用途，这将是一个巨大的未开发资源。“药物重定位”可以跳过临床I期，相比于新药物大大地缩减研究成本和投入时间。到目前为止，从已知的药物中发现新的适应症，成功重定位的药物已经有100多种。如何从已知药物中发现对于前列腺癌有治疗效果的药物是本文探讨的问题。

信息熵是信息论中用于度量信息量的一个概念。对于基因来说，它的信息熵越大，所包含的信息量就越大，具有的生物学意义就越大，计算出基因的信息熵再求其互信息值来提取特征基因是常用的特征基因提取方法。但是用互信息方法提取上调基因和下调基因的过程中，对阈值的选取要求比较高，阈值选取不当，会导致上调基因和下调基因选取不当。

Fold Change (简称FC)是差异表达倍数，对于基因表达数据来说，计算其logFC值来筛选特征基因，如果logFC值为正则为上调基因；如果为负则为下调基因。若只用FC值来提取特征基因，计算所有基因的FC值再取对数则计算量大且时间成本高。但如果将两个方法结合，先利用信息熵提取出一些致病信息量大的特征基因，再通过计算它们的FC值再筛选出上调基因和下调基因，不仅省去选取阈值的步骤，而且只需计算信息熵大的基因的FC值，大大缩减计算步骤和计算量，耗时低，结果也更精确。如何在计算信息熵的基础上结合FC值提取出前列腺癌的特征表达基因是本文探讨的问题。本文首先从TCGA数据库中获取前列腺癌与癌旁的基因表达数据，利用R软件将数据进行标准化预处理；然后利用互信息算法算出每条基因的信息熵，从大到小排列，将与前列腺肿瘤密切相关即信息熵值最大的前5000个特征基因筛选出来，再利用计算其logFC值，筛选出上调基因和下调基因；最后通过cmap数据库 [5] 分析，检索出具有与肿瘤基因相反的基因标签的药物。gliclazide (格列齐特)作为一种治疗中型糖尿病的药物，与胰岛素合用治疗胰岛素依赖型糖尿病，可减少胰岛素用量，经分析比对得到的负相关分值最高，表明对于前列腺癌可能具有较好的治疗效果。Luteolin (木犀草素)、phenoxybenzamine (竹林胺)、naloxone (盐酸纳洛酮)等化合物也具有较高的负相关分值，表明极可能对前列腺癌有治疗效果。

2. 数据和方法

2.1. 基因表达数据

TCGA是美国国家癌症研究所(National Cancer Institute)和美国人类基因组研究所(National Human Genome Research Institute)共同监督的一个项目，旨在应用高通量的基因组分析技术，以帮助人们对癌症有个更好的认知，从而提高对于癌症的预防、诊断和治疗能力。作为目前最大的癌症基因信息数据库，TCGA数据库的全面不止体现在众多的癌型上，还体现在多组学数据，主要收录基因表达数据、mRNA表达数据、拷贝数变异、DNA甲基化、SNP等等，是癌症研究者十分重要的数据来源。TCGA数据库最大的优势就是包含各种人类癌症(包括亚型在内的肿瘤)的临床数据，针对每个癌型都有规范的，大样本量的临床数据，比在GEO数据库中获取的基因表达数据更丰富。TCGA覆盖了人体六十个组织或器官的三十八种癌型以及亚型，四十个projects，近四万个患者，而且在不断更新中，更新速度也比其他数据库快，对于癌症的研究有十分丰富的资源。

本文的前列腺基因表达数据来自TCGA数据库，在TCGA数据库中，样品名称格式为TCGA-19-2619-10A，其癌症样本和正常样本的分类标准与最后的10A有关，编号从01~09是癌症样本；10~29是癌旁样本。按照此分类标准利用R软件进行分类筛选，获得前列腺癌与癌旁的基因表达数据，其中包括488个患病样本和12个健康样本，包括60,482条基因(https://cancergenome.nih.gov/)。

2.2. 特征基因提取

通过TCGA数据库得到的基因表达数据为原始数据，要首先对数据进行预处理，利用标准分数(z-score)将得到的前列腺基因表达数据进行标准化处理，得到标准化的数据，其标准化公式为：

$z_{ij}=\frac{g_{ij}-mean\left(g_i\right)}{std\left(g_i\right)}$ (1)

其中g_ij表示基因i在样本j中的表达值， $mean\left(g_i\right)$ 表示基因i在所有样本中的平均值， $std\left(g_i\right)$ 表示基因i在所有样本中的标准差 [6] 。

传统的特征基因提取方法通常把所有基因FC (差异表达倍数)的值计算出来，根据FC值的正负再筛选出上调基因和下调基因，FC值为正则为上调基因，FC值为负则为下调基因。通过此方法提取特征基因的计算量大，要把六万多条的基因每个都计算一边，耗时多，最后得到的上调基因和下调基因的数量过多，不具有针对性，再拿到cmap数据库中去检索，得到的结果也不精确，影响最后的结论。也有研究者采用只求出基因的信息熵和互信息值的方法，对于复杂的基因关系，熵和互信息的方法能有效抓住基因与基因之间的关联性，利用熵和互信息值能有效的提取出复杂疾病的致病基因 [7] [8] [9] 。熵是对不确

定性的度量，在信息论中，熵是用来衡量一个随机变量出现的期望值。设基因变量 $X=\left[x_1,x_2,\cdots ,x_n\right]$ 是一

个基因表达模式，基因变量X的熵表示该模式所包含的信息量公式为：

$H\left(X\right)=-{\displaystyle \underset{x\in X}{\sum }p\left(x\right){\log }_2^{p\left(x\right)}}$ (2)

利用公式将得到的数据中60,482条基因的信息熵计算出来，取熵值较大即含有较多致病信息的部分基因作为特征基因，再计算这些特征基因之间的互信息值，不仅编程过程复杂且计算量大，而且采取此方法计算得到的互信息矩阵还需要再选取适当的阈值来筛选出可以拿到cmap数据库中去比对的基因。而阈值的选取是一个复杂的过程，选的阈值过高或过低都会影响最后的检索结果，所以选取不当会影响结果的精确性，使得检索到的药物不具有代表性和针对性。

针对传统特征提取方法耗时多，计算量大，步骤复杂的情况，本文采取信息熵和差异表达倍数相结合的方法来提取特征基因。首先计算出每条基因的信息熵，然后根据信息熵的值降序排列，取信息熵大即生物学意义高，含有较多致病信息的前5000个基因作为特征基因，然后再计算这5000个特征基因的基因差异表达倍数值(Fold Change)。前列腺癌的基因表达数据包含60,482条基因，如果用上述传统的方法，需要计算60,482个基因的FC值，采取本方法只需计算信息熵大的5000个基因的FC值，计算量小，节省了计算时间，FC值的计算公式为：

$\log F{C}_i={\log }_2\left(\frac{\frac{1}{s}{\displaystyle \underset{k\in T}{\sum }{z}_{ik}}}{\frac{1}{M}{\displaystyle \underset{k\in N}{\sum }{z}_{ik}}}\right)$ (3)

其中z_ik表示基因i在样本k中的标准化表达值，T是癌症病人样本的集合，N是正常样本的集合，S和M分别为癌症病人样本和正常样本的个数(本文中S = 488，N = 12)，最终将求得的 $\log FC$ 值以0为分界点，把 $\log FC$ 值大于0的所有基因作为上调基因集合，小于0的所有基因筛选为下调基因集合。

2.3. Cmap分析

Connectivity Map (简称cmap)是一个基因表达谱数据库，其利用小分子药物处理人类细胞后的基因表现差异，建立一个小分子药物、基因表现与疾病相互关连的生物应用数据库。以基因表达谱为所建立之基因、疾病与药物的关联性，协助学者们在药物开发领域上，快速利用基因表达谱的数据比对出与疾病高关联性的药物、推论大部分药物分子的主要化学结构，并能够归纳出药物分子可能作用的机制方向。近年来的研究趋势也显示出利用cmap基因表达谱数据库应用在疾病治疗与药物开发的领域上，可提供越来越精确的方向。目前cmap第二版已经发展成收录了1309种药物表达谱，有超过7000笔的基因表达谱资料。每一种药物分子会以不同浓度(10 nm、100 nm、1 µM与10 µM)处理在不同的细胞株(breast、prostate、leukemia与melanoma cell line)，并处理不同的时间点(6与12小时)。理论上讲，与疾病和药物相关的任何基因表达数据都可以在cmap数据库中进行高效率的查询比对，从数据库揭示药物、基因和疾病三者之间潜在的联系 [10] 。

通过R软件筛选出前5000个信息熵大的特征基因，并利用公式计算其FC值，根据正负将特征基因分成上调基因和下调基因，最后得到上调基因2778个和下调基因2222个。接下来在上调基因中选取FC值较大和下调基因中FC值较小的部分基因去检索cmap数据库。即将得到的上调基因中降序排列，提取前200个 $\log FC$ 值最大的上调基因作为检索标签，将下调基因升序排列，提取前200个 $\log FC$ 值最小的下调基因作为检索标签，将它们存为.grp文件，检索cmap数据库 [11] 。将前列腺癌基因表达标签与药物处理基因标签进行统计比较 [12] 。依据表达谱的相似性给每个前列腺癌-药物配对计算一个分值，如果分值为负数，则表明这种药物与癌症基因有相反的基因标签，对癌症基因具有抑制作用，即可能对前列腺癌具有较好的治疗效果 [13] [14] 。所以在检索的过程中，删除试验次数较少的药物( $n<4$ )，关注药物得分Mean分值为负值的药物 [15] 。

3. 结果与展望

3.1. 结果分析

Cmap的分析结果如表1所示，根据表中数据可以看出负相关分值最高的是gliclazide (格列齐特)，分值为−0.69，格列齐特是一种治疗中型糖尿病的药物。经分析和比对，得到的负相关分值最高，表明其对于前列腺癌可能具有较好的治疗效果；表中还可以看出排在后面的Luteolin (木犀草素)、phenoxybenzamine (竹林胺)、alpha-estradiol ( $\alpha$ -雌二醇)，naloxone (盐酸纳洛酮)等化合物也具有较高的负相关分值，表明极可能对前列腺癌有治疗效果。其中表中的alpha-estradiol ( $\alpha$ -雌二醇，别名雌二醇)是经皮肤吸收的雌激素治疗剂，目前已经被用来治疗晚期前列腺癌。通过表1可以看出排在它上面的药物最后的检索分值的负相关性均高于它，所以这几种药物很可能对治疗前列腺癌有治疗效果。

Table 1. Candidate prostate cancer drugs screened from the connectivity map database

表1. Connectivity map数据库筛选出的候选抗前列腺癌药物

注：Mean表示药物检索得分值，n为药物在cmap数据库中重复试验的次数，enrichment为前列腺癌症基因标签与药物基因标签相似的聚合度。

根据在cmap数据库中的结果比对显示，gliclazide (格列齐特)极有可能对前列腺癌有治疗作用，作为治疗糖尿病的降糖药物已经成为二型糖尿病一线治疗的标准，并在糖尿病临床治疗中积累了大量成功经验，具有很好的耐受性，多数患者依从性良好。如果有实例进一步验证该药对前列腺癌有治疗作用，那么它作为已经成功上市的药物，则节省了研发新药的时间，对其稍加改善，就可以投入到前列腺癌的治疗中。

3.2. 结论与展望

在提取特征基因过程中，通过表2可以看出re表示在求得信息熵之后，取前5000个熵最大的基因，然后根据公式求得的每个基因的logFC值，此方法只需分别计算这5000个基因的差异表达倍数即可。简单，易操作，程序编写也容易。但是如果用信息熵和互信息的方法来获取特征基因，如表3所示，在同样求得信息熵之后，取前5000个熵最大的基因，根据互信息的定义，要计算出这5000个基因两两之间的互信息值，最后得到的5000 × 5000的矩阵，编程步骤繁琐，计算量大，而且还要对求得的矩阵选取合适的阈值来筛选特征基因，阈值的选取既不能太高，也不能太低，选取不当会导致特征基因不具有代表性，影响最后到cmap数据库检索的结果。

Table 2. The code of using R software to calculate FC value

表2. 利用R软件求FC值的代码

Table 3. The code of mutual information value using R software

表3. 利用R软件求互信息值的代码

本文通过信息熵和Fold Change相结合的新方法提取前列腺癌中的特征基因，利用cmap数据库将基因与药物进行比对打分，最后得到与治疗前列腺癌有关的药物gliclazide (格列齐特)、Luteolin (木犀草素)、phenoxybenzamine (竹林胺)等。

与传统的特征基因提取方法比较，本文采取的方法耗时短，在利用R软件计算过程中，节省了时间，基于信息熵和Fold Change算法提取特征基因，比只利用信息熵和互信息与只利用FC值的方法提取特征基因用时少、成本低和准确性高，但数据分析结果还需要临床试验的进一步验证，希望有条件的实验室能完成这一工作。为药物重定位提供了新的途径，推动生物医药产业的发展。

基金项目

国家自然科学基金资助项目(61703290)，辽宁省科技厅自然科学基金资助项目(20180550133)；辽宁省教育厅科学技，术研究项目(LQN201710)。