Homology modelling

1) Preparing the project directory

Create one directory for each project of homology modeling. An example project is included in biskit: biskit/test/Mod/project/

Starting point (and only input) is a fasta file with the primary sequence of the target protein (60 characters per line). The default name for this file is target.fasta. If you use default file names and default parameters, you can run the following scripts without any options -- as long as you are in the project folder. Refer to the full documentation for any non-standard options, file names, etc. Advanced options of each script are available with: script.py -help

Ultra-quick instructions -- for a standard run w/o cross-validation, call in this order:

search_sequences.py
search_templates.py
clean_templates.py
align.py
model.py

2) Search homologous sequences

>>> search_sequences.py

3) Search homologous 3D structures (Templates)

>>> search_templates.py

4) Clean template PDB files

>>> clean_templates.py

5) Optional: Set up cross-validation

>>> setup_validation.py

6a) Built multiple alignment (non-parallel)

For main project:

>>> align.py

Repeat in each cross-validation folder (project_root/validation/????/)

>>> align.py

6b) Or built multiple alignment both for the main project and for the cross-validation (parallelised)

>>> pvm  ## start pvm if it is not running
>>> align_parallel.py

7a) Build models (non-parallel)

For main project:

>>> model.py

Repeat in each cross-validation folder (project_root/validation/????/)

>>> model.py

7b) Or build models both for the main project and for the cross-validation (parallelised)

>>> pvm  ## start pvm if it is not running
>>> model_parallel.py

8) Optional: Evaluation of cross validation

Example: benchmark.py -d list-of-folders

>>>  benchmark.py -d  ~/biskit/test/Mod/project/validation/*

Example: analyse.py -d list-of-folders

>>>  analyse.py -d  ~/biskit/test/Mod/project/

A) Prepare a project folder for an homology modeling

Put a file with the target sequence into an empty directory.

B) Search for homolgous sequences: search_sequences.py

Syntax:  search_sequences.py [-q |target.fasta| -o |outFolder| -log |logFile|
               -db |database| -e |e-value-cutoff| -aln |n_alignments|
               -psi
               -... additional options for blastall (see SequenceSearcher.py)]

Options: -q      fasta file with query sequence (default: ./target.fasta)
         -o      output folder for results (default: .)
         -log    log file (default: STDOUT)
         -db     sequence data base
         -e      E-value cutoff for sequence search
         -aln    number of alignments to be returned
         -simcut similarity threshold for blastclust (score < 3 or % identity)
         -simlen length threshold for clustering
         -ncpu   number of CPUs for clustering
         -psi    use PSI Blast instead, experimental!!

Result: A folder "sequences" in the project folder containing two files:

• sequences/all.fasta: the list of all the homologous sequences found
• sequences/nr.fasta: a list of non redundant homologous sequences (1/cluster)

You can edit the search result by removing (or adding) sequences from nr.fasta.

C) Search of homolgous 3D structures (Templates): search_templates.py

Syntax: search_templates.py [-q |target.fasta| -o |outFolder| -log |logFile|
                            -db |database| -e |e-value-cutoff|
                            -aln |n_alignments| -psi
                             -... additional options for blastall
                             (see SequenceSearcher.py)]

Options: -q    fasta file with query sequence (default: ./target.fasta)
         -o      output folder for results (default: .)
         -log    log file (default: STDOUT)
         -db     sequence data base
         -e      E-value cutoff for sequence search
         -aln    number of alignments to be returned
         -simcut similarity threshold for blastclust (score < 3 or % identity)
         -simlen length threshold for clustering
         -ncpu   number of CPUs for clustering
         -psi    use PSI Blast instead, experimental (currently not working)!!

Result: A folder "templates" in the project folder containing:

• templates/all/: directory with all homologous structures (PDB format).
• templates/all.fasta: the sequences of all template structures
• templates/blast.out: the blast alignment
• templates/cluster_result.out: clustering of all template sequences
• templates/nr.fasta: non-redundant template sequences (1/cluster)
• templates/templates.fasta: non-redundant template sequences taken directly from the PDB-files
• templates/nr/: directory with the non-redundant template PDBs.
• templates/nr/chain_index.txt: the PDBs (plus chain IDs) which will be cleaned for modeller and T-Coffee in the next step.

Edit chain_index.txt and templates.fasta to remove or add templates or force using different chains. If adding templates, make sure the new entry in chain_index.txt points to an existing file (best put into templates/nr/).

D) Cleaning of the PDB files: clean_templates.py

Needs two input files (see Search of homologous structures):

• templates/nr/*.pdb
• templates/nr/chain_index.txt

This step prepares the PDB files in the directory templates/nr for the multiple alignement step (T-Coffee) and the building of the models (Modeller).

The script changes non standard residues to standard ones, removes hydrogen atoms, removes atoms with multiple configurations etc...

Syntax:  clean_templates.py [-o |output_folder| -i |chainIndex| -log |logFile|]

Options: -o       output folder for results (default: .)
         -i       chain index file for templates
         -log   log file (default: STDOUT)

Result: A folder "t_coffee" in the folder /templates containing:

• templates/t_coffee/*.alpha: CA input files for the structure alignment
• templates/modeller/*.pdb: input PDBs for Modeller

E) Setting up the cross-validation: setup_validation.py

This step creates a new homology modelling project for each template structure. The aim is to re-model all known template structures (pretending we don't know their structure). setup_validation.py links in the cleaned structures from the previous steps.

Syntax:  setup_validation.py [ -o |project folder(s)| ]

Options:  -o  one or several project folders (default: current)

Result: A folder "validation" containing:

• validation/|code|_pdb: a sub-folder for each template structure

  Each sub-folder will contain (taking validation/1EX7 as example):

  • reference.pdb: known structure of 1EX7
  • target.fasta: the primary sequence of 1EX7
  • t_coffee_template_files: the PDB files used as templates for re-building 1EX7
  • sequences: link to the sequences directory
  • templates/modeller: see D
  • templates/t_coffee: see D
  • templates/templates.fasta: see C

F) Multiple alignment with or without parallelisation: align.py and align_parallel.py

Build multiple alignments with T-Coffee. (A) one-by-one with by running align.py in the project folder and in each single cross-validation sub-folder. OR (B) by running the parallised align_parallel.py in the project folder.

Syntax: align.py [ -o |outFolder| -log |logFile| -h |host_computer| ]

Options: -o      output folder for results (default: .)
         -log    log file (default: STDOUT)
         -h      different computer for calculation  (default: local computer)
                -> must be accessible w/o password via ssh, check!

Syntax: align_parallel.py [ -d |list of folders| -h |host| -pdb |pdbFolder|
                       -ft |fastaTemplates| -fs |fastaSequences|
                       -fta |fastaTarget| -fe |ferror|]

Options: -d            [str], list of project directories [/validation/*]
         -h            int, number of hosts to be used [10]
         -a            first add hosts to pvm [yes]
         -pdb  str, pdbFolder for the pdb *.alpha
         -ft           str, path to 'templates.fasta'
         -fs           str, path to 'nr.fasta'
         -fta          str, path to 'target.fasta'
         -fe           str, path to the error file

G) Building the models: model.py and model_parallel.py

As before, the Modeller step can be run in two different ways: (A) Either one-by-one non-parallel (model.py) or (B) parallelised (model_parallel.py).

Syntax:  model.py [ -o |outFolder| -log |logFile| -h |host_computer| ]

Options:  -o    project folder (default: .)
          -log  log file (default: STDOUT)
          -h    host computer for calculation  (default: local computer)
                -> must be accessible w/o password via ssh, check!

Syntax:  model_parallel.py -d |list of folders| -h |host|
                         [-fta |fastaTarget| -pir |f_pir|
                          -tf |template_folder| -sm |starting_model|
                          -em |ending_model| -fe |ferror|]

Options:   -d   [str], list of project directory (full path)
           -h   int, number of hosts to be used
           -fta  str, path to find 'target.fasta'
           -pir  str, alignment filename
           -tf  str, directories for input atom files
           -sm   int, index of the first model
           -em   int, index of the last model
           -fe   str, filename to output errors from the Slave

H) Evaluation step by cross validation: benchmark.py

The first evaluation step is to determine how well we can reproduce a template structure by modeling them from the remaining templates. As the structure of each template is available we can then compare them with the dummy structures that we produce using Biskit/Mod.

First execute benchmark.py in each sub-directory of the validation folder.:

Syntax:   benchmark.py -d |list of folders|

Options:  -d   [str], list of project directory

To run all at once:

>>> benchmark.py -d validation/

Result: A folder validation/benchmark containing:

• validation/????/benchmark/Fitted_??.pdb: Each dummy model superimposed on its known structure. The fitted structure is called fitted_00.pdb refering to the dummy model number 00 etc. Currently, the iterative superimposition method implemented in biskit/Biskit/rmsFit.py is used (max. 10 iterations).
• validation/????/benchmark/rmsd_aa.out: gives the all-atom rmsd of the different dummy models. (1) without iterative fitting, (2) with iterative fitting and (3) the percentage of atoms that has been removed during the iterative fitting.
• validation/????/benchmark/rmsd_ca.out: same as above, but only for C-alpha atoms
• validation/????/benchmark/rmsd_res_??: gives the C-alpha rmsd for each residue in the current dummy model.
• validation/????/benchmark/PDBModels.list: pickled PYTHON list of PDBModels. Same as modeller/PDBModels.list but now each model contains the benchmark information in the atom and residue profiles: 'rmsd_aa', 'rmsd_ca', 'rmsd_res'. See PDBModel.profile()!

I) Evaluate the dummy models: analyse.py

Here, the aim is to give an accurate way to evaluate the models created by the module. Three sorts of results are provided. The first, global, the second local, and the last one is a simple visualisation of the result directly on the 3D structure of the model (score RMSD).:

Syntax:  analyse.py -d |list of folders| [ -s ]

Options:  -d   [str], list of project directories
          -s   show 'final.pdb' in PyMol

Result: A folder "analyse" containing:

• analyse/global_results.out: for each template, PDB code, the heavy atom rmsd(rmsd_aa) with and without iterative fitting, the C-alpha rmsd (rmsd_ca) with and without iterative fitting, the rmsd_aa with the percentage of outliers during the iterative fit, Moreover, the percentage of identities (mean) of the target with its templates, the best Modeller score and finally the number of templates which served to remodel each template.
• analyse/local_results.out: contains for each residue of the target the rmsd profile given by each template using the multiple alignment. Then, a score, which is the mean rmsd. This score will served as final evaluation for the models of our target.
• modeller/final.pdb: the "best" model where the B-factor (Temperature factor) are replaced by the score rmsd for visualisation, e.g., in VMD or Pymol.