I'm new to Rosetta and learning demos.
I tried rosetta-3.5/rosetta_demos/public/calculate_e6ap_ubch7_ddgbind
to calculate ddG by mutation, but the average difference between dg_wt
and dg_mut converges to zero, of all 5 mutations in the demo.
I used the command below to calculate scores, using a relaxed structure
calculated before. The mutation_script.xml is attached, same as provided
in the demo except removing <Add mover_name=relax/>.
rosetta_scripts.linuxgccrelease
-database $ROSETTA3_DB
-s relaxed_1C4Z_0001.pdb
-ignore_unrecognized_res
-parser:protocol mutation_script.xml
-out:path:pdb output_files/
-out:path:score output_files/
-nstruct 1000
I noticed that no residues around the mutation moved by <PackRotamersMover>,
but others moved, therefore I suspect "jump=3" option of <DDG> filter,
though I don't understand the jump number meaning well.
Where should I find a description of the jump number?
Or, is there anything I forgot to calculate ddG?
Thank you,
| Attachment | Size |
|---|---|
 mutation_script_xml.txt | 705 bytes |
Jumps are simply rigid body translations and rotations. They're used by the internal foldtree to keep track of connections between chains. (And other non-chemically connected pieces of the pose.) jump=3 means the third jump of the pose. Usually this is the jump that specifies the rigid body orientation of the fourth chain. (The first chain doesn't need a rigid body degree of freedom, the first jump connects the second chain, the second jump the third chain, etc.) What the jump corresponds to can be changed by changing the foldtree, so the third jump is not always top the fourth chain, so you have to be a little careful of jump designations, especially when using foldtree manipulations.
What values are you observing that are converging to zero? One "gotcha" in Rosetta scripts is that the scores of filters in the scorefile, as well as those printed at the end of the logfile, are those for the filters applied to the pose as it is at the end of the run not when the filter was applied during the middle of the run. So what happens is that it runs through the script, does the packing, calculates the ddG, does the mutation, calculates the ddG again, and then on the final, mutated structure calculates the ddG twice more, once for the dg_wt filter and once for the dg_mut filter. It's the values for the final two calculations (done on the same structure) which are reported in the scorefile and the end of the logfile.
Thank you, rmoretii.
I'm still confused by jump.
You said "the third jump is not always top the fourth chain", then how should I confirm the value is right? In the demo example, there are two proteins E6AP (chains ABC) and UbcH7 (chain D), and they appear in order A,B,C,D in the input pdb. Therefore, jump=3 is used to calculate ddG between (ABC) & (D), isn't it?
I averaged dg_wt and dg_mut values in score.sc and their difference (dg_mut - dg_wt). They were same as the values found at the end of output pdb files. The difference is not necessarily zero for each pose, though it converged to zero (+/-0.001).
What the jump refers to depends on how the foldtree is set up. For a four chain protein (A,B,C,D), the default foldtree is currently set up to have jump 1 connect chain B to residue 1 by a rigid body transform, jump 2 connect chain C to residue 1, and jump 3 connect chain D to residue 1. So using ddG with jump=3 will calculate (ABC)-(D), jump=2 will calculate (ABD)-(C), etc.
That's what should happen now if your PDB is set up normally and there's no further foldtree manipulation. If you have chain D coming before C in the PDB, it will be chain C rather than chain D that's the last chain, so it will get jump 3 instead of chain D. Also, there are various foldtree manipulations protocols do. You may end up with a linear foldtree where jump 1 connects the end of A to the start of B, jump 2 connects the end of B to the start of C, etc, so a ddG calcuation with jump=2 would calculate (AB)-(CD) instead.
Working with the foldtrees and jumps in RosettaScripts is currently a weak point, though I believe (though am not certain) Rosetta 3.5 should have a PoseInfo filter which should print out basic pose information to the tracer when applied, including the current foldtree, which might help debug things.
--
The values of dg_mut (and of dg_wt) are supposed to be identical between the score.sc file and the end of the PDB file - it's the same calculation reported in different locations. Because you've defined dg_mut and dg_wt to be identical, they're also going to be identical to each other in those locations - though as the ddG filter involves a stochastic packing step, there will be some small numerical variation. As mentioned, the values found in the scorefile (also the ones at the end of the PDB file) are for the filters as applied to the final pose - that is, the structure which is output in the PDB. The value of "dg_wt" in the scorefile is that for the filter applied to the structure with the mutation, not that prior to the mutation.
If you want the value of the filter as it was in the middle of the protocol, you need to look through the logfile to find where the filter was called in the middle of the protocol (when the pose still had the wild type residue), and get the value from what was printed there.
Thanks again, rmoretii.
I find a culprit for the zero ddG.
Now I understand that jump target depends on the chain order in a pdb file,
and that it divide chains in pdb to 2 groups. If I want to calculate ddG for
(AB)-(CD), use jump=2 provided chains appear in order A,B,C,D in a pdb file.
I'll try PoseInfo afterward.
As you mentioned, I extracted dg_wt and dg_mut values from log files.
This time (dg_mut - dg_wt) converges a certain value, from -6 to +6.
Those are expected values. Therefore, the jump is not the culprit.
Is there a way to output a caluculated ddG value in the middle of the protocol
using RosettaScript, as in the case? I think the example script (xml) in
rosetta_demos/public/calculate_e6ap_ubch7_ddgbindshould be corrected.If you want to calculate (AB)-(CD), jump=2 probably *won't* work for you given current defaults. That's because jump 2 currently is set up to connect chain C with the ABD complex. If you do want to do that sort of ddG calculation, you need to do foldtree manipulation beforehand to set up the foldtree appropriately. For this particular example, the AtomTree mover (https://www.rosettacommons.org/manuals/archive/rosetta3.4_user_guide/Mov...) with the simple_ft=1 option would probably help you.
I unfortunately don't know of a way in the release version to store a filter value from the middle of the protocol. Your best bet is to grab the value from the tracer output.
thanks for pointing out.
I'll check the foldtree/atomtree documents, and probably ask again in another thread.
For the filter scores, I take them from a log file for the moment, with a hope for a future release.