缺乏运动降低成年人胰岛素敏感性的加权基因共表达网络分析*

刘会平 刘湘茹 李 智

【提 要】 目的 通过加权基因共表达网络分析(weighted gene co-expression network analysis，WGCNA)分析久坐影响成年人胰岛素敏感性的关键基因及其相关功能和机制。方法 针对下载于NCBI的人类基因表达汇编(gene expression omnibus，GEO)数据库的数据集GSE9103，利用WGCNA进行模块化分析，并通过相关性分析和网络可视化确定模块内关键hub基因及其调控网络，再对关键hub基因进行基因本体GO分析和基因集富集GSEA分析。结果 通过构建共表达网络，最终得到5个基因模块，其中yellow模块与久坐呈正相关，而NDUFA9为其枢纽基因；tan模块与久坐呈负相关，而MYBPC2为其枢纽基因。进一步对GO和GSEA分析发现，yellow模块基因包括NDUFA9主要发挥线粒体形成、有氧呼吸、脂肪酸降解等生物学功能并参与物质代谢、三羧酸循环、阿尔茨海默氏症、帕金森病相关通路；而tan模块基因包括MYBPC2则与骨骼肌收缩、糖原代谢、心脏发育相关。结论 本研究通过构建基因共表达网络筛选出与久坐生活方式显著相关的2个关键基因模块和2个枢纽基因；这2个模块基因分别通过线粒体功能、物质代谢和骨骼肌收缩、糖原代谢而影响成年人对胰岛素的敏感性。

胰岛素敏感性降低是导致肥胖、代谢综合症、糖尿病、高脂血症等多种慢性疾病的重要临床问题，随着年龄的增长，其发生风险显著增加。线粒体功能障碍表现为线粒体还原酶活性、蛋白质合成和表达以及DNA(mtDNA)丰度的下降，多见于上述人群[1]。已有报道显示，在老年、2型糖尿病患者和肥胖者以及2型糖尿病患者的后代等群体内，胰岛素敏感性与肌肉线粒体功能密切相关。肌肉线粒体功能障碍是产生胰岛素敏感性下降甚至胰岛素抵抗的重要影响因素[2]。也有研究表明，久坐的生活方式引起显著的胰岛素抵抗，定期运动可提高胰岛素敏感性，均与年龄无关，但有研究发现老年人较年轻人具有显著的线粒体功能障碍[3]。因此，有必要进一步针对久坐的生活方式引起胰岛素敏感性下降的具体机制进行深入研究。

本研究利用加权基因共表达网络分析(WGCNA)，从大量基因中筛选与久坐的生活方式密切相关的枢纽基因，结合对基因模块的基因本体分析和基因集富集分析，进行功能注释，探讨久坐的生活方式较定期运动更可能引起胰岛素抵抗的相关分子机制与调控网络。

资料与方法

1.数据集和研究对象

本研究分析的数据集GSE9103下载于NCBI的GEO数据库(http：//www.ncbi.nlm.nih.gov/geo)。该研究纳入的40名健康人中，老年人(59～76岁)20人，年轻人(18～30岁)20人。久坐组研究对象(10名老年人、10名年轻人)的运动时间每周少于2次，每天少于30分钟；而运动组研究对象(10名老年人、10名年轻人)在过去至少4年的时间内每周运动6天且每天骑车或跑步至少1小时。研究对象均需填写业余时间活动调查表，并经过病史、体检、分级跑步机测试、血化验(如血脂、激素水平、血糖检测)的筛选。有代谢性疾病或心血管疾病病史、血糖>99mg/dl、体重指数>28kg/m2、应用影响研究测量结果的药物、贫血、妊娠、药物/物质成瘾滥用者被排除。基因组数据的获得方法如下：研究对象经过一次分级跑步机测试，此后7天内严格控制饮食与运动，并在测试后第630小时利用局麻下经皮细针穿刺取腓肠肌活组织(约350～400mg)进行基因组分析。详细信息参见已发表论文[3]。

2.基因数据的获得

提取研究对象腓肠肌活组织中的总RNA并应用Human Genome U133 Plus 2.0 进行基因表达微阵列分析，获得该数据集series matrix数据，以方差最大的探针为基本单位进行后续分析。

3.基因共表达网络构建与关键枢纽基因的筛选

利用R软件的WGCNA工具包，构建基因共表达网络并筛选出与表型高度相关的基因模块。再通过计算模块中每个基因与该模块(module membership，MM)的相关系数，筛选出核心基因(GS 绝对值>0.2，MM绝对值>0.8)，纳入string数据库(https：//string-db.org/)构建基因相互作用网络，再利用Cytoscape和cytoHubba MCC法进行可视化分析与前10位关键作用基因的鉴定筛选。

4.基因GO与GSEA分析

Metascape 平台基因注释分析数据库可对输入基因进行基因本体(GO)分析。分别提交tan和yellow模块包含的基因，输入物种和分析物种均选择为“H.sapiens”，设置P<0.01，进行 GO 分析，保存结果，并按照每个条目涉及的靶点数目进行排序，筛选排名靠前的生物过程和通路，采用 GraphPad Prism 7.0 软件对结果可视化。

应用Gene Set Enrichment analysis 3.0版软件，以hubgenes的表达值为表型，先后以MSigDB在线基因集数据库的KEGG/ GO gene sets作为参考基因集，对基因表达谱进行基因集富集分析(GSEA)。随机组合次数设置为1000次。以P<0.05，且FDRq<0.25为具有统计学意义。

5.统计分析

WGCNA[4]基于两个假设：表达模式近似的基因可能接受共调控或处于同一通路；表达网络符合无尺度分布。由此根据表达相似性，将基因划分为不同模块并找出枢纽基因。进行WGCNA之前，首先进行数据预处理，构建N×M矩阵(N=样本数，M=基因数)。通过层次聚类方法，去除离群样本。基于pickSoftThreshold函数筛选软阈值构建邻接矩阵和拓扑矩阵，以使基因分布符合无尺度网络。根据相异度对基因进行聚类，使用动态剪切法将聚类树剪切为不同模块，以每个模块的第一主成分基因为特征向量基因。根据表型与特征向量基因的相关性挖掘相关性P值最小且具有统计学意义的模块。通过计算模块所包含基因的内部连接度筛选出枢纽基因。

结 果

1.基因模块与核心基因鉴定

共5000个基因用于WGCNA。以邻接矩阵权重参数(软阈值)β=16、每个基因网络模块最少的基因数目为30作为标准，最终获得black、blue、greenyellow、tan和yellow 5个模块，以|MM|>0.9和|GS|>0.2为参数筛选出其核心基因分别为GORASP2、LINC00491、S100A4、MYBPC2和NDUFA9。再根据各个模块的特征向量，分别计算模块与 2 个表型(久坐和定期运动)的相关性。如图1所示，yellow模块(400个基因)与久坐呈正相关(COR=0.53)，而tan模块(30个基因)与久坐呈负相关(COR=-0.62)。

图1 与定期运动相比，久坐者共筛选出5000基因与久坐显著相关的基因模块与核心基因

2.yellow模块枢纽基因网络分析

将yellow模块的枢纽基因上传string数据库，利用cytoscape进行网络可视化，并利用cytohubba插件的MCC方法筛选出前10位基因，包括ATP、MRPL、NDUF家族等在共表达网络和蛋白相互作用网络皆重要的关键作用基因(图2)。

图2 分析明确久坐相关的基因与yellow模块基因显著相关并进行网络可视化筛选yellow模块的前10位关键基因

3.yellow模块枢纽基因的功能分析

如图3和表1所示，yellow模块基因行GO分析显示其具有线粒体形成、有氧呼吸、脂肪酸代谢等生物学功能如mitochondrial envelope/matrix/respiratory chain complex III，aerobic respiration，fatty acid degradation/beta-oxidation using acyl-CoA dehydrogen等。进一步KEGG分析则富集了KEGG_OXIDATIVE_PHOSPHORYLATION、KEGG_CITRATE_CYCLE_TCA_CYCLE、KEGG_PORPHYRIN_AND_CHLOROPHYLL_METABOLISM、KEGG_ALZHEIMER_DISEASE、KEGG_PARKINSONS _ DISEASE等与糖脂代谢、三羧酸循环、阿尔茨海默氏症、帕金森病密切相关的信号通路。

图3 Metascape平台对yellow模块基因进行GO分析明确其生物学功能

表1 久坐者骨骼肌中显著上调且与NDUFA9共表达的yellow模块基因的GSEA分析

4.tan模块枢纽基因网络分析

tan模块中，与MYBPC2密切相关的蛋白包括多种能够调节骨骼肌功能、减少细胞凋亡、抑制炎症、调节血脂代谢和抑制血管新生等的蛋白质，如：GLRX、PHKA1、NIBAN1、CYP27A1、TNNI2、CD38、TNNT3、ATP2A1、PPP1R3D、CASQ1、MTHFD2、ANKRD23、TNFRSF1A等(图4)。

5.tan模块枢纽基因的功能分析

tan模块基因行GO分析富集了调节骨骼肌收缩、糖原/蛋白代谢、凋亡、氧化还原酶活性及心脏发育等生物学功能(图5)。KEGG分析按本研究标准未能成功富集信号通路。

图4 分析明确久坐相关的基因与tan模块基因显著相关并进行网络可视化筛选tan模块的关键节点基因

图5 Metascape平台对tan模块基因进行GO分析明确其生物学功能

讨 论

据报道，yellow模块的代表基因如ATP、NDUF和MRPL家族蛋白分别具有调节线粒体呼吸链功能和线粒体内蛋白质合成的生物学功能，是影响线粒体形成和功能的重要因素。其中，ATP5F1B、ATP5MC3和ATP5ME编码的蛋白为线粒体合成酶的重要亚基，对线粒体合成ATP具有重要的调节作用。已有研究表明，上述ATP家族蛋白表达升高与肥胖、高血压[5]、血脂代谢异常[6]、胰岛素抵抗[7]、2型糖尿病、肾积水、恶性肿瘤侵袭特性和转移耐药等临床问题密切相关[8]。NDUFA9 和NDUFS2编码的蛋白质分别是NADH脱氢酶复合体的疏水蛋白部分和核心亚基，NADH脱氢酶复合体是位于线粒体内膜的电子传输链中的第一个酶复合物。因此，NDUF家族蛋白对调节线粒体呼吸链功能具有重要作用。NDUF家族蛋白可引起新生儿Leigh综合征、恶性肿瘤的转移潜能、视神经病变、早发精神分裂症等临床问题[9-10]。MRPL15和MRPL41编码的蛋白质是线粒体核糖体蛋白，有助于线粒体内蛋白质的合成。二者的临床报道很少，目前发现其表达增加与ER阳性的乳腺癌相关[11]，而相对较多的基础研究主要发现其与细胞周期抑制和凋亡相关[12]。本研究利用骨骼肌标本的基因组学分析发现，与定期运动相比，久坐的生活方式导致上述基因增加，并发现其通过发挥促进线粒体的形成、呼吸链、脂肪酸降解、细胞色素C释放、辅酶A水解等功能调节糖脂代谢、三羧酸循环、阿尔茨海默氏症、帕金森病相关通路。因此，上述基因及其调节的相关通路，可能是久坐患者发生胰岛素抵抗的重要机制，其对胰岛素抵抗以及相关疾病的影响值得深入研究。

另一方面，久坐者的MYBPC2基因显著低表达。该基因编码的是快速型蛋白肌球蛋白结合蛋白C，与肥厚性心肌病的发生、发展相关[13]。本研究发现的MYBPC2网络互作蛋白中，GLRX编码蛋白催化谷胱甘肽还原反应，在抗氧化防御系统中发挥重要作用[14]。PHKA1编码磷酸化酶激酶，异常可引起糖原贮积病[15]。CYP27A1编码细胞色素P450超家族蛋白，催化涉及药物代谢和胆固醇，类固醇和其他脂质合成的许多反应，对于调节胆固醇的体内平衡非常重要[16]。TNNI2编码快缩肌蛋白，参与组成肌钙蛋白复合体，调节骨骼肌与血管平滑肌功能，还具有抑制肿瘤发生发展的作用[17]。TNNT3、ATP2A1、ANKRD23及CASQ1均参与骨骼肌收缩与功能的调节[18-20]。PPP1R3D编码丝苏氨酸蛋白磷酸化酶，而TNFRSF1A编码肿瘤坏死因子超家族蛋白，二者均具有广泛的调节作用，参与多种细胞生物学活动如细胞生存、凋亡和炎症等[21]。本研究发现tan模块中MYBPC2及其网络互作蛋白在久坐者体内较规律运动者显著减低，可能为久坐者的骨骼肌收缩、糖原代谢等活动减低的适应性反应，但同时可能是引起胰岛素抵抗的重要机制。

本研究通过对基因数据集GSE9103进行基因组学WGCNA分析，针对定期运动者与久坐者的基因表达谱，构建了与久坐密切相关的基因共表达网络，并对核心枢纽基因及其功能机制进行了多重分析，从基因水平上揭示了缺乏运动导致胰岛素抵抗的发生机制，对于探讨胰岛素抵抗以及相关疾病的调控机制和靶点具有重要意义。