NAMDAnalyzer

NAMDAnalyzer is a python based API providing several analysis routines for NAMD generated files.

Installation:

Unix and Windows

To install it within your python distribution, use

make [build] (for openMP version) or [build_cuda] (for CUDA accelerated version - recommended)
make install

or

python [setup.py, cuda_setup.py]
python setup.py install

Start the interpreter:

Initialization of the ipython console can be done using the following command:

ipython -i <path to NAMDAnalyzer.py> -- <list of files to be loaded>

Quick start

The program is organized on a master class contained in NAMDAnalyzer.Dataset.

Open log file and plot data

To analyze log file, the following can be used:

import NAMDAnalyzer as nda

data = nda.Dataset('20190326_fss_tip4p_prelax.out') #_Import log file

#_Another log file can be append to the imported one using data.logData.appendLOG() method

#_To plot data series using keywords given in data.logData.etitle, removing first 500 frames of minimization
data.logData.plotDataSeries('TEMP KINETIC TOTAL', begin=501)

#_To plot data distribution
data.logData.plotDataDistribution('KINETIC', binSize=20)

#_Data can be fitted with any model using
data.logData.plotDataDistribution('KINETIC', fit=True, model=your_model_function, p0=init_parameters)

Outputs of previous code:

_images/log_dataSeries.png _images/log_dataDist.png

Load trajectories, selection, and analysis

import NAMDAnalyzer as nda

d = nda.Dataset('psfFile.psf', 'dcdFile.dcd')

#_Trajectories can be append to already loaded ones using either d.appendDCD('dcdFile.dcd')
#_or d.appendCoordinates('pdbFile.pdb').


#_Coordinates and cell dimensions can be accessed in various ways

#_For several frames
sel = data.selection('name OH2 and within 3 of protein and resid 40:80 frame 0:100:2')

#_For one frame
sel = data.selection('bound to name OH2 and within 3 of protein and resid 40:80 frame 10')

sel.coordinates()

# or

d.dcdData[sel]

#_Coordinates can also be accessed directly using
d.dcdData[:4000,10:100:2,0] #_To get x coordinate of first 4000 atoms, in every two frame from 10 to 100

#_For cell dimensions
d.cellDims[10:100:2] #_To get corresponding cell dimensions



#_To compute RMSD per atom for molecules aligned in all frames
d.getRMSDperAtom(selection='protein and segname V1', align=True, frames=slice(0, None))

#_To compute and plot RMSD per atom for molecules aligned in all frames
d.plotRMSDperAtom(selection='protein and segname V1', align=True, frames=slice(0, None))



#_To compute radial pair distribution function for water around a protein region
from NAMDAnalyzer.dataAnalysis.RadialDensity import RadialNumberDensity

rdf = RadialNumberDensity( 'name OH2 and within 3 of protein and resid 40:80',
                                   'name OH2 and within 3 of protein and resid 40:80',
                                   dr=0.1, maxR=15, frames=range(0,1000,5) )

rdf.compDensity()
rdf.plotDensity()

#_To compute radial pair distribution density for water around each residue
from NAMDAnalyzer.dataAnalysis.RadialDensity import ResidueWiseWaterDensity

rdf = ResidueWiseWaterDensity( 'protein', dr=0.1, maxR=15, frames=range(0,1000,5) )

rdf.compDensity()
rdf.plotDensity()


#_This can be linked to residue wise residence time
from NAMDAnalyzer.dataAnalysis.ResidenceTime import ResidenceTime

rt = ResidenceTime(data, 'name OH2 and within 3 of protein and segname V5')
rt.compResidueWiseResidenceTime(dt=25)
rt.plotResidueWiseResidenceTime()



#_To plot averaged distances between a residue and the rest of the protein using a parallel plot
d.plotAveragedDistances_parallelPlot('protein and resid 53', 'protein', maxDist=10, step=2)

#_To plot the same distances but using a chord diagram
cd = d.plotAveragedDistances_chordDiagram('protein and resid 53', 'protein', maxDist=10, step=2)
cd.show()

Outputs of previous code:

_images/ubq_rmsdPerAtom.png _images/radialDistWater.png _images/residueWise_radial3D.png _images/averagedDistances_parallel.png _images/averagedDistances_chord.png

Analysis of rotations

import NAMDAnalyzer as nda
from NAMDAnalyzer.dataAnalysis.Rotations import Rotations

d = nda.Dataset('psfFile.psf', 'dcdFile.dcd')


#_To analyze O-H1 water vectors for O being within 3 angstrom of protein region
rot = Rotations(d, 'name OH2 and within 3 of protein and resid 40:80',
                   'name H1 and bound to name OH2 and within 3 of protein and resid 40:80',
                axis='z', nbrTimeOri=20)

rot.compRotationalRelaxation()
rot.compOrientationalProb()

rot.plotRotationalRelaxation()
rot.plotOrientationalProb()

#_To analyze water dipole moment orientation relative to protein surface
from NAMDAnalyzer.dataAnalysis.Rotations import WaterAtProtSurface

surf = WaterAtProtSurface(d, protSel='protein and segname V2 V3 V4', minR=1, maxR=5, maxN=5, watVec='D')
surf.compOrientations()
surf.generateVolMap()
surf.writeVolMap('myFileName', frame=0)

surf.plotOrientations()

#_A PDB file and a DX volumetric map file are generated and can be directly imported into VMD
#_for visualization exactly the same way it is done with APBS.

Outputs of previous code:

_images/rotRelaxation.png _images/rotOrientationProb.png _images/water_orient.png _images/water_orientations_hist.png

Analysis of hydrogen bonds

import NAMDAnalyzer as nda
from NAMDAnalyzer.dataAnalysis.HydrogenBonds import HydrogenBonds

d = nda.Dataset('psfFile.psf', 'dcdFile.dcd')

#_To analyze hydrogen bonds auto-correlation
#_The 'hydrogens' argument is optional, if None, they are obtained from hydrogens bound to donors
#_maxTime is tha maximum number of frame, maxR is the maximum distance for acceptor, hydrogen distance
#_step is the frame increment from origin to maxTime, minAngle is the minimum angle to accept hydrogen bond
#_between acceptor-hydrogen and donor-hydrogen vectors

hb = HydrogenBonds(d, donors='name OH2', acceptors='name OH2', hydrogens=None, maxTime=50
                    nbrTimeOri=20, step=1, maxR=2.5, minAngle=130)

#_For continuous auto-correlation (default if 'continuous' not given)
hb.compAutoCorrel(continuous=1)

#_For intermittent auto-correlation
hb.compAutoCorrel(continuous=0)

#_To plot the result
hb.plotAutoCorrel(corrType='continuous')
hb.plotAutoCorrel('intermittent')

Outputs of previous code:

_images/hbContinuous.png _images/hbIntermittent.png

Mean-squared displacement and neutron backscattering

import NAMDAnalyzer as nda
from NAMDAnalyzer.dataAnalysis.backscatteringDataConvert import BackScatData

d = nda.Dataset('psfFile.psf', 'dcdFile.dcd')


#_Defines some q-values for incoherent scattering function
qVals = [0.2, 0.4, 0.6, 0.8, 1, 1.2, 1.4, 1.6, 1.8]

bs = BackScatData(d)


#_To compute MSD for non exchangeable hydrogens in protein for increasing time steps,
#_without center of mass motion
msd = []

for frame in range(0, 200, 5):
    bs.compMSD(frameNbr=frame, selection='protNonExchH', alignCOM=True)
    msd.append( bs.MSD )

import matplotlib.pyplot as plt

times = np.arange(0, 200, 5) * d.timestep * d.dcdFreq[0:200:5] * 1e9
msd   = np.array(msd)

plt.plot(times, msd[:,0])
plt.xlabel('Time [ns]')
plt.ylabel('MSD [$\AA^{2}$]')

plt.show()


#_To compute and plot incoherent intermediate function, EISF and inoherent scattering
#_function for water hydrogens with 200 time steps

bs.compScatteringFunc(qVals, nbrTimeOri=50, selection='waterH', alignCOM=True, nbrTS=200)

bs.plotIntermediateFunc()
bs.plotEISF()
bs.plotScatteringFunc()

Outputs of previous code:

_images/msd.png _images/bs_interFunc.png _images/bs_EISF.png _images/bs_scatFunc.png

Indices and tables