Preparing the non-periodic model

Note

The scripts for this step of the tutorial can be found in p450_enzyme/2_qmmm.

In this stage we prepare a model suitable for QM/MM calculations based on an MD snapshot generated from the solvation workflow in the previous stage.

Cutting a solvent shell

Unlike the initial MD runs in the solvation workflow, the QM/MM calculation is non-periodic and does not require a cubic water box. We can reduce the system size and therefore also the computational expense of the calculation by cutting out a solvent shell around the protein.

The script to cut the solvent shell is cut_solvent_shell.py. Firstly, we read in a pre-prepared snapshot in PQR format (to include charges):

# Read in pre-prepared snapshot PQR file (including charges)
p450 = Fragment(coords='p450_full_snapshot_sample.pqr')

Alternatively, if you have successfully run the long solvation workflow example in the previous stage, you can read in a snapshot from that instead, e.g.:

p450 = Fragment(coords='../1_solvation/_snapshots/_snapshot_19.pqr')

Note

A PDB file may also be used instead of a PQR file, if you also provide a list of charges, i.e.:

p450 = Fragment(coords='snapshot.pdb', charges=[1.00,...])

Once the Fragment has been read in, we identify a list of non-solvent atoms using an atom selection command:

# Identify non-solvent atoms to cut shell around
indices_non_solvent = p450.selectByCommonNames('non-solvent')

Next, we can cut the shell itself, using a special selection command defined for this purpose:

# Cut the shell
indices_p450_w_water_shell = p450.selectByShell(around=indices_non_solvent,
      padding=7.0, unit='angstrom', boundary='incl', convex=True)

The selectByShell method is controlled by the following options:

  • around specifies the system to cut around, defined by the previous atom selection.

  • padding defines the size of the shell to cut, in units defined by unit.

  • boundary - if 'incl' (default) will keep whole residues in the shell, if any part of them is within the padding distance. Alternatively, excl will remove whole residues unless all of the atoms are within the padding distance.

  • convex defines the selection algorithm. If True (default), the convex hull algorithm will be used which is the most efficient option, but results in a slightly different (smoother) shell than an explicit distance criterion (convex=False).

Once the shell selection has been made, a new Fragment object is created using the selection:

# Create a new Fragment object from the shell
p450_w_water_shell = p450.getSelected(indices_p450_w_water_shell)

Finally, the new fragment needs to be saved to disk for use in the following QM/MM calculation:

# Save the new Fragment for use in the QM/MM calculation
p450_w_water_shell.save('p450_w_water_shell.pqr')
p450_w_water_shell.save('p450_w_water_shell.pdb')

Before moving on to the next step, open p450_w_water_shell.pdb in your visualiser and see how it differs from the system used in the solvation workflow.