大家好！今天小编和大家分析一篇23年3月发表在Scientific Reports（IF:4.6）杂志的文章《Bioinformatics analysis of rheumatoid arthritis tissues identifies genes and potential drugs that are expressed specifically》。
研究涉及骨科相关疾病的坏死性凋亡机制，因为坏死性凋亡是一种独特的调节性细胞死亡模式。然而，坏死性凋亡相关基因在类风湿性关节炎（RA）中的作用尚未得到很好的描述。作者分别从基因表达综合（GEO）、京都基因和基因组百科全书（KEGG）数据库和基因组富集分析（GSEA）下载了RA相关数据信息和坏死性凋亡相关基因。鉴定了113个与RA相关坏死性凋亡相关的基因，这些基因与细胞因子介导的信号通路、坏死性凋亡和程序性坏死密切相关。随后，通过LASSO、RF和SVM-RFE三种机器学习算法，在48个坏死性凋亡相关差异基因中鉴定出FAS、MAPK8和TNFSF10为关键基因，这些关键基因对区分RA患者和健康对照具有良好的诊断能力。根据功能富集分析，这些基因可能调节多种途径，如B细胞受体信号传导、T细胞受体信号通路、趋化因子信号通路和细胞因子-细胞因子受体相互作用，并在RA中发挥相应作用。此外，预测了48种针对关键基因的靶向药物和31种化学结构式。在ceRNA网络中，关键基因与复杂的调控关系相关。根据CIBERSORT分析，FAS、MAPK8和TNFSF10可能与RA患者免疫微环境的改变有关。研究发展了诊断有效性，并为RA的机制提供了见解。在临床应用之前，需要进一步的研究来测试其对RA的诊断价值。
背景
慢性自身免疫性疾病类风湿关节炎(RA)可引起肢体畸形和肢体残疾，但其在疾病和疾病进展方面存在个体差异。成年人群中RA患病率约为1%，主要影响中年至老年妇女。RA除了是世界上导致残疾的主要原因外，还对公共卫生和经济产生重大影响。RA的临床特征以滑膜细胞浸润和存活为特征，导致关节软骨损伤、滑膜增厚、骨侵蚀和局部炎症细胞浸润，导致慢性炎症。RA的发病机制涉及肿瘤坏死因子α、白细胞介素6和细胞间基质等炎性细胞因子激活破骨细胞，导致RA患者和软骨的破坏。针对RA发病机制的临床治疗主要有生物制剂、非甾体类抗炎药和糖皮质激素以及手术治疗，但对于RA滑膜仍缺乏有效的治疗靶点。因此，寻找RA的潜在治疗靶点是至关重要的。
坏死性凋亡是一种受调控的细胞死亡的新模式，发现于1985年，并于2005年命名。许多因素会引发坏死性凋亡，包括：肿瘤坏死因子(TNF)家族成员和其他致病微生物。半胱氨酸蛋白酶8 (caspase 8)的活性抑制减少细胞凋亡但激活坏死性凋亡。RIP1和RIP3结合形成一种称为“坏死体”的复合体可诱导RIP3自磷酸化。导致细胞膜破裂并最终破坏细胞。最近的研究表明，坏死性凋亡机制在骨科相关疾病(包括骨关节炎)的发病机制中起着重要作用。在骨关节炎软骨中，坏死性凋亡标志物RIPK1和RIPK3的表达在肿瘤坏死因子ɑ等炎症因子刺激后上调，并抑制坏死性凋亡，显著减少关节软骨破坏，延缓骨关节炎进展，有效证实坏死凋亡在骨关节炎进展中起重要作用。关于RA和坏死性凋亡的报道很少，其作用机制仍不清楚。为了从坏死性凋亡的角度进一步了解RA的发病机制，发现潜在的靶点很重要。
本研究利用基因表达综合数据库(Gene Expression Omnibus, GEO)对免疫浸润进行系统分析。因此，评估坏死性凋亡在RA发展中的作用。同时，研究坏死性凋亡在RA中的潜在功能机制，关键基因与浸润免疫细胞的关系，以更好地了解RA进展中潜在的分子免疫过程。
方法：
数据下载和整理：为了开展这项研究，使用关键词“类风湿关节炎”从Gene Expression Omnibus (GEO, http://www.ncbi.nlms.nih.gov/geo)检索RA数据。GEO为该研究提供了所有数据，因此不需要知情同意和伦理批准。
差异表达分析：对坏死性凋亡相关基因进行差异分析。
功能富集分析：通过研究基因本体（GO）富集分析和KEGG富集分析，了解NRDEGs在RA中的生物学功能和意义。
确定诊断RA的最佳关键基因：NRDEG使用三种机器学习算法（LASSO，SVM-RFE和RF）求和。使用最小绝对收缩和选择算子（LASSO）和支持向量机（SVM）对RA的诊断标志物进行分类。
基因组富集分析：对关键基因进行GSEA分析。
免疫细胞浸润分析
关键基因药物预测与ceRNA网络的构建：筛选关键基因相关靶向药物；进行mRNA-miRNA-lncRNA的ceRNA网络构建。
关键基因的验证。
研究结果
坏死性凋亡相关差异基因的鉴定
在实验数据集（GSE55235和GSE55457）中，包括23例RA患者滑膜组织和20例健康对照滑膜组织。差异分析结果显示，113个坏死性凋亡相关基因中，共有48个发现RA滑膜组织样本与健康对照样本之间存在差异，其中36个基因上调，12个基因下调。聚类热图显示了坏死性凋亡差异基因（NRDEGs）在RA样本和健康对照样本中的表达模式（图1A），NRDEG之间的相关性如图所示（图1B）。FAS与MAPK8、BIRC、JAK2和STAT1呈正相关，MAPK8与PYCARD、CASP1和TNFSF10呈负相关。IFNA1和IFNA8仅与FAD和PLA2G4C呈正相关。
NRDEG的富集分析
富集了RA分化基因和NRDEGs（包括GO富集分析和KEGG分析），以了解NRDEGs在RA中的生物学功能和相关途径。GO富集分析结果显示，类风湿性关节炎中T细胞活化、质膜外侧和酰胺结合均显著富集，NRDEG与细胞因子介导的信号通路、膜筏和细胞因子受体结合密切相关。其中，“膜筏”的分子功能普遍富集（图A，C）。此外，KEGG结果显示，“EB病毒感染”和“人T细胞白血病病毒1感染”在RA中显著富集，而NRDEGs与“坏死性凋亡”和“NOD样受体信号通路”等通路密切相关（图B，D）。综上所述，富集分析证据表明，NRDEGs可能通过调节“细胞因子介导的信号通路”、“膜筏”、“坏死性凋亡”和“NOD样受体信号通路”在RA的发展中发挥重要作用。
NRDEG中的3个基因被确定为RA的关键基因
鉴于RA患者和健康个体之间的明显差异，评估了NRDEG的诊断潜力。通过三种机器学习算法（LASSO、SVM-RFE和RF）筛选训练数据集中有意义的关键基因，以区分RA患者。在LASSO算法中，从48个NRDEG中选择了12个基因（图3A，B）。同时，在SVM-RFE算法中，识别了19个NRDEG（最大精度= 0.935，最小RMSE = 0.065）（图3C，D）。RF算法用于确定NRDEG的基因重要性,鉴定了8个重要性大于0.9位的基因（图3E-F）。最后，结合了三种机器学习算法的结果，其中三个基因被确定为关键基因（FAS，MAPK8和TNFSF10）。为了研究三个关键基因是否可用于区分RA样本和健康对照，绘制了ROC曲线。如图3H所示，关键基因的AUC均大于0.85。此外，还绘制了验证数据集的ROC曲线，以进一步验证关键基因的诊断潜力。结果表明，验证数据集中3个关键基因的AUC值均大于0.8，表明这3个关键基因在区分RA样本和正常对照样本方面具有一定的准确性和特异性（图3)上述证据表明，关键基因在区分RA患者与健康对照样本方面具有较高的准确性和特异性。
关键基因的GSEA分析
由于不知道RA中关键基因的潜在功能，因此进行了GSEA分析以了解它们在RA中的作用和意义。与关键基因富集相关的前六种途径如图4所示。综合分析表明，关键基因主要参与免疫应答（B细胞受体信号通路、T细胞受体信号通路、趋化因子信号通路和细胞因子-细胞因子受体相互作用）和多种疾病通路（病毒性心肌炎、张力性关节炎、扩张型心肌病和肥厚型心肌病）。
免疫细胞浸润分析
除了发病机制复杂外，RA还受免疫系统的影响很大。为了研究RA患者和健康对照之间免疫微环境的差异，采用了CIBERSORT算法。如图5A所示RA样本中浆细胞、T细胞CD8、T细胞滤泡辅助细胞和巨噬细胞M1的比例高于健康对照样本，而T细胞CD4记忆静息、NK细胞活化、单核细胞和肥大细胞活化在健康对照样本中表达更多。此外，还分析了特征基因高低表达与免疫微环境的关系。FAS和MAPK8与它们呈显著正相关，而TNFSF10在免疫细胞T细胞滤泡辅助细胞、浆细胞、CD8T细胞、CD4记忆激活T细胞和γδ细胞中与它们呈负相关。FAS和MAPK8之间存在显着的负相关，并且在TNFSF10和树突状细胞休息时发现了显着的正相关（图5B-D）。基于这些证据，FAS，MAPK8和TNFSF10可能强烈影响RA患者的免疫微环境。
关键基因靶向药物预测与ceRNA网络构建
DGIdb数据库可能靶向药物的发现，Cytoscape软件在图6中显示了结果。靶向签名基因的药物48种，其中FAS靶向19种，靶向MAPK8的39种，TNFSF10未发现相关靶向药物。值得注意的是，在预测的靶向药物中，TANZISERTIB，CC-401和BENTAMAPIMOD是已知的MAPK8抑制剂。此外，通过DrugBank数据库检索了相关靶标药物的结构式，如图7所示。在预测的48种靶向药物中共检索到31个化学结构式。下一步是使用starBase和Miranda数据库构建ceRNA网络。如图8所示网络由238个节点（3个关键基因、87个microRNA和148个lncRNA）和264个边缘组成。调节FAS的lncRNA的数量包括59个lncRNA，它们竞争结合has-miR-650，has-miR-625-5p，has-miR-766-3p和has-miR-338-3p。在ceRNA网络中，20个lncRNA可以与has-miR-93-3p，has-miR-1238-3p和has-miR-214-3p结合，以调节关键基因MAPK8。has-miR-1238-3p和has-miR-214-3p调节关键基因MAPK8。关于TNFSF10，共有26个miRNA参与调控。
验证浓缩特征基因的表达
GSE77298数据集用于验证RA样品中关键基因的表达。实验数据与3个关键基因（FAS、MAPK8和TNFSF10）表达趋势的比较是一致的。其中，RA患者MAPK8（p=0.0018）的表达低于健康对照样本，而FAS（p=0.0096）和TNFSF10（p=0.0018）在RA样本中的表达量更高（图9）。
