try ai
科普
编辑
分享
反馈
  • 计算化学生物学

计算化学生物学

SciencePedia玻尔百科
核心要点
  • 玻恩-奥本海默近似是计算生物学的基石,它通过分离快速运动的电子和缓慢运动的原子核的运动来简化模拟。
  • 模拟通过使用确定性的分子动力学或概率性的蒙特卡洛等方法来探索分子体系的势能面,以预测分子行为、结构和动力学。
  • 混合QM/MM方法通过用量子力学处理体系中小的、反应性的核心区域,而用计算成本较低的经典物理学处理较大的环境,从而能够研究化学反应。
  • 计算方法在现代药物设计中不可或缺,它采用QSAR、虚拟筛选和基于片段的发现等技术来识别和优化潜在的治疗分子。
  • 整合建模通过将计算系综与来自核磁共振(NMR)和X射线小角散射(SAXS)等技术的实验数据相结合,创建出更准确、更动态的分子体系图像。

引言

观察和预测分子行为的能力是现代生物学和医学的基础。然而,生命的基本组成部分——蛋白质、DNA和小分子——其错综复杂的舞蹈发生在时间和空间尺度上,快到、小到不可能直接观察。计算化学生物学应运而生,为探索这个看不见的世界提供了一台虚拟显微镜。它所面对的核心问题是巨大的:支配分子的量子力学规则极其复杂,对于任何大于几个原子的体系都无法求解。本文通过解释使模拟生物学成为可能的巧妙近似和强大算法,来弥合这一差距。在第一章“原理与机制”中,我们将深入探讨基础理论,如玻恩-奥本海默近似和势能面的概念,并探索分子动力学和蒙特卡洛等核心模拟方法。接下来,“应用与跨学科联系”一章将展示这些工具如何应用于实际问题,从理性药物设计到复杂实验数据的解释,彰显计算对生命科学的变革性影响。

原理与机制

为了在分子水平上模拟生命错综复杂的舞蹈,我们首先必须面对一个相当不便的真相:世界是量子力学的。一个蛋白质、一条DNA链、一个水分子——每一个都是由原子核和电子组成的沸腾集合,受制于众所周知的奇异量子力学规则。对于任何大于几个原子的体系,直接、全面的模拟都是一场计算上的噩梦,即使是我们最强大的超级计算机也无法企及。那么,我们该如何着手呢?计算化学生物学的艺术在于一系列绝妙的近似和巧妙的视角,它们使我们能够构建一个易于处理但又在物理上忠实于现实的虚拟世界。我们的旅程始于这些思想中最重要的一个。

两种时间尺度的故事:量子世界与经典世界

想象一下,试图拍摄一只悬停在乌龟旁的蜂鸟。蜂鸟的翅膀是一片模糊,在乌龟迈出蹒跚一步的时间里,它已经振动了数千次。原子的世界与此非常相似。电子是蜂鸟;原子核是乌龟。一个电子的质量至少比最轻的原子核(单个质子)小1836倍。由于分子中的所有粒子都感受到相当的电力,牛顿第二定律(F=maF=maF=ma)告诉我们,轻巧的电子比笨重的原子核加速和移动得快得多得多。

这种巨大的时间尺度差异是关键。电子完成其轨道并重新调整其构象的速度如此之快,以至于从它们的角度看,原子核基本上是固定不动的。反过来,当原子核缓慢地移动时,它们感受到的不是单个电子的瞬时位置,而是来自一个经过时间平均、弥散开的负电荷云的力。

这一洞见被形式化为​​玻恩-奥本海默近似(BOA)​​,它几乎是所有计算化学生物学赖以建立的基石。BOA允许我们将电子的运动与原子核的运动解耦。这是两者之间的一项“君子协定”:

  1. 首先,我们将原子核固定在某个几何构型上。然后,我们求解电子在这些固定原子核的静电场中运动的纯粹量子力学问题。这为我们提供了该特定核排布下的电子基态能量。
  2. 然后,我们对原子核的每一种可能排布重复此过程。

通过这样做,我们创建了一个映射:对于任何给定的核位置集合,我们都有一个相应的能量。这个映射就是我们新的、简化的宇宙。我们用一个完整、纠缠的量子问题换来了一个简单得多的问题:经典原子核在一个预先计算好的能量景观上运动。从本质上讲,我们已经搭建好了化学戏剧即将上演的舞台。虽然BOA对于大多数基态生物过程非常稳健,但在某些情况下它可能会彻底失效,例如在光化学反应期间,或在称为锥形交叉点的特殊几何构型附近,此时电子态之间过于接近,电子无法再瞬时调整。

生命的舞台:势能面

玻恩-奥本海默近似的结果是一个概念上的奇迹:​​势能面(PES)​​。想象在一个与原子坐标数量一样多维度的空间中,有一个广阔起伏的景观。这个景观上任意一点的高度代表了当原子以该特定构型排列时体系的势能。这个景观就是所有分子行动的舞台。

PES的地理形态就是分子的化学性质。

  • ​​谷和盆地:​​ 景观中的深谷对应于低能量的稳定状态。这些是蛋白质的折叠结构、DNA的双螺旋,或是紧密契合在结合口袋中的药物分子。
  • ​​山隘:​​ 连接一个谷到另一个谷的最低能量路径必须翻越一个“山隘”。这些被称为​​鞍点​​的通道代表了化学反应或构象变化的过渡态。这个山隘的高度就是决定过程速率的活化能垒。
  • ​​曲率:​​ 景观的形状告诉我们分子的动力学信息。一个谷壁的陡峭程度——即其曲率——决定了分子振动的频率。一个狭窄、陡峭的峡谷对应于刚性的高频振动(如共价键的伸缩),而一个宽阔、平浅的盆地则对应于柔性的低频运动。

我们可以用微积分的工具来数学地描述这个景观。PES上任意一点的梯度(斜率)的负值给出了作用在每个原子上的力,将体系推向“下坡”方向。二阶导数被封装在一个称为​​海森矩阵(Hessian)​​的矩阵中,描述了局部曲率。在稳定谷底,海森矩阵是正定的,意味着景观在所有方向上都向上弯曲,像一个碗。在过渡态,它在除一个方向外的所有方向上都向上弯曲,而这个方向指向反应路径。

原子之舞:模拟运动

舞台已经搭好,我们如何让演员——原子——动起来呢?有两种主要的哲学思想,两种不同的原子之舞编排。

分子动力学:发条宇宙

第一种方法,​​分子动力学(MD)​​,将体系视为一个微型发条机器。原子被当作由弹簧连接的经典球体,根据牛顿运动定律运动。在每一刻,我们通过寻找势能面的下坡斜率来计算每个原子上的力。一个微小的推动将体系带到一个新的构型,然后我们重新计算力并重复此过程。通过将数百万个这样的小步骤串联起来,我们生成了一条轨迹——一部原子运动的电影。

但这有一个难题:“时间步长的暴政”。PES不是平滑的;它具有各种尺度的特征。我们体系中最硬的弹簧是共价键,尤其是那些涉及轻氢原子的共价键。这些键以大约10飞秒(10−1410^{-14}10−14 s)的周期振动。为了精确捕捉这种狂热的运动,我们的积分时间步长 Δt\Delta tΔt 必须更小,通常约为1飞秒。这个严苛的限制意味着,即使是一毫秒长的模拟——在生物学意义上只是一眨眼的功夫——也需要惊人的 101210^{12}1012 个计算步骤。为了克服这个问题,人们使用了一些巧妙的技术,例如通过算法冻结这些快速振动(使用像SHAKE这样的约束),或者人为地增加氢原子的质量来减慢它们的速度,这可以允许使用更大、更高效的时间步长。

蒙特卡洛:醉汉的漫步

第二种哲学思想,​​蒙特卡洛(MC)模拟​​,采用了一种不同的、统计学的方法。我们不是计算力并遵循确定性路径,而是通过“醉汉漫步”的方式来探索能量景观。我们从某个构型开始,并提出一个小的、随机的移动。然后我们参考能量景观来决定是否接受这个新步骤。这个决定由巧妙的 ​​Metropolis算法​​ 支配:

  1. 如果提议的移动将体系带到更低的能量(下坡),我们总是接受它。
  2. 如果提议的移动将体系带到更高的能量(上坡),我们仍可能接受它。接受这个能量上不利的移动的概率取决于温度和能量惩罚的大小。

这个“有时接受上坡移动”的规则是该方法的精髓所在。它防止模拟陷在它找到的第一个谷中,并允许它探索整个景观,甚至跨越能垒。经过许多步之后,这种随机行走并非漫无目的;它保证会根据​​玻尔兹曼分布​​所描述的热力学概率来抽样构象。这个分布告诉我们,在给定温度下,一个体系将以与其能量成指数递减的概率占据各个状态,P(i)∝exp⁡(−Ei/kBT)P(i) \propto \exp(-E_i / k_B T)P(i)∝exp(−Ei​/kB​T)。通过从这种行走中收集统计数据,我们可以通过用每个构象的玻尔兹曼概率对其贡献进行加权,来计算宏观的、实验可测量的平均值,例如肽中螺旋结构的平均含量。

搭建舞台:从量子纯粹性到经典实用主义

我们之前讨论PES时,仿佛它是已知的,但它实际上从何而来?创建这个能量景观是模拟中最关键、计算成本最高的部分。

经典方法:力场与溶剂

对于绝大多数大规模生物分子模拟而言,在每一步都求解完整的电子量子问题是不可能的。取而代之,我们使用一个巧妙的近似:​​经典力场​​。力场是一套模仿真实PES的经验性函数和参数。它将能量描述为简单项的总和:用于键长和键角的谐振子弹簧、用于扭转旋转的周期函数,以及用于非键相互作用的对相互作用项——用于短程排斥和长程吸引的Lennard-Jones势,以及用于静电相互作用的库仑定律。

体系能量的很大一部分来自其环境,对生物分子而言,这几乎总是水。对溶剂建模是一个关键的选择。在​​显式溶剂​​模型中,我们在模拟盒子中填充数千个独立的、可移动的水分子,每个水分子都与溶质及其他水分子相互作用。这种方法非常准确,但计算成本极高。另一种选择是​​隐式溶剂​​模型,我们用具有水的介电性质的连续介质来代替离散的水分子。这是一种更粗糙但效率高得多的近似,用原子细节换取了计算速度。

量子-经典混合方法:QM/MM

但是,如果我们想模拟一个化学反应,其中共价键正在断裂和形成,该怎么办呢?经典弹簧不再是有效的描述。对于这些情况,我们需要​​混合量子力学/分子力学(QM/MM)​​这种“两全其美”的方法。在这里,我们将体系划分:

  • ​​QM区域​​:体系中小的、化学活性核心(例如,活性位点中的底物和关键酶残基)用完整的、动态的量子力学计算来处理。
  • ​​MM区域​​:体系中庞大的其余部分(蛋白质的主体、溶剂)用经典力场来处理。

这两个区域通过静电作用进行最重要的沟通。QM区域的电子云被MM原子的电荷极化,反过来,MM原子也感受到由QM区域的原子核和电子云产生的电场。QM/MM是一种强大的多尺度技术,它使我们能够将计算火力精确地集中在最需要的地方。

终极奖赏:预测真实世界的可观测量

这套庞大的计算机器的目的是产生可以与真实世界实验相比较的数字。我们想要预测结构、计算反应速率,并确定药物与其靶点结合的紧密程度。

这种预测能力的一个绝佳例子是​​同源建模​​。如果我们有一个新蛋白质的氨基酸序列,我们如何预测其三维结构?我们可以在数据库中搜索进化上相关的蛋白质(同源物)的已知结构。其核心原理是,进化对蛋白质结构——其稳定的能量景观——的保守性远超其对序列的保守性。许多不同的序列可以成功地折叠成相同的低能形状。因此,一个微弱但统计上显著的共同祖先迹象,通常通过捕捉整个蛋白质家族进化模式的复杂基于谱的方法发现,是选择模板的更好指南,远胜于一个简单但可能纯属巧合的高比例相同氨基酸。

也许我们能计算的最重要的量是​​自由能(GGG)​​。它是化学的真正货币,支配着所有过程的自发性和平衡。与势能(UUU)不同,自由能包含了熵(SSS)的影响——熵是衡量无序程度或体系可用微观选项数量的指标(G=U+PV−TSG = U + PV - TSG=U+PV−TS)。计算药物与蛋白质结合的自由能是计算机辅助药物设计的“圣杯”。

直接模拟结合事件通常太慢。取而代之,我们使用一个深刻的热力学“会计技巧”:​​热力学循环​​。因为自由能是状态函数——其变化仅取决于起点和终点,而与所走的路径无关——我们可以构建一个封闭的变换环路。如果物理结合过程是循环中“困难”的一环,我们可以通过计算一系列更容易、非物理的“炼金术”环节的自由能变化来计算其 ΔG\Delta GΔG,这些环节共同完成了这个循环(例如,让一个配体在溶剂中“消失”,然后在结合位点“重现”)。由于围绕一个闭环的总 ΔG\Delta GΔG 必须为零,所以难以获得的值可以从容易获得的值中求出。

更值得注意的是,非平衡统计力学最近的突破,如​​Crooks涨落关系​​,使我们能够从快速、不可逆的过程中确定平衡自由能差。通过反复将配体从其结合口袋中拉出,并每次测量所需的功,我们可以分析功值的分布,从而精确地提取出平衡结合自由能。这是一个惊人的结果,是摩擦的耗散世界与热力学平衡的永恒领域之间的深刻联系,展示了基本原理在揭开分子世界秘密方面的持续力量。

应用与跨学科联系

在我们完成了对驱动计算化学生物学世界的基本原理和机制的探索之后,人们可能会想坐下来欣赏理论机器的优雅。但科学的核心不是一项观赏性运动。这些原理的真正美妙之处不在于其抽象的表述,而在于其解决实际问题、连接不同领域并引导我们走向新发现的力量。这正是理论联系实际的地方——或者,更贴切地说,是算法与活性位点相遇的地方。在本章中,我们将探讨我们所学的概念如何应用于整个科学领域,改变我们设计药物、解释实验以及理解生命分子之舞本质的方式。

分子识别的语言

在我们模拟生物过程之前,我们必须首先学会说它的语言:能量的语言。分子间的相互作用受制于力的微妙 interplay,我们的首要任务是将这些力转化为计算机能理解的定量能量函数。

在这场分子戏剧中,一个主导角色是静电力。你可能从基础化学中记得,像氨基酸这样的分子的电荷取决于其环境的酸度——即pHpHpH值。这不是一个次要的细节;它是决定分子行为的核心特征。通过应用简单的Henderson-Hasselbalch关系式,我们可以计算出在给定生理pHpHpH下蛋白质各部分的平均电荷。这个简单的计算是为模拟“参数化”分子的第一步。一个始于酸碱化学教科书练习的计算,最终上升为决定宏观结构的关键因素。对于像天然无序蛋白(IDPs)这样巨大而柔性的链,这些依赖于pH的电荷分布支配着整个构象系综。高的净正电荷或负电荷将导致链膨胀,因为同种电荷相互排斥,使蛋白质保持在伸展的动态状态。相反,在等电点——即净电荷为零的特定pHpHpH值——这些排斥力消失,可能使蛋白质塌陷甚至聚集。一个蛋白质的命运,其功能或功能障碍,可能取决于这种微妙的静电平衡。

能量语言中的另一个关键术语是去溶剂化。细胞中的分子并非处于真空中;它们畅游在水的海洋里。水是一种高极性溶剂,是带电和极性基团的绝佳稳定环境。当药物与蛋白质结合,或两个蛋白质相互对接时,它们必须将这些水分子推开。这样做的能量代价是什么?我们可以用源自Max Born的一个优美而简单的想法来估计它。通过将离子建模为带电球体,我们可以计算其在不同环境中的静电自能。将一个带电基团从水的高介电常数舒适环境(εw≈80\varepsilon_{\mathrm{w}} \approx 80εw​≈80)移动到蛋白质内部极性低得多的低介电常数环境(εeff≈4−10\varepsilon_{\mathrm{eff}} \approx 4-10εeff​≈4−10),会产生显著的能量惩罚。这种“去溶剂化惩罚”是分子识别中的一个主要驱动力。与其说是药物和其靶点之间有强大的吸引力,不如说是整个体系通过最小化那些不幸地从其水性家园中被“流放”的极性基团数量而获益。这就是疏水效应的本质,一种并非源于吸引力,而是源于与溶剂排斥的力。

寻找“正确”的形状

一旦我们有了能量函数,下一个巨大挑战便随之而来:找到与最低能量相对应的正确的原子三维排布,即构象。这就是“搜索问题”。

对于一个与蛋白质结合的小分子,这不仅涉及探索其位置和方向,还涉及其内部分子柔性。一个典型的类药分子可能有几个可旋转的键。如果我们以30∘30^{\circ}30∘的增量(每个键12个状态)对每个键的旋转进行采样,数字会以惊人的速度增长。对于一个只有7个这样键的分子,可能的构象数是12712^7127,将近3600万。检查每一个构象的能量在计算上是无法承受的。而且这甚至还没有考虑蛋白质靶标的柔性!这种“组合爆炸”或“维度灾难”是为什么暴力搜索不可能的原因。它解释了为什么计算生物学领域充满了巧妙的启发式搜索算法——遗传算法、模拟退火、蒙特卡洛方法——这些算法旨在智能地导航这些浩瀚无垠的搜索空间,而无需访问每一个点。

当我们的搜索算法提出一个结构——例如,一个蛋白质的预测模型——我们如何判断其质量?我们如何知道是否找到了一个好的答案?像均方根偏差(RMSD)这样简单的度量标准可能会产生误导。一个预测可能总体RMSD很低,但却搞错了蛋白质的基本“折叠”方式。学术界已经开发出更复杂的工具,如全局距离测试(GDT_TS)分数。GDT_TS不是进行单一比较,而是提出一系列问题:蛋白质残基中有多大比例在其正确位置的1 A˚1\,\text{\AA}1A˚范围内?有多大比例在2 A˚2\,\text{\AA}2A˚范围内?4 A˚4\,\text{\AA}4A˚范围内?8 A˚8\,\text{\AA}8A˚范围内?通过对这些多项选择题的结果进行平均,我们得到了一个更稳健、更有意义的预测质量度量。高的GDT_TS分数告诉我们,我们已经在多个分辨率上正确地捕捉到了蛋白质的拓扑结构,从局部细节到全局架构。这种严格的自我评估使得该领域能够取得真正的进步。

有目的地设计分子

有了这些计算能量和搜索结构的工具,我们可以从理解自然转向设计自然。这就是理性药物设计和虚拟筛选的领域。

在这个领域最强大的思想之一是定量构效关系(QSAR)。QSAR的目标是找到分子结构与其生物活性之间的数学关联。但是你如何将“结构”放入一个方程式中呢?诀窍是将复杂的化学特征提炼成数值“描述符”。例如,Kier-Hall连接性指数是一个巧妙编码分子碳骨架支化程度的数字。通过比较像正戊烷这样的线性分子与像新戊烷这样的高度支化分子的指数,我们可以看到该指数如何量化这一拓扑特征。通过为一组已知药物计算数十或数百个此类描述符,我们可以使用统计方法建立一个模型:活性=f(描述符1,描述符2,… )\text{活性} = f(\text{描述符}_1, \text{描述符}_2, \dots)活性=f(描述符1​,描述符2​,…)。然后,这个模型可以用来预测新的、未经测试的分子的活性,让化学家能够优先决定在实验室中合成和测试哪些分子。

一种更基于物理的药物设计方法是基于片段的先导化合物发现(FBLD)。FBLD不是试图一次性找到一个大的、高亲和力的分子,而是采取一种更耐心的策略。首先,我们找到两个或多个非常小的“片段”,它们微弱但高效地结合在蛋白质靶标上相邻的口袋里。然后,我们设计一个化学连接子将它们缝合成一个单一的、更大的配体。希望连接后的分子的亲和力会大于其各部分之和。然而,热力学告诉我们没有免费的午餐。将两个自由漂浮的片段连接成一个分子的行为本身就带有熵罚;我们限制了它们独立翻滚和平移的自由。此外,连接子本身可能会因扭曲以将片段置于其最佳结合姿态而引入应力。一次成功的FBLD努力,是组合片段的有利结合能足够强大,以克服这些固有的惩罚。

无论我们是从片段设计分子还是筛选现有库,我们都需要一种衡量成功的方法。想象一下,使用对接程序筛选一个包含一百万个化合物的库,其中只有100个是真正有活性的。我们如何知道我们的方法是否有效?我们不能等到合成所有一百万个化合物。相反,我们使用像富集因子(EF)这样的指标。如果我们查看计算排序列表的前1%(10,000个化合物),并发现我们捕获了100个活性物中的80个,那么我们的方法表现得非常出色。EF通过比较我们顶部部分中的命中率与随机机会预期的命中率来量化这一点。例如,一个16的EF值意味着我们的方法比猜测好16倍。这样的结果提供了强有力的验证,也许可以证明像使用多个受体构象来模拟柔性这样复杂方法的高计算成本是合理的。

理论与实验的对话

也许计算化学生物学中最激动人心的前沿是其与实验技术之间深刻且不断发展的伙伴关系。这就是整合建模的世界,在这里计算和实验不再在各自的领域中运作,而是进行着持续、富有成效的对话。

这种对话发生的最直接方式之一是通过使用实验约束。一个仅由理论能量函数引导的分子模拟,可能会漫游到物理上不现实的构象空间区域。实验数据,例如通过核磁共振(NMR)测量的两个质子之间的距离,可以充当引导绳。我们可以将这些数据作为惩罚项直接纳入我们的能量函数中。该势可以设计成具有“平底”,意味着只要距离保持在实验观察到的范围内,模拟就不会受到惩罚,但如果偏离这个容差范围,就会受到二次能量惩罚。在模拟退火等模拟方法中,算法在不断地做决策。一个违反约束的提议移动被接受的可能性较小,但并非不可能。这种概率性接受允许系统逃离局部最小值,同时仍然受到实验数据这只“磁手”的温和引导。

这种整合将我们引向最终、最微妙,也许也是最深刻的应用。许多关键的生物分子不是静态、刚性的实体。它们是动态、柔性的机器,以一个巨大的、相互转换的构象系综形式存在。对这类分子的溶液进行的实验,无论是X射线小角散射(SAXS)还是核磁共振(NMR),测量的都不是单个结构,而是整个群体的平均属性。这里存在一个关键的区别:可观测量的平均值与平均结构的可观测量不同。数学上,这是因为将结构与可观测量联系起来的函数是非线性的:⟨O(x)⟩≠O(⟨x⟩)\langle O(\mathbf{x}) \rangle \neq O(\langle \mathbf{x} \rangle)⟨O(x)⟩=O(⟨x⟩)。

例如,NMR交叉峰(NOE)的强度对两个质子之间的距离 rrr 极其敏感,其标度关系为 ⟨r−6⟩\langle r^{-6} \rangle⟨r−6⟩。由于这种反六次方依赖性,平均值在很大程度上由最紧密接触的构象主导,即使这些构象是瞬时的且稀疏分布的。一个基于单一平均结构的模型会完全忽略这一点,无法重现数据。同样,SAXS曲线取决于原子间距离的正弦,这是另一个非线性关系。为了正确地建模这些系统,我们的计算方法必须改变其目标。我们不应寻求单一的最佳结构,而必须生成一个构象系综,当其性质被适当地平均后,能够共同重现实验数据。这迫使我们拥抱分子世界的动态、统计现实,远远超越简单的“锁与钥”范式,进入生命本身丰富而复杂的舞蹈之中。