What is CS-ROSETTA?

To date, interpretation of isotropic chemical shifts in structural terms is largely based on empirical correlations gained from the mining of protein chemical shifts deposited in the BMRB, in conjunction with the known corresponding 3D structures. Chemical-Shift-ROSETTA (CS-ROSETTA) is a robust protocol to exploit this relation for de novo protein structure generation, using as input parameters the 13Ca, 13Cb, 13C', 15N, 1Ha and 1HN NMR chemical shifts. These shifts are generally available at the early stage of the traditional NMR structure determination procedure, prior to the collection and analysis of structural restraints. CS-ROSETTA approach, as shown below, utilizes SPARTA-based selection of protein fragments from the PDB, in conjunction with a regular ROSETTA Monte Carlo assembly and relaxation method. Evaluation of 16 proteins, varying in size from 56 to 129 residues yielded full atom models that have 0.7-1.8 angstrom root-mean-square deviations for the backbone atoms relative to the experimentally determined X-ray or NMR structures. The strategy also has been successfully applied in a blind manner to eight structural genomics targets with molecular weights up to 16 kDa, whose conventional NMR structure determination was conducted in parallel. This protocol potentially provides a new direction for high-throughput NMR structure determination, in particular in structural genomics.

Contents

Download and installation

Preparation of input chemical shift table
    Identification and exclusion of flexible tails and loops

How to use CS-ROSETTA
    Fragment selection
    Protein structure generation
    Evaluation of CS-ROSETTA models

   setenv rosettaDir /home/software/ROSETTA
   setenv csrosettaDir /home/software/CSROSETTA

   if (-e $baseDir/com/csrosettaInit.com) then
      source $baseDir/com/csrosettaInit.com
   endif

Preparation of input chemical shift table

CS-ROSETTA utilizes the protein backbone 15N, 1HN, 1HA, 13CA, 13CB and 13C' chemical shifts as inputs and search a structural database for best matched fragments. The chemical shifts need to be in TALOS format, as defined at http://spin.niddk.nih.gov/NMRPipe/talos/#preparing shifts. If starting from a standard BMRB format chemical shift file, a C-shell script $baseDir/com/brmb2talos.com can be used to generate a TALOS-format chemical shift file using:

   brmb2talos.com bmrbCS.str > inCS.tab

   DATA SEQUENCE MQYKLVINGK TLKGETTTKA VDAETAEKAF KQYANDNGVD GVWTYDDATK
   DATA SEQUENCE TFTVTE

   VARS RESID RESNAME ATOMNAME SHIFT 
   FORMAT %4d %1s %4s %8.3f
   1 M CA 54.519
   1 M CB 29.320
   1 M HA 4.189
   2 Q C 174.318
   2 Q CA 55.632
   2 Q CB 30.865
   2 Q HN 8.347
   2 Q HA 5.109
   2 Q N 123.775
   ...

Residues from the disordered tails and loops are identified from the input chemical shifts and the criteria in the paper:
(1) S2 < 0.7 from the RCI analysis, and/or
(2) no "good" predictions using the program TALOS

Positively identified flexible residues in the N- and C-terminal tails should be excluded from the structure prediction by preparing a 'truncated' input chemical shift table file excluding data (both sequence and chemical shifts) from those residues. For long flexible loops, it is advantageous to exclude all terms from the ROSETTA full-atom energy and the chemical shifts rescoring (see details in Evaluation of CS-ROSETTA models)

Top

How to use CS-ROSETTA

CS-ROSETTA has been designed to work in a black-box manner, and is supported by multiple C-shell/C++/python scripts/programs. To use CS-ROSETTA for protein structure prediction, users only require to: (1) run master script runCSRjob.com (provided with an input experimental chemical shift file) to select fragments from the structure database and to prepare a starting package for the subsequent ROSETTA structure generation, (2) run runRosetta.com (prepared in the previous step) to perform ROSETTA structure prediction and generate full-atom models, (3) run runCSrescore.com to add a chemical shift term to the CS-ROSETTA full-atom "energies". A flow diagram of the CS-ROSETTA protocol implemented in this package looks as shown below:

1. Fragment selection - runCSRjob.com

The script runCSRjob.com is the master script to generate MFR fragments and inputs for ROSETTA. To use this script, an input chemical shift filename MUST be specified in a command line such as:

   runCSRjob.com inCS.tab

   set OUTPUT_DIR = rosetta
   set MFR_EXCL = "PDB1_ PDB2A PDB3X"
   set PDB_NAME = t000
   set CHAIN = _
   set TAG = aa
   set N_MODEL = 5000 

The output of running this script will be the MFR-selected 9-residue and 3-residue fragments in ROSETTA format (with file names of $OUTPUT_DIR/$TAG$PDB_NAME$CHAIN09_05.200_v1_3 and $OUTPUT_DIR/$TAG$PDB_NAME$CHAIN03_05.200_v1_3, respectively), a FASTA sequence file ($OUTPUT_DIR/$PDB_NAME$CHAIN.fasta), a ROSETTA paths definition file ($OUTPUT_DIR/paths.txt) and a script file ($OUTPUT_DIR/runRosetta.com) for starting ROSETTA structure generation.

It is recommended for users to read through and understand this script, which would facilitate use of the MFR program and other supporting scripts, and solve any problems during the MFR search (for example, by checking the intermediate and log files). This script (1) Prepares input for the MFR program, (2) Runs the MFR fragment search, (3) Convert fragments to ROSETTA format and (4) Prepares inputs and scripts for running ROSETTA structure generation. In details:

(1) Prepare MFR inputs

The current MFR program in NMRPipe requires a TALOS format chemical shift table file as input (see Section of "Preparing of chemical shift table"). A "reference PDB" coordinate file is needed to provide the MFR program with the protein's amino acid sequence. For this purpose, runCSRjob.com creates a dummy "reference PDB" coordinate file from its FASTA sequence using:

   fasta2pdb protein.fasta > dummy_ref.pdb

(2) Run MFR fragment search

The script runCSRjob.com then conducts a MFR fragment search following a command line of:

   mfr.com -cs inCS.tab -ref ref.pdb -out out.tab -excl $MFR_EXCL >& log

where $MFR_EXCL is a list of protein names (without .pdb suffix) that users like to exclude from the MFR structure database prior to the fragment search. MFR fragment selection from searching the structural database is a time-consuming job (1-3 hours); the 'log' file can be used to check the progress.

(3) Convert fragments to ROSETTA format

The MFR-selected 9-residues and 3-residues fragment candidates, generated after running the MFR fragment search, have file names of the types frag9.$PDB_NAME.mfr.tab and frag3.$PDB_NAME.mfr.tab, respectively. Script mfr2rosetta.com can be utilized to convert these fragments to standard ROSETTA format:

   mfr2rosetta.com -mfr frag9.$PDB_NAME.mfr.tab -segLength 9 > frag9.$PDB_NAME.rosetta.tab
   mfr2rosetta.com -mfr frag3.$PDB_NAME.mfr.tab -segLength 3 > frag3.$PDB_NAME.rosetta.tab

(4) Prepare ROSETTA running package

runCSRjob.com generates at the last step a new directory (defined by variable $OUTPUT_DIR) and prepares the following required inputs and script for running ROSETTA structure generation:

• 9-residues and 3-residues fragment files $TAG$PDB_NAME$CHAIN09_05.200_v1_3 and $TAG$PDB_NAME$CHAIN03_05.200_v1_3

• FASTA sequence file $PDB_NAME$CHAIN.fasta

• ROSETTA paths definition file paths.txt

• starting script file runRosetta.com for ROSETTA structure generation

2. Protein structure generation - runRosetta.com

Script runRosetta.com generated by running the script runCSRjob.com can be used to start the standard ROSETTA structure generation (fragment assembly and full atom relaxation). This script includes solely a standard ROSETTA command line defining the inputs names, the parameters for fragment assembly and full-atom relaxation

   ${rosetta} $TAG $PDB_NAME $CHAIN -silent -output_silent -increase_cycles 10 -new_centroid_packing 
   -abrelax -output_chi_silent -stringent_relax -vary_omega -omega_weight 0.5 -farlx -ex1 -ex2 
   -termini -short_range_hb_weight 0.50 -long_range_hb_weight 1.0 -no_filters -rg_reweight 0.5
   -rsd_wt_helix 0.5 -rsd_wt_loop 0.5 -output_all -accept_all -do_farlx_checkpointing -relax_score_filter 
   -record_irms_before_relax -acceptance_rate 1.0 -filter1a 10000 -filter1b 10000 -nstruct $N_MODEL

Running script runRosetta.com in the working directory (where the selected fragments, FASTA sequence and paths.txt are stored) will start the standard ROSETTA structure prediction. The output of the ROSETTA run is a so-called silent output file with name $TAG$PDB_NAME.out, which includes the scores and full-atom descriptions for all accepted models.

The ROSETTA is intensively computationally demanding (in the order of 5-10 minutes per model). Therefore, use of computer clusters is highly recommended. In order to run multiple ROSETTA jobs in parallel for a given project, users can simply run the same script runRosetta.com for each 'parallel' job in the same working directory; ROSETTA will output the results from each job to a single silent output file. (Note: Please replace the option "-output_silent_gz" with "-output_silent" in the previously distributed runRosetta.com script, or download the newest distribution, in order to fix the problem during parallel ROSETTA jobs)

3. Evaluation of CS-ROSETTA models - runCSrescore.com

A ROSETTA silent output file $TAG$PDB_NAME.out contains a header in its first two lines (first line: the sequence information, second line: definition of scores and values), and for each model, the scores and values defined in the header (first line) and the residue-specific description of the full-atom model. An example ROSETTA silent output file is shown below:

SEQUENCE: MQYKLVINGKTLKGETTTKAVDAETAEKAFKQYANDNGVDGVWTYDDATKTFTVTE
SCORE: score env pair vdw hs ss sheet cb rsigma hb_srbb hb_lrbb rg co contact rama bk_tot fa_atr fa_rep fa_sol h2o_sol hbsc fa_dun fa_intra fa_pair fa_plane fa_prob fa_h2o h2o_hb gsolt sasa omega_sc description
SCORE: -116.07 -12.24 -1.62 0.31 -3.25 -72.74 0.34 15.99 -28.93 -18.75 -29.62 11.42 17.69 0.00 -13.19 -134.45 -165.57 7.34 84.55 0.00 -3.13 22.82 0.07 -3.11 0.00 -15.74 0.00 0.00 45.19 3958.98 7.55 S_0001_9019
1 E -81.616 95.109 -180.850 0.000 0.000 0.000 -60.979 -52.601 -70.896 0.000 S_0001_9019
2 E -103.775 140.428 181.156 0.662 1.438 3.457 -61.314 178.939 7.568 0.000 S_0001_9019
3 E -128.914 145.607 180.081 -1.325 0.472 6.551 -62.403 85.133 0.000 0.000 S_0001_9019
...

   extract_lowscore_decoys.py silent_file.out N > new_silent_file.out

If residues in flexible loops are positively identified (e.g., by RCI analysis), the energy terms for those residues should be excluded from the ROSETTA full-atom energy terms. Therefore, full-atom energies of models in a ROSETTA silent output file could be rescored by a script runRosettaRescore.com:

   runRosettaRescore.com silent_file.out ref.pdb

   ${rosetta} -extract_segment -segments seg.txt -n ref.pdb -s silent_file.out -fullatom -fa_input -all -termini

   6 60 80 99 120 130

The script runCSrescore.com cen be used to apply the chemical shift based rescoring for ROSETTA full-atom models, using the following command line:

   runCSrescore.com silent_file.out inCS.tab

1. extracts all-atom PDB models from the ROSETTA silent output file silent_file.out
2. for each set of PDB coordinates: runs SPARTA chemical shift predictions; calculates the chi2 values between the SPARTA-predicted chemical shifts and the input experimental chemical shifts (stored in inCS.tab); stores the chi2 values as a table file CS_chi2.txt
3. collects the raw all-atom energies from silent_file.out and generates a table file name.rawscore.txt, apply corrections to raw energies using the chemical shift scores collected in (2) and generates a table file name.rescore.txt with the rescored energy information.
4. calculates the RMSD deviation relative to the model with the lowest re-scored full-atom energy and generates a final output table file name.rms.rescore.txt.

Details of each step include:

(1) Full-atom PDB coordinates extraction

To rescore the ROSETTA full-atom models, the full-atoms PDB coordinates are required to be 'extracted' using the full-atom description encoded in the silent output file, which can be done by using script extract_pdb.com and a command line such as:

   extract_pdb.com silent_file.out

This script actually runs the following ROSETTA command to generate full-atom PDB coordinates from the information in silent output file 'silent_file.out':

   rosetta -extract -s silent_file.out -fa_input -all -termini -write_atoms_only

Note that the generated PDB coordinates will be stored in the directory defined in paths.txt. The extracted PDB coordinate files will be named according to their 'description' labels defined in the silent output.

   OUTPUT PATHS:
   movie ./
   pdb ./

(2) Calculate predicted chemical shifts and "chemical shifts scores" for ROSETTA full-atom models

The ROSETTA full atom models are evaluated by comparing the initial chemical shift inputs with chemical shifts predicted for the models by SPARTA. SPARTA takes the standard PDB coordinates for proteins and predicts the backbone 15N, 1HN, 1H , 13C , 13C and 13C' chemical shifts using the following command line:

   sparta -in inPDB.pdb -ref inCS.tab

where inPDB.pdb contains the input PDB coordinates file, inCS.tab is the input experimental chemical shifts (the same input as used for the MFR fragment search). The output is a table file with predicted chemical shifts, and the chi2 value between the predicted and experimental chemical shifts is calculated.

The script runCS_rescore.com runs SPARTA chemical shift prediction for each set of full-atom PDB coordinates generated above, and stores all outputs in a default directory "./pred".

(3) Collect energy score and apply correction using chemical shifts chi2 value

For each full-atom model, its name and raw ROSETTA full-atom energy score are extracted from the ROSETTA silent output file silent_file.out and stored in a file "name.rawscore.txt". Next, the energy scores will be 'corrected' by using the corresponding chemical shifts chi2 values (stored in file "CS_chi2.txt"). The new file with the model names and re-scored energy will be "name.rescore.txt", which for example contains the following contents:

   # name raw_energy chi2 rescore_energy
   S_0001_9019 -116.07 251.854 -53.1065
   S_0002_6198 -124.99 140.1 -89.965
   S_0003_5600 -123.80 115.019 -95.0452
   S_0004_2308 -118.78 148.107 -81.7533
   S_0005_3837 -111.47 75.3417 -92.6346
   S_0006_9327 -119.05 110.179 -91.5052
   ...

(4) Calculate RMSD values to the lowest energy model

The C_alpha RMSD values between each model and the model with the lowest re-scored energy are calculated using the script pdbrms. pdbrms is a program written in C++ and can be used to calculate the RMSD values between one PDB coordinate and a set of protein PDB coordinates. For example, to calculate the C-alpha RMSD between ref.pdb and 1.pdb, 2.pdb, 3.pdb, the following command line can be used:

   pdbrms ref.pdb 1.pdb 2.pdb 3.pdb

The script runCSrescore.com will identify the model with the lowest rescored energy, and calculate the C_alpha RMSD values to this model for all models in the ./output directory. The final output file "name.rms.rescore.txt" will contain the model names, C-alpha RMSD values and re-scored energies:

   # name rmsd rescored_energy
   S_0001_9019 1.883 -53.1065
   S_0002_6198 1.384 -89.965
   S_0003_5600 0.984 -95.0452
   S_0004_2308 2.788 -81.7533
   S_0005_3837 2.718 -92.6346
   S_0006_9327 1.943 -91.5052
   ...

By default, the ROSETTA running directory (./rosetta) contains the following contents upon finishing the above CS-ROSETTA protein structure generation:

Top

How to select CS-ROSETTA models

Criteria for convergence and accepting models

After finishing CS-ROSETTA structure generation, users have to decide whether the ROSETTA models are acceptable. For this purpose, it is convenient to plot the "landscape" of (re-scored) ROSETTA full-atom energies of all models with respect to their C_alpha RMSD values relative to the lowest-energy model, using the data stored in a file "name.rms.rescore.txt".

1. If the low energy models cluster within less than C_alpha RMSD of about 2 angstrom from the model with the lowest (re-scored) energy (see the following example plot from the structure prediction of protein GB3), the structure prediction is deemed successful and the 10 lowest energy models are accepted.
    
2. If no clustering around the low energy model is observed (see the following example plot generated for protein nsp1), the structure prediction has not converged and the low energy models can not be accepted at this stage.
   

Number of models required

By using the current method implemented in CS-ROSETTA package, 5,000 to 20,000 predicted CS-ROSETTA models are generally required to obtain the convergence. For small proteins (<= 90-100 amino acids), 1,000 to 5,000 predicted CS-ROSETTA models often sufficient. ROSETTA takes about 5-10 minutes to calculate one all-atom model on a single 2.4GHz CPU.

Top

FAQs

Why no (or no sufficient, i.e., <200) fragment candidates obtained after MFR fragment search using runCSRjob.com? Answer Current MFR program (implemented in NMRPipe) search for fragment candidates with a default threshold of 3.0 for chemical shift matching score (-csThresh 3.0), which is good for most cases. If the experimental chemical shifts involving in a given target fragment are unusual, there will be no (or no sufficient) fragment candidates could be selected with csThresh<3.0. When this happened, please first check and make sure the experimental chemical shifts involving in the problematic target fragments, as well as the referencing, are correct. (Script runCSRjob.com will list the number of selected fragment candidates for each target fragment, and the "problematic target fragment" will be listed with a warning message) Users can "force" the MFR program to generate fragment candidates by using larger csThresh value, for example 5.0, by: copy script runCSRjob.com to current directory replace the line "mfr.com -cs ${CS_NAME} -ref ${REF_NAME} -out ${PDB_NAME}.mfr.tab -excl ${MFR_EXCL}>& ${PDB_NAME}.mfr.log" with "mfr.com -cs ${CS_NAME} -ref ${REF_NAME} -out ${PDB_NAME}.mfr.tab -excl ${MFR_EXCL} -csThresh 5.0>& ${PDB_NAME}.mfr.log" run runCSRjob inCS.tab to get fragment candidates. Note that those fragment candidates, which could only be obtained by MFR fragment search with a larger csThresh value, are commonly not accurate. Use of "H", instead of "HN", as atom name for amide proton in the experimental chemical shift talble is another reason of no fragment candidates (a bug of current MFR module)

How to exclude some proteins from the structural database during the fragment search? Answer find the exact names of the proteins you like to exclude in the file $CSROSETTA/PDBH/resolution.tab copy script runCSRjob.com to current directory find the line set MFR_EXCL = "NONE_", and repalce NONE_ with the list of protein names run runCSRjob inCS.tab to get fragment candidates.

Can I use incomplete chemical shifts data for CS-ROSETTA structure generation? Answer Current CS-ROSETTA protocol was initially evaluated using ~20 proteins with (nearly) complete chemical shift assignments (15N, 13Ca, 13Cb, 13C', 1Ha and 1HN). However, CS-ROSETTA protocol is still good for proteins with certain degree of 'incomplete' chemical shift data, for example, two successful proteins shown in our paper don't have 13C' chemical shifts. The impact of the chemical shift data completeness on the performance of CS-ROSETTA generation is various, please consult us in order to better use your chemical shift data with different completeness.