Difference between revisions of "Homology Based Structure Predictions Hemochromatosis"
(→3dbxA) |
(→MSA1) |
||
| Line 231: | Line 231: | ||
=== MSA1 === |
=== MSA1 === |
||
| − | <br style="clear:both;"> |
||
<figure id="MSA1basedModel"> |
<figure id="MSA1basedModel"> |
||
[[File:hemochromatosis_modeller_MSA1.gif|200px|thumb|left|model built by using MSA1 superimposed with 1A6Z]] |
[[File:hemochromatosis_modeller_MSA1.gif|200px|thumb|left|model built by using MSA1 superimposed with 1A6Z]] |
||
</figure> |
</figure> |
||
| + | <br> |
||
| − | |||
| + | The model resembles the original protein, and the MSA shows a common gap on position 19-42 of the target sequence. This is the gap in the 3DBX model simple alignment, on which the best model has been built (based on TMAlign score). |
||
| − | |||
| + | <br style="clear:both;"> |
||
=== MSA2 === |
=== MSA2 === |
||
Revision as of 11:51, 4 June 2012
Hemochromatosis>>Task 4: Homology based structure predictions
Contents
Riddle of the task
After endless battles and deadly traps you have finally reached the tomb's final chamber. As you enter it you notice that there is no sign of the treasures that were promised by the old map you found months ago. Suddenly you hear a loud noise behind you and a solid wall of stone blocks the only entrance into the room. You are trapped! You look around and notice something on the walls. On the left wall are four runes in an ancient language. Luckily its the same language as the notes on the map you deciphered and they translate into four single letters:
- C N O I
On the opposite wall you can see four simple symbols:
- a triangle
- a square
- a circle
- and a diamond (dt. Raute)
After further investigation you notice that the four symbols can be pushed into the wall, but you don't know what would happen and which one(s) to push.
What do you do?
Some hints:
- You only have to push one button and only once.
Short Task Description
Detailed description: Homology based structure predictions
In this task we want to assess the quality of 3D models built with homology information. For this we have employed several model building methods:
- Modeller
- SwissModel
- I-Tasser
- 3D-Jigsaw
The models were then evaluated by eye, TM-Score, TM-Align, and SAP.
Regarding the extra task:
- Extra diligence task: define a radius of 6 Angstrom around the catalytic centre and calculate the all atom RMSD in that region
As there is no defined active site for our protein, we couldn't do it.
Protocol
A protocol with a description of the data acquisition and other scripts used for this task is available here.
PDB templates
In order to find templates for our models we performed several searches for homologs with COMA and HHPred. We also reused the sequences from Task 2. However none of these methods yielded homologs with a sequence identity above 40% (except 1a6z which is HFE itself) that could be mapped to a PDB structure. The best results for COMA and HHPred are listed in <xr id="coma_t"/> and <xr id="hhpred_t"/> respectively. Therefore we could not generate models with a sequence identity above 80% and the 40%-80% range was limited to its lower bound. In addition to those shown below, we also used 2iad_B (P01921, 21.10% identity) from task 2.
<figtable id="coma_t">
| PDB ID | e-Value | Identities | Positives |
|---|---|---|---|
| 1a6z_A | 1.00E-63 | 100% | 100% |
| 1t7v_A | 1.20E-63 | 34% | 62% |
| 3nwm_A | 1.60E-57 | 30% | 57% |
| 1frt_A | 4.90E-65 | 28% | 66% |
| 2wy3_A | 2.80E-63 | 26% | 66% |
| 3ov6_A | 2.70E-55 | 20% | 67% |
| 1u58_A | 1.00E-54 | 19% | 59% |
| 3d2u_A | 7.50E-60 | 16% | 70% |
| 3dbx_A | 9.70E-59 | 16% | 71% |
| 3it8_D | 3.40E-52 | 15% | 65% |
</figtable>
<figtable id="hhpred_t">
| PDB ID | e-Value | Identities | Similarity |
|---|---|---|---|
| 1a6z_A | 1.80E-69 | 100% | 1.623 |
| 1k5n_A | 2.80E-68 | 40% | 0.725 |
| 1s7q_A | 5.80E-78 | 37% | 0.655 |
| 1t7v_A | 7.40E-68 | 36% | 0.702 |
| 3p73_A | 1.10E-69 | 35% | 0.638 |
| 3bev_A | 1.00E-69 | 34% | 0.684 |
| 2yf1_A | 2.40E-74 | 32% | 0.617 |
| 1zs8_A | 7.30E-68 | 30% | 0.553 |
| 2wy3_A | 3.30E-68 | 29% | 0.496 |
| 1cd1_A | 1.60E-67 | 21% | 0.394 |
</figtable>
Models
For our models we selected 1k5nA, 1zs8A, 2iadB, and 3dbxA as templates. We chose them to have a wide variety of sequence identities. The MSA models were created with a subset of these 4 templates: MSA1 contains all four templates, MSA2 the lower three (1zs8A, 2iadB, and 3dbxA), and MSA3 only 2iadB and 3dbxA. For the Modeller models we used both alignment methods (simple and 2d) to create the single template models. For the evaluation we compared the models to the "native" (complexed with beta-2-microglobulin only) and complex (complexed with beta-2-microglobulin and transferrin receptor) structure of HFE, 1a6zA and 1de4A respectively. Both structures contain only 275 of HFE's 348 residues (namely 23-297), thus excluding the signal peptide and the transmembrane and cytoplasmic region.
Modeller
In the following segment the models were evaluated that were built with modeller. The modelbuilding was based on different alignments, created with modeller itself and valuated against our HFE protein (1A6Z and 1DE4).
The presented pictures show the resulting models (in green) superimposed with the 1A6Z protein (red). The yellow lines indicate which positions PYMol had aligned to superimpose them.
1k5nA
<figure id="1K5NbasedModels">
</figure>
Both models are nearly the same, showing a good superimposition with chain A from 1A6Z. Also the model generated using the 2d alignment of modeller does not improve TM-score, TMAlign score and RMSD together. This is caused by the fact that the structure incorporating alignment method generates nearly the same alignment as the "normal" one (from modeller). Also the 2dalign() call aligns the first occuring amino acid way before the following residues. This might be an error and tried to avoided by manual correction.
1zs8A
<figure id="1ZS8basedModels">
</figure>
On the first look both structures look the same but as Table TODO shows, the 2d alignment based model has a Higher TM-Score, TMAlign score and RMSD score. This is caused by the fact that the 2d alignment (which incorporates structure) introduces more gaps but aligns the structure a bit better, resulting in the improvement of the model (against the normal-modeller-alignment based one).
2iadB
<figure id="2IADbasedModels">
</figure>
As indicated by the alignment lines and also detectable in the ones without the alignment lines these models have not a good predicted structure as the ones presented before. This also reflects the low TMAlign and high RMSD scores (compared to the other models). But it is also noticable that the beta sheets and the two helices have nearly the right orientation. Therefore homology (0.5 TMScore) is nearly achieved.
The 2d alignment introduces more gaps, but results in a more accurate model through this (0.49 TMAlign score against 0.40).
3dbxA
<figure id="3DBXbasedModels">
</figure>
Both of the 3DBX based models resemble the original HFE protein mostly, also indicated by their TMAlign scores. However, in this case the 2d alignment does not improve but worsen the prediction (against the simple alignment). The small slide in the 2d alignment of the first residues (agains the normal alignment) seems to cause this.
MSA1
<figure id="MSA1basedModel">
</figure>
The model resembles the original protein, and the MSA shows a common gap on position 19-42 of the target sequence. This is the gap in the 3DBX model simple alignment, on which the best model has been built (based on TMAlign score).
MSA2
<figure id="MSA2basedModel">
</figure>
MSA3
<figure id="MSA3basedModel">
</figure>
Evaluation
<figtable id="modeller_scores_native">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| 1K5N_2d | 272 | 0.1686 | 0.0846 | 0.0487 | 0.83298 | 2.200 (over 272 atoms) |
| 1K5N_simple | 272 | 0.1649 | 0.0846 | 0.0506 | 0.83358 | 2.325 (over 272 atoms) |
| 1ZS8_2d | 272 | 0.1698 | 0.0827 | 0.0432 | 0.85494 | 1.642 (over 272 atoms) |
| 1ZS8_simple | 272 | 0.1550 | 0.0790 | 0.0423 | 0.79841 | 2.302 (over 271 atoms) |
| 2IAD_2d | 272 | 0.1725 | 0.0836 | 0.0460 | 0.49103 | 2.166 (over 269 atoms) |
| 2IAD_simple | 272 | 0.1337 | 0.0607 | 0.0349 | 0.40162 | 3.705 (over 272 atoms) |
| 3DBX_2d | 272 | 0.1742 | 0.0892 | 0.0496 | 0.81698 | 2.374 (over 267 atoms) |
| 3DBX_simple | 272 | 0.1684 | 0.0882 | 0.0496 | 0.86512 | 1.524 (over 267 atoms) |
| MSA1 | 272 | 0.1680 | 0.0855 | 0.0496 | 0.83014 | 2.366 (over 270 atoms) |
| MSA2 | 272 | 0.3218 | 0.1811 | 0.0855 | 0.72682 | 1.889 (over 270 atoms) |
| MSA3 | 272 | 0.1530 | 0.0708 | 0.0358 | 0.40265 | 3.679 (over 272 atoms) |
</figtable>
<figtable id="modeller_scores_complex">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| 1K5N_2d | 272 | 0.1674 | 0.0800 | 0.0469 | 0.85238 | 1.986 (over 272 atoms) |
| 1K5N_simple | 272 | 0.1638 | 0.0800 | 0.0478 | 0.83836 | 2.142 (over 272 atoms) |
| 1ZS8_2d | 272 | 0.1715 | 0.0827 | 0.0450 | 0.86026 | 1.720 (over 272 atoms) |
| 1ZS8_simple | 272 | 0.1550 | 0.0790 | 0.0414 | 0.82124 | 2.059 (over 271 atoms) |
| 2IAD_2d | 272 | 0.1706 | 0.0836 | 0.0441 | 0.48549 | 2.279 (over 269 atoms) |
| 2IAD_simple | 272 | 0.1318 | 0.0653 | 0.0358 | 0.40088 | 3.655 (over 272 atoms) |
| 3DBX_2d | 272 | 0.1737 | 0.0873 | 0.0487 | 0.83662 | 2.088 (over 269 atoms) |
| 3DBX_simple | 272 | 0.1687 | 0.0901 | 0.0515 | 0.86877 | 1.492 (over 269 atoms) |
| MSA1 | 272 | 0.1698 | 0.0873 | 0.0524 | 0.84313 | 2.150 (over 270 atoms) |
| MSA2 | 272 | 0.3223 | 0.1783 | 0.0846 | 0.70465 | 2.011 (over 270 atoms) |
| MSA3 | 272 | 0.1574 | 0.0754 | 0.0377 | 0.40026 | 3.637 (over 272 atoms) |
</figtable>
One can see in both figures and tables (tm-align scores) that all models built via modeller (except the one based on 2IAD) lead to fairly good models. Also it seems that using multiple templates does not mean that the resulting models get better. This is indicated e. g. between the tm-align scores of MSA2 and MSA3 and the predicted structure of their models. It is plausible that the aligned sequence of 2IAD in the MSAs leads in case of MSA3 to a disruption of the modelpositions, but in case of MSA2 or MSA1 this one template is not enough anymore to alter the resulting model greatly. Unfortunately the sequences we used have a big overlap in the MSA, it would be informative how modelling performs when different fragments were modelled from different templates.
SwissModel
For SwissModel we only used the four single templates 1k5nA, 1zs8A, 2iadB, and 3dbxA. The alignments for each model can be found here.
1k5nA
<figtable id="swiss_1k5nA_stats">
</figtable>
The alignment from SwissModel for 1k5nA spans the residues 26 to 299 of HFE. It has a sequence identity of 38.71% and contains only three single gaps in 1k5nA's sequence, none of which break up a secondary structure. The estimated QMEAN Z-score is -1.92. <xr id="swiss_1k5nA_stats"/> shows the summary of the quality assessment provided by SwissModel. The error rates per residue fluctuate rapidly with peaks of over 4 anstrom. The anolea estimation shows two long regions with unfavorable scores from residue 69 to 110 and 128 to 206.
1zs8A
<figtable id="swiss_1zs8A_stats">
</figtable>
Similar to 1k5nA the alignment spans from 27 to 298 with a sequence identity of 29.56%. It contains three single gaps and one gaps of 3 residues within 1zs8A that don't split a secondary structure, but also one bigger gap of 7 residues that breaks up a helix. The QMEAN Z-score of -2.88 is worse than that for 1k5nA which isn't surprising given the drop of 9% sequence identity. The error rates are also worse with a maximum of 8 angstrom and they almost never go below 1 angstrom. The anolea graph exhibits the same unfavorable regions as in 1k5nA, though the first region is even worse and the second a bit shorter.
2iadB
<figtable id="swiss_2iadB_stats">
</figtable>
In contrast to the other three templates, the alignment for 2iadB is rather short (111 to 299). It contains one single gap and two 2 residue gaps, neither of them breaking a secondary structure in 2iadB. This and the sequence identity of only 21.76% might be the cause for the very bad QMEAN Z-score of -3.42. The error rates don't go as high as for 1zs8A, but they never drop below 2 angstrom for the first half of the alignment and improve only slightly in the second half. This is also reflected in the anolea distribution as there are almost no favorable regions during the first half. The second half (starting around 222) gets much better which correlates with the regions from the other templates.
3dbxA
<figtable id="swiss_3dbxA_stats">
</figtable>
3dbxA has the longest alignment of all templates with a length of 275 residues (24 to 298 in HFE), but it also has the lowest sequence identity with only 18.93 percent. The alignment has only two single gaps in 3bdxA's sequence, but one of them splits a sheet and the other one a helix. Nevertheless it has a QMEAN Z-score of -2.8 which is even slightly better than for 1zs8A. The error rates average around 3 angstrom for the first half of the alignment and slightly improve for the second. The anolea graph again shows the two unfavorable regions from 1k5nA and 1zs8A, though this time the second one is almost neutral and the first is slightly longer.
Evaluation
TM-Score again seems to have problems to correctly align both sequences and therefore provides no meaningful data (TM-Score, GDT-TS, and GDT-HA in <xr id="swiss_scores_native"/>). When comparing the TM-Score from TM-Align to the corresponding scores for the Modeller results (cf. <xr id="modeller_scores_native"/> and <xr id="modeller_scores_complex"/>)you can see that SwissModel does neither perform better nor worse than Modeller on a general basis.
Surprisingly 3dbxA is provides the best model (TM-Align 0.89) despite its low sequence identity and outperforms those from Modeller. Especially the weighted RMSD of 1.1 is quite good compared to the other models. In the case of 1zs8A SwissModel performs better than the simple alignment based Modeller model, but worse than the 2d alignment one. Regarding the weighted RMSD SwissModel outperforms Modeller for both alignment methods. 2iadB again performs worst (cf. Modeller results) as it suffers from the short alignment with HFE. Even with TM-Align the TM-Score barely reaches 0.5 which is the minimum threshold for similar folds. Despite its high sequence identity 1k5nA performs about as well as 1zs8A regarding the TM-Score, but has a much worse weighted RMSD.
<figtable id="swiss_scores_native">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| 1k5nA | 250 | 0.1626 | 0.0809 | 0.0478 | 0.84456 | 2.121 (over 272 atoms) |
| 1zs8A | 249 | 0.1449 | 0.0662 | 0.0377 | 0.83755 | 1.514 (over 271 atoms) |
| 2iadB | 165 | 0.1218 | 0.0680 | 0.0450 | 0.50849 | 2.805 (over 187 atoms) |
| 3dbxA | 252 | 0.1684 | 0.0836 | 0.0478 | 0.89308 | 1.111 (over 272 atoms) |
</figtable>
<figtable id="swiss_scores_complex">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| 1k5nA | 250 | 0.1611 | 0.0781 | 0.0469 | 0.85087 | 2.009 (over 272 atoms) |
| 1zs8A | 249 | 0.1450 | 0.0689 | 0.0377 | 0.83904 | 1.501 (over 271 atoms) |
| 2iadB | 165 | 0.1201 | 0.0671 | 0.0432 | 0.49068 | 3.172 (over 187 atoms) |
| 3dbxA | 252 | 0.1679 | 0.0818 | 0.0460 | 0.88762 | 1.203 (over 272 atoms) |
</figtable>
I-Tasser
For the I-Tasser models we used the option to provide a template. As templates we have used so far 1k5nA and 2iadB. The runtimes were extremely high with 4 (2iadB) to 5 (1k5nA) days. As there is a limit to only one submission per account and considering the long runtime, we won't be able to process more templates until tuesday. It should also be noted that I-Tasser incorporated 1a6z and 1de4 into the modelling process, both of which are pdb entries for HFE.
2iadB
<figure id="itasser_2iadB_1">
</figure>
<figure id="itasser_swiss_2iadB_1">
</figure>
I-Tasser provided 5 models for 2iadB. The best one with a C-score of -1.46 and an estimated TM-Score of 0.53±0.15 was selected for the evaluation.
Compared to all other models this one looks the worst (see <xr id="itasser_2iadB_1"/>). A whole helix (upper front in the figure) is not predicted. This is especially weird as this helix is quite good modelled in the SwissModel for 2iadB (<xr id="itasser_swiss_2iadB_1"/>). Additionally the majority of the sheets is missing.
On the other hand I-Tasser correctly predicted a part of the transmembrane helix around 307 to 330 (lower right corner in the figure) which is not included in the pdb file for 1a6z. I-Tasser also predicted a helix at the beginning of HFE in the signal peptide region, but this region is also not contained in the pdb file nor is it specified as a helix in uniprot, but it was also predicted as helix by PsiPred and ReProf in Task 3.
1k5nA
<figure id="itasser_1k5nA_1">
</figure>
The best model for 1k5nA provided by I-Tasser has a C-score of -0.75 and an estimation of 0.62±0.14 for the TM-Score.
The model fits the reference structure (1a6zA) quite well (see <xr id="itasser_1k5nA_1"/>). The long helices are well aligned, but the sheets seem to be quite shifted. I-Tasser also predicted a helix in the signal peptide region (background, center), just like in the previous 2iadB model.
Evaluation
The scores for both models are shown in <xr id="itasser_scores_native"/> and <xr id="itasser_scores_complex"/>. TM-Align calculates a TM-Score of around 0.79 for the best model for 2iadB. Although this is much higher than the other methods achieved for 2iadB it shows that 2iadB seems to be a quite bad template for HFE as even the incorporation of the pdb structures of HFE (1a6z, 1de4) in the modelling process didn't raise the score above many of the other models based on low identity templates. The same is true for 1k5nA. Overall the I-Tasser results are pretty underwhelming given that it is considered one of the best methods out there.
<figtable id="itasser_scores_native">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| 1k5nA_1 | 272 | 0.1682 | 0.0836 | 0.0478 | 0.84911 | 2.305 (over 272 atoms) |
| 2iadB_1 | 272 | 0.1699 | 0.0873 | 0.0506 | 0.78679 | 2.333 (over 270 atoms) |
</figtable>
<figtable id="itasser_scores_complex">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| 1k5nA_1 | 272 | 0.1666 | 0.0836 | 0.0496 | 0.85045 | 2.217 (over 272 atoms) |
| 2iadB_1 | 272 | 0.1692 | 0.0855 | 0.0478 | 0.79875 | 2.161 (over 272 atoms) |
</figtable>
3D-Jigsaw
The first attempt with 3D-Jigsaw contained the following models as those seemed to be about the best:
- 1zs8A (Modeller, 2d align)
- 3dbxA (Modeller, simple align)
- 3dbxA (SwissModel)
- MSA1 (Modeller)
- MSA2 (Modeller)
3D-Jigsaw failed to generate new models with these.
The second attempt contained all 15 models from Modeller and SwissModel combined. This time 3D-Jigsaw was able to generate 5 new models.
Models
<figure id="3dJigsawAll">
</figure>
All models seem to be almost identically (see <xr id="3dJigsawAll"/>). They also match the pdb structure (1a6zA) quite well. The only differences for the models are around the coiled regions not included in the pdb structure (residues 1-22 and 298-348). This is no surprise as all previous models also had the problem that the templates lacked these regions.
Evaluation
It is not surprising that all 5 models have about the same scores (cf. <xr id="jigsaw_scores_native"/> and <xr id="jigsaw_scores_complex"/>) as they almost look alike. Although they can compete with the best models from the other methods none of them are really astonishing. The TM-Score (TM-Align) and weighted RMSD are both lower than several of the provided models (e.g. 3dbxA).
<figtable id="jigsaw_scores_native">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| model_1 | 272 | 0.1647 | 0.0864 | 0.0533 | 0.83038 | 2.359 (over 272 atoms) |
| model_2 | 272 | 0.1648 | 0.0873 | 0.0542 | 0.83125 | 2.351 (over 272 atoms) |
| model_3 | 272 | 0.1648 | 0.0873 | 0.0542 | 0.83127 | 2.350 (over 272 atoms) |
| model_4 | 272 | 0.1648 | 0.0873 | 0.0542 | 0.83127 | 2.350 (over 272 atoms) |
| model_5 | 272 | 0.1648 | 0.0873 | 0.0542 | 0.83127 | 2.350 (over 272 atoms) |
</figtable>
<figtable id="jigsaw_scores_complex">
| Model | Common residues | TM-Score | GDT-TS | GDT-HA | TM-Align | Weighted RMSD |
|---|---|---|---|---|---|---|
| model_1 | 272 | 0.1637 | 0.0836 | 0.0506 | 0.83550 | 2.166 (over 272 atoms) |
| model_2 | 272 | 0.1637 | 0.0836 | 0.0506 | 0.83622 | 2.158 (over 272 atoms) |
| model_3 | 272 | 0.1637 | 0.0836 | 0.0506 | 0.83627 | 2.157 (over 272 atoms) |
| model_4 | 272 | 0.1637 | 0.0836 | 0.0506 | 0.83625 | 2.158 (over 272 atoms) |
| model_5 | 272 | 0.1637 | 0.0836 | 0.0506 | 0.83626 | 2.157 (over 272 atoms) |
</figtable>
Conclusion
TODO: make text :P
- HFE structure more conserved than sequence
- low identity =/= bad model
- alignment important
- correlation
