PDBReporter options

[Saurabh Belsare]

Hi,

I'm simulating a solvated protein. Is there any way to configure PDBReporter to write the trajectory only for the protein atoms, and not write out the waters?

Thank you.

[Jason Swails]

The only way you can do this is if you write your own PDB reporter class that does this.

Alternatively, if it doesn't have to be a PDB file that you generate (but rather, HDF5, NetCDF, or DCD files will suffice), then you can use the OpenMM reporters that are included with MDTraj instead. These offer an "atomSubset" argument that does exactly this.

[Mark Williamson]

This is snippet example of what Jason is referring to:

Code: Select all

import mdtraj.reporters

....<snip>...

# Get the index of protein atoms only, 
# to enable omission of water and counter ions in the production trajectory

traj = mdtraj.load('minimisation.pdb')
top, bonds = traj.top.to_dataframe()

protein_indices = np.where(top.resName[(top.resName != "HOH" ) & (top.resName != "CL")])[0]

....<snip>...

simulation.reporters.append(mdtraj.reporters.HDF5Reporter('production.h5', 5000, atomSubset=protein_indices))

[Saurabh Belsare]

Thanks for the code!

I could save these binary formats. However, I want to further analyze the trajectories that are generated using my own scripts. Is there any of these binary formats that could then easily be converted to PDB for processing after the simulation is completed?

Thank you.

[Jason Swails]

All of them . The "mdconvert" program that is part of MDTraj will do this for you. Unless you use the HDF5 trajectory format, though, I think you will need to create a stripped-down topology (because no program will be able to generate a PDB file from a DCD file containing frames with 100 atoms and a topology file containing 1000 atoms). So the HDF5 reporter is probably your best bet here.

Upsides of using a binary reporter compared to a PDB reporter are twofold:

1. The files themselves are a lot smaller and can be written a lot faster, meaning that your simulations will run faster

2. The unit cell dimensions can be saved every step instead of just once at the beginning

[Saurabh Belsare]

Thanks! This is definitely the way I should be going, instead of writing the whole system and then removing the waters. However, when I'm trying to set up the simulation as outlined in http://mdtraj.org/latest/examples/openmm.html, I am unable to use PME for the nonbondedMethod. It gives the following error:

ValueError: Requested periodic boundary conditions for a Topology that does not specify periodic box dimensions

I was able to use PME when I loaded the same pdb structure using PDBFile(). Is there some other way I need to load the periodic box dimensions when I use the mdtraj.load() to load a pdb file?

Thanks.