2

我有一些使用 MDAnalysis 的工作代码,它将质量时间序列的残差中心保存在一个数组中,但我想知道是否有更 Pythonic 或整体高效/快速的方式(理解、数组操作......)来做到这一点。

import MDAnalysis as mda
mdau = mda.Universe(pdb, xtc)
arr = np.empty(( len(mdau.select_atoms("protein").residues), len(mdau.trajectory), 3 ))
# 288 protein residues, 1250 frames and 3 xyz-coordinates per center of mass; this array shape is important
for ts in mdau.trajectory:
    for num, res in enumerate(mdau.select_atoms("protein").residues):
        arr[num, ts.frame] = res.atoms.center_of_mass()

我正在使用的 .pdb 和 .xtc 文件可以在以下链接中下载: https ://submission.gpcrmd.org/dynadb/files/Dynamics/11579_dyn_169.pdb https://submission.gpcrmd.org/dynadb/files/动态/11576_trj_169.xtc

4

1 回答 1

3

您可以对代码进行一些更改:

  • 在 for 循环之外选择您的蛋白质原子,而不是在每次迭代时
  • 通过将compound='residues'参数传递给center_of_mass方法,对残差的质心计算进行矢量化
  • 使用ag.n_residuesandu.trajectory.n_frames属性

以下是使用这些建议的代码更新:

import numpy as np
import MDAnalysis as mda

u = mda.Universe('11579_dyn_169.pdb', '11576_trj_169.xtc')
protein = u.select_atoms("protein")
arr = np.empty((protein.n_residues, u.trajectory.n_frames, 3))

for ts in u.trajectory:
    arr[:, ts.frame] = protein.center_of_mass(compound='residues')

于 2022-01-19T14:57:58.700 回答