<p>OK, I tried to answer half of your question. I made a test answer showing how to get the delta-phi , delta-psi between two structure aligned with <strong>PyMOL</strong>, here my code , not sure it is the fastest/cleanest way  but only one I could figure out.
If you align more structures to a reference you will end up with a raw_alignment that looks like:</p>
<pre><code>[('3hr8', 1988), ('1xmv', 1764), ('1xp8', 1776)]
[('3hr8', 1997), ('1xmv', 1773), ('1xp8', 1785)]
[('3hr8', 2001), ('1xmv', 1777), ('1xp8', 1789)]
[('3hr8', 2017), ('1xp8', 1808)]
[('3hr8', 2028), ('1xmv', 1804), ('1xp8', 1817)]
[('3hr8', 2049), ('1xp8', 1831)]
[('3hr8', 2058), ('1xmv', 1829), ('1xp8', 1839)]
[('3hr8', 2093), ('1xmv', 1860)]
[('3hr8', 2134), ('1xmv', 1906)]
[('3hr8', 2141), ('1xmv', 1914)]
[('3hr8', 2205), ('1xmv', 1972)]
[('3hr8', 2217), ('1xmv', 1978)]
[('3hr8', 2229), ('1xmv', 1990)]
[('3hr8', 2295), ('1xmv', 2021)]
[('3hr8', 2345), ('1xp8', 2081)]
[('3hr8', 2352), ('1xp8', 2089)]
[('3hr8', 2361), ('1xmv', 2085), ('1xp8', 2094)]
[('3hr8', 2372), ('1xmv', 2099), ('1xp8', 2106)]
[('3hr8', 2380), ('1xp8', 2114)]
[('3hr8', 2389), ('1xmv', 2116), ('1xp8', 2119)]
[('3hr8', 2397), ('1xmv', 2124)]
</code></pre>
<p>so you'll have to use <code>cmd.iterate</code> over it selecting all the possible couples of reference-target structures and get their CA atoms coupled-indexes extract the corresponding resi (residue number/index) and use <code>phi_psi</code> command to get the two dihedral of interest, like it is done below:</p>
<pre><code>import pymol


from pymol import cmd

from pymol import stored

import time


pymol.finish_launching()


# get two structures
# fetch from remote
cmd.fetch('3hr8 ')
cmd.fetch('3hr8 1xp8')


# load from locale
# cmd.load('3hr8.cif' , '3hr8')
# cmd.load('1xp8.cif', '1xp8')

time.sleep(3)

cmd.remove('hetatm')


### Super and get raw alignment ###
cmd.do('super 1xp8, 3hr8, object = aln')

print('super ended !!!!')


### get raw alignment ###
raw_aln = cmd.get_raw_alignment('aln')

cmd.do('set seq_view , 1')

# print residue pairs (atom index)
for idx1, idx2 in raw_aln:
    print('%s`%d -&gt; %s`%d' % tuple(idx1 + idx2))


## To print residue numbers instead of atom indices:
def myfunc(resi,model,index,name):

    if name == 'CA':
        
        idx2resi[model, index] = resi
        
    else :
        
        idx2resi[model, index] = None
        
    
idx2resi = {}
cmd.iterate('aln', 'myfunc(resi,model,index,name)', space={'idx2resi': idx2resi , 'myfunc': myfunc})


# print(len(idx2resi), len(raw_aln))

errors = []

# print residue pairs (residue number)
for idx1, idx2 in raw_aln:
    
    # print(idx1 , idx2)

    if idx2resi[idx1] and idx2resi[idx2] :
        
            a = idx1[0]+'...'+idx2resi[idx1]
            
            b = idx2[0]+'...'+idx2resi[idx2]
            
            a_die = cmd.phi_psi(f'{idx1[0]} and resi {idx2resi[idx1]}' )
            
            if a_die :
                
                # print(a_die , ' OK')
                pass
                
            else :
                
                a_die[idx1] = (0,000000 , 0,00000000)
                
                # print('\n\n error set to 0 !!!!!!!!' , a_die)
                
                errors.append((idx1 , idx2resi[idx1]))
            
            b_die = cmd.phi_psi(f'{idx2[0]} and resi {idx2resi[idx2]}' )
            
            if b_die :
                
                # print(b_die , ' OK')
                pass
            
            else :
                
                b_die[idx2] = (0,000000 , 0,00000000)
                
                # print('\n\n error set to 0 !!!!!!!!' , b_die)
                
                errors.append((idx2, idx2resi[idx2]))
            
            # print(a, b)
            
            # print(a_die, b_die)

            stored.a_die_resi = str()
            
            stored.b_die_resi = str()
            
            
            cmd.iterate(f&quot;{idx1[0]} and resi {idx2resi[idx1]}&quot; , &quot;stored.a_die_resi = resn&quot; )
            
            cmd.iterate(f&quot;{idx2[0]} and resi {idx2resi[idx2]}&quot; , &quot;stored.b_die_resi = resn&quot; )
        
            # print('%s -&gt; %s' % (a , b ), ' phi_psi :' , a_die)
            
            s = (f&quot;{a} {stored.a_die_resi} ---&gt; phi_psi: {a_die[idx1]}&quot;
                 
            f&quot;   {b} {stored.b_die_resi} ---&gt; phi_psi: {b_die[idx2]}&quot;
            
            f&quot;   DELTA : {abs(a_die[idx1][0] -b_die[idx2][0])}  , {abs(a_die[idx1][1] -b_die[idx2][1])}&quot;)
            
            print(s)
            


print('\n\n Errors : \n ', errors)
</code></pre>
<p>If for same reason (gap in structure) one of the residues doesnt return phi_psi its value gets assigned to 0 and is appended in the error list, probably would be better to remove it from the results, but took already a lot of time to answer this ... so ...).</p>
<p>Last bit of output:</p>
<pre><code>3hr8...313 VAL ---&gt; phi_psi: (-64.3411865234375, -44.5438232421875)   1xp8...319 ILE ---&gt; phi_psi: (-64.2035903930664, -27.23161506652832)   DELTA : 0.13759613037109375  , 17.31220817565918
3hr8...314 GLN ---&gt; phi_psi: (-61.05774688720703, -39.42118453979492)   1xp8...320 ALA ---&gt; phi_psi: (-75.51972198486328, -28.401968002319336)   DELTA : 14.46197509765625  , 11.019216537475586
3hr8...315 PHE ---&gt; phi_psi: (-56.17056655883789, -43.013572692871094)   1xp8...321 TYR ---&gt; phi_psi: (-71.33622741699219, -41.687713623046875)   DELTA : 15.165660858154297  , 1.3258590698242188
3hr8...316 LEU ---&gt; phi_psi: (-62.212162017822266, -40.03855895996094)   1xp8...322 ILE ---&gt; phi_psi: (-69.78271484375, -35.62084197998047)   DELTA : 7.570552825927734  , 4.417716979980469
3hr8...317 LYS ---&gt; phi_psi: (-62.3416633605957, -42.122344970703125)   1xp8...323 ALA ---&gt; phi_psi: (-59.84389114379883, -19.88322639465332)   DELTA : 2.497772216796875  , 22.239118576049805
3hr8...318 ASP ---&gt; phi_psi: (-77.59044647216797, -13.573332786560059)   1xp8...324 GLU ---&gt; phi_psi: (-95.1947250366211, -7.630105972290039)   DELTA : 17.604278564453125  , 5.9432268142700195


 Errors : 
  [(('1xp8', 1386), '223')]
</code></pre>
<p>here a nice pics of <strong>PyMOL</strong> showing the correspondence between aligned/superposed atoms :</p>
<p><a href="https://i.sstatic.net/Gg75e8QE.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/Gg75e8QE.png" alt="enter image description here" /></a></p>
<p>try playing with the alignment algorithm and different selections for target and reference and see if you get any different results:</p>
<p><code>cmd.super( 'mobile' , 'target' , object='aln') </code></p>
<p>vs</p>
<p><code>cmd.align( 'mobile' , 'target' , object='aln') </code></p>
<pre><code>cmd.select('target' , '/3hr8//A//CA')
cmd.select('mobile' , '/1xp8//A//CA')
cmd.super( 'mobile' , 'target'  , cycles=0, transform = 0, object='aln') 
</code></pre>
<p>vs</p>
<pre><code>cmd.select('target' , '3hr8')
cmd.select('mobile' , '1xp8')
cmd.super( 'mobile' , 'target'  , cycles=0, transform = 0, object='aln') 
</code></pre>
<p>and cycles parameter :</p>
<blockquote>
<p>cycles = int: maximum number of outlier rejection cycles {default: 5}</p>
</blockquote>

CC BY-SA 4.0
2024-05-24T17:40:37.303Z
2024-05-24T19:38:17.613Z
2024-05-24T19:38:17.613Z
pippo1980
<p>I aligned multiple protein structures from the same protein family to a reference structure in pymol and now I want to calculate the Difference in Phi psi values for each 'aligned' residue pair between Reference Structure and each Target Structure, I know about extracting the phi psi values using phi_psi command and I have the sequence alignment but I don't know how to extract Delta Phi Psi values for every aligned residue between Reference and Target.</p>
<p>Here the code for fetching and generate a multiple alignment using PyMOL <a href="https://pymolwiki.org/index.php/Super" rel="nofollow noreferrer">Super</a>:</p>
<blockquote>
<p>super aligns two selections. It does a sequence-independent (unlike align) structure-based dynamic programming alignment followed by a series of refinement cycles intended to improve the fit by eliminating pairing with high relative variability (just like align). super is more robust than align for proteins with low sequence similarity.</p>
</blockquote>
<p>between PDBs <a href="https://www.rcsb.org/structure/3hr8" rel="nofollow noreferrer">3HR8</a> (reference) , <a href="https://www.rcsb.org/structure/1XMV" rel="nofollow noreferrer">1XMV</a> and <a href="https://www.rcsb.org/structure/1XP8" rel="nofollow noreferrer">1XP8</a> :</p>
<pre><code>import pymol
from pymol import cmd
import time


list_pdbs = ['1xp8' , '3hr8', '1xmv']

pymol.finish_launching()

time.sleep(5)

for _i in list_pdbs :
    try :
        
        cmd.fetch(_i)
        
    except :
        
        print('cant load ', _i)
              
cmd.remove('hetatm')

sel = 'name CA'
ref = list_pdbs[1]
method = 'super'

cmd.do(f&quot;extra_fit {sel} , {ref} , {method}, object = aln_super_all&quot;)

cmd.do('set seq_view , 1')

cmd.select('aln_super_all')

print(cmd.get_type('aln_super_all'))

</code></pre>

CC BY-SA 4.0
2024-02-02T15:52:04.670Z
Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking.
Community
CC BY-SA 4.0
2024-05-22T17:02:59.887Z
and at list provide the first half of the data i.e. the proteins lets say one reference and two other from the family and how you aligned them
pippo1980
<pymol>
How to Calculate Delta Phi Psi between aligned pair of residue between Reference and Target Protein Structure in Pymol
# How to Calculate Delta Phi Psi between aligned pair of residue between Reference and Target Protein Structure in Pymol

(Asked by: username_0 on 2024-02-02)

I aligned multiple protein structures from the same protein family to a reference structure in pymol and now I want to calculate the Difference in Phi psi values for each 'aligned' residue pair between Reference Structure and each Target Structure, I know about extracting the phi psi values using phi_psi command and I have the sequence alignment but I don't know how to extract Delta Phi Psi values for every aligned residue between Reference and Target.
Here the code for fetching and generate a multiple alignment using PyMOL Super:

super aligns two selections. It does a sequence-independent (unlike align) structure-based dynamic programming alignment followed by a series of refinement cycles intended to improve the fit by eliminating pairing with high relative variability (just like align). super is more robust than align for proteins with low sequence similarity.

between PDBs 3HR8 (reference) , 1XMV and 1XP8 :
import pymol
from pymol import cmd
import time


list_pdbs = ['1xp8' , '3hr8', '1xmv']

pymol.finish_launching()

time.sleep(5)

for _i in list_pdbs :
    try :
        
        cmd.fetch(_i)
        
    except :
        
        print('cant load ', _i)
              
cmd.remove('hetatm')

sel = 'name CA'
ref = list_pdbs[1]
method = 'super'

cmd.do(f"extra_fit {sel} , {ref} , {method}, object = aln_super_all")

cmd.do('set seq_view , 1')

cmd.select('aln_super_all')

print(cmd.get_type('aln_super_all'))

**Tags:** pymol
> **Comments:**
> **username_1**: Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking.  
> **username_2**: and at list provide the first half of the data i.e. the proteins lets say one reference and two other from the family and how you aligned them  

---

## 1 answers

### Answer by username_2 on 2024-05-24. Score: 0

OK, I tried to answer half of your question. I made a test answer showing how to get the delta-phi , delta-psi between two structure aligned with PyMOL, here my code , not sure it is the fastest/cleanest way  but only one I could figure out.
If you align more structures to a reference you will end up with a raw_alignment that looks like:
[('3hr8', 1988), ('1xmv', 1764), ('1xp8', 1776)]
[('3hr8', 1997), ('1xmv', 1773), ('1xp8', 1785)]
[('3hr8', 2001), ('1xmv', 1777), ('1xp8', 1789)]
[('3hr8', 2017), ('1xp8', 1808)]
[('3hr8', 2028), ('1xmv', 1804), ('1xp8', 1817)]
[('3hr8', 2049), ('1xp8', 1831)]
[('3hr8', 2058), ('1xmv', 1829), ('1xp8', 1839)]
[('3hr8', 2093), ('1xmv', 1860)]
[('3hr8', 2134), ('1xmv', 1906)]
[('3hr8', 2141), ('1xmv', 1914)]
[('3hr8', 2205), ('1xmv', 1972)]
[('3hr8', 2217), ('1xmv', 1978)]
[('3hr8', 2229), ('1xmv', 1990)]
[('3hr8', 2295), ('1xmv', 2021)]
[('3hr8', 2345), ('1xp8', 2081)]
[('3hr8', 2352), ('1xp8', 2089)]
[('3hr8', 2361), ('1xmv', 2085), ('1xp8', 2094)]
[('3hr8', 2372), ('1xmv', 2099), ('1xp8', 2106)]
[('3hr8', 2380), ('1xp8', 2114)]
[('3hr8', 2389), ('1xmv', 2116), ('1xp8', 2119)]
[('3hr8', 2397), ('1xmv', 2124)]

so you'll have to use cmd.iterate over it selecting all the possible couples of reference-target structures and get their CA atoms coupled-indexes extract the corresponding resi (residue number/index) and use phi_psi command to get the two dihedral of interest, like it is done below:
import pymol


from pymol import cmd

from pymol import stored

import time


pymol.finish_launching()


# get two structures
# fetch from remote
cmd.fetch('3hr8 ')
cmd.fetch('3hr8 1xp8')


# load from locale
# cmd.load('3hr8.cif' , '3hr8')
# cmd.load('1xp8.cif', '1xp8')

time.sleep(3)

cmd.remove('hetatm')


### Super and get raw alignment ###
cmd.do('super 1xp8, 3hr8, object = aln')

print('super ended !!!!')


### get raw alignment ###
raw_aln = cmd.get_raw_alignment('aln')

cmd.do('set seq_view , 1')

# print residue pairs (atom index)
for idx1, idx2 in raw_aln:
    print('%s`%d -> %s`%d' % tuple(idx1 + idx2))


## To print residue numbers instead of atom indices:
def myfunc(resi,model,index,name):

    if name == 'CA':
        
        idx2resi[model, index] = resi
        
    else :
        
        idx2resi[model, index] = None
        
    
idx2resi = {}
cmd.iterate('aln', 'myfunc(resi,model,index,name)', space={'idx2resi': idx2resi , 'myfunc': myfunc})


# print(len(idx2resi), len(raw_aln))

errors = []

# print residue pairs (residue number)
for idx1, idx2 in raw_aln:
    
    # print(idx1 , idx2)

    if idx2resi[idx1] and idx2resi[idx2] :
        
            a = idx1[0]+'...'+idx2resi[idx1]
            
            b = idx2[0]+'...'+idx2resi[idx2]
            
            a_die = cmd.phi_psi(f'{idx1[0]} and resi {idx2resi[idx1]}' )
            
            if a_die :
                
                # print(a_die , ' OK')
                pass
                
            else :
                
                a_die[idx1] = (0,000000 , 0,00000000)
                
                # print('\n\n error set to 0 !!!!!!!!' , a_die)
                
                errors.append((idx1 , idx2resi[idx1]))
            
            b_die = cmd.phi_psi(f'{idx2[0]} and resi {idx2resi[idx2]}' )
            
            if b_die :
                
                # print(b_die , ' OK')
                pass
            
            else :
                
                b_die[idx2] = (0,000000 , 0,00000000)
                
                # print('\n\n error set to 0 !!!!!!!!' , b_die)
                
                errors.append((idx2, idx2resi[idx2]))
            
            # print(a, b)
            
            # print(a_die, b_die)

            stored.a_die_resi = str()
            
            stored.b_die_resi = str()
            
            
            cmd.iterate(f"{idx1[0]} and resi {idx2resi[idx1]}" , "stored.a_die_resi = resn" )
            
            cmd.iterate(f"{idx2[0]} and resi {idx2resi[idx2]}" , "stored.b_die_resi = resn" )
        
            # print('%s -> %s' % (a , b ), ' phi_psi :' , a_die)
            
            s = (f"{a} {stored.a_die_resi} ---> phi_psi: {a_die[idx1]}"
                 
            f"   {b} {stored.b_die_resi} ---> phi_psi: {b_die[idx2]}"
            
            f"   DELTA : {abs(a_die[idx1][0] -b_die[idx2][0])}  , {abs(a_die[idx1][1] -b_die[idx2][1])}")
            
            print(s)
            


print('\n\n Errors : \n ', errors)

If for same reason (gap in structure) one of the residues doesnt return phi_psi its value gets assigned to 0 and is appended in the error list, probably would be better to remove it from the results, but took already a lot of time to answer this ... so ...).
Last bit of output:
3hr8...313 VAL ---> phi_psi: (-64.3411865234375, -44.5438232421875)   1xp8...319 ILE ---> phi_psi: (-64.2035903930664, -27.23161506652832)   DELTA : 0.13759613037109375  , 17.31220817565918
3hr8...314 GLN ---> phi_psi: (-61.05774688720703, -39.42118453979492)   1xp8...320 ALA ---> phi_psi: (-75.51972198486328, -28.401968002319336)   DELTA : 14.46197509765625  , 11.019216537475586
3hr8...315 PHE ---> phi_psi: (-56.17056655883789, -43.013572692871094)   1xp8...321 TYR ---> phi_psi: (-71.33622741699219, -41.687713623046875)   DELTA : 15.165660858154297  , 1.3258590698242188
3hr8...316 LEU ---> phi_psi: (-62.212162017822266, -40.03855895996094)   1xp8...322 ILE ---> phi_psi: (-69.78271484375, -35.62084197998047)   DELTA : 7.570552825927734  , 4.417716979980469
3hr8...317 LYS ---> phi_psi: (-62.3416633605957, -42.122344970703125)   1xp8...323 ALA ---> phi_psi: (-59.84389114379883, -19.88322639465332)   DELTA : 2.497772216796875  , 22.239118576049805
3hr8...318 ASP ---> phi_psi: (-77.59044647216797, -13.573332786560059)   1xp8...324 GLU ---> phi_psi: (-95.1947250366211, -7.630105972290039)   DELTA : 17.604278564453125  , 5.9432268142700195


 Errors : 
  [(('1xp8', 1386), '223')]

here a nice pics of PyMOL showing the correspondence between aligned/superposed atoms :

try playing with the alignment algorithm and different selections for target and reference and see if you get any different results:
cmd.super( 'mobile' , 'target' , object='aln') 
vs
cmd.align( 'mobile' , 'target' , object='aln') 
cmd.select('target' , '/3hr8//A//CA')
cmd.select('mobile' , '/1xp8//A//CA')
cmd.super( 'mobile' , 'target'  , cycles=0, transform = 0, object='aln') 

vs
cmd.select('target' , '3hr8')
cmd.select('mobile' , '1xp8')
cmd.super( 'mobile' , 'target'  , cycles=0, transform = 0, object='aln') 

and cycles parameter :

cycles = int: maximum number of outlier rejection cycles {default: 5}


---

