You are here

Why does the sum of the per-residue total energy not match the whole score?

3 posts / 0 new
Last post
Why does the sum of the per-residue total energy not match the whole score?
#1

At the bottom of a scored PDB file is a POSE_ENERGIES_TABLE. The sum of total, which is already a dot product with the weights, does not match the total score.

Here is an example in Python3, either with or without PyRosetta —This is a question about the C++ not the Python bindinds, hence my posting in Rosetta General

## with PyRosetta
import pandas as pd
pose = pyrosetta.pose_from_pdb('1ubq.r.pdb')
pyrosetta.rosetta.core.pose.remove_nonprotein_residues(pose) # there are no ligands
print(pose.num_chains(), pose.num_jump()) # there is 1 chain and 0 jumps.
pose.energies().clear_energies()
scorefxn = pyrosetta.get_fa_scorefxn()
single = scorefxn(pose)
scores = pd.DataFrame(pose.energies().residue_total_energies_array())
summed = sum(scores.total_score)
print(f'sum={summed}, call={single}, from_energy={pose.energies().total_energy()}')
# sum=-265.05077538177227, call=-307.92166266135933, from_energy=-307.92166266135933
pose.dump_scored_pdb('scored.pdb', scorefxn)

## agnostic
import re
from io import StringIO
with open('scored.pdb') as fh:
    block = re.search('#BEGIN_POSE_ENERGIES_TABLE\s+([\w\W]+)#END_POSE_ENERGIES_TABLE\s+',
                      fh.read())\
              .group(1)
totals = pd.read_csv(StringIO(block), sep=' ').total
print(f'expected={totals[1]}, summed={sum(totals[2:])}')
# expected=-307.922, summed=-265.05079000000006

What exactly causes this difference?

There are no non-zero terms absent from the per residue table. The above is a single chain protein with no ligands.

While I cannot figure out where `total_energies_` intercepts with the `residue_total_energies_` attribute, in the ScoreFunction.cc file there is:

pose.energies().total_energy() = pose.energies().total_energies().dot( weights_ );

Which means that residue_total_energies_ and total_energies_ are calculated separately somewhere. And beyond that, I cannot figure out.

Personally, I am interested in this as I was tring to wrap my head around the fact `scorefxn.get_sub_score` in PyRosetta is also non additive (example below) and I believe this lies at the heart of it.

trues = pyrosetta.rosetta.core.select.residue_selector.TrueResidueSelector().apply(pose)
all_sub = scorefxn.get_sub_score(pose, trues)
ResidueSpanSelector = pyrosetta.rosetta.core.select.residue_selector.ResidueSpanSelector
fore = ResidueSpanSelector(1, 38).apply(pose)
aft = ResidueSpanSelector(39, 76).apply(pose)
scorefxn.get_sub_score(pose, fore) + scorefxn.get_sub_score(pose, aft), all_sub
# -118.6, -214.1

 

Category: 
Post Situation: 
Fri, 2021-07-09 08:36
matteoferla

Because you have to decompose hydrogen bonds first. Use the PerResidueEnergyMetric (a SimpleMetric) or do it in a script with this code: 

 

    //Sum residue energies, making sure that we have pair-energies to do the summation properly

    weights =  local_scorefxn.weights() //Create the EnergyMap 

    emopts = core.scoring.methods.EnergyMethodOptions( local_scorefxn.energy_method_options() ) 
    emopts.hbond_options().decompose_bb_hb_into_pair_energies( True )
    local_scorefxn.set_energy_method_options( emopts )

 

Fri, 2021-07-09 09:04
jadolfbr

Thank you so much!
Your code makes the values identical in favour of the better value as hoped.

(I just googled the command `decompose_bb_hb_into_pair_energies` and found it is in the documentation —ops)

Mon, 2021-07-12 02:23
matteoferla