PyRosetta (6): Movers for Backbone Folding

!pip install pyrosettacolabsetup
import pyrosettacolabsetup; pyrosettacolabsetup.install_pyrosetta()
import pyrosetta; pyrosetta.init()
from pyrosetta import *
from pyrosetta.teaching import *
init()

Introduction

One of the most basic operations in protein structure and design algorithms is manipulation of the protein conformation.

In Rosetta, these manipulations are organized into mover(s).

A Mover object simply changes the conformation of a given pose.

It can be simple, like a single φ or ψ angle change, or complex, like an entire refinement protocol.

Movers(RosettaScripts)

MoveMap class

Most Movers require a MoveMap object to specify which degrees of freedom are fixed and which are free to change.

pyrosetta.rosetta.core.kinematics.MoveMap()

Example1:

my_movemap = MoveMap()
my_movemap.set_bb(True) # Allowing a mover to change the backbone structure

Example2:

movemap.set_bb(False) # Prohibit a mover to change the backbone structure
movemap.set_bb(50, True) # Allow a mover to chage the backbone structure of residue #50
movemap.set_bb(51, True) # Allow a mover to chage the backbone structure of residue #51
# Only the backbone structure at residue #50 and #51 can be modified by a mover

Basic Backbone Folding Process

Making a trial move

Apply a move (change) to the pose (structure).

Scoring the move

Compare the energy scores before and after the move.

(Therefore, the pose should be cloned as a backup before the move.)

And, compare the energy score after the move and the lowest energy score in history.

If the score after this move is the lowest, store the _pose after this move _as the best pose, and renew the lowest energy score in history.

Deciding whether or not to accept the move

For the decision step, we need to make a subroutine that either accepts or rejects the new conformation based on the Metropolis criterion.

The difference between the energy scores after (t+1) and before (t) a move is $ \Delta E = E_{t+1} - E_{t} $.

The Metropolis criterion has a probability of accepting a move as $ P = \exp( -\Delta E / kT ) $.

When $ ΔE ≥ 0 $, the Metropolis criterion probability of accepting the move is $ P = \exp( -\Delta E / kT ) $.

When $ ΔE < 0 $, the Metropolis criterion probability of accepting the move is $ P = 1 $.

Use $ kT = 1 \text{\ Rosetta \ Energy \ Unit \ (REU)} $ .

NOTICE: Search “Boltzmann function” and “Simulated annealing” for more information about $ kT $.

If the move is accepted, the pose after the move will be used for the next round of trial move.

If not, the pose before the _move _will be used for the next round of trial move.

Iterations

For each iteration, the 3 steps above are executed.

The final output of this program should be the lowest energy conformation that is achieved at any point during the simulation.

NOTICE: Search “greedy algorithm” and “Monte Carlo” for more information about the algorithm design.

Basic backbone movers

SmallMover & ShearMover classes

Mover	Description
SmallMover	Makes “small-move-style” torsion moves (no propagation minimization)
ShearMover	Makes “shear-style” torsion moves that minimize downstream propagation

For convenience, the SmallMover and ShearMover can do multiple rounds of perturbation.

They also check that the new φ/ψ combinations are within an allowable region of the Ramachandran plot by using a Metropolis acceptance criterion based on the rama score component change. (The rama score is a statistical score from Simons et al. 1999, parametrized by bins of φ/ψ space.)

Because they use the Metropolis criterion, we must also supply kT.

kT = 1.0 # temperature
n_moves = 1 # iterations

movemap = MoveMap()
movemap.set_bb(True)

small_mover = SmallMover(movemap, kT, n_moves)
shear_mover = ShearMover(movemap, kT, n_moves)

BackrubMover

Mover	Description
BackrubMover	Makes local rotations around two backbone atoms

Minimization Mover MinMoverclass

The MinMover carries out a gradient-based minimization to find the nearest local minimum in the energy function, such as that used in one step of the Monte-Carlo-plus-Minimization algorithm of Li & Scheraga.

min_mover = MinMover()

MonteCarlo class

Introduction

The MonteCarlo object is an encapsulated object that creates a whole MonteCarlo simulation.

That is, it can decide whether to accept or reject a trial conformation, and it keeps track of the lowest-energy conformation and other statistics about the search.

Having the Monte Carlo operations packaged together is convenient, especially if we want multiple Monte Carlo loops to nest within each other or to operate on different parts of the protein.

# Assume that we have instantiated a mover `my_mover`
# Instantiate a `MonteCarlo` object `mc`
mc = MonteCarlo(my_pose, scorefxn, kT)
# my_pose: a Pose object
# scorefxn: a ScoreFunction object
# kT: a parameter - temperature; int or float 

# apply the mover
my_mover.apply(my_pose) 

# deside whether or not to accept the move
flag = mc.boltzmann(my_pose)
# flag is boolean
# True: accept
# False: reject

# In practice, these 3 steps should be written into a loop.

A Monte Carlo Simulation process is (in pseudo code):

def monte_carlo_sampling(initial_state, kT, N_iteration):
    current_state = initial state
    for i in range(N_iteration):
        new_state = random_move(current_state)
        #fe_delta: delta (difference of) free energies
        fe_delta = free_energy(new_state)-free_energy(current_state) 
        if fe_delta < 0 or uniform(1) < exp(-fe_delta / kT): 
            # exp(-fe_delta / kT): boltzmann distribution
            current_state = new_state

MonteCarlo Methods

mc.boltzmann(my_pose) # See above
mc.show_scores()
mc.show_counters()
mc.show_state()

TrialMoverclass

A TrialMover combines a specified Mover with a MonteCarlo object.

Each time a TrialMover is called, it performs a trial move and tests that move’ s acceptance with the `MonteCarlo object.

It is designed to test the effects of a Mover.

trial_mover = TrialMover(small_mover, mc)

for i in range(10):
    trial_mover.apply(test)
    
print(trial_mover.num_accepts())
print(trial_mover.acceptance_rate())

# After the trial, information about the trial can also be visited by
# mc.show_state()
# so that different movers can be compared

SequenceMoverclass and RepeatMoverclass

seq_mover = SequenceMover()
seq_mover.add_mover(small_mover)
seq_mover.add_mover(shear_mover)
seq_mover.add_mover(min_mover)

n_repeats = 3
repeat_mover = RepeatMover(my_mover, n_repeats)