Phytochelatins are a group of polypeptides which have gained a lot of momentum in the field of synthetic
biology and metabolic engineering as they hold the potential to be employed as detoxification modules
more information see: Design).
As phytochelatins are an important compartiment of our project, we wanted to understand its ion
capabilities on a deeper level. Thus our aim was to characterize the polypeptide and to explore
possibilities to improve its manganese ion chelation.
As a result, we generated two altered sequences.
In the original sequence (1), we have a lot of glutamic acid - cysteine pairs. We assumed that the
acid pulls ions through electronegativity and that the thiol groups of the cysteine then bind the heavy
One way in which we tried to improve this binding was through exchanging the cysteines with histidines,
we assumed that the amino acid would bind to metal ions via its imidazole ring.
In the second altered sequence, we replaced cysteine with aspartic acid, since its electronegative
properties (in addition to Glutamic acid) could bind heavy metal attraction.
In the following, those three generated sequences will be referred to as:
Origin, for the Cysteine including sequence.
His, for the Histidin including sequence.
Asp, for the Aspartic acid including sequence.
For all sequences, we needed to determine a three dimensional structure and test its stability. After
we planned to analyze them further by for example determining and comparing the binding energies via
That is a technique to study protein binding-unbinding processes. As well it is possible to estimate the
binding free energy between those states.
We built models as specified in the following:
Figure 1: Workflow for the bioinformatic analysis of phytochelatin and its
Protein structure prediction
To model the three dimensional structures of our versions of the phytochelatin from their sequences, we
evaluated several web tools listed here:
Those tools all employ different methods based on either homology modelling, threading or ab initio, to
predict the structures based on the sequence.
Homology modelling methods can be used if there exist (enough) homolog proteins With similar sequence
known structure, as they compare the query peptide sequence against a database. The structures of close
matches can be used as templates. Threading is often used to find template structures in a parallelized
If there only exists little or no knowledge about homologous proteins, ab-initio
modelling might be necessary. In this case, a free energy minimization is performed while iterating over
different possible structures for smaller parts of the sequence.
In our case, there exist phytochelatins with known structures and similar sequences, but with much less
repeats than observed in our phytochelatin. After extensive research, we decided to work with a
threading method “I-TASSER”, as it aligns fragments of the query sequence with fragments of a database.
Furthermore, the structures obtained with “I-TASSER” scored highest in terms of confidence scores in
comparison to all other tools we tested.
Briefly, I-TASSER uses a procedure using four different steps: (8)
Template identification by LOMETS using Threading (3)
Fragment structure reassembly by replica-exchange Monte Carlo simulations (4) using ab
the threading results from step 1
Atomic level structure refinement by REMO (5) and FG-MD (6) using molecular dynamic simulations,
H-bond network optimization
Structure-based function interpretation by COFACTOR (7) by evaluating similar template proteins
about their functions and structures
I-TASSER generates a confidence score to measure the quality of a model-result. Thus, for all of our
structures we chose the model with the highest confidence score (C-Score). Models with a C-Score of -1,5
expected to have a correct fold. Our models have the following C-Scores:
Origin: -1,28 (Figure 2)
His: -1,12 (Figure 3)
Asp: -0,75 (Figure 4)
As these are acceptable confidence scores, we used the Origin, His and Asp models
for our next step.
Figure 2: I-TASSER result for the Origin sequence.
Figure 3: I-TASSER result for the His sequence.
Figure 4: I-TASSER result for the Asp sequence.
Molecular dynamics simulation
Next, we needed to test whether the predicted structures are plausible in a Molecular Dynamics (MD)
simulation. The method simulates the movement of a set of atoms over time by calculating and analyzing
forces between them. The structure being stable for a defined period of time is a good indicator for its
stability and thus plausibility in vivo.
To achieve this, we used GROMACS 2018.6 a MD package capable of simulating Molecular Mechanical (MM) and
Quantum Mechanical (QM) forces applied to all the particles of a molecule. (9)
Using MD simulations it can be assumed whether a molecule can withstand the forces of a defined system
as well calculate the way it will change when exposed to a specific system.
Those systems´ forces are defined via force fields. A force field is a set of parameters for each atom
or group in the simulation that is used to calculate the energy of the system. Since several different
fields for various molecule types exist it is imperative to find the most suitable for the problem at
Although, we can presume that our Origin sequence will be stable, since we got it from an already
iGEM-Part, we commenced our MD simulation with the Origin structure as well in order to prove this
presumption and to receive reference values for the other sequences. Thus, we started by following the
“Lysozyme in Water Tutorial” (10) with the exceptions mentioned here:
In our first approach we used the AMBER99-SB Force Field, because AMBER is pretty well known and was
recommended to us by Dr.
Thiel . Furthermore, we found articles on how to extend the
AMBER99-SB Force Field by values of heavy metals (in focus: Manganese).
Generally, it is challenging to come across force fields with these properties already built in so we
decided that the aforementioned approach is presumably the most promising. Accordingly, we employed the
following parameters provided as .mdp files:
In brief, the parameters mostly differ by the following aspects:
The most important thing here was to set the coulomb for electrostatic, and van der Waals for
intramolecular forces, effects to 1.4 (instead of 1.0). (13)
The simulation time was risen up to 50ns (instead of 1ns) for a more reliable test of
Figure 5: 50ns Videos of our Phytochelatin-sequences.
As a result, the Origin and our His variant molecules sustained in the AMBER Force Field for 50 nm,
our Asp prototype fell apart after ~48ns. This indicates that the Origin and His variant may be stable
vivo, an important indication for a putative application of the His variant in the wet lab.
Concerning the Asp variant, no force field is proven to model 100% of real conditions and a stability
ns are, compared to publications which require 20 ns, quite strong. Thus, we decided to place the Asp
variant within the OPLS-AA/L force field. The OPLS-AA/L force field is an all atom force field which was
used in the referred tutorial as well. (11) In this force field, the Asp variant remained very stable,
making it a second promising candidate for our wet lab work.
Of note, the His modification showed quite a lot of gyration in comparison to the other two variants,
could eventually lead to the protein being torn apart on the long run.
Binding Site Prediction:
In addition to our MD simulation, we used another tool to elucidate the probable binding sites of metal
(Zn2+, Cu2+, Fe2+, Fe3+, Ca2+, Mg2+, Mn2+, Na+, K+) for Origin, His and Asp using IonCom. IonCom is a
webtool from zhanglab to predict possible binding sites to metals. It combines methods like ab initio
training with template-based transferals (2).
Contrary to our suspicion of additional sulfide-bonds from the thiol-groups of the Cysteine in the
sequence, we found Glutamic Acid (E15 & E19) to bind to Mn2+. Additionally, the Origin sequence does not
have active sites for Fe2+, Fe3+, Na+, K+)). In contrast, our His and Asp are both binding to all 7
different metal types and therefore appear to be more versatile than the Origin sequence.
To conclude, we suggest testing and using both of our phytochelatin variants additionally to the
sequence in the wet lab, as they appear to be stable and more versatile in terms of ion binding.
Despite not being able to execute all the planned tasks to evaluate phytochelatin alternatives and
them to the origin sequence, we can conclude that both of our variants are acceptable, stable
since they sustained for 50 ns in the forcefield. Asp, although torn apart in the AMBER99-SB force
performed extremely stable in the OPLS force field.
Both phytochelatin variants have the potential to improve the binding of phytochelatin to metals, but
considering the higher radius of gyration, we suggest to focus on the Asp variant.
As an outlook, the performance of Umbrella sampling to estimate the free energy between bound and
structures would present us with a putatively higher variant amongst our alternatives. Furthermore,
considering quantum mechanical aspects in our force field would potentially increase the yield of our
simulation, as recommended by our collaboration partner from Waterloo.