脊髓损伤是神经外科一类重症疾病，可分为创伤性和非创伤性脊髓损伤，前者较为常见，常由外部物理因素引起^[1], 后者常由肿瘤、缺血或先天性疾病引起^[2]。脊髓损伤通过神经元损害、肌肉废用等途径引起骨骼肌萎缩^[3-4]。目前已有学者对脊髓组织损伤后继发的线粒体功能障碍、氧化应激激活、肌肉蛋白质合成与降解等过程进行了相关研究^[5-6], 但其病理、生理的复杂性限制了该病的治疗进展。本研究对脊髓损伤后肌肉萎缩患者的肌肉组织进行生物信息学分析，并探讨该病的发病机制，以期为潜在治疗方法提供理论依据，现报告如下。

1.   资料与方法

1.1   基因芯片的筛选

选择GEO数据库中GPL570[HG-U133_Plus_2]Affymetrix Human Genome U133 Plus 2.0 Array平台检测，编号为GSE21497的芯片进行分析。此芯片从招募的10位脊髓损伤肌肉萎缩的患者中获取骨骼肌样本，其中包括脊髓损伤患者男9例和女1例(6例为四肢瘫痪, 4例为截瘫)，平均年龄44岁。入选标准: 18岁以上患者; 格拉斯哥昏迷评分大于13分者; 无肌肉挤压伤、低氧损伤、系统性败血症、全身性炎症者; 自身免疫性疾病及恶性肿瘤者。分别在脊髓损伤后第2、5天使用5mm Berstrom活检针从10位患者的股外侧肌获取肌肉样本，取下后立即置入液氮中冻存，储存于-80℃待测。以患者第2天的股外侧肌样本的基因表达谱数据为对照组，以第5天样本的基因表达谱数据作为实验组。

1.2   芯片数据处理与差异基因提取

使用R软件中affyPLM、RColorBrewer包进行芯片数据分析，得出相对对数表达图，观察各组样本数据表达强度。然后调用R软件中affyPLM、affy包通过RMA算法对芯片进行背景校正。利用R软件中limma包，设置差异基因筛选条件: logFC>1或logFC < (-1), 且adj. P. val < 0.05。将提取的差异基因导入OmicShare(https://www.omicshare.com/)在线网站制作出火山图、热图。

1.3   差异基因富集分析

将筛选出的差异基因导入Metascape在线网站，进行GO分析、KEGG信号通路分析，以柱状图表现富集结果，获取基因相关的生物学过程(BP)、细胞组分(CC)、分子功能(MF)。

1.4   构建蛋白质互作用网络

将筛选出的差异基因导入Metascape在线网站，进行蛋白互作网络构建。利用MCODE算法检测出相互作用结果紧密的蛋白质模块中的关键差异基因。

2.   结果

2.1   脊髓损伤后肌肉萎缩患者肌肉组织中的差异基因

使用R软件得出相对对数表达图，每个样品的中心均较接近纵坐标0, 结果表明样本平行试验一致性较强，见图 1A。样本背景校正后，设置差异倍数大于2及校正后P值小于0.05, 筛选出294个差异基因。2组样本中共有8个上调表达差异基因, 286个下调表达差异基因，见表 1。将差异基因导入OmicShare在线网站制作热图和火山图，见图 1B、图 1C。

图  1  芯片质量评估与差异基因表达图

A: 脊髓损伤肌肉萎缩患者基因芯片中样本的相对对数表达图; B: 患者第2天肌肉样本与第5天肌肉样本间差异最显著的前100个基因聚类热图，红色为表达水平高的基因，绿色为表达水平低的基因; C: 火山图，红色的点为患者肌肉样本中第5天相对于第2天上调表达差异基因，绿色的点为下调表达差异基因，黑色的点为无表达差异基因。

下载: 全尺寸图片 幻灯片

表  1  脊髓损伤肌肉萎缩患者第2、5天肌肉样本的差异基因

上调(8个基因)

LOC646903、GADD45A、PRRT3-AS1、MSTN、MAP3K7CL、RCSD1、E2F8、LSMEM1

下调(286个基因)

TNFRSF18、PRB4、CA9、SEMA5B、GCK、KLRG1、LINC00868、KIAA0040、C16orf71、FAM223B、LOC283194、CSF3、LILRB1、LOC101928111、CA14、LOC101929897、LGALS9、SHISA3、CCDC64B、RTCA、C11orf84、RASSF4、ASPRV1、PORCN、ANKRD34A、MEDAG、LOC101060424、TGM7、BEX4、ST3GAL1、KIF3B、PEX11A、PPARGC1A、F3、XAGE-4、FOXP1-IT1、PARG、HOXC11、CLEC2L、RETSAT、TNFSF13、MYL3、PROCR、ANKRD2、TP53INP2、NEDD4L、NID1、MFAP5、PRKCD、DAXX、FCN3、MYH7、MAP1B、HSD17B11、BBS10、PTPRU、ANGPTL7、SLC27A4、AOX1、FHL5、FABP7、METTL7A、KDR、FAM110D、NDUFC1、RC3H1、KCNB1、SMAP2、CDC37L1、RNF25、NMNAT1、SMAD3、ARID5B、CHN2、SSPN、WBP1L、CYP4B1、CPAMD8、ASAH1、LMCD1、UGP2、ZNF438、C1orf162、PSD、TFPT、ACSM5、COL6A3、RPN1、FAM166B、PDE7A、PIK3CB、DIO2、TSPYL2、SGMS1、MLXIPL、LIX1L、SLC7A2、NMRK2、NCOA3、FEM1C、SMS、MOCOS、VLDLR、FOXO3、LIFR、DDIT4、GPT2、HOMER1、GPX7、MTFP1、EPHB1、TMEM45A、CHAC1、ALOX5AP、SLC19A2、GOT1、PTPRG、SORT1、ITGAM、RASSF1、SLC16A9、LAMB1、P2RY14、MYOM2、PLIN3、ANG、AKR1B10、PRCC、ZBTB16、PIBF1、PKD2、MFGE8、ELF1、LZTS3、HCLS1、THBS2、SLC38A2、METRNL、ECI2、ACSL1、ACTN2、ZNF385B、HSF2、LINC00312、FIBIN、SPSB1、PDK4、SAMD13、CA3、ITGB1BP2、PRKAG3、MT1M、SIRT1、SRPX、NPTX2、ENKD1、C1orf198、LTBP2、YAP1、GPCPD1、FKBP5、FSD2、CPM、CLIC4、FLJ36848、KIAA1161、CNST、ANKRD27、LPIN1、MYOM3、DNAJA4、BLM、MMP3、SLC25A1、RAB23、PRELP、CASQ2、IL1R2、PANK1、PIK3R1、COL5A2、INHBB、ITGA10、ADAMTS5、COL4A1、VIT、RNASE4、CLIC6、ACACB、EBF1、APCDD1、COQ10A、DNAJB5、TSPAN5、BZW1、MYLK4、FBN1、RND3、FOSL2、FAM65B、ROBO1、ANO5、ENPP5、DHRS7C、AGTR1、DHCR24、ANKMY2、GLUL、S100A9、PCYOX1、CD163、IGFBP6、WIPI1、FABP3、RBP4、NFIL3、LINC00948、FAM46C、MIOS、SLC38A1、C5orf27、TNKS1BP1、CNN3、CSPP1、PLAT、C20orf166、SOD2、TMEM208、VCL、SLC2A4、ARHGEF3、GADD45G、LYVE1、SAMHD1、SMAD7、MYL2、GALNT15、SIMC1、SNAI2、FOXO1、HSPA2、SH2D1B、C1QC、TIMP4、VSIG4、MRC1、GGT7、ARG2、RASSF5、SMAD9、IRX3、C1QA、CALML6、SERPINA5、MIR646HG、MYH2、LOC101927391、ABRA、OLFML2A、KLF5、LINC01405、DAAM2、TCAP、IGFBP3、MGST1、RASD1、LINC00844、MAOA、CIDEC、F13A1、GYS1、SLC25A33、C1QB、FBP2、MYOZ2、MAFB、C12orf75、LRRN3、LOC643659、PLTP、CTGF、ECHDC3、FAM179A、PEBP4、TNNC1、TNNI1

下载: 导出CSV 

| 显示表格

2.2   差异基因的GO、KEGG分析

GO富集中, CC中差异基因主要集中于肌纤维组分、肌质、异构SMAD蛋白复合体、肌球蛋白复合物、肌膜、细胞器外膜。分子MF中差异基因主要集中于肌肉的结构组分、肌动蛋白结合、生长因子结合、整联蛋白结合、碳水化合物结合、β-连环蛋白结合、热休克蛋白结合、碳酸脱水酶活性、脂肪酸连接酶活性、阴离子跨膜转运蛋白活性。BP中差异基因主要集中于肌肉器官发育过程、细胞衰亡的正向调控、细胞增殖、分化的负向调控、能量平衡、脂肪细胞分化、细胞对外部刺激的反应、对氮化合物的反应、对酸性化学物质的反应、活性氧代谢过程、辅因子生物合成、突触修剪等。见图 2。

图  2  差异基因的GO富集和KEGG图

A: GO富集于CC图; B: GO富集于MF图; C: GO富集于BP图; D: KEGG分析图。

下载: 全尺寸图片 幻灯片

KEGG信号通路主要包括: 补体系统信号通路、FoxO信号通路、黏着斑信号通路、氮素代谢信号通路、精氨酸和脯氨酸代谢信号通路、过氧化物酶体增殖物激活受体信号通路、谷胱甘肽代谢信号通路、铁死亡信号代谢通路、烟酸和烟酰胺代谢信号通路、半乳糖代谢信号通路等。

2.3   蛋白质相互作用分析及关键差异基因获取

蛋白质互作用网络中, 4个蛋白相互作用较丰富模块中的关键基因共有18个: TNNI1、TNNC1、TCAP、ACTN2、MYL3、MYL2、DNAJA4、DNAJB5、GYS1、UGP2、AKR1B10、GOT1、C1QA、C1QB、C1QC、CASQ2、DAXX、SLC2A4, 见图 3。

图  3  蛋白质互相作用网络图

A: 差异基因的蛋白互作网络图; B: 4个蛋白质互作丰富模块中的关键差异基因。

下载: 全尺寸图片 幻灯片

3.   讨论

脊髓损伤导致肌肉萎缩，亦使全身多系统继发严重功能障碍^[6]。研究^[7-8]表明，在脊髓损伤后6周内，骨骼肌平均横截面积比对照组低18%~46%, 且身体总脂肪率增高。肌肉萎缩可加重胰岛素抵抗并促进糖尿病发展^[9], 同时下肢肌肉收缩力量减弱，引起循环血量不足，导致运动耐力降低、呼吸功能减弱等^[10-11]。目前对于该病的研究仍以动物实验模型为主，本研究选择GEO数据库中编号为GSE21497的患者骨骼肌样本进行研究，更具客观意义。

通过分析筛选出294个差异基因，经GO分析富集于CC中的差异基因主要存在于肌质、肌膜、肌球蛋白复合物中; MF中差异基因主要参与肌肉构成以及β-连环蛋白、热休克蛋白等物质的结合过程; BP中差异基因主要参与肌肉发育、突触修剪、活性氧代谢等过程。研究^[5]发现，肌肉蛋白质的合成和降解失衡、氧化应激反应、自噬等过程在该病的进程中起重要作用。CHEN Z等^[12]发现，损伤的肌肉组织中FoxM1转录因子表达量下降，引起Wnt通路中β-连环蛋白信号过度激活，导致肌肉萎缩和再生不良。通过促进FoxM1转录可显著上调肌肉干细胞中Apc的表达，降低β-连环蛋白水平，拮抗肌肉损伤。CRISTOFANI R等^[13]发现，小热休克蛋白B8分子通过介导错误折叠蛋白的自噬可缓解肌肉萎缩。

KEGG信号通路中差异基因主要集中在补体系统、FoxO信号通路、谷胱甘肽代谢等通路。研究^[14]发现，肌肉损伤后使用重组趋化因子CCL5可促使补体系统C3a和C3aR信号表达增强，促进肌肉再生。FoxOs通路参与体内众多重要生物过程，如自噬、ROS解毒、DNA修复等^[15-16]。学者^[17]发现，抑制FoxO1、FoxO3、FoxO4基因的激活，可缓解自噬、蛋白降解导致的神经性肌肉萎缩。NINFALI C等^[18]学者提出, ZEB1可抑制FoxO3转录活性，从而抑制肌肉萎缩基因Fbxo32和Trim63启动子的转录。ABDULLAH M等^[19]通过建立犬肌肉营养不良模型发现，损伤后组织中胶原蛋白纤维化与肌肉萎缩在相关，而脯氨酸/精氨酸代谢途径可能是炎症相关胶原蛋白合成、肌纤维化的核心过程。

通过PPI分析及MCODE算法筛选出18个关键差异基因，其中部分基因已被报道与脊髓损伤后肌肉萎缩过程紧密相关。肌钙蛋白在骨骼肌的收缩过程中起重要作用，脊髓损伤后肌肉组织内源性肌浆网中Ca²⁺含量降低约10%，与慢速骨骼肌型肌钙蛋白I 1 (TNNI1)的表达相关，通过上调骨骼肌中TNNI1的表达，可增强脊髓损伤者康复运动的抗疲劳性，改善疾病预后^[20-22]。肌肉特异性糖原合酶1(GYS1)缺失可导致葡萄糖代谢功能和运动能力受损^[23], 促进GYS1上调可提高肌肉收缩能力并减缓萎缩。补体级联反应在神经损伤区域被激活，保护神经免受感染、缓解损害，其具体机制尚未阐明^[24-25]。研究^[26]证实, C1q、C3、C5等补体在神经损伤后可调控神经再生。此外，肌源性转录因子MYOD1可协同E-box调控肌联蛋白帽(TCAP), 促进对骨骼肌新陈代谢^[27]。其余关键差异基因与肌肉萎缩的研究较匮乏，仍需进一步实验探索。

目前，治疗脊髓损伤后肌肉萎缩的方法包括手术、药物、基因工程等^[28-31]。早期手术减压可缓解组织继发性缺氧、缺血，降低细胞损伤程度^[32-33]。河豚毒素、苯妥英钠等钠通道阻滞剂可维持离子稳态平衡，缓解轴突病变^[34-35]。促红细胞生成素上调SDF-1α与G蛋白偶联的CXCR4受体等趋化因子的表达，可募集骨髓间充质干细胞，增强组织抗凋亡能力，改善神经功能等^[36-37]。甲泼尼龙等激素可上调抗炎因子、降低组织氧化应激反应、减轻水肿等^[38]。学者^[39-41]提出，激素治疗不仅在运动或神经功能恢复上与对照组无明显差异，还可能会加重感染、胃肠道出血等症状。胚胎干细胞移植和多能干细胞诱导可调节炎症反应，改变微环境，分泌营养因子，使神经再生^[42-44]。

综上所述，本研究所筛选出的差异基因在脊髓损伤病程中涉及微环境变化、神经元变性坏死、肌肉组织纤维化、能量获取、线粒体功能障碍等重要过程。临床可对上述基因的相关机制进行深入研究，以为组织损伤后神经肌肉系统的生长、代谢、凋亡等方面的治疗提供理论依据，从而改善患者预后。