Comprehensive Protein-Ligand Interaction Analysis: Fragment Molecular Orbital Calculation on the Complexes of Human Protease Renin and its Inhibitors

Abstract

To elucidate the molecular recognition mechanism between renin and its inhibitor, we analyze intermolecular interaction energies based on fragment molecular orbital (FMO) calculations for twenty different complexes of various inhibitors. We discuss a relationship between the experimental activity value of inhibitors and the calculated binding energy, and clustering analyses for inter-fragment interaction energies (IFIEs) between the inhibitor and the amino acid residues in renin. We estimated the sum of IFIEs as binding energy between an inhibitor and renin, and found that the calculated binding energies have a relatively strong correlation (R² = 0.73) with the experimental IC₅₀ values of each inhibitor. The high-activity correlation between the calculated and experimental values can lead to predicting the effects of drugs and the activity value of new compounds. In addition, we carried out a detailed interaction energy analysis between inhibitors and the amino acid residues in renin, and performed clustering of inhibitors not only by their structure/binding mode but also by the characteristics of interaction, such as energy values, energy patterns, and interacting amino acid residues. As a result, we found that the difference in interaction due to a slight difference in structure, such as the addition/replacement of a single atom/functional group, can be related to the difference in IC₅₀ values. Consequently, the inhibitors in the finally classified clusters tend to show the same order of IC₅₀ value. These results indicate that the structure and the activity value of inhibitor are related to each other through the interaction between an inhibitor and relevant amino acid residues, and suggest that it is certainly possible to predict the IC₅₀ values of inhibitors. Therefore, we consider that our FMO-IFIE analysis would be a useful method to contribute toward innovative drug discovery.

1. Introduction

The protein is one of the most important carriers on life activities. Then, the biological functions such as enzyme reactions, drug actions and expression of taste and smell, which are exhibited by proteins are essentially due to the electron-electron interaction between active substances such as organic compounds or nanomaterials and protein or between proteins. Thus, in order to understand and elucidate the functions of proteins, it is necessary to know the electronic structure of a protein or that of a complex of active substance and protein as accurately as possible. For that purpose, an all-electron quantum chemical calculation that can accurately describe the behavior of electrons is essential.

We have focused our attention on the fragment molecular orbital (FMO) calculation [1-3] for the complexes of proteins with ligands to elucidate the function of proteins by analyzing the interaction between a protein and a ligand [4-7]. In this study, we focused on a system with multiple co-crystal structures of the complex and investigated whether activity could be predicted using FMO calculations based on X-ray structures. We selected a human protease, renin, as a target protein. Renin is an aspartic protease which exists in human kidney and contains 340 amino acid residues [8]. The active site amino acid residues of renin are Asp38 and Asp226. In Figure 1, we illustrate the structure of renin and the active site amino acid residues, Asp38 and Asp226.　Importantly, renin is the first component of the renin-angiotensin system in hypertension. In this paper, we perform the FMO calculations at FMO2-MP2/6-31G* level on the twenty complexes of renin and its inhibitors. For the FMO calculations, we selected twenty different complexes. Figure 2 shows the selected inhibitors with their PDB-ID and the activity value on 50% inhibitory concentration (IC₅₀) obtained from ChEMBL [9]. Here, complexes of 2G24 (No.19) and 2IKO (No.20) PDB entries have the same inhibitor.

Figure 1. Structure of renin (left side) and its active site amino acid residues Asp38 and Asp226 (right side)

Color code: green C, red O, blue N, white H, and yellow S.

Figure 2. Data set of inhibitors with IC₅₀ value

The primary purpose of the present study is to predict the activity values of individual inhibitors on the basis of an FMO calculation. Additional aims are to clarify relationships among the structures of each inhibitor, the activity values, IC₅₀, of each inhibitor, and the calculated energies such as interaction energy and binding energy for the complexes and to elucidate the molecular recognition mechanism of renin and its inhibitors. Through a series of work, our goal is to obtain knowledge on predicting the effects of drugs and the degree of side effects in the future.

In the next section, we describe the procedure of the biomolecular quantum chemical calculation including an FMO calculation to obtain the interaction and binding energies between renin and its inhibitors. In Section 3, we show our calculation results and discuss a relationship between the binding energy of complexes and IC₅₀ values, detailed analysis of interaction energy between inhibitors and the amino acid residues in renin, and clustering of inhibitors. In the last section, we conclude the paper by summarizing the whole picture of this system.

2. Computational Method

2.1 Molecular modeling

As mentioned above, we selected twenty different renin-inhibitor complexes for the FMO calculation. The X-ray crystallographic structures of complexes were obtained from Protein Data Bank (PDB) [10]. In this study, all crystallization water was removed. For reconstruction of renin, since residues surrounding missing residues tend to have a high temperature factor (B factor), two residues before and after each missing residue were deleted. The missing atoms and missing residues were complemented by homology modeling with SWISS MODEL [11-15]. Then, The N- and C-terminus of renin were in their native charged forms, NH₃⁺ and COO^–, respectively. For renin, hydrogen atoms were added by using AMBER11 software [16]. For inhibitors, hydrogen atoms were added by using GaussView. The charge of inhibitors was assigned by AM1-BCC of Antechamber implemented in AMBER11 [17]. In order to neutralize the system, counterions were placed. Subsequently, we carried out the energy minimization by a molecular mechanics calculation using AMBER11 with the ff99 and GAFF force fields for renin and inhibitors, respectively [18, 19]. The complemented residues/atoms, cysteine with disulfide bond, and added hydrogen atoms were optimized concurrently. Since the protonation states of the active site residues Asp38 and Asp226 are still under discussion and have not been clarified, we have, as an initial step in our research, assumed that Asp38 and Asp226 are not protonated in the present calculations. Furthermore, in this work, we performed the calculations with all inhibitors are electrically neutral because of following two reasons: (1) the charge states of the inhibitors are unclear, and (2) slight structural changes in inhibitors or amino acid residues have a significant effect on the electrostatic interaction between the amino acid residues in renin and the charged inhibitors.

2.2 FMO calculation

After energy minimization, counterions were removed and the resulting complexes were subsequently subjected to an FMO calculation. All FMO calculations were carried out at the FMO2-MP2/6-31G* level with the Cholesky decomposition approximation [20] by using the MIZUHO/ABINIT-MP program [21] on the K computer. The K computer was a Japanese supercomputer developed by RIKEN and Fujitsu. It was completed in 2011 and was one of the fastest supercomputers in the world at the time. The K computer was used for research in science, medicine, and the environment. The name "K" comes from the Japanese word "kei" (京), meaning 10 quadrillion, which refers to its target performance of 10 quadrillion calculations per second. In our FMO calculations, each amino acid residue in renin was treated as a single fragment. On the other hand, the inhibitors were divided into some fragments. Each system was divided into 338 to 340 fragments in the present calculation. All FMO calculations were performed under vacuum conditions. We obtained the inter-fragment interaction energy (IFIE) between inhibitor fragments and the fragments of protein (amino acid residues), and evaluated the sum of these IFIEs, called IFIE sum, as binding energy between an inhibitor and renin. We also estimated the sum of IFIEs between inhibitor fragments and a fragment of protein as interaction energy between an inhibitor and an amino acid residue in renin. In addition, in order to discuss the details of interaction, we carried out the pair interaction energy decomposition analysis (PIEDA) calculation [22, 23] and evaluated the energy components. Furthermore, we examined the characteristic indications of inhibitor binding with the IFIEs between each amino acid residue in renin and inhibitors by using visualized cluster analysis (VISCANA) method [24] implemented in BioStationViewer. The FMO data were registered in the FMO database (FMODB) [25, 26], and their entry codes (FMODB IDs) are listed in Table 1.

3. Results and Discussion

3.1 Relationship between the binding energy and IC₅₀

We obtained the binding energy between an inhibitor and renin by the FMO calculation mentioned above. Table 1 lists IC₅₀ and pIC₅₀ (= −log₁₀IC₅₀) values for each inhibitor, and the IFIE sum (the calculated binding energy) for each complex with PDB and FMODB IDs [26]. In Figure 3, we show the relationship between the IFIE sum and pIC₅₀. As shown in Figure 3, the IFIE sum can be correlated to pIC₅₀ of the inhibitors (R² = 0.7284; R = −0.8535). This result indicates that it is possible to predict the activity value such as IC₅₀ of the inhibitors by the IFIE sum.

Table 1. IC₅₀ and pIC₅₀ values for each inhibitor, and IFIE sum for each complex with PDB and FMODB IDs [26]

Figure 3. Relationship between IFIE sum and pIC₅₀ for renin and its inhibitor complexes

3.2 IFIE analysis including PIEDA calculation

We also obtained the interaction energy between an inhibitor and each amino acid residue in renin by the procedure mentioned in Section 2.2. For a detailed discussion on the interactions between inhibitors and renin, we performed the IFIE analysis including PIEDA calculation. In PIEDA calculation, the IFIE is divided into four components, electrostatic interaction (ES), exchange repulsion (EX), charge transfer with higher order mixed terms (CT + mix), and dispersion interaction (DI). The interaction energies between the inhibitors and the selected amino acid residues in renin for each complex and their energy components are shown in Figure S1 and S2 in Supplementary information, respectively. Here, as an example, the interaction energies and their energy components of 2V0Z and 2V16 are shown in Figure 4 and 5, respectively. All inhibitors have strong attractive interactions with the active site amino acid residues, Asp38 and Asp226. In addition to Asp38 and Asp226, each inhibitor has relatively strong attractive interactions (|energy| > 5 kcal/mol) with the particular amino acid residues (Table 2). As shown in Table 2, we also found that many inhibitors strongly interact with the identical amino acid residues such as Gln19, Tyr20, Tyr83, Ser84, Thr85, Phe119, Phe124, and Ala229. Furthermore, from the results of PIEDA (Figure 5 and Figure S2), we found that the main component of interaction energies between inhibitors and the active site amino acid residues, Asp38 and Asp226, and that between inhibitors and Ser84/Thr85 are ES. Moreover, it is found that the main component of interaction energy between inhibitors and the aromatic amino acid residues such as Tyr83, Phe119, and Phe124 is DI. It will be considered that the DI such as the van der Waals interaction is the main source of the CH/π interaction.

Figure 4. Interaction energy between inhibitor and the selected amino acid residues in renin (IFIE) of 2V0Z and 2V16

Figure 5. Energy components of interaction energy (PIEDA) of 2V0Z and 2V16

Table 2. Strongly interacting amino acid residues with inhibitor for each complex

The amino acid residues listed in the table show an attractive interaction of more than 5 kcal/mol.

3.3 Clustering of inhibitors including VISCANA method

We also focused on the binding mode of inhibitor and renin to examine the relationship between the structure and IC₅₀ value for each inhibitor. The inhibitors having similar structures tend to show similar binding mode with renin. First, we classified the inhibitors from a viewpoint of the structure of inhibitor and the binding mode. The resulting classification is as follows: (a) 2V0Z and 2V16, (b) 3O9L and 3OAD, (c) 3VYD, 3VYE, and 3VYF, (d) 3Q4B, 3Q5H, and 3GW5, (e) 4GJA, 4GJC, and 4GJD, (f) 2G1R, 2G24, and 2IIKO, (2G24 and 2IKO have the same inhibitor). The inhibitors of 3D91, 3Q3T, 4GJ7, and 2BKT have the individual structures and the binding modes, respectively. Figure 6 illustrates the structures and the binding modes of each complex.

In classification (a) (Figure 6), both of 2V0Z and 2V16 have two hydrogen bonding interactions: (1) between the hydrogen atom of hydroxyl group in the inhibitor and the oxygen atom of side chain carboxyl group of Asp38, and (2) between the hydrogen atom of primary amine in the inhibitor and the oxygen atom of side chain carboxyl group of Asp226.

In classification (b) (Figure 7), for both of 3OAD and 3O9L, a hydrogen bonding interaction is found between the hydrogen atom of secondary amine of piperidine ring in the inhibitor and the oxygen atom of side chain carboxyl group of Asp226.

In classification (c) (Figure 8), for all of three complexes, 3VYD, 3YE, and 3VYF, there are two hydrogen bonding interactions: (1) between the hydrogen atom of secondary amine of piperidine ring in the inhibitor and the oxygen atom of side chain carboxyl group of Asp226, and (2) between the oxygen atom of carbonyl group in the inhibitor and the hydrogen atom of main chain of Thr85. For both of 3VYE and 3VYF, there is an additional hydrogen bonding interaction between the oxygen atom of carbonyl group in the inhibitor and the hydrogen atom of side chain hydroxyl group of Ser84.

In classification (d) (Figure 9), for three complexes, 3Q4B, 3Q5H, and 3GW5, the hydrogen atom of secondary amine in each inhibitor exhibits a hydrogen bonding interaction with the oxygen atom of side chain carboxyl group of Asp38 (3Q4B and 3GW5) or Asp226 (3Q5H). For both of 3Q4B and 3Q5H, respectively there are three additional hydrogen bonding interactions: (1) between the oxygen atom of carbonyl group in the inhibitor and the hydrogen atom of side chain hydroxyl group of Ser84, and (2) between the oxygen atom of tetrahydropyran ring in the inhibitor and the hydrogen atom of side chain hydroxyl group of Thr85, and (3) between the oxygen atom of carbonyl group in the inhibitor and the hydrogen atom of main chain of Tyr20.

In classification (e) (Figure 10), for three complexes, 4GJA, 4GJC, and 4GJD, respectively have two hydrogen bonding interactions: (1) between the hydrogen atom of secondary amine of piperidine ring in each inhibitor and the oxygen atom of side chain carboxyl group of Asp38 (4GJA and 4GJC) or Asp226 (4GJD), (2) between the oxygen atom of carbonyl group in the inhibitor and the hydrogen atom of main chain of Ser84.

In classification (f) (Figure 11), for all of three complexes, 2G1R, 2G24, and 2IKO, there are two hydrogen bonding interactions between the hydrogen atom of primary amine in the inhibitor and the oxygen atoms of side chain carboxyl group of Asp38 and Asp226, respectively. For 2G1R, there is an additional hydrogen bonding interaction between the oxygen atom of carbonyl group in the inhibitor and the hydrogen atom of main chain of Tyr20.

In the other case (Figure 12), for 3D91, two hydrogen bonding interactions are found: (1) between the hydrogen atom of hydroxyl group in the inhibitor and the oxygen atom of side chain carboxyl group of Asp38, and (2) between the oxygen atom of hydroxyl group in the inhibitor and the hydrogen atom of side chain hydroxyl group of Ser84. For 3Q3T, two hydrogen atoms of primary amine in the inhibitor exhibit two hydrogen bonding interactions with the oxygen atom of side chain carboxyl group of Asp38 and Asp226, respectively. For 4GJ7, there are two hydrogen bonding interactions: (1) between the hydrogen atom of secondary amine of pyrrolidine ring in the inhibitor and the oxygen atom of side chain carboxyl group of Asp226, and (2) between the oxygen atom in the inhibitor and the hydrogen atom of main chain of Tyr20. For 2BKT, there is a hydrogen bonding interaction between the hydrogen atom of secondary amine of piperazine ring in the inhibitor and the oxygen atom of side chain carboxyl group of Asp226.

Figure 6. Structures and the binding modes of 2V0Z and 2V16

Light blue broken lines indicate the hydrogen bonding interaction.

Figure 7. Structures and the binding modes of 3O9L and 3OAD

Light blue broken lines indicate the hydrogen bonding interaction.

Figure 8. Structures and the binding modes of 3VYD, 3VYE, and 3VYF

Light blue broken lines indicate the hydrogen bonding interaction.

Figure 9. Structures and the binding modes of 3Q4B, 3Q5H, and 3GW5

Light blue broken lines indicate the hydrogen bonding interaction.

Figure 10. Structures and the binding modes of 4GJA, 4GJC, and 4GJD

Light blue broken lines indicate the hydrogen bonding interaction.

Figure 11. Structures and the binding modes of 2G1R, 2G24, and 2IKO

Light blue broken lines indicate the hydrogen bonding interaction.

Figure 12. Structures and the binding modes of 3D91, 3Q3T, 4GJ7, and 2BKT

Light blue broken lines indicate the hydrogen bonding interaction.

However, it should be noted that the inhibitors with similar structure/binding mode do not necessarily show similar IC₅₀ values. Thereupon, we performed VISCANA calculation and classified the inhibitors not only from a viewpoint of each structure/binding mode but also in terms of the interaction pattern between an inhibitor and amino acid residues in renin. The result of VISCANA for twenty complexes is shown in Figure 13. We found that the IFIEs between the inhibitor and amino acid residues located more than 10 Å away from the inhibitor were negligibly weak. Therefore, Figure 13 shows the IFIEs between the inhibitor and residues located within 10 Å of the inhibitor. In Figures 14 and 15, we show the IFIE sum for each complex (Figure 14) and the selected interaction energies with some amino acid residues (Figure 15(a)-(d)) from Figure 13. The arrangement of the complexes on the abscissa in Figures 14 and 15 is the same as the order on the ordinate in Figure 13. In Figure 14, it is found that the inhibitors of 2BKT, 2G24, and 2IKO weakly bind to renin. In Figure 15(a)-(d), we display the characteristic interaction energies, such as interactions between inhibitor and Asp38/Asp226, that between inhibitor and Ser41, those between inhibitor and Tyr83/Ser84/Thr85, and that between inhibitor and Ala229, respectively. All inhibitors interact very strongly with the active site amino acid residues Asp38 or Asp226 (Figure15 (a)). As seen in Figure 15(b), for many complexes, the interaction between inhibitor and Ser41 is repulsive. However, the interactions for 3VYD, 3VYE, 3VYF, 2V0Z, and 2V16 are strongly attractive, and those for 3O9L and 3OAD are very weakly attractive (Figure 15(b)). As shown in the interaction energies between inhibitor and Tyr83/Ser84/Thr85 (Figure 15(c)), we found that the energy values and pattern for 3Q4B and those for 3Q5H are similar. This similarity was also found in 4GJA and 4GJC, and 3VYE and 3VYF (Figure 15(c)). For 3O9L and 3OAD, it was observed that the energy values are slightly different, but energy patterns are similar (Figure 15(c)). In Figure 15(d), the interaction energies between inhibitor and Ala229 for 2BKT, 2G24, 2IKO, 3VYE, and 3VYF are very weak, but those for 3Q4B and 3Q5H are strongly attractive. It is notable that, for example, for the group of 3VYD, 3VYE, and 3VYF, although the structure of inhibitors and binding modes are similar to each other, the values and pattern of interaction energies for some particular amino acid residues of 3VYD are different from those of 3VYE and 3VYF (Figures 13, Figure 15(c) and 15(d)). Interestingly, from Figure 16, also for the group of 3Q4B, 3Q5H, and 3GW5 and the group of 4GJA, 4GJC, and 4GJD, we found the same feature mentioned above as for the group of 3VYD, 3VYE, and 3VYF. That is, for three inhibitors in the group, the structure of inhibitors and binding modes are similar to each other, and two of the three inhibitors show the similar energy pattern for Tyr83/Ser85/Thr86.

Figure 13. VISCANA for twenty complexes

Abscissa and ordinate represent amino acid residues in renin located within 10 Å of the inhibitor

and values of interaction energies for each complex, respectively. Color code: red, Attractive

interaction; blue, Repulsive interaction. Colored lines on the left side of PDB IDs correspond to colored enclosed lines in Figure 16.

Figure 14. IFIE sum for each complex

Figure 15. Interaction energies between inhibitor and the selected amino acid residues for each complex

(a) Asp38/Asp226, (b) Ser41, (c) Tyr83/Ser84/Thr85, and (d) Ala229.

In Figure 16, we show the resulting cluster of inhibitors considering the structure of inhibitors, the binding mode of inhibitor and renin, and the interaction pattern between inhibitor and the particular amino acid residues in renin. Blue and red circles in Figure 16 indicate the functional group interacting with Asp38 and Asp226, respectively.

Complexes with inhibitors having attractive interaction with Ser41 are enclosed by orange line. Among them, complexes are divided by similarity of structure of the inhibitors, similarity of binding form, and presence of interactions with Tyr83, Ser84, and Thr85 into cluster-1 (3VYD, 3VYE and 3VYF: enclosed by green line) and cluster-2 (2V0Z and 2V16: enclosed by purple line). Cluster-1 and -2 are the same as classification (c) and (a) mentioned in section 3.3, respectively. Moreover, cluster-1 is classified by intensity of interaction with Ala229 and interaction pattern with Tyr83, Ser84, and Thr85 into 3VYD and 3VYE, 3VYF (enclosed by red line).

Cluster of 3O9L and 3OAD (enclosed by yellow line) has very weak attractive interaction with Ser41, similar structure of the inhibitors, similar binding mode, and similar pattern for almost amino acid residues, each other. This cluster is the same as classification (b) in section 3.3.

Cluster of 4GJA, 4GJC, and 4GJD (enclosed by magenta line) is classified by similarity of structure of the inhibitors, similarity of binding form, and presence of interactions with Tyr83, Ser84, Thr85, and Ala229. This cluster is the same as classification (e) in section 3.3. Furthermore, 4GJA and 4GJC (enclosed by light blue line) are classified by having a similar interaction pattern for Tyr83, Ser84, and Thr85, and relatively strong repulsive interaction with Ser41.

Cluster of 3Q4B, 3Q5H, and 3GW5 (enclosed by purple line) possesses relative strong attractive interaction with Ala229, similar structure of the inhibitors, and similar binding mode. This cluster is the same as classification (d) in section 3.3. Farther, 3Q4B and 3Q5H (enclosed by light green line) show a similar interaction pattern for Tyr83, Ser84, and Thr85, each other, and strongly interact with Ala229.

Cluster of 2G1R, 2G24, and 2IKO (enclosed by dark red line) has similar structure of the inhibitors and similar binding mode. This cluster is the same as classification (f) in section 3.3. Strong attractive interaction with Ala229 is found in 2G1R, whereas interaction with Ala229 is very weak in 2G24 and 2IKO. In addition to 2G24 and 2IKO, 2BKT also shows very weak binding of the inhibitor and renin.

Inhibitors, 3D91, 4GJ7, and 3Q3T (enclosed by black line), have no similarity to other inhibitors.

We found that the inhibitors in the finally classified clusters such as (3VYE, 3VYF), (2V0Z, 2V16), (3O9L, 3OAD), (4GJA, 4GJC) and (3Q4B, 3Q5H) tend to show the same order of IC₅₀ value, respectively. This clustering result is also consistent with the result of VISCANA (Figure 13). That is, we can discuss the correlation between the structure and the activity values such as IC₅₀ of each inhibitor by this clustering. In other words, our computational results suggest that differences in interaction due to slight differences in structure, such as the addition/replacement of a single atom/functional group, can be related to the differences in activity value. Furthermore, our results also indicate that VISCANA method can extract the characteristic interactions and identify the important amino acid residues involved in the interaction and the activity value.

In the present work, we performed the FMO calculations with all inhibitors being electrically neutral. However, most of inhibitors interact with the active site amino acid residues, Asp38 and/or Asp226, in renin at their amine, so, it is likely that each inhibitor actually has different charge states.

Although we need to carry out the FMO calculation including the charge states of inhibitor in the future, there are some difficult problems to be solved in the calculation: a quantitative comparison of interaction energies of complexes with ligands having different charges, and a high-precision structural optimization such as QM calculation, etc. The slight structural changes in inhibitors or amino acid residues have a significant effect on the electrostatic interaction between the amino acid residues in renin and the charged inhibitors. Then we performed a trial calculation with the protonated ligands by an Auto-FMO protocol [27], and found that the complexes with different protonation states showed significant differences in the interactions, particularly the electrostatic interaction. Therefore, the computation including the charge states of inhibitor is in the preliminary stage at the moment.

Figure 16. Clustering of inhibitors considering the structure of inhibitors, the binding mode of inhibitor and renin, and the interaction pattern between the inhibitor and the particular amino acid residues in renin

Blue and red circles indicate the functional group interacting with Asp38 and Asp226, respectively.

4. Conclusion

We performed FMO calculations for twenty different renin-inhibitor complexes to explore the relationship between the calculated interaction energies such as IFIEs and IFIE sum, and IC₅₀ (pIC₅₀) values of inhibitors. In the present work, we dealt with all inhibitors as electrically neutral. The binding energy (IFIE sum) between the inhibitor and renin showed a relatively high correlation with pIC₅₀ for each inhibitor. This result suggests that the FMO computations can enable us to evaluate and to predict the activity values of individual inhibitors.

Additionally, we carried out the detailed IFIE analysis including PIEDA calculation, and the clustering of inhibitors including VISCANA method. Consequently, we found that many inhibitors interact with some particular amino acid residues including the active site Asp38 and Asp226, and identified the main component of interaction energies between inhibitor and each amino acid residue. Moreover, we were able to classify the inhibitors not only by their structures/binding mode but also by the pattern of interaction energy and the interacting amino acid residues. The inhibitors in the consequently classified, same cluster were found to be more likely to show the same order of IC₅₀ values. These results indicate that the structure and the activity value of inhibitor are to be related to each other through the interaction between an inhibitor and relevant amino acid residues.

At the moment, the reason of the strong correlation between the binding energy and IC₅₀ value and the role of the particular amino acid residues interacting with inhibitor are not clear. However, as shown in the present and previous works [4-7], the IFIE analysis is a useful method for obtaining detailed information on protein-ligand binding. The interactions for the selected amino acid residues in protein can lead to much useful information for developing more effective and safer agents than current ones. We consider that our FMO computations give us a clue to elucidate the molecular recognition mechanism of renin and its inhibitors. The goal of our ongoing study is to make good use of our computational methods to in silico drug design and to contribute toward the innovative drug discovery.

Now, there are some issues to be solved in the future: (1) identification of the amino acid residues involved in the molecular recognition mechanism of renin and its inhibitors, (2) protonation states of the active site amino acid residues, Asp38 and Asp226, (3) clustering of inhibitors including PIEDA calculation, (4) roles of surrounding water molecules (solvent effects) in protein-ligand binding, (5) charge states of the inhibitors. For problem (2), we intend to carry out FMO simulations considering all possible protonation states of the active site amino acid residues, Asp38 and Asp226. We believe that these FMO calculations may allow the prediction of the protonation state of the active site amino acid residues. For problem (4), we are preparing to perform FMO calculations that take water molecules into account. The strong correlation observed between interaction energies calculated by FMO calculations and activity values, pIC₅₀, suggests that the binding of the inhibitors to renin is enthalpy-driven. The consideration of crystallographic water molecules and the investigation of hydration structures using molecular dynamics (MD) simulations are subjects for future study. FMO calculations that include water molecules are expected to provide insights into interaction networks mediated by water molecules. In this paper, we used the renin-inhibitor complexes directly in our calculations, as PDB structures are available for all of them. However, because the number of available co-crystal structures for typical targets is limited, actual ligand design often relies on a single structure. Therefore, in practical ligand design, we consider that it is important to perform IFIE analysis based on FMO calculations by first conducting docking simulations using a selected structure, followed by structural optimization. Further, we are planning similar FMO calculations toward the complexes of a human protease and its inhibitor. These additional biomolecular simulations may enable us to predict not only the effects of drugs but also the degree of side effects. We are carrying out a collaborative research supported by the Basis for Supporting Innovative Drug Discovery and Life Science Research (BINDS) program as a part of Platform Project for Supporting Drug Discovery and Life Science Research from Japan Agency for Medical Research and Development (AMED).

Supplementary Information

The interaction energy between an inhibitor and each amino acid residue in renin (IFIE) and their energy components (PIEDA) are shown in Figure S1 and S2, respectively.

Acknowledgements

This research was performed as activities of the FMO drug design consortium (FMODD). This work used computational resources of the K computer (project ID: hp190119). MIZUHO/ABINIT-MP software package on PIEDA calculation and MIZUHO/BioStation Viewer 3.0 were provided by Mizuho Information and Research Institute, Inc.. This research was partially supported by Platform Project for Supporting Drug Discovery and Life Science Research (BINDS program) from AMED under Grant Number JP20am0101113.

References

• [1]  Kitaura, K.; Ikeo, E.; Asada, T.; Nakano, T.; Uebayasi, M. Fragment Molecular Orbital Method: an Approximate Computational Method for Large Molecules. Chem. Phys. Lett. 1999, 313, 701-706.doi:10.1016/S0009-2614(99)00874-X
• [2]  Fedorov, D.G.; Kitaura, K. The Fragment Molecular Orbital Method: PRACTICAL APPLICATION TO LARGE MOLECULAR SYSTEM; CRC Press: Boca Raton, 2009.
• [3]  Tanaka, S.; Mochizuki, Y.; Komeiji, Y.; Okiyama, Y.; Fukuzawa, K. Electron-correlated Fragment-molecular-orbital Calculations for Biomolecular and Nano Systems. Phys. Chem. Chem. Phys. 2014, 16, 10310-10344. doi:10.1039/C4CP00316
• [4]  Yagi, Y.; Hattori, Y.; Naoshima, Y. Large-Scale Quantum Chemical Calculation on the Complexes of HIV-1 Protease and Inhibitors: A Relationship of Calculation Energies with Drug Effects of Inhibitors. Trans. Jpn. Soc. Simulat. Tech. 2012, 4, 41-50. doi:10.11308/tjsst.4.41
• [5]  Yagi, Y.; Tanaka, T.; Imagawa, A.; Moriya, Y.; Mori, Y.; et al. Large-Scale Biomolecular Chemical Computations toward the Prediction of Burkholderia cepacia Lipase Enantioselectivity. J. Adv. Simulat. Sci. Eng. 2014, 1, 141-160. doi:10.15748/jasse.1.141
• [6]  Yagi, Y.; Kimura, T.; Kamezawa, M.; Naoshima Y. Biomolecular Chemical Simulations toward Elucidation of the Enantioselectivity and Reactivity of Lipases in Organic Synthesis. CBIJ. 2018, 18, 21-31. doi:10.1273/cbij.18.21
• [7]  Sheng, Y.; Watanabe, H.; Maruyama, K.; Watanabe, C.; Okiyama, Y.; et al. Towards Good Correlation between Fragment Molecular Orbital Interaction Energies and Experimental IC₅₀ for Ligand Binding: A Case Study of p38 MAP Kinase. Comput. Struct. Biotech. J. 2018, 16, 421-434. doi:10.1016/j.csbj.2018.10.003
• [8]  Imai, T.; Miyazaki, H.; Hirose, S.; Hori, H.; Hayashi, T.; et al. Cloning and Sequence Analysis of cDNA for Human Renin Precursor. Proc. Natl. Acad. Sci. USA. 1983, 80, 7405-7409. doi:10.1073/pnas.80.24.7405
• [9]  Mendez, D.; Gaulton, A.; Bento, A.P.; Chambers, J.; Veij, M.D.; et al. ChEMBL: towards direct deposition of bioassay data. Nucleic Acids Res. 2019, 47(D1), D930-D940. doi:10.1093/nar/gky1075
• [10]  wwPDB consortium. Protein Data Bank: the single global archive for 3D macromolecular structure data. Nucleic Acids Res. 2019, 47(D1), D520–D528. doi:/10.1093/nar/gky949
• [11]  Waterhouse, A.; Bertoni, M.; Bienert, S.; Studer, G.; Tauriello, G.; et al. SWISS-MODEL: homology modelling of protein structures and complexes. Nucleic Acids Res. 2018, 46(W1), W296-W303. doi:10.1093/nar/gky427
• [12]  Bienert, S.; Waterhouse, A.; de Beer, T.A.P.; Tauriello, G.; Studer, G.; et al. The SWISS-MODEL Repository - new features and functionality. Nucleic Acids Res. 2017, 45(D1), D313-D319. doi:10.1093/nar/gkw1132
• [13]  Guex, N.; Peitsch, M.C.; Schwede, T. Automated comparative protein structure modeling with SWISS-MODEL and Swiss-PdbViewer: A historical perspective. Electrophoresis 2009, 30 S162-S173. doi:10.1002/elps.200900140
• [14]  Benkert, P.; Biasini, M.; Schwede, T. Toward the estimation of the absolute quality of individual protein structure models. Bioinformatics 2011, 27, 343-350. doi:10.1093/bioinformatics/btq662
• [15]  Bertoni, M.; Kiefer, F.; Biasini, M.; Bordoli, L.; Schwede, T. Modeling protein quaternary structure of homo- and hetero-oligomers beyond binary interactions by homology. Sci. Rep. 2017, 7, 10480. doi:10.1038/s41598-017-09654-8
• [16]  Case, D.A.; Darden, T.A.; Cheatham, III, T.E.; Simmerling, C.L.; Wang, J.; et al. AMBER 11; University of California, San Francisco, 2010.
• [17]  Jakalian, A.; Jack, D.B.; Bayly, C.I. Fast, efficient generation of high‐quality atomic charges. AM1‐BCC model: II. Parameterization and validation. J. Comput. Chem. 2002, 23, 1623-1641. doi:10.1002/jcc.10128
• [18]  Wang, J.; Cieplak, P.; Kollman, P.A. How well does a restrained electrostatic potential (RESP) model perform in calculating conformational energies of organic and biological molecules?. J. Comput. Chem. 2000, 21, 1049-1074. doi:10.1002/1096-987X(200009)21:12<1049::AID-JCC3>3.0.CO;2-F
• [19]  Wang, J.; Wolf, R.M.; Caldwell, J.W.; Kollman, P.A.; Case, D.A. Development and testing of a general amber force field. J. Comput. Chem. 2004, 25, 1157-1174. doi:10.1002/jcc.20035
• [20]  Okiyama, Y.; Nakano, T.; Yamashita, K.; Mochizuki, Y.; Taguchi, N.; et al. Acceleration of fragment molecular orbital calculations with Cholesky decomposition approach. Chem. Phys. Lett. 2010, 490, 84-89. doi:10.1016/j.cplett.2010.03.001
• [21]  Mizuho Research & Technologies, Ltd. MIZUHO/ABINIT-MP and MIZUHO/BioStation Viewer 3.0; 2013. https://www.mizuho-rt.co.jp/solution/research/life/macromolecule/biostation/index.html (accessed June 05, 2022, in Japanese).
• [22]  Fedorov, D.G.; Kitaura, K. Pair interaction energy decomposition analysis. J. Comp. Chem. 2007, 28, 222-237. doi:10.1002/jcc.20496
• [23]  Tsukamoto, T.; Kato, K.; Kato, A.; Nakano, T.; Mochizuki, Y.; et al. Implementation of Pair Interaction Energy Decomposition Analysis and its Applications to Protein-Ligand Systems. J. Comput. Chem. Jpn. 2015, 14, 1-9. doi:10.2477/jccj.2014-0039
• [24]  Amari, S.; Aizawa, M.; Zhang, J.; Fukuzawa, K.; Mochizuki, Y.; et al. VISCANA: Visualized Cluster Analysis of Protein−Ligand Interaction Based on the ab Initio Fragment Molecular Orbital Method for Virtual Ligand Screening. J. Chem. Inf. Model. 2006, 46(1), 221-230. doi:10.1021/ci050262q
• [25]  Takaya, D.; Watanabe, C.; Nagase, S.; Kamisaka, K.; Okiyama, Y.; et al. FMODB: The World’s First Database of Quantum Mechanical Calculations for Biomacromolecules Based on the Fragment Molecular Orbital Method. J. Chem. Inf. Model. 2021, 61(2), 777-794. doi: 10.1021/acs.jcim.0c01062
• [26]  FMO database, https://drugdesign.riken.jp/FMODB/ (accessed Oct 10, 2025)
• [27]  Watanabe, C.; Watanabe, H.; Okiyama, Y.; Takaya, D.; Fukuzawa, K.; et al. Development of an automated fragment molecular orbital (FMO) calculation protocol toward construction of quantum mechanical calculation database for large biomolecules. CBIJ. 2019, 19, 5-18. doi:10.1273/cbij.19.5