Part 2: Simulating Kinases with CHARMM and Amber. Retaining Water.

3 min readNov 4, 2023

Welcome back! From Part 1, we prepared our kinase by:

Adding Hydrogens to Carbon and Hetereo-atoms
Minimizing the hydrogens so they are in the most favorable state.

One important thing to note is when the PDB file is provided to you, you might see water trapped within the structure.

These waters are part of the protein structure and perhaps important to it’s shape. We want to retain these waters and not remove them as a solvent which is a key point.

For this part, we are going to take the waters, add hydrogens and minimize them in charmm. I grabbed all the water atoms out and placed them in their own separate PDB file available in the repository under step 2.

Cheminformatics-Teaching-Material/molecular_dynamics/simulating_kinases/step_2 at main ·…

Resources, Code, and Other things I use to teach Cheminformatics. …

github.com

Let’s go ahead and setup our variables and the water segid is W00A which can be found in the PDB File:

ATOM   2692  OH2 TIP3 A 507     25.464  -5.700  -9.555  1.00 28.74      W00A
ATOM   2693  OH2 TIP3 A 508     -4.363  12.732  -3.532  1.00 26.53      W00A
ATOM   2694  OH2 TIP3 A 509     17.162 -11.054 -26.296  1.00 14.34      W00A
ATOM   2695  OH2 TIP3 A 510      0.915 -16.168 -28.743  1.00 15.09      W00A
ATOM   2696  OH2 TIP3 A 511     18.802 -21.089 -16.846  1.00 35.84      W00A
ATOM   2697  OH2 TIP3 A 512     29.821  -7.016 -29.030  1.00 28.65      W00A
ATOM   2698  OH2 TIP3 A 513     12.611 -18.802 -28.010  1.00 29.04      W00A
ATOM   2699  OH2 TIP3 A 514    -13.959   5.032   3.437  1.00 33.48      W00A
ATOM   2700  OH2 TIP3 A 515     15.195 -19.143 -24.569  1.00 30.33      W00A
ATOM   2701  OH2 TIP3 A 516      9.107 -21.531 -21.490  1.00 19.31      W00A
ATOM   2702  OH2 TIP3 A 517      5.754   3.176 -12.816  1.00 21.11      W00A

Read it in in CHARMM:

! Inputs and Variables
!-------------------------
set in  =  water_xray
set segi = W00A
set out =  results/step2_h_xwat

Next we want to pass in the topology (the atom types) and parameter (bonding and non-bonding forces) files and read them into CHARMM:

set toppar = ../toppar 
set top = @TOPPAR/top_all22_prot.rtf 
set par = @TOPPAR/par_all22_prot.prm
set watstr = @TOPPAR/toppar_water_ions.str

! Read Topology and Parameter Files 
!----------------------------------
open read card unit 10 name @top 
read rtf card unit 10 
close unit 10 

open read card unit 10 name @par 
read param card flex unit 10
close unit 10

In FORTRAN, you assign a file a number or pointer to give an address when opening or closing the file. This is gives the rest of the code a reference number to look back too. In more modern languages this gets taken care of under the hood. Since we open the file and then close it after we can use the same number 10.

Some numbers are reserved so you might want to check the CHARMM documentation.

Read in the water stream file which are parametes for TIP3 (a method in how we model water) PDB and generate the molecules in CHARMM:

! Read Water Input
!------------------
stream @watstr


! Read Sequence Information
!---------------------------
open unit 10 read form name @{in}.pdb 
read sequ pdb unit 10 
generate @{segi} setup warn first none last none noangle nodihedral 
close unit 10

Read the coordinates and have a look at how CHARMM places the coordinates:

! Read Coordinate Information 
!-----------------------------
open read card unit 10 name @{in}.pdb 
read coor pdb unit 10 resid 
close unit 10 

! Print Heavy Atoms with Unknown Coordinates 
!--------------------------------------------
coor print sele ( .not. INIT ) .and. ( .not. hydrogen ) end

Build the hydrogens:


HBuild

And then write your output files:

!----------------
! Write Files 
!-----------------
open unit 10 write card name @{out}.pdb 
write coor pdb unit 10 
close unit 10 

open unit 10 write card name @{out}.psf 
write psf card unit 10 
close unit 10 

open unit 10 write card name @{out}.crd 
write coor card unit 10 
close unit 10 

stop