Copying angles (pose.set_phi()) from a substructure with the same sequence to pose doesn’t produce the exact same substructure?

[Anonymous]

Dear Rosetta Gurus,

I’m writing an extension to the Abinitio Protocol, the code I will be describing is called exactly after setup_fold() in AbrelaxApplication. Now let me describe the input and what I’m trying to do:

I use as input two poses, let’s call them “model” and “prediction”, both are sequencially identical.
“prediction” is the extended pose created in AbrelaxApplication::setup_fold(), “model” is a pose my class loads exterally from a pdb file of an idealized protein in a Centroid fashion.
I’m trying to copy the structure of sequentially continuous chunks/parts from “model” into “prediction”, I do this by iterating over the residues of the chunk while calling
prediction.set_phi(index, model.phi(index))/set_psi()/set_omega(), you’ll find the complete function at the end of this post. After finishing this, I dump the structure of “prediction” into a pdb file.
What I expect to have, is an extended pose with my chunks having the exact same structure as “model”
Surprisingly, I find that the substructures that I’ve copied into “prediction” are not exactly the same as the ones in “model”, most alignements in pymol gives an RMSD of ~1.

Why did this happen?, and how can I copy the exact structures? (Which is very important for my case)

My theory is that it should have something to do either with the initialization in setup_fold() or the way the atomtrees work. Hopefully there is a work around.

Thanks for your help!
yours,
Mabrouk

Appendix – code of copy function:

//frag in “model in our case and pose is “prediction”
//bb_fragment is just an extension to pose with no real change

void
protocols::abinitio::bb_lib::use_fragment(bb_fragment* frag,core::pose::Pose& pose)
{
core::Size tgtlength = pose.total_residue();
core::Size srclength = frag->total_residue();

core::Size srsbegin, tgt_index;
srsbegin = frag->get_begin();
for( core::Size ires = 1; ires <= srclength ; ires++){
tgt_index = srsbegin + ires ;

pose.set_phi( tgt_index, frag->phi(ires) );
pose.set_psi( tgt_index, frag->psi(ires) );
pose.set_omega( tgt_index, frag->omega(ires) );

}
}

[Anonymous]

Can you manually check the actual values for torsions/angles/distances (say in pymol, or preferably something that offers higher precision)? From a casual glance your torsions should be correct, but you haven’t done anything to the bond lengths or bond angles, so unless they are all perfectly ideal you will accumulate many small errors (which propagate fast) from nonideal lengths/angles. I see in your post that the source and destination should both have preidealized bond geometries, but it’s worth checking.

[Anonymous]

bb_fragment is your own class, right (I can’t find it or bb_lib in 3.4, which I assume you’re using?)

Also, you’re not a developer, correct? It’s 100% awesome for users to hack in their own code, but if you’re a core developer, I’m honor-bound to point out that raw pointers (bb_fragment*) are forbidden in favor of our smart pointer system.

[Anonymous]

Is it possible you’ve got a degrees/radians mismatch? This seems unlikely given that you only get small errors.

[Anonymous]

@smlewis:
1/ After setting the angles using
pose.set_phi( tgt_index, frag->phi(ires) );
I check the angles like this:
tr.Info<< "phi n." << ires << "value in frag"<< frag->phi(ires) << " value in pose " << pose.phi(tgt_index) << std::endl;
The angles in the debug output are always the same (3 numbers after the comma)
I will check the output in pymol to see if there is a difference.
2/ Yeah bb_lib and bb_fragment are my classes. I’m not core developer. Actually I’ve tried using the smart pointer system, but I got some weird bug and after a night debugging, I choose to ignore the whole thing and pray that the segmentation won’t fault
3/ About the degrees/radian mismatch, how do I check which system is used?

Thanks for your help

[Anonymous]

And by the way I’m using Rosetta 3.4

[Anonymous]

This is really strange, I agree with you that your code ought to work fine as written.

What sort of chunks are you using? Do the chunks overlap? Do the chunks cover the whole new pose (the originally extended one)? What happens if you just assign the entire pose as one big fragment and overwrite _all_ the backbone torsions as in one loop then dump out the new PDB? Is it a normally connected monomer, or are there gaps/multiple chains/ligands/other weirdness?

[Anonymous]

The structures I’m using are proteins from the pdbs, the chunk are continuous parts from the idealized native structure, the chunks cover only a part of the structure. The Protein is normally connected with no gaps or any weirdness.
I’ll try overwriting all the angles first thing tomorrow and get back to you. Although I think it will work with no problems this way, because I believe the source of the problem should be in the update function in atomtree, which is may be influenced somehow by the extended structure when updating the angles.
By the way, everything work well when I use pose.copy_segment(srclength, *frag, frag->get_begin(), 1) instead of use_fragment(..). Of course in this case the coordinate are the one which are copied and not the angles.

[Anonymous]

Yahoo!!
I’ve solved the problem, the chunks I was using didn’t have the phi value for the first residue nor the psi and omega value for the last atoms, so when they were copied they were copied as zeros which brought some problems.
Also it seems that I copied in the wrong place (1 residue more). I’ve checked beforhand that possibility but somehow I’ve missed it… I noticed it again when I copied the whole protein.
Anyway thanks alot for the help. Now everything works perfectly

[Anonymous]

If this :

tr.Info<< "phi n." << ires << "value in frag"<< frag->phi(ires) << " value in pose " << pose.phi(tgt_index) << std::endl; checks out, then you don’t have a degrees/radians mismatch, it would show up as huge discrepancies here. For the record, pose.set_whatever takes degrees, as is documented in src/core/pose/Pose.hh. I can’t tell what your fragment class returns, of course. Where is the dump_pdb operation relative to this code? If the sequences match, you could try doing replace_residue operations instead of setting torsion angles; that will bring along most of the bond angle/length changes, if there are any. pose.replace_residue(core:Size position_to_replace, Residue residue_to_replace_with, true)

[Anonymous]

So, first, I checked the angles in pymol, there is a small differance, which proves what you said, a small difference is propagating:
dihedrals in model
(-0.67404101871493505, -0.97408505506333432), (-1.1714136099171257, -0.51712636175411353), (-1.374361118089741, -0.61138579104245194), (-1.1430177074798937, -0.59640137577712082), (-1.3291127811232226, -0.79210875054576879), (-1.1163555164463814, -0.54730221287611069),
dihedral in prediction
(-0.67357784493794992, -0.97412315669481264), (-1.1706627791225646, -0.51791671099084347), (-1.3738396224678908, -0.61242785515246234), (-1.1429102575423693, -0.59667634675645487), (-1.3289145856270514, -0.79196007390643552), (-1.1158190186632198, -0.54751335827445746),

Second, my bb_fragment class is an exact copy of pose, it inherits from it everything and adds just small variable (as name of the fragment and stuff that has nothing todo with the physics)

Third, after doing this for multiple chunks I call pose.dump_pdb(“results/tmpresult.pdb”);
So my code overall is straightforward:

for fragments frag
use_fragment(&frag, pose);
endfor
pose.dump_pdb(“results/tmpresult.pdb”);

I can’t use replace residue, because I plan in the future to copy structures with a different sequence to “prediction”. This is just a simplified test case I made to find the source of the problem, all my code is disabled except for this part.

[Anonymous]

And shouldn’t the bong length be constant in an idealized structure?

[Anonymous]

Absolutely, but it’s something to check. The subtle changes in your torsions are enough to cause the problem, though.

[Anonymous]

Glad to hear it, congrats on getting it working. I’m thrilled to hear there are folks out there extending Rosetta – we made the libraries modular for the purpose but the learning curve is so steep we aren’t sure if anyone does it or not.

[Author]

Posts