自由能计算揭示苏氨酸对抗冻蛋白与冰晶结合能力的影响

崔韶丽，张维佳，邵学广，2，蔡文生

（1.南开大学分析科学研究中心,有机新物质创造前沿科学中心,化学学院,天津市生物传感与分子识别重点实验室，2.药物化学生物学国家重点实验室,天津300071）

抗冻蛋白（Antifreeze proteins，AFPs），又称为冰结合蛋白，是一类能控制冰晶生长和抑制冰晶重结晶的蛋白质，可以保护生物体在零度以下的环境中生存［1，2］.AFPs主要存在于昆虫、微生物、植物和鱼类等生物体内［3，4］.通常，昆虫AFPs和微生物AFPs被认为是超活性抗冻蛋白，而植物AFPs和鱼类AFPs则被认为是中等活性抗冻蛋白［5，6］.AFPs由于具有优异的抗冻性能被应用于组织和器官及食品工业的低温保存领域［7，8］，因此AFPs抗冻机制［5，6，9］和其中的AFPs与冰-水界面相互作用［10，11］吸引了不少研究者的兴趣.

近年来，研究者们提出锚定包合物模型来描述AFPs与冰-水界面的结合［12，13］.该模型认为AFPs通过携带自己的“冰”（即锚定包合水，Anchored clathrate water，ACW）与冰晶结合，即在冰结合位点（Icebinding site，IBS）通过疏水作用使液态水分子按晶格方式排列，并通过氢键锚定液态水分子晶格，锚定的液态水分子晶格又通过匹配特定的冰面使AFPs与冰晶结合.Garnham等［12］认为该模型适用于任何的AFP，不论其抗冻活性水平如何.超活性昆虫抗冻蛋白（Hyperactive insect AFPs，hypAFPs），如TenebriomolitorAFP（TmAFP），Choristoneura fumiferanaAFP等的IBS具有平整的结构且含有较多的苏氨酸残基（Thr），IBS基序为Thr-x-Thr（TxT，x是非保守的氨基酸，主要是疏水氨基酸）.通过Thr残基的两亲性基团锚定与冰结合［14，15］，吸附在冰晶的基面和至少一面与基面正交的晶面.苏氨酸重复距离和基面氧原子重复距离形成完美匹配［图1（A）］［13，16］.Scotter等［17］发现能够与基面结合是AFPs具有超活性的基础.Bar等［18］对TmAFP的突变研究表明IBS上Thr被丝氨酸或缬氨酸取代，其活性会降低很多，说明Thr在蛋白与冰晶相互识别和结合过程中起到很重要的作用［19］.中等活性黑麦草抗冻蛋白（Lolium perenneAFP，LpAFP）与TmAFP的尺寸相似［图1（B）］［20］，其IBS也具有疏水且平整的结构，与冰晶基面结合，却具有中等的抗冻活性.然而，如果突变野生型LpAFP的IBS使其具有较多的Thr和具有完美Thr-x-Thr结构，其与冰结合的能力和活性水平还需要进一步明确.

Fig.1 Schematic representation of Tm AFP(PDB ID:1EZG)aligned to the basal plane of ice(A),crystal structures of Tm AFP and Lp AFP(PDB ID:3ULT)(B)

分子动力学（Molecular dynamics，MD）模拟是在原子水平上揭示分子之间相互作用机制的一种常用而有效的方法［21~26］.如Hudait等［27］利用MD模拟研究了超活性抗冻蛋白TmAFP冰结合位点处的氢键和疏水基团对TmAFP-冰结合自由能的贡献，结果表明，AFPs与冰结合时氢键和疏水基团同等重要，Mochizuki等［28］观察到抗冻糖蛋白通过疏水基团与冰可逆地结合，Hudait等［13］还计算了TmAFP与冰面的结合自由能，发现其与冰晶基面结合能力最强.本文通过突变LpAFP的IBS上的部分非苏氨酸残基，将其IBS上苏氨酸含量逐步增加，其中包括一个使每个β片段均具有TxT基序的突变，利用MD模拟LpAFP及其突变体与冰晶的结合，并进一步计算结合过程的自由能曲线来研究了Thr含量对其吸附冰晶能力和抗冻活性的影响.

1 理论和方法

1.1 模型建立

从蛋白质晶体结构数据库中获得LpAFP的晶体结构（PDB ID：3ULT），利用CHARMM-GUI［29］建立4个LpAFP突变体的结构，每个突变体中突变的残基如表1所示，突变的目标氨基酸均为苏氨酸.目的是逐步提高苏氨酸的含量.所有冰的模型都是利用GenIce［30］软件生成.构建的带有周期性边界条件的AFP-冰-水模型如图2（A）所示，AFP的IBS距离冰晶基面为1 nm，且IBS上两排氨基酸平行于冰晶基面的a1轴［图2（B）］，每个体系中加入3个Na+作为抗衡离子，为了避免其对水结冰动力学的影响，将其约束在盒子的边缘［31］.为了防止冰晶的融化，将冰晶约束在初始的位置，力常数均为1000 kJ·mol‒1·nm‒2.所建模型的细节见表1.

Table 1 Details of the molecular assemblies investigated in this study

Fig.2 Lp AFP/mutant⁃ice⁃water assemblies for MD(A),schematic diagram of the a⁃axis of ice(B)and Lp AFP/mutant⁃ice⁃water assemblies for free⁃energy calculations(C)

1.2 动力学模拟

所有MD模拟均采用GROMACS 2020.2软件［32］，使用CHARMM27力场［33］参数描述AFPs，使用TIP4P/ICE模型［34］中的参数描述水分子，根据文献［26，35］报道，CHARMM27力场和TIP4P/ICE水模型相结合可用于描述AFPs的吸附过程.使用V-rescale［36］和Parrinello-Rahman［37］控制温度和压强分别为268 K和1×105Pa.采用LINCS算法［38］将所有涉及氢原子的共价键限制在其平衡长度.用2 fs时间步长对运动方程进行积分，范德华和短程静电相互作用的截断半径为1.2 nm，长程静电相互作用采用粒子网格埃瓦尔德（Particle Mesh Ewald，PME）方法［39］计算.首先采用最速下降算法，使每个体系的初始能量最小化50000步.然后，在正则系综（NVT）中进行了100 ps模拟，在等温等压系综（NPT）中进行了100 ps模拟，最后在NPT系综中进行了400 ns的MD模拟.使用VMD 1.9.3［40］进行MD轨迹的可视化和分析.

1.3 自由能计算

自适应偏置力（ABF）［41，42］方法是一种基于热力学积分的重要采样方法，通过施加在由集合变量构成的反应坐标方向上的偏置力来增强该方向上的采样，具有严谨的收敛性，但是在使用中存在很多局限，如该方法需要实现集合变量间的相互独立，并与体系中存在的约束正交等大大限制了其在复杂体系中的应用.因此，我们［43］实现了扩展自适应偏置力（eABF）的计算方法基本解决了经典ABF方法的缺点.本文采用本课题组最新发展的WTM-eABF方法［44，45］来计算AFPs与冰结合过程的自由能曲线，该方法结合了Well-tempering metadynamics［46，47］和eABF［43，48］“填谷”和“削峰”的优点，极大地提高了采样效率，并已得到了广泛应用［44，49］.从平衡轨迹中获得AFP吸附时的构象来构建如图2（C）所示的模型，用于利用WTM-eABF方法计算AFP与冰晶基面可逆吸附的自由能变化.为了防止冰的生长造成吸附过程的不可逆，将温度控制在280 K（高于熔点10 K）［28］.另一方面，为了防止冰晶的融化，将冰晶约束在其初始的位置，力常数为1000 kJ·mol‒1·nm‒2.选择AFPs的质心与冰中氧原子的质心之间的距离在y轴（即沿冰生长的方向）上的投影作为吸附过程的反应坐标d（其取值范围为1.9 nm≤d≤3.3 nm）.用来计算自由能的每个体系模拟时间为800 ns，总的模拟时间为2.4μs.

2 结果与讨论

2.1 Lp AFP及其突变体IBS周围液态水分子的分布

液态水分子的分布对AFPs与冰的结合具有较大的影响，为了分析LpAFP及其突变体IBS周围液态水分子的分布情况，计算了水分子相对于IBS氨基酸的径向分布函数，如图3（A）所示.可见，LpAFP和不同突变体对峰的位置和形状无明显的影响，但峰的高度变化明显.苏氨酸越多其周围液态水分子分布就越多.这一结果表明提高IBS上苏氨酸的含量有助于结合更多的水分子，与已报道的苏氨酸的甲基形成疏水水笼，羟基形成氢键锚定水分子的结论吻合［13，15，27］，进而可能有利于抗冻蛋白与生长的冰晶相结合.

Fig.3 Radial distribution function of water(oxygen atoms)around the IBS of Lp AFP and its mutants obtained from MD simulations(A),hydrogen bond lifetime correlation functions between any water molecules and amino acid residues at the IBS for wild⁃type Lp AFP and its mutants adsorbed on ice crystals(B)

2.2 Lp AFP及其突变体的ACW稳定存在时间的比较

ACW的稳定存在主要依靠AFPs与冰晶之间形成的氢键［50］.为了研究LpAFP及其突变体形成的ACW的稳定性差异，计算了AFP吸附在冰晶上时，IBS的残基与所有水分子（包含新生长的冰）之间的氢键寿命自相关函数C（t）［51］，定义如下：

当冰结合位点的氨基酸和所有水分子（包含新生长的冰）之间存在氢键时，h（t）为1；当不存在氢键时，h（t）为0.这种自相关函数可以反映氢键能够维持的时间以及随时间的衰减速度，结果如图3（B）所示.可见，随着IBS上苏氨酸数量的增多，C（t）衰减速度减慢，其中突变体4的C（t）衰减速度最慢，说明其与冰晶形成氢键的稳定性最强.结果表明，苏氨酸的含量越高，越有利于氢键的形成.

2.3 自由能的比较

根据前面的分析，如果IBS上具有完美重复的TxT基序（Mutant 4），推测其与冰晶的结合力最强.为了进一步证明此结论，选择了Mutant 2，Mutant 4和LpAFP为例，计算了它们与冰结合过程的自由能曲线，并分析了吸附时IBS与冰晶表面形成的氢键网络.图4（A）给出了描述LpAFP及其突变体2和4与冰晶结合过程的自由能曲线.可见，IBS上苏氨酸含量越高，AFPs与冰晶的结合能力越强，且具有规整TxT基序时与冰晶的结合能力显著增强.推测这是由于苏氨酸上甲基形成了疏水水笼［图4（A）插图］，羟基形成氢键锚定液态水分子，且IBS上苏氨酸重复距离与冰晶基面重复氧原子之间距离相匹配［13，16］，使得粘合AFP与冰面的氢键网络的稳定性明显增强.Chakraborty等［52］计算了野生型AFPIII及其几个突变体与冰的结合自由能，发现结合自由能的大小与实验观察到的热滞活性高低呈线性相关，因此推测图4（A）中抗冻活性高低顺序为Mutant 4>Mutant 2>LpAFP.图4（B）~（D）为自由能曲线随模拟时间的演化，表明了在800 ns内计算已经收敛.

Fig.4 Free⁃energy profiles characterizing the adsorption process for Lp AFP,Mutant 2 and Mutant 4(A),free⁃energy profiles at different simulation times for the Lp AFP(B),Mutant 2(C)and Mutant 4(D)

2.4 Lp AFP及其突变体吸附在冰面时形成的氢键网络的比较

为了更直观地展示突变体与冰晶结合能力的提高，对IBS周围0.5 nm之内水分子与IBS的氨基酸、通道水分子（Channel water，CW）形成的氢键网络进行了分析，如图5所示.其中，CW是ACW的一部分，参与AFPs与冰结合时形成的氢键网络.可见，野生型的LpAFP虽然可以吸附在冰晶的基面，但是CW个数较少，与冰结合时形成的氢键网络不完整，所以其结合能力并不强［图5（A）］.突变体2形成的氢键网络已较为完整［图5（B）］，突变体4则可以形成紧密完整的氢键网络［图5（C）］，CW个数较多且从模拟轨迹上可以看出其稳定存在的时间变长，CW的稳定存在为突变体与冰晶的强结合提供了基础［50］.总之，苏氨酸含量越高，其IBS与冰面之间的氢键网络也越稳定，且具有规整TxT重复基序会使形成的氢键网络最稳定.

Fig.5 Hydrogen bond structures of water molecules within 0.5 nm below the IBS of Lp AFP(A),Mutant 2(B)and Mutant 4(C)adsorbed on ice crystals

3 结 论

利用分子动力学模拟结合自由能计算的方法，比较了中等活性LpAFP的IBS苏氨酸含量对其与冰晶结合能力的影响.结果表明，IBS上苏氨酸含量越高，LpAFP与冰晶的结合在能量上越有利；当突变IBS使其具有与hypAFP相同的TxT基序时，其与冰的结合能力最强.进一步探究了其与冰晶结合能力增强的原因，结果表明，IBS周围液态水分子分布随苏氨酸含量的增加而增强.进一步分析了ACW稳定存在的时间和IBS与冰面之间的氢键网络都表明AFPs与冰晶的结合能力随苏氨酸含量的增加而更稳定.由此得出，增加苏氨酸的含量是提高中等活性抗冻蛋白抗冻活性的有效方法.研究结果为冰生长抑制剂的设计提供了一套计算方法和理论指导.