我试图使用TPRParser
来解析从生成的TPR
文件GROMACs
。不幸的是,它抛出了一个我从未见过的错误:
Traceback (most recent call last):
line 43, in <module>
topo = mda.topology.TPRParser.TPRParser(TPR).parser(tpr_resid_from_one=True)
AttributeError: 'TPRParser' object has no attribute 'parser'
有人能帮我删除这个错误吗?
这是我完整的代码和数据文件可在这里
import argparse
# define and parse command line arguments
parser = argparse.ArgumentParser()
parser.add_argument('input', metavar='INPUT', help='XTC input file')
parser.add_argument('--output', metavar='FILENAME', help='H5MD output file for saving trajectory data')
parser.add_argument('-v', '--verbose', action='store_true', help='be verbose')
args = parser.parse_args()
import MDAnalysis as mda
import numpy as np
import matplotlib.pyplot as plt
import os, sys
from collections import OrderedDict
XTC = args.input
TPR = os.path.splitext(XTC)[0] + '.tpr'
if not os.path.isfile(TPR):
print("Error: TPR file {0} does not exists".format(TPR))
sys.exit(1)
# the following list was generated manually using
# grep "^ *1[A-Z]" -B1 lp400start.gro | grep -v "^ *1[A-Z]" | cut -c 1-5
protein_seq_length = [ 38, 38, 137, 137, 252, 252, 576, 576, 1092, 1092, 1016, 1016, 1285, 1285 ]
protein_resid_range = []
i = 0
for n in protein_seq_length:
protein_resid_range += [ slice(i, i + n) ]
i += n
print("Expecting {0:d} residues in {1:d} proteins".format(i, len(protein_resid_range)))
# build reduced universe to match XTC
# (ignore warnings that no coordinates are found for the TPR)
topo = mda.topology.TPRParser.TPRParser(TPR).parser(tpr_resid_from_one=True)
u0 = mda.Universe(topo)
u = mda.Merge(u0.select_atoms("not resname W and not resname WF and not resname ION"))
u.load_new(XTC)
# this selection does not work, it mixes the atoms of both EPHA proteins
# FIXME delete this output
print(u.select_atoms("segid seg_0_EPHA"))
# select proteins through range of residue IDs
for idx in protein_resid_range:
res = u.residues[idx]
print("Protein of type {0:s} is composed of {1:d} residues and {2:n} atoms:n{3}n").format(res[0].segment.segid, len(res), len(res.atoms), res)
# segments
protein_segments = u.segments
# build the fragments
# (a dictionary with the key as the protein name -- I'm using an
# OrderedDict so that the order is the same as in the TPR)
protein_fragments = OrderedDict((seg.segid[1:], seg.atoms.fragments) for seq in protein_segments)
# this doesn't seem to be helpful either FIXME delete this output
for f in protein_fragments:
print(f)
# analyze trajectory
timeseries = []
for ts in u.trajectory[1:4]:
coms = []
for name, proteins in protein_fragments.items():
# loop over all the different proteins
# unwrap to get the true COM under PDB (double check!!)
coms.extend([p.center_of_mass(unwrap=True) for p in proteins])
timeseries.append(coms)
timeseries = np.array(timeseries, dtype=float)
if args.output:
with opne(args.output, 'w') as outfile:
# make clear distinction in frames
outfile.write('# Array shape: {0}n'.format(timeseries.shape))
outfile.write('# New frame: {0}n'.format(u.trajectory.time))
# Iterating through a n-dim array
for data in timeseries:
np.savetxt(outfile, data, fmt='%-7.2f')
您不需要显式地使用TPRParser。而不是
### do NOT use this code
topo = mda.topology.TPRParser.TPRParser(TPR).parser(tpr_resid_from_one=True)
u0 = mda.Universe(topo)
使用
u0 = mda.Universe(TPR, tpr_resid_from_one=True)
Universe()
类自动使用正确的解析器(基于文件扩展名)。几乎没有必要显式地使用解析器。请使用Universe()
。