Hello all,

I'm interested in applying RMSD constraints to a pose using pyrosetta.

I currenly have a native pose which is **deformed manually**, such the the backbone RMSD with respect to the native is 1 Angstroms.

`pose.residue(r).set_xyz(a, py.rosetta.numeric.xyzVector_double_t(new_coords[0], new_coords[1], new_coords[2]))`

However, when I relax (cartesian), I notice the structure relaxing back to the native, ie the RMSD decreases to ~0 Angstroms.

Im keen to apply constraints on the **deformed coordinates before relaxing** such that the RMSD of 1 Angstroms is maintained.

Currently, I have been applying coordinate constraints on the backbone of the deformed structure using a harmonic potential with a std dev of 0.5 Angstroms. But this is not ideal, as the backbone RMSD can still change.

There are parallels to this, for example using PLUMED for molecular dynamics simulations. I wish to know if I can implement this somehow in Rosetta? I would be really grateful to any insights. Feel free to contact me for any details/data.

Is there a general solution to the application of constraints to much more complex collective variables (CVs) in Rosetta?

Thanks,

Akshay

First off, a bit of terminology clarification. When Rosetta says "constraints", it's referring to what is normally termed "restraints" in molecular dynamics parlance -- that is, an adjustment to the score function which biases conformations according to a particular specified functional form. Rosetta doesn't particularly have a specific term for what MD typically calls "constraints" (that is, turning off sampling for particular degrees of freedom). Because Rosetta generally uses Monte Carlo sampling, if you don't want to sample a particular degree of freedom ... you just don't. (Though it's sometimes a bit of effort to get a sampling scheme that samples the particular DOFs you are interested in.)

The other thing to keep in mind is that, as an adjustment to the energy function, there's actually two parts to a Rosetta constraint. There's the specification of the functional form itself (the constraint file or equivalent), which gets added to the Pose, and then there's the constraint score terms which get enabled in the ScoreFunction. If you don't have the constraints enabled in the scorefunction you're using for sampling, the constraints won't affect the sampling. -- This is used for various protocols to vary how heavily constraints are obeyed, without having to re-define them in the Pose. Indeed, the FastRelax protocol will ramp the constraint term in the ScoreFunction during its multiple stages unless you set the options to turn that off. (This could be an issue you're seeing. If you didn't explicitly turn off constraint ramping, the final structure will be optimized without constraints.)

Theoretically, the Rosetta constraint system can address any particular functional form you desire, so long as it can be expressed relative to atomic coordinates. That said, in practice doing more complex functions would require adding additional C++ code and making a new constraint type. (There may be a way to write a new constraint term with PyRosetta, though hooking it into the constraint file reader might be a trick.) Much easier is standard sums of individual atom pair, coordinate, angle and dihedral constraints.

The upshot is that if you want to keep the backbone from moving at all, you probably want to look into MoveMaps. The FastRelax protocol only changes the backbone during minimization, and MoveMaps are how minimization determine which degrees of freedom get to be minimized. It should be easy to turn off all backbone DOFs and only allow minimization of the sidechains. Alternatively, you can just skip FastRelax and use a protocol which is specific to sidechain degrees of freedom. The MinPacker is probably what I would recommed. It doesn't use the same protocol as FastRelax, but if you have a fixed backbone it should give you relatively equivalent results. (It packs and independently minimizes each sidechain rotamer in context. It doesn't do any ramping repulsive, but if you have a fixed backbone you shouldn't need the "breathing" motion that enables to optimize and collapse things.)

If you want to allow some backbone movement, then FastRelax with constraints is a decent way to go. However, since they're restraints instead of MD constraints, you may need to play around with the strength of the constraints in order to find a setting which keeps things in the range you want. You can do this either by adjusting the coordinate_constraint scoreterm in the score function, or you can adjust the sd term of the constraints themselves (sd is measured in Angstrom - smaller values are stronger constraints).

There isn't an RMSD specific constraint, but if you want to hit an RMSD of 0, then the sum of harmonic coordinate constraints to your reference structure is equivalent. (As a square root of zero means a zero inside the radical, which implies the sum of the squared deviations is zero, which is equivalent to a sum of squared/harmonic coodinate constraints set to zero.) RMSD constraints to non-zero values get to be a bit harder, as in order to use them in minimization you would need to also define the gradient of the overall scoreterm with respect to each of the individual DOFs. Probably doable, but it would take some coding, and hasn't been anything people have necessarily been interested in. -- Which prompts me to encourage you to consider what it is you *really* want to model, biochemically speaking. Rosetta and MD take slightly different approaches to how they represent and treat problems, so in certain situations you could be taking an MD hammer to a Rosetta "helical ring-shank nail" (a.k.a screw).

Hello,

Thanks a lot for your elaborate post. I broadly understand the points you've made. I've created a custom relax script without ramping down the coordinate constraints and have benchmarked the process. Plus I've tweaked the movemap movers accordingly. All of this seems to work for me as of now. I needed to correct the bond lengths and angles as I had manually deformed the positions, so I did not freeze any backbone atom positions but instead used heavier constraint penalties (and std dev) than usual.

In my use case, Im deforming the atoms based on the direction of a set of Normal modes whose magnitude is determined by the RMSD. In reality, I needed to apply constraints on the RMSD projected along the Normal mode vector, which is an added layer of complexity in defining the constraints (in pyrosetta). Another approach Im considering is applying one-sided constraints (for each backbone x,y,z directions) during relaxation such that the direction of allowed change is aligned to the normal mode deformation. Im guessing this should do the job for now.

Thanks again,

Akshay,