首页 > 其他分享 >MindSponge分子动力学模拟——增强采样(2024.11)

MindSponge分子动力学模拟——增强采样(2024.11)

时间:2024-11-01 15:09:13浏览次数:3  
标签:采样 phi 2024.11 0.00 psi MindSponge ATOM np import

技术背景

关于增强采样(Enhanced Sampling)算法的具体原理,这里暂不做具体介绍,感兴趣的童鞋可以直接参考下这篇综述文章:Enhanced sampling in molecular dynamics。大致的作用就是,通过统计力学的方法,使得目标分子的CV(Collective Variables)具有一个尽可能大的采样子空间,并且可以将其还原回真实的自由能面。常用的增强采样算法,有早期的伞形采样,到后来大家常用的MetaDynamics以及高老师的温度积分增强采样算法(ITS)。在MindSponge中已经实现了MetaDynamics算法和ITS算法,本文我们使用MetaDynamics算法来做一个演示。

准备工作

使用MindSponge,可以参考本系列文章先了解一下MindSponge的安装和基本使用方法。这里我们用一个比较简单的多肽体系进行一个测试,相应的pdb文件为:

CRYST1    0.000    0.000    0.000  90.00  90.00  90.00 P 1           1
ATOM      1  H1  ACE A   1      -1.838  -6.570  -0.492  0.00  0.00      
     
ATOM      2  CH3 ACE A   1      -0.764  -6.587  -0.283  0.00  0.00      
     
ATOM      3  H2  ACE A   1      -0.392  -7.533  -0.746  0.00  0.00      
     
ATOM      4  H3  ACE A   1      -0.592  -6.446   0.740  0.00  0.00      
     
ATOM      5  C   ACE A   1      -0.006  -5.404  -0.828  0.00  0.00      
     
ATOM      6  O   ACE A   1      -0.544  -4.619  -1.673  0.00  0.00      
     
ATOM      7  N   ALA A   2       1.278  -5.323  -0.423  0.00  0.00      
     
ATOM      8  H   ALA A   2       1.622  -5.845   0.368  0.00  0.00      
     
ATOM      9  CA  ALA A   2       2.284  -4.164  -0.399  0.00  0.00      
     
ATOM     10  HA  ALA A   2       2.098  -3.653   0.505  0.00  0.00      
     
ATOM     11  CB  ALA A   2       3.651  -4.787  -0.566  0.00  0.00      
     
ATOM     12  HB1 ALA A   2       4.274  -4.031  -0.972  0.00  0.00      
     
ATOM     13  HB2 ALA A   2       3.977  -5.106   0.419  0.00  0.00      
     
ATOM     14  HB3 ALA A   2       3.697  -5.612  -1.274  0.00  0.00      
     
ATOM     15  C   ALA A   2       1.995  -3.152  -1.576  0.00  0.00      
     
ATOM     16  O   ALA A   2       1.544  -2.065  -1.221  0.00  0.00      
     
ATOM     17  N   NME A   3       2.255  -3.614  -2.845  0.00  0.00      
     
ATOM     18  H   NME A   3       2.788  -4.485  -2.929  0.00  0.00      
     
ATOM     19  CH3 NME A   3       1.991  -2.802  -4.055  0.00  0.00      
     
ATOM     20 HH31 NME A   3       2.561  -1.891  -3.988  0.00  0.00      
     
ATOM     21 HH32 NME A   3       1.897  -3.419  -4.937  0.00  0.00      
     
ATOM     22 HH33 NME A   3       0.985  -2.388  -3.930  0.00  0.00      

END

该构象可以保存成pdb文件然后直接用VMD进行可视化:

其他的力场文件我们直接使用MindSponge已经自带的amber.ff14sb力场即可。

普通MD案例

在测试增强采样算法之前,我们可以先跑一段普通的分子模拟测试一下:

from mindspore import nn, context
import mindspore as ms
# 固定随机种子,确保朗之万控温的随机数可以被复现
ms.set_seed(0)
# 这里使用的是一个未编译版本的mindsponge,所以要把sponge所在路径添加到系统路径中
import sys
sys.path.insert(0, '../..')
from sponge import ForceField, Sponge, set_global_units, Protein, UpdaterMD, WithEnergyCell
from sponge.callback import RunInfo, WriteH5MD
from sponge.colvar import Torsion
from sponge.function import PI

context.set_context(mode=context.GRAPH_MODE, device_target='GPU')
set_global_units('A', 'kcal/mol')

system = Protein('mol.pdb', template=['protein0.yaml'], rebuild_hydrogen=False)
energy = ForceField(system, parameters=['AMBER.FF14SB'], use_pme=False)
min_opt = nn.Adam(system.trainable_params(), 1e-03)

dihedrals = Torsion([[4, 6, 8, 14], [6, 8, 14, 16]])
run_info = RunInfo(20)

opt = UpdaterMD(system=system,
                time_step=1e-3,
                integrator='velocity_verlet',
                temperature=300,
                thermostat='langevin')

sim = WithEnergyCell(system, energy)
md = Sponge(sim, optimizer=opt, metrics={'dihedrals': dihedrals})
cb_h5md = WriteH5MD(system, 'test_meta.h5md', save_freq=10, write_image=False)
md.run(5000, callbacks=[run_info, cb_h5md])

运行之后会在本地保存一个h5md格式的轨迹文件,可以在vscode中使用h5web拓展工具打开:

找到dihedral的value,用matrix的形式查看,然后可以导出csv格式(如:path_sink.csv):

导出csv之后,可以使用如下的python脚本进行绘图:

import numpy as np
import matplotlib.pyplot as plt

def gaussian2(x1, x2, sigma1=1.0, sigma2=1.0, A=0.5):
    return np.sum(A*np.exp(-0.5*(x1**2/sigma1**2+x2**2/sigma2**2))/np.pi/sigma1/sigma2, axis=-1)

def potential_energy(position, psi, phi, sigma1, sigma2):
    # (A, )
    psi_, phi_ = position[:, 0], position[:, 1]
    # (A, R)
    delta_psi = psi_[:, None] - psi[None]
    delta_phi = phi_[:, None] - phi[None]
    # (A, )
    Z = -np.log(gaussian2(delta_psi, delta_phi, sigma1=sigma1, sigma2=sigma2, A=2.0)+1)
    return Z

data = np.genfromtxt('./path_sink.csv', delimiter=',')
phi = data[:, 0]
psi = data[:, 1]

num_grids = 100
num_levels = 10
psi_grids = np.linspace(-np.pi, np.pi, num_grids)
phi_grids = np.linspace(-np.pi, np.pi, num_grids)
grids = np.array(np.meshgrid(psi_grids, phi_grids)).T.reshape((-1, 2))

Z = potential_energy(grids, phi, psi, 1.0, 1.0).reshape((psi_grids.shape[0], phi_grids.shape[0])).T
X,Y = np.meshgrid(psi_grids, phi_grids)
levels = np.linspace(np.min(Z), np.max(Z), num_levels)

plt.figure()
plt.title("Biased MD Traj")
plt.xlabel(r'$\phi$')
plt.ylabel(r'$\psi$')
fc = plt.contourf(X, Y, Z, cmap='Greens', levels=levels)
plt.colorbar(fc)

plt.xlim(-np.pi, np.pi)
plt.ylim(-np.pi, np.pi)
plt.plot(phi, psi, 'o', alpha=0.4, color='red')
plt.savefig('meta.png')

画出来的轨迹效果如下图所示:

可以发现,在不加任何的Bias的时候,整个轨迹还是比较集中的在一个区域,一般该区域就是对应了一个能量极小值点附近的区域。

MetaDynamics案例

从上一个章节的结果中可以看出,常规的分子模拟采样方法很容易陷入到一个局部区域中,这样使得我们很难可以观测到其他采样子空间所对应的构象和相关信息,因此我们可以引入一个MetaDynamics增强采样算法,做一个有偏估计:

from mindspore import nn, context
import mindspore as ms
ms.set_seed(0)
import sys
sys.path.insert(0, '../..')
from sponge import ForceField, Sponge, set_global_units, Protein, UpdaterMD, WithEnergyCell
from sponge.callback import RunInfo, WriteH5MD
from sponge.colvar import Torsion
from sponge.function import PI
from sponge.sampling import Metadynamics

context.set_context(mode=context.GRAPH_MODE, device_target='GPU')
set_global_units('A', 'kcal/mol')

system = Protein('mol.pdb', template=['protein0.yaml'], rebuild_hydrogen=False)
energy = ForceField(system, parameters=['AMBER.FF14SB'], use_pme=False)
min_opt = nn.Adam(system.trainable_params(), 1e-03)

# 这里定义的CV是phi和psi角,是一个2维的CV
dihedrals = Torsion([[4, 6, 8, 14], [6, 8, 14, 16]])
run_info = RunInfo(20)

# 配置Meta的参数,主要是高斯波包的高度、宽度、更新频率、CV范围、CV格点数等等
metad = Metadynamics(
    colvar=dihedrals,
    update_pace=10,
    height=2.5,
    sigma=0.4,
    grid_min=-PI,
    grid_max=PI,
    grid_bin=50,
    temperature=300,
    bias_factor=100,
    use_cutoff=True,
)

opt = UpdaterMD(system=system,
                time_step=1e-3,
                integrator='velocity_verlet',
                temperature=300,
                thermostat='langevin')

sim = WithEnergyCell(system, energy, bias=metad)
md = Sponge(sim, optimizer=opt, metrics={'dihedrals': dihedrals})
cb_h5md = WriteH5MD(system, 'test_meta.h5md', save_freq=10, write_image=False)
md.run(5000, callbacks=[run_info, cb_h5md])

用类似的方法,可以计算得增强采样之后的CV轨迹:

可以看到,采样子空间在Meta的作用下已经扩展到几乎整个CV空间,这使得我们可以更快的去分析整个采样空间各处的自由能的相对大小。

总结概要

本文介绍了在MindSponge中进行分子动力学模拟以及增强采样的实现方法。通过使用MetaDynamics增强采样算法,我们可以将分子模拟的采样子空间,从某个能量极小值区域,扩大到尽可能大的采样子空间。

版权声明

本文首发链接为:https://www.cnblogs.com/dechinphy/p/mindsponge-meta.html

作者ID:DechinPhy

更多原著文章:https://www.cnblogs.com/dechinphy/

请博主喝咖啡:https://www.cnblogs.com/dechinphy/gallery/image/379634.html

标签:采样,phi,2024.11,0.00,psi,MindSponge,ATOM,np,import
From: https://www.cnblogs.com/dechinphy/p/18518508/mindsponge-meta

相关文章

  • 电流采样电阻
    AP5165和电感(L)、电流采样电阻(RS)形成一个自振荡的连续电感电流模式的降压型恒流LED控制器。VIN上电时,电感(L)和电流采样电阻(RS)的初始电流为零,LED输出电流也为零。这时候,CS比较器的输出为高,内部功率开关导通,SW的电位为低。电流通过电感(L)、电流采样电阻(RS)、LED和内部功率开关......
  • 小波下采样,即插即用
    论文:Haarwaveletdownsampling:AsimplebuteffectivedownsamplingmoduleGitHub地址:https://github.com/apple1986/HWD论文地址:                    https://www.sciencedirect.com/science/article/pii/S0031320323005174 ......
  • 蓝桥杯备赛&学习计划 2024.11—2025.4
    学习计划概览2024年11月到12月-巩固基础,学习基本算法。2025年1月到2月-学习中级算法和数据结构。2025年3月-进阶算法学习和刷题练习。2025年4月-复习重点知识,专注于比赛准备。详细周计划2024年11月:基础知识&基础算法第1-2周:复习基本控制结构(循环、条件语......
  • 点云学习笔记2——使用VoxelGrid滤波器进行点云降采样(c++)
    #include<iostream>#include<pcl/point_cloud.h>#include<pcl/io/pcd_io.h>#include<pcl/point_types.h>#include<pcl/filters/voxel_grid.h>#include<pcl/common/common_headers.h>#include<pcl/io/pcd_io.h>#inclu......
  • 点云学习笔记4——点云滤波降采样后进行4PCS粗配准【四点一致集配准算法(4-Point Congr
    #include<iostream>#include<pcl/point_cloud.h>#include<pcl/point_types.h>#include<pcl/filters/voxel_grid.h>#include<pcl/common/common_headers.h>#include<pcl/io/pcd_io.h>#include<pcl/visualization/cloud_vi......
  • 重采样方法(交叉验证法)——基于glm与LOOCV法(Weekly数据集分析)
    Chapter5:Exercise7读取数据集Weekly数据集通常指的是在统计、数据分析或机器学习领域中,一个以周为单位进行记录的数据集合。以下是对Weekly数据集的一个详细介绍:一、数据来源与背景Weekly数据集可能来源于多个领域,如金融、经济、市场营销等,这些领域通常需要按周跟踪......
  • 提高ADC采样精度:C语言中的滤波与取平均值技巧
    在嵌入式系统中,ADC(模数转换器)是常用的组件,用于将模拟信号转换为数字信号。然而,由于噪声和其他干扰因素,ADC采样值可能会波动,导致读数不稳定。为了提高ADC读数的准确性,常用的方法是进行滤波和取平均值。本文将详细介绍如何在C语言中实现ADC采样值的滤波和取平均值,并提供详细的代......
  • 2024.10.27~2024.11.3
    2024.10.27这么说吧,csp-s打的不好,是时候做出些调整了约法n章:1.在NOIP之前把ybt刷完,保守估计一天5道题2.一道题若超出一个半小时内没有A就换下一道题,并在博客中记录此题并整理思路,有时间补完3.模拟赛我的得分要有以下两种评估:切题得分和难题高分暴力得分4.禁用一个月B站,休息......
  • 实验干货|电流型霍尔传感器采样设计03-信号调理
    在前两篇博客中,将霍尔输出的电流信号转换成了有正有负的电压信号,但是DSP需要采集0~3V的电压信号,因此需要对信号缩放并抬升至全部为正的信号。常见的方法是,通过比例放大(缩小)电路对信号进行放缩,通过加法电路抬升基准电平。这里分为两步,首先设计基准电平。设计基准电平DSP的A......
  • 过采样与欠采样技术原理图解:基于二维数据的常见方法效果对比
    在现实场景中,收集一个每个类别样本数量完全相同的数据集是十分困难的。实际数据往往是不平衡的,这对于分类模型的训练可能会造成问题。当模型在这样一个不平衡数据集上训练时,由于某个类别的样本数量远多于其他类别,模型通常会更擅长预测样本量较大的类别,而在预测小类别时表现不......