基于基因共表达网络分析的扩张和限制型心肌病分子特征比较

时间：2024-08-31

黎桑王光斌李政文曾伟饶妮妮,3*

1(电子科技大学信息生物中心，成都 611731)2(电子科技大学生命科学与技术学院，成都 610054)3(电子科技大学广东电子信息工程研究院，广东东莞 523808)

引言

心血管疾病一直都是威胁全人类健康的复杂疾病，而心肌病在心血管疾病的发病中占据了重要的位置，并且有着逐年上升的趋势。其中，原发性心肌病的病因尚未明确，笔者所研究的DCM和RCM都属于原发性心肌病[1]。

DCM的发病率占心肌病的80%以上，并且其中有58.3%死于该病[2]。因此，识别DCM的生物标志物、探究DCM的发生发展机制以及筛选DCM的治疗靶点等问题都十分重要。迄今为止，已找出超过40个与DCM致病相关的基因，如表1所示[3]。但是，这些基因对于理解DCM的发生发展机制还远远不够。

表1 与DCM致病相关的基因[3]Tab.1 Genes associated with DCM disease[3]

RCM虽然发病率不高，仅占心肌病发病率的4.5%，但是它却是最为严重的心肌病，其预后性较差，并且在临床上常常容易被误诊。相关文献[4]报道，RCM可能被误诊为DCM、病毒性心肌炎等。在DCM和晚期肥厚型心肌病(hypertrophic cardiomyopathy，HCM)中，也会有在RCM中出现的舒张功能障碍的情况。可见，RCM的精准诊断一直是临床上存在的一大难题。因此，通过生物信息学的方法筛选出RCM的生物标志物，将能够在一定程度上改善这一问题。通过构建RCM相关的差异表达基因的共表达网络，可以加深对RCM的发生发展机制的理解。

目前，国内外对于RCM致病基因的文献报道较少，但TNNI3[5]、TNNT2[6]、ACTC[6]、MYH7[7]、DES[8]、TPM1[9]、MYL3[9]、MYL2[9]等都被证明与RCM相关。

为此，本研究首先识别两类心肌病相关的差异表达基因；接着，通过构建差异表达基因的共表达网络，分别筛选与两类心肌病相关的基因标志物，并对其功能进行分析。在此基础上，对两类心肌病进行比较，希望在分子层面上加以区分，减少误诊的几率，为精准医疗提供一些参考。

1 材料与方法

1.1 RNA-Seq数据

本研究数据来自于EMBL-EBI数据库(NCBI上编号为GEO-GSE71613)，是全基因高通量转录组测序(RNA-Seq)的数据，为Fastq格式。在DCM或RCM患者做心脏器官移植手术时，取出病变组织并进行活体检查，收集RNA样本，通过高通量转录组测序得到该数据；共有8组样本，其中4个正常样本(有8个Fastq文件)、2个DCM样本和2个RCM样本。数据下载网址为http://www.ebi.ac.uk/arrayexpress/experiments/E-GEOD-71613/。

1.2 DCM/RCM相关差异表达基因的筛选

采用标准化的RNA-Seq数据处理方法来筛选DCM/RCM相关的差异表达基因，包括质量控制(FastQC)、参考序列比对(Tophat)、计数数据的转换(HTSeq)、归一化(EDASeq)和差异表达(DESeq)等步骤。

本研究采用了FastQC进行质量控制，去除一些重复和低质量的reads。在参考序列比对中，基因参考序列是hg38(UCSC数据库)。在Linux系统下运行tophat调用bowtie2，就能完成参考序列比对，得到格式为SAM或者BAM的文件。计数数据的转换由以下两步完成：先将Map后的BAM文件通过Samtools软件[10]转化成排序过的SAM文件，再将排序过的SAM文件经过HTSeq软件[11]进行定量分析得到计数数据。归一化是用来降低因测序的深度不同所产生的影响，以及潜在的其他重复性实验所带来的影响，采用了R语言的软件包EDASeq[12]来完成。最后，采用DESeq2软件包进行差异表达基因的筛选。在众多的差异表达分析方法中，DESeq2最适合于样本数较少的数据[13]。因此，本研究选用DESeq2。

1.3 差异表达基因的共表达网络构建方法及生物属性分析

基因共表达网络的构建是利用R语言软件包WGCNA实现的，其中选取皮尔逊相关系数作为基因共表达网络的参考相关系数。假设两个基因的表达谱变量分别为X和Y，这两者之间的皮尔逊相关系数的计算公式如下：

(1)

式中：E表示求数学期望；cov表示求协方差；ρX,Y为相关系数，其值在(-1,1)之间。

用Matlab和R语言的WGCNA软件包提取共表达网络的关系对，将所得的关系对输入软件Pajek中，绘制出差异表达基因的共表达网络图。

进一步，需要验证差异表达基因的共表达网络属于生物网络，才能进行后续的生物功能分析。生物网路应该具备3个属性：一是小世界的特性。构建与差异表达基因共表达网络同规模和平均度数相同的随机网络，如果其特征路长L1和聚集系数C1与差异表达基因共表达网络的两个对应参数L和C相比，满足L≥L1并且C≫C1，则所构建的共表达网络具备小世界的特性[15]。二是无尺度分布的特性。三是生物网络具有一定的鲁棒性和适应性。如果网络的度分布符合幂率分布，则具有第二、三两种属性。

1.4 基于差异表达基因的共表达网络识别疾病相关重要基因

重要节点在很大程度上涵盖了网络的重要信息，根据相关文献[16]，取度排序(由大到小)前5%的差异表达基因作为网络中的重要节点(即重要基因)。

2 结果

2.1 DCM和RCM相关的差异表达基因

对于DCM，筛选出451个差异表达基因，其中有上调基因218个、下调基因233个；对于RCM，筛选出1 326个差异表达基因，其中有上调基因594个、下调基因732个。

2.2 构建的共表达网络及生物属性分析

2.2.1DCM共表达网络及生物属性

根据本文第1.3节中所介绍的皮尔逊相关网络的硬阈值的选取条件[14]，取硬阈值为0.98，得到用DCM相关差异表达基因构建的共表达网络，其中网络的特征参数为LDCM=6.737 13，CDCM=0.590 259 02；随后，构建了与DCM的共表达网络同规模及平均度数相同的随机网络，其特征路长和聚集系数分别为L1=2.924 35和C1=0.020 325 89，满足LDCM≥L1和CDCM≫C1。因此，所构建的DCM共表达网络具备小世界的特性。DCM共表达网络的度分布如图1所示，符合幂律分布，因此所构建的共表达网络具有一定的鲁棒性和适应性。由此可知，用DCM相关差异表达基因构建的共表达网络属于生物学网络。

图1 DCM共表达网络的度分布。(a)DCM共表达网络中各节点的度的大小分布；(b)DCM共表达网络中各节点的度的概率分布Fig.1 Degree distribution of DCM co-expression network. (a) Degree number of each node in DCM co-expression network; (b) Probability distribution of degree of each node in DCM co-expression network

2.2.2RCM共表达网络及生物属性

同理，取硬阈值为0.98，得到用RCM相关差异表达基因构建的共表达网络，其中网络的特征参数为LRCM=6.737 13，CRCM=0.590 259 02；同理，构造了与RCM共表达网络同规模及平均度数相同的随机网络，得到随机网络的特征路长L2和聚集系数C2为L2=2.743 34，C2=0.015 283 91，满足LRCM≥L2和CRCM≫C2。因此，RCM共表达网络具备小世界的特性。RCM共表达网络的度分布如图2所示，符合幂律分布，该网络具有无尺度分布特性和一定的鲁棒性及适应性。由此可知，用RCM相关差异表达基因构建的共表达网络属于生物学网络。

图2 RCM共表达网络的度分布。(a)RCM共表达网络中各节点的度的大小分布；(b)RCM共表达网络中各节点的度的概率分布Fig.2 Degree distribution of RCM co-expression network. (a) Degree number of each node in RCM co-expression network; (b) Probability distribution of degree of each node in RCM co-expression network

2.3 识别的DCM/RCM相关重要基因

2.3.1DCM的重要节点基因

筛选了与DCM相关的21个重要基因(基因标志物)，如表2所示。其中，5个基因有文献支持，与DCM的发生发展密切相关，分别为NRG1[17-19]、MAB21L2[20-21]、POSTN[22]、IGFBP3[23]和NPR3[24]。由此可知，本研究的部分预测结果得到了现有生物证据的支持，提示这些预测结果的可靠性。

表2DCM相关的21个重要基因及其生物功能

Tab.221DCM-relatedkeygenesandtheirbiologicalfunctions

序号GENE NAME基因所对应的生物功能1SYT12调解钙离子所依赖的突触传递运输的膜2COMP 编码细胞外基质的非胶原蛋白;钙离子结合和蛋白质绑定3FOSB调节细胞增殖、分化和转换4NRG1 编码出膜糖蛋白,能够介导细胞间信号传导;在器官系统的生长发育中起作用5ART4 NAD(P)+蛋白精氨酸ADP-核糖基转移酶的活性;与蛋白质代谢通路相关6MAB21L2 TGF—β信号通路(转化生长因子信号通路)的下游靶向7FBXL16 与免疫系统通路、 I类MHC介导的抗原加工和表达通路相关8PLXDC2受体活性相关9SMYD2与基因表达与染色质构型通路相关10SNAP25 基因产物是突触前等离子体膜蛋白,参与神经递质释放的调节11CP 编码结合大部分铜等离子体的金属蛋白,参与铁(II)转铁蛋白铁(III)转铁蛋白的过氧化反应12COL8A1 基因产物是一个短链胶原蛋白和角膜内皮细胞的基底膜的主要组成部分13IGFN1蛋白编码基因14NPR3 编码3个利钠肽受体之一;利钠肽受体是一种调节血容量和压力、肺动脉高压、心脏功能以及一些代谢和生长过程的小肽15POSTN 编码分泌细胞外基质蛋白,具有组织发展和再生的功能,包括伤口愈合和心肌梗死后心室重塑16DIRAS3生长抑制相关;肿瘤抑制功能17IGFBP3 胰岛素样生长因子结合蛋白,促进细胞的分化、增殖以及个体的生长发育18RGS1 通过绑定激活减弱G的蛋白质信号活动,提高从 GTP转化到 GDP的速率19B2M与分解代谢过度的低蛋白血症相关20BMX 编码蛋白与Stat通路、调节的几种类型癌细胞的分化和致瘤性的通路相关21OMD与新陈代谢、艾滋病毒生命周期通路相关

2.3.2RCM的重要节点基因

筛选了与RCM相关的65个重要基因(基因标志物)，如表3所示。其中，2个基因有文献支持，与RCM的发生发展密切相关，分别为TNNI3[5]和DES[8]。由此可知，本研究的部分预测结果得到了现有生物证据的支持，提示这些预测结果的可靠性。

表3RCM相关的65个重要节点基因及其生物功能

Tab.365RCM-relatedkeygenesandtheirbiologicalfunctions

序号GENE NAME基因所对应的生物功能1PON2 编码蛋白在人体组织、膜结合中广泛表达,可以作为细胞抗氧化剂,保护细胞免受氧化应激;该基因缺失与血管类疾病、糖尿病等变型疾病相关2FN3KRP 转移酶活性、激酶活性相关,蛋白质代谢、γ羧化作用等通路相关3TRAK2 受体结合和GABA受体结合相关;与新陈代谢通路相关4KMT2A 与转录因子活性、特异的DNA序列结合、同类蛋白结合相关5F2R 属于G蛋白耦合的受体家族,是凝血因子II受体,能调节血栓性反应6ARHGAP33磷酸肌醇绑定的功能位点,与胞内运输相关7LRRN3蛋白质复杂的绑定相关8FIBIN蛋白编码基因9TNNI3 在心肌组织表达,该基因突变会导致家族性的RCM10CELF2 调控pre-mRNA可变的剪接,可能与信使RNA编校和翻译相关11SHPK磷酸转移酶活性相关12CREB3L1 转录因子活性,特异的DNA序列结合,染色质绑定13SPECC1 位于细胞核,在睾丸和一些癌症细胞系中高表达14ANKZF1蛋白编码基因15EPHX2编码蛋白位于胞质和过氧化物酶体 16PTCH1 肿瘤抑制基因,基因的突变与基底细胞痣综合症、食管鳞状细胞癌、毛上皮瘤、膀胱移行细胞癌以及前脑无裂畸形相关17ADK 催化ATP到腺苷的转移,腺苷对心血管有很大的影响18ZNF599核酸绑定有关;与基因的表达通路相关19UBN2转录因子活性、特异的DNA序列结合20ADAMTS17促进乳腺癌细胞的生长和存活21UBXN7转录因子绑定和泛素绑定22GAB1 分支管腺增生介质;在细胞的生长反应、运输和凋亡中起重要作用23TRIB2 与Wnt / Hedgehog / Notch、DNA损伤等通路相关24EI24 抑制细胞生长,通过细胞凋亡蛋白酶9和线粒体通路来诱导凋亡细胞死亡25PATZ1 染色质建模和转录调控相关,与雄激素受体信号通路、雄激素受体的活性通路相关26CKM 转移酶活性、激酶活性,与新陈代谢等通路相关27DES 肌间线蛋白基因,与同类蛋白结合、细胞骨架的结构组成相关 28CYS1细胞器生物起源和保持的信号通路有关29PDGFD 与生长因子活性、血小板源生长因子受体结合相关30XAF1与干扰素γ信号通路、免疫系统通路相关31AKAP6 离子通道绑定和绑定蛋白激酶相关,与依赖cAMP的蛋白激酶激活、DAG/IP3信号通路相关

续表

3 讨论

3.1 基于生物信息学工具的DCM相关基因标志物的功能分析

3.1.1DCM相关基因标志物的GO分析

将DCM相关的21个重要节点基因用GO分析的系统工具PANTHER进行富集分析，得出以下结果：在生物学过程中，重要节点基因主要富集在细胞生理(54.5%)、代谢(31.8%)、多细胞生物(31.8%)、组织或生物起源细胞组件(27.3%)等过程中；在分子功能中，主要富集在蛋白结合(40.0%)、催化活性(27.3%)等方面；在细胞组成成分上，主要富集在细胞外基质(extracellular matrix，ECM,占27.3%)、膜(18.2%)、细胞组分(18.2%)等成分上。

3.1.2DCM相关基因标志物的信号通路分析

同样，由PANTHER可得出DCM相关的信号通路分析结果，包括整合素信号通路、p53通路、5-HT受体信号通路、β肾上腺素能受体信号通路、G蛋白的受体激活通路等。

Xiong等学者通过将大鼠的MDM4基因剔除，发现它最终发展成DCM，平均存活了153天；并且发现大鼠的MDM4基因突变是通过p53信号通路导致心肌细胞的凋亡，因此推测出p53通路可能在DCM的发生发展中有着重要的作用，也提示p53通路可能与DCM的发生发展机制相关[25]。如图3所示，p53基因在细胞生长周期中起着负调节的作用，会间接促进细胞的凋亡。p53蛋白激活其靶基因IGFBP3(该基因是本研究所筛选出的DCM的基因标志物)，IGFBP3会抑制IGF(胰岛素样生长因子)，最终导致细胞凋亡。p53蛋白激活其靶基因IGFBP3，而IGFBP3间接作用于IGF-1/mTOR通路。在动物实验中，Zhang等学者将小鼠的心脏进行特异性的mTOR敲除，之后可导致小鼠发生致死性DCM[26]。因此，可推测出p53/IGF-1/mTOR通路可能与DCM的致病相关，但具体的作用机制尚未明确。由此可知，本研究的部分预测结果得到了现有生物证据的支持，提示这些预测结果的可靠性。

图3 p53通路与DCM的关系(其中，黑色箭头实线表示直接作用，黑色箭头虚线表示间接作用，黑色直线加上小短线代表抑制作用，蓝色箭头虚线代表可能导致DCM)Fig.3 Relationship between p53 pathway and DCM (The black solid line with arrow represents direct regulation, the black dotted line with arrow indicates indirect regulation, the black solid line with short line represents inhibition and the blue dotted line with arrow indicates that DCM were possibly induced)

3.2 基于生物信息学工具的RCM相关基因标志物的功能分析

3.2.1RCM相关基因标志物的GO分析

RCM的65个基因标志物用系统工具PANTHER进行GO功能富集分析，得出如下结果：在生物学过程中，主要富集在细胞生理(38.5%)、代谢(36.9%)、应激(12.3%)等过程中；在分子功能中，主要富集在催化活性(32.3%)、蛋白结合(23.1%)等方面；在细胞组成成分上，主要富集在细胞组分(18.5%)、细胞的细胞器(13.8%)、细胞外基质(7.7%)等成分上。

3.2.2RCM相关基因标志物的信号通路分析

同样，由PANTHER可得出DCM相关的信号通路的结果，包括Wnt信号通路、Notch信号通路、EGF受体信号通路、血管生成通路、血液凝固通路等路径。

RCM的基因标志物FRZB和FBXW11参与了Wnt信号通路，而Wnt信号通路在心脏的发育以及心肌细胞的分化过程中起着非常重要的作用。Cohen等学者通过生物实验，对小鼠的胚胎进行Wnt5a和Wnt11的同时敲除，导致了经典Wnt信号的激活，最终造成心脏发育缺陷[27]。越来越多证据表明，Wnt信号通路参与了心脏重构和心衰的进程[28]，虽然RCM的致病通路在国内外文献中很少涉及，但是研究发现RCM的基因标志物能关联到此通路。因此猜测，在RCM的病程中，RCM很容易发展为心衰(与RCM预后性差符合)，并且可能与Wnt通路相关。

RCM的基因标志物NOTCH4介导了Notch信号通路。NOTCH4编码Notch蛋白，间接作用于MAPK信号通路。在所发现的MAPK信号通路中，JNK MAPK和p38 MAPK信号通路在应激反应中发挥重要作用，如炎症反应、细胞凋亡等。文献[29]表明，p38 MAPK信号通路在RCM中有着重要的作用。除此之外，在心脏发育中，Notch信号通路也发挥着重要的调控作用，参与了心肌的分化过程，并且在血管的损伤修复中起着重要的作用。

3.3 DCM与RCM的比较

3.3.1DCM与RCM的基因标志物比较

DCM与RCM只有一个相同的基因标志物——NPR3，该基因已被证实可以保护心肌细胞[24]。

3.3.2GO分析结果的比较

DCM和RCM有着相同的GO功能富集分析条目。在生物学过程中，两者都主要与细胞生理和代谢过程相关；在分子功能方面，两者都具有蛋白结合和催化活性两类功能；在细胞组成成分上，两者都富集在细胞组分和细胞外基质上。

然而，DCM与/RCM也有着不相同的GO功能富集分析条目。在生物学过程中，DCM基因标志物富集在多细胞生物、组织或生物起源细胞组件等过程中，而心肌细胞凋亡与DCM的发生发展息息相关，这与其基因标志物富集在多细胞生物等过程相符；RCM的基因标志物则还富集在应激过程中，说明RCM的发生发展可能与应激性相关。在细胞的组成成分上，DCM的基因标志物富集在膜成分上，而RCM的基因标志物则富集在细胞的细胞器上。基因标志物在细胞中位置的分布，在一定程度上反映了其调控蛋白的功能。

有学者对细胞中的钠钾钙镁离子及ATP与正常人细胞中的含量及其转运功能进行研究，发现细胞中离子的异常变化及其转运的障碍在DCM病程的发生发展中起着重要的作用[30]，而细胞的膜正是与离子的转运直接相关联的。这个证据间接地证明了上述分析的可靠度。

3.3.2信号通路分析结果的比较

DCM与RCM与相同的信号通路相关，分别为EGF受体信号通路和异三聚体G蛋白信号通路。但是，介导同一通路的基因标志物却是不同的，如DCM的基因标志物NRG1和 RCM的基因标志物GAB1分别介导EGF受体信号通路。

然而，两类心肌疾病相关的生物信号通路也存在较大差异。DCM的基因标志物参与P53通路、炎症介导信号通路、β肾上腺素能受体信号通路以及G蛋白受体激活通路等，DCM病程的发生发展机制主要与炎症通路、细胞凋亡以及Ca离子等密切相关。RCM的基因标志物参与Wnt信号通路、Notch信号通路等，这些通路都与心衰有着非常密切的关系。由此推测，RCM的预后性差可能与基因标志物介导这些通路相关，使得RCM更容易向心衰发展。

3.4 本研究的局限与进一步的工作

在本研究中，分别识别出了DCM和RCM相关的基因标志物，且部分基因标志物有文献支持，如DCM致病相关基因MAB21L2和POSTN、RCM致病相关基因TNNI3和DES。虽然现有相关文献和生物实验证实了上述基因在疾病的发生发展进程中起到了重要作用，但是其作为靶标的特异性并未得到证实。因此，下一步的工作需要对这些可能的基因标志物进行特异性的验证，以找到这两类心肌疾病的治疗靶点。

此外，本研究发现了与两类心肌疾病相关的一些新的基因，如在DCM中找出的SNAP25、COL8A1、B2M等基因，在RCM中找出的TSC1、FBXW11、NOTCH4等基因。虽然这些新发现的基因与疾病的直接关系在生物实验和现有的文献中并未得到证实，但是GO和KEGG的分析表明，这些基因都与心血管疾病相关，其所介导的通路也与这两类心肌疾病息息相关，因而具有作为两类心肌疾病标志物的潜力，值得开展进一步的生物实验来验证它们的功能。

最后，笔者对DCM和RCM的发生发展机制进行了一些新的阐释。由于现有的相关文献及动物实验证据太少，因此对于完全理解DCM和RCM复杂的发病机制还存在较大距离，还有待后续的深入研究。