2型糖尿病患者心肌间充质细胞基因表达变化及相关环境化学物的筛选

刘 辉，段佳琪，任晓梅，李熙蕾，陈 静，李 鹏，张荣强

1陕西中医药大学公共卫生学院，陕西咸阳 712046 2兰州大学公共卫生学院营养与食品卫生学研究所，兰州 730000

糖尿病(diabetes mellitus，DM)是一种常见的内分泌代谢异常综合征，主要表现为血糖升高，严重危害人类健康。2017年，18～99岁的DM患者全球有4.51亿，预计2045年将升至6.93亿[1]；我国20～79岁居民中共有DM患者约1.14亿，预计2045年将升至1.20亿[2]。继肿瘤、心血管之后，DM已经发展成为第3种对人类健康威胁程度很大的慢性疾病[3]，成为全世界重要的健康问题。90%的DM患者为2型糖尿病(type 2 diabetes mellitus，T2DM)[4]，T2DM的危险因素、发病机制、预防策略研究成为该领域的重要内容之一。近年来，越来越多<30岁的居民被诊断为T2DM，久坐和肥胖是其重要的危险因素[3，5]，然而其发病机制仍然未明。研究表明，DM主要由遗传因素与环境因素导致的胰岛素分泌绝对或相对减少而造成[6]。积极探索遗传和环境因素在T2DM发病中的作用机制，受到该领域研究者越来越多的重视。

组学和生物信息学技术为糖尿病的机制研究提供了新思路和新方法[7]，并取得一定研究成果。李中辉等[8]通过分析芯片数据集GSE7014与GSE29221，发现CDK9、TXN以及NDUFS8可能是T2DM潜在治疗靶点。任晓梅等[9]研究糖尿病肾病(diabetic kidney diseases，DKD)基因芯片数据，认为炎症和应激反应与DKD发病有一定关联。Leopold等[10]发现诸多代谢因素中，DM是对心血管影响最大因素之一；Lamb等[11]认为是否患DM与心血管疾病(cardiovascular disease，CVD)的发生存在密切关系，但糖尿病与CVD的具体关系还需深入研究。心肌间充质细胞(cardiac mesenchymal cells，CMSC)对心肌细胞的成型及CVD发病均具有重要作用。

为了探讨遗传和环境因素对T2DM患者CMSC细胞基因表达的影响及作用机制，本研究利用来源于公共数据库(Gene Expression Omnibus，GEO)的数据集GSE106177[12]，首先筛选T2DM患者CMSC的差异表达基因，通过富集分析诠释其生物学功能，建立PPI网络、microRNA-Gene交互作用网络，同时分析敏感差异表达基因与环境化学物的交互作用，综合探讨了遗传和环境因素对T2DM患者心肌细胞基因表达特征的影响机制。

资料和方法

差异表达的基因数据以“Diabetes”为关键词在美国国立生物技术信息中心(National Center for Biotechnology Information，NCBI)的GEO数据库进行检索，获得Spallotta F等提交的基因测序数据GSE106177(https：//www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE106177.)，该数据集于2017年10月25日提交，并于2019年1月23日进行更新，检测平台为Illumina NextSeq 500 (Homo sapiens)。GSE106177共包含14例CMSC样本，其中7例来自T2DM供者，7例来自非2型糖尿病(non-type 2 diabetes mellitus，NT2DM)供者；该数据共检测了36 317个基因的表达。

数据预处理本数据集属于RNA测序数据，提供了基因外显子水平的表达测量，有多个探针被定位到同一基因的情况，此时需要计算该基因的平均值；还可能有无用或空白数据，此时需要剔除此类数据，然后进行标准化处理以确保每个样品的表达分布是同质的。本研究采用Network Analyst[13- 14]在线分析软件对原始数据进行log2CPM (log2-counts per million)标准化处理，以保证后续分析结果准确、可靠。

差异表达基因分析将GSE106177数据集导入Network Analyst(https：//www.networkanalyst.ca/)，采用Limma法筛选两组的差异表达基因(DEGs，adjustedPvalue<0.05，Log2fold change>2)。

主成分分析将高维的测序数据GSE106177投射到三维空间中进行主成分分析(principal component analysis，PCA)，实现“降维”，观察各组样本间的总体分布状况，鉴别是否存在奇异样本。采用Network Analyst对数据预处理后GSE106177数据集行PCA分析。 Network Analyst的具体运行参数严格依据其操作说明进行分析(https：//www.networkanalyst.ca/)。

生物功能和KEGG信号通路富集分析将筛选出的差异表达基因上传至在线分析软件String 11.0(https：//string-db.org/cgi/input.pl)，根据需要调节主要参数(Confidence>0.05)，对两组DEGs进行生物功能及其涉及的KEGG通路进行富集分析。

PPI网络的建立将两组DEGs相应的蛋白名称上传至Network Analyst，为降低假阳性率，调节适当的过滤参数(Degree cutoff >5)，利用Cytoscape 3.7.1进行可视化，建立差异表达基因的PPI网络。

microRNA-Gene交互作用分析利用Network Analyst和Cytoscape工具，分析差异表达基因与microRNA的交互作用，调节过滤参数(Degree cutoff>5)，绘制相应的交互作用网络图，以了解差异表达基因上游的microRNA。

差异表达基因与环境化学物交互作用分析除上述遗传因素外，本研究将PPI网络中的核心节点基因、差异表达最明显的两个DEGs(上、下调各选取1个)上传至CTD(http：//ctdbase.org/)，挖掘可能影响此3个重要差异基因的环境化学物，以解释环境因素对T2DM患者CMSC基因表达的影响机制。

结 果

测序数据预处理将原始测序数据GSE106177上传至Network Analyst中，选择log2CPM对数据进行标准化，结果显示各样本基因表达数据均数基本为一条直线，各样本基因表达数据分布密度基本一致，数据质量良好，可进一步分析(图1)。

主成分分析基于原始数据集，本研究构建了各样本基因表达模式PCA图，结果显示：PC1轴、PC2轴、PC3轴的数据方差贡献率分别是26.5%、13.7%、9.2%(图2)。

差异表达基因分析

根据筛选条件，经差异分析工具共发现135个DGEs(AdjustedP<0.05、log2fold change>2.0)，其中上调基因58个(42.96%)，下调基因77个(57.04%)；差异表达基因分析火山图见图3A；差异最明显的Top20 DEGs见表1；上调和下调最明显的8个基因在两组的表达状况分别见图3B～I。

NT2DM：非2型糖尿病；T2DM：2型糖尿病
NT2DM：non-type 2 diabetes mellitus；T2DM：type 2 diabetes mellitus
图1经数据预处理后各样本基因表达数据的箱式图(A)和分布密度曲线(B)
Fig1Box plot (A) and distribution density curve (B) of gene expression data of each sample after data preprocessing

PCA：主成分分析
PCA：principal component analysis
图2T2DM组与NT2DM组基因表达模式PCA图
Fig2PCA of gene expression patterns of T2DM and NT2DM groups

生物功能富集分析生物功能富集分析结果显示，两组DEGs主要参与多细胞生物发育、解剖结构发展、系统发育等生物学过程；涉及高密度脂蛋白颗粒受体结合、蛋白质结合等分子功能及质膜部分、质膜的内在成分等细胞组分(表2)。

KEGG信号通路富集分析两组的差异表达基因主要涉及肝细胞癌、库欣综合征、胃癌、黑色素瘤以及胆固醇代谢等信号通路(表3)。

差异表达基因的PPI网络结果建立DEGs的PPI并可视化后，结果显示，在PPI网络中，UBC与其他24个蛋白(ACTN2、HIST1H4E、APBB1IP等)存在相互作用关系，是最为重要的中心节点。删除UBC后，可看出整个网络结构的稳定性降低，且结构变得异常涣散(图4)。

差异表达基因与microRNA交互作用关系差异表达基因的microRNA-Gene交互作用网络(Degree cutoff >5)中共有7个microRNA与DEGs有交互作用，分别为hsa-mir- 8485、hsa-mir- 190a- 3p、hsa-mir- 5011- 5p、hsa-mir- 1277- 5p、hsa-mir- 26b- 5p、hsa-mir- 124- 3p、hsa-mir- 335- 5p，其中hsa-mir- 8485与9个DEGs具有靶向作用关系(图5)。

T2DM重要相关基因与环境化学物的交互作用分析将PPI网络中的核心节点基因、差异表达最明显的2个DEGs(上、下调各选取1个)上传至CTD(http：//ctdbase.org/)，结果显示，UBC与Bisphenol A、Cadmium等9种环境化学物具有交互作用；DNER与Arsenic、Bisphenol A等8种环境化学物具有交互作用；CNTN1与Arsenite、Atrazine等6种环境化学物具有交互作用。构建上述3个基因所涉及环境化学物的Venn图，3者均与Bisphenol A(双酚A)存在交互作用(表4、图6)。

表1 T2DM与NT2DM组表达差异前20位的基因Table 1 Top 20 differentially expressed genes in T2DM and NT2DM groups

T2DM：2型糖尿病；NT2DM：非2型糖尿病
T2DM：type 2 diabetes mellitus；NT2DM：non-type 2 diabetes mellitus

图3差异表达基因的火山图(A)及前8个基因的表达模式(B～I)
Fig3Volcanic plot of differentially expressed genes(A)and expression patterns of the first 8 differentially expressed genes(B-I)

表2 差异表达基因的GO功能富集结果Table 2 GO function enrichment results of differentially expressed genes

续表2

表3 KEGG信号通路富集分析结果Table 3 KEGG signaling pathway enrichment results of differentially expressed genes

讨 论

T2DM近年来在我国的患病率不断上升，从1980年的0.67%升高到2013年的11%[15- 16]，33年增长了16倍之多。有研究报道，T2DM的发展会使心血管疾病发生风险显著增高，约为非糖尿病者的2至2.5倍[17]，其所引起的心脑血管疾病等并发症也严重影响着居民的身体健康状况，是我国DM患者死亡原因中占比较重的一部分[16]。本研究基于生物信息学技术，对T2DM患者CMSC测序数据进行深入分析，旨在了解其基因变化情况、功能和信号途径富集情况、PPI交互作用网络、microRNA-Gene交互作用等遗传特征，并通过分析DM重要相关基因与环境化学物的交互作用，了解环境因素对T2DM患者CMSC的影响。

本研究结果显示，T2DM患者CMSC发生明显改变，其主要差异表达基因主要参与多细胞生物发育、解剖结构发展、系统发育等生物学过程，提示T2DM患者CMSC对心肌细胞的功能和发育发挥重要作用，其基因表达的改变可能影响T2DM患者的心肌功能。

图4差异表达基因的PPI网络图
Fig4PPI network of differentially expressed genes

蓝色菱形代表microRNAs，红色椭圆代表上调DEGs，绿色椭圆代表下调DEGs
Blue diamonds represent microRNAs，red ovals represent up-regulated DEGs，and green ovals represent down-regulated DEGs
图5差异表达基因的microRNA-Gene交互作用网络
Fig5MicroRNA-Gene interaction network of differentially expressed genes

表4 T2DM重要相关基因与环境化学物的交互作用分析Table 4 Analysis of interactions between important DEGs related to T2DM and environmental chemicals

图63个T2DM重要相关基因所涉及环境化学物的Venn图
Fig6Venn diagram of environmental chemicals involved in three key T2DM related genes

本研究还发现，两组差异表达基因的PPI拓扑网络中存在1个典型的子网络，其核心节点是UBC。UBC全称为ubiquitin C，人UBC基因定位于12号染色体，基因ID为7316[18]。Ryu等[19]发现UBC基因有助于维持对细胞应激和耐热性的反应。Marinovic等[20]在胰岛素缺乏的大鼠体内证实UBC反应的肌肉特异性，通过电泳迁移率变动分析和体内基因组足迹证明胰岛素缺乏增加了Sp1(一种通用的转录激活因子)与UBC启动子中富含GC元件的结合，得出了糖皮质激素通过涉及肌肉特有的Sp1机制增加UBC转录的结论。位于UBC末端水解肽键的UCHL1在脑、睾丸/卵巢和胰岛中富集表达[21- 22]。Chu等[23]通过蛋白质组学筛选，认为UCHL1是用棕榈酸盐处理的MIN6β细胞中最上调的蛋白质，同时Lopez-Avalos等[22]研究证明UCHL1在胰腺β细胞中的重要性。综合本研究结果，提示UBC可能通过影响胰岛的正常功能而参与DM的发生发展。

microRNA-Gene交互作用结果发现，差异表达基因受7个microRNA的调控。有研究报道mir- 8485(mir-NID1)可以与TDP- 43基因结合，抑制NRXN1的表达，影响人体健康[24]。mir- 190在肝细胞癌和大多数胰腺癌组织和细胞系中表达上调[25- 26]。在CVD中，mir- 26涉及控制心肌梗死后内皮细胞生长和血管生成相关的关键信号通路[27]，且研究发现mir- 26b- 5p表达上调，心肌能力得到改善，N端脑钠肽(N-terminal-pro-brain natriuretic peptide，NT-proBNP)水平显著降低[28]。mir- 26b的过表达促进了人前脂肪细胞中的脂肪形成，从而影响肥胖[29]，且mir- 26b- 5p在T2DM和健康对照者中的表达存在显著差异[30]。急性心肌梗死患者外周血中mir- 124水平显著升高，mir- 124比其他诊断标志物有助于早期检测，可被认为是诊断和评估急性心肌梗死发病和发展的潜在生物标志物[31- 32]。另有结果表明，血清miR- 124在急性缺血性卒中受到抑制，从而促进神经炎症和脑损伤[33]。研究人员发现，通过抑制mir- 124表达可以介导姜黄素对糖尿病肾病中足细胞的有益作用[34]。Kay等[35]研究表明，mir- 335- 3p/5p是心脏分化过程的调节因子。上述结果表明，mir- 26、mir- 124、mir- 335- 5p与心脏之间存在关联，mir- 26b可能与DM的发生与发展相关，与本研究分析结果基本一致。

T2DM的危险因素不仅仅包括遗传与环境，还有遗传与环境交互作用[36]。本研究在上述遗传因素分析的基础上，探索了T2DM重要关联基因的环境化学物交互作用。通过CTD分析T2DM重要关联基因与环境化学物的相互作用发现，UBC所关联的环境化学物数目最多，其次为DNER，最后为CNTN1。3个基因均与双酚A这种环境化学物存在一定关系，这提示双酚A可能是T2DM患者心肌功能紊乱的危险因素。环境化学物暴露能够通过一系列方式扰乱基因组的功能[37]， Huang等[38]对Embase和Medline (PubMed) 数据库中16项研究双酚A和T2DM相关性的文献进行了Meta分析，结果显示双酚A暴露水平与T2DM患病风险呈正相关。Lang等[39]对1455名18至74岁的美国成年人研究发现，较高的尿双酚A浓度与糖尿病有关联。目前已有多项研究表明，双酚A还可以明显提高心脏代谢紊乱的风险，明显增加CVD的发病风险[40- 41]。本研究结果与上述结果基本一致，并为上述研究提供了一定的基因组学基础。综上，本研究结果显示，T2DM患者CMSC基因的表达模式发生明显改变，UBC可能在上述过程中发挥重要生物学作用。双酚A暴露也可能影响T2DM患者心肌细胞的发育和正常功能，上述结论尚待在今后实验中进一步深入研究证实。