Remove Water

The second stage is to mark the waters that overlap the protein and lipid. For this we use VMD's atom selection capabilities. We are trying to build a system with hexagonal boundary conditions, so we chop three times, rotating the system 120 degrees each time to chop a different set of atoms. We also mark waters that are inside the hydrophobic core of the membrane, waters that are too far from the membrane, and waters that overlap the protein or membrane.

Once we hav marked the waters that we don't want in the final system, we run psfgen again, deleting the unwanted waters from the structure. We will use the resulting pdb and psf files for all of our subsequent simulations.

To do:

- Try modifying the values in the script to mark a different range of waters. Note that if you increase the width or height of your system, you may have to change the size of the periodic cell in your simulations.

- Use VMD to see which waters you are deleting. Rerun the first half of the script by hand to mark waters again. In the Graphics menu, change the coloring method to Beta. Waters that will be removed should be colored differently.

Run this step!     Continue

remove_water.vmd

# STEP 6: Remove Water

mol delete all
mol load psf pope_gram_wat.psf pdb pope_gram_wat.pdb

set all [atomselect top all]
$all set beta 0

# tag waters outside of the periodic cell (30A hexagon)
for { set i 0 } { $i < 3 } { incr i } {
  set sel [atomselect top \
    "water and same residue as (not hydrogen and (x < -14.99 or x > 14.99))"]
  $sel set beta 1
  $all move [transaxis z 120]
}

# tag waters in the lipid tail region
set sel [atomselect top \
  "water and same residue as (not hydrogen and (z > -15 and z < 15))"]
$sel set beta 1

# tag waters outside of the periodic cell (64A)
set sel [atomselect top \
  "water and same residue as (not hydrogen and (z < -31.99 or z > 31.99))"]
$sel set beta 1

# Select the marked waters.  All of our atom selections have been
# "by residue", so we only need to select the oxygen atoms to get
# the water molecules that have been marked.
set badwater [atomselect top "name OH2 and beta > 0"]

# Read the structure and coordinates into psfgen
resetpsf
readpsf pope_gram_wat.psf
coordpdb pope_gram_wat.pdb

# Delete the marked waters from the structure
foreach segid [$badwater get segid] resid [$badwater get resid] {
  delatom $segid $resid
}

# Write out the correctly solvated structure
writepsf grama.psf
writepdb grama.pdb

# Read it back into VMD for display
mol delete all
mol load psf grama.psf pdb grama.pdb