4.2.2. 本地联系人分析 MDAnalysis.analysis.contacts

此模块包含用于分析本机联系人的类 Q 在弹道上。构象的天然接触是存在于参考结构和构象中的接触。参考结构中的触点始终定义为比距离近 radius 。构象的自然接触比例可以用不同的方法来计算。此模块支持下面列出的3种不同指标，以及自定义指标。

  1. 硬切 ：把原子算作接触 i 和 j 必须至少与参考结构中的位置一样接近。

  2. 软切削 ：原子对 i 和 j 如果距离为0，则指定的软势为1；如果距离与参考中的距离相同，则指定为1/2；如果距离较大，则指定为0。有关势和参数的准确定义，请参阅函数 soft_cut_q() 。

  3. 半径切割 ：把原子算作接触 i 和 j 不能比某一距离更远 radius 。

“母语接触者的比例” Q(t) 是介于0和1之间的数字，其计算方法为给定时间范围内的本地接触总数除以参考结构中的接触总数。

4.2.2.1. 接触分析示例

4.2.2.1.1. 一维接触分析

作为一个例子，我们分析了当ADK酶打开时盐桥的打开(“解压缩”)；这是MDAnalysis中的一个示例轨迹。**

import numpy as np
import matplotlib.pyplot as plt
import MDAnalysis as mda
from MDAnalysis.analysis import contacts
from MDAnalysis.tests.datafiles import PSF,DCD
# example trajectory (transition of AdK from closed to open)
u = mda.Universe(PSF,DCD)
# crude definition of salt bridges as contacts between NH/NZ in ARG/LYS and
# OE*/OD* in ASP/GLU. You might want to think a little bit harder about the
# problem before using this for real work.
sel_basic = "(resname ARG LYS) and (name NH* NZ)"
sel_acidic = "(resname ASP GLU) and (name OE* OD*)"
# reference groups (first frame of the trajectory, but you could also use a
# separate PDB, eg crystal structure)
acidic = u.select_atoms(sel_acidic)
basic = u.select_atoms(sel_basic)
# set up analysis of native contacts ("salt bridges"); salt bridges have a
# distance <6 A
ca1 = contacts.Contacts(u, select=(sel_acidic, sel_basic),
                        refgroup=(acidic, basic), radius=6.0)
# iterate through trajectory and perform analysis of "native contacts" Q
ca1.run()
# print number of averave contacts
average_contacts = np.mean(ca1.results.timeseries[:, 1])
print('average contacts = {}'.format(average_contacts))
# plot time series q(t)
fig, ax = plt.subplots()
ax.plot(ca1.results.timeseries[:, 0], ca1.results.timeseries[:, 1])
ax.set(xlabel='frame', ylabel='fraction of native contacts',
       title='Native Contacts, average = {:.2f}'.format(average_contacts))
fig.show()

第一张图显示，当ADK打开时，当酶打开时，大约20%存在于关闭状态的盐桥消失。他们以一种循序渐进的方式打开(电影更清楚地说明了这一点 AdK_zipper_cartoon.avi )。

注意事项

针对不同模拟建议的截止距离

  • 对于全原子模拟，截止值=4.5？

  • 对于粗粒度模拟，截止值为6.0？

4.2.2.1.2. 二维接触分析(Q1-Q2)

分析ADK在闭合构象和开放构象之间的单尺度转变，并绘制出在Q1-Q2上投影的轨迹 [Franklin2007] **

import MDAnalysis as mda
from MDAnalysis.analysis import contacts
from MDAnalysisTests.datafiles import PSF, DCD
u = mda.Universe(PSF, DCD)
q1q2 = contacts.q1q2(u, 'name CA', radius=8)
q1q2.run()

f, ax = plt.subplots(1, 2, figsize=plt.figaspect(0.5))
ax[0].plot(q1q2.results.timeseries[:, 0], q1q2.results.timeseries[:, 1],
           label='q1')
ax[0].plot(q1q2.results.timeseries[:, 0], q1q2.results.timeseries[:, 2],
           label='q2')
ax[0].legend(loc='best')
ax[1].plot(q1q2.results.timeseries[:, 1],
           q1q2.results.timeseries[:, 2], '.-')
f.show()

将生成的路径与 MinActionPath result for AdK [Franklin2007] 。

4.2.2.1.3. 编写您自己的联系人分析

这个 Contacts 类被设计为可扩展以用于您自己的分析。作为一个例子，我们将分析ADK的酸性基团和碱性基团何时彼此接触；这意味着在参考文献中形成的至少一个接触比2.5？更近。

为此，我们定义了一个新函数来确定是否有任何接触距离大于2.5；该函数必须实现 Contacts ：：

def is_any_closer(r, r0, dist=2.5):
    return np.any(r < dist)

前两个参数 r 和 r0 由以下公司提供 Contacts 当它呼唤时 is_any_closer() 而其他参数可以作为关键字args使用 kwargs 中的参数 Contacts 。

接下来，我们将创建 Contacts 类，并使用 is_any_closer() 函数作为参数 method 并运行分析：

# crude definition of salt bridges as contacts between NH/NZ in ARG/LYS and
# OE*/OD* in ASP/GLU. You might want to think a little bit harder about the
# problem before using this for real work.
sel_basic = "(resname ARG LYS) and (name NH* NZ)"
sel_acidic = "(resname ASP GLU) and (name OE* OD*)"

# reference groups (first frame of the trajectory, but you could also use a
# separate PDB, eg crystal structure)
acidic = u.select_atoms(sel_acidic)
basic = u.select_atoms(sel_basic)

nc = contacts.Contacts(u, select=(sel_acidic, sel_basic),
                       method=is_any_closer,
                       refgroup=(acidic, basic), kwargs={'dist': 2.5})
nc.run()

bound = nc.results.timeseries[:, 1]
frames = nc.results.timeseries[:, 0]

f, ax = plt.subplots()

ax.plot(frames, bound, '.')
ax.set(xlabel='frame', ylabel='is Bound',
       ylim=(-0.1, 1.1))

f.show()

4.2.2.2. 功能

MDAnalysis.analysis.contacts.hard_cut_q(r, cutoff)[源代码]

计算本地联系人的分数 Q 为了一次硬切断。

截止值可以是浮点型，也可以是 ndarray 与形状相同的 r 。

参数:

  • r (ndarray) -- 距离矩阵

  • cutoff (ndarray | float) -- 截断值以计算距离。可以是大小与距离相同的ndarray的浮点

返回:

Q --接触部分

返回类型:

float

MDAnalysis.analysis.contacts.soft_cut_q(r, r0, beta=5.0, lambda_constant=1.8)[源代码]

计算本地联系人的分数 Q 对于软性切断

原生联系函数定义为 [Best2013]

\[Q(r，r_0)=\frac{1}{1+e^{\beta(r-\lambda r_0)}}\]

不同模拟类型的合理值为

  • 所有原子 ： lambda_constant = 1.8 (无单位)

  • 粗粒化 ： lambda_constant = 1.5 (无单位)

参数:

  • r (array) -- 时间t处的接触距离

  • r0 (array) -- 时间t=0时的接触距离，参考距离

  • beta (float (default 5.0 Angstrom)) -- 开关功能的柔软性

  • lambda_constant (float (default 1.8, unitless)) -- 参考距离公差

返回:

Q --本地联系人的比例

返回类型:

float

引用

[Best2013] (1,2,3)

Robert B. Best, Gerhard Hummer, and William A. Eaton. Native contacts determine protein folding mechanisms in atomistic simulations. Proceedings of the National Academy of Sciences, 110(44):17874–17879, 2013. doi:10.1073/pnas.1311599110.

MDAnalysis.analysis.contacts.radius_cut_q(r, r0, radius)[源代码]

计算本地联系人 Q 基于单个距离半径。

参数:

  • r (ndarray) -- 原子间距离阵列

  • r0 (ndarray) -- 未用来填满 Contacts API

  • radius (float) -- 原子间形成接触的距离

返回:

Q --接触部分

返回类型:

float

引用

[Franklin2007] (1,2,3,4,5,6,7,8,9)

Joel Franklin, Patrice Koehl, Sebastian Doniach, and Marc Delarue. MinActionPath: maximum likelihood trajectory for large-scale structural transitions in a coarse-grained locally harmonic energy landscape. Nucleic Acids Research, 35(suppl_2):W477–W482, 2007. doi:10.1093/nar/gkm342.