I've been using the ddg_monomer application with high resolution protocol flags and the alternative mutation file format (to obtain more than one mutation) on the constant domain of a TCR (alpha and beta chain).
I'm interested in 3 ddg values:
- With a single mutation on the alpha chain (mut1)
- With a single mutation on the beta chain (mut2)
- With the pair of mutations (mut1+2)
However, in one set of mutation I'm introducing two cysteines (1 on the alpha chain + 1 on the beta chain) - to yield a disulfide bond between the two chains. (These mutations are based on the observation of this specific disulfide bond at the same positions in another specie…)
So I was expecting to obtain ddg values > 0 for mut1 and mut2 and for mut3 an increase in stability (lower ddg) due to the disulfide bond.
However, I obtain a higher ddg value for the mut1+2.
ddG: description total fa_atr fa_rep fa_sol fa_intra_rep fa_elec pro_close hbond_sr_bb hbond_lr_bb hbond_bb_sc hbond_sc dslf_fa13 coordinate_constraint an
ddG: T48C 1.089 1.728 0.654 -2.855 0.003 2.586 0.007 -0.026 -0.836 0.397 1.041 0.008 0.000 0.000
ddG: S147C 1.842 -1.355 1.629 -0.665 0.015 0.092 -0.007 -0.009 0.100 0.465 1.443 -0.002 0.000 0.000
ddG: T48CS147C 2.450 1.173 0.950 -5.977 0.006 6.373 0.012 0.042 -0.595 1.295 1.387 0.005 0.000 0.000
I see in the output pdb's that the disulfide bond is not made. Is there a way to help ddg_monomer along with the formation of the disulfide bond?
I then tried to force the disulfide bond (after the ddf_monomer analysis) using
the -in:fix_disulf with a txt-file containing the positions of the cysteines:
score_jd2.linuxgccrelease -database /services/tools/rosetta/2016.20/main/database -in:file:silent mut_T48CS147C.out -in:file:fullatom -out:pdb -overwrite -score:weights talaris2014 -out:file:scorefile score_dis.sc -in:fix_disulf dis.txt
Followed by minimization of the structure with the forced disulfide bond:
relax.default.linuxgccrelease -in:file:s mut_T48CS147C_round_1_0001.pdb -score:weight talaris2014 -in:fix_disulf dis.txt
And then the score was calculated again:
score_jd2.linuxgccrelease -database /services/tools/rosetta/2016.20/main/database -in:file:s mut_T48CS147C_round_1_0001_0001.pdb -in:file:fullatom -out:pdb -overwrite -score:weights talaris2014 -out:file:scorefile score_cc_single.sc
This seemed to "fix" the problem - I now have a disulfide bond. But then I obtained a non-comparable "total score"-value, as this was relaxed. So therefore I relaxed the pdb file, with the mut1+2 without the forced disulfide bond, and obtained an ever lower "total score". Thus I once again ended up with the problem that the introduction of the disulfide bond was not favorable (as I expected it to be).
Are there anything you can suggest to solve this ddg_monomer+disulfide problem for me?
The ddg_monomer protocol really isn't set up for considering the introduction of disulfides. So you're correct that that's probably not the best way to go about things.
I think that the approach you take is a reasonable one: introduce the cysteines by some method (ddg_monomer is as good as any), then impose the disulfide and relax. As you've surmised, you need to treat all the structures you're comparing similarly. That is, you need to put each of the single mutants and the non-cysteine parent molecule through the same relax process as the disulfide structure. (Comparing the two-cys-non-disulfide protein to the disulfide protein is only potentially interesting
As to why the disulfide isn't being predicted as favorable, I'm not quite sure. The first thing is that the Rosetta disulfide energy hasn't really been set up to predict disulfide formation - it should predict the geometry of the disulfide once formed, but it really hasn't been optimized to do an energetic prediction with/without the disulfide. (Generally, the formation of a covalent bond between the sulfurs is favorable enough that if you have oxidizing conditions you will make a disulfide for any two cysteines in appropriate geometry, so predicting if two cysteines will/will not form a disulfide isn't all that interesting a question in most cases.)
The other issue you might run into is due to backbone rearragements due to disulfide formation. The energy in forming a disulfide bond is great enough that it can substantially torque the backbone around. If there's a great enough change in backbone conformation between the disulfide and non-disulfide forms, Rosetta may have difficulty predicting that backbone change (particularly if you use protocols which make a near-fixed-backbone assumption, as ddg_monomer does). So you may be seeing a high energy for the disulfide containing protein because the backbone is not quite at the right spot.