109

Biochem. J. (1992) 288, 109-116 (Printed in Great Britain)

Molecular-dynamics investigation of molecular flexibility in ligand binding Boryeu MAO Upjohn Research Laboratories, Kalamazoo, MI 49001, U.S.A.

The molecular flexibility of an inhibitor in ligand-binding process has been investigated by the mass-weighted moleculardynamics simulation, a computational method adopted from the standard molecular-dynamics simulation and one by which the conformational space of a biomolecular system over potential energy barriers can be sampled effectively. The bimolecular complex of the aspartyl proteinase from Rhizopus chinensis, rhizopuspepsin, and an octapeptide inhibitor was previously studied in a mass-weighted molecular-dynamics simulation; the study has been extended for investigating the molecular flexibility in ligand binding. A series of mass-weighted molecular-dynamics simulations was carried out in which libration of the inhibitor dihedral angles was parametrically controlled, and threshold values of dihedral angle libration amplitudes were observed from monitoring the sampling of the enzyme binding pocket by the inhibitor in the simulations. The computational results are consistent with the general notion of molecular-flexibility requirement for ligand binding; the freedom of dihedral rotations of side-chain groups was found to be particularly important for ligand binding. Thus the critical degree of molecular flexibility which would contribute to effective enzyme inhibition can be obtained precisely from the modified molecular-dynamics simulations; the procedure described herein represents a first step toward providing quantitative measures of such a molecular-flexibility index for inhibitor molecules that have been otherwise targeted for optimal protein-ligand interactions.

INTRODUCTION The description of structural and dynamical properties of biological macromolecules is essential for understanding many biochemical processes. By a semi-empirical potential-energy function that describes inter- and intra-molecular interactions in these biomolecular systems, it has been possible to simulate their structure and dynamics on high-speed digital computers [1-4]; molecular properties that are dependent on local potentialenergy surfaces in the conformational space can be obtained with a reasonable precision and without an excessive computing cost, given that potential-energy functions are generally accurate. For many problems in structural biochemistry and structural biology, however, it is necessary to sample the molecular conformational space extensively, and because of the large number of degrees of freedom in a biomacromolecule, it is often difficult to assess whether molecular structures in the multidimensional conformational space have been sampled sufficiently. For molecular-dynamics simulations, the difficulty is due to the time scale of numerical integration and potential-energy barriers that exist in the conformational space; both of these two types of limitations can be substantially removed in the so-called 'massweighted molecular-dynamics' simulation (MWMD) [5]. In this variant of the molecular-dynamics method, atomic particles in the molecular system are assigned masses that are higher than their respective physical values by a uniform multiplicative weighting factor, cl (co > 1); at the same time, internal degrees of freedom such as bond lengths and valence angles are strengthened with harmonic constraining potentials for maintaining the structural integrity, whereas dihedral rotations and non-bonded interactions remain unchanged. On account of its numerical stability, this approach was shown previously to have an improved overall conformational sampling in comparison with standard and conventional high-temperature molecular-dynamics simulations [6], and in contrast with sampling local regions in

the conformational space [7,8], potential energy barriers in the 'soft' degrees of freedom of biomacromolecular systems can be overcome much more readily in this methodology. The improved conformational sampling in systems ranging from short linear oligopeptides [5.6] and constrained polypeptides [9] to an inhibitor-enzyme complex [10] provided many interesting structural characteristics of these biomolecular systems. In the MWMD simulation of the complex of the aspartyl proteinase rhizopuspepsin (from the fungus Rhizopus chinensis) and its octapeptide inhibitor [10], the dihedral angle space of individual polypeptide backbone 1D and T angles of the inhibitor was found to be extensively sampled; moreover, the kinetically energized ligand molecule was able to negotiate its path from the crystallographically determined position in the complex toward the opening of the binding pocket and eventually became dissociated from the enzyme. These results suggest that molecular flexibility of the ligand molecule is essential in the ligand dissociation process, and imply that the molecular flexibility is also important in the association process. In the present paper the requirement of molecular flexibility of a ligand molecule in the binding process has been investigated in an extension of the previous MWMD study of the complex of rhizopuspepsin and the octapeptide inhibitor D-His-Pro-Phe-His-PheT[CH2NH]Phe-Val-Tyr. Specifically, rotational freedom of and T dihedral angles of the polypeptide backbone and X angles of the amino acid side chains of the inhibitor are parametrically controlled during MWMD simulations, and the effect of such conformational restriction on the ability of the ligand molecule to become dissociated from rhizopuspepsin is monitored. It is found that, when the molecular flexibility of the octapeptide inhibitor is reduced by conformational restrictions on dihedral rotations of the polypeptide, the rigidified molecule could not reach the opening of the binding pocket to become dissociated. Moreover, threshold values of the molecular flexibility required for the octapeptide ligand were quantitatively determined.

Abbreviations used: MWMD, mass-weighted molecular-dynamics; RMS, root Vol. 288

(I

mean square.

B. Mao

110 (a)

1992

III

Molecular flexibility in ligand binding Results from MWMD simulations in which atoms in the enzyme are allowed to move are also presented; the implication of molecular flexibility in the study of protein-ligand interactions is. discussed. METHODS All molecular-dynamics simulations were carried out within the macromolecular program package CHARMM [11]. In MWMD simulations, two modifications are added to the otherwise standard molecular-dynamics procedure: the uniform weighting of all atomic masses by a numerical factor, and the installation of harmonic constraining potentials on rigid internal co-ordinates [5]; both features of the MWMD method are implemented with standard procedures and options of CHARMM. Previous applications of this methodology indicated that the underlying potentials for unconstrained components such as dihedral rotations are not affected by these modifications [5]; physical characteristics of MWMD simulations was also described. For MWMD studies of molecular flexibility in enzyme-ligand complexes, a starting structure of rhizopuspepsin complexed with the octameric polypeptide inhibitor D-His-Pro-Phe-HisPheT[CH2-NH]Phe-Val-Tyr was prepared from the crystallographically determined atomic co-ordinates of the complex [12] as described previously [10]; this structure, with the 'flap' of rhizopuspepsin (residues Gly70-Gly85 inclusive) poised in the crystallographically observed 'closed' position, is referred to as the 'closed-flap complex'. From the closed-flap complex, a second structure is then generated in which atoms in the flap are shifted toward the solvent, while the inhibitor is kept in its original position in the closed-flap complex; the position shift vector for each atom in the flap is the vector difference between its position in the complex and its position in the structure of ligand-free rhizopuspepsin [13]. This structure, referred to as the 'open-flap complex', henceforth, represents a hypothetical penultimate complex of the ligand-free enzyme and the inhibitor after the latter has entered the binding site but before the flap closes down to form the final complex. The comparison of crystallographically determined structures of rhizopuspepsin with and without the bound inhibitor, and the comparison of the closed-flap and open-flap complexes are shown in Fig. 1. For the reduced peptide bond between residues Phe5m and Phe6 of the inhibitor, the atom type of the nitrogen atom in Phe6 is changed to that of a terminal nitrogen atom in CHARMM [11]. This in turn changes the bond type of the reduced peptide bond to one without the partial double-bond character found in normal peptide bonds; the dihedral rotation force constant for this bond type is 2510 J (600.0 cal)/mol rad2. These changes are reasonably accurate within the first order of approximation in that, in a normal molecular-dynamics simulation of the complex starting from its crystal structure [10], this nitrogen atom has a position shift from its starting position and a root-mean-square (RMS) position fluctuation comparable with those of backbone atoms of the four-residue fragment from His4 to Val' of the inhibitor [0.029 nm (0.29 A) and 0.04 nm (0.40 A) respectively, compared with 0.028 nm (0.28 A) and 0.039 nm (0.39 A) for the backbone average]; the RMS amplitude of the

libration of the reduced peptide bond is 150 compared with about 9° for a normal peptide bond. In the previous study of rhizopuspepsin-inhibitor complex [10], the (F and T dihedral angles of the inhibitor polypeptide backbone, as well as the X angles of the amino acid residue side chains, have the complete freedom to rotate (whereas bond lengths, bond angles, amide planes and improper torsional angles are constrained as a part of the standard MWMD protocols). In order to investigate the role of the molecular flexibility of the inhibitor in the ligand-binding process, the degree of rotational freedom of all proper dihedral angles of the octapeptide (i.e., backbone (D and I, and sidechain X dihedral angles) is controlled parametrically by constraining potentials placed on these internal co-ordinates; in the present study, constraining potentials for all proper dihedral rotations have the same force constant, and this too is implemented using standard procedures in CHARMM. In a separate set of MWMD simulations, only (D and T angles are constrained; i.e., the rigidity of only the polypeptide backbone is controlled, whereas side-chain groups are free to rotate. As in the previous study [10], the mass-weighting factor, w, is 30.0, atomic velocities are assigned from the Gaussian distribution at 600 K, and rigid degrees of freedoms are constrained with harmonic potentials which maintain the geometry of the covalent and other planar and chiral structures (force constants of the harmonic potential being 4184 MJ/mol nm2 (10000 kcal/mol A2) for bond lengths, 16736 kJ (4000 kcal)/mol rad2 for valence angles, and 1674 kJ (400 kcal)/mol rad2 for amide planes and improper dihedral rotations). In most simulations reported here, protein atoms (in either the open-flap or the closed-flap complex) are held fixed in space as a first-order approximation. For comparisons, several simulations are repeated with mobile protein atoms. The effects of mobile protein atoms on the molecular flexibility of the inhibitor are simulated as follows. Amino acid residues of rhizopuspepsin in each of the two complexes are divided into two groups: a core region and a shell region; those residues with one or more atoms located within the 1.5 nm (15 A) distance of any inhibitor atom in the complex form the core region, and the remaining residues constitute the shell region. Atoms in the shell region are held fixed, whereas each atom in the core region is allowed to move within a three-dimensional harmonic potential [force constant 209.2 MJ/mol nm2 (500 kcal/mol A2) and centred at its fixed position in either complex]. Protein atoms in the core region are also mass-weighted, with w equal to 30.0. Typically, in any one of these simulations, the centres of atomic-position distributions for core atoms are located essentially at their respective fixed crystallographic positions [shifted by an average of only 0.00 14 nm (0.0 14 A)], and the average of RMS position fluctuations for core atoms is about 0.027 nm (0.27 A). In simulations with fixed protein atoms, there are 92 atoms from the inhibitor for molecular-dynamics integration, and in simulations with mobile core atoms, there are 1676 atoms for integration; the number of atom pairs for non-bonded interactions increases from 6728 pairs to 96061 pairs. The total central processor time for a 20 ps simulation increases 6-fold overall, from 45 min (on an IBM 3090J processor) to 320 min. For characterizing the molecular flexibility of the inhibitor molecule, the time series for the values of a backbone dihedral -

-

Fig. 1. Comparison of structures of the rhizopuspepsin with and without the bound polypeptide inhibitor D-His-Pro-Phe-His-PheVICH2-NHIPhe-Val-Tyr (a) Crystallographically determined structures, in C' tracing, of the complex (continuous line [12]) and the free enzyme (broken line [13]). The C" atom of Gly"' in the flap has the largest position shift, of 0.11 nm (1.10 A). A 0.5 nm (5 A) bar is shown on the upper-right corner for reference. (b) A close-up view of the inhibitor-binding site from (a). The C' atom of Phe5 of the inhibitor, and C` atoms of the catalytic residues (Asp35 and Asp218) are marked with small squares. (c) The inhibitor-binding site of the closed-flap complex (continuous lines) and the open-flap complex (broken lines). Note that positions of the two terminal residues of the inhibitor have been generated. See the Methods section for description of the closed-flap and the open-flap complexes.

Vol. 288

B. Mao

112 Table 1. Fluctuations of 14 backbone dihedral angles of the octapeptide inhibitor D-His-Pro-Phe-His-Phe'lICH2-NHJPhe-Val-Tyr during MWMD simulations of rhizopuspepsin-inhibitor complexes The fluctuations are indicated by the widths (in degrees) of the distributions of numerical values of the dihedral angles. Polypeptide backbone dihedral angles (D and T, except (i Of D-His5 and 02 of Pro2) are labelled sequentially from the DHis' residue on the Nterminus to the Val8 residue on the C-terminus. As a reference, dihedral angle fluctuations in the MWMD simulation in which proper dihedral angles were not constrained [10] were listed in Column 1. Columns 2 and 3 show dihedral fluctuations of the inhibitor in 200 ps simulations of the closed-flap complex and the open-flap complex respectively, in which the dihedral angles were constrained at 1674 kJ (400 kcal)/mol rad2. Averages from the 14 dihedral angles are given in parentheses. -

Value (0) Dihedral angle

Unconstrained

Closed-flap

Open-flap

2 3 4 5 6 7 8 9 10 11 12 13 14

323 359 356 302 309 321 305 296 312 360 358 352 360 360 (334)

86.03 80.35 80.40 78.12 71.38 78.12 70.21 77.01 78.05 75.96 83.37 96.05 89.99 125.09 (83.58)

98.89 82.63 79.10 82.26 75.22 74.54 68.29 82.32 82.63 76.82 90.85 80.10 94.93 121.19 (84.98)

Average...

angle (.1 or T) is computed from Cartesian co-ordinate sets saved during a dynamics simulation. From these values, the width of the distribution is obtained (Table 1); the collection of these widths for all backbone dihedral angles is a measure of the polypeptide backbone flexibility in that a larger average width reflects a larger sampled region in the multidimensional space and a larger degree of freedom for torsional motion, and hence a higher molecular flexibility (see the Appendix for a further discussion).

As a reference for the reduced molecular flexibility of a conformationally constrained polypeptide, standard moleculardynamics simulations (without mass weighting and constraining potentials for covalent structures) were also carried out for a hexapeptide (Cys-Tyr-Phe-Gln-Asn-Cys), in the linear and the disulphide-linked cylic forms [9]. These standard simulations at 300K were performed by using the protocol described previously [5,6]. RESULTS

In the previous mass-weighted molecular-dynamics study of the rhizopuspepsin-inhibitor complex [10], proper dihedral angles of the polypeptide inhibitor (D and T angles of the backbone and X angles of side-chain groups) were shown to have nearly complete freedom for rotation that is not hindered by the presence of atoms in the ligand binding pocket of the enzyme (Table 1, column 1). When these dihedral angles are constrained with harmonic potentials [force constant 1674 kJ (400 kcal)/ mol * rad2], the molecular flexibility for dihedral angle fluctuation is reduced: the average width of dihedral distributions for the 15 backbone dihedral angles of the inhibitor is reduced from 3340 to 83.60 and 85.00 in the closed-flap and the open-flap complexes respectively (Table 1, columns 2 and 3). Owing to this reduced molecular flexibility, the inhibitor molecule remains in the binding site of the closed-flap (Fig. 2) and the open-flap complexes (Fig. 3a), whereas, within the first 200 ps of the previous simulation, the inhibitor molecule with torsional freedom was able to locate the opening of the binding pocket and negotiate the path to the solvent space and become dissociated from the complex [10]. As the force constant of the harmonic potentials for constraining dihedral rotations is lowered from 1674 kJ (400 kcal)/mol rad2, the octapeptide inhibitor regains flexibility; Table 2 summarizes the MWMD simulations that were performed to characterize the molecular flexibility of the ligand in the rhizopuspepsin-inhibitor complex. The minimal requirement of molecular flexibility for the ligand is the threshold value of the dihedral libration amplitude of the energized molecule for escaping the enzyme 268° for the closed-flap complex binding pocket, i.e. (simulations A-C), and 1710 for the open-flap complex (simulations E-G). The larger torsional freedom required of the inhibitor for escaping the closed-flap complex (2680 as against 1710) is as expected from the fact that the binding pocket (and thus the opening of this pocket) in this complex is smaller than

Fig. 2. Position distribution of Cx of Phe5 of the inhibitor in the binding pocket of the closed-flap complex in MWMD simulation The thick lines represent the polypeptide fragments from rhizopuspepsin that constitute the binding pocket; this is a slightly magnified view from Fig. l(c). The positions of the C0 atom of Phe5 were obtained from the MWMD simulation of the closed-flap complex, with a force constant of 1674 kJ (400 kcal)/mol rad2 for dihedral angles 1, I and X (simulation A in Table 2); positions at 0.2 ps intervals from 20 to 200 ps are displayed.

1992

Molecular flexibility in ligand binding

113

(a)

(b)

Fig. 3. Position distribution of C' of Phe5 of the inhibitor in the binding pocket of the open-flap complex in MWMD simulations Note that the flap to the left of the inhibitor is shifted toward the solvent relative to its position in the closed-flap complex in Fig. 2. Position distributions of the C' atom are obtained from the MWMD simulations of the open-flap complex (with a force constant of 1674 kJ (400 kcal)/mol * rad2 for dihedral angles b, T and X)- (a) Protein atoms are fixed (simulation E in Table 2). Atomic positions are displayed at 0.2 ps intervals from 20 to 200 ps. (b) Protein atoms are allowed to move (simulation E' in Table 2). Atomic positions are displayed at 0.2 ps intervals from 20 to 60 ps.

that of the open-flap complex; a more flexible structure (with larger torsional librations) is required of the inhibitor in order for it to negotiate the path through the opening of the binding

pocket. For comparison, the molecular flexibility of the linear and the disulphide-linked cyclic forms of a hexapeptide (Cys-Tyr-PheGln-Asn-Cys) at room temperature is computed from standard molecular-dynamics simulations of these two molecules, and the results are shown in Table 3. Owing to the constraint of the intramolecular disulphide bond in the cyclic form of the hexapeptide, the molecular flexibility is reduced, as indicated by the smaller widths of the distributions of numerical values of individual dihedral angles. This comparison indicates that intramolecular linkages generally reduce the conformational flexibility of a molecule. The requirement of molecular flexibility can also be studied when only the backbone dihedral angles, 1D and T, are constrained. In the closed-flap complex, the additional molecular flexibility acquired from rotational freedom of side-chain groups does not allow the inhibitor to escape (simulation D in Table 2). The open-flap complex, however, has a larger binding pocket Vol. 288

and a larger opening of the pocket compared with the closed-flap complex, and the flexibility of side-chain groups allows the inhibitor to leave the binding site more readily, even when fluctuations of the polypeptide backbone are reduced to an average of 36.80 (simulation I). The flexibility of the backbone needs to be further reduced before it, too, becomes too rigid for the inhibitor to leave the binding pocket (simulation J in Table 2 and Fig. 4). The effects of mobile protein atoms on molecular flexibility were also studied in MWMD simulations in which enzyme atoms in the binding pocket are allowed to move within atom-centred constraining potentials (simulations D', E', I' and J' in Table 2). The comparison of these simulations with corresponding ones with fixed protein atoms indicates that the dihedral fluctuations of the inhibitor are somewhat reduced, owing to more frequent collisions with mobile protein atoms (i.e., the size of the ligand binding pocket is effectively smaller, owing to the motion of mobile protein atoms); the requirement of molecular flexibility of the inhibitor for dissociation is generally not affected by mobile protein atoms. Fig. 3(b) shows the position distribution of Cm of Phe5 of the inhibitor obtained from simulation E'.

B. Mao

114 Table 2. MWMD simulations for investigating molecular flexibility in ligand binding For each simulation, the force constant for dihedral angle constraining potentials (in kcal/mol *rad2; note 1 kcal 4.184 kJ), the average backbone dihedral angle fluctuation (0) (e.g. those in the last row of Table 1) and S.D. (in parentheses) are listed; also indicated for each simulation is whether the inhibitor remains bound (b) in the binding pocket. Simulations A-D' are based on the closed-flap complex, and E-J' the open-flap complex. Results from simulations in which all proper dihedral angles (1', ' and X) are constrained are listed under (a); the results from constraining only backbone angles SD and ' are listed under (b). In simulations D', E', I' and J', protein atoms are allowed to move (see the text for details); these simulations were carried out for 60 ps, whereas all other simulations (in which protein atoms are fixed) were carried out for 200 ps.

(a) SD, ' and X constrained

(b) I and ' constrained

Angle Force Angle Simul- Force ation constant fluctuation Bound? constant fluctuation Bound? A B C D D' E E' F G H I I' J J'

400.0 25.0 16.0

83.6 (13.3) 233.9 (47.9) 268.0 (54.6)

b b -

400.0 82.4 (11.4) 400.0 71.2 (8.6) 400.0 400.0 80.0 50.0 2000.0

85.0 (12.7) 70.0 (8.9) 142.2(19.9) 170.7 (29.7) 39.8 (5.7)

b b

b b b -

b

2000.0 2000.0 6000.0 6000.0

36.8 (4.7) 32.8 (4.7) 25.1 (3.1) 20.0 (1.9)

-

b b

DISCUSSION In the previous MWMD study of the rhizopuspepsin-inhibitor complex, the backbone dihedral angles ((F and I) and side-chain dihedral angles (X) have the complete freedom for rotation, reflecting the molecular flexibility of a normal linear polypeptide molecule. The 3600 range of each backbone dihedral angle ((D or T) is essentially fully accessible [10], indicating that the molecular

Table 3. Fluctuations of backbone dihedral angles of the hexapeptide CysTyr-Phe-Gln-Asn-Cys during molecular-dynamics simulations at 300K Polypeptide backbone dihedral angles (D and T, except IDl of Cys1) are labelled sequentially from the Cys' residue on the N-terminus to the Cys6 residue on the C-terminus. Columns 1 and 2 are dihedral fluctuations from 100 ps simulations of the linear molecule and the disulphide-linked cyclic molecule respectively. Averages from the 11 dihedral angles are given in parentheses. Angle fluctuation (0)

Angle

(1) Linear

(2) Cyclic

1 2 3 4 5 6 7 8 9 10 11

304.05 139.37 270.58 252.71 172.03 149.44 268.08 97.65 355.50 171.32 284.90

90.80 95.97 149.27 99.82 118.21 113.17 115.49 74.97 184.60 124.10 355.69

Average... (224.15)

(138.37)

flexibility is not significantly hindered by intermolecular interactions within the ligand binding pocket of the enzyme. Within the first 200 ps of the MWMD simulation, the kinetically energized inhibitor molecule was able to locate the opening of the binding pocket to the solvent by assuming a proper conformation (or conformations) to fit through the opening of the binding pocket and became dissociated from the enzyme. These results suggested that molecular flexibility of the inhibitor molecule is essential for the observed dissociation process; a rigid molecule would not have a sufficient conformational flexibility in the negotiation of a path through the opening of the binding pocket. The results also implied that the molecular flexibility of a ligand is essential for- the association process as well, since a conformationally restricted molecule (e.g. due to intramolecular linkages) might not have proper conformations in its somewhat limited conformation 'repertoire' that would allow it to fit

Fig. 4. Position distribution of C" of Phe' of the inhibitor in the binding pocket of the open-flap complex in MWMD simulation Protein atoms are fixed, and only dihedral angles SD and T are constrained [force constant 25.1 M J (6000 kcal)/mol rad']; side-chain dihedral angles X are not constrained (simulation J in Table 2). Cx positions from the MWMD simulation are displayed at 0.2 ps intervals from 20 to 200 ps.

1992

115

Molecular flexibility in ligand binding through the opening from either the binding site (namely the dissociation process) or the solvent space (namely the association or binding process). Taking that study as a reference, one could then study the role of molecular flexibility in ligand binding explicitly by carrying out MWMD simulations in which the rotational degrees of freedom of the inhibitox are parametrically controlled by constraining potentials. When the polypeptide inhibitor becomes too rigid (i.e., the fluctuations of dihedral angles of the molecule fall below some threshold value), it loses the molecular flexibility for negotiating the 'reaction path' in entering or leaving the binding pocket of its target enzyme. From such MWMD simulations, these minimal requirements for molecular flexibility have been determined: fluctuations in dihedral angles should be greater than 2680 and 171° respectively for the closed-flap and the openflap complexes. For the open-flap complex, the larger binding site and its larger opening allow a more rigid molecule to enter or leave compared with the closed-flap complex (i.e. an average torsional libration of 1710 as against 2680). Whereas the physical ligand-binding process may be expected to involve a complex sequence of conformational changes in both the inhibitor and the enzyme, the closed-flap and the open-flap complexes utilized in the present study are based on observed rhizopuspepsin structures and represent two limiting cases for studying ligand-binding processes. The actual requirement of molecular flexibility is likely to lie between the two limits determined in the present study. Moreover, the binding process could involve large-scale interdomain motions on the part of the enzyme, which have been observed recently in other members of the aspartyl-proteinase family [14] (but only imperceptibly, if at all, in rhizopuspepsin; see, e.g., Figs. 1 and 2). If such concerted skewing motions are present in an enzyme, which are expected to take place on a time scale longer than that of the encounter process, the molecularflexibility requirement for the inhibitor binding would be further complicated (and expected to be made different) by the skew of globular domains. In the theoretical simulation study reported here, the conformational restraints on the polypeptide flexibility were achieved by placing constraining potentials on proper dihedral angles of a polypeptide molecule; chemical modifications, on the other hand, have been successful in introducing 'real' conformational constraints into polypeptide molecules [15,16,17]. As shown in Table 3, a typical side-chain cyclization (e.g. disulphide linkage between two Cys residues) considerably reduces the molecular flexibility, and other peptidomimetic chemical modifications, such as incorporation of multi-ring moieties, are expected to reduce the molecular flexibility even further. The results from the present studies suggest that a minimum molecular flexibility must be retained in a ligand molecule or it will become too inflexible to undergo conformation changes required in the normal ligand-binding process. To the extent that the inhibitor molecule samples the interior of the binding pocket in the present simulations, the results on molecular flexibility requirement for reaching the binding pocket opening are expected to remain essentially unchanged when explicit solvent molecules are included; the eventual escape of the inhibitor from the binding pocket is also expected to be unaffected by the fluid state of solvent molecules that might be present in the simulations. Although the importance of molecular flexibility in guest-host molecular recognition is generally recognized (see, e.g., Roberts [18]), the dynamic roles of molecular flexibility in protein-ligand .interactions within the binding site have not been accounted for explicitly in current theoretical and computational methodologies [19-21], mostly owing to the lack of suitable techniques for sampling the molecular conformational space. For example, in Vol. 288

the exploration of the three-dimensional structure of an enzyme to search for ligand structures that allow optimal interactions with the binding pocket, discrepancies between experimental and predicted ligand-binding activities were observed, and it has been suggested that inappropriate binding geometry could have been concluded for some target ligands [22]; it is possible that this type of situation had resulted not only from insufficient surveying of the enzyme binding pocket, but also from not taking an explicit account of the required molecular flexibility of the inhibitor: the inhibitor may actually bind in an unintended position or conformation, as it cannot negotiate its conformation for reaching the anticipated (and desirable) binding site, or the rates of dissociation and association are affected by reduced molecular flexibility. A related note is that a ligand molecule optimized strictly for its complementarity with the enzyme binding pocket may not have sufficient dynamical flexibility within the binding site that is critical for its entering the binding site and reaching the final binding conformation as noted here; moreover, the different molecular-flexibility requirements for the inhibitor in the simulations of open-flap and closed-flap complexes indicate that the conformation of the binding pocket of the enzyme in the ligand-free state must also be considered. The methodology described here represents a first step toward incorporating such information into structure-based molecular studies. In summary, the MWMD method has been shown to be an effective procedure for sampling conformation space, especially for regions in the conformational space that are separated by potential energy barriers [5]. The present study has taken advantage of the capability of MWMD for sampling flexible degrees of freedom in a biomolecular system and investigated the molecular flexibility in ligand binding. From studying the dissociation of an inhibitor from the binding site of its target enzyme, it is possible to infer the molecular flexibility required for the association process, and a quantitative measure of the required molecular flexibility for the ligand binding process can be determined precisely. This information could complement direct simulations of the binding process such as the Brownian dynamics simulation [23]. In MWMD simulations, the higher region of the potential-energy surface is sampled, owing to the larger kinetic energy from mass weighting [5]; as a result, the opening of the binding pocket experienced by the inhibitor molecule in the actual binding process is expected to be somewhat smaller than that indicated in these simulations. The measure of molecular flexibility derived from the MWMD simulations thus represents a lower bound of the required flexibility for ligand binding. Finally, the results presented here (namely Table 2) further demonstrate the robustness of the MWMD simulation approach for conformational studies of biomolecular systems. REFERENCES 1. Scheraga, H. A. (1981) Biopolymers 20, 1877-1899 2. McCammon, J. A. & Harvey, S. C. (1987) Dynamics of Proteins and Nucleic Acids, Cambridge University Press, Cambridge 3. Brooks, C. L., III, Karplus, M. & Pettitt, B. M. (1988). Proteins: A Theoretical Perspective of Dynamics, Structure, and Thermodynamics, John Wiley and Sons, New York 4. Karplus, M. & Petsko, G. A. (1990) Nature (London) 347, 631-639 5. Mao, B. (1991) Biophys. J. 60, 611-622 6. Mao, B. & Friedman, A. R. (1990) Biophys. J. 58, 803-805 7. Bennett, C. H. (1975) J. Comput. Phys. 19, 267-279 8. Pomes, R. & McCammon, J. A. (1990) Chem. Phys. Lett. 166, 425-428 9. Mao, B., Maggiora, G. M. & Chou, K. C. (1991) Biopolymers 31, 1077-1,086 10. Mao, B. (1991) Biophys. J. 60, 966-973 11. Brooks, B. R., Bruccoleri, R. E., Olafson, B. D., States, D. J., Swaminathan, S. & Karplus, M. (1983) J. Comp. Chem. 4, 187-217

B. Mao

116 12. Suguna, K., Padlan, E. A., Smith, C. W., Carlson, W. D. & Davies, D. R. (1987) Proc. Natl. Acad. Sci. U.S.A. 84, 7009-7013 13. Suguna, K., Bott, R. R., Padlan, E. A., Subramanian, E., Sheriff, S., Cohen, G. H. & Davies, D. R. (1987) J. Mol. Biol. 196, 877-900 14. Sali, B., Veerapandian, B., Cooper, J. B., Moss, D. S., Hoffman, T. & Blundell, T. L. (1992) Proteins 2, 158-170 15. Al-Obeidi, F., Hadley, M. E., Pettitt, B. M. & Hruby, V. J. (1989) J. Am. Chem. Soc. 111, 3413-3416 16. Saragovi, H. U., Fitzpatrick, D., Raktabutr, A., Nakanishi, H., Kahn, M. & Greene, M. I. (1991) Science 253, 792-795 17. Hruby, V. J., Kazmierski, W., Kawasaki, A. M. & Matsunaga, T. 0. (1991) in Peptide Pharmaceuticals (Ward, D., ed.), pp. 135-184, Elsevier, New York

18. Roberts, G. C. K. (1991) in Host-Guest Molecular Interactions: From Chemistry to Biology (Chadwick, D. J. & Widdows, K., eds.), pp. 169-182, John Wiley and Sons, New York 19. DesJarlais, R. L., Sheridan, R. P., Dixon, J. S., Kuntz, I. D. & Venkataraghvan, R. (1986) J. Med. Chem. 29, 2149-2153 20. Howard, A. E. & Kollman, P. A. (1988) J. Med. Chem. 31,1669-1675 21. Cohen, N. C., Blaney, J. M., Humblet, C., Gund, P. & Barry, C. D. (1990) J. Med. Chem. 33, 884-894 22. Weber, A. E., Halgren, T. A., Doyle, J. J., Lynch, R. J., Siegl, P. K. S., Parsons, W. H., Greenlee, W. J. & Patchett, A. A. (1991) J. Med. Chem. 34, 2692-2701 23. Northrup, S. H., Allison, S. A. & McCammon, J. A. (1984) J. Chem. Phys. 80, 1517-1524

APPENDIX This discussion addresses the question of obtaining a quantitative measure of the sampling of the conformational space of a polypeptide molecule. In our example, the molecular backbone conformation is specified by the 14 backbone (b,T) angles (Table 1 of the main paper). Appendix Fig. lA shows that the 180

-018

180

-180

U5(o) Fig. 1A. Distributions of numerical values of the backbone IF angles of residues Phe5 and D-His' Points from simulation A in Table 2 are drawn in small filled dots, and those from simulation B in small crosses. For simulation A, the width of distribution of 's values is -780, and -88° for 1; for simulation B, the widths are 1830 and 2330 respectively. -

-

numerical values of these dihedral angles in a dynamics trajectory are normal distributions, not only for single dihedral angles, but also in the two-dimensional space of a pair of them. The sampling of the conformational space is more precisely described by an approximately equiaxial elliptical distribution. For each of the

two distributions in Appendix Fig. IA, a 2 x 2 matrix A can be constructed, of which the two elements in the first row are the products AT1 AP1 and AT1 A'T5, each summed over the trajectory, and in the second row are the two sums of products AT5 *T1 and A'5 AT5; AT1 and AT5 are the values of the two dihedral angles about their respective averages. Similar to the treatment of a moment of inertia matrix in classical mechanics, eigenvectors and eigenvalues from diagonalization ofthe matrices are then the major axes, and their lengths, of the ellipses that characterize the distributions. For both simulations A and B the axes are nearly coincident with 5 and 1axes; the axial lengths of the ellipse are 19.8° (T,5) and 23.80 (T1) for simulation A, and correspondingly 54.4° and 72.00 for simulation B; they are approximately a quarter of the full widths, since the distributions are skewed toward the origin. That these distributions have approximately equal axial lengths in the two perpendicular directions indicates that the motions of the two dihedral angles are not correlated in any fashion. These quantitative features hold for other pairs of dihedral angles; in fact, the procedure can be carried out to the multi-dimensional space spanned by 14dimensional vectors whose components are the numerical values of the 14 dihedral angles in Table 1 of the main paper. For both simulations A and B, the 14-dimensional distributions are again approximately spherical, with an average of 21.9° (variance 4.20) for the axial lengths of the distribution for simulation A, and an average of 62.60 (variance 22.50) for the 14 axial lengths of the distribution for simulation B; these average axial lengths are again about a quarter of the corresponding widths in Table 2 of the main paper. Although the above variances indicate that the distributions are only approximately spherical, and that there is some degree of correlation among motions of the dihedral angles, ellipsoid analysis shows an approximately 3-fold increase in the sampling of each dimension from simulation A to B, the same amount of increase reported in Table 2 of the main paper according to the increase of the widths for the two simulations. Thus the widths of the distributions of numerical values listed in Table 2 of the main paper are reasonably accurate measures of the sampling of the polypeptide conformational space in the present context; more detailed analyses of the ellipsoidal representation (e.g. the total volume and the characteristics of the ellipsoidal axes) are beyond the immediate aim of the present paper.

Received 9 March 1992/7 May 1992; accepted 29 May 1992

1992

Molecular-dynamics investigation of molecular flexibility in ligand binding.

The molecular flexibility of an inhibitor in ligand-binding process has been investigated by the mass-weighted molecular-dynamics simulation, a comput...
1MB Sizes 0 Downloads 0 Views