pyrosetta.rosetta.ObjexxFCL is a Python module wrapped by PyBind11 around a C++ wrapper for Fortran, so I expect this matryoshka to be a problem. The documentation says these are being slowly phased out and it is uncommon finding a method that requires a ObjexxFCL FArray as opposed to a vector. However, I wanted to use one...
I was looking into the class DecoySetEvaluation in toolbox so I wanted an initialised pyrosetta.rosetta.ObjexxFCL.FArray1_double_t or pyrosetta.rosetta.ObjexxFCL.FArray2_double_t. But no constructor is defined and there seem to be no class methods. Here is a full example:
# init
import pyrosetta
pyrosetta.init(extra_options='-no_optH false -mute all -ignore_unrecognized_res true -load_PDB_components false')
# load a multimodel PDB (& split it) as a proxy for a MCMC mover .dump_poses output
original = pyrosetta.toolbox.rcsb.pose_from_rcsb('1L2Y')
mini = pyrosetta.rosetta.protocols.grafting.return_region(original, 1,20)
n = [pyrosetta.rosetta.protocols.grafting.return_region(original, 1+i,20+i) for i in range(0, original.total_residue(),20)]
# n is a list. Instancemethods use vectors
poses = pyrosetta.rosetta.utility.vector1_core_pose_Pose()
altposes = pyrosetta.rosetta.utility.vector1_std_shared_ptr_const_core_pose_Pose_t()
poses.extend(n) # used by a few methods —such as MCMC movers .dump_poses()
altposes.extend(n) # for: pyrosetta.rosetta.core.io.pdb.dump_multimodel_pdb(altposes, 'test.pdb')
# analyse
dse = pyrosetta.rosetta.protocols.toolbox.DecoySetEvaluation()
dse.reserve(len(poses))
for pose in poses:
dse.push_back(pose)
# this works fine:
r = pyrosetta.rosetta.utility.vector1_double()
dse.rmsf(r)
print(r)
# this requires an FArray1
weights = xxx(pyrosetta.rosetta.ObjexxFCL.FArray1_double_t, dimension=len(poses))
xxx.fill(weights, 1)
dse.center_all(weights)
In Fortran an array is declared with dimensions
real, dimension(5) :: myvarname
There is a Dimension class in the ObjexxFCL module that can be instantiated but it does nothing that I can see, but might be a way in.
So How does one instantiate a FArray?
If it is not possible to instantiate a FArray, are there any methods that return one that could be co-opted?
And also, is it possible to convert a FArray2 to pose?
# this requires a FArray2, but would seem to fill it as a pose that is not a pose
notpose = xxx(pyrosetta.rosetta.ObjexxFCL.FArray2_double_t,...)
dse.compute_average_structure()
pose = xxx.FArray2_to_pose(notpose)
This is not an insurmountable problem for this example as I have been using some GROMACS utilities that compute the RSMF and the average structure. But it would be nice to know as I have seen ObjexxFCL arrays requirements in a few places...
I re-encountered this problem just now —weights for something else. Even though I do not know what spits out an instance of an ObjexxFCL, I had figured out how to get the RMSF working (root-mean square fluctuation aka. RMSD from averaged decoy).
But actually, RMSF does not require an averaged decoy. It is just that DecoySetEvaluation has the strangest way of initialising.
The RMSF is:
dse.ref_pose is the first pose, not the average pose, which is not accessible, but not important as the RMSF seems fine.
I've never heard of or used the DecoySetEvaluation. I would just use SimpleMetrics; RunSimpleMetrics mover to get RMSD, etc.
Yes, I normally use
pyrosetta.rosetta.core.scoring.CA_rmsd& co. for RMSD of a pose vs. reference. But given an ensemble (with no reference) I would want it against the averaged positions —either per residue or per pose/decoy in the ensemble (above was per residue). That is not possible with the SimpleMetric classes, right?I simply searched RMSF in the documentation and found DecoySetEvaluation. However, I really should admit that I have not actually really used DecoySetEvaluation for anything as outputting as PDB file or string and reloading the data in ProDy is way easier and I trust it more —even if it sounds like a bonkers thing to do. I mean say the above is:
Same goes for clustering.
My original problem was instantiating one of the ObjexxFCL object and I reincountered it for weights for something (forgot already) —but I am guessing aren't actually meant to be used pythonically anyway. The RMSF thing was solely parenthetical as going through ProDy is fine.