How to create an ObjexxFCL FArray?

[Anonymous]

    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…

[Anonymous]

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.

# given a list of poses or vector1_core_pose_Pose

# trp-cage miniprotein (38 models)

original = pyrosetta.toolbox.rcsb.pose_from_rcsb('1L2Y')

mini = pyrosetta.rosetta.protocols.grafting.return_region(original, 1,20)

nmr = pyrosetta.rosetta.utility.vector1_core_pose_Pose()

n = [pyrosetta.rosetta.protocols.grafting.return_region(original, 1+i,20+i) for i in range(0, original.total_residue(),20)]

nmr.extend(n)



# init the decoyset

dse = pyrosetta.rosetta.protocols.toolbox.DecoySetEvaluation()  # no arguments accepted

dse.reserve(len(nmr))

for pose in nmr:

dse.push_back(pose)

print(len(nmr), dse.n_decoys())

# get the rmsf

results = pyrosetta.rosetta.utility.vector1_double()

dse.rmsf(results)

print(results)

The RMSF is: 

vector1_double[1.61104, 0.384238, 0.193109, 0.191801, 0.237876, 0.153467, 0.128347, 0.248295, 0.319172, 0.212154, 0.257973, 0.296829, 0.305155, 0.424018, 0.561687, 0.331984, 0.275279, 0.324276, 0.325702, 1.13208]

dse.ref_pose is the first pose, not the average pose, which is not accessible, but not important as the RMSF seems fine.

ref = dse.ref_pose()

pyrosetta.rosetta.core.scoring.all_atom_rmsd(ref, nmr[1])

 

[Anonymous]

I’ve never heard of or used the DecoySetEvaluation.  I would just use SimpleMetrics; RunSimpleMetrics mover to get RMSD, etc. 

[Anonymous]

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:

import prody as pdy

ensemble = pdy.Ensemble('trpcage')

ag = pdy.parsePDB('1L2Y.pdb')

ensemble.setCoords( ag.getCoords() )

ensemble.addCoordset( ag.getCoordsets() )

ensemble.iterpose() # ref_coords is now the averaged.

ensemble.getRMSFs() # this is per atom though and needs reshaping

ensemble.getRMSDs() # RMSDs is per conformation against the average

# etc.

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.

[Author]

Posts