Flex ddG Tutorial

The most up-to-date version of this tutorial is available on GitHub.

Introduction

In this activity, you will utilize the Flex ddG [KB2018] protocol within Rosetta to computationally model and predict changes in binding free energies upon mutation (interface ΔΔG).

This protocol uses the "backrub" protocol [CS2008] implemented in Rosetta to sample conformational diversity.

This tutorial was originally written against weekly release "Rosetta 2017.52", released on Wednesday, January 3, 2018, and that release is still the reference version for the published benchmark numbers. More recent releases also work: the tutorial has been verified to run unmodified on Rosetta 2022.45 (release-334). If you use a newer Rosetta, note that predicted ΔΔG values will not reproduce the 2017.52 numbers to the last decimal place, since the underlying score function implementation has continued to change.

The flex ddG protocol is outlined below (Fig 1. from [KB2018]):

Example runs

Example 1: Run Flex ddG on a specific set of mutations

1. From within your downloaded copy of this tutorial, open run_example_1.py in your editor of choice.

2. Find the rosetta_scripts_path at the top of run_example_1.py and check that it is set to the appropriate location of your compiled Rosetta rosetta_scripts binary.

3. Run python run_example_1.py. The full command line call to each instance of Rosetta will be displayed, and will look something like this:

   /home/user/rosetta/source/bin/rosetta_scripts -s /home/user/flex_ddG_tutorial/inputs/1JTG/1JTG_AB.pdb -parser:protocol /home/user/flex_ddG_tutorial/ddG-backrub.xml -parser:script_vars chainstomove=B mutate_resfile_relpath=/home/user/flex_ddG_tutorial/inputs/1JTG/nataa_mutations.resfile number_backrub_trials=10 max_minimization_iter=5 abs_score_convergence_thresh=200.0 backrub_trajectory_stride=5 -restore_talaris_behavior -in:file:fullatom -ignore_unrecognized_res -ignore_zero_occupancy false -ex1 -ex2

   Important flags explained:

   • -ex1 -ex2 tell Rosetta's side chain packing algorithm to sample extra subrotamers for chi1 and chi2 angles (Packer documentation)
   • mutate_resfile_relpath= is an input file that tells Rosetta which protein positions to mutate in the ΔΔG calculation. The resfile must start with NATAA. See the Rosetta documentation for more information on resfiles.
   • number_backrub_trials= is the number of backrub sampling steps. 35,000 steps gives good performance for a variety of inputs, although it is likely similar performance could be achieved for certain cases with less sampling.
   • max_minimization_iter= is the maximum number of minimization gradient descent steps to take. 5000 is the normal, benchmarked value.
   • abs_score_convergence_thresh= is the maximum allowed change in total model score after minimization (comparing initial score to score after minimization). If this change exceeds this threshold, then another minimization cycle will be started.
   • backrub_trajectory_stride= after every N backrub steps, finish the flex ddG calculations. This allows protocol performance to be judged at intermediate number of backrub sampling steps without running a separate trajectory. In general, this can be set to the same number as number_backrub_trials, but it may be useful for benchmarking purposes to set this at an intermediate value.
   • chainstomove= This argument is passed as a Rosetta Script variable to the InterfaceDdGMover as the "chain_name" argument. It defines one side of the interface by specifying all of the chains that make up that side. Can be a single chain as in "B" or multiple chains, as in "L,H".

4. Output will be saved in a new directory named output

Example 2: Run Flex ddG for single site saturation mutagenesis

This example covers the commonly desired use case of evaluating the energies of all possible mutations at a single residue site in the interface. To do this in Rosetta, it is necessary to create a resfile for each possible amino acid mutation, and run the flex ddG protocol with each of these resfile as inputs. In this example, run_example_2_saturation.py is a modified version of the first example script that has been modified to automatically create resfiles for all 20 possible canonical amino acid mutations, and then run flex ddG on those resfiles. You can also create the resfiles yourself manually before running the protocol.

1. From within your downloaded copy of this tutorial, open run_example_2_saturation.py in your editor of choice.
2. Find the rosetta_scripts_path at the top of run_example_2_saturation.py and check that it is set to the appropriate location of your compiled Rosetta rosetta_scripts binary.
3. Run python run_example_2_saturation.py. The full command line call to each instance of Rosetta will be displayed.
4. Output will be saved in a new directory named output_saturation

Analysis

In normal usage, you would run the flex ddG protocol 35+ times (at 35,000 backrub steps each run), and average the resulting ΔΔG predictions for best performance. For the purposes of making this tutorial run quickly on an average laptop, we will generate fewer output models for many fewer backrub and minimization steps.

Python analysis

These Python packages are required in order to run the analysis, and can be installed via pip: pip install numpy pandas.

Run the analysis script for example 1 as follows:

python analyze_flex_ddG.py output

Or for example 2:

python analyze_flex_ddG.py output_saturation

The script will print to the terminal (in separate table blocks) the wild type interface binding ΔG score (wt_dG), the mutant interface ΔG (mut_dG), and the ΔΔG of binding post-mutation. These scores are also written to a .csv file in analysis_output. Scores are reported separately for each checkpoint of the backrub trajectory (with the tutorial's default settings, 5 backrub steps and 10 backrub steps). For the mutant ΔΔG, the ΔΔG score is also calculated and reweighted with the fitted GAM model [KB2018].

The backrub_steps column is labelled using the backrub_trajectory_stride that each run was launched with, which the script reads back out of the run's own ddG.db3. You therefore do not need to edit the analysis script to match your run. If a database does not record the stride, the script warns and falls back to a default of 5; pass --stride N to override it. The stride affects only the checkpoint labels, never any energy.

Extract structures

If you are interested in viewing or using the generated backrub, wildtype minimized, or mutant minimized structures, you can extract them from the struct.db3 file in the output. A convenience wrapper script is provided to do this, and can be run as follows:

python3 extract_structures.py output

The script will recursively find all output struct.db3 files, run Rosetta to output PDBs, and rename the PDBs to more informative names.

Note that extraction uses the score_jd2 binary rather than rosetta_scripts. This is a separate executable, built alongside rosetta_scripts by the standard Rosetta build. Check the score_jd2_path variable near the top of extract_structures.py and set it to the location of your compiled binary. As above, the extracted PDBs are named using the stride read back out of each struct.db3, and --stride N overrides it.

Additional reading

The Rosetta documentation wiki can provide additional context for how to adapt this Rosetta Scripts protocol to your specific use case. In particular, the page on RosettaScripts and the section of that page that explains XML variable substitution might prove helpful.

References

[KB2018]

(1, 2, 3) Kyle A. Barlow, Shane Ó Conchúir, Samuel Thompson, Pooja Suresh, James E. Lucas, Markus Heinonen, and Tanja Kortemme. Flex ddG: Rosetta Ensemble-Based Estimation of Changes in Protein–Protein Binding Affinity upon Mutation. J. Phys. Chem. B, February 2018. doi: 10.1021/acs.jpcb.7b11367. URL: https://pubs.acs.org/doi/pdf/10.1021/acs.jpcb.7b11367

[KB2017]

Kyle A. Barlow, Shane Ó Conchúir, Samuel Thompson, Pooja Suresh, James E. Lucas, Markus Heinonen, and Tanja Kortemme. Flex ddG: Rosetta Ensemble-Based Estimation of Changes in Protein–Protein Binding Affinity upon Mutation. bioRxiv Preprint, November 2017. URL: https://www.biorxiv.org/content/early/2017/11/17/221689

[CS2008]

Smith, C. A.; Kortemme, T. Backrub-Like Backbone Simulation Recapitulates Natural Protein Conformational Variability and Improves Mutant Side-Chain Prediction. Journal of Molecular Biology 2008. DOI:10.1016/j.jmb.2008.05.023.