ramachandranPlot3D

From p3d
Jump to: navigation, search

The frequency of a given phi,psi angle combination can only be obtained by using a non-redundant data set of proteins. Such a pre-compiled set of protein chains was taken directly from the PISCES web site [1]. The chosen set was cullpdb_pc30_res2.0_R0.25_d090317_chains5053, and contained 5053 unique protein chains. Files were downloaded and a symmetric link was created in order to add the Dunbrack chain identifier to the pdb file name (i.e. 1OZH.pdb.gz was linked to 1OZHA.pdb.gz). As a result the loading/parsing of the file was drastically accelerated since the p3d.protein module loaded only the one chain into memory for investigation.

ramachandranPlot3D.py

usage

ramachandranPlot3D.py <bin size for angles> <pdb file(s)>

The output is e.g.

-90.0 	 -30.0 	 9 
-90.0 	 -40.0 	 13
-60.0 	 -50.0 	 3
60.0 	 10.0 	 1
-90.0 	 140.0 	 2
-70.0 	 130.0 	 1
-80.0 	 -40.0 	 36
-130.0 	 120.0 	 1
-130.0 	 70.0 	 1

columns are phi, psi and occurence.

result

The 3D Ramchandran plot of a non-redundant protein data set using 5 deg bin size. Shown the phi and psi angles in degrees on the x and y axis, respectively. The z axis illustrated as colour code shows the frequency of the observed phi psi combination. Plots were created using using ProFit.

Ramchandran plot of observed Phi Psi combinations in non redundant protein data set

code

#!/usr/bin/env python
'''Ramachandran Plot in 3D by Ch. Fufezan 2009

usage: ramachandranPlot3D.py <bin size for angles> <pdb file>'''

import sys, os, time
from p3d import protein as protein
from collections import defaultdict as ddict

if (__name__ == '__main__'):
	if (len(sys.argv) < 3):
		print (__doc__)
		sys.exit(1)
	binSize = float(sys.argv[1])
	if 180 % binSize != 0: sys.exit(1)
	ramachandranPlot3D = ddict(int)
	pdbs = sys.argv[2:]
	StartTime = time.time()
	for k,entry in enumerate(pdbs):
		print('Analysing',entry,end='...')
		try:
			pdb = protein.Protein(entry,DunbrackNaming=True,BSPTree=False)
			''' 
			query Set alphas to have one atom per residue from both chains
			'''
			for i,alpha in enumerate(pdb.query('alpha and model 1')):
				PhiPsi = alpha.calcPhiPsi(allowAlternativeConfs=False)
				ramachandranPlot3D[((PhiPsi[0][0]//binSize)*binSize,(PhiPsi[0][1]//binSize)*binSize)] += 1
				'''
				calcPhiPsi returns a list but allowAlternativeConfs= is set to False
				the list has only one element and thus we can directly access 
				phi and psi with PhiPsi[0][0] and PhiPsi[0][1] respectively ...
				'''
			print('done - ',k)
		except:
			print('FAILED - ',k)
	for coords,value in ramachandranPlot3D.items():
		print(coords[0],'\t',coords[1],'\t',value)