间充质干细胞成骨分化与成脂分化关键基因的生物信息学分析

王碧雯,骆亚莉,b*,李 研,马 玉,周啸天

(甘肃中医药大学a.甘肃省高校重大疾病分子医学与中医药防治研究省级重点实验室;b.基础医学院病理教研室,中国甘肃 兰州 730000)

间充质干细胞(mesenchymal stem cells,MSCs)有多向分化和自我更新的优良特性,这些优良特性使其具有修复和再生受损组织的治疗潜力[1～2]。MSCs直接或间接参与骨缺损[3]、软骨缺损[4]、心肌梗死[5]、神经损伤[6]等组织损伤的修复,这使其在临床应用中具有良好的前景[7]。目前,针对MSCs成骨和成脂分化相关转录因子、蛋白质的研究多见于体外或体内实验,从生物信息学角度对成骨成脂分化关键基因展开分析的综合研究少见,本文通过生物信息学技术在众多相关基因中筛选出关键基因,并试图揭示成骨成脂分化可能存在的关联点和作用机制,这对于进一步研究MSCs异常成骨成脂分化的发生机制具有重要意义,同时可为骨质疏松症等骨与脂肪生成不平衡的相关疾病提供治疗思路。

1 材料与方法

1.1 文献检索及基因统计

以“MSC”“osteogenic differentiation”“adipogenic differentiation”为关键词,通过PubMed检索公开发表的与MSCs成骨分化和成脂分化相关的研究文献,纳入文献包括与MSCs成骨分化和MSCs成脂分化有关的基因或基因产物的临床研究、论著、随机对照试验、综述及系统评价文献,排除无法获取摘要的文献、重复文献及其他与研究目的不符的文献,剔除重复、错误的基因或基因产物,最终摘选出MSCs成骨分化与成脂分化相关的差异表达基因(differentially expressed genes,DEGs)或基因产物。

1.2 生物功能和通路富集分析

将1.1所得有关基因分别导入在线软件DAVID(https://david.ncifcrf.gov/tools.jsp),Select I-dentifier选择 Official gene symbol,Species选择Homo sapiens,Background设置为Homo sapiens,设定P＜0.01,对相关基因进行基因本体论(Gene Ontology,GO)分析与京都基因和基因组数据库(K-yoto Encyclopedia of Genes and Genomes,KEGG)通路分析。

1.3 基因所表达蛋白质相互作用网络构建和关键基因筛选

将1.1所得成骨分化有关基因和成脂分化有关基因分别上传至STRING 11.0(http://STRING-db.org/)软件进行分析,物种设置“Homo sapiens”,构建由导入基因表达产物所组成的相互作用网络,以TSV格式导出,并将源文件导入Cytoscape软件进行可视化分析。利用Cytoscape软件中的插件Centiscape计算出网络中每个节点的拓扑特性,继而根据节点度数(degree)和节点中介性(betweenness centrality)筛选出关键节点,即MSCs成骨分化与成脂分化的关键基因表达产物。

1.4 成骨分化与成脂分化共同关键基因筛选

将1.3筛选所得的关键基因表达产物导入Venny 2.1.0(http://bioinfogp.cnb.csic.es/tools/venny/index.html)在线工具,对成骨分化和成脂分化的关键基因进行交互比较,得到成骨分化和成脂分化共同关键基因。

2 结果

2.1 文献检索及入选情况

经PubMed数据库检索,共获得成骨分化相关文献2 083篇,成脂分化相关文献1 028篇,按照纳入标准进一步筛选出符合要求的成骨分化文献279篇,成脂分化文献164篇;最终筛选出成骨分化基因及表达产物217个,成脂分化基因及表达产物227个。

2.2 GO分类和KEGG通路分析

利用DAVID数据库对MSCs成骨分化相关的DEGs进行GO分类,得到P＜0.01的分类281种,其中生物过程(biological process,BP)条目221个,细胞组分(cellular component,CC)条目19个,分子功能(molecular function,MF)条目41个。在生物过程层面,这些成骨分化DEGs主要与RNA聚合酶Ⅱ启动子转录的正调控和负调控、成骨细胞分化等有关(排名前10位的条目见图1A);在细胞组分层面,这些基因大多参与细胞体、黏着斑、质膜外侧、早期内体膜等的组成(排名前10位的条目见图1B);其分子功能的变化主要集中在转录因子结合、生长因子活性、转录激活子活性以及RNA聚合酶Ⅱ核心启动子近端区域序列特异性结合等(排名前10位的条目见图1C)。进一步对上述基因展开KEGG通路分析,结果显示其主要涉及调节干细胞多能性的信号通路、转化生长因子β(transforming growth factor-β,TGF-β)信号通路、Hippo信号通路等13条信号通路(P＜0.01),排名前10位的通路信息见图1D。

图1 成骨分化差异表达基因的GO分类和KEGG通路分析(前10位)Fig.1 Enrichment analysis of GO and KEGG pathway of DEGs in osteogenic differentiation(top 10)

利用DAVID数据库对MSCs成脂分化相关的DEGs进行GO分类,得到P＜0.01的分类317种,其中BP条目238个、CC条目24个、MF条目55个。在生物过程层面,这些成脂分化DEGs主要与转录正调控、DNA模板化、RNA聚合酶Ⅱ启动子转录的正调控和负调控等有关(排名前10位的条目见图2A);在细胞组分层面,这些基因大多参与胞质溶胶、SMAD蛋白复合物、细胞外区域等的组成(排名前10位的条目见图2B);其分子功能的变化主要集中在蛋白质结合、转录因子结合、转录激活子活性、RNA聚合酶Ⅱ核心启动子近端区域序列特异性结合等(排名前10位的条目见图2C)。成脂分化相关DEGs的KEGG通路分析结果显示,其主要涉及癌症通路、调节干细胞多能性的信号通路、AMPK信号通路等47条信号通路(P＜0.01),排名前10位的通路信息见图2D。

图2 成脂分化差异表达基因的GO分类和KEGG通路分析(前10位)Fig.2 Enrichment analysis of GO and KEGG pathway of DEGs in adipogenic differentiation(top 10)

2.3 差异基因所表达蛋白质的相互作用网络

将筛选出的217个成骨分化基因、227个成脂分化基因分别输入STRING数据库,得到成骨分化和成脂分化的网络数据,将所得网络数据分别导入Cytoscape软件构建可视化网络图(图3)。网络数据显示,成骨分化的PPI(protein-protein interaction)网络图共193个节点、2 271条边、平均节点度数23.5;成脂分化的PPI网络图共206个节点、2 126条边、平均节点度数20.6,这些节点之间绝大部分存在着密切的相互作用关系,且有多个节点处于网络的中心,属于整个通路的核心位置。

图3 PPI网络图(A)成骨分化相关蛋白质的PPI网络图;(B)成脂分化相关蛋白质的PPI网络图。Fig.3 PPI network diagram(A)PPI network diagram of osteogenic differentiation related proteins;(B)PPI network diagram of the adipogenic differentiation related proteins.

2.4 关键基因筛选结果

将STRING中所得PPI网络数据导入Cytoscape软件进行分析。节点的度数越高,表示与其相互作用的基因越多,其在定向分化过程中起的作用越大。节点的中介性值越大,说明该基因对于整个网络的调控影响越大[8]。在基因调控网络中,关键基因的调控变化在整个网络中起着至关重要的作用,而节点度数与节点中介性值越高的节点所具有的影响力越大,因此,对于关键基因的筛选可依据节点度数与节点中介性的值。

通过Cytoscape软件对PPI网络数据进行分析发现,在成骨分化网络中,节点度数与节点中介性同时大于或等于平均值(x)的节点有42个,节点度数大于或等于平均值+标准差(+s)的节点有30个,而节点中介性大于或等于平均值+标准差(+s)的节点有15个。具体分析这3种方法筛选出的节点后,本研究筛选关键基因的标准设定为节点度数与节点中介性同时大于或等于平均值+标准差(+s)[9],结果共筛选出14个成骨分化关键基因(表 1,图 4A)。

表1 成骨分化关键基因的筛选Table 1 Screening of key genes for osteogenic differentiation

在成脂分化网络中,节点度数与节点中介性同时大于或等于平均值()的节点有45个,节点度数大于或等于平均值+标准差(+s)的节点有35个,而节点中介性大于或等于平均值+标准差(+s)的节点有25个。以节点度数与节点中介性同时大于或等于平均值+标准差(+s)为筛选标准,结果显示共筛选出20个成脂分化关键基因(表2,图4B)。

图4 关键基因所编码蛋白质的PPI网络图(A)成骨分化关键基因对应的PPI网络图;(B)成脂分化关键基因对应的PPI网络图。Fig.4 PPI network diagram for screening key genes(A)PPI network diagram of key genes for osteogenic differentiation;(B)PPI network diagram of key genes for adipogenic differentiation.

表2 成脂分化关键基因的筛选Table 2 Screening of key genes for Adipogenic differentiation

2.5 成骨分化与成脂分化的共同关键基因

为了评估成骨分化关键基因与成脂分化关键基因之间的不同交互作用,利用Venny 2.1.0对成骨分化和成脂分化中的关键基因进行交互分析比较。结果见图5,在14个成骨分化关键基因和20个成脂分化关键基因中,JUNB、SMAD2、SMAD3、MAPK13(mitogen-activated protein kinase 13)、BMP4(bone morphogenetic protein 4)、SRC 和 RUNX2等7个基因共同参与调控成骨分化与成脂分化。

图5 成骨分化与成脂分化差异表达关键基因的韦恩图7个成骨分化关键基因分别是HSP90、TP53、VEGFA、BMP2、CD44、NOTCH1和COL1A1;13个成脂分化关键基因分别是CCND1、IGF1、IGF1R、SOX2、SIRT1、MAPK8、PPARG、GAPDH、FOXO1、HDAC1、CDH2、TNF 和 SREBF1;7 个成骨分化和成脂分化共有的关键基因分别是 JUNB、SMAD2、SMAD3、MAPK13、BMP4、SRC 和 RUNX2。Fig.5 Venny diagram of the intersection of key genes differentially expressed between osteogenic differentiation and adipogenic differentiationThe 7 key genes for osteogenic differentiation are HSP90,TP53,VEGFA,BMP2,CD44,NOTCH1,COL1A1.The 13 key genes for adipogenic differentiation are CCND1,IGF1,IGF1R,SOX2,SIRT1,MAPK8,PPARG,GAPDH,FOXO1,HDAC1,CDH2,TNF,SREBF1.The 7 key genes shared by osteogenic differentiation and adipogenic differentiation are JUNB,SMAD2,SMAD3,MAPK13,BMP4,SRC,RUNX2.

3 讨论

MSCs的一个重要特征是具有多系分化潜能[10],这一特性使其具有广泛的临床应用前景。分化潜能的稳定保持是发挥MSCs治疗作用的基础。目前,针对MSCs成骨分化和成脂分化的研究以体外实验和体内实验为主,虽然生物信息学分析技术在疾病和基础研究的基因筛查中已被广泛使用,但对于MSCs成骨分化和成脂分化相关基因的系统生物学信息分析较少。本研究旨在通过生物信息学方法筛选MSCs发生成骨成脂分化的关键基因,试图从生物信息角度发现成骨分化与成脂分化关联的新思路,为MSCs异常成骨成脂分化的机制研究提供可探讨的靶点基因,为临床骨质疏松症等成骨脂肪生成失衡相关疾病提供可能的治疗靶点。

文中通过文献挖掘摘选出217个与MSCs成骨分化有关的DEGs或基因产物和227个成脂分化有关的DEGs或基因产物。GO分析显示,成骨分化相关基因或基因产物主要与RNA聚合酶Ⅱ启动子转录的正调控和负调控、成骨细胞分化等过程有关,参与细胞体的组成,介导的分子功能有转录因子结合、生长因子活性等(图1);成脂分化相关基因或基因产物主要与转录正调控、DNA模板化、RNA聚合酶Ⅱ启动子转录的正调控和负调控等过程有关,参与胞质溶胶、SMAD蛋白复合物等的组成,介导的分子功能主要集中在蛋白质结合、转录因子结合等(图2)。KEGG通路富集分析结果表明,成骨分化差异基因主要涉及调节干细胞多能性的信号通路、TGF-β信号通路、Hippo信号通路等(图1D);成脂分化差异基因则主要涉及癌症的途径、调节干细胞多能性的信号通路、AMPK信号通路等(图2D)。当前,成骨分化的分子途径研究以BMP通路和Wnt/β-catenin通路最为经典[11],TGF-β信号通路、Hippo信号通路作为成骨分化的重要途径也已被报道[12],我们的结果提示这两条通路或许可以得到更高的关注。MSCs的多能性与其分化特性紧密联系[13],故推测干细胞多能性信号通路参与其成骨分化和成脂分化的发生发展过程。研究报道,脂肪源性干细胞参与促进肿瘤的发生、发展和复发过程,与癌症具有显著相关性,肿瘤微环境对周围组织脂肪微环境的调控起着重要作用[14],或许这可以解释成脂分化与癌症途径的密切相关性。另外,AMPK作为脂代谢的调控分子,可以通过磷酸化脂质合成过程中的一系列酶对脂肪代谢进行转录调控,这可能与AMPK信号通路在MSCs成脂分化差异基因的通路分析结果中具有显著性有关。

MSCs成骨与成脂分化之间存在一种理论上的反向关系,即MSCs向成骨方向分化的同时会抑制其成脂方向的分化[15]。骨质疏松症是骨与脂肪形成不平衡的代表性疾病之一,常以成骨生成减少、脂肪生成增多为特点[16]。MSCs成骨或成脂分化方向的不同取决于不同的信号通路及其调控的转录因子。许多信号通路遵循成骨和成脂分化之间的反向平衡,即成骨/抗成脂作用(或反之)。这些信号通路包括Wnt/β-catenin信号通路、HH(Hedgehog)信号和Nell-1信号等[15]。成纤维细胞生长因子-2(fibroblast growth factor-2,FGF-2)[17]、TGF-β1[18]和Notch信号通路[19]同样支持这种反向平衡关系。但也有一些信号通路例外,它们具有促成骨/促成脂的双向性。BMP信号通路是参与诱导成骨分化的中枢信号通路之一,在骨形成中起着至关重要的作用[20]。但已有研究发现,BMP通过SMAD1/5/8和MAPK激活可诱导脂肪形成,这可能是由于诱导的剂量和结合受体类型不同[21]。从剂量来看,较低浓度的BMP-2可诱导脂肪生成,而较高浓度的BMP-2则有利于成骨分化[22]。从受体类型来看,BMPR-IA信号转导一般可诱导成脂效应,而BMPR-1B信号通路则可诱导成骨作用[23]。同样,胰岛素样生长因子-Ⅰ (insulin like growth factor 1,IGF-Ⅰ)被发现既有促成骨又有促成脂的作用[24]。

文中筛选出了14个成骨分化关键基因和20个成脂分化关键基因(图4),其中7个基因(JUNB、SMAD2、SMAD3、MAPK13、BMP4、SRC、RUNX2)为两者共有(图5),推测这7个基因在成骨成脂双向调节的通路中起关键作用。我们通过文献分析发现,SMAD2、SMAD3、MAPK13、BMP4 和 JUNB 均与成骨或成脂双向调控的BMP信号通路有关[16,25]。典型BMP信号通路依赖两个丝氨酸-苏氨酸激酶细胞表面BMP受体(BMPRs)结合而产生作用。BMPR-Ⅱ在与BMP配体结合后启动信号传递,随后招募、磷酸化和激活参与MSCs分化的BMPRIA和BMPR-IB,进而磷酸化并激活调节型SMAD1/5/8,调节型SMAD结合物进入细胞核激活或抑制基因的表达[23,25～26]。在非典型BMP信号通路中,MAPK、胞外信号调节激酶(extracellular signalregulated kinase,ERK)和c-Jun N端激酶(c-Jun N-terminal kinase,JNK)等下游BMP信号被激活,影响脂肪分化转录因子过氧化物酶体增殖物激活受体(peroxisome proliferator-activated receptor gamma,PPARγ),从而发生成脂分化[25～26]。

此外,在筛选出的成骨分化关键基因(图4A)中,我们发现血管内皮生长因子(vascular endothelial growth factor,VEGF)作为成骨分化相关因子鲜有提及,但有文献证实其在成骨分化以及软骨分化过程中具有重要的作用[27～28]。研究发现,VEGFA可以调节成骨分化早期骨发育时周围血管的调节因子[28]。这提示VEGFA可作为骨生成障碍疾病骨发育早期的可能治疗靶点。肿瘤蛋白p53(tumor protein p53,TP53)作为肿瘤抑制因子已被人们所熟知,Tataria等[29]发现p53抑癌基因的缺失会促进MSCs成骨分化,这进一步证实TP53在骨肿瘤(例如骨肉瘤)发生过程中对MSCs的成骨分化调控具有重要意义,可作为骨肿瘤的可能治疗靶点。在筛选出的成脂分化关键基因(图4B)中,肿瘤坏死因子(tumor necrosis factor,TNF)除了是促炎性因子,同时也是脂肪组织代谢的重要调节因子[30～31],这提示TNF可能作为炎症相关疾病中MSCs成脂分化的重要靶点发挥作用。

综上所述,本文通过对成骨分化和成脂分化相关基因的生物信息学筛选,得到了调节两者平衡的关键基因,为治疗骨与脂肪生成失衡相关疾病(例如骨质疏松症)探索了可能机制。对于筛选的成骨成脂分化关键基因,本课题组拟进行进一步的生物学实验验证。