Task alignments

Most prediction methods are based on comparisons to related proteins. Therefore, the search for related sequences and the alignment to other proteins is a prerequisite for most of the analyses in this practical. Hence we will investigate the recall and alignment quality of different alignment methods.

Theoretical background talks

The introductory talks should given an overview of

• pairwise alignments and high-throuput profile searches (e.g. Fasta, Blast, PSI-Blast, HHsearch)
• multiple alignments (e.g. ClustalW, Probcons, Mafft, Muscle, T-Coffee, Cobalt) and MSA editors (e.g. Jalview)

with special attention to advantages and limitations of theses methods.

You may want to look in my (Arianes) talk again. Here is the link.

Where to run the analyses

• You can run the analyses on your own computers.
• You can also use the student computer pool: i12k-biolab0n.informatik.tu-muenchen.de, where n goes from 1 to 9 (or more?). The file server does not have blast etc. installed!

Sequence searches

For every native protein sequence for every disease employ different tools for database searching: The methods to employ (minimally) are:

• Searches of the non-redundant sequence database big_80:
  • Blast
  • PSI-Blast using standard parameters with all combinations of
    • 2 iterations
    • 10 iterations
    • default E-value cutoff (0.002)
    • E-value cutoff 10E-10
  • HHblits (HHsearch) using standard parameters, since there is no big_80 for HHblits, search against Uniprot
• Data can be found in /mnt/project/pracstrucfunc13/data/ or /mnt/project/rost_db/data/

Note: Save intermediate files, e.g. a3m and hhm for the alignments and HMMs generated by HHblits, or checkpoint files for PSI-Blast. We will reuse this later.

For evaluating the differences of the search methods:

• compare the result lists (e.g. how much overlap, distribution of %identity and E-values)
• validate the result lists -- e.g. (you do no need to do all!)
  • using CATH, SCOP or COPS (/mnt/project/pracstrucfunc12/data/COPS/) to check whether found pdb entries fall into the same fold class
  • using GO to check whether sequences have common GO classifications
  • any other ideas how you could validate that the hits are really related?

Note: When we do comparisons, the data needs to be comparable. Therefore:

• Check the outcome of your simple blast search. If there are many significant hits, increase the number of reported hits (-v, -b or max_target_seqs depending on blast version and output format) until no more relevant hits are found. Use that parameter also for the PSI-Blast searches and use a similar setting for HHblits (Think about why we ask you to do this.)
• And of course: If your PSI-Blast and HHblits searches hit the limit (and your blast search didn't), also increase the number of reported hits!

CAVE: If your data set gets large, the PSI-Blast searches will take a while.

Note: There are a few catches that arise from the differences in how the tools operate:

• HHblits searches against the clustered Uniprot version. In the output the cluster representatives are listed together with the cluster members.
  • If you compare the representatives against a PSI-Blast result for big_80, you will get more hits for big_80.
  • If you compare the representatives plus the cluster members against big_80, you will get fewer hits for big_80.
• Come up with a way to generate comparable results. (There is also a complete database "big" which you can use for searching -- reusing the profiles from your big_80 search. -- Think about why we don't ask you to start out with a search against big.)
• big_80 is generated with CD-HIT, which prefers long sequences over shorter ones. Hence the number of pdb hits in your big_80 search is going to be low. Likewise, the Uniprot database for hhblits does not contain pdb structures. So, if you want to do the quality check using structure data, come up with a way to generate comparable results.
• To get all hits in pdb with HHblits (not just clustered hits), you can also use the pdb_full database.

Multiple sequence alignments

For calculating multiple sequence alignments, create a dataset of diverse sequences. Generate groups of 10 (20 for the third group) sequences where

• one contains only sequences with low sequence identity (<30%) (also low mutual similarity!)
• one contains only sequences with high sequence identity (>60%)
• one contains sequences covering the whole range of sequence identity.

Ideally there should be at least two sequences with pdb-structures in each group. You can use the structure of your target sequence as a second structure for 3D-Coffee.

The alignment methods to use on each of these groups are:

• ClustalW
• Muscle
• T-Coffee with
  • default parameters ("t_coffee your_sequences.fasta)
  • use of 3D-Coffee

Note:

• ClustalW should be on your path on the student machines, there is a version of T-Coffee on /mnt/opt/T-Coffee/bin/. If you include that in your path, you also have muscle.
• You may also use MAFFT if you like. Probably you should then leave out one of the other methods. ;-)

Compare your alignments (qualitatively! You do not need to run statistics!). Things to look for are:

• How many conserved columns?
• How many gaps?
• Are functionally important residues conserved?
• Are there gaps in secondary structure elements?
• Where do functionally important residues stand out most?

Points for discussion:

• Observe how the sequence identity in the groups of sequences influences the alignments.
  • Do all methods cope with low similarity?
  • Are residues that are aligned in the high similarity group still aligned when the low similarity sequences are added?
• Does the incorporation of structural information (3D Coffee) help?
  • Does it make a difference how many structures you include?
• Overall, what would be your criteria for a good alignment?
• Based on your experience, which method would you like to use in the future?