Introduction to MDTraj

Start by loading up a trajectory from disk. MDTraj will automatically parse the file extension and use the appropriate loader.

>>> import mdtraj as md
>>> t = md.load('trajectory.xtc', top='trajectory.pdb')
>>> print t
<mdtraj.Trajectory with 100 frames, 22 atoms at 0x109f0a3d0>

To load files that don’t contain topology information, like Gromacs XTC files, we need to supply something with the top keyword argument that describes the topology, for example a PDB file.

If I’m interested in only a subset of the frames of the trajectory, I can slice it

>>> # lets take a look at the first ten frames
>>> print t[1:10]
<mdtraj.Trajectory with 9 frames, 22 atoms at 0x109b91310>

>>> # or maybe the last frame?
>>> print t[-1]
<mdtraj.Trajectory with 1 frames, 22 atoms at 0x109b97810>

There’s a lot of information in the trajectory object. The most obvious is the Cartesian coordinates. They’re stored as a numpy array under xyz. All of the distances in the Trajectory are stored in nanometers. The time unit is picoseconds. Angles are stored in degrees (not radians).

>>> print t.xyz.shape
(100, 22, 3)
>>> print np.mean(t.xyz)
0.89365752249053032

>>> # the simulation time (in picoseconds) of th first 10 frames
>>> print t.time[0:10]
array([ 0.002,  0.004,  0.006,  0.008,  0.01 ,  0.012,  0.014,  0.016,
      0.018,  0.02 ], dtype=float32)

>>> # or the unitcell lengths in the last frame? (in nanometers of course)
>>> t.unitcell_lengths[-1]
array([ 2.,  2.,  2.], dtype=float32)

Saving the trajectory back to disk is easy.

>>> # the hdf5 format stores the topology inside the file for convenience
>>> t[::2].save('halftraj.h5')

>>> # the format will be parsed based on the extension, or you can call the
>>> # format-specific save methods
>>> t[0:10].save_dcd('first-ten-frames.dcd')

The trajectory contains a reference to a topology object, which can come in handy. For example, if you want to save a copy of your trajectory with only alpha carbons present, you can do that pretty easily.

>>> atoms_to_keep = [a.index for a in t.topology.atoms if a.name == 'CA']
>>> t.restrict_atoms(atoms_to_keep)  # this acts inplace on the trajectory
>>> t.save('CA-only.h5')