骨关节炎的关键通路与免疫细胞浸润模式及中药预测分析

李玲玲 黄海量

山东中医药大学康复医学院 山东济南 250355

骨关节炎(osteoarthritis, OA)是以关节软骨退变、继发性骨赘形成、软骨下骨重塑和滑膜炎症为特征的慢性骨关节病，临床主要表现为关节周围疼痛、肿胀、功能受限和关节畸形等[1]，严重影响着患者的健康和生活质量。研究显示，OA的发病率女性高于男性，农村高于城市[2]。随着人口预期寿命的延长及社会进入老龄化，其患病率和致残率也逐渐升高。预计到2030年将有近4亿人患有OA[3]，越来越多的患者承受肉体乃至精神上的痛苦与折磨，同时也给许多家庭和社会带来了巨大的精神压力和沉重的经济负担[4]。

有研究认为，OA的发生可能与年龄、性别、关节损伤、肥胖和代谢等多种因素有关，其病理生理机制及其复杂，多个生物学过程和多条信号通路参与了其损伤进程，实验研究难以系统、整体分析信号通路和免疫浸润模式[5]。基于此，本研究利用生物信息学的方法进行通路富集分析和免疫浸润模式分析，并预测对OA有潜在治疗作用的中药，为OA发病机制和药物治疗提供了一定参考和借鉴。

1 资料与方法

1.1数据来源 以GEO中的GSE55235基因表达芯片为研究对象，共包括OA滑膜组织10例，正常滑膜组织10例。

1.2筛选差异表达基因 使用GEO2R在线分析工具(https://www.ncbi.nlm.nih.gov/geo/geo2r/)，设定校正后P<0.05且|logFC|≥2，筛选OA滑膜组织和正常滑膜组织之间的差异表达基因(Differentially expressed genes, DEGs)。

1.3GO和KEGG分析 采用DAVID在线数据库(https://david.ncifcrf.gov/)对筛选出的DEGs进行基因本体功能富集分析(GO)和京都基因与基因组百科全书通路富集分析(KEGG)，选择物种为“智人”，设定P<0.05。对GO富集前10的生物过程(biological process, BP)、细胞组成(cellular component, CC)和分子功能(molecular function, MF)三个部分，以及KEGG富集前20的信号通路使用微生信平台(http://www.bioinformatics.com.cn/)对结果进行可视化。

1.4PPI网络构建与核心基因的筛选 将上述DEGs上传到STRING数据库(https://string-db.org/)，设置物种为“智人”，蛋白最低相互作用分值为0.40，并隐藏离散的节点。将节点信息以“TSV”格式保存并导入Cytoscape 3.8.2，构建蛋白-蛋白相互作用网络(protein-protein interaction, PPI)，并利用cytoHubba插件按照度值筛选核心靶点。

1.5免疫细胞浸润模式分析 利用CIBERSORT算法对GSE55235基因表达芯片中人类免疫细胞亚型的表达矩阵进行去卷积分析，计算22种免疫细胞的相对比例，设定P<0.05过滤样本。使用R语言barplot函数和pheatmap包绘制正常组与OA组中免疫细胞组成的柱状图和热图，使用corrplot包对OA组中的免疫细胞进行相关性分析并绘制共表达热图，使用vioplot包分析正常组和OA组免疫细胞浸润差异并绘制小提琴图。

1.6核心基因的中药预测 将筛选出的DEGs与Coremine Medical数据库(https://www.coremine.com/)相互映射筛选对OA有潜在治疗作用的中药。按照口服生物利用度≥30%且类药性≥0.18从TCMSP数据库(https://tcmsp-e.com/)筛选其活性成分和靶点信息，使用UniProt(https://www.uniprot.org/)对靶点名称进行标准化处理，Cytoscape 3.8.2软件用于构建“中药-活性成分-作用靶点”网络图。

2 结果

2.1DEGs 共筛选出235个DEGs，其中上调DEGs101个和下调DEGs134个，火山图见图1。

图1 DEGs的火山图

2.2GO和KEGG富集分析 DEGs共富集到194个GO条目(P<0.05)，其中生物过程(BP)143个，细胞组成(CC)21个和分子功能(MF)30。分别选取BP、CC和MF富集前10条目绘制条形图结果，见图2。BP涉及免疫反应、炎症反应、细胞粘附和凋亡过程等；CC涉及质膜、胞外区和细胞外间隙等；MF涉及蛋白质同源二聚化活性、转录因子活性、序列特异性DNA结合和相同蛋白结合等。同时，KEGG富集分析得到37条信号通路(P<0.05)，富集前20条通路的气泡图如图3所示。与OA相关的信号通路主要有TNF信号通路、破骨细胞分化、MAPK信号通路、NF-kappa B信号通路、趋化因子信号通路、NOD样受体信号通路、Toll样受体信号通路等。

图2 DEGs的GO富集分析

图3 DEGs的KEGG富集分析

2.3PPI网络和核心基因 采用STRING在线数据库构建了DEGs的PPI，见图4，包括203个蛋白节点和751次相互作用，平均度值为7.4。利用Cytoscape软件中cytoHubba插件进行网络拓扑学分析，以度值排前10的DEGs作为PPI中的核心基因，见图5。

图4 DEGs的PPI网络

图5 核心差异表达基因

2.4OA免疫细胞浸润分布组间差异 通过CIBERSORT反卷积方法计算并最终获得14个可信样本，横坐标轴左边4个为对照组，其余10个为OA组。OA组与正常组免疫细胞浸润差异分析结果见图6， OA组织中B细胞和肥大细胞的浸润程度与对照组存在显著差异，见图7(P<0.05)。

图6 OA组和对照组22种免疫细胞的浸润组成

图7 OA组与正常组免疫细胞浸润差异

2.5骨关节炎免疫细胞相关性分析 免疫细胞占比的相关性分析结果,见图8，蓝色代表负相关，红色代表正相关。结果显示，M1型巨噬细胞与活化的CD4+记忆性T细胞(r=0.72)、静息的树突状细胞呈强正相关(r=0.7)，静息的CD4+记忆性T细胞与M1型巨噬细胞、静息的树突状细胞呈较强正相关(r均=0.69)。活化的肥大细胞与静息的肥大细胞呈强负相关(r=-0.7)，初始B细胞和记忆B细胞呈较强负相关(r=-0.65)。

图8 OA免疫细胞相关性分析

2.6中药预测结果 通过Coremine Medical数据库预测具有潜在治疗OA作用的中药，见表1，其中郁金、姜黄和黄芩是频次前3的中药，故选取3味中药构建中药-活性成分-作用靶点网络，见图9。图中椭圆为郁金、姜黄和黄芩的有效活性成分，中间黄色长方形为三者的共同作用靶点，绿色长方形为非共同作用靶点。郁金、姜黄和黄芩的有效活性成分共同参与ESR1、F10、HSP90AA1、PRKACA、PGR、PTGS1、CHRM3、CHRM1、NCOA2、ADRA1A、CHRM2、CHRNA7、NR3C2、PPARG、AKR1B1等基因的调控。

表1 具有潜在治疗OA作用的中药

图9 中药-活性成分-作用靶点网络

3 讨论

随着高通量基因组学技术的发展，从微阵列获得的基因表达谱已被广泛用于筛查与各种疾病进展相关的基因和通路[6]。本研究利用GSE55235基因表达芯片数据筛选出OA患者和健康对照组之间的DEGs共235，其中IL6、IL1β、MMP-9、PTGS2、VEGFA、CXCL8、JUN、MYC等为核心差异基因。IL6和IL1β是最为常见的促炎细胞因子，有研究显示细胞因子IL-6、IL-1β可促进软骨细胞的炎症和凋亡，抑制聚集蛋白聚糖的表达，刺激基质金属蛋白酶(MMPs)的表达，从而促进OA疾病的发生发展[7]。MMPs是一个锌依赖性的蛋白水解酶家族，主要通过分解胶原和蛋白多糖等参与细胞外基质的降解，在OA关节软骨基质及软骨细胞破坏的病理过程中起着重要作用[8]。MMP-9又可以通过增加促炎介质如TNF-α 和IL-1的表达加剧局部炎症，导致软骨加速退化和软骨下骨过度生长等[9]。此外，MMP-9和PTGS2 (Cox-2)在激活时可以通过降解和裂解胶原蛋白II和软骨特异性蛋白多糖参与OA软骨的降解[10]。血管内皮生长因子(vascular endothelial growth factor，VEGF)的表达在晚期OA升高，VEGF可上调软骨细胞中MMP-1、MMP-3、MMP-13的表达，下调关键细胞外基质成分(包括聚集蛋白聚糖和II型胶原蛋白)的表达，加速软骨退变的进程[11]。VEGF还可以诱导促炎细胞因子IL-6、CXC趋化因子IL-8的产生，VEGF除了促血管生成活性外，还具有强大的促炎症作用[12]。JUN是激活蛋白1(AP-1)的一个亚基，后者参与软骨细胞中MMP-13的表达，Gao X等[13]研究认为，JUN参与了破骨细胞的分化。MYC通过灭活MAPK信号通路抑制细胞凋亡并促进软骨细胞增殖，防止OA的疾病进展[14]。

为进一步深入了解OA相关的生物途径和信号通路，对筛选出的235个DEGs进行GO和KEGG富集分析。GO富集分析结果表明，参与OA的生物过程主要涉及免疫反应、炎症反应、细胞粘附和凋亡过程等，这些生物过程无疑与OA疾病的发生发展密切相关。KEGG 分析结果表明，与OA相关的信号通路主要有TNF信号通路、破骨细胞分化、MAPK信号通路、NF-kappa B信号通路、趋化因子信号通路、NOD样受体信号通路、Toll样受体信号通路等。TNF信号通路可以介导KOA炎症反应，促进炎性细胞因子的产生，导致滑膜增生和炎症。有研究发现，NF-κB和MAPK均可以通过协调多层信号网络导致基质降解酶、促炎细胞因子和炎症介质的产生[15-16]。OA患者破骨细胞诱导的骨吸收和成骨细胞诱导的重塑之间的平衡会逐渐失调，增强的破骨细胞介导的骨吸收也是OA的潜在病理机制[17]。Toll样受体和NOD样受体能各自通过信号转导产生级联反应以刺激不同的效应分子的表达，共同作用于炎症过程[18]。

免疫细胞浸润结果显示，OA组织中B细胞和肥大细胞的浸润程度与对照组存在显著差异(P<0.05)。肥大细胞是免疫系统的重要传感器和效应细胞，可以诱导炎症介质的释放和先天免疫细胞的募集[19]，并通过抑制成骨细胞活性(如 IL-1、TNF)和/或促进破骨细胞生成(如组胺、TNF、IL-6)来调节或诱导骨代谢[20]。免疫细胞相关性分析显示，M1型巨噬细胞与活化的CD4+记忆性T细胞(r=0.72)和静息的树突状细胞呈强正相关(r=0.7)，静息的CD4+记忆性T细胞与M1型巨噬细胞、静息的树突状细胞呈较强正相关(r=0.69)。巨噬细胞可以根据环境刺激获得不同的表型，M1型巨噬细胞主要促进炎症的发生，可以分泌IL-1β、IL-18、IL-6等炎症相关细胞因子，进而导致软骨退化和骨赘形成，而M2型巨噬细胞主要负责抗炎和组织修复[21]，所以能使M1转变为M2表型的药物也将是治疗OA的有效方法。OA属中医“痹症”、“骨痹”、“痿症”范畴，证属本虚标实、本痿标痹[22-23]。通过Coremine Medical数据库预测了对OA有潜在治疗作用的中药，其中郁金、姜黄和黄芪具有较高频次，并构建了中药-活性成分-作用靶点网络。有研究显示黄芩的有效成分黄芩素可以通过保护软骨下骨、抑制血管生成和滑膜增殖来缓解OA[24]。郁金和姜黄均为姜黄属植物，其中根茎为姜黄、块根为郁金，具有抗炎、抗氧化、抗癌等多重生物活性[25]。

综上所述，本研究以GEO中的GSE55235基因表达芯片为研究对象，利用生物信息学的方法筛出235个DEGs，通过富集分析，发现OA疾病的发生发展与免疫反应、炎症反应、细胞粘附、凋亡等生物过程和TNF信号通路、破骨细胞分化、MAPK信号通路、NF-kappa B等信号通路密切相关。其次，通过CIBERSORT算法和Coremine Medical数据库分析了免疫细胞浸润模式，并预测了具有潜在治疗OA作用的中药，为OA发病机制和药物治疗提供了新的见解和思路。