You are here

Different Pose Energy Values before and after Dumping

9 posts / 0 new
Last post
Different Pose Energy Values before and after Dumping
#1
Im developing a Pyrosetta MC de-novo protein fold prediction method and at the end of the script I 'd want to print the energy values and RMSD from the native structure of my simulations and dump their poses as pdb files. so after doing that i noticed that if I import (pose_from_pdb) these poses and rescore them (using the same energy function I used in the script to score their energies) I get different energy values. Here i list the final part of my script which scores and dumps the poses: print_counter=0 f=open("output.dat","w") while(print_counter"smaller-than"decoy_number):{ #decoy_counter= number of simulations I am running in parallel print_pose.assign(pose_vector[print_counter]) #pose_vector[] = vector where the poses are stored print_rmsd_pose=rmsd_vector[print_counter] #rmsd_vector[]= vector where their RMSD values(previously calculated) are stored scorefnx_low(print_pose) J=print_pose.energies().total_energy() f.write(str(print_counter)+" "+str(J)+" "+str(print_rmsd_pose)+"\n") print_pose.dump_pdb("output_pose"+str(print_counter)+".pdb") print_counter+=1 } f.close() Here I list the commands I use to import the poses from pdb and rescore them: to_centroid = SwitchResidueTypeSetMover( 'centroid' ) scorefxn_low3 = create_score_function( 'score3' ) Pose=pose_from_pdb('pose_something.pdb') to_centroid.apply(Pose) scorefxn_low3(Pose) The Energy values obtained doing this are quite different(e.g. energy value on file = 77.384237, energy value on screen = 98.28273). Anyone has any idea about how this can be possible? Thanks in andvance Alfredo
Post Situation: 
Wed, 2013-02-06 09:35
oppopomoz

With fullatom poses, this is a normal behavior (although not of this _magnitude_). It is caused by the fact that the PDB spec is only carries 3 decimal digits (a coordinate can't be more precise than, say, 34.224), but Rosetta uses way more precision. Changes of this size would still be concerning.

With centroid...who knows? I wouldn't guess the centroid scorefunction is that much more sensitive but it might be. Are you using this one:? Is it the same in both cases?

# score3 from rosetta++, used in stage 4 of the ClassicAbinitio protocol.
env 1.0
pair 1.0
cbeta 1.0
vdw 1.0
rg 3.0
cenpack 1.0
hs_pair 1.0
ss_pair 1.0
rsigma 1.0
sheet 1.0

STRAND_STRAND_WEIGHTS 0 6

If you've got constraints in your protocol, I can _very_ readily see them causing this problem if they aren't reloaded.

Are the two poses that score differently structurally identical? You've done rather a lot of operations (a fullatom->centroid conversion). There is a function somewhere that compares pose coordinates I can find if you like (or just diff the ATOM lines of the PDBs).

Wed, 2013-02-06 09:49
smlewis

Hey, no fair, you changed the original thread title! I liked that title.

Wed, 2013-02-06 09:49
smlewis

Tnx for answering.

I just found out that this difference appears with the fragment insertion i perform during the folding simulation so as i insert a fragment my simulations wud end up presenting this difference as output. this doesnt happen performing another kinds of movers.

here my mover:

movemap = MoveMap()
movemap.set_bb( True )

long_frag_mover = ClassicFragmentMover( fragset_long , movemap )

scorefnx_low = create_score_function( 'score3' )

folding_mover = SequenceMover()
folding_mover.add_mover( long_frag_mover )
mc = MonteCarlo( test_pose , scorefnx_low , kT )
trial = TrialMover( folding_mover , mc )
folding = RepeatMover( trial , cycles )
mc.reset( test_pose )

folding.apply( test_pose )

Here more details:
I only use centroid poses and i score them using the score3 u just mentioned, no constraints. The two poses seem to be identical cos when i calculate the RMSD between them i get 0.0A (using CA_rmsd(pose-x,pose-y) which calculates the RMSD through the backbone only i guess, am i right?). if u can find me that function it wud be cool! tnx a lot!!!!

ps deeply sorry for the change of title, i'll make up for it ;)

Thu, 2013-02-07 09:10
oppopomoz

The function core.scoring.CA_rmsd() will only compute the rmsd over the C alpha atoms, not the entire backbone. There are other functions in the score.scoring namespace that will do other rmsd calculations, like bb_rmsd(), bb_rmsd_including_O(), all_atom_rmsd(), etc.

The function Steven is probably thinking about is core.pose.compare_atom_coordinates( pose1, pose2, n_dec_places=3).

By the way, if you want to dump poses but keep all the internal coordinate precision, the best way to do that is to use a binary silent file. See the core.io.silent.SilentFileData class for the io file interface, and the core.io.silent.BinaryProteinSilentStruct class for the per-structure "binary" format representation that goes into the SilentFileData container class.

Thu, 2013-02-07 11:16
rmoretti

"The two poses seem to be identical cos when i calculate the RMSD between them i get 0.0A (using CA_rmsd(pose-x,pose-y) which calculates the RMSD through the backbone only i guess, am i right?)"

This calculates only alpha carbon RMSD. It's unlikely for two centroid poses to have their CA atoms match if the other backbone atoms do not match. You could have a differently placed terminus on one that causes a clash.

pose coordinates:

The function in C++ is core::pose::compare_atom_coordinates, here's the header entry:

///@brief this function compares pose atom coordinates for equality; it is not the == operator because it does not compare all pose data.
///@author Steven Lewis
///@param[in] lhs one pose to compare
///@param[in] rhs one pose to compare
///@param[in] n_dec_places number of decimal places to compare for the coordinates (remember == doesn't work for float); defaults to 3 which is PDB accuracy
bool compare_atom_coordinates(
Pose const & lhs,
Pose const & rhs,
Size const n_dec_places = 3);

I assume it works in python.

scorefunctions:

Unfortunately, like Pose, ScoreFunction does not have an == operator to ensure you are using precisely the same scorefunction in each case...consider the ScoreFunction.show() operation to double check?

fragments:
no idea why this would be related to fragments.

Thu, 2013-02-07 11:36
smlewis

Thanks a lot for the answers. I indentified where the error appears:

After inserting the fragments if i rescore the pose and print the energy contributions i obtein something like :

total_energy vdw 67.615
total_energy cenpack 10.133
total_energy pair -3.288
total_energy env -7.894
total_energy cbeta 33.045
total_energy rg 13.846
total_energy hs_pair 0.000
total_energy ss_pair 0.000
total_energy rsigma 0.000
total_energy sheet 3.776

after dumping , reloading , rescoring the pose and printing its energy contributions, i obtein:

total_energy vdw 67.615
total_energy cenpack 10.133
total_energy pair -3.288
total_energy env -7.894
total_energy cbeta 33.045
total_energy rg 13.846
total_energy hs_pair 0.000
total_energy ss_pair 0.000
total_energy rsigma 0.000
total_energy sheet 0.000

as u can see the sheet contribution becomes 0.000 corrisponding to the difference of energy which appears

HERE I LIST THE COMMAND I RAN:

to_centroid = SwitchResidueTypeSetMover( 'centroid' )
scorefxn_low = create_score_function( 'score3' )
pose1=pose_from_pdb('pose.pdb')
to_centroind.apply(pose1)
scorefxn_low(pose1)
->72.8439 #i just made up this number
fragset = ConstantLengthFragSet(3)
fragset.read_fragment_file(“test_in_3mer.frag”)
movemap = MoveMap()
movemap.set_bb(True)
frag_mover = ClassicFragmentMover(fragset, movemap)
frag_mover.apply(pose1) #once or more times
->67.8439 #i just made up this number

pose1.energies().show()

total_energy vdw 67.615
total_energy cenpack 10.133
total_energy pair -3.288
total_energy env -7.894
total_energy cbeta 33.045
total_energy rg 13.846
total_energy hs_pair 0.000
total_energy ss_pair 0.000
total_energy rsigma 0.000
total_energy sheet 3.776

pose1.dump_pdb(pose.pdb)

pose2=pose_from_pdb('pose.pdb')

to_centroind.apply(pose2)
scorefxn_low(pose2)
->64.1028

pose2.energies().show()

total_energy vdw 67.615
total_energy cenpack 10.133
total_energy pair -3.288
total_energy env -7.894
total_energy cbeta 33.045
total_energy rg 13.846
total_energy hs_pair 0.000
total_energy ss_pair 0.000
total_energy rsigma 0.000
total_energy sheet 0.000

hope someone knows why this happens

Thanks a lot in andvance

Mon, 2013-02-11 07:43
oppopomoz

I *think* I know what is happening. Looks like the fragment insertion assigns secondary structure to those residues of a pose it just replaced. When you reload the pose, the information is not there. So here is code that 'should' work to get ss back. Note that here, the energies from the ss may be different since were going to be assigning them to the whole pose, and the method of ss identification is different from that of what the fragments use.

from rosetta.core.scoring.dssp import Dssp

sec_struct = Dssp(pose)
sec_struct.insert_ss_into_pose( pose );

This should do the trick. You could run this before fragment insertion to keep the energies fairly tight together.

Mon, 2013-02-11 10:09
jadolfbr

last asnwer has been particolarly useful, Thanks a lot.

when I use the command "pose=pose_from_pdb('randompdb.pdb')" if i print its structure with "print pose.secstruct()" I obtein:

"LLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLL"

whatever "randompdb.pdb" I import. this seems pretty strange to me, does anybody know why this happens? btw i can get back the right structure using :
sec_struct = Dssp(pose)
sec_struct.insert_ss_into_pose( pose )

so this loss of secondary structure might explain why I got different energy values for the same pose before and after pdb-dumping it I guess.

Tue, 2013-02-12 09:19
oppopomoz