INDEPENDENT MODELING PROJECT
of
C687: Computing Methods in Biochemistry

Katianna Pihakari

*	PURPOSE	Model the 'hinge-bending' motion of a bacteriophage M6I T4 lysozyme.
*	MOTIVATION	The T4 lysozyme consists of two domains that are connected through a 'waist' region. The active site cleft is formed between the two domains. The two domains undergo a hinge-bending motion perpendicular to an axis passing through the waist. This motion allows reactants to enter and products to exit the active site.
*	SYSTEM	The mutant M6I, Met6 replaced by isoleucine, is studied. Five conformers differing in the hinge-bending angle have been observed in the two crystals formed by this mutant.

---

CONTENTS :

*	Introduction
*	Methods Extended Multiple Time Step Integrators Random Coil Transformation Umbrella Sampling
*	Studied System
*	Results In vacuo In water solution
*	Conclusions
*	References
*	List of Figures

---

Figure 1. The bacteriophage M6I T4 lysozyme. The conformer is modelled in water solvent and has a hinge-bending angle of 95 degrees. The secondary structure of the molecule is shown on the left while the solvent accessible Connolly surface is shown on the right.

INTRODUCTION

The bacteriophage T4 lysozyme forms two domains that are connected through a 'waist' region as seen in Fig.1. The thin waist of the molecule allows it to undergo a hinge-bending motion perpendicular to an axis passing through the waist. The opening and closing of the active site cleft allows reactants to enter and products to exit the active site.

In the present study the mutant M6I (Met6 residue replaced by an isoleucine) of the T4 lysozyme is modelled. The mutation enhances the hinge-bending motion but does not affect the catalytic activity or thermal stability of the molecule[1]. However, the mutant is temperature sensitive. This 164 amino acid monomeric enzyme consists of 11 alpha-helices out of which 7 has a negatively charged residue close to its N-termini. The total charge of the M6I T4 lysozyme is +8. The Engleman-Steitz hydrophobicity of the mutant as a function of the residue number is shown in Fig.2. The mutant conformers form two kinds of crystals: trigonal isomorphous crystals and orthorhombic crystals with four independent molecules in the asymmetric unit. Each of the four M6I lysozyme molecules in the orthorhombic crystal has a different hinge-bending angle. The difference between the smallest and the largest hinge-bending angle in these crystal conformers is 30 degrees. The mutant forms are distinctly different from the wild type crystal form[2]. The M6I T4 lysozyme crystals have an unusually high amount of solvent.

Figure 2. The Engleman-Steitz hydrophobicity of the M6I mutant of T4 lysozyme as a function of the number of the residue.

METHODS

All the calculations presented in this paper are done by using the extended system multiple time step molecular dynamics (MD) simulation methods as implemented in the computational package PINY_MD. The visualization is done with a commercial program Insight.

The extended multiple time step integrators[3,4,5] enabled the use of large time steps, 6 fs, during the simulation runs.

The name of the method comes from the fact that one large integration step is divided into N segments. In each segment only a portion of the forces on each atom are calculated. Thus, the computably inexpensive strong intramolecular forces are calculated often and the expensive long range interactions only once per large time step.

The present study concentrates on a large motion and slow time scale domain/breathing motion of a large molecule. This type of a calculation is time consuming. Thus, perhaps the most crucial method that enabled this study to be performed was the random coil variable transformations.

A large molecule can be considered to consist of a large number of subunits. A subunit can consist, for example, of a residue. However, the size of a subunit is related to the persistance length of the molecule, i.e. the length scale on which one given unit can flex/bend enough to interact with any other unit. To a first approximation, let's assume the subunits are connected by harmonic bonds.

Thus, consider a model potential

with the normal modes
.
We can define new masses for the modes by writing
.
They all have the same frequency and evolve in the same time scale.

The new masses can now be used to define conjugate momenta to u's. This transformation that can be used in general, {p_r,r}->{p_u,u}, is non-canonical and thus changes the dynamics. The particular transformation used in the present study is called Staging or Levy flight. It has been shown to work well in random coil calculations[6].

The probability of conformers as a function of a studied quantity such as the hinge-bending angle might have local/global minimum between the local/global maxima. While optimizing a structure of a molecule the structure might thus 'get stuck' in one of the probability maxima without ever visiting the other possible conformers.

The idea behind umbrella sampling is that a biasing potential is added to the real potential. The probability associated with the biasing potential has an umbrella shaped when plotted as a function of the studied quantity. This potential is used to 'drag' the molecule through the space and thus all the possible conformers get visited.

The weighted histogram method is employed to exactly round the bias.

STUDIED SYSTEM

Figure 3. The bacteriophage M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees.

The bacteriophage M6I T4 lysozyme, was studied both in gas, i.e. in vacuo, and in the water solution phase. The structure of the molecule can be seen in Fig.3. As already mentioned all the calculations were done by using the computational package PINY_MD. The used force field is CHARMM22[7]. For both the phases, the temperature was kept at a constant value of 300K and the periodic boundary conditions (PBC) were implemented in all the three dimenstions.

In the gas phase the number of particles, volume and the temperature of the system were kept constant. Thus the ensemble used was NVT. The temperature was held constant with Nose-Hoover chain (NHC)[8] thermostatting. One lysozyme molecule that consists of 2645 atoms was modelled in a cell. Six umbrella sampling windows were used for calculations.

The water solution phase lysozyme was studied using 10138 water molecules. Eight Cl- atoms were added into the solvent to neutralize the system. The water molecules were treated with an all-atom model with bonds and bend of the molecule constrained. The ensemble of the solution phase system was NPT, i.e. the number of particles, the pressure and the temperature were kept constant. Tne umbrella sampling windows were considered.

The calculations of both the systems were started with well equilibrated configurations. Each of the umbrella sampling windows of the in vacuo lysozyme was modelled for 120 ps and the solution lysozyme for 240 ps. In order to study the effect of the length of the calculation to the properties of the molecule, the sequence of the in vacuo conformers was repeated once giving a total of two whole simulations and the in solution conformers two times. The calculations were began with the final conformers of the previous calculations.

RESULTS

IN VACUO

The use of the umbrella sampling method showed that the most probable hinge-bending angle of the M6I T4 lysozyme molecule in vacuo at 300K is approximately 70 degrees. This hinge-bending angle is formed by the center-of-masses of the three domains forming the molecule. The probability distribution function of the hinge-bending angle is shown in Fig.4. A conformer with a hinge-bending angle of 70 degrees was chosen to be analyzed more thoroughly. As mentioned earlier two simulations of 120 ps were performed so that the latter run started from the configuration that the first one ended with. As can be seen from Fig.4, the plot of the probability is different. The difference defines the error of the simulation. In the present study the error can be seen to be small, a couple of degrees for the most probable conformer.

Figure 4. The true probability of the conformer of the M6I T4 lysozyme in vacuo as a function of the hinge-bending angle.

As can be seen from the ribbon representation in Fig.5, the in vacuo conformer with the angle of 70 degrees has a closed active site cleft. The closing of the cleft in vacuo is partly due to the hydrogen bonds formed inside and across the cleft. Due to the lack of water, intermolecular hydrogen bond formation is impossible and the bonds need to be formed between the acceptors and donors of the lysozyme.

Figure 5. The ribbon representation of M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees.

The Fig.6 shows the solvent accessible Connolly surface of the studied conformer of the lysozyme. It can be seen that no large hydrophilic or hydrophobic regions can be found on the surface of the present conformer.

Figure 6. The solvent accessible surface of the in vacuo M6I T4 lysozyme with a hinge-bending angle of 70 degrees. The surface is colored in accordance to the hydrophobicity of the residues so that blue is hydrophilic and red is hydrophobic.

In the backbone unit of a peptide there are two defined rotation angles: phi which is the rotation around the N-C(alpha) bond and psi which is the rotation around the C(alpha)-C bond. If phi and psi and plotted against each others the plot obtained is a Ramachandran plot. The Ramachandran plot for the discussed form of lysozyme is shown in Fig.7. It shows that most of the secondary structure consists of alpha-helices. This can be concluded also by studying the residues involved in the different secondary structure elements: more than 100 residues of the total 164 are involved in alpha-helices. Most of the remaining residues are involved in forming a beta-sheet on the wall of the active site cleft.

Figure 7. The Ramachandran plot of the M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees.

With the exception of glycine, all the amino acids have four different groups attached to C(alpha) and thus they are chiral molecules. More specifically, the nature has chosen them all to be of so called L-form[9]. They can, however, temporarily break this chirality. In the present study the standard deviations of the backbone structures were calculated. The standard deviation of 1.5 was chosen as a limit were the 'chirality is broken'. In Fig.8 backbone atoms of the residues that had a standard deviation of C(alpha) chirality between 1.5 and 2.0 are shown in red and more than 2.0 are shown in blue. By comparing the M6I lysozyme conformer in vacuo with a different hinge-bending angle, it can be noted that the C(alpha)-chirality seems to be broken in the regions that undergo the largest movements when the hinge-bending angle is changed from 85 degrees to 105 degrees[10].

Figure 8. M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees. The backbone atoms of the residues where C(alpha)-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.

Comparing the in vacuo conformer that has a hinge-bending angle of 70 degrees to the crystal structures~\cite{brook} it can be seen that no good alignment can be obtained. The best of the alignments with all the five crystal structures in Brookhaven Protein Databank entry 150L is shown in Fig.9. The crystal structure is shown in orange in Fig.9 and the in vacuo conformer in blue. The root-mean-square (RMS) deviation calculated with Insight in aligned positions was 2.7 even in the best case. It can also be seen in Fig.9 that the hinge-bending angle of the shown crystal structure is larger than of the in vacuo simulation. This is expected because the crystal structure has water molecules in the active site cleft that keeps the cleft open while the active site cleft of the in vacuo conformer is empty.

Figure 9. The best obtained alignment of the in vacuo M6I T4 lysozyme conformer with a crystal structure (Broohakhaven Protein Databank entry 150L). The crystal structure is shown in orange and the in vacuo conformer in blue.

M6I T4 lysozyme in water solution is closer to the interest of many people than the in vacuo study. As already mentioned, the solvent was not only water molecules but 8 Cl- ions were added to neutralize the system.

The probability distribution functions of the hinge-bending angles, given by the umbrella sampling study, shows that in solution phase the most probable hinge-bending angle is 95 degrees. The difference with the in vacuo conformer is 15 degrees and will be discussed later on. As mentioned already earlier three production simulations were done so that the following started from the configuration the previous one ended with. Each of the simulations were 240 ps which gives a total time of 2400 ps for each calculation. Fig.10 shows the probabilities of the conformers of all the three production runs as a function of the hinge-bending angle. It can be seen that the probability curve differs with a simulation more in solution than in gas phase. Thus the calculational error is bigger in the solution phase. The qualitative results remain the same in all the three solution phase simulations.

Figure 10. The true probability of the conformer of the M6I T4 lysozyme in water solution as a function of the hinge-bending angle.

The ribbon representation of the solution conformer with the hinge-bending angle of 95 degrees is shown in Fig.11. All the conformers shown during this discussion are obtained during the first production simulation. As can be seen in Fig.11 the active site cleft is noticeably larger than in the gas phase conformer shown in Fig.5. This allows water molecules fill the cleft, the reactants enter and the products leave the cleft. In the in vacuo conformer the cleft was reasoned to close due to the hydrogen bonds across the cleft. In the present, solution phase, system the hydrogen bonds are formed between the water molecules and the lysozyme. Thus the water molecules exist in the cleft constantly and keep it open.

Figure 11. The ribbon representation of M6I T4 lysozyme in water solution with the hinge-bending angle of 95 degrees.

The solvent accessible Connolly surface of the molecule is shown in Fig.12. In addition to the in vacuo conformer shown in Fig.6, no clear hyprophilicity trend on the surface can be found in the present conformer either. The solvent accessible active site cleft is, however, clearly seen in Fig.12. In order to see the orientation of the molecule a ribbon of the molecule is shown in the upper picture.

Figure 12. The solvent accessible surface of the M6I T4 lysozyme in solution with a hinge-bending angle of 95 degrees. The surface is colored in accordance to the hydrophobicity of the residues so that blue is hydrophilic and red is hydrophobic.

The Ramachandran plot of the discussed solution conformer in Fig.13 shows that the large amount of residues included in alpha-helices, as discussed for the in vacuo conformer, has not decreased. Also in the present conformer more than 100 residues are involved in the alpha-helix elements. Compared to the in vacuo conformer it can be seen that the solution conformer is more ordered which can be seen as more clear peaks in the Ramachandran plot of the solution phase system. The better ordering of the solution conformer was expected.

Figure 13. The Ramachandran plot of the M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees.

Comparing the Fig.14 to the corresponding in vacuo conformer figure, Fig.8, it can be seen that no similarities are found in the broken C(alpha) chiralities of these two phases. In addition, comparing the corresponding figures of the three solution production runs, in Fig.14-16, no similarities can be found even among the runs of the same system.

Figure 14. M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.

Figure 15. M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees. The second production run. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.

Figure 16. M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees. The third production run. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.

In order to study the change in the structure of the solution conformer through the three production runs, the conformers were superimposed and aligned from both ends. The alignments that can be seen in were done one end at a time. In Fig.17 the first run is shown with red, the second with green and the third with blue. It can be seen that no remarkable changes occured in the structure when the simulation was continued from 240 ps to 720 ps.

Figure 17. The ribbon representations of the solution conformers for all the three production runs. First run is red, second green and third blue. On the left the residues 15-60 are aligned and on the right the residues 80-160.

Comparing the studied solution phase conformer to the same crystal structures mentioned earlier, a better alignment for the solution phase than for the in vacuo conformer is found. The best alignment of the solution phase conformer with a crystal structure is shown in Fig.18. In this figure the crystal structure is shown in orange and the solution phase conformer in blue. The RMS deviation calculated by Insight for the shown alignment is 1.6. This value is approximately half of the RMS deviation for the best in vacuo crystal conformer alignment. It can be seen from the figure that the solution phase conformer has a larger hinge-bending angle than the crystal structure. This is caused by the packing forces that exist in the crystal structure but not in the solution phase conformer.

Figure 18. The best obtained alignment of the solution phase M6I T4 lysozyme conformer with a crystal structure (Broohakhaven Protein Databank entry 150L). The crystal structure is shown in orange and the in vacuo conformer in blue.

---

CONCLUSIONS

This study shows that the hinge-bending motion of M6I mutant of T4 lysozyme can be succesfully modeled with the extended system multiple time step molecular dynamics by using umbrella sampling methods and random coil transformations. The mutant was studied both in vacuo and in water solution.

As was expected, the hinge-bending angle of both the in vacuo and in the solution phase conformers differs from the hinge-bending angles of the five crystal structures (Brookhaven Protein Databank entry 150L). In the in vacuo system the active site cleft was shown to collapse and in the solution phase it was proved to be more open than in the crystal conformers.

When the simulation time was increased from 240 ps to 720 ps in the solution phase study, no remarkable change in the structure of the molecule was observed. In addition, no differences was found in the quantitative results.

REFERENCES

[1] Faber and Matthews, A Mutant T4 Lysozyme Displays Five Different Crystal Conformations, Nature 348, 15 Nov 1990.

[2] Weaver and Matthews, Structure of Bacteriophage T4 Lysozyme refined at 1.7 A Resolution, J.Mol.Biol. 193, 189 (1987).

[3] Martyna, Tobias, Klein, Constant Pressure Molecular Dynamics Algorithms, J. Chem. Phys. 101, 4177 (1994).

[4] Martyna, Tuckerman, Tobias, Klein, Explicit reversible integrators for extended system dynamics, Mol. Phys. 87,1117 (1995).

[5] Tuckerman, Berne, Martyna, Reversible multiple time scale molecular dynamics, J. Chem. Phys. 97, 1990 (1992).

[6] Tuckerman, Berne, Martyna, Klein, J. Chem. Phys. (1993).

[7] MacKerell, Jr., A. D. et al., All-atom empirical potential for molecular modeling and dynamics Studies of proteins, J. Phys. Chem. B 102, 3568 (1998).

[8] Martyna, Klein, Tuckermann, Nose-Hoover chains: The canonical ensemble via continuos dynamics, J. Chem. Phys. 97, 2635 (1992).

[9] Branden and Tooze, Introduction to Protein Structure, Garland Publishing, NY, USA, 1991.

[10] Pihakari, Samuelsson, Martyna, PNAS, 1999, submitted.

LIST OF FIGURES

Figure 1.	The bacteriophage M6I T4 lysozyme. The conformer is modelled in water solvent with a hinge-bending angle of 90 degrees. The secondary structure of the molecule is shown on the left while the solvent accessible Connolly surface is shown on the right.
Figure 2.	The Engleman-Steitz hydrophobicity of the M6I mutant of T4 lysozyme as a function of the number of the residue.
Figure 3.	The bacteriophage M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees.
Figure 4.	The true probability of the conformer of the M6I T4 lysozyme in vacuo as a function of the hinge-bending angle.
Figure 5.	The ribbon representation of M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees.
Figure 6.	The solvent accessible surface of the in vacuo M6I T4 lysozyme with a hinge-bending angle of 70 degrees. The surface is colored in accordance to the hydrophobicity of the residues so that blue is hydrophilic and red is hydrophobic.
Figure 7.	The Ramachandran plot of the M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees.
Figure 8.	M6I T4 lysozyme in vacuo with the hinge-bending angle of 70 degrees. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.
Figure 9.	The best obtained alignment of the in vacuo M6I T4 lysozyme conformer with a crystal structure (Broohakhaven Protein Databank entry 150L). The crystal structure is shown in orange and the in vacuo conformer in blue.
Figure 10.	The true probability of the conformer of the M6I T4 lysozyme in water solution as a function of the hinge-bending angle.
Figure 11.	The ribbon representation of M6I T4 lysozyme in water solution with the hinge-bending angle of 95 degrees.
Figure 12.	The solvent accessible surface of the M6I T4 lysozyme in solution with a hinge-bending angle of 95 degrees. The surface is colored in accordance to the hydrophobicity of the residues so that blue is hydrophilic and red is hydrophobic.
Figure 13.	The Ramachandran plot of the M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees.
Figure 14.	M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.
Figure 15.	M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees. The second production run. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.
Figure 16.	M6I T4 lysozyme in solution with the hinge-bending angle of 95 degrees. The third production run. The backbone atoms of the residues where CA-chirality is defined to be broken shown as balls. If the standard deviation is between 1.5 and 2.0 atoms are shown in red and if the standard deviation is more than 2.0 the atoms are shown in blue.
Figure 17.	The ribbon representations of the solution conformers for all the three production runs. First run is red, second green and third blue. On the left the residues 15-60 are aligned and on the right the residues 80-160.
Figure 18.	The best obtained alignment of the solution phase M6I T4 lysozyme conformer with a crystal structure (Broohakhaven Protein Databank entry 150L). The crystal structure is shown in orange and the in vacuo conformer in blue.