Skip to main content
Bounty Ended with 25 reputation awarded by CommunityBot

For modelling, I wrote some code for setting spans with the archetypical angles for different SS: herehere. I do not recall whence the data came (maybe Wikipedia?), but it may be handy anyway.

For modelling, I wrote some code for setting spans with the archetypical angles for different SS: here. I do not recall whence the data came (maybe Wikipedia?), but it may be handy anyway.

For modelling, I wrote some code for setting spans with the archetypical angles for different SS: here. I do not recall whence the data came (maybe Wikipedia?), but it may be handy anyway.

Source Link
Matteo Ferla
  • 4.3k
  • 6
  • 19

I do not know if there is a function in Biopython.PDB for dihedral space transformations, but due to the lack of answers I would say no. Biopython.PDB is nice and clean, but I would use something with more molecular mechanics power... or should I say potential (:drum::laughing:).

Personally, I use PyRosetta, but it is not open source (code available to academics), so others may opt for OpenForceField. PyRosetta requires a licence (free for academia), which is different from the Rosetta proper licence.

The code below uses a helper module I wrote, so I'll use the installer therein (which is shorthand for the official pip wheel, which provides different options etc., while the official shorthand is actually pyrosettacolabsetup)

pip3 install pyrosetta-help
install_pyrosetta -u 👾👾👾 -p 👾👾👾

Ligands?

First lets make sure there's nothing interesting ligand wise.

NB. RCSB has a new shiny API which may be better... Not tried it: currently using European.

import requests

def retrieve_nonpolymer_data(entry_id: str) -> dict:
    """cf. https://www.ebi.ac.uk/pdbe/api/doc/"""
    response: requests.Response = requests.get(f' https://www.ebi.ac.uk/pdbe/api/pdb/entry/molecules/{entry_id.lower()}')
    response.raise_for_status()
    return response.json()[entry_id.lower()]

info = retrieve_nonpolymer_data('1enh')
chem_resns = [c for entity in info for c in entity.get('chem_comp_ids', [])]
print(chem_resns)  # ['HOH']

If is were something that is not water, I made a list here, which is based off Peter Curran's repo for the HotSpots algorithm, which allows the elimination of worthless ligands like glycol.

I mention this because the route differs if there is a ligand.

Init and read

import pyrosetta
import pyrosetta_help as ph
from types import ModuleType

logger = ph.configure_logger()
pyrosetta.init(extra_options=ph.make_option_string(no_optH=False,
                                                ex1=None,
                                                ex2=None,
                                                #mute='all',
                                                ignore_unrecognized_res=True,
                                                load_PDB_components=False,
                                                ignore_waters=True)
                               )

# shorthand for sanity
prc: ModuleType = pyrosetta.rosetta.core

# this route strips ligands
pose: Pyrosetta.Pose = pyrosetta.toolbox.pose_from_rcsb('1enh')

Note ignore_unrecognized_res, load_PDB_components and ignore_waters.

Angles

This will not happen to an non-Rosetta user, but I get tripped out by the fact that ψ & φ data is not stored in residue, but in the pose.

for residue in pose.residues:  #: prc.conformation.Residue
    # residue.get_phi() # ---> Error!
    pass

For the residue index in pose, Rosetta was written in Fortran so is one-indexed.

import numpy as np
import numpy.typing as npt
phi:npt.NDArray[float] = np.fromiter(map(pose.phi, range(1, pose.total_residue()+1)), dtype=float)
psi:npt.NDArray[float] = np.fromiter(map(pose.psi, range(1, pose.total_residue()+1)), dtype=float)

Change φ/ψ angles

Time for some fun in dihedral space, say:

r = 2
phi = 180
pose.set_phi(r, phi)

For modelling, I wrote some code for setting spans with the archetypical angles for different SS: here. I do not recall whence the data came (maybe Wikipedia?), but it may be handy anyway.

Save

To save a PDB file:

pose.dump_pdb('👾👾👾.pdb')

There is also pose.dump_scored_pdb, but that is only for friends of J. Willard Gibbs.

Minimisation

It is said in the comments that energy minimisation is not required. But here is a snippet that minimises the sidechains only

scorefxn = pyrosetta.get_fa_scorefxn()  # ref2015 will do.
before: float = scorefxn(pose)
cycles = 5
movemap = pyrosetta.MoveMap()
movemap.set_bb(False)  # backbone
movemap.set_chi(True)  # sidechain
relax = pyrosetta.rosetta.protocols.relax.FastRelax(scorefxn, cycles)
relax.set_movemap(movemap)
relax.apply(pose)
after: float = scorefxn(pose)
print(after - before)  # Willard says "Remember free energy is a potential so it should be negative"

If the BB tweaks are minimal, there's also the option of scoring/minimising against the electron density etc.