sidechain
- Sidechain ModellingΒΆ
Tools and algorithms to model sidechains given backbone coordinates. The full module is heavily based on SCWRL4 [krivov2009] . The according paper describes the modelling of sidechains using two different rotamer models. A rigid model, as well as a flexible model. Both models are implemented in ProMod3 and can be applied in flexible ways.
The following code fragment shows an example of a basic sidechain reconstruction algorithm using the functionality in the module. Note, that this code will crash as soon as you have structures containing all the weirdness the PDB throws at us. In this case, you should use the full fletched sidechain reconstruction pipelines available in the modelling module.
from ost import io,mol
from promod3 import sidechain
# load a protein
prot = io.LoadPDB('data/1CRN.pdb')
# load rotamer library
library = sidechain.LoadBBDepLib()
# we need a rotamer constructor to create any rotamers or
# frame residues
rot_constructor = sidechain.SCWRL4RotamerConstructor(False)
# create new entity from protein only containing the amino acids
prot = mol.CreateEntityFromView(prot.Select("peptide=true"), True)
# gather some data, the rotamer ids and backbone torsion angles
torsion_angles = list()
rotamer_ids = list()
for r in prot.residues:
rotamer_ids.append(sidechain.TLCToRotID(r.GetName()))
phi_handle = r.GetPhiTorsion()
psi_handle = r.GetPsiTorsion()
# set typical default values for an alpha helix...
phi = -1.0472
psi = -0.7854
if phi_handle.IsValid():
phi = phi_handle.GetAngle()
if psi_handle.IsValid():
psi = psi_handle.GetAngle()
torsion_angles.append((phi,psi))
# first build a frame representing the rigid parts including
# cystein sidechains
frame_residues = list()
for i,r in enumerate(prot.residues):
frame_residue = rot_constructor.ConstructBackboneFrameResidue(
r, rotamer_ids[i], i,
torsion_angles[i][0],
torsion_angles[i][1],
i == 0,
i == (len(rotamer_ids)-1))
frame_residues.append(frame_residue)
frame = sidechain.Frame(frame_residues)
# let's build up rotamer groups
rotamer_groups = list()
aa_with_rotamers = list()
for i,r in enumerate(prot.residues):
if r.GetName() == "ALA" or r.GetName() == "GLY":
continue
rot_group = rot_constructor.ConstructFRMRotamerGroup(
r, rotamer_ids[i], i, library,
torsion_angles[i][0], torsion_angles[i][1],
i == 0, i == (len(rotamer_ids)-1), 0.98)
# calculate pairwise energies towards the rigid frame
# the values will be stored and used by the solving algorithm
# later on when considering self energies
rot_group.SetFrameEnergy(frame)
# remove super unlikely rotamer in rotamer group
# e.g. those who clash with the frame
rot_group.ApplySelfEnergyThresh()
# finally add it and keep track of the indices
rotamer_groups.append(rot_group)
aa_with_rotamers.append(i)
# buildup a graph given the rotamer groups
graph = sidechain.RotamerGraph.CreateFromFRMList(rotamer_groups)
# and get a solution out of it using a minimal width tree approach
solution = graph.TreeSolve()[0]
# let's finally apply the solution to the residues
for (i, rot_group, sol) in zip(aa_with_rotamers,
rotamer_groups,
solution):
rot_group[sol].ApplyOnResidue(prot.residues[i])
io.SavePDB(prot, "example_reconstruction.pdb")
Contents: