Typing characteristics of metabolism-related genes in osteoporosis
骨质疏松症中的代谢相关基因的类型特征
发表期刊：Front Pharmacol
发表日期：2022 Sep 15
影响因子：5.988
DOI: 10.3389/fphar.2022.999157
一、研究背景
骨质疏松症是一种全身性的代谢性骨病，发病率不断上升。它的特点是骨量减少，骨的微观结构恶化，骨脆性和骨折的风险增加。骨质疏松症主要表现为骨质疏松性骨折，包括身体不同部位的骨折，尤其是髋部骨折。
人类包括多种生物，如动物，主要通过碳水化合物、脂肪和蛋白质来摄取潜在的能量，但也有涉及几种酶的特定生化级联进行分解。在生理上，许多代谢途径中的基因和路径对许多细胞的代谢功能至关重要。因此，这些基因和途径的失调或失衡会导致细胞功能紊乱和各种代谢性疾病。尽管有大量证据表明骨病和代谢紊乱之间可能存在联系，但这些代谢网络中的具体参与者和分子相互作用仍然不清楚。同时，组学技术促进了对各种临床疾病的分子变化的开发和探索。同时，利用基因组学、转录组学、蛋白质组学和代谢组学等多组学技术，是促进临床疾病正确治疗的关键。
二、材料与方法
1、数据来源
1） 以 "osteoporosis "和 "Homo sapiens "为关键词：GSE56814（42个正常和31个骨质疏松症样本）、GSE56815（40个正常和30个骨质疏松症样本）、GSE35959（9个正常和5个骨质疏松症样本）和GSE7429（10个正常和10个骨质疏松症样本）；GSE35959和GSE7429微阵列数据被用作外部验证数据集
2） 以前的研究已经确定了2752个代谢相关的基因，编码所有已知的人类代谢酶和转运器
2 、分析流程
1） 聚类分析：使用集合相似性网络融合和共识聚类算法（SNF-CC）来观察骨质疏松症患者的基因表达模式和聚类识别；根据聚类结果绘制多组学热图，并用于从多组学数据中识别、验证和可视化分子疾病亚型
2） 新陈代谢相关的骨质疏松症亚型的差异性基因表达谱：R语言limma包进行亚型之间的基因比较分析，并使用Venn diagram方法得出共同的代表性DEGs
3） 通路富集分析和蛋白质图谱的可视化：基因富集；使用Proteomaps将差异分析中的蛋白质数据可视化
4） 筛选与骨质疏松症代谢相关的特征基因：利用排名靠前的mrmr基因构建了一个SVM分类器，并应用增量特征选择（IFS）方法来确定作为特征基因的最佳数量
5） 特征基因表达的数据验证：验证特征基因的表达情况，并通过绘制小提琴图将数据可视化
6） 免疫相关物分析：xCell
三、实验结果
01 - 骨质疏松症亚型的差异性基因表达谱
作者从Gene Expression Omnibus（GEO）数据库下载了骨质疏松症患者的基因表达和临床数据，包括GSE56814、GSE56815、GSE35959和GSE7429。首先，对GSE56814和GSE56815数据集的表达谱进行主成分分析，揭示基因表达数据，如图2A,B所示。其次，通过整合相似性网络融合和共识聚类算法（SNF-CC），在骨质疏松症患者队列中计算CPI（蓝线）和差距统计（红线），确定了最佳的聚类数量。其中，当k=3时，共识矩阵热图仍能保持清晰明确的边界，说明样本的聚类是稳定和稳健的（图2C）。将这三种骨质疏松症亚型定义为C1、C2和C3。最后，通过主成分分析确定了三种骨质疏松症亚型，如图2D所示。
此外，作者用代谢基因热图绘制了代谢途径，显示高BMD（中位绝对偏差）组与C1亚型聚在一起，低BMD组与C3亚型聚在一起（图3A），绘制了不同BMDS之间三种亚型的百分比条形图。其中，C1亚型在高BMD骨质疏松症患者中占较高比例，C3亚型在低BMD骨质疏松症患者中占较高比例，C2亚型在不同BMD骨质疏松症患者中占类似比例。如图3B所示，使用R语言limma包进行亚型之间的基因差异分析，并绘制维恩图以确定每个亚型的代表性基因，其中1501个在三个亚型之间共享。同时，对于不同的亚型C1、C2和C3，分别确定了1351、8和158个差异表达的基因（图3C），并绘制了火山图。
02 - 代谢相关差异基因的多组学富集图谱
为了可视化不同亚组的特征基因的富集信息的差异，并利用其他生物信息学工具确认富集分析中得到的结果，作者将三个亚组的特征基因分别上传到网络上的互动软件proteomaps。在创建的蛋白图可视化矩形中，整个区域被划分为代表顶级类别的彩色编码多边形，顶级类别区域又被细分为疾病区域亚类、功能相关蛋白共享共同区域亚类和相关基因共享共同区域亚类，如图4所示。C1亚型在代谢方面关于生物合成、氨基酸代谢、脂质和类固醇代谢的富集，主要蛋白是K00718，主要基因是FUT。C2亚型在中央碳代谢、糖酵解途径、K00002蛋白、基因AKR1A1方面代谢丰富。C3亚型代谢富集于多糖的生物合成、糖代谢途径，主要蛋白为k07968、k01197、k03909，主要基因为B4GALT3、HYAL2、TFPI。
03 - 骨质疏松症中代谢相关基因的筛选和差异性表达
作者可视化了不同骨质疏松症亚型的差异基因，如图5A所示。根据筛选出的前160个mRMR基因，通过构建100个SVM分类器并应用增量特征选择（IFS）方法来确定作为特征基因的最佳数量，如图5B所示。当基因数量为10个时，MCC峰值为0.386，这10个基因被认为是与骨质疏松症的代谢有关的特征基因，包括GPR31、GATM、DDB2、ARMCX1、RPS6、BTBD3、ADAMTSL4、COQ6、B3GNT2和CD9。随后，根据10个候选基因在不同的骨质疏松症亚型中的表达情况绘制了可视化的小提琴图。GPR31、DDB2、RPS6、ADAMTSL4和COQ6在C1和C2亚型中表达，C1>C2。GATM、ARMCX1、BTBD3、B3GNT2和CD9在C1和C2亚型中表达，C1<C2。在C2和C3亚型中GATM和ARMCX1的表达，C2>C3。ADAMTSL4和COQ6在C2和C3亚型中表达，C2<C3。有趣的是，在C1与C3的比较中没有统计学意义，如图5C所示。
04 - 特征基因表达数据的验证
根据GSE35959和GSE7429数据集中GPR31、GATM、DDB2、ARMCX1、RPS6、BTBD3、ADAMTSL4、COQ6、B3GNT2和CD9基因表达情况，我们绘制了小提琴图进行数据可视化的表达验证。其中，在综合数据中GPR31的表达，高BMD<低BMD，如图6A所示。GPR31在GSE7429数据中表达，高BMD>低BMD，DDB2在GSE7429数据中表达，高BMD>低BMD，如图6B所示。在GSE35959数据中GPR31和DDB2表达，正常组＜骨质疏松组，而ARMCX1和ADAMTSL4在GSE35959数据中表达，正常组＞骨质疏松组，如图6C所示。
05 - 免疫相关因素分析
最后，用xCell来推断64种免疫和基质细胞类型的GPR31、GATM、DDB2、ARMCX1、RPS6、BTBD3、ADAMTSL4、COQ6、B3GNT2和CD9基因。如相关图所示，这些细胞分别由16、15、14、10、17、13、27、19和24个代表细胞组成。每个基因的代表细胞主要包括NKT（自然杀伤性T细胞）、iDC（不成熟树突状细胞）；CLP（普通淋巴祖细胞）、cDC（Xonventional树突状细胞）；Th1细胞（1型T辅助细胞）、浆细胞、CD4 Tem（CD4效应记忆T细胞）、类开关记忆B细胞；肥大细胞、巨噬细胞M2。CD8 T细胞，B细胞；mv内皮细胞（微血管内皮细胞）；类开关记忆B细胞，NKT（自然杀伤T细胞）；Th1细胞（1型T辅助细胞），CD4 Tem（CD4效应记忆T细胞）；成纤维细胞，嗜碱性细胞，MSC（间充质干细胞）；血小板，巨核细胞，HSC（造血干细胞），如图7所示。
四、结论
基于对骨质疏松症患者基因表达的聚类分析和机器学习，作者确定了不同的代谢相关亚型和骨质疏松症的特征基因，这将有助于为骨质疏松症的代谢相关发病机制提供新思路，为后续研究提供新方向。
