作者：孙娜 史传道 张荣强 张帆 李星慧 张章

单位：陕西中医药大学公共卫生学院流行病与卫生统计学教研室

摘要：目的  利用基因芯片公共数据库中RNA标本芯片数据，分析绝经后骨密度降低程度不同的骨质疏松症患者外周血基因表达谱及分子网络变化。方法  从基因芯片公共数据库(Gene Expression Omnibus，GEO)中下载美国奥马哈20名绝经后妇女的外周血RNA标本基因芯片数据，利用基因云、基因大数据分析(Gene-Cloud of Biotechnology Informs，GCBI)实验平台、GenClip2.0、Sytoscape 3.5.1等软件，筛选差异表达基因。结果  差异基因主要涉及MAPK调节、细胞内信号转导的调控、细胞运动、免疫效应等生物学功能。基因网络中以MAPK3、ARRB2、STK11、USP9X、MYO1C、GTF2H1等基因degree值较大，基因MAPK3、ARRB2、STK11、MYO1C、GTF2H1在不同骨密度人群中表达差异有统计学意义(P＜0.05)，其中基因MAPK3、ARRB2、STK11和MYO1C诊断价值有统计学意义(P＜0.05)。结论  MAPK3、ARRB2、STK11和MYO1C基因为蛋白网络核心节点。MAKP3和ARRB2基因与骨质疏松的发生关系密切。

关键词：骨质疏松；基因表达；生物信息学

骨质疏松症(osteoporosis，OP)是一种主要以骨量降低和骨组织微结构破坏为特征的全身性骨骼疾病，现已成为当今全球仅次于心血管疾病、受害人群最多、最具危害性且尚未能治愈的慢性疾病。绝经后骨质疏松症(postmenopausal osteoporosis，PMOP)为Ⅰ型原发性骨质疏松，是一类由雌激素缺乏引起的，以骨骼强度减弱为特点的骨骼疾病，其特征为骨密度(bone mineral density，BMD)下降、骨脆性增加。随着老龄化社会的到来，中国已进入骨质疏松高发期，严重威胁着人类的身心健康，给家庭、社会与经济带来了沉重的负担。BMD是诊断骨质疏松和预测骨质疏松骨折的重要手段。近年来，针对该病的研究越来越多，但对于基因谱及生物通路网研究还不充分，探明绝经后骨质疏松分子机制，成为亟待解决的问题。本研究通过对基因芯片公共数据库(Gene Expression Omnibus，GEO)中收录的一组绝经后妇女外周血基因芯片数据资料进行生物信息学分析，探讨PMOP患者分子网络改变及其可能涉及的分子生物学功能，为预测该病提供一定的理论参考依据。

资料与方法

资料来源

以“Postmenopausal Osteoporosis"为关键词在NCBI数据库GEO中检索,得到由XIAO等在2008年6月提交、2016年7月更新的基因芯片数据集GSE7429(https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE7429)。该数据集来自于美国奥马哈20名绝经后妇女，按照脊椎或髋关节Z评分标准，将其分为低BMD组(Z评分 < －0.84)和高BMD组(Z评分＞0.84)的外周血组织RNA标本，由Creighton University Medical Center公司利用Affymetrix Human Genome U133A Array基因表达芯片平台检测获得。此数据样本病例的纳入标准和排除标准均有严格详细规范。该数据集包含22 283个数据资料。]]>

分析方法

数据处理及差异基因(DEGs)筛选：将数据集GSE7429导入基因云、基因大数据分析(Gene-Cloud of Biotechnology Informs，GCBI)实验软件平台，进行(均数=0，标准差=1)标准化处理，并通过(P＜0.01，q＜0.01)过滤后行两独立样本t检验，筛选出低BMD和高BMD患者的差异表达基因。

差异基因的GO富集分析：将筛选出的差异基因列表上传至人类基因功能和网络分析软件GenCliP 2.0，调节各类分析的主要参数，分析两组数据差异表达基因的生物学通路、分子功能、生物学过程预测分析。Fold Enrichment代表富集倍数。

http://www.networkanalyst.ca/)软件进行分析，根据具体分析结果适当调节可信度和附加节点参数，分析其所编码的蛋白质之间的相互作用联系图，找出处在关键节点的蛋白质。]]>

关键基因的诊断能力评价：采用GraphPad Prism软件对两组关键基因绘制受试者工作特征曲线(receiver operating characteristic curve，ROC)，AUC值代表ROC下方面积大小。通常，AUC的值介于0.5到1.0之间，越接近于1，指标的诊断效果越好，取值在0.7~0.9间，诊断有一定的准确性。

结果

DEGs分析结果

采用GCBI实验平台对标准化数据进行分析，并将差异表达基因绘制聚类热图，可以看出差异表达基因主要聚类为2类，能够很好将低BMD和高BMD人群区分开(图 1)。

从火山图可以反映出总体基因的表达情况(图 2)。图中每个点代表一个基因，颜色用以区分基因是否差异表达，橙色的点代表差异表达基因，蓝色的点代表没有差异表达的基因。

差异基因GO富集

按基因差异表达倍数(Fold-Change)大小排序，将差异较大的前250个基因导入GenCliP 2.0在线软件，对差异基因所涉及的生物学功能进行分析，结果表明这些靶基因主要涉及MAPK调节、细胞内信号转导的调控、细胞运动、免疫效应等生物学功能(图 3)。

蛋白质—蛋白质网络分析

采用NetworkAnalyst在线软件分析这些差异基因编码的蛋白质间的作用，整个网络以MAPK3(degree=113)、ARRB2(degree=42)、STK11(degree=38)、USP9X(degree=31)、MYO1C(degree=29)、GTF2H1(degree=28)等基因degree值较大，为此蛋白相互作用网络中心节点，删除这些节点蛋白后，网络结构涣散(图 4)。

关键基因表达比较

进一步对6个关键基因做两独立样本t检验，与低BMD组相比，基因MAPK3、ARRB2、STK11、MYO1C在高BMD组高表达，而USP9X和GTF2H1为低表达，MAPK3、ARRB2、STK11、MYO1C和GTF2H1基因在两组间表达差异有统计学意义(P＜0.05)，而基因USP9X表达在两组间差异无统计学意义(P＞0.05)(表 1)。

差异基因诊断能力评价

对5种关键基因分别制作ROC曲线(图 5)，根据AUC值大小评价指标的诊断价值，依次为基因MAPK3、ARRB2、STK11、MYO1C和GTF2H1，其中基因GTF2H1诊断能力无统计学意义(表 2)。

讨论

随着人类寿命的延长和人口结构老龄化的加剧，骨质疏松症已成为世界普遍关注的公共卫生问题，其在世界常见病、多发病中居第7位，且每年骨质疏松症的患病率呈上升趋势。如何预防和治疗骨质疏松症，已成为全球医药卫生事业亟待解决的问题。因此，本课题组采用功能基因组学预测骨质疏松症的发生。

本研究通过对20名不同水平BMD的骨质疏松症患者外周血RNA标本基因芯片数据进行分析，并对差异基因的生物学功能分析，结果表明这些靶基因主要涉及MAPK调节、细胞内信号转导的调控、细胞运动、免疫效应等生物学功能。基因网络中以MAPK3、ARRB2、STK11、USP9X、MYO1C、GTF2H1等基因degree值较大，其中基因MAPK3、ARRB2、STK11、MYO1C、GTF2H1在不同BMD人群中表达差异有统计学意义，基因MAPK3、ARRB2、STK11和MYO1C诊断价值有统计学意义。

雌激素通过作用于雌激素受体促进成骨细胞增值和破骨细胞凋亡，发挥抗骨质疏松的作用。MAPK信号通路是存在于真核细胞内的一类丝/苏氨酸蛋白激酶，可以被细胞外信号或刺激激活。雌激素通过激活MAPK通路促进骨形成，对抗骨质疏松的发展。绝经后妇女由于体内雌激素水平下降，从而更易导致骨质疏松发生，本研究结果与此相符。有研究表明，ARRB2基因通过与其他受体相互作用抑制炎症反应。骨质疏松是导致骨折发生的高危因素，因此高BMD组骨质疏松患者该基因表达增高。MYO1C基因是一种被广泛表达的运动蛋白，它将肌动蛋白细胞骨架与细胞膜连接起来，参与细胞核DNA的转录、胞外分泌和膜贩运等过程。STK11蛋白激酶作为一种抑癌因子在肿瘤发生过程中的多个方面起作用，但有关MYO1C和STK11在骨质疏松方面的研究还很少，这为我们后期的研究提供了一定的思路。

综上，本研究采用生物信息学研究方法对不同BMD骨质疏松症患者基因芯片数据进行分析，差异基因主要涉及MAPK调节、细胞内信号转导的调控、细胞运动、免疫效应等生物学功能，MAKP3和ARRB2基因与骨质疏松的发生关系密切，STK11和MYO1C基因在骨质疏松上的研究还有待于进一步探索。