癌症是威胁生命的主要疾病之一,其发病率在全球范围内逐年上升
[1]。近年来,肿瘤治疗相关研究热度不断上升,急需探索新的肿瘤治疗模式
[2-4]。化疗作为目前临床应用最广泛的治疗手段,仍存在明显的局限性
[5-6],包括由非特异性受体亲和性引发的全身毒性、药物溶解性差、生物利用度欠佳、易产生耐药性以及治疗周期长等
[7-9]。这些缺陷也推动了新疗法开发,特别是针对肿瘤血管系统的抗血管生成疗法(饥饿细胞疗法)
[10-11]。饥饿细胞疗法通过抑制或阻断内皮依赖性血管生成,切断肿瘤的营养来源和能量供应,最终使肿瘤细胞因缺氧而死亡
[12-14]。这种疗法的靶点通常是细胞代谢过程中的关键分子、信号通路或酶
[15-16]。
血管内皮生长因子受体2(vascular endothelial growth factor receptor 2,VEGFR2)是一种在内皮细胞中表达的酪氨酸激酶受体
[17]。VEGF,VEGFR2信号通路,通过在肿瘤血管生成过程中提供氧气和营养物质来刺激肿瘤生长。临床实验结果显示,在多种实体瘤中抑制VEGFR2,可诱导肿瘤缺氧和坏死
[18]。因此,VEGFR2成为癌症治疗的一个有吸引力的靶点
[19-20]。目前,美国食品药品监督管理局(Food and Drug Administration,FDA)已经批准抗VEGFR2药物用于治疗多种肿瘤
[21-22],VEGFR2药物包括舒尼替尼
[23]、索拉非尼
[24]、阿帕替尼
[25]、伦伐替尼
[26]及替沃扎尼
[27]等。同时,FDA还对很多候选药物(如卢西他尼
[28]、伐他拉尼
[29]等),进行临床试验。虽然VEGFR2药物的相关研究取得一些进展,但目前的VEGFR2抑制剂可引起高血压、蛋白尿、甲状腺功能减退、白质脑病综合征和动脉血栓形成等
[30],严重制约其长期临床使用
[31-32]。因此,开发全新骨架结构VEGFR2抑制剂具有重要意义。计算机辅助药物设计(computer aided drug design, CADD)通过整合多学科技术,实现靶向聚焦的理性筛选
[33],革新了药物发现领域
[34]。整合分子对接、药效团建模、分子动力学(molecular dynamics, MD)模拟和结合自由能计算等策略,在药物发现领域发挥了重要作用
[35]。例如,通过同源建模方法构建VEGFR2蛋白、使用药效团筛选技术研发的用于治疗软组织肉瘤的帕唑帕尼
[36];基于药效团筛选技术发现的用于治疗乳腺癌、前列腺癌的三嗪-吲哚衍生物
[35];通过分子对接与药效团筛选技术发现的用于治疗晚期肾细胞癌的阿昔替尼
[37];通过分子对接技术发现的用于治疗卵巢癌的奥拉帕利
[38];通过已知抑制剂进行药效团建模、使用分子对接技术发现的用于治疗类风湿性关节炎的乌帕替尼
[39];通过分子对接、分子动力学和AI筛选技术研发的治疗人类免疫缺陷病毒(HIV)感染的洛匹那韦
[40]。这些研究结果显示,计算机辅助药物设计在药物发现中展现出巨大潜力。文中通过结构虚拟筛选、分子动力学模拟和结合自由能计算的多层次策略,并结合体外验证实验,发现了具有新型分子骨架的VEGFR2抑制剂。
1 实验方法与测试
1.1 材料与软件
小分子化合物(上海陶术生物科技有限公司);人血管内皮生长因子受体2(VEGFR2)酶联免疫分析(ELISA)试剂盒(上海科艾博生物技术有限公司)。
使用GNINA 1.3
[41]软件进行分子对接分析,使用GROMACS 2024.4
[42]软件进行分子动力学模拟,使用gmx_MMPBSA
[43]工具计算MM-PBSA结合自由能
[44]。
1.2 分子对接虚拟筛选
VEGFR2蛋白晶体结构(PDB ID: 3WZE)从蛋白质结构数据库(
https://www.rcsb.org)下载,结构分辨率为1.90 Å。文中筛选的配体小分子来自ZINC、 PubChem、 ChemDiv数据库(共计222万个分子)。
GNINA 1.3软件关键参数的设置:以索拉非尼分子与VEGFR2结合时的结合位点为参考定义对接盒子范围,搜索穷尽度(exhaustiveness)为16,采用CNN打分模式,以CNN Pose Score大于0.8作为初筛阈值。人工目视筛选标准化原则:配体与目标蛋白之间形成氢键、卤键、盐桥及π-π堆积等非共价键相互作用,并且二者形成的构象无空间冲突。
1.3 分子动力学模拟
通过GROMACS 2024.4软件对经分子对接虚拟筛选方法得到的化合物进行分子动力学模拟。使用Amberff19SB
[45]、GAFF2
[46]力场分别建立蛋白质和配体相互作用的拓扑参数。每个蛋白配体复合物体系都溶解在一个立方体水盒中,水分子使用TIP3P模型
[47],通过添加Cl
-维持每个体系保持电中性。首先,使用最速下降算法对复合物体系进行10 000步的能量最小化模拟;随后,分别进行500 ps的NVT、NPT平衡模拟;最后,进行100 ns的无约束分子动力学模拟。对分子动力学模拟轨迹进行均方根偏差(
d)、均方根波动(
srm)、回转半径(
rg)和溶剂可及表面积(
s)分析,进而评估复合物体系是否达到平衡。
1.4 结合自由能的计算
MM-PBSA是一种速度和精度兼备的计算结合自由能(Gb)的方法,其基本原理是计算两个溶剂化分子在结合和游离状态的Gb差异。从复合物体系平衡阶段的最后10 ns轨迹中提取1 000帧,间隔时间为10 ps,计算复合物体系的Gb。对于一个蛋白质-配体复合物体系,Gb越低,复合物的结合亲和力越强。
1.5 药代动力学预测
通过在线预测网站SwissADME (
http://www.swissadme.ch/index.php),对筛选出的化合物分子进行药代动力学预测,包括相对分子质量、脂水分配系数(lg
P)、拓扑极性表面积(
st)、水溶性(log
S)、分子柔性(可旋转键数,
Nb)、氢键受体数目(
Na)和氢键供体数目(
Nd)。
1.6 抗肿瘤活性测试
VEGFR2在胃癌细胞中常常呈高度表达状态,促进肿瘤中血管生成与肿瘤进展
[48]。选用人胃癌细胞株MGC-803(武汉普诺赛生命科技有限公司)进行体外实验。通过MTT比色法评价筛选出化合物的抗肿瘤活性,再通过细胞实验验证虚拟筛选结果。
将MGC-803置于含有10%胎牛血清、w=1.0%青霉素-链霉素双抗溶液DMEM(Dulbecco’s modified eagle medium)的完全培养基,于37 ℃、φ(CO2)=5%的培养箱中培养。将处于对数生长期的MGC-803接种到96孔板(密度为2 000个/孔)培养过夜。然后加入溶解在100 μL培养基中不同浓度化合物溶液,并置于培养箱中孵育72 h。再将10 μL ρ(MTT)=5 μg/mL的磷酸盐缓冲液(PBS)加入每个孔,并置于培养箱中孵育4 h。最后,用BioTek微孔板检测仪(美国Agilent Technologies公司)在570 nm处测定各孔溶液的光密度,并计算ρ(IC50)。每组实验平行重复3次。
1.7 酶活性测试
收集对数生长期的MGC-803,以1×106 个/mL密度接种至6孔板,置于37 ℃、φ(CO2)=5%的培养箱。细胞贴壁后,更换新鲜培养基,再向每孔分别加入20 μmol/L化合物溶液处理24 h。然后除去含化合物溶液的培养基,用pH=7.4 PBS溶液清洗细胞2~3次,用w=0.25%胰蛋白酶进行消化。收集细胞,离心分离,用PBS溶液稀释细胞悬液,计数,采用超声细胞破碎仪裂解细胞。在4 ℃、3 000 r/min条件下离心分离20 min,收集上层清液,置于-80 ℃冻存备用。
取出VEGFR2酶联免疫分析试剂盒,室温下平衡20 min。同时,样品在室温下平衡60 min。向标准品孔中分别加入50 μL ρ=0,0.625,1.25,2.5,5.0,10.0,20.0 ng/mL酶标准品。向样品孔中分别加入10,40 μL待测样品(化合物溶液)稀释液,再加入100 μL酶标试剂。使用封板膜封板,于37 ℃温育60 min。揭掉封板膜,弃去液体,快速甩干。每孔加满洗涤液,静置30 s,弃去洗涤液,重复上述操作5次。向上述孔中先加入50 μL显色剂A,再加入50 μL显色剂B,轻轻震荡混匀,置于恒温培养箱中于37 ℃避光孵育15 min。最后,向每孔加入50 μL终止液,终止反应。以空白孔调零,测定λ=450 nm处各孔溶液的光密度(D)。每组实验平行重复3次。根据标准曲线计算待测样品的浓度。
2 结果与分析
2.1 分子对接筛选结果
通过GNINA 1.3软件进行分子对接筛选。将VEGFR2蛋白中索拉非尼分子结合的位点定义为活性口袋,该口袋主要由保守的
β折叠构成,重要氨基酸残基为GLU71、CYS105、ASP182、PHE183。以索拉非尼分子作为对接筛选的阳性对照,以数据库中的222万个化合物为分子对接筛选的配体。通过分子对接打分(CNN pose score)与人工目视筛选相结合方法,共筛选出22个候选化合物(
表1)。
分析22个候选化合物与VEGFR2蛋白的结合模式。阳性对照组(索拉非尼分子)与VEGFR2蛋白中CYS105、ASP182形成氢键,与GLU71形成两个氢键,与PHE183形成π-π堆积。化合物1与VEGFR2蛋白中GLU71、ASP182形成氢键,与PHE183形成π-π堆积。其他化合物(2~22)均与GLU71、CYS105、ASP182形成氢键,其中,化合物5、化合物18与GLU71残基形成两个氢键。
2.2 MD模拟与MM-PBSA结合自由能计算结果
首先,使用GROMACS 2024.4软件对由对接筛选技术得到的22个化合物进行100 ns的MD模拟。然后,使用gmx_MMPBSA工具进行MM-PBSA结合自由能计算。最终,筛选出8个复合物(化合物1~化合物8)(
图1)。MD模拟结果显示,8个化合物与VEGFR2蛋白结构域之间形成的相互作用与索拉非尼相似,如VEGFR2蛋白中PHE183与配体形成π-π堆积,GLU71、CYS105、ASP182与配体形成氢键,这些作用有助于复合物稳定。同时,在整个100 ns模拟期间,这些相互作用保持不变,进一步证实,这些氨基酸残基对于维持VEGFR2位点活性和稳定小分子化合物有重要作用。
8个化合物与VEGFR2形成的8个复合物的均方根偏差(
d)均稳定在0.2 nm左右,表明复合物达到了平衡状态。其中,化合物1与VEGFR2蛋白结合形成复合物的分子动力学模拟结果见
图2。均方根被动(
srm)可表示分子中各个原子的自由程度。对8组氨基酸残基的轨迹进行
srm分析,发现只有位于蛋白质外边缘的氨基酸残基表现出柔性,而结合口袋中的氨基酸残基相对稳定。对活性口袋中关键残基的
srm进行分析,结果显示,GLU71、ASP182及ASP183的
srm分别为0.05,0.04,0.04 nm,均较低,显示VEGFR2与小分子配体结合,使这些关键氨基酸残基失去柔性。回转半径(
rg)是衡量分子紧凑性的重要参数,可评估蛋白质与配体结合后的构象及其结合紧密性。8个复合物的
rg值平均稳定在2 nm左右,表明在100 ns的模拟过程中,这些复合物均保持着高度紧密构象。表面积(
s)反映了生物分子表面可被溶剂接触的程度,其变化可评估蛋白质的折叠状态。8个复合物的
s显示,在100 ns内蛋白趋于稳定。
MM-PBSA是计算小分子抑制剂与靶蛋白之间
Gb的一种有效、可靠的方法。一般,
Gb越低,配合物越稳定,配体的活性和效力也越高。文中选取8个复合物结合过程最后10 ns的动力学轨迹进行MM-PBSA计算(
表2),
表2中,
Gb分项包括范德华作用能(
Ew)、静电能(
Ee)、溶剂化极化能(
Es)、非极性溶剂化能(
En)、气相自由能(
Gg)和溶剂化自由能(
Gs)。由
表2可知,
Ew,
Ee,
En,
Gg均有助于蛋白与配体的结合,其中,
Ew,
Gg的作用尤为显著;所有复合物总的
Gb均为负值,表明这些复合物在热力学上是稳定的,且其结合过程可以自发进行。
2.3 药代动力学预测结果
通过SwissADME工具对索拉非尼和筛选得到的8个候选化合物进行药代动力学性质预测(
表3)。所有化合物的相对分子质量(
M)均低于500 g/mol,符合Lipinski规则中对口服药物相对分子质量的要求。脂水分配系数lg
P=0~5,表明8个候选化合物均具有适宜的亲脂性,且亲脂性可促进细胞膜通过性。拓扑极性表面积(
st)为60~130 Å
2,有利于胃肠道吸收,但候选化合物6的
st稍高(131.15 Å
2),可能影响人体对它的吸收。水溶性(log
S)预测结果显示,8个候选化合物的水溶性均较低(log
s<-4)。其中,化合物8的水溶性较差,可能限制它的生物利用度。8个候选化合物的可旋转键数(
Nb)为5~9,表明8个候选化合物分子具有适度的柔性。
2.4 体外生物活性实验
抗肿瘤活性实验结果显示,化合物1~化合物5均可抑制MGC-803细胞的活性(
表4),显示这5个化合物能进入细胞并与靶蛋白结合,进而发挥抑制作用。其中,化合物1的抑制性最强,其
c(IC
50)=(5.54±0.25) μmol/L,而化合物7、化合物8的抑制性较差。
体外酶活性实验结果显示,与对照药物索拉非尼(抑制率(
ri)为(70.35±9.43)%)相比,化合物5、化合物6、化合物8均能较好地抑制VEGFR2酶活性(
表5)。对
ri进行单因素方差分析,发现
ri的差异具有统计学意义(
p<0.01)。化合物5、化合物6和化合物8的
ri分别为(79.32±7.86)%,(72.06±14.29)%,(76.86±8.87)%,显示这3个化合物能够较好地与VEGFR2的活性位点结合,进而发挥抑制作用。尽管化合物5(
Gb=-33.40 kJ/mol)比索拉非尼(
Gb=-77.00 kJ/mol)有更高的
Gb,但由于化合物5的结构应变能比索拉非尼更低,导致化合物5的
ri (79.32%)高于索拉非尼(70.35%)。
3 结论
采用分子对接、分子动力学模拟以及结合自由能计算等技术与方法,筛选出潜在的VEGFR2抑制剂。首先通过GNINA 1.3软件对VEGFR2靶点进行基于结构的虚拟筛选,从222万个小分子中筛选出22个候选化合物。这些候选化合物与VEGFR2蛋白的结合模式和相互作用方式与参考药物索拉非尼的相似。进一步,分子动力学模拟筛选出8个与VEGFR2紧密结合的化合物。结合形成的8个复合物在100 ns内达到平衡,均方根偏差稳定在0.2 nm左右;关键结合氨基酸残基的均方根波动均在0.05 nm左右,证实筛选出的8个化合物与VEGFR2可稳定结合。此外,回转半径和溶剂可及表面积计算结果也显示,复合物构象结合紧密。MM-PBSA结合自由能显示,化合物2的结合自由能与索拉非尼的相当。由于熵效应的计算量大,在计算过程中未考虑阳性对照(索拉非尼)与筛选出化合物的熵效应,但是通过分析候选化合物相对于索拉非尼的相对结合能,可帮助预测反应活性。
体外活性实验结果显示,筛选出的8个化合物与索拉非尼分子有近似的半抑制浓度c(IC50),其中,化合物1显示出最好的VEGFR2酶活性抑制作用,其c(IC50)=(5.54±0.25) μmol/L。化合物5、化合物6及化合物8对VEGFR2酶活性抑制率高于索拉非尼,其中,化合物5的抑制率最高((79.32±7.86)%)。该研究结果为VEGFR2靶向抑制剂的开发提供了新分子骨架,并为后续的先导化合物优化提供了重要参考。