Optimize Conformations

In this part we will optimize the geometries of conformations generated in the previous step using NAMD.

Configuration

Let’s find the location of NAMD executable:

In [1]: from prody.utilities import which

In [2]: namd2 = which('namd2')

In [3]: namd2
Out[3]: '/usr/local/bin/namd2'

We will need a force field file for energy minimization. VMD ships with CHARMM force field files. We can write a tcl script to find and write their location as follows:

In [4]: tcl_cmd = '''package require readcharmmpar
   ...: package require readcharmmtop
   ...: global env
   ...: set outfile [open charmmdir.txt w]
   ...: puts $outfile $env(CHARMMPARDIR)
   ...: puts $outfile $env(CHARMMTOPDIR)
   ...: close $outfile
   ...: exit'''
   ...: 

In [5]: with open('where_is_charmmpar.tcl', 'w') as inp:
   ...:    inp.write(tcl_cmd)
   ...: 

This can be run in vmd from ipython as below:

In [6]: !vmd -dispdev text -e where_is_charmmpar.tcl
/usr/local/lib/vmd/vmd_LINUXAMD64: /usr/lib/x86_64-linux-gnu/mesa/libGL.so.1: no version information available (required by /usr/local/lib/vmd/vmd_LINUXAMD64)
Info) VMD for LINUXAMD64, version 1.9.3 (November 30, 2016)
Info) http://www.ks.uiuc.edu/Research/vmd/                         
Info) Email questions and bug reports to vmd@ks.uiuc.edu           
Info) Please include this reference in published work using VMD:   
Info)    Humphrey, W., Dalke, A. and Schulten, K., `VMD - Visual   
Info)    Molecular Dynamics', J. Molec. Graphics 1996, 14.1, 33-38.
Info) -------------------------------------------------------------
Info) Multithreading available, 16 CPUs detected.
Info)   CPU features: SSE2 AVX 
Info) Free system memory: 58GB (91%)
Info) No CUDA accelerator devices available.
Info) Dynamically loaded 2 plugins in directory:

Info) /usr/local/lib/vmd/plugins/LINUXAMD64/molfile

1.3

1.2

file4

Info) VMD for LINUXAMD64, version 1.9.3 (November 30, 2016)
Info) Exiting normally.
vmd > 

We then read the output file to get the parameter directory:

In [7]: inp = open('charmmdir.txt', 'r')

In [8]: lines = inp.readlines()

In [9]: inp.close()

In [10]: import os

In [11]: par = os.path.join(lines[0].strip(), 'par_all27_prot_lipid_na.inp')

In [12]: top = os.path.join(lines[1].strip(), 'top_all27_prot_lipid_na.inp')

In [13]: par
Out[13]: '/usr/local/lib/vmd/plugins/noarch/tcl/readcharmmpar1.3/par_all27_prot_lipid_na.inp'

In [14]: top
Out[14]: '/usr/local/lib/vmd/plugins/noarch/tcl/readcharmmtop1.2/top_all27_prot_lipid_na.inp'

To configure this computer

Let’s make a folder for writing optimization input and output files:

In [15]: mkdir -p p38_optimize

We will write an NAMD configuration file for each conformation based on min.conf:

In [16]: import glob

In [17]: conf = open('conformational_sampling_files/min.conf').read()

In [18]: for pdb in glob.glob(os.path.join('p38_ensemble', '*.pdb')):
   ....:     fn = os.path.splitext(os.path.split(pdb)[1])[0]
   ....:     pdb = os.path.join('..', pdb)
   ....:     out = open(os.path.join('p38_optimize', fn + '.conf'), 'w')
   ....:     out.write(conf.format(
   ....:         out=fn, pdb=pdb,
   ....:         par=par))
   ....:     out.close()
   ....: 

Optimization

Now we will run NAMD to optimize each of these conformations. We make a list of commands that we want to execute:

In [19]: os.chdir('p38_optimize')  # we will run commands in this folder

In [20]: cmds = []

In [21]: for conf in glob.glob('*.conf'):
   ....:     fn = os.path.splitext(conf)[0]
   ....:     cmds.append('namd2 ' + conf + ' > ' + fn + '.log')
   ....: 

In [22]: cmds[:2]
Out[22]: ['namd2 p38_28.conf > p38_28.log', 'namd2 p38_13.conf > p38_13.log']

We will run these commands using multiprocessing module. We will allocate 3 processors for the job:

In [23]: from multiprocessing import Pool

In [24]: pool = Pool(3) # number of CPUs to use

In [25]: signals = pool.map(os.system, cmds)

signals will collect the output from execution of NAMD. If everything goes right, we should have only 0s.

In [26]: set(signals)
Out[26]: {0}

All NAMD output should be in p38_optimize folder. We go back to origional folder as follows:

In [27]: os.chdir('..')