研究不同立场对B-DNA到A-DNA构型转化过程产生的结果,生物化学论文

研究不同立场对B-DNA到A-DNA构型转化过程产生的结果,生物化学论文Abstract:TheaimofthepresentworkistoinvestigateandcomparetheeffectofthelatestCHARMMandAMBERforcefields〔includingbsc1andOL15〕ontheB-DNAtoA-DNAconversionthroughexploringthefree-energychangesoftheconversionprocess.Theextendedadaptivebiasingforce〔eABF〕methodwasutilizedtoperformthefree-energycalculations.Theresultsshowedthatthefree-energyprofilescharacterizingthetransitiondiffersignificantlyforthesetwoforcefields.TheAMBERforcefieldperformsbetterthantheCHARMMforcefieldinaqueoussolution.Thestructureneartheglobalminimumofthefree-energyprofilebytheAMBERforcefieldpresentsB-form,inagreementwiththeexperimentalresults,whilethemoststablestructurebytheCHARMMforcefieldlocatesbetweenA-andB-form.DeepanalysisoftheradialdistributionfunctionsofthecounterionsaroundDNArevealsthatthedistributionofcounterionsinminorgrooveusingtheCHARMMforcefieldislowerthanthatusingtheAMBERforcefield.Therefore,fortheCHARMMforcefield,therepulsionofphosphatesbackbonecouldnotbeproperlycounteractedbycounterions,asaresult,theminorgroovebecomeswider,causingaslightconformationalchangetowardsA-form.Keyword:B-DNA;A-DNA;Conformationalconversion;CHARMMforcefield;AMBERforcefield;Free-energycalculation;DNA作为生命信息传递的载体，主要存在B-DNA,A-DNA和Z-DNA3种形式[1].通常，在水活度较高的环境〔如细胞〕中，DNA主要以B构型的形式存在。A-DNA则主要存在于水活度较低的环境〔如高盐或高乙醇溶液〕中，当DNA与RNA配对时，同样会呈现A构型[2].B构型与A构型DNA的构造如此图1所示，两者均呈现右旋状态，但A构型小沟较宽，螺旋较短，呈现更严密的螺旋状态。DNA构造的改变与生命活动息息相关，分子动力学〔MD〕模拟被广泛应用于研究DNA构造改变的机理[3,4].而分子动力学模拟的准确性极大程度上依靠于分子力场的精到准确程度。CHARMM[5,6],AMBER[7,8],GROMOS[9]及OPLS[10]等力场应用于生物大分子体系已较为成熟。对于DNA,应用较广泛的力场包括CHARMM和AMBER力场，华而不实AMBER力场中的bsc1[7]和OL15[8]力场作为描绘叙述DNA的最新版力场被广泛应用。随着力场的不断优化与发展，其参数的准确性不断提高，但缺乏仍然存在，因而发现并改良现有力场的缺乏是优化力场参数的目的。已有研究结果表示清楚，CHARMM力场描绘叙述DNA,其构造倾向于A构型，AMBER则稳定于B构型[11].因而，不同力场的应用范围存在差异，通常，在水环境中AMBER力场能够更好地描绘叙述B构型，但在水活度较低的环境中，CHARMM力场对A-DNA的描绘叙述则愈加准确，与实验结果相吻合[12,13].固然不同力场下的构造信息研究得比拟具体，但对该现象的本质解释并不特别清楚。同时，随着力场参数的发展，不断有新的力场被提出并广泛应用[7,8].为了深切进入研究最新版CHARMM和AMBER力场对DNA构型转变的影响，本文采用分子模拟结合自由能计算的方式方法[14,15]比拟了不同力场下水溶液中B-DNA到A-DNA转化的自由能变化。通过自由能计算能够实现对B构型到A构型转变经过的分析，结果表示清楚，CHARMM力场下的DNA小沟较宽，与实验结果存在差异，采用径向分布函数对该现象进行分析，发现不同力场下DNA周围的离子分布存在较大区别，不同的离子分布是造成DNA稳定构造存在差异的主要因素。1、理论和方式方法1.1模型建立B-DNA和A-DNA的初始构造由AMBER提供的NucleicAcidBuilder插件构建，序列分别为GCGCGC和ATATAT.将构建好的DNA模型置于带有周期性边界条件的水立方盒子中，每个体系参加10个K+作为抗衡离子，盒子大小约为6.1nm6.0nm6.4nm,每个体系的原子数均为20000个左右。用于模拟自然的生理环境。AT序列〔PDBcode:4J2I〕以及GC序列〔PDBcode:1QC1〕的晶体构造参考自蛋白质晶体构造数据库[16,17].1.2动力学模拟采用NAMD2.12软件进行MD模拟[18],分别使用CHARMM36力场、AMBER-OL15力场和AMBER-bsc1力场参数描绘叙述DNA;使用TIP3P模型[19]中的参数描绘叙述水分子。采用SKAKE/RATTLE算法[20,21]将非水分子中含有氢原子的共价键的长度限制在其平衡值，采用SETTLE算法[22]保持水分子的刚性。采用恒温恒压的Langevin动力学方式方法和Langevin活塞方式方法[23]将温度和压力分别控制在300K和1.01105Pa.范德华截断半径为1.2nm,长程静电互相作用采用粒子网格埃瓦尔德〔PME〕方式方法计算。对运动方程积分的时间步长为2fs,每个体系在自由能计算之前均经历了2000步的能量最小化经过。从自由能计算结果中选取最稳定构造进行10ns的平衡模拟，应用该轨迹进行DNA构造的参数分析及抗衡离子和水分子的径向分布分析。采用VMD软件[24]进行轨迹分析。应用CURVES+软件[25]进行DNA构造参数分析，分析经过中DNA两末端分别除去两对碱基，选取中间残基进行分析。1.3自由能计算采用扩展自适应偏置力〔eABF〕方式方法计算了转化经过的自由能变化[15],以RMSD〔RMSDB-RMSDA〕作为该转化经过的反响坐标[26].选取DNA上所有重原子作为RMSD计算的参考原子，标准B-DNA与标准A-DNA的RMSD差值为0.3nm,-0.3RMSD0.3.当RMSD0时，DNA的构造偏向于B构型；当RMSD0时，DNA构造偏向于A构型。为了提高计算效率，反响途径均被分为2个窗口：-0.3RMSD0和0RMSD0.3.每个窗口的模拟时间为150ns,每个体系的模拟时间为300ns,总模拟时间为1.2s.2、结果与讨论2.1自由能曲线图2〔A〕给出了ATATAT序列B-DNA到A-DNA转变的自由能变化。对于AMBRER力场〔包括bsc1力场和OL15力场〕,最稳定的构造为B-DNA,bsc1力场与OL15力场对应的RMSD值分别为-0.14和-0.15nm,如此图3〔A〕和〔C〕所示，其构造非常类似。对于CHARMM力场，存在一个范围较广的低能区，由提取能量最低点的构造[图3〔B〕]可见，DNA两末端的氢键毁坏较为严重，氢键的毁坏使末端碱基运动较为灵敏，这是导致其稳定构造RMSD值范围较广的原因之一。图2〔B〕给出了GCGCGC序列的B-DNA到A-DNA转变的自由能变化。可见，在水环境中不同力场下DNA的最稳定构造存在差异。如此图4〔A〕和〔C〕所示，对于AMBER力场无论是bsc1力场还是OL15力场，B-DNA到A-DNA转化的经过中最稳定的构造为B-DNA,对应的RMSD值分别为-0.12和-0.13nm.对于CHARMM力场，最稳定的构造则偏离标准的B-DNA,表现为介于B构型与A构型中间的构造[图4〔B〕],对应的RMSD值为-0.02nm,与实验上表现出来的标准B-DNA存在差异。以上结果显示，AMBER力场在描绘叙述水溶液中B-DNA到A-DNA转化的最稳定构造愈加准确，该结果与旧版力场的模拟结果类似[12,13],且bsc1力场与OL15力场结果几乎一样，因而，下面只选取bsc1力场进行分析。2.2DNA构造参数分析为了研究在不同力场下DNA全局能量最小点的稳定构造，提取出全局能量最低点的代表构造进行平衡模拟并分析其构造参数。将不同序列DNA的构造参数与晶体构造进行了比拟〔表S1,见本文支持信息〕,结果表示清楚，对于〔AT〕6序列，AMBER力场下的绝大部分构造参数能够更接近实验结果。而对于〔GC〕6序列，CHARMM和AMBER力场对不同参数的精到准确程度不同。华而不实不同力场差距较大的是大沟和小沟的宽度，为了进一步探究其原因，进行了深切进入分析。为了进一步探究DNA在不同力场下的构造差异，对DNA大沟和小沟的宽度进行了统计分析。图5描绘叙述了不同力场下不同碱基序列的大沟宽度分布情况，发现对于〔AT〕6序列，bsc1力场下的大沟宽度分布范围为1.1~1.4nm,与晶体构造基本吻合；但CHARMM力场则出现较窄的大沟宽度，与晶体构造存在差异[图5〔A〕].对于〔GC〕6序列，bsc1力场下的大沟宽度相较于晶体构造略宽，CHARMM力场下的大沟宽度与晶体构造较为类似[图5〔B〕].图6描绘叙述了DNA小沟宽度的分布情况，〔AT〕6序列晶体构造的小沟宽度分布范围在0.3~0.6nm之间，bsc1力场从一定程度上能够反映晶体构造的小沟宽度，但是仍然存在小沟宽度偏大的现象；但对于CHARMM力场，小沟宽度分布在0.6~1.3nm之间，该结果完全偏离晶体构造提供的信息[图6〔A〕].对于〔GC〕6序列，bsc1力场能从一定程度上反映晶体构造的小沟宽度，但CHARMM力场下的小沟宽度在分布上明显右移，呈现较宽的小沟宽度[图6〔B〕].可见，bsc1力场对DNA小沟宽度的描绘叙述更为准确，而CHARMM力场下不同序列的小沟宽度均大于晶体构造。由此可见，AMBER力场〔bsc1〕参数下DNA大沟和小沟的宽度明显比CHARMM力场愈加精到准确。2.3抗衡离子对DNA构造的影响为了探究不同力场下大沟和小沟宽度存在较大差异的原因，对DNA周围抗衡离子的分布进行了探究。图7和图8分别描绘叙述了〔GC〕6序列与〔AT〕6序列在不同力场下DNA骨架、大沟和小沟周围的离子分布情况。由图7〔A〕,〔C〕和图8〔A〕,〔C〕可见，无论是〔GC〕6序列还是〔AT〕6序列，bsc1力场下骨架及小沟周围的离子分布密度均明显高于CHARMM力场。bsc1力场下严密的离子分布更好地中和了磷酸骨架上的负电荷，降低了2条磷酸骨架之间的静电排挤作用，因而与CHARMM力场相比，bsc1力场下的DNA构造小沟的宽度较窄。可见，CHARMM力场不能准确描绘叙述小沟宽度的主要因素是，在CHARMM力场下阳离子在小沟的分布较少，导致带负电荷的磷酸骨架静电排挤较大，小沟变宽。由图7〔B〕和图8〔B〕可见，对于一样序列，bsc1力场和CHARMM力场下离子在DNA大沟周围的分布密度差异不同不大，然而〔GC〕6序列与〔AT〕6序列之间却存在较大差异。根据文献[27]报道，对于〔GC〕6序列的DNA,离子在大沟的分布密度高于〔AT〕6序列，这在模拟中也得到了很好的验证。同时，小沟周围的离子分布对于不同的碱基序列也差异不同较大，由图7〔C〕和图8〔C〕可见，对于〔AT〕6序列的DNA,离子在小沟分布的倾向明显高于〔GC〕6序列，离子在小沟周围的严密分布更好地中和了磷酸骨架的负电荷，降低了静电排挤作用，使小沟变窄，这很好地解释了实验上〔AT〕6序列DNA的小沟宽度小于〔GC〕6序列的现象。可见，离子在DNA周围的分布情况是影响DNA大沟和小沟宽度的重要因素。2.4水分子对DNA构造的影响DNA周围水分子的分布情况同样对DNA的构造具有较大的影响。对于〔GC〕6序列，不同力场下DNA周围水分子的分布相差不大〔图9〕;但对于〔AT〕6序列，CHARMM力场下DNA的大沟和小沟周围的水分子分布较多〔图10〕.为了探究水分子对DNA构造的影响，我们提取了平衡模拟中DNA的典型构造，由图3〔B〕可见，DNA末端碱基氢键被水分子毁坏较为严重。因而，对于CHARMM力场，〔AT〕6序列末端碱基的氢键易被周围的水分子毁坏，这也是〔AT〕6序列大沟和小沟周围水分子分布较多的原因。对于碱基之间存在3条氢键的〔GC〕6序列，相比于构成2条氢键的〔AT〕6序列愈加稳定，CHARMM力场与bsc1力场下DNA周围水分子的分布情况则没有较大区别。3、结论利用分子动力学模拟结合自由能计算的方式方法比拟了CHARMM和AMBER力场〔包括bsc1力场和OL15力场〕对水溶液中B-DNA到A-DNA转化经过的影响，分别得到了不同力场下最稳定的DNA构造。结果表示清楚，在水环境中不同力场下DNA的最稳定构造存在差异，CHARMM力场下DNA最稳定构造的小沟较宽，介于B构型与A构型之间；AMBER力场下的DNA小沟较窄，稳定于B构型。进一步探究了不同力场下DNA稳定构造存在差异的原因，首先分析了DNA周围离子的分布情况。结果表示清楚，CHARMM力场下DNA小沟周围的离子密度明显低于AMBER力场，不能很好地抵消磷酸骨架的排挤作用，是CHARMM力场小沟较宽且趋向A构型的主要原因。进一步分析了水分子的分布对DNA构造的影响，表示清楚CHARMM力场下的〔AT〕6序列末端碱基之间的氢键被水分子毁坏较为严重，与AMBER力场相比，CHARMM力场下的氢键互相作用较弱。在模拟生理环境的水溶液中AMBER力场能够更好地描绘叙述实验现象。研究结果为分子模拟实验设计提供了理论指导，并对进一步的参数优化提供了帮助。以下为参考文献[1]SaengerW.,DefiningTermsfortheNucleicAcids,Springer,NewYork,1984[2]IvanovV.I.,MinchenkovaL.E.,MinyatE.E.,Frank-KamenetskiiM.D.,SchyolkinaA.K.,J.Mol.Biol.,1974,87〔4〕,817-833[3]MukherjeeA.,LaveryR.,BagchiB.,HynesJ.T.,J.Am.Chem.Soc.,2008,130〔30〕,9747-9755[4]DumontE.,WibowoM.,Roca-SanjunD.,GaravelliM.,AssfeldX.,MonariA.,J.Phys.Chem.Lett.,2021,6〔4〕,576-580[5]FoloppeN.,MacKerellA.D.Jr.,J.Comput.Chem.,2000,21〔2〕,86-104[6]MacKerellA.D.Jr.,BashfordD.,BellottM.,DunbrackR.L.Jr.,EvanseckJ.D.,FieldM.J.,FischerS.,GaoJ.,GuoH.,HaS.,J.Phys.Chem.B,1998,102〔18〕,3586-3616[7]IvaniI.,DansP.D.,NoyA.,PrezA.,FaustinoI.,HospitalA.,WaltherJ.,AndrioP.,Go1iR.,BalaceanuA.,Nat.Methods,2021,13〔1〕,55[8]ZgarbovM.,LuqueF.J.,ponerJ.I.,CheathamIIIT.E.,OtyepkaM.,JureckaP.,J.Chem.TheoryComput.,2020,9〔5〕,2339-2354[9]SoaresT.A.,HnenbergerP.H.,KastenholzM.A.,KrutlerV.,LenzT.,LinsR.D.,OostenbrinkC.,vanGunsterenW.F.,J.Comput.Chem.,2005,26〔7〕,725-737[10]RobertsonM.J.,Tirado-RivesJ.,JorgensenW.L.,J.Chem.TheoryComput.,2021,11〔7〕,3499-3509[11]SongC.,XiaY.Y.,ZhaoM.W.,LiuX.D.,LiF.,JiY.J.,HuangB.D.,YinY.Y.,J.Mol.Model.,2006,12〔3〕,249-254[12]CheathamT.E.,CrowleyM.F.,FoxT.,KollmanP.A.,Proc.Natl.Acad.Sci.USA,1997,94〔18〕,9626-9630[13]YangL.,PettittB.M.,J.Phys.Chem.,1996,100〔7〕,2564-2566[14]ShenW.,ShaoX.G.,CaiW.S.,Chem.J.ChineseUniversities,2021,37〔10〕,1809-1816〔沈文，邵学广，蔡文生。高等学校化学学报，2021,37〔10〕,1809-1816〕[15]FuH.H.,ShaoX.G.,ChipotC.,CaiW.S.,J.Chem.TheoryComput.,2021,12〔8〕,3506-3513[16]Acosta-ReyesF.J.,SubiranaJ.A.,PousJ.,Snchez-GiraldoR.,CondomN.,BaldiniR.,MalininaL.,CamposJ.L.,Biopolymers,2021,103〔3〕,123-133[17]TimsitY.,MorasD.,EMBOJ.,1994,13〔12〕,2737[18