骨质疏松症(osteoporosis，OP)是一种骨量降低、骨微结构退化，骨脆性增加、骨强度减低的全身代谢性骨病[1]。患者发生骨折的风险较常人高，累及脊柱、髓关节和腕关节，在一定程度上使之丧失运动能力，预后较差，导致患者生活质量严重下降[2]。OP在中医学理论中属于“骨痹”的范畴。肾虚是OP的主要原因，肾中精气的盛衰决定骨骼的强弱。补肾中药及复方可以通过多种信号通路调节成骨细胞和破骨细胞的增生、分化，使之达到动态平衡，发挥治疗作用[3]。仙灵骨葆胶囊由淫羊藿、续断、补骨脂、知母、丹参和地黄组成，具有滋补肝肾、活血通络、强筋健骨等功效[4]。动物实验及临床研究均证实仙灵骨葆胶囊对OP、股骨头坏死、骨关节炎等疾病有显著治疗效果，已成为临床防治骨关节疾病的首选药物[5-9]。
组织蛋白酶K(cathepsin K，CatK)是半胱氨酸蛋白酶家族的成员之一，与OP的关系密切，是近些年备受关注的靶点蛋白[10-11]。CatK在人类破骨细胞中大量生成，主要参与骨胶原基质蛋白的降解，其过量表达可导致OP。通过对CatK结构功能的抑制可以有效减少其对骨有机质的降解[12-13]，而对骨形成无明显影响，因此其抑制剂的研究对治疗OP具有十分重要的意义。分子对接是基于配体与受体作用的锁匙原理，可以有效地确定与靶受体活性部位空间和电性相匹配的小分子化合物；而药效团模型则是对活性分子重要药效元素特征的整合，目前均广泛地应用于药物的设计和筛选研究[14-15]。本文以分子对接和药效团模拟技术为主要研究方法，对仙灵骨葆胶囊中所含的抗OP活性化合物进行了初步的筛选和讨论分析。
材料与方法 受体蛋白准备从RCSB蛋白质数据库http://www.rcsb.org中下载获得CatK不同结合位点蛋白-小分子抑制剂晶体复合物结构文件，其中活性位点与别构位点PDB ID分别为4x6i和5j94。采用AutoDock Vina 1.1.2分子模拟软件[8]对不同蛋白分子分别进行去水、加氢、加电荷和结构优化等前期处理。
小分子化合物库的构建通过文献调研和数据库搜索，对仙灵骨葆胶囊复方中单味药淫羊藿、续断、丹参、知母、补骨脂和地黄的化学成分进行整理分析(表 1)；利用OpenBabel v 2.3.1软件对所有小分子化合物进行三维结构生成、加氢和电荷处理操作，并以MMFF94分子力场进行优化以确保稳定性。去除重复的化合物，获得共5 509个小分子的化合物库，并保存为sdf文件。
| 复方组成 | 植物拉丁名 | 代表性化学成分种类 |
| 淫羊藿 | Epimedium，brevicornu | 黄酮苷、生物碱类 |
| 续断 | Dipsacus asperoides | 环烯醚萜苷类 |
| 丹参 | Salvia miltiorrhiza | 丹参酮、水溶性酚酸类 |
| 知母 | Anemarrhena aspho-deloides | 甾体皂苷、木质素和黄酮类 |
| 补骨脂 | Psoralea corylifolia | 香豆素类、黄酮和单萜酚类 |
| 地黄 | Rehmannia glutinosa | 环烯醚萜苷类 |
首先，利用AutoDock Vina分别将4x6i和5j94蛋白晶体结构中的配体3Y1、1XF重新对接(Re-Docking)到相应的结合位点，结合口袋分别由不同晶体复合物中配体小分子所在坐标为中心10 Å范围以内的所有氨基酸残基构成；同时在对接的过程中，小分子构象生成数目设置为20，Gird Box盒子大小为18.75×18.75×18.75 Å，其他参数使用默认值。根据对接后配体分子与原晶体结构中配体分子构象均方根偏差(root mean square deviation, RMSD)是否小于2.5 Å，判断对接参数设置的合理性及方法的适用性。其中，为了避免使用原晶体结构配体小分子构象重对接可能带来的“欺骗性”嫌疑，将其分离提取出后采用OpenBabel程序重新生成配体的3D结构作为对接的输入，这样生成的配体的构象便与晶体构象无关。化合物库的筛选采用相同的方法参数，整个对接过程使用Python脚本命令自动完成预处理靶点蛋白与所有配体小分子的对接计算模拟，并对不同构象对接后的打分值进行收集、整理和分析。同时，基于已经优化准备好的蛋白晶体复合物结构，利用Ligandscout 4.0软件生成不同结合位点的初始药效团模型(默认参数设置)；再根据受体复合物结合腔内关键氨基酸残基与小分子相互作用模式，选择其中比较重要的药效特征元素，对药效团模型进行优化。
虚拟筛选与结果聚类分析化合物库的筛选以分子对接模拟所得打分值(docking score)靠近前10%的化合物，进一步采用药效团模型进行复筛。然后剔除其中不符合Lipinski规则的小分子化合物，再根据与药效团模型的匹配打分(pharmacophore-fit score)，保留前30%的化合物。基于分子骨架指纹相似性的比较，并以群组间平均距离为度量标准，并对筛选结果进行聚类分析。
结果 分子对接方法准确性验证将CatK靶点不同结合位点蛋白晶体结构4x6i与5j94中配体3Y1、1XF提取出来，然后分别进行重对接处理；计算对接前后复合物体系中小分子构象的差异(图 1)，其RMSD值分别为1.32Å和0.987Å。
|
| 图 1 不同结合位点蛋白-小分子复合物重对接后最佳叠合构象及其打分 Figure 1 Optimal conformation and score for the re-docking of protein-ligand complex in different binding sites |
根据蛋白-小分子结合模式，对不同结合位点蛋白受体与小分子配体化合物的相互作用进行初步地分析。在活性结合位点，3Y1分别与4x6i蛋白活性位点Tyr67、Gly66以及Asn161发生疏水和氢键相互作用；而在别构位点，1XF与5j94蛋白之间则主要通过静电相互作用与Asn199，再经氢键作用分别与氨基酸Lys122、Arg123和Ala 1发生结合(图 2)。
|
| 图 2 不同结合位点蛋白-配体相互作用模式分析 Figure 2 Interaction analysis of protein and ligand for different binding-site modes |
根据晶体结构中配体小分子与蛋白结合位点关键氨基酸之间的相互作用，对不同位点蛋白-小分子结合模式中比较重要元素特征进行提取。活性结合口袋中3Y1母体结构芳环、羰基和氰基分别与4x6i蛋白分子中Tyr67、Gly66以及Asn161氨基酸之间发生疏水和氢键相互作用，可以作为主要的药效特征元素；而别构位点元素特征则主要有5j94中Asn199、Lys122和Ala 1氨基酸与配体小分子1XF的离子和氢键相互作用(图 3)。
|
| 图 3 基于蛋白-小分子复合物结构的药效团模型构建及其主要的化学元素特征 Figure 3 Building of pharmacophore modeling based on the protein-ligand complex structure and their main chemical features |
基于不同结合位点4x6i和5j94蛋白复合物晶体结构，采用分子对接和药效团模型组合方法对仙灵骨葆胶囊小分子化合物库进行筛选，最终分别得到98个和75个化合物。同时，采用基于分子骨架指纹相似性比较的方法，提取出其中比较有代表性骨架结构的化合物。聚类后样本距离矩阵的分析结果中对角线上颜色较深的方形区域分别代表一种结构骨架相似的化合物类群，例如：其中点A是预测结合在4x6i活性位点一个黄酮苷类化合物，点B则是筛选得到与5j94别构位点相结合的三萜类化合物(图 4)。
|
| 图 4 基于骨架指纹相似性比较的化合物聚类分析 Figure 4 Clustering analysis of the screened compounds according to the similarity of molecular skeleton fingerprints A:4X6i活性位点黄酮苷类化合物；B:5j94别构位点三萜类化合物 |
根据筛选化合物的聚类分析结果，结合Docking Score和Pharmacophore-Fit Score打分值进行综合评估，分别从活性位点和别构位点所筛选化合物中得到综合评分最佳的9个化合物(图 5, 6)。其中，针对不同位点结合蛋白最终预测筛选得到Top 9小分子化合物的对接与药效团模型匹配分值(表 2)。
|
| 图 5 基于4x6i蛋白结构结合位点筛选最靠前的9个活性化合物 Figure 5 Top 9 compounds that screening based on the active binding-site structure of 4x6i |
|
| 图 6 基于5j94蛋白别构位点筛选的9个活性化合物 Figure 6 Nine potential inhibitors for allosteric binding-site of 5j94 |
| 4x6iReaxys ID | DScore | PScore | 5j94Reaxys ID | DScore | PScore | |
| 30293362 | -8.4 | 46.12 | 15504189 | -6.6 | 58.18 | |
| 19950807 | -7.5 | 46.48 | 2955826 | -6.5 | 57.34 | |
| 28370395 | -7.3 | 46.82 | 20304885 | -6.3 | 58.32 | |
| 2227585 | -7.2 | 45.97 | 8130930 | -6.2 | 57.91 | |
| 4887090 | -7.0 | 46.35 | 2163967 | -6.1 | 57.84 | |
| 64198 | -7.0 | 46.66 | 2184909 | -6.1 | 57.79 | |
| 63067 | -7.0 | 46.07 | 8177094 | -6.1 | 58.42 | |
| 7100695 | -6.8 | 46.54 | 4210044 | -6.0 | 57.74 | |
| 4720754 | -6.5 | 47.87 | 269701 | -5.8 | 57.60 | |
| DScore:分子对接模拟打分值；PScore:药效团模型匹配打分值 | ||||||
中药复方化学成分复杂，呈现多成分-多靶点的综合相互作用与药效特征，从整体上对疾病网络进行调理，具有不良反应小及耐受性高的特点。分子对接和药效团模型分析作为常用的分子计算模拟技术，可以初步评价靶点蛋白与配体小分子相互作用的好坏，实现对中药活性成分的快速筛选，阐明方剂药效作用物质基础的关键。本文主要以OP疾病发生、发展和治疗中具有关键作用的CatK为靶点，采用分子对接和药效团模型分析的组合虚拟筛选技术与方法，对仙灵骨葆胶囊复方中药的药效物质进行研究。
研究结果发现，由靶点CatK不同结合位点蛋白晶体结构出发进行分子对接模型的构建，其对接前后复合物中小分子构象RMSD值均＜2.5 Å，表明该对接方法较好地重复出了原受-配体结合的构象，表明其参数设置合理、程序适用。同时，对图 2和图 3进行比较分析，可以看出分子对接过程中蛋白-小分子主要相互作用模式与药效团模型所提取出来的特征元素基本保持一致，说明采用Ligandscout软件构建基于晶体复合物结构的药效团模型能够很好地反映出蛋白-小分子相互作用特征，具有一定的合理性。
通过对仙灵骨葆胶囊复方中药小分子化合物库的筛选和骨架聚类分析，并结合分子对接和药效团模型匹配打分值进行综合评估，以4x6i和5j94蛋白晶体复合物中原抑制剂配体小分子重对接后的打分值进行比较(3Y1=-6.5，1XF=-4.5)，可以发现不同结合位点的Top 9化合物均显示出更好的对接亲和力；由此，可以初步推测其可能具有比各自蛋白晶体结构中配体小分子更为理想的抑制活性。
综上所述，将分子对接技术应用于探索药物治疗OP是可行的，同时也给中药复方靶点研究提供参考，为从分子水平研究中药作用机制提供了一定的依据及治疗OP作用机制研究提供一定的参考。
| [1] | Rachner TD, Khosla S, Hofbauer LC. Osteoporosis: now and the future[J]. Lancet, 2011, 377: 1276–1287. DOI:10.1016/S0140-6736(10)62349-5 |
| [2] | 夏泰. OP的危害及其防治[J]. 中华老年医学杂志, 2008, 27: 719–720. |
| [3] | An J, Yang H, Zhang Q, et al. Natural products for treatment of osteoporosis: the effects and mechanisms on promoting osteoblastmediated bone formation[J]. Life Sci, 2016, 147: 46–58. DOI:10.1016/j.lfs.2016.01.024 |
| [4] | 倪力刚, 王伟, 李春雯, 等. 仙灵骨葆胶囊对OP性髋部骨折的研究近况[J]. 中国骨质疏松杂志, 2011, 17: 1014–1018. DOI:10.3969/j.issn.1006-7108.2011.11.019 |
| [5] | 侯静怡, 祝乃强, 王学海, 等. 仙灵骨葆胶囊联合阿仑膦酸钠治疗原发性OP的临床应用效果探讨[J]. 双足与保健, 2018, 15: 165–166. |
| [6] | 李现林. 仙灵骨葆治疗中老年膝骨性关节炎的疗效观察[J]. 中国中医骨伤科杂志, 2007, 3: 41–42. |
| [7] | 唐卡毅. 仙灵骨葆防治骨质疏松症的研究进展[J]. 重庆医学, 2009, 38: 1535–1537. DOI:10.3969/j.issn.1671-8348.2009.12.063 |
| [8] | 杨海俊, 高祖玲, 陈莉, 等. 仙灵骨葆胶囊治疗肾性骨病肝肾亏虚证临床研究[J]. 河南中医, 2018, 38: 1442–1444. |
| [9] | 郑勇, 高建美, 龚其海. 仙灵骨葆的药理作用和临床应用研究进展[J]. 中国新药与临床杂志, 2018, 1: 1–5. |
| [10] | 何伟涛, 刘康, 孙金谞, 等. 组织蛋白酶K与OP治疗的研究进展[J]. 中国慢性病预防与控制, 2008, 14: 670–673. |
| [11] | 王敏, 张月琴. 以组织蛋白酶K为靶点的新药筛选[J]. 中国抗生素杂志, 2005, 30: 373–379. DOI:10.3969/j.issn.1001-8689.2005.06.017 |
| [12] | 李明, 吕厚辰, 张里程, 等. 组织蛋白酶K抑制剂Odanacatib治疗OP的研究进展[J]. 解放军医学院学报, 2015, 3: 285–287. DOI:10.3969/j.issn.2095-5227.2015.03.023 |
| [13] | 晁爱军. 组织蛋白酶K抑制剂Odanacatib研究进展[J]. 天津医药, 2014, 42: 1149–1152. DOI:10.3969/j.issn.0253-9896.2014.11.029 |
| [14] | 宋向岗, 周威, 陈超, 等. 基于分子对接方法的川芍治疗脑缺血的物质基础及分子机制研究[J]. 中国中药杂志, 2015, 40: 2195–2198. |
| [15] | 任洁, 魏静. 分子对接技术在中药研究中的应用[J]. 中国中医药信息杂志, 2014, 21: 123–125. |
| (收稿日期：2019-05-06) |