Sonia Gupta

The flexibility of the mutated (T44A & I6A) B1 domain of protein G was probed using Energy Minimization and Molecular Dynamics (at room temperature) with a particular emphasis on the solvent exposed guest site T53. This was an attempt to study the changes in the local flexibility of my protein at room temperature. For this, a particularly chosen, solvent exposed guest site T53 was mutated to a different amino acid with a different length of side chain (T53A) to observe its effect on the local flexibility. This could also give us information about the changes in solvent accessibility of position 53 as a result of the mutations. During these calculations, the nearest neighbors of the surrounding host environment were also taken into account. The results were analyzed in terms of the changes observed in the hydrogen bond distances between these neighbors and the guest site amongst many other factors. An attempt has also been made to calculate the value of S² (order parameter) for the backbone N-H of residue 53 in both the structures. The order parameter reflects the amplitude of motion of the backbone N-H. All these calculations were done under the assumption that only the residues lying close to the guest site would get affected significantly after mutation. So the investigations were done in a sphere of a defined radius, centered on the guest site, freezing everything outside that sphere. All the calculations were done under vacuum. It was assumed that similar results would be obtained in vivo. (Vacuum was considered only to get results in a shorter period of time.)

The host protein taken up for this study was the immunoglobulin G (IgG) binding, 56 amino acid B1 Domain of Protein G (GB1) from group G Streptococcus, having two mutations (I6A and T44A). The structure of the B1 domain has been determined both by NMR spectroscopy and high-resolution X-ray crystallography^{1, 2, 3}. The small size, remarkable stability and the absence of any prosthetic groups or disulfide bridges made it a good candidate for our study. Moreover, quite a lot is already known about this protein, and particularly the B1 domain, from various ongoing studies.

Streptococcal and Staphylococcal strains produce several cell surface proteins that bind to the constant region of immunoglobulin allowing the bacteria to evade the host immune system⁴. Protein G is a large, multidomain cell surface protein of group G Streptococcus that is thought to evade the host defenses through its protein binding properties. It contains several highly conserved immunoglobulin G (IgG) binding domains^{2, 5-6}. The B1 domain (56 amino acid residues, for sequence see Sequence 1), which is involved in binding to the F_c region of human IgG⁷, is one of the smallest stable folded globular domains known and has been studied extensively by biophysical methods.

The domain (Figure 1) consists of two beta-hairpins (residues 1-20 and residues 42-56) that are associated to form a four-stranded mixed antiparallel/parallel beta-sheet. A single alpha-helix (residues 22-37) lies across one face of the sheet. There are no proline residues or disulfide bonds in the B1 domain. There is an extensive network of hydrogen bonds in this domain (Figure 2).

Figure 1: Ribbon representation of the overall secondary structure of the B1 domain of protein G (Gronenborn et al. 1991)

Figure 2: The complete hydrogen bonding pattern within the B1 domain. The arrows are depicting hydrogen bonds (Gronenborn et al. 1993)

This complete hydrogen bonding pattern of B1 domain⁸ showed that the guest site T53 is in close interactions with I6, T44 and T55 positions which thus seemed to comprise the surrounding host environment (Figure 3) for the guest site. To provide an optimal host environment, in order to maximize the solvent accessibility of the guest site and to maintain a reasonably stable protein, a combination of mutations (I6A and T44A) was made in addition. T53 was specifically chosen⁸ as the guest site because (a) it is centrally located on the solvent exposed face of the beta-strand. This provides a homogeneous environment to study (b) its side chain extends away from the alpha-helix, which avoids any guest-helix interaction (c) its carbonyl and amide functionalities are fully involved in hydrogen bonding.

Figure 3: Illustration of the nearest residues (circled in light) of the surrounding host environment for the guest site, position 53 (circled in dark). The arrows are depicting hydrogen bonds (Gronenborn et al. 1993)

The host B1 domain undergoes reversal thermal denaturation at ~87 ^oC^{1, 9} yet its stability at physiologically relevant temperatures is not unusually high⁹. It is most stable at a pH of 5.2⁸.

All molecular modeling studies were done using the InsightII software on the Silicon Graphics O2 workstations. The structure of the wild type was obtained from the Protein Data Bank (PDB accession no. 2GB1). Residues I6 & T44 were changed to Alanines. Forcefield cff91 in the Biopolymer module was used to fix the potentials.

After trying spheres of varying radii, centered on T53, a sphere of 8Å radius was chosen to be the optimum (having 20-21 residues). The output file of residue listing was read everytime to see how many residues were there within that particular radius. Then, using the Constraints option in Discover, everything outside this sphere was frozen for all calculations, since we were only interested in investigating the local flexibility around the point of mutation, i.e., T53. The dielectric constant was set to 4, having distance dependence turned on, and the cutoff (under the Variables option in Discover) was set to a 100 since we wanted to have all the Van Der Waals interactions ON during minimization. The gradient type was set to Conjugate, for 6000 iterations with a maximum derivative of 0.1 (since we were already so close to the real structure to begin with). Under these specified constraints, these structures were energy minimized using Discover, with only the charges being ON. Each minimization took around 3 hours. Minimization generated three files- .mdf, .inp and .car, all having the object name as their prefix. These files were used further to run molecular dynamics on the energy-minimized structures.

For running dynamics the .inp file had to be altered appropriately for my system. The dynamics script was appended to the minimization script already present in the .inp file. Keeping the charges turned ON, and the Morse terms & the cross terms turned OFF, the duration of the total dynamics run was set to 5nanoseconds for the T53 mutant (GBT) and 4nanoseconds for the T53A mutant (GBA). A 10picoseconds-time period was allowe d for the system to rise to a temperature of 300K from 0K, and after which it was let to equilibrate at 300K, in both cases. A time step of 1femtosecond was used and conformers were collected after every 5picoseconds. This enabled us to collect 1000 different frames (or conformers) for GBT and 800 for GBA, each having generated after 5000 iterations. All these dynamics runs were performed under vacuum.

where a is the angle between the N-H bond vector at time t = 0 and at time t . The brackets signify the final average taken over all the conformers being considered for this calculation. This required the usage of a special script (see Appendix), Marty & I especially wrote, which took the x, y & z coordinates as the input. Various time scales were tried to come up with the one that gave the most complete result.

The total energy (in kcal) was graphed for all the conformers using the Analysis module. The backbone N-H bond length of residue 53 (T/A) was also plotted for all the conformers using the Analysis module. Both the graphs were edited using the Graph option in the same module. Ramachandran and Hydrophobicity Plots were graphed for all the 56 residues using the Homology module in splatter. These were edited using the Graph option iconized on the left-hand side panel of the screen.

A movie was also made which showed the movement of the protein, keeping the conformers fixed at the backbone N of residue 53. A script was especially written (see Appendix) to automate this process.

The secondary structure of B1 domain of the wild type (WT) protein G is as shown in Figure 1.

Figure 1. Structure of the B1 domain of the wild type protein G made using Kabsch_Sander’s method in the secondary render option. (red = alpha-helix, yellow = 4 antiparallel beta-strands, blue = turns, green = random coils).

The structures for the B1 domains of GBT and GBA were obtained from the .car file generated after minimization followed by molecular dynamics runs of 5ns and 4ns, respectively (Figure 2). All the subsequent analyses of the results were done on these structures. The structure of the protein showing the fixed sphere of residues can be seen in Figure 3.

Figure (a)
Figure (b)

Figure 2. Post dynamics ribbon structures of the two mutants (a) GBT (b) GBA having two additonal mutations, I6A & T44A common to both. Mutation at position 53 has been highlighted in red for (a) and in green for (b). It is to be noted here that the two structures, obtained after the dynamics runs, look slightly different in a few places. That may be because of the mutation at position 53.

Figure 3. CPK structure showing the residues that were fixed in blue, the ones within the 8Å sphere of interest in yellow, and the guest site 53 (T/A) in red. The ones in yellow include residues 3-9, 14-16, 27, 30, 39, 42-47, 49-56. These same residues were held fixed for both mutants in order to be consistent.

Ramachandran plots (Figure 4) were drawn to confirm if the value for the backbone phi and psi angles had been violated or not while doing the calculations for all the 56 residues in both mutants.

Figure (a)
Figure (b)

Figure 4. Ramachandran plots for (a) GBT (b) GBA. The backbone phi angles have been graphed against the backbone psi angles over the range of 56 residues. The Homology module was used for plotting.

The hydrophobicity plots (Figure 5) reflect the hydrophobicity of all the residues in both the mutants. The values look acceptable; those buried inside are fairly hydrophobic and those exposed are hydrophilic.

Figure (a)
Figure (b)

Figure 5. Hydrophobicity plots for (a) GBT (b) GBA. The extent of hydrophobicity for each of the 56 residues has been depicted.The Homology module was used to generate this.

The solvation energy calculations (Figure 6) have been depicted in the form of a bar graph. It can be seen that the total solvent accessible area of the A53 mutant is the highest and that of the T53 is the lowest. The solvent exposed surface area of residue 53 increases as we move from the WT to GBT. This makes sense as the threonine in GBT bears a hydroxyl group, which should be more exposed to the solvent than the methyl group of alanine in GBA.

Figure (a)
Figure (b)

Figure 6. Bar graphs for (a) the total solvent accessible surface area and, (b) the exposed surface area of residue 53. Both graphs show a relative comparison of values for the WT, GBT and GBA. (Plotted using Microsoft Excel).

The total energy spanned over the range of conformers (Figure 7) shows a random variation. An average value for GBT is seen to be ~115 kcal and that for GBA is ~185 kcal.

Figure (a)
Figure (b)

Figure 7. Total energies (kcal) of (a) GBT (b) GBA graphed for all the conformers. The Analysis module was used to make these.

The backbone N-H bond length of residue 53 also shows a random variation (Figure 8) over the range of conformers. The maximum range spanned by it is ~0.94-1.08Å for both the mutants.

Figure (a)
Figure (b)

Figure 8. Graphs showing the variation in the backbone N-H bond length of residue 53 for (a) GBT and (b) GBA over all conformers. Majority of the conformers in GBT is showing an average ideal value of ~1.0115Å with a standard deviation of 0.024 and in GBA is showing an average of ~1.0114Å with a standard deviation of 0.007;. (The ideal N-H bond length is known to be 1.02Å).

The consistence in the NH bond length in both mutants shows that all the conformers remain somewhat the same in this respect. So, this could mean that this length of ~1.011Å is probably ideal for NH bonds in both the mutants, irrespective of the side chain present (threonine in GBT & alanine in GBA).

The distance of the pseudo atom (which is the geometrical center of the alpha-helix) also shows a very random pattern over a range of conformers (Figure 9).

Figure (a)
Figure (b)

Figure 9. Graphs showing the variation in the distance between the geometrical centroid of the alpha-helix and the backbone N of residue 53 for (a) GBT (b) GBT. Average values observed for (a) and (b) are ~11.25 & ~10.55Å respectively. The same in the WT is seen to be 10.25Å.

The large deviation of this distance in GBT leads us to propose that the helix moves closer when an alanine side chain is introduced in position 53, in place of a threonine. This is what would have been expected in general too, since an alanine is much smaller in size than a threonine.

Some important hydrogen bond distances, distance of residue 53 from a couple of residues in the alpha-helix and a few in the adjacent beta-strands have also been calculated (Figure 10) to see how much the mutation at position 53 affects the rigidity of this protein.

Figure 10. Comparisons of some important distances have been tabulated here, to comment on the effect of mutation on a residue in one beta-strand, on the neighboring strands and the helix.

(1) Distances in the WT:

Within Beta strands-

H-Bond distances::

GB1:THR53:O-GB1:THR44:NH 2.22

GB1:THR53:NH-GB1:THR44:O 2.61

GB1:VAL54:NH-GB1:ILE6:O 2.86

GB1:ILE6:NH-GB1:PHE52:O 1.98

Distances::

GB1:THR53:C-GB1:ILE6:C 4.47

Between Beta Strand & Helix-

GB1:THR53:C-GB1:PHE30:C 9.41

GB1:THR53:C-GB1:GLU27:C 9.93

(2) Distances in the post dynamics GBT:

Wthin Beta strands-

H-Bond distances::

GB1:THR53:O-GB1:ALA44:NH 3.13

GB1:THR53:NH-GB1:ALA44:O 2.78

GB1:VAL54:NH-GB1:ALA6:O 3.91

GB1:ALA6:NH-GB1:PHE52:O 2.57

Distances::

GB1:THR53:C-GB1:ALA6:C 4.95

Between Beta Strand & Helix-

GB1:THR53:C-GB1:PHE30:C 10.43

GB1:THR53:C-GB1:GLU27:C 11.43

(3) Distances in the post dynamics GBA:

Wthin Beta strands-

H-Bond distances::

GB1:ALA53:O-GB1:ALA44:NH 3.92

GB1:ALA53:NH-GB1:ALA44:O 2.79

GB1:VAL54:NH-GB1:ALA6:O 3.22

GB1:ALA6:NH-GB1:PHE52:O 2.48

Distances::

GB1:ALA53:C-GB1:ALA6:C 4.65

Between Beta Strand & Helix-

GB1:ALA53:C-GB1:PHE30:C 10.17

GB1:ALA53:C-GB1:GLU27:C 10.15

Cluster graphs for both the mutants have also been made to draw an estimate of the number of conformers, which are somewhat alike in the whole ensemble of conformers.

Figure (a)
Figure (b)

Figure 11. Cluster graphs for both (a) GBT and (b) GBA have been made using the Analysis module. Different colors have been used for different ranges of the RMS deviations. Blue patches are seen to become dominant in both the graphs, both of which are within the range of 1.5.

The dominant blue patches in the cluster graphs drawn for both show a good conservation of a similar kind of geometry over a range of conformers. This shows that the structures obtained after dynamics had a fairly favorable geometry, which was thus conserved. These conserved structures were within the range of a RMS deviation of 1.5, which seemed to be quite good.

Cones were made (Figure 12) using the backbone N-H bonds of residue 53 of all the conformers. This was done to see how this NH bond was changing its position over all the conformers. The changing NH bond lengths could also be nicely visualized in the cone. (A graphical representation has also been shown for this bond length in Figure 8).

Figure (a)
Figure (b)

Figure 12. Cones of the backbone NH bonds of residue 53 for (a) GBT (N in purple; H in red), and (b) GBA (N in red; H in green). It is to be noted that the dynamics run durations for (a) and (b) were 5ns and 4ns, respectively. N was kept fixed at the very bottom of the cone. The varying heights of the spikes, each representing a NH bond, also showed how this bond length was changing with time. This was generated using a special perl script (see Appendix).

The S² values for residues 52, 53 and 54 for both mutants (Figure 13) were also graphed to see the change in rigidity of the backbone NH bonds with time. These graphs showed some very good results. In all of them, S² is seen to be varying a little in the beginning and is then seen to be leveling off after a very long time. This leveling off signifies that after a reasonably long time the flexibility of this bond is becoming almost constant. This supports the cluster graph results quite well. Both show that a majority of the conformers are almost the same.

Figure 13. Graphs for the S² values have been drawn over the total time of the dynamics run. (Note: The closer the value is to 1, the more rigid the bond is). These have been drawn using Microsoft Excel.

These S² graphs for GBT were drawn first. In all of them, the value of S² was dying to a plateau value after ~4ns. This meant that the structure of the protein was not changing much after that. So, since a dynamics run of only 4ns seemed to be enough to see any significant changes in the protein, the run time for GBA was reduced to 4ns.

The overall results indicate that the protein has relaxed in terms of the close packing. So it can be said that the hydrogen bonding pattern of this protein is somewhat sensitive to mutations. But within a particular mutant the structures of the conformers remain quite conserved. So, the protein remains folded throughout the dynamics run

Similar calculations will be done on a third mutant, T53M.