3

I want to get all the dihedral angles of a protein in Pymol (phi, psi, chi1, chi2, chi3, chi4) but I only manage to find a function that can shows me the phi and psi.

For instance:

 PyMOL>phi_psi 1a11
 SER-2:    (   67.5,  172.8 )
 GLU-3:    (  -59.6,  -19.4 )
 LYS-4:    (  -66.4,  -61.7 )
 MET-5:    (  -64.1,  -17.9 )
 SER-6:    (  -78.3,  -33.7 )
 THR-7:    (  -84.0,  -18.1 )
 ALA-8:    (  -85.7,  -40.8 )
 ILE-9:    (  -75.1,  -30.8 )
 SER-10:   (  -77.6,  -47.0 )
 VAL-11:   (  -61.3,  -27.4 )
 LEU-12:   (  -60.7,  -47.5 )
 LEU-13:   (  -71.1,  -38.6 )
 ALA-14:   (  -46.2,  -50.7 )
 GLN-15:   (  -69.1,  -47.4 )
 ALA-16:   (  -41.9,  -52.6 )
 VAL-17:   (  -82.6,  -23.7 )
 PHE-18:   (  -53.4,  -63.4 )
 LEU-19:   (  -61.2,  -30.4 )
 LEU-20:   (  -61.1,  -32.3 )
 LEU-21:   (  -80.6,  -60.1 )
 THR-22:   (  -45.9,  -34.4 )
 SER-23:   (  -74.5,  -47.8 )
 GLN-24:   (  -83.5,   11.0 )

It's missing the chiral angles. Does anyone know how to get all the dihedral angles?

Many thanks!

pippo1980
  • 2,181
  • 3
  • 14
  • 30
Daniel Bonetti
  • 2,306
  • 2
  • 24
  • 33

4 Answers4

3

You can get arbitrary dihedral angles with get_dihedral. Create four selections, each with a single atom and then use it like this:

get_dihedral s1, s2, s3, s4

It's exposed to the Python API as cmd.get_dihedral(). I suggest writing a Python script that uses this function along with cmd.iterate() to loop over residues. Create a dict so that on each residue you can look up a list of atom quadruples that define the chi-angles.

Praxeolitic
  • 22,455
  • 16
  • 75
  • 126
  • Iterate loops ove atoms see : https://pymolwiki.org/index.php/Iterate The iterate command executes a Python expression for all atoms in a selection – pippo1980 Jul 18 '23 at 18:38
1

You can easily do it in R. This is the link containing information on how to calculate the main chain and side chain Torsion/Dihedral Angles: http://thegrantlab.org/bio3d/html/torsion.pdb.html

But first you have to install the Bio3D package for R: http://thegrantlab.org/bio3d/download

After installing the package, load it by typing library(bio3d) at the R console prompt.

>library(bio3d)

This R script answers your question:

#returns the file path of the current working directory.

getwd()   

#sets the working directory to where you want.

setwd("home/R/Rscripts") 

#fetches the pdb file from the protein data bank and saves to dataframe 'pb'

pb <- read.pdb("insert PDB ID")  

#trim to protein only

pb.prot <- trim.pdb(pb, "protein") 

#calculates the torsion angles of the protein and save to dataframe 'tor'

tor <- torsion.pdb(pb.prot)  

#to get the number of rows and columns of 'tor'

dim(tor$tbl) 

#identify each row by their chain, residue ID and residue Number obtained from your PDB entry

res_label <- paste(pb.prot$atom$chain[pb.prot$calpha], pb.prot$atom$resid[pb.prot$calpha], pb.prot$atom$resno[pb.prot$calpha], sep="-") 

rownames(tor$tbl) <- res_label

#creates a table of the torsion angle

torsion <- tor$tbl  

#For example, to look at the angles for VAL, residue 223 from chain A

tor$tbl["A-VAL-223",]

#writes out the table to a file

write.table(torsion, file = "torsion_angles.txt", quote = F, sep = "\t")    

Your output file which is saved in your working directory will contain a table of the chain-resID-resNo and their corresponding phi, psi, chi1, chi2, chi3, chi4, and chi5 values. Goodluck!

Cave
  • 201
  • 1
  • 4
  • 14
0
#install bio3d library and call
library(bio3d)

#returns the file path of the current working directory.
 getwd()   

 #sets the working directory to where you want.
 setwd("home/R/Rscripts") 


 #fetches the pdb file from the protein data bank and saves to dataframe 'pb'
 pb <- read.pdb("insert PDB ID")  

 #trim to protein only
 pb.prot <- trim.pdb(pb, "protein")

 #calculates the torsion angles of the protein and save to dataframe 'tor'
 tor <- torsion.pdb(pb.prot) 

 #to get the number of rows and columns of 'tor'
 dim(tor$tbl) 

 #identify each row by their chain, residue ID and residue Number obtained from your PDB entry
 res_label <- paste(pb.prot$atom$chain[pb.prot$calpha], pb.prot$atom$resid[pb.prot$calpha], pb.prot$atom$resno[pb.prot$calpha], sep="-") 

 rownames(tor$tbl) <- res_label

 #creates a table of the torsion angle
 torsion <- tor$tbl  

 #For example, to look at the angles for VAL, residue 223 from chain A
 tor$tbl["A-GLY-65",]

 #convert "double" into a datatype
 dataframe_df=as.data.frame.matrix(torsion)

 #write dataframe to a .csv file
 write.csv(dataframe_df, file="name.csv", row.names=TRUE,col.names=TRUE)
Suraj Rao
  • 29,388
  • 11
  • 94
  • 103
0

copying from Calculate Psi/Phi and/or dihaedral Pymol - Pymol API

the /model/chain/segi/resn/resi/atom_name atom selection was something I wasnt able to find

until I stumbles into Selection Macros :

Selection Macros allow to represent a long atom selection phrase such as:

PyMOL> select pept and segi lig and chain B and resi 142 and name CA

in a more compact form:

PyMOL> select /pept/lig/B/142/CA

I tried to get phy , psi and omega ; double and triple check my code to see If I got it right :

#!/usr/bin/env python3

import pymol

from pymol import cmd, stored

pymol.finish_launching()

# pymol.finish_launching(['pymol', '-q', '-W' , '1200' , '-H' , '900' ])

cmd.fab('AAAGGAAGGVVPCAGCA', 'peptide')

# cmd.fab('AG', 'peptide')


# cmd.save('peptide.pdb' , 'peptide')

# cmd.hide('sticks' , 'peptide')

cmd.hide('everything' , 'peptide')

cmd.show('lines' , 'peptide')

a = cmd.select('selected' , 'peptide and resi 1-5 and name CA')

print('\n a --->   ' , a, '\n\n')

stored.dict = {}


cmd.iterate('selected' , "stored.dict[int(resi)] = [model, segi, chain, resn, int(resi)]")

print('\n stored.dict -------> ', stored.dict,'\n\n')

for i in stored.dict.keys() :
    
    print( i , type(i))
    

print('stored.dict.keys ', stored.dict.keys() , type(stored.dict.keys()))


def dihedral_res( res_id , selected):
    
    if int(res_id) - 1 in [i for i in selected.keys()] :
        
        print(selected[res_id][3], '   can calculate phy', '   ____________' , res_id , int(res_id))
        
        i = selected[res_id - 1]
        
        # print(i)
        
        s1 = "/{}/{}/{}/{}`{}/C".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        # print(s1)

        i = selected[res_id]

        s2 = "/{}/{}/{}/{}`{}/N".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        s3 = "/{}/{}/{}/{}`{}/CA".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        s4 = "/{}/{}/{}/{}`{}/C".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
       
        cmd.select('dihedral_phy_'+str(res_id) ,  s1 +' '+  s2 +' '+  s3 +' '+ s4)
        
        try:
            
            dihedral_val_phy = cmd.get_dihedral(s1, s2, s3, s4, state=0)
                
        except:
            
            dihedral_val_phy = None

    else : 
        
        print(selected[res_id][3], '   cannot calculate phy' , '____________' , res_id , int(res_id))
        
        dihedral_val_phy = None
        
    if int(res_id) + 1 in [i for i in selected.keys()] :
        
        print(selected[res_id][3], '   can calculate psi', '   ____________' , res_id , int(res_id))
        
        i = selected[res_id ]
        
        # print(i)
        
        s1 = "/{}/{}/{}/{}`{}/N".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        # print(s1)

        s2 = "/{}/{}/{}/{}`{}/CA".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        s3 = "/{}/{}/{}/{}`{}/C".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        i = selected[res_id + 1]
        
        s4 = "/{}/{}/{}/{}`{}/N".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
       
        cmd.select('dihedral_psi_'+str(res_id) ,  s1 +' '+  s2 +' '+  s3 +' '+ s4)
        
        try:
            
            dihedral_val_psi = cmd.get_dihedral(s1, s2, s3, s4, state=0)
                
        except:
            
            dihedral_val_psi = None

    else : 
        
        print(selected[res_id][3], '   cannot calculate psi' , '____________' , res_id , int(res_id))
        
        dihedral_val_psi = None    
        
    if int(res_id) - 1 in [i for i in selected.keys()] :
        
        print(selected[res_id][3], '   can calculate omega', '   ____________' , res_id , int(res_id))
        
        i = selected[res_id - 1]
        
        # print(i)
        
        s1 = "/{}/{}/{}/{}`{}/CA".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        s2 = "/{}/{}/{}/{}`{}/C".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        # print(s1)

        i = selected[res_id]

        s3 = "/{}/{}/{}/{}`{}/N".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
        s4 = "/{}/{}/{}/{}`{}/CA".format(i[0] , i[1] , i [2] , i[3] , i[4])
        
       
        cmd.select('dihedral_omega_'+str(res_id) ,  s1 +' '+  s2 +' '+  s3 +' '+ s4)
        
        try:
            
            dihedral_val_omega = cmd.get_dihedral(s1, s2, s3, s4, state=0)
                
        except:
            
            dihedral_val_omega = None

    else : 
        
        print(selected[res_id][3], '   cannot calculate omega' , '____________' , res_id , int(res_id))
        
        dihedral_val_omega = None
    
    return dihedral_val_phy , dihedral_val_psi , dihedral_val_omega
   
for i in stored.dict.keys() :

    print( dihedral_res( i, stored.dict) )
    

snippet of printed output for my code calculating phy, psi, omega for first 5 residues of a fake peptide AAAGGAAGGVVPCAGCA . Hope I got the directions right.

ALA    can calculate phy    ____________ 3 3
ALA    can calculate psi    ____________ 3 3
ALA    can calculate omega    ____________ 3 3
(180.0, -179.39898681640625, -180.0)

here the result hihlighting the residues involved in one of such dihedrals:

enter image description here

for the other dihedrals , I assume could works the same way. I bet there is a better way to create a function that returns the dihedrals but not an expert in Python or protein molecular graphics. Just for the sake of completeness, enclosed a table with the naming convention of the different aminoacid side chain dihedrals, hope the table is good [taken from http://www.mlb.co.jp/linux/science/garlic/doc/commands/dihedrals.html :

enter image description here

pippo1980
  • 2,181
  • 3
  • 14
  • 30