基于生物信息学及分子动力学探讨牛膝影响铁死亡治疗骨质疏松机制

肖剑伟，蔡旭，郭粉莲，黄新民，汪荣盛

（1.深圳市福田区风湿病专科医院风湿免疫科，广东 深圳 518000；2.上海市光华中西医结合医院风湿免疫科，上海 200052）

骨质疏松症（Osteoporosis,OP）被定义为一种全身性骨骼疾病，其特征是骨量、强度和微结构的系统性损害，从而增加了脆性骨折的发生率[1]。调查显示，我国绝经后的女性，OP的患病率为48%[2]。骨质疏松症的治疗方式包括生活方式措施和药物治疗。生活方式包括改变饮食结构、补充矿物质和锻炼等。而药物治疗最常用的治疗骨质疏松症的药物包括双膦酸盐、阿仑膦酸盐、唑来膦酸、地诺舒单抗和甲状旁腺激素（PTH）衍生的药物以及激素替代药物等。然而，由于部分药物存在较多不良反应限制了药物的长期使用。

铁死亡是最近发现的一种非凋亡的与脂质过氧化相关的细胞死亡的形式[3]，大量研究显示其参与了神经系统、癌症等疾病的发生发展过程[4,5]。对铁死亡通路相关靶点的调控有望成为治疗疾病的一种新策略。

中药中含有丰富的天然化合物，其在化学结构和生物活性方面也表现出丰富的多样性，并且毒性更小。因此，使用中药治疗骨质疏松成为当前研究热点之一。在大量治疗骨质疏松的中药中，牛膝味苦，性平，归肝肾经，具有强筋骨、补肝肾等功效，十分适合治疗以“肾虚”“脾虚”为病机的骨质疏松。研究显示，用牛膝提取物对去卵巢骨质疏松大鼠治疗后，治疗组的骨密度（Bone mineral density,BMD）较对照组均有所提升，提示牛膝的活性成分能提升去势大鼠的BMD[6]。而王艳等[7]的研究发现，用牛膝治疗去卵巢大鼠能够抑制其骨量丢失，从而达到预防骨质疏松的目的。以上结果显示，牛膝具有的抗OP的疗效，然而其具体机制尚需进一步明确。

本研究以生物信息学为基础，筛选出与OP发病及铁死亡相关的关键基因并进行试验验证，同时结合分子对接和分子动力学探讨牛膝通过影响铁死亡治疗OP的作用机制。

1 材料与方法

1.1 数据检索及数据差异基因分析

以“Osteoporosis”为关键词，在GEO数据库（www.ncbi.nlm.nih.gov/geo） 中检索，检索到GSE56116、GSE2208和GSE7158共3个数据集，见表1。使用R软件（V 4.04）的limma包对编码蛋白的mRNA进行差异分析，取adjust P＜0.05及|LogFC|值大于0.5。

表1 GEO数据集Table 1 Details of the data set

1.2 预测牛膝作用靶点

通过中药系统药理学技术平台（TCMSP）[8]检索牛膝化学成分，将口服生物利用度（oral bioavailability,OB）≥40%及类药性（drug likeness,DL）≥0.18作为活性小分子化合物的筛选条件，将具有较高活性的化合物进行筛选提取并将小分子化合物作用靶点转化为基因名并使用Cytoscape 3.7.1软件绘制“药物化合物靶基因”相互作用网络图。

1.3 铁死亡相关基因及靶点获取

通过FerrDb数据库（www.zhounan.org/ferrdb）[9]，检索铁死亡相关基因。将上述结果与1.1、1.2步骤得到的靶点取交集，获得牛膝影响OP铁死亡途径的潜在治疗靶点并绘制韦恩图。

1.4 试验验证

1.4.1 材料试验采用小鼠单核细胞RAW264.7细胞株购自深圳市豪地华拓生物科技有限公司。ALOX12一抗（Abcam，ab168384）和二抗、GPX4一抗和二抗（Abcam，ab125066）；内参-actin（CST，#4970）。琼脂糖凝胶（Solarbio，P1200）。BioTraceTMNT硝酸纤维素膜（Millipore，IPVH00010）。DMEM培养基、胎牛血清和胰蛋白酶（GIBCO，12430054、10100147、12605036）；黄芩素（Sigma-Aldrich，491-67-8）、MTT试剂（Sigma-Aldrich，CT02）、二甲基亚砜（DMSO）（Sigma-Aldrich，C6164）、核因子-B受体活化因子配基（receptor activator of NF-B ligand,RANKL）（Sigma-Aldrich，207621-35-0）、巨噬细胞集落刺激因子（macrophagecolonystimulating factor,M-CSF）（Sigma-Aldrich，SRP3221）和青链霉素（sigma-Aldrich,P4333）。

1.4.2 蛋白质印迹法（Western blot,WB）检测潜在治疗靶点根据上述步骤得到的差异基因ALOX12，通过WB验证其在OP患者中的表达情况。选取福田区风湿病专科医院风湿免疫科符合《中国老年骨质疏松症诊疗指南(2018)》OP诊断标准，排除合并有恶性肿瘤、近6个月内服用过激素类药物以及严重肝肾功能障碍的患者，共纳入2020年4~6月确诊的OP患者4例，其中男性2例、女性2例，年龄52~65岁，平均（58.8±5.4）岁。按照OP患者性别、年龄比例选择4名深圳市福田区风湿病专科区院体检健康者作为对照组，其中男性2名、女性2名，年龄50~63岁，平均（57.3±5.9）岁。2组之间年龄差异无统计学意义（t=3.757，P=0.624）。本研究通过医院伦理委员会审查（伦理审查编号FS202001012），研究对象均自愿参与，签署知情协议书。ALOX12一抗和二抗（Abcam，ab168384）；内参-actin（CST，#4970）。琼脂糖凝胶（Solarbio，P1200）。BioTraceTMNT硝酸纤维素膜为美国（Millipore，IPVH00010）产品。

抽取静脉血5 mL，样本加入5 mL已经混有苯甲基磺酰氟（PMSF）的蛋白裂解液（每1 mL裂解液加入10L PMSF）进行蛋白定量。取适量的蛋白样品上样，进行SDS变性10%聚丙烯酰胺凝胶电泳（SDS-PAGE），室温摇动封闭1 h。将膜按一抗蛋白和-actin蛋白的印迹位置剪开，置于含对应一抗的Blotto溶液中，4℃摇动作用过夜。TBST溶液中摇动漂洗后加入二抗室温作用，摇动漂洗后将膜置于Western LightningTMChemiluminescence Reagent显色剂中1 min。常规曝光后进行显影、定影处理。以-actin为内参，使用Image J（V 1.8.0）软件分析各条带灰度值与内参的灰度值比值作为目的蛋白的相对表达量。

1.4.3 WB检测黄芩素对RAW264.7细胞GPX4、ALOX12表达的影响RAW264.7细胞在含有100U/mL青霉素和100g/mL链霉素、10%FBS、10%M-CSF、10%RANKL的DMEM中培养，并置于温度为37℃、CO2浓度为5%的培养箱中。黄芩素用微量的DMSO溶液溶解，然后用DMEM培养基将黄芩素溶液稀释至试验所需浓度，使DMSO的最终浓度不超0.1%。

取对数生长期的RAW264.7细胞，常规消化收集细胞，计数并调整细胞密度至1 104个/孔，接种到96孔板中，当细胞达到80%~90%融合度时，根据既往文献黄芩素体外使用配制浓度报告[10]，将试验空白对照组及不同浓度的黄芩素组（5mol/L、10mol/L和20mol/L）处理48 h。

收集上述各组细胞用裂解液提取细胞蛋白，进行蛋白定量。取适量的蛋白样品上样，进行SDS变性10%聚丙烯酰胺凝胶电泳（SDS-PAGE），室温摇动封闭1 h。将膜按一抗蛋白和-actin蛋白的印迹位置剪开，置于含对应一抗的Blotto溶液中，4℃摇动作用过夜。TBST溶液中摇动漂洗后加入二抗室温作用，摇动漂洗后将膜置于WesternLightningTMChemiluminescence Reagent显色剂中1 min。常规曝光后进行显影、定影处理。以-actin为内参，使用Image J软件分析计算GPX4及ALOX12蛋白的相对表达。

1.5 分子对接

明确ALOX12在OP患者中的表达情况后，通过RCSB网站检索得到靶点的晶体结构，Pubchem数据库中获取中药活性小分子结构，使用AutoDock Vina对蛋白及小分子进行半柔性对接，选取结合能最小的结果作为最佳配体。

1.6 分子动力学

使用Gromacs进行分子动力学模拟[11]。通过ATB网站（atb.uq.edu.au/）[12]生成配体小分子位置限制文件。通过Gromacs内置命令将小分子及靶蛋白转换为拓扑结构文件。使用Gromos54a7_atb力场及SPC水分子模型，创建大小为1.2的十二面体水盒子并加入6个钠离子来中和系统电荷。采用最速下降法使系统能量最小化。使用V-rescale进行步长为2 fs，总时间长为100 ps的NVT系综平衡，然后使用Berendsen算法在1 atm压力下进行100 ps的NPT系综平衡。对蛋白质中的重原子（非氢原子）施加位置限制后进行总时长为20 ns的动力学模拟，步长设定为2 fs。根据模拟得到的结果在去除体系的周期性边界后提取20 ns时的结构作为最终复合物构象。使用Gromacs内置命令计算复合物的均方根偏差（rootmeansquaredeviation,RMSD）、溶剂可及表面（solvent accessible surface area,SASA）、蛋白回旋半径（Gyrate,GR）和氢键。使用g_mmpbsa软件处理GROMACS生成的结构（.tpr）和轨迹（.xtc）文件，并计算小分子与蛋白的结合自由能[13]。使用VMD[14]、PyMOL[15]、LigPlot+软件[16]对数据进行可视化处理。

1.7 统计分析

使用R软件（V 4.04）进行统计学分析。临床样本的蛋白表达水平使用t检验进行比较。P＜0.05具有统计学意义。

2 结果

2.1 差异分析结果

结果显示与正常组对比，GSE56116有760个mRNA存在表达差异，其中上调295个，下调465个；GSE2208有114个mRNA存在表达差异，其中上调78个，下调36个；GSE7158有276个mRNA存在表达差异，其中上调152个，下调124个。

2.2 牛膝作用靶点

经过筛选后得到牛膝有效成分9个。9个有效成分的作用靶点共106个，将其转化为官方基因名（Gene Symbol）。使用Cytoscape 3.7.1将药物、化合物和靶基因绘制“药物化合物靶基因”相互作用网络图，见图1。该互作网络包含了116个节点和222条边。

图1 “药物 化合物 靶基因”相互作用网络图Fig.1” Drug compound target gene” interaction network

2.3 牛膝治疗OP靶点

根据1.3步骤获得铁死亡相关基因97个，与2.1、2.2步骤所得结果取交集（图2），得到潜在通过影响铁死亡治疗OP的的靶基因1个（ALOX12）。

图2 GSE56116、GSE22085、GSE7158、铁死亡相关基因及牛膝作用靶基因交集的韦恩图Fig.2 Venn diagram of intersection genes among the GSE56116,GSE22085,GSE7158,ferroptosis-associated genes and target genes of Achyranthes bidentata

2.4 试验结果

2.4.1 蛋白质印迹法结果以-actin为内参，在OP患者全血中，与对照组相比，ALOX12蛋白表达明显上升（正常组：0.63±0.17，OP组：0.99±0.34，t=19.38，P＜0.05），见图3。

图3 WB结果显示OP组患者的ALOX12蛋白表达水平升高Fig.3 WB results showed the levels of ALOX12 protein expression increased in the OP group

2.4.2 黄芩素干预对RAW264.7细胞GPX4、ALOX12表达结果与对照组比较，不同浓度的黄芩素组RAW264.7细胞中ALOX12的蛋白表达水平明显降低，差异有统计学意义，见图4。而GPX4蛋白表达无明显差异，表明黄芩素能抑制RAW264.7细胞中ALOX12的表达，而对GPX4蛋白表达无影响。

图4 不同浓度黄芩素处理后对RAW264.7细胞GPX4、ALOX12表达的影响（**<0.001）Fig.4 Effects of different concentrations of baicalein on the expression of GPX4 and ALOX12 in RAW264.7 cells（**<0.001）

2.5 分子对接结果

根据2.2及2.3步骤筛选得到小分子化合物黄芩素（baicalein）作用的靶点及铁死亡相关基因ALOX12。分子对接结果显示，黄芩素与ALOX12的Arg189、His365形成氢键，与Leu366、Thr396、Glu369、Asn400、Thr662、Leu193和Val190形成疏水键，提示黄芩素与ALOX12有较强的结合能力，可能通过与其结合后抑制ALOX12活性，见图5 A。

图5 分子对接结果Fig.5 Molecular docking results

2.6 分子动力学模拟结果

黄芩素与ALOX12的复合体在整个动力学运行过程中逐渐稳定，RMSD波动值小于0.1 nm（图6 A）。复合物的溶剂可及表面呈逐渐缩小趋势（图6 B），复合物体系的蛋白回旋半径（Rg）在模拟的过程中也逐渐缩小（图6 C），提示蛋白与小分子结合后蛋白逐渐收敛并使得结构趋于稳定。复合物体系氢键数量随着模拟时间变化仍保持2个以上的连接（图6 D），同时，与模拟前结构，活性小分子对比位置无明显变化，提示两者结合稳定（图5 B）。

图6 蛋白小分子复合物体系的RMSD、SASA、RG和氢键数量Fig.6 RMSD,SASA,RG and number of hydrogen bonds of protein-small molecule complex systems

结合能计算结果显示，黄芩素与ALOX12的总结合自由能为（107.893kJ/mol±16.671kJ/mol）。结果显示，黄芩素与ALOX12主要通过范德华力（154.481kJ/mol±10.798 kJ/mol）和静电势能（68.847 kJ/mol±11.810 kJ/mol），而极性溶剂化作用则对结合起抑制作用（101.368 kJ/mol±19.941 kJ/mol）。

3 讨论

在骨质疏松症发病的各种原因中，骨稳态的破坏是其发病的关键即由成骨细胞介导的骨形成减少和由破骨细胞驱动的骨破坏增加，而铁死亡在这个过程中发挥着重要的作用。铁死亡的特征包括细胞内铁过载和活性氧类（Reactive oxygen species,ROS）的大量堆积。而既往研究发现，铁过载和氧化应激是骨质疏松发病的重要因素之一，高度提示铁死亡参与了骨质疏松的发病过程。骨质疏松症患者存在明显的铁过载[17]，其通过氧化应激损伤促进骨质疏松症的发生和发展[18]。研究显示，铁过载可影响成骨细胞的分化和活性[19]，并且可以诱导成骨细胞的前体细胞骨髓间充质干细胞的自噬和凋亡[20,21]，而使用降低铁过载的药物可以改善骨质疏松情况[22]。在氧化应激的状态下，过载的ROS作为诱导铁死亡模式启动的重要因素。ROS在大多数情况可以被机体内抗氧化酶清除。铁超载可通过Fenton反应导致ROS形成增多[23]。ROS参与了成骨细胞的分化，过多ROS生成使正常的骨重建过程发生异常，使得骨形成减少而骨吸收增加，从而导致骨质疏松[24,25]，降低ROS水平能够抑制破骨细胞分化[26]。另一个报道显示[27]，在去卵巢小鼠中，降低ROS水平可以抑制破骨细胞的产生从而阻止雌激素缺乏引起的骨质疏松。试验显示在卵巢切除术后小鼠中，Loureirin B可通过抑制ROS活性来抑制破骨细胞的产生从而防止小鼠的骨丢失，进一步证实了ROS在骨质疏松发病中的作用[28]。YANG等[29]发现内皮细胞分泌的外泌体是细胞间通讯的重要介质，其通过抑制铁蛋白吞噬作用依赖的铁死亡，外泌体逆转了糖皮质激素诱导的成骨细胞对成骨的抑制作用。以上多个研究结果表明，铁离子和氧化应激在骨质疏松过程中具有重要的作用。而由铁介导的氧化应激损伤导致的铁死亡是近来新发现的细胞死亡形式，由铁离子介导的铁死亡可能导致骨质疏松发生及病情进展。调控铁死亡相关的靶点有望成为治疗骨质疏松潜在手段。

铁死亡作为一种由铁过载和脂质过氧化驱动的细胞死亡方式，其实质是细胞内脂质氧化物出现代谢障碍，在铁离子催化作用下发生代谢异常，导致脂质过氧化从而诱导细胞死亡。铁死亡有明显区别于细胞凋亡的形态学特征，如线粒体膜破裂增加。脂质过氧化作为一种自由基驱动的反应，主要影响细胞膜中的不饱和脂肪酸。脂质合成代谢失调会导致大量多不饱和脂肪酸积累，同时会消耗细胞内的谷胱甘肽（GSH），导致更多的过氧化脂质产生，从而提高了细胞对铁死亡的敏感性。其中，多不饱和脂肪酸（PUFA）的氧化是由脂氧合酶（ALOXs）途径介导的。LOXs是一种非血红蛋白含铁酶，主要功能是过氧化不饱和脂肪酸[30]。人类LOX家族共有6个成员（ALXOE3、ALOX5、ALOX12B、ALXO12、ALOX15B和ALOX15）。在铁死亡途径中，铁可以增加LOXs的活性，随后LOXs通过介导PUFA过氧化产生花生四烯酸等，促进铁死亡。其中，ALOX12可以直接过氧化细胞膜上的磷脂形成磷脂过氧化物，然后在亚铁离子的催化下，脂质过氧化物与其他脂质发生链式氧化反应，导致铁死亡的发生[31]。多个研究显示，ALOX12的多态性与骨质疏松发病呈正相关[32,33]，提示ALOX12与骨质疏松存在密切关系。一项研究显示，黄芩素能够与铁结合，抑制血小板中ALOX12的表达[34]。Dai等[35]发现，黄芩素能够抑制ALOX12，减少铁死亡，从而改善CCl4诱导的小鼠急性肝损伤。Bo Chu等[36]通过CRISPR/Cas9技术来生成敲除ALOX12的U2OS细胞的亚克隆，敲除了ALOX12基因的U2OS细胞铁死亡受到抑制。由此可推断，ALOX12是黄芩素的作用靶点之一，通过靶向ALOX12抑制其活性，可以影响铁死亡过程，见图7。

图7 牛膝治疗OP的机制图Fig.7 Mechanism diagram of Achyranthes bidentata for the treatment of OP

祖国传统医学根据OP的表现将其归属于“骨痹”、“骨痿”和“骨枯”等病范畴。肾藏精，精血养骨髓，髓枯则筋痿，肾虚能够诱发骨质疏松。从西医学角度来看，骨髓干细胞能够分化为成骨细胞及破骨细胞。铁死亡会导致炎症水平升高，产生的各种炎症因子如肿瘤坏死因子-等，能够抑制骨髓间充质干细胞成骨分化的生物效应，使骨形成减少导致骨量丢失，这与中医肾虚导致的“骨枯”具有异曲同工之处。姚新苗认为OP多属本虚标实，本虚责于脾肾而标实责于瘀血，在运用补肾壮骨、健脾益气治疗的同时应当辅以活血化瘀之药，常选用赤芍、川芎、丹参和牛膝等[37]。陈雯等归纳中药治疗骨质疏松用药频次，活血化瘀药作为重要治法，在活血兼有养血功能的中药处方中，当归使用频次为44次，骨碎补使用频次为28次，牛膝使用频次为26次，排在活血化瘀药物的前三位[38]，提示牛膝在治疗OP中有着广泛的应用。现代研究证实，牛膝及含有牛膝的复方药物治疗骨质疏松疗效明确，如含有牛膝的左归丸、无敌丹胶囊具有滋阴补肾、填精益髓的作用，在临床已被证实能够改善OP患者症状[39,40]。牛膝含有大量的活性小分子成分，其中作为有效成分的黄芩素具有抗炎抗氧化作用[41]。研究显示黄芩素可以抑制RANKL诱导的破骨细胞的形成以及抑制破骨细胞的骨吸收功能[10]，提示黄芩素是治疗骨质疏松的潜在药物分子。

本研究显示，在OP患者中，ALOX12蛋白表达较健康组明显升高，提示其可能参与了OP的发病。作为牛膝有效成分之一的黄芩素，通过分子对接及分子动力学结果显示，其能通过多个氢键及疏水键与ALOX12紧密对接并形成稳定结构。在动力学模拟过程中，整个复合物体系能够快速收敛并达到稳态，说明黄芩素能够通过影响ALOX12蛋白的结构稳定性从而发挥抑制其活性的作用。试验结果显示，不同浓度的黄芩素可以抑制RANKL诱导的RAW264.7细胞ALOX12蛋白的表达，而不影响GPX4蛋白的表达。以上结果提示ALOX12可能作为牛膝治疗OP的潜在靶点，抑制ALOX12的活性，达到治疗骨质疏松的目的。

综上所述，本研究以生物信息学、网络药理学、分子对接及分子动力学结合试验研究验证牛膝有效活性成分黄芩素可能通过靶向ALOX12影响OP铁死亡途径从而达到治疗OP的作用，初步阐明了牛膝有效活性成分治疗OP的潜在作用机制，为进一步探讨牛膝通过铁死亡途径治疗OP提供了理论依据。然而，本研究的结果未通过动物模型研究证实，需要进一步通过使用体内外试验来验证。