Inverse Kinematics and Protein Loop Closure

1. Inverse Kinematics and Protein Loop Closure Presenter: Chittaranjan Tripathy February 21, 2008 Figures are taken from the references unless otherwise stated. Special thanks to John MacMaster for allowing me to use some of his slides and figures.

2. y0 End Effector revolute joint 2 revolute joint 1 link 1 x0 Kinematics Kinematics is the branch of mechanics that studies the constrained motion of rigid objects connected by joints, without regard to forces. Degrees of Freedom The #DOF is the # of independent position variables which would have to be specified in order to locate all parts of the mechanism. Forward kinematics Where will the end effector move when we change its DOFs? i.e. F(link distances and joint angles) = the pose of the end-effector = “position and orientation” (pose) of end-effector Inverse kinematics Given the position and orientation of the end-effector, calculateall possible sets of joint angles which could be used to attain this given position and orientation. i.e., R(the pose of the end-effector) = link distances and joint angles The inverse kinematics problem for manipulators with six revolute joints has been studied for over 40 years.

3. Four Solutions of the PUMA 560 (more possible too) Inverse Kinematics The general 6R manipulator: is a robotic arm with seven rigid links connected by six revolute joints and has 6 DOF. The inverse kinematics problem for the general 6R manipulator has been studied for over 40 years. Its analytic solution was found in the late 1990’s. From: Tolani, Goswami, and Badler, Graphical Models62, 353–388 (2000)

4. algebraic geometric Solution strategies Closed form Numerical Heuristic Inverse Kinematics is Hard • Attach local frames to each joint using Denavit-Hartenberg notation. • Link transformation (moving one coordinate frame to another) can be concatenated by matrix multiplication of the type: • Where each transformation ii+1T is a function of the joint variables (combination of a translation and rotation matrix). • 0NT relates frame {N} with frame {0}. • Solving the above matrix system is hard in general. • - Can we find just one solution? • - Can we find all solutions? A three-link planar arm Differ in approach only Widely used

5. Residue #1 Residue #2 Residue #3 Residue #4 R1 H O R3 H O R5 + Ca N C’ Ca N C’ Ca N H3 C’ Ca N C’ Ca N C’ O R2 H O R4 H O Protein Backbone: A Kinematic Chain Residue #1 is the N-terminal amino acid while residue #m is the C-terminal residue. Φ: dihedral C(i-1) - N(i) - CA(i) - C(i) φ: dihedral N(i) - CA(i) - C(i) - N(i+1) C’ / N - Cα / C’ N - Cα / \ C’ C’ Φ = ? Φ = ?

6. Protein Loop Closure: An Inverse Kinematics Problem Given the positions of the BB atoms of the stationary anchor, assign values to the DOFs (Φ, φ)’s of the kinematic chain modeling the fragment so that the backbone atoms of the mobile terminusassume their target pose in the stationary anchor. Further: Which are the loops that are in good agreement with experimental data (electron density map in X-ray crystallography, NOEs and RDCs in NMR spectroscopy)

7. Why Loop Closure is an Interesting Problem? Why are Loops Important? • Loops on the surface of a protein are often flexible. • Loops play important roles in binding, recognition, and active sites of proteins and enzymes. • Often difficult to characterize by X-ray crystallography as they often introduce disorder in the protein crystal. • So we see reported structures having well-defined secondary structure elements but the loops are missing! Why Loop Closure is Interesting? • IK is interesting in its own right; extensively studied in CS, Mechanical Engineering and Robotics, and in Structural biology. • Unlike Secondary structure elements in a protein, loops do not have a stereotypical structure, therefore pose a difficult challenge to compute them efficiently (with provable guarantees) from experimental data. • Often less experimental data is available for loops, and flexibility of loops may lead to larger experimental error or data that is difficult to interpret.

8. Approaches to Solve Protein Loop Closure • Ab Initio Methods • Database (Loop Library) approach. • A Greedy Heuristic: Cyclic Coordinates Descent (CCD) Algorithm • Analytical Solution to Tripeptide (6 DOF) Loop Closure.

9. Modeling Missing Loops: Ab Initio Methods • Sample from discrete set of conformational parameters, such as from the (Φ, Ψ) map, and then refine through • Monte Carlo searches with simulated annealing. • Genetic algorithms. • Dynamic programming. • And a few others… • Robotics inspired probabilistic sampling method. (Kavraki et. al. 1996) • Sample loop conformations ignoring constraints and later enforce the constraints through gradient descent. • Need to solve an inverse kinematics (IK) problem.

10. Assemble Database Modeling Missing Loops: Database Method • Search for candidate loops that satisfy geometric constraints in homologous proteins available in a structural database (e.g. Loop Library). • Drawback: Limited loop diversity. Your search may not hit a single loop in the library! • A Hybrid approach: Assemble long missing loops from small fragments sampled from a loop library (Kolodny et. al. 2005). • Capability: Can model loops up to 15 residues long. Accuracy decreases with the length of the loop. Anchors Closed Loop

11. Hybrid Approach: Ab Initio + Database Approach: • C-space is discretized by the loops/fragments in the loop library. • Represent candidate loops as a sequence of rigid building blocks (fragments) concatenated without any DOF. • Choose the fragments from the loop library (database). Issues and Technicalities: • Combinatorial Explosion due to a huge search tree. • Grow the tree from both the ends (bi-directional search). • What is an optimal length of a fragment (coarseness of sampling)? • Typically 4/5/6. • How to determine the position of a new fragment? • Best superimpose end three Cα atoms of the already grown chain with the first three Cα atoms of the new fragments. • How to eliminate bad loops from the ensemble? • Eliminate loops that don’t agree with the experimental data. • Bad geometry, other steric and energetic parameters.

12. Unidirectional Construction for short loops • Fragments are chosen from the loop library. • Loops are closed approximately (10 A tolerance) • 3 residue overlap (ex. uses 2) between successive fragments for alignment, and 1 residue overlap to ensure final closure (C-ter). • Total #Chains • L = #fragments in the Library, • f = #res in a fragment (paper uses 5, ex. uses 4) • l = #res in the loop • 3? #res overlap • Do you see a complete L-ary tree of depth N? • - Each path from root->leaf encode a chain. • Note: Only a small fraction of N ensure loop closure. • Works well for loops with length L < 9. 2 res overlap 1 res overlap

13. Bidirectional Construction for Long loops • Fragments are chosen from the loop library. • Loops are closed approximately (10A tolerance) • 3 residue overlap between successive fragments for alignment, and 2 residue overlap at both ends of the middle fragment to ensure final closure. • Grow the loop from both ends (N-ter and C-ter), and let them meet somewhere in the middle. • Mark all positions that are the coordinates of the last two Ca atoms of first half-loop (from N-ter), and allow for end point error tolerance (10A voxel). • Enumerate all half-loops from the C-ter end and store those that have end two Ca fall in one of the previously marked voxels. These marked points are called “valid”. • Regenerate the first half-loops corresponding to the valid marked points and assemble the loops. • Reduces the time complexity by a factor of 2 in the exponent.

14. Results E. Coli TEM1 Beta-Lactamase Galactose Oxidase Modeling 8 residue long loops. Native conformation is shown as a dashed line, and the Cα trace of the top 5 template loops are shown using solid lines. The cRMS for 1BTL is [0.56-1.41A0], and for 1GOF it is [4.5-5.76A0]. Bottom Line: We get an ensemble of loops that closes the gap geometrically. How about satisfying experimental data? We rely on a generate-and-test framework.

15. The Algorithm (A Greedy Heuristic): • Generate (many) initial loop conformation(s) by sampling the values of the Degree of Freedoms (i.e. (Φ, Ψ) pairs) uniformly at random in [-π, +π]. • Fix the loop at one end • Repeat until the closure criterion is satisfied • For each DOF of loop picked in some order • Minimize closure distance for DOF Closure Distance = sqrt(Sum of squared distances of N, Cα and C atoms of final residue from their target positions) Closure criterion:Closure distance ≤ cutoff distance ε Order: • Sequential ordering of DOFs from N to C terminus. • Random permutation of DOFs Cyclic Coordinate Descent (CCD) CCD for Robot Inverse Kinematics (Wang & Chen ’91) CCD for Protein loops (Canutescu & Dunbrack ’03)

16. Cyclic Coordinate Descent… Why CCD? • Simple to implement and extremely fast! • CCD algorithm is an optimization based algorithm to solve IK problems. Here the problem is recast as a minimization problem. • Numerically stable. • (Some) External constraints on DOFs can be integrated with predictable behavior. • Linear time complexity in the number of DOFs. Drawbacks: • It is not guaranteed to return all solutions. • It may miss out a solution even if it exists. • Some of the initial fragments don’t even close while doing a CCD on them. • Not friendly towards integrating additional constraints from certain type of experimental data.

17. Fixed Target Moving C-ter (before) Moving C-ter (after) Working of CCD Algorithm Goal: Find the optimal dihedral rotation that minimizes:

18. Minimize To minimize S, set Working of CCD Algorithm Since the angle of rotation α (an extremum of S) is given by (do 2nd derivative test)

19. Multiply last two terms by: Define and Rewriting is minimum when Working of CCD Algorithm A Cute way of deriving α: Advantage of this formulation: We have sin and cos defined explicitly. Use atan2(y, x) to return θ in correct quadrant instead of doing a second derivative test.

20. A B C Results: CCD Algorithm • Two Implementations: • No constraints on the dihedrals • Bias using Ramachandran probability Map • Using the map: • For a new proposed angle φnew by CCD, propose a new ψnew using CCD. Using the map compare to (φold, ψold). • If Prob(φnew, ψnew)/Prob(φold, ψold) ≥ 1 • Accept the new position (φnew, ψnew) • Else Accept the new position (φnew, ψnew) with probability Prob(φnew, ψnew)/Prob(φold, ψold). Cα trace of the lowest RMS loop generated from 5000 trials of the CCD Ramachandran Map method for loops of 4, 8, 12 amino acids, compared with X-ray (dark) structures. (A) Loop 1EJ0A 74-77, (B) 1CTQ 144-151, and (C) 1EGU 508-519. Bottom Line: We get an ensemble of loops that closes the gap geometrically. How about satisfying experimental data? Again a generate-and-test framework!

21. Tri-peptide loop closure 16-deg polynomial in one variable Loop Conformations What is Achieved here: • Solve tri-peptide loop closure (six-torsion loop closure) analytically. • Torsion angles need not be consecutive (intervening fragments must be rigid). • Can be useful in sampling longer loops when combined with an existing loop construction algorithm. • Can be used to implement a set of local moves for Monte Carlo minimization. Tri-Peptide Loop Closure The Problem:finding the ensemble of possible backbone structures of a chain segment (with 6 DOFs) of a protein molecule that is geometrically consistent with preceding and following parts of the chain whose structures are given.

22. Cα2 Rigid Body Rigid Body C2 C1 N2 N3 Cα1 Cα3 N1 C3 Fixed in Space Tri-Peptide Loop Closure variables: τi (i=1,2,3) constraints: θi (i=1,2,3) • A possible motion involving the six dihedrals can be represented by τi (i=1,2,3). • The constraints (bond angles) θi (i=1,2,3) remain fixed => the motion is coupled.

23. Representation in Different Ref. Frames In the frame of three fixed Cα atoms. The same configurations in the original frame of fixed atoms N, Cα1, Cα3,C3.

24. Rep. in Different Ref. Frames

25. Cα2 Rigid Body C2 C1 N2 N3 Cα1 Cα3 N1 C3 Fixed in Space Choosing the Dihedrals

26. xyz Formulating Loop Closure Equations • Doing a complex derivation, we finally arrive at: • A degree 16 polynomial equation in a single variable. • Upper-bound on #solutions (= loop conformations) is 16. • The authors found 10 real solutions (at most) choosing suitable peptide torsion angles and bond angles.

27. Acknowledgements Prof. Bruce Donald John MacMaster Ed Triplett Michael Zeng

28. Thank You…

29. References • A. A. Canutescu and R. L. Dunbrack Jr. Cyclic coordinate descent: A robotics algorithm for protein loop closure. Protein Science, 12:963-972, 2003. • I. Z. Emiris, E. D. Fritzilas, and D. Manocha. Algebraic algorithms for determining structure in biological Chemistry. International Journal of Quantum Chemistry, Spec. Issue on Symbolic Methods, 2005. • E. Coutsias, C. Seok, and M. Jacobson, and K. Dill. (2004). A Kinematic View of Loop Closure. Journal of Computational Chemistry, 25, 510-528. • R. Kolodny, L. Guibas, M. Levitt and P. Koehl. Inverse kinematics in biology: the protein loop closure problem. International Journal of Robotics Research, 24, 151-163 (2005). • L. Wang, R. Mettu, and B. R. Donald. A Polynomial-Time Algorithm for De Novo Protein Backbone Structure Determination from NMR Data. Journal of Computational Biology, 13(7):1276-1288, 2006. • J. J. Craig. Introduction to Robotics: Mechanics and Control. 2nd Edition, Boston, MA: Addison-Wesley 1989, 450pp. • Shehu A, Clementi C, Kavraki LE. Modeling protein conformational ensembles: from missing loops to equilibrium fluctuations. Proteins. 2006 Oct 1;65(1):164-179.