张天一,徐 苗,丁 欢,田菲菲
(西南交通大学生命科学与工程学院,中国四川成都610031)
宫颈癌是全球高发病率和高致死率的妇科恶性肿瘤,严重影响女性的身心健康和生活质量[1]。高危型人乳头瘤病毒(high-risk human papillomavirus,hrHPV)的持续感染是导致宫颈上皮细胞癌变的主要因素[2~4],其中HPV16和HPV18是引起宫颈癌的最常见的两种亚型,在宫颈癌和宫颈上皮瘤样病变中的检出率达到70%,分别引起宫颈腺癌和宫颈鳞癌[5~6]。此外,研究表明HPV感染也与肛门癌、阴茎癌和口腔癌等癌症的发生密切相关[7~10]。目前,市场上已有3种针对特定亚型HPV的预防性HPV疫苗[11~12],但这些疫苗较为昂贵,在发展中国家尚未普及[13~14]。而且,现阶段对已感染HPV的患者仍缺乏有效的治疗方法,每年有约27万人死于宫颈癌,五年的存活率仍维持在60%左右[7,15]。小分子药物因具有较好的空间分散性、良好的药物代谢动力学性质和非免疫原性,与传统细胞毒性药物不同,成为治疗宫颈癌的先导化合物的研究热点[16]。
高危型乳头瘤病毒的早期蛋白E6和E7可改变细胞内信号通路,破坏细胞周期调节的检查点,帮助病毒及其宿主细胞逃避免疫反应,因而在HPV的致癌特性中起关键作用[17~18]。hrHPV E6最显著的特性之一是与泛素连接蛋白E6AP、细胞周期调控蛋白p53形成三元复合物,并通过泛素化降解途径诱导p53的降解,阻止细胞凋亡,允许受损DNA的复制和变异的积累,这是hrHPV感染促使细胞癌变的重要因素[19~22]。通过模拟E6AP的LxxLL基序,一种被命名为pep11**的结合肽和小分子化合物木犀草素及其衍生物被发现能够特异性结合E6,竞争性抑制E6AP与E6的结合,从而增加HPV16阳性细胞的p53蛋白水平并引起细胞凋亡,表明E6的LxxLL结合区具有作为药物靶点的可行性[23~24]。根据近期获得的HPV16 E6/E6AP/p53三元复合物的晶体结构[25],我们发现E6AP与p53都和E6紧密结合,但E6AP与p53两者之间几乎没有直接作用,因此可以推测,抑制E6与p53的结合同样能够达到提高p53水平和恢复细胞凋亡能力的目的,但以E6的p53结合区域作为药物靶点的研究有待填补。
通过查阅文献,我们发现一些具有芳香性的化合物[26~29],如白藜芦醇、紫檀芪等多酚或黄酮化合物,能够诱导HPV阳性细胞在G2期积累,而细胞中野生型p53的存在对G2期的阻滞是必需的[30],提示这些化合物能够解除hrHPV E6蛋白对p53的抑制作用,有望开发为治疗宫颈癌的药物。Sangthong等[29]合成的几种蒽-9,10-二酮衍生物能够干扰HPV E6表达并诱导细胞凋亡,其中以苄胺为取代基的4-(苄基氨基)-9,10-二氧代-4a,9,9a,10-四氢蒽-1-基4-乙基苯磺酸盐的IC50值为0.3 μmol/L,显著高于其他取代基。Sangthong等认为这是由于芳香环的π-π作用相比静电作用可导致更高的凋亡能力,而支链结构可能降低了化合物进入靶标的能力,从而导致细胞凋亡能力丧失。遗憾的是,该团队并未更深入阐述其分子机制。我们的文献分析发现,富含芳香环、无侧链且具有宫颈癌细胞抑制作用的小分子化合物RITA(reactivation of p53 and induction of tumor cell apoptosis)和姜黄素(curcumin)符合上述特点,且两者分子大小相近、结构相似并对称,可能具有相似的作用方式。RITA是从美国国家癌症研究所化合物库中筛选出的,能够通过抑制HDM-2(human double minute-2)对p53蛋白的泛素化降解发挥抗癌作用的小分子化合物[31]。但是,有实验研究发现,RITA是通过抑制hrHPV E6蛋白介导的泛素-蛋白酶体途径对p53蛋白的降解,诱导p53蛋白的积累,从而恢复HPV16和HPV18阳性细胞的周期调控和凋亡功能,而非通过抑制HDM-2对p53的泛素化降解的方式产生抗癌作用[32~34]。因此,我们推测在HPV阳性细胞中RITA能够直接作用于E6蛋白,阻断E6蛋白与p53的结合,抑制复合物形成,避免p53蛋白被E6/E6AP介导的泛素-蛋白酶体途径降解,从而发挥肿瘤抑制功能。姜黄素是从姜科、天南星科植物的根茎中提取的一种二酮类化合物,具有降血脂、抗肿瘤、抗炎、利胆、抗氧化等作用。有文献报道,姜黄素可能是以HPV E6蛋白为靶点诱导细胞自噬并抑制HPV的复制[35~36]。虞艳霞等[37]的实验结果表明姜黄素能够对宫颈癌的放射治疗起到增敏作用。另外,相关研究表明姜黄素在与顺铂、紫杉醇等化疗药的联合使用中也呈现出对宫颈癌细胞的生长抑制作用[38~39]。这些结果显示姜黄素与RITA有着近似的作用效果。考虑到两者在分子结构方面具有一些相似的特征,我们猜想:RITA和姜黄素在抑制HPV阳性细胞p53蛋白的降解方面可能具有相同的分子机制。基于此,我们采用分子对接和分子动力学模拟的方法,分别对RITA、姜黄素与HPV16/18 E6蛋白p53界面的结合模式进行研究,以探究RITA和姜黄素抑制宫颈癌细胞的可能的分子机制,并为宫颈癌的药物设计提供理论指导。
使用量子化学波函数分析程序Multiwfn[40]计算并分析RITA与姜黄素分子的表面静电势分布,同时分析配体在E6蛋白复合物中的相对朝向以及次级键的成键倾向。使用可视化软件VMD[41]对计算结果制图。
在 NCBI数据库(https://www.ncbi.nlm.nih.gov/)获取HPV18 E6蛋白的氨基酸序列(GI:1258167973),NCBI BLAST显示其与HPV16 E6蛋白氨基酸序列的一致性为59.64%。随后,分别使用Modeller 9.20 软件、在线服务器 SWISS-MODEL(https://www.swissmodel.expasy.org/)和在线服务器I-TASSER(https://zhanglab.ccmb.med.umich.edu/I-TASSER/)[42~44]进行同源建模,并使用SAVES(http://servicesn.mbi.ucla.edu/SAVES/)对模型质量进行评价。
在Chemical Book(www.chemicalbook.com)化合物数据库网站下载RITA(CAS#:213261-59-7)和姜黄素(CAS#:458-37-7)的分子文件作为配体。从PDB数据库(www.rcsb.org)下载HPV16 E6/E6AP/p53三元复合物的晶体结构文件(PDB ID:4XR8),取其中的E6蛋白结构(H链)作为对接的受体分子。使用分子对接软件AutoDock 4.2.6[45]进行对接操作。对接采用遗传算法,其余设置选择默认参数。根据E6蛋白晶体结构和文献[25]可知,HPV16 E6蛋白与p53的作用区域大致可分为3个亚界面(图1),我们分别以这3个亚界面作为对接区域。HPV18 E6蛋白对应HPV16 E6三个亚界面的残基选择对接区域。
选择HPV16 E6与RITA和姜黄素对接结果中符合条件的构象,使用the PRODRG Server(http://davapc1.bioch.dundee.ac.uk/cgi-bin/prodrg/submit.html)构建蛋白质配体复合物文件,用GROMACS 5.1分子模拟套件进行分子动力学模拟[46]。采用适用于蛋白质分子动力学模拟的Gromacs54a7力场,添加SPC水模型,构建蛋白质和配体复合物的原子文件与拓扑文件,盒子形状设置为六面体,盒子最近处距离蛋白质边缘的距离设置为1.2 nm,加入0.15 mmol/L的氯离子和钠离子以平衡电荷。使用最陡梯度法对该体系进行能量最小化。对能量最小化后的体系依次进行200 ps的热浴(NVT,number of particles,volume,and temperature)和1 000 ps的压浴(NPT,number of particles,pressure,and temperature)系综模拟,之后进行30 ns的位置非限制性分子动力学模拟。模拟过程的步长时间设置为2 fs,温度设置为310 K以模拟人体温度。
使用PLIP在线服务器(http://plip.biotec.tudresden.de/plip-web/plip/index)[47]对分子对接和分子动力学模拟的结果进行配体和受体结合模式的分析,使用UCSF Chimera程序[48]显示蛋白质配体复合物构象。
图1 HPV16 E6蛋白上与p53蛋白作用的3个亚界面Fig.1 Three sub-interfaces of HPV16 E6 protein interacting with p53 protein
分子表面静电势的计算结果经VMD软件可视化处理后如图2所示:RITA和姜黄素表面静电势的分布大体相似,分子长轴两端的静电势绝对值较大,芳香环结构表面的静电势绝对值较小;不同的是,具有3个芳香环的RITA仅在两端具有较大的静电势分布,而姜黄素分子的中部由于存在两个酮羰基,所以也具有较大的静电势。相似的静电势分布意味着配体在复合物中具有近似的相对朝向及氢键位点,即比较一致的构象和结合模式,显示了RITA和姜黄素能够产生相同药理作用的结构基础。
Modeller 9.20、SWISS-MODEL 和 I-TASSER构建的HPV18 E6蛋白模型以SAVES在线服务器进行模型质量评价。Verify 3D评价结果显示,SWISS-MODEL建模结果中82.69%的氨基酸残基的3D-1D分数高于0.2,建模质量较好,而另外两者未通过评价(图3A)。PROCHECK评价结果显示,Modeller建模结果最好(图3A),但拉氏图显示SWISS-MODEL建模结果只有ARG152在不允许区(图3B),远离选择作为对接区域的3个亚界面,对于本文的研究结果影响甚微。ERRAT评价方面,虽然SWISS-MODEL的分数83.783 8低于I-TASSER的92.666 7,但通过查看每个残基分数发现不合理残基主要位于远离结合对接区域的C末端(图3C),说明SWISS-MODEL模型中3个亚界面的残基评价比较合理。综合考虑各项评价,我们认为SWISS-MODEL建模结果更适用于本文研究。
2.3.1 HPV16 E6蛋白的分子对接结果及分析
通过查看HPV16 E6蛋白上与p53蛋白作用的3个亚界面的对接结果中能量最低的构象,对RITA和姜黄素与HPV16 E6蛋白的结合模式进行比较和分析,结果如图4所示。
亚界面Ⅰ为HPV16 E6蛋白的6~18位残基区域。对接结果显示,在此区域RITA更倾向于与E6蛋白的ARG8、ARG10、GLN14残基形成氢键,而姜黄素则更多通过疏水作用与E6蛋白的GLN14、THR17、GLU18 残基结合。
亚界面Ⅱ为HPV16 E6的41~49位残基区域。对接结果中的最低能量构象显示,RITA与E6蛋白的ILE23、ALA46、PHE47残基形成疏水键,与THR17、ASP49残基形成氢键;姜黄素与HPV16 E6蛋白的 LYS11、PRO13、PHE47残基形成疏水键,与LYS11、THR17、ASP49形成氢键。两个小分子的两端都与E6蛋白THR17、ASP49残基形成氢键,形成了比较一致的配体结合构象。
亚界面Ⅲ为HPV16 E6的97~115位残基区域,结果显示RITA与PRO109、PRO112形成疏水键,与 SER97、ILE101、LYS115 形成氢键;姜黄素未形成疏水键,而是通过甲氧基和羟基与LEU99、ILE101、LYS108、PRO109、LYS115 形成较多氢键。
2.3.2 HPV18 E6蛋白的分子对接结果及分析
通过查看HPV18 E6蛋白上与p53蛋白作用的3个亚界面的对接结果中能量最低的构象,对RITA和姜黄素与HPV18 E6蛋白的结合模式进行比较和分析,结果如图5所示。
亚界面Ⅰ为HPV18 E6蛋白的6~18位残基区域。对接结果显示,在此位点RITA和姜黄素都与E6蛋白10~20位残基形成丰富的氢键,且都与TYR12形成疏水键,但RITA也能与ASP16、THR19和GLU20形成疏水键。
图2 RITA和姜黄素分子表面的静电势分布Fig.2 Electrostatic potential distribution on the surfaces of RITA and curcumin
图3 HPV18 E6同源建模结果的SAVES评价(A)三种方法建模结果的评价总览;(B)SWISS-MODEL建模结果的主拉氏图。用于评价氨基酸残基的空间合理性,允许区为红色,最大允许区为黄色,不允许区为白色;(C)SWISS-MODEL建模结果的ERRAT残基评价。通过重原子间非键作用评价氨基酸残基分布,绿色为正确,黄色为警告,红色为错误。Fig.3 SAVES evaluation of homology modeling for HPV18 E6(A)Overview of evaluation of three homology modeling methods;(B)Main Laplace diagram of SWISS-MODEL modeling results.To evaluate the spatial rationality of amino acid residues,the allowable area is red,the maximum allowable area is yellow,and the unallowable area is white;(C)ERRAT residue evaluation of SWISS-MODEL modeling results.To evaluate the distribution of amino acid residues by non-bonding between heavy atoms,green is correct,yellow is warning,and red is wrong.
图4 HPV16 E6蛋白对接构象及结合模式左图:E6蛋白亲疏水表面的配体叠加构象,其中RITA显示为黄色,姜黄素显示为绿色;右图:RITA、姜黄素分别与E6蛋白的相互作用。Fig.4 Protein docking conformation for HPV16 E6 and binding modeLeft panel:RITA is shown in yellow in the ligand-superimposed conformation of the protein-hydrophobic surface,and curcumin is shown in green;Right panel:Interaction between RITA/curcumin and E6 protein.
亚界面Ⅱ为HPV18 E6的41~49位残基区域。对接结果显示,RITA的噻吩结构与E6蛋白的 PHE45、PHE49 形成 π-π 堆积,与 GLN26、GLU46形成氢键,且与PHE49形成疏水键;姜黄素一端的芳香环结构与HPV18 E6蛋白的PHE45残基形成π-π堆积(这与RITA的结合情况相似),而另一端的芳环结构仅与E6蛋白的LYS110、ASN113形成氢键,与LYS110形成疏水键。我们分析认为这可能是由于姜黄素的分子结构中部缺乏芳香环,所以与PHE49形成疏水键而非π-π相互作用。
亚界面Ⅲ为HPV18 E6的97~115位残基区域,结果显示 RITA 与 PHE49、LEU102、PRO111形成疏水键,与 LEU101、ILE103、LYS110、ASN113、LYS117形成氢键;姜黄素与TYR99、PRO111形成疏水键,与 PHE49、LYS110、ASN113、LYS117形成氢键,且与TYR99形成π-π堆积。
根据Martinez-Zapien等[25]解析的HPV16 E6/E6AP/p53晶体结构,HPV16 E6蛋白的PHE47、ILE23、HIS24和TYR43残基为其与p53结合贡献疏水作用,ASP49与p53的HIS115以氢键形式相互作用,F47R或D49A突变会干扰三元复合物的形成并抑制p53的降解。这同我们的RITA和姜黄素与HPV16 E6蛋白亚界面Ⅱ的对接结果的成键情况非常吻合,而且不论是HPV16还是HPV18的E6蛋白,RITA和姜黄素在亚界面Ⅱ均有更好的构象重叠,与亚界面Ⅱ残基的结合方式也更一致,因此我们选择亚界面Ⅱ对接结果中的蛋白质配体复合物作为初始构象进行分子动力学模拟研究。
图5 HPV18 E6蛋白对接构象及结合模式左图:E6蛋白亲疏水表面的配体叠加构象,其中RITA显示为黄色,姜黄素显示为绿色;右图:RITA、姜黄素分别与E6蛋白的相互作用。Fig.5 Docking conformation for HPV18 E6 and binding modeLeft panel:RITA is shown in yellow in the ligand-superimposed conformation of the protein-hydrophobic surface,and curcumin is shown in green;Right panel:Interaction between RITA/curcumin and E6 protein.
在30 ns分子动力学模拟过程中,无论是RITA还是姜黄素都没有从HPV16/18 E6蛋白的41~49位的疏水沟中脱离,表明该区域的结合稳定。于分子动力学模拟25 ns之后选择复合物体系能量较低时刻的构象进行分析,与对接结果相比,RITA和姜黄素与HPV16/18 E6蛋白结合模式中,疏水作用的比重显著增大(图6),推测是E6蛋白表面41~49位残基形成的疏水沟与溶液环境共同作用的结果。从30 ns的动态模拟中多次截取瞬时构象进行观察,发现HPV16 E6蛋白的PHE47残基的侧链苯环与配体的芳香环以疏水键或π-π堆积的方式发生相互作用,而HPV18 E6的PHE47残基侧链包埋于蛋白质内,但PHE45和PHE49残基的侧链苯环暴露于蛋白质表面,代偿PHE47的侧链苯环与配体发生疏水作用和π-π相互作用。在HPV18 E6蛋白复合物的30 ns分子动力学模拟中,RITA与E6蛋白的相对构象变动较小,而姜黄素沿疏水沟方向发生较大的位移,这与HPV16 E6蛋白复合物的模拟结果相差较大。对比分析HPV16与HPV18的E6蛋白41~49位残基的构成及配体结合成键情况,推测出现这种差异是由于HPV18 E6亚界面Ⅱ表面存在PHE45和PHE49侧链的两个苯环,芳香环结构更为密集的RITA恰好能与两个苯环同时形成疏水键和π-π相互作用,而姜黄素的芳香环结构位于分子两端,当一端与PHE45(或PHE49)侧链形成疏水键或π-π相互作用时,由于空间分布和构象能量的原因,另一端的芳香环与PHE49(或PHE45)侧链相距较远,而且由Multiwfn计算结果可知,姜黄素分子中部有较高的表面静电势分布,在41~49疏水沟中难以发生相互作用并稳定构象;HPV16 E6在此区域只有PHE47侧链苯环位于蛋白质表面,RITA分子中更为紧凑的3个芳环结构使其与PHE47的相互作用具有更高的自由度,并易于在多种构象间相互转化。这也能够解释在Zhao等[34]的实验中为什么RITA对HeLa细胞(HPV18)的抑制效果高于CaSki细胞(HPV16)。
图6 分子动力学模拟25 ns之后复合物较低能量的构象及结合模式(A)HPV16 E6-RITA复合物;(B)HPV16 E6-姜黄素复合物;(C)HPV18 E6-RITA复合物;(D)HPV18 E6-姜黄素复合物。Fig.6 The low-energy binding patterns and interactions of compounds by molecular dynamics simulation after 25 ns(A)HPV16 E6-RITA complex;(B)HPV16 E6-curcumin complex;(C)HPV18 E6-RITA complex;(D)HPV18 E6-curcumin complex.
基于已测得的晶体结构和相关文献[25]报道,我们推测小分子抑制剂RITA和姜黄素可能通过作用于hrHPV E6蛋白的p53结合界面,抑制E6蛋白对抑癌蛋白p53的降解。本研究通过同源建模构建HPV18 E6蛋白模型,使用AutoDock软件将RITA和姜黄素与HPV16/18 E6蛋白的3个亚界面分别进行分子对接,比较两个小分子在这3个区域的对接模式,结果显示RITA和姜黄素与HPV16 E6蛋白或HPV18 E6蛋白均在亚界面Ⅱ有着相似的结合模式和较为一致的构象重叠。为更进一步研究体液环境中RITA和姜黄素与亚界面Ⅱ的结合模式,使用GROMACS 5.1进行30 ns的分子动力学模拟,结果表明RITA和姜黄素均能与HPV16/18 E6蛋白的41~49位残基形成较为稳固的疏水作用,E6蛋白亚界面Ⅱ的PHE残基的侧链苯环是提供疏水作用和π-π相互作用的重要结构。结合文献报道的HPV16 E6 F47R/G定点突变导致其诱导p53降解能力降低的实验结果[49~51],可以设想,抑制剂通过竞争性结合HPV16 E6蛋白PHE47可能获得同样结果,同理HPV18 E6蛋白的PHE45和PHE49可作为抑制剂结合靶点。这一结果或可作为宫颈癌特异性治疗药物分子设计的基础,而且受体的亚界面Ⅱ不仅是E6蛋白与p53蛋白相互作用的区域,本身也是暴露在E6蛋白表面最大的疏水区域,具有作为药物靶点的可能性。此外,从配体方面来看,小分子抑制剂的芳香性基团对结合贡献了较多的疏水作用和π-π相互作用,可作为设计与筛选更多小分子抑制剂的结构特征。对此我们也做了虚拟筛选方面的尝试,且筛选出的化合物比较符合这一特征规律,结合能低的主要是具有3~4个芳香环的多酚化合物,与RITA和姜黄素结构比较相似。HPV的E6蛋白是具有多种生物活性的低相对分子质量蛋白质,其结构极易变化,在我们的分子动力学模拟过程中结合了不同配体的E6蛋白的构象发生了显著变化,这种蛋白质构象的变化可能会破坏三元复合物的形成并使p53免于降解。
虽然RITA和姜黄素在体外和体内实验都显示了对宫颈癌的治疗作用,但由于耐药性和生物利用度低等特点尚未能够应用于临床[52-53]。很多研究证明,RITA、姜黄素等小分子对癌细胞的抑制除了通过恢复p53途径外还有多种其他作用途径,如 IRE1α/XBP1、p38 和 JNK/SAPK 途径[54~56],而且对于RITA、姜黄素激活p53途径的分子机制和结合模式也存在争议,仍需更多研究。本研究进行分子对接和分子动力学模拟,旨在通过比较两者与HPV16/18 E6蛋白的结合模式,寻找潜在的药物靶点,为宫颈癌的研究和新药设计提供更多参考。在后续研究中,我们将以本研究的结论作为理论基础,以靶点结构和小分子抑制剂特征为出发点,进行更多探索,以期发现更有治疗价值的活性小分子。