Prosculpt

Protein design and sculpting using Rosetta and Deep learning methods (RFDiff and Alphafold2)

Description

The script prosculpt_run.py runs a pipeline to automate the processes of generating protein structures with RFdiffusion, sequence generation with ProteinMPNN, and folding and evaluation with AF2 and Rosetta. Specifically, the script uses motif scaffolding to generate structures (see RFdiffusion github repository).

The main steps are as follows:

1. The main inputs are a protein PDB file and a yaml file that specifies how to generate the new structure around the original protein (which parts to keep, etc.).
2. The RFdiffusion module generates new protein structures (backbone atoms only) based on the options provided.
3. Generated PDB files are preprocessed using helper scripts and ProteinMPNN generates sequences (aminoacid residues) for these structures.
4. Sequences are prepared for the AF2, which folds them into structures.
5. Additional evaluation parameters are calculated for these structures using a separate script with Rosetta.
6. The final output is a CSV file containing AF2 structure path, evaluation parameters etc.

A working Collab can be found here

Requirements

The script requires the prosculpt package. It assumes that RFdiffusion, proteinMPNN, and AF2 are installed and that the correct paths are provided in the installation.yaml config file. Additionall biopython, hydra-core, pandas and scipy and pyRosetta are required.

For running tests, Python version must be >= 3.7 (it needs the capture_output arg).

Installation

The easiest way to install prosculpt is using a container manager like apptainer (or singularity). To do so, follow these instructions. Installation time is quick, but SIF files can take several minutes to download.

Begin by cloning prosculpt:

git clone https://github.com/ajasja/prosculpt.git
cd prosculpt

Create a new environment called prosculpt (or any name of your choice) using python or any venv manager and activate it:

python -m venv prosculpt_venv/
source prosculpt_venv/bin/activate

Install it and its dependencies

pip install .

Install pyrosetta into the environment

pip install pyrosetta --find-links https://west.rosettacommons.org/pyrosetta/quarterly/release.cxx11thread.serialization

At a location of your choice download the appropriate SIF files, executing the following commands:

singularity pull rfdiff.sif docker://rosettacommons/rfdiffusion
singularity pull proteinmpnn.sif docker://rosettacommons/proteinmpnn
singularity pull pymol.sif docker://jysgro/pymol:3.1.0_amd_arm

Depending on the prediction model, follow the appropriate instructions (Multiple instructions can be followed and the prediction model selected as desired in the configuration)

For using Boltz2, install boltz2 into the environment

pip install boltz[cuda] -U

it might be necessary to install compilers for boltz2 to run correctly. The exact command for this is system-dependent. If using conda, it can be done by

conda install -c conda-forge compilers

For Alphafold 3 use, request the model parameters following the instructions on the Alphafold 3 github. Then download the sif file using

singularity pull af3_2.sif docker://kosinskilab/alphafold3

For using Colabfold (Note: Colabfold is slower and less convenient than using Boltz2 or AF3 and is not recommended)

singularity pull colabfold.sif docker://unionbio/colabfold:w4KMVR7WrKDlCbdQ1BYrjQ-test

Copy or rename the file config/installation.yaml.singularity_template as config/installation.yaml. Then replace the default values with the ones corresponding to your system. Mainly:

• Replace the run command of the sif files with the appropriate commands and paths for your system.
• Replace the prosculpt_python_path with the python path of your prosculpt env (You can get it by running "which python" with the env active)
• Replace the slurm options with those that correspond with your cluster (these can then be overwritten by individual job settings)
• If multiple users will be using the same installation, it is good to set the boltz2 cache parameter to the current user so that only one set of parameters are downloaded.

If using Boltz2 or Colabfold it is necessary to run a job to allow them to download the parameters. In a slurm cluster you can run:

Colabfold:

python run_tests.py -t unconditional

Boltz2:

python run_tests.py -t unconditional_boltz

Usage

To run Prosculpt on Slurm, use the slurm_runner.py passing it an input yaml file with the job options. Note that slurm_runner no longer supports passing arguments inline with the exception of ++output_folder.

python slurm_runner.py job_parameters.yaml

By default, slurm's own output/error log files (sbatch's -o/-e) are forced into a logs/ subdirectory of the job's output_dir, regardless of any slurm.output/slurm.error set in the job yaml or in config/installation.yaml - so a job's slurm log always ends up next to its own results instead of wherever a job/site config happened to point it. Pass --allow-custom-log-path to opt back into respecting slurm.output/slurm.error instead.

python slurm_runner.py job_parameters.yaml --allow-custom-log-path

If you want to run Prosculpt locally, it can be called as

python prosculpt_run.py -cd config_directory -cn config_file_name

(note that the config file name must be written without the yaml extension)

Examples

For usage examples see the /Examples or /Minimal_Prosculpt_script folders

Skipping RfDiff

Sometimes you may want to skip RfDiff and only pass the pdb to ProteinMPNN and AF with a few residues to change.

Add skipRfDiff: True and designable_residues: [A8, A9, A10, A13, A85] to the yaml config, where designable_residues contains the residues you want to change. If there are non-redesigned chains in the model that you want included in the final model, add the letter chain (without any numbers) to the designable residues.

In order to actually produce diverse results, --sampling_temp and --backbone_noise passed to ProteinMPNN are increased to 0.3 and 1, respectively.

Working with natural proteins

To be modelled correctly, natural proteins require structure prediction softwares to work with multiple sequence alignments. In order to obtain good results when working with natural proteins, the option "use_a3m" needs to be set to True. This option requires the user to input an MSA file for each chain in the reference structure for which the user wants alignment information to be used. To obtain these MSAs just run each chain in the reference structure separately through an MSA server and use the A3M file it produces. The easiest way is to run the sequence for your input structure through Alphafold. Each file needs to contain within the name "Chain_X" where X is the chain ID of the corresponding chain in the reference structure file. These alignment files need to be placed together in a directory and the "a3m_dir" option needs to be set to this directory in the yaml config file or passed through command line. For clarification, look at the binder example in the Examples directory.

Prosculpt generates new "partial alignment" files for each design where natural proteins are modelled using MSAs and designed proteins (or designed portions of hybrid proteins) are not.

Backbone filtering

For backbone filtering, custom scripts called plugins are used. A plugin is a python script containing a function called "filter_backbone" which takes a pdb file and any set of arguments and returns true if the pdb passes the filter and false if it doesn't. The content of the function can be defined as needed by the user. To call filters add the following code to your configuration yaml. Examples can be found in the plugins folder and usage example can be found in the examples folder

rfdiff_backbone_filters:
    - filter_name: [filter_name]
      filter_script: [filter script file path]
      delete_failed: [true/false]

Metrics and scoring script

Prosculpt calculates the following metrics

Metric	Meaning
RMSD	Full model RMSD to RFDiff backbone
pLDDT	full model pLDDT
RMSD_sculpted	RMSD of any residues not present in original structure to RFDiff backbone
RMSD_fixed_chains	RMSD of residues in non-redesigned chains to RFDiff backbone
RMSD_motif	RMSD of non-redesigned residues in chains containing redesigned residues
pae	Predicted alignment error

In addition, a scoring script can be added to calculate extra metrics. The scoring script is a python file called "scoring_script.py" in the scripts folder that outputs a CSV called "rosetta_scores.csv", with a "model_path" column. This CSV file gets merged with the prosculpt metrics for the final output. A default scoring script is included in the scripts folder, which calculates radius of gyration, charge and sap.

Inputs:

Most input parameters are documented in the examples, as well as the run.yaml config file. However, here's some additional info about them:

• contig: explained in detail in the RFdiffusion github repository.
  • Please note: The whole argument must be enclosed in double quotes and square brackets!
• output_dir: output directory. Along the pipeline, each module will create its own subdirectory. The AF2 models are in the end renamed and copied to a subdirectory called final_pdbs
• num_designs_rfdiff: number of structures generated by RFdiffusion
• num_seq_per_target_mpnn: number of sequences generated per one RFdiffusion structure (each sequence is by default folded 5 separate times by AF2)

A note on outputs: when running on a cluster for each task a separate final_outputs.csv will be created. Run merge_csv.py along with the argument --output_dir to merge all csv files into one file.

Passing additional arguments to ProteinMPNN

Add parameters as required by ProteinMPNN to the pass_to_mpnn group in run.yaml (including trailing --).

Passing additional arguments to AlphaFold

Add parameters as required by AlphaFold to the pass_to_af group in run.yaml (including trailing --). Switches should be passed with an empty string value, like this --templates: "".

Examples for contigs

Contigs are the most important input (for now, guiding potentials are another input to be explored in the future). Here are a few guidelines to ease the start of your projects.

• Ranges starting with letters will be taken from the input.pdb
• Ranges without letters will generate that many residues
• /0 (mind the space!) will create a chain break
• For connecting two parts, you could use a contig such as this: '[A1-37/30-60/A42-70]'. If the original PDB has another chain B, it will not be given to RFdiffusion in this case.
• For connecting two parts with respect to another structure: '[C33-60/4-7/B1-30/0 C61-120]'. Chain B is given to RFdiffusion but no new structures are generated on it. Note: in the generated structure PDB, the chains will be labeled as A (connected A and B) and B (C in original).
• For connecting multiple chains: '[E10-68/5-15/C4-72/5-15/D78-145/5-15/A1-60/]'
• For leveraging symmetry: '[10/A29-41/10/0 10/B29-41/10/0 10/C29-41/10]'. In this case, the additional symmetry parameter must be: symmmetry=C3 (three subunits). Also, you cannot use contigs with interval lengths to be generated e.g. '[10-20/A29-41/10-20]'.
• Important: if you wish to force a number of newly generated residues, pass it as range: [A1-7/3-3/A11-12/1-1/A14-84/1-1/A86-96] (and not [A1-7/3/A11-12/1/A14-84/1/A86-96])
• More info is available in the RFdiff repo.

Automatic restart

In case of errors in any of the sub programs (RfDiffusion, ProteinMPNN, AlphaFold) prosculpt can automatically restart, so that the previous steps are not lost. To do so, pass auto_restart: n, where n ... number of allowed restarts, to the .yaml config.

Caveats/Troubleshooting

Binder design

• In the input pdb, target must start with chain B, and use subsequent chain letters (C, D ...) if it has multiple chains.
  • Otherwise, errors arise when passing into ProteinMPNN.
• Targets should not have missing residues s. t. distance between two residues would be larger than chain_break_cutoff_A
  • Otherwise, rechain.py will start a new chain at that point, causing errors when passing into ProteinMPNN.
  • Note that even without missing residues, calculated distance after the RfDiffusion step may sometimes be just slightly above 2 Å around prolines, so adjust accordingly.
• The chains_to_design parameter should not be present at all in the .yaml config file!
  • Otherwise, target will be fixed incorrectly and thus its aminoacids changed.