Fragment Molecular Orbital Based Interaction Analyses on COVID-19 Main Protease - Inhibitor N3 Complex (PDB ID: 6LU7).

The worldwide spread of COVID-19 (new coronavirus found in 2019) is an emergent issue to be tackled. In fact, a great amount of works in various fields have been made in a rather short period. Here, we report a fragment molecular orbital (FMO) based interaction analysis on a complex between the SARS-CoV-2 main protease (Mpro) and its peptide-like inhibitor N3 (PDB ID: 6LU7). The target inhibitor molecule was segmented into five fragments in order to capture site specific interactions with amino acid residues of the protease. The interaction energies were decomposed into several contributions, and then the characteristics of hydrogen bonding and dispersion stabilization were made clear. Furthermore, the hydration effect was incorporated by the Poisson-Boltzmann (PB) scheme. From the present FMO study, His41, His163, His164, and Glu166 were found to be the most important amino acid residues of Mpro in interacting with the inhibitor, mainly due to hydrogen bonding. A guideline for optimizations of the inhibitor molecule was suggested as well based on the FMO analysis.

Introduction

The 2019 novel coronavirus (SARS-CoV-2 or its disease name COVID-19) caused the pneumonia outbreak in Wuhan, China, in late December 2019 and has rapidly spread around the world.[1] As of March 19, 2020, around 190 000 people were infected and over 7700 have succumbed to the epidemic (as of May 25, 5 200 000 and 337 000, respectively). The World Health Organization (WHO) has declared this novel coronavirus outbreak a global health emergency, and subsequently has described the outbreak as pandemic on March 11. At present, there are no targeted therapeutics, and effective treatment options remain very limited. Although a large number of researchers around the world are engaged in developing antiviral drugs against COVID-19, it is also well-known that new drug discovery and development is a long, costly, and rigorous scientific process.

Many researchers investigating the genetic and functional data of SARS-CoV-2 compare it with other coronaviruses (CoVs) to design proper infection control strategy and seek potential drugs that can prevent and/or cure this serious epidemic. CoVs infect humans and other animal species, causing a variety of highly prevalent and severe diseases, including Severe Acute Respiratory Syndrome (SARS) and Middle East Respiratory Syndrome (MERS). CoVs are positive-sense, single-stranded RNA viruses featuring the largest viral RNA genomes known to date. The SARS-CoV-2 genome is comprised of approximately 30 000 nucleotides. As for the development of relevant drugs, viral proteases are well-known common targets in dealing with human viruses such as the HIV virus and hepatitis C virus. Protease inhibitors would thus be remarkably effective in blocking the replication of coronavirus, including the SARS and the MERS, providing a promising foundation for the development of anticoronaviral therapeutics. Based on its functions, the main protease (Mpro) or chymotrypsin-like protease (3CLpro) is suggested to be a potential drug target to combat COVID-19, which is highly conserved among coronaviruses. Sequence alignment revealed that the Mpro of SARS-CoV-2 shares 96% similarity with that of previous SARS-CoV.[2] Studies for identifying the inhibitors of COVID-19 Mpro have then quickly begun for discovering and developing drugs against the disease. For example, on the basis of homology models for SARS-CoV-2 Mpro and docking simulations of 1903 small drug compounds, Xu et al.[3] selected four promising drug candidates according to the docking score, 3D similarity of binding mode to known ligands and the binding free energy evaluation. They thus suggested that nelfinavir can be a potential inhibitor against the SARS-CoV-2 Mpro.

The X-ray crystal structure of SARS-CoV-2 Mpro was released on February 5, 2020 (and revised a week later) by Liu et al.[4] It is of dimer-like structure in complex with a peptide-like inhibitor called N3,[5] and very similar to that of previous SARS-CoV[6,7] reflecting the 96% sequence homology between them. We can thus anticipate the development of dependable structure-based design of inhibitors against SARS-CoV-2 Mpro. In order to obtain reliable strategies for this aim, it is now desirable to know accurate information on the manner of ligand-protein interactions in the pharmacophore.

As illustrated in a recent book of in silico drug design,[8] the fragment molecular orbital (FMO) method[9−12] provides an efficient tool for performing ab initio quantum-chemical calculations for biomolecular systems and accurately analyzing their intermolecular interactions in terms of the interfragment interaction energies (IFIEs). Namely, the set of IFIEs can be used as numerical indices to capture the nature of interactions in a given target system. In fact, the FMO method has been successfully applied to the comprehensive analyses of protein–ligand[13−16] and protein–protein[17−19] interactions, including viral molecular systems. The aim of the present work is to apply the FMO method to the analysis of the interactions between SARS-CoV-2 Mpro and its inhibitor N3 to elucidate the essence of its pharmacophore structure and associated molecular recognition.

In the following sections, after addressing the present methodology by means of the FMO-IFIE, we show the calculated results on the interactions between the amino-acid residues of Mpro and the fragments of N3. With the aid of the advanced analyses based on the Pair Interaction Energy Decomposition Analysis (PIEDA)[20] and the incorporation of the solvent effect[21,22] as well, we will provide some useful information toward a rational drug design for the inhibitors against SARS-CoV-2 Mpro.

Method of Calculation

In this study, we adopted the crystal structure of the complex[4] between COVID-19 Mpro and inhibitor N3 and processed it for the FMO based interaction analyses. The setup of calculations is described as follows.

Processing of the Crystal Structure

First, the crystal structure of the MPro:N3 complex (PDB ID: 6LU7; Version 2, 2.16 Å resolution, 4) was downloaded. As stated in the Introduction, although Mpro usually functions as a homodimer, we calculated only the monomeric unit stored in the PDB file since the dimer interface has no direct interaction with N3 ligand. We further processed the molecular structure with the MOE program.[23] All of the 84 crystal water molecules were preserved and subjected to subsequent FMO calculation. Figure illustrates the peptide-like N3 inhibitor as the ligand to bind to the protein. Note that this ligand has a covalent bond to Cys145 as the result of the Michael addition. We performed standard modeling operations such as hydrogen attachment and pKa adaptation (pH 7). The charged N- and C-terminuses were assigned to Ser1 and Gln306, respectively. Molecular mechanics based energy minimization (with AMBER10:EHT force field[24] in MOE[23]) was carried out for all the attached hydrogen atoms, the ligand, and amino acid residues within 4.5 Å of the ligand, where a tether mask of 1.0 Å was used for the crucial pharmacophore region. The entire structure of the processed protein–ligand complex and the close-up view of pharmacophore are shown in Figure . Figure is a schematic illustration of the ligand–residue interactions, based on geometrical relations judged by MOE.[23] This illustration may be informative as a guide in discussing the FMO results.

Figure 1

Structure of inhibitor N3 (ligand).

Figure 2

Structure of complex between COVID-19 main protease and inhibitor. Left: overall view; right: selected residues in pharmacophore.

Figure 3

Schematic situation of interactions in pharmacophore based on geometrical relations drawn by MOE.[23]

FMO Calculations

The FMO method[9] was originated by Kitaura in 1999, at the two-body expansion (FMO2) for energy aswhere I and J are fragment indices and the second term in the right-hand side is defined as the increment from the dimer[10−12]Namely, the FMO2 energy is given as the sum of monomer and dimer energies including effective many-body contributions. For a target system, the two-body FMO calculation consists of the monomer stage (in which a mutual consistency of electrostatic potential (ESP) is required) and the dimer stage, and these computations are highly parallelized for efficient processing. By extracting the ESP part, eq can be rewritten asThe first term in this equation is the fictitiously isolated monomer energy, whereas the second term corresponds to pair interaction energy (PIE)[11] or interfragment interaction energy (IFIE),[12,25] which is quite useful for interaction analyses[10,12] (see refs (13−19) if necessary). We therefore used IFIE for the interaction analysis in this paper.

In the FMO calculation of a protein,[10−12] the protein is usually fragmented at the Cα carbon atom with a sp3 hybridization[26] as exemplified in Figure , where there is no hydrogen-capping for the bond segmentation via the technique of bond detachment atom (BDA). In this study, we regarded an amino acid residue as a fragment to maintain a one-to-one correspondence between a residue and a fragment. It should be noted that care is necessary when the >C=O part in a peptide bond is involved in interactions because the formal assignment of reside is shifted as in Figure (27) (it will be addressed later). In the present case, as the N3 inhibitor is too large to be regarded as a fragment, we divided it into 5 fragments as shown in Figure . This fragmentation pattern was similar to that often used to divide peptides. To be specific, each fragment in Figure was designed so as to mimic an amino acid residue fragment of a peptide. We chose this pattern because the overall molecular structure of N3 is analogous to that of a peptide. The fragmentation of both Mpro and N3 was done with the FMOe program,[28] and the total number of fragments became 395 (=306 (protein) + 5 (ligand) + 84 (water)).

Figure 4

Example of fragmentation for mini protein consisting of five residues (total number of fragments is thus five). The bonding electron pair at segmentation point is expressed as a colon symbol. Care is necessary for fac that the segmentation is not made at the peptide bond (refer to texts).

Figure 5

Fragmentation of inhibitor moiety. The bonding electron pair at segmentation point is expressed as a colon symbol with yellow mask.

Our FMO program ABINIT-MP[12] was used for the present study. The electron correlation effect was incorporated by the second-order Møller–Plesset perturbation (MP2) calculation[29] which is efficiently executable with ABINIT-MP.[30−32] A partial renormalization[33] was utilized to reduce the trend of overstabilizations by MP2. The 6-31G* basis set[29] was adopted as a standard choice in FMO calculations.[9,12]

The decomposition analysis of PIE (called as PIEDA[20]) can grasp the nature of interactions in detail. We thus adopted the PIEDA option in ABINIT-MP.[34] In PIEDA, an interaction energy, PIE or IFIE, between fragments I and J is expressed aswhere ES, EX, CT, and DI mean “ElectroStatic”, “EXchange repulsion”, “Charge-Transfer and mix terms”, and DIspersion”, respectively. The DI energy is usually evaluated by the MP2-correlated calculation. Based on our experiences, contributions from ES and CT in eq take large values for the hydrogen bonding, while the so-called CH/π[35] and π/π interactions have large DI portions. As the hydrogen-bonding is potentially influenced by the hydration effect, the Poisson–Boltzmann (PB) model[21,22] was employed to incorporate the hydration effect. The PB grand iteration was iterated until the reaction field of hydration converged to 1.0 × 10–5 au in energy. Once converged, PIEDA was performed to obtain the respective contributions including the additional “Solvation” term. To identity the type of interaction, we introduced the following index with quantities of eq (inspired by ref (36))where the “(Solv)” term in the denominator is valid only for the case of PB solvation. The quantity of eq is multiplied by 100% (not necessarily normalized) for discussion in the next section. Note that a large value of FES+CT suggests the existence of hydrogen bonding.

Additional Model without N3-Cys145 Bond

We additionally considered a noncovalent bonding model in which the N3 ligand is not connected to Cys145. This model was generated from the 6LU7 structure through the following steps with the MOE program;[23] (I) The covalent bond between N3 and Cys145 was cleaved to restore the α,β-unsaturated carbonyl reactant of the Michael addition reaction, (II) constrained optimization (tether = 1) was performed for the region around N3 (4.5 Å) by fixing the main chain of residues, and (III) selected optimization for N3, Cys145, and two connected residues. Note that the point of fragmentation in N3 was shifted due to the presence of double bond by which the Michael addition takes place. Finally, the FMO calculation was carried out without PB.

Results and Discussion

Since the ligand moiety was divided into five fragments as shown in Figure , we describe the interaction analyses in the order of fragment index. We will directly cite the evaluated values in the tables for detailed discussions. The corresponding bar graph representations of PIEDA are given in Figures S1–S5 (in the Supporting Information) for graphical understanding.

Fragment 1

Table compiles the results of IFIE and decomposed contributions (PIEDA) interacting with Fragment 1 of the ligand, where the listing threshold is set as 2.0 kcal/mol The distance between the main chain >C=O of Thr190 and the N–H part of Fragment 1 is as close as 1.99 Å as illustrated in Figure , suggesting that they are forming a typical hydrogen bond. As addressed already, there is a concern of assignment shift[27] that this carbonyl part belongs to theAla191 fragment in the segmentation of FMO scheme.[10−12] The IFIE of Ala191 fragment is −11.58 kcal/mol. The ES term of the IFIE is −9.15 kcal/mol, and the CT term is −4.52 kcal/mol, indicating an electrostatic interaction with charge transfer. The FES+CT index is 69.5% certainly. These features indicate the existence of hydrogen bond.

Table 1

IFIE and PIEDA Results for Fragment 1

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	DI (kcal/mol)	FES+CT (%)
Leu167	2.53	–4.45	–1.80	1.90	–2.23	–2.32	63.5
Gln189	2.67	–3.00	–1.23	0.36	–0.82	–1.32	60.8
Ala191	1.99	–11.58	–9.15	8.10	–4.52	–6.01	69.5
Gln192	2.50	–2.06	–1.32	0.94	0.40	–2.08	30.7

Figure 6

Close-up view of interacting residues for Fragment 1.

The IFIE and PIEDA values under the PB condition are given in Table . For the hydrogen bonding on Ala191 (actually Thr190 residue), the hydration effect increments the stabilization. The sums of IFIE values (of all residues in the protease) for Fragment 1 with and without hydration are −26.20 and −24.21 kcal/mol, respectively. The interaction between Fragment 1 and protease can be characterized as hydrogen bonding type, consistent with the situation in Figure .[23]

Table 2

IFIE and PIEDA Results under Hydration for Fragment 1

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	Solv (kcal/mol)	DI (kcal/mol)	FES+CT (%)
Leu167	2.53	–4.46	–1.00	1.90	–2.18	–0.90	–2.28	50.0
Gln189	2.67	–3.37	–1.41	0.36	–0.78	–0.21	–1.33	58.7
Ala191	1.99	–12.49	–9.52	8.02	–4.52	–0.46	–6.01	68.5
Gln192	2.50	–1.85	–1.47	0.92	0.39	0.44	–2.11	39.3

Fragment 2

Fragment 2 has a Val-like side chain. Table lists the IFIE and PIEDA results, and Figure shows a couple of important geometrical features. There is a hydrogen bond between the main chain >C=O part of Glu166 and the N–H part of Fragment 2, where the distance is found to be 1.91 Å. The assignment shift due to the fragmentation[27] appears again, and thus the IFIE of Leu167 fragment should be of interest. The corresponding value is −14.50 kcal/mol, and its FES+CT index is 87.3%. Another notable stabilization by hydrogen bonding (−10.69 kcal/mol) is associated with the side chain of Gln189, where the associated distance is 2.44 Å as indicated in Figure . Note that the importance of Glu166 and Gln189 is consistently suggested in Figure .

Table 3

IFIE and PIEDA Results for Fragment 2

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	DI (kcal/mol)	FES+CT (%)
Glu166	2.75	–3.92	–2.25	0.16	–0.50	–1.34	67.2
Leu167	1.91	–14.50	–16.91	7.03	–1.87	–2.74	87.3
Pro168	3.22	2.13	2.11	0.05	0.34	–0.38	118.4
Gln189	2.44	–10.69	–6.24	5.45	–4.65	–5.25	67.5
Thr190	3.26	–3.25	–4.60	–0.01	1.79	–0.42	87.0
Ala191	3.81	–3.57	–3.37	0.01	0.00	–0.20	94.4

Figure 7

Close-up view of interacting residues for Fragment 2.

The PB results of IFIE and PIEDA are compiled in Table . It is notable that the hydration enhances the stabilization for Leu167 fragment (actually Glu166 residue) but provides destabilization for Gln189 although their amounts are small. The sums of IFIE values for Fragment 2 with and without hydration are −26.92 and −28.41 kcal/mol, respectively. In other words, the PB hydration leads to a certain electrostatic shielding for Fragment 2.

Table 4

IFIE and PIEDA Results under Hydration for Fragment 2

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	Solv (kcal/mol)	DI (kcal/mol)	FES+CT (%)
Glu166	2.75	–3.15	–2.10	0.16	–0.60	0.69	–1.31	81.3
Leu167	1.91	–14.75	–17.06	7.08	–1.82	–0.21	–2.73	86.5
Pro168	3.22	1.68	2.05	0.05	0.35	–0.39	–0.37	146.3
Gln189	2.44	–8.51	–5.84	5.51	–4.63	1.63	–5.19	74.6
Thr190	3.26	–3.80	–4.49	–0.01	1.76	–0.65	–0.41	72.0
Ala191	3.81	–4.26	–3.47	0.01	–0.02	–0.58	–0.20	81.7

Fragment 3

As shown in both Figures and 8, the main chain >N–H part of Glu166 and the >C=O part of Fragment 3 have a hydrogen bonding, where the distance of interest is 1.88 Å and the corresponding IFIE value is −9.18 kcal/mol found in Table . There is another hydrogen bonding with the side chain of Gln189 (distance of 1.90 Å), whose IFIE value is −10.96 kcal/mol. These two interactions have sizable FES+CT values, as expected. Interestingly, several residues such as His41 and Met165 provide vital stabilizations with small FES+CT values, suggesting that dispersion type interactions play an additional role in binding of Fragment 3. Both are hydrophobic interactions with the Leu-like side chain of Fragment 3; they are CH/π interaction with the imidazole ring of His41 (see Figure S6) and hydrophobic interaction with the thioether chain of Met165.

Figure 8

Close-up view of interacting residues for Fragment 3.

Table 5

IFIE and PIEDA Results for Fragment 3

residue name	distance (Å)	IFIE (kcal/mol)	ESa (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	DI (kcal/mol)	FES+CT (%)
His41	2.45	–3.51	–0.80	2.91	–1.36	–4.25	33.7
Met49	2.47	–2.28	–1.07	1.78	–0.49	–2.49	38.5
His163	5.77	–7.20	–7.20	0.00	0.00	0.00	100.0
Met165	2.15	–5.53	–2.70	7.72	–4.47	–6.07	54.2
Glu166	1.88	–9.18	–13.81	8.12	–0.53	–2.96	82.9
Leu167	3.69	–3.11	–2.54	0.03	–0.19	–0.41	86.9
Arg188	2.95	–2.29	–0.95	0.84	–0.59	–1.58	49.4
Gln189	1.90	–10.96	–13.31	10.73	–3.09	–5.29	75.6

The dimer-ES approximation[37] was adopted due to a long distance from this fragment.

Table lists the PB results. The hydrogen bonding stabilizations of Glu166 and Gln189 are increased by the hydration. The sums of IFIE values for Fragment 3 with and without hydration are −49.45 and −47.50 kcal/mol, respectively. These values are about twice larger than the corresponding values of Fragment 1 and Fragment 2.

Table 6

IFIE and PIEDA Results under Hydration for Fragment 3

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	Solv (kcal/mol)	DI (kcal/mol)	FES+CT (%)
His41	2.45	–3.65	–0.94	2.90	–1.35	–0.02	–4.24	35.0
Met49	2.47	–2.23	–1.07	1.79	–0.56	0.10	–2.49	40.5
His163	5.77	–7.05	–7.17	0.00	0.00	0.12	0.00	101.7
Met165	2.15	–5.15	–2.34	7.73	–4.47	0.00	–6.07	52.9
Glu166	1.88	–11.03	–14.68	8.08	–0.43	–1.00	–3.00	79.1
Leu167	3.69	–3.00	–2.41	0.03	–0.18	–0.04	–0.40	85.5
Arg188	2.95	–2.29	–0.92	0.83	–0.60	0.00	–1.60	48.7
Gln189	1.90	–11.80	–14.05	10.89	–3.04	–0.21	–5.39	75.3

Fragment 4

Figure suggests that His164 has a hydrogen bonding with the ligand, and this bonding is actually formed between its main chain >C=O part and the >N–H part of Fragment 4 as illustrated in Figure (distance of 1.87 Å). Due to a concern of assignment shift,[27] the corresponding stabilization is found as the IFIE value of −15.78 kcal/mol for Met165 fragment (see Table ). Notably, the neighbored residue His163 has a large stabilization of −37.25 kcal/mol as IFIE, where the protonated side chain ring of His 163 plays a vital role in interactions with the >C=O part of the ring in Fragment 4 (close distance of 1.74 Å). The FES+CT index of His163 is as large as 89.2%, just indicating the hydrogen bonding. Glu166 provides a stabilization of −8.73 kcal/mol, and a relatively small FES+CT value (57.0%) implies a contribution from dispersion. Phe140 and Asn142 have dispersion type stabilizations as well.

Figure 9

Close-up view of interacting residues for Fragment 4.

Table 7

IFIE and PIEDA Results for Fragment 4

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	DI (kcal/mol)	FES+CT (%)
His41	3.46	–2.36	–0.92	0.03	–0.57	–0.89	62.6
Phe140	2.73	–2.27	–0.43	1.10	–1.01	–1.93	42.7
Leu141	2.69	–7.35	–5.44	0.88	–0.90	–1.89	77.0
Asn142	2.37	–3.66	1.28	2.24	–2.39	–4.79	18.8
His163	1.74	–37.25	–41.11	14.93	–5.46	–5.61	89.2
Met165	1.87	–15.78	–18.53	10.13	–3.65	–3.73	85.6
Glu166	2.22	–8.73	–4.48	4.93	–3.30	–5.88	57.0
Leu167	4.19	–3.85	–3.58	0.01	–0.11	–0.16	95.8

The PB results are listed in Table . The hydration leads to a stabilization for Met165 (actually His164 fragment). In contrast, Glu166 is rather destabilized via the electrostatic shielding (positive value of “Solv” term). Because the charged carboxylic acid in the side chain of Glu166 is exposed to the solvent, opposite-signed and sizable charges should be induced on its solvent surface. The interaction with these surface charges effectively suppresses the “ES” component between Fragment 4 and Glu166 as the “Solv” one. These can also be interpreted from more microscopic changes in electron density distribution using natural population analysis (NPA)[38,39] (refer to Figure S7). One carboxyl O atom of Glu166 and H of Fragment 4 that form a hydrogen bonding are depolarized with the magnitudes of 0.008e and 0.018e in atomic charge by the hydration, respectively. This weakens the “ES” interaction slightly than in vacuo. Another carboxyl O atom exposed to the solvent is rather polarized with the magnitude of 0.110e while inducing the solvent charges reactively. The presence of such a solvent-charge induction increases the significance of the “Solv” term because it incorporates the electrostatic shielding. The sums of IFIE values for Fragment 4 with and without hydration are −79.58 and −84.85 kcal/mol, respectively. These values are outstanding in comparison with other four Fragments. Namely, this Fragment 4 should be a crucial moiety of the ligand (inhibitor N3).

Table 8

IFIE and PIEDA Results under Hydration for Fragment 4

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	Solv (kcal/mol)	DI (kcal/mol)	FES+CT (%)
His41	3.46	–1.84	–0.78	0.03	–0.57	0.36	–0.89	71.8
Phe140	2.73	–2.22	–0.27	1.11	–0.95	–0.18	–1.93	36.6
Leu141	2.69	–6.14	–5.33	0.91	–0.89	1.04	–1.87	88.2
Asn142	2.37	–4.89	0.10	2.22	–2.72	0.33	–4.82	36.8
His163	1.74	–35.33	–40.09	14.96	–5.54	0.79	–5.45	90.7
Met165	1.87	–17.09	–19.10	10.04	–3.51	–0.78	–3.75	83.3
Glu166	2.22	–1.81	–3.25	4.94	–3.13	5.16	–5.51	94.8
Leu167	4.19	–4.17	–3.67	0.01	–0.10	–0.24	–0.17	90.2

Fragment 5

According to Figure , His41 can interact with Fragment 5. This suggestion is confirmed by the calculated results in Table . A hydrogen bonding forms (distance of 2.16 Å in Figure ) between the >N–H part of side chain ring of His41 and the >C=O part of Fragment 5. The IFIE value is −10.57 kcal/mol, and consistently the FES+CT value is 83.4%. CH/π interaction between phenyl ring of Fragment 5 and Thr25 is also observed (Figure S6), whose main energy component is DI (−3.17 kcal/mol).

Table 9

IFIE and PIEDA Results for Fragment 5

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	DI (kcal/mol)	FES+CT (%)
Thr25	2.60	–4.58	–1.81	2.11	–1.72	–3.17	52.7
Leu27	2.40	–2.79	–1.48	2.24	–1.29	–2.27	55.0
His41	2.16	–10.57	–8.97	2.29	–1.77	–2.13	83.4
Asn142	3.32	–4.61	–3.20	0.06	–0.51	–0.97	79.3
Gly143	2.38	–2.34	–2.40	1.10	1.44	–2.49	27.8
His163	5.21	2.80	2.84	0.00	0.00	–0.03	101.1

Figure 10

Close-up view of interacting residues for Fragment 5.

The PB hydration increases the stabilization for His41 as shown in Table . The sums of IFIE values for Fragment 5 with and without hydration are −17.77 and −21.32 kcal/mol, respectively, indicating an electrostatic shielding. The role of Fragment 5 is observed to be rather minor relative to other four fragments.

Table 10

IFIE and PIEDA Results under Hydration for Fragment 5

residue name	distance (Å)	IFIE (kcal/mol)	ES (kcal/mol)	EX (kcal/mol)	CT (kcal/mol)	Solv (kcal/mol)	DI (kcal/mol)	FES+CT (%)
Thr25	2.60	–4.04	–1.60	2.11	–1.87	0.49	–3.16	56.5
Leu27	2.40	–2.76	–0.79	2.27	–1.40	–0.61	–2.23	43.5
His41	2.16	–11.12	–9.23	2.30	–1.76	–0.27	–2.16	81.9
Asn142	3.32	–4.41	–3.28	0.07	–0.38	0.16	–0.96	82.1
Gly143	2.38	–2.17	–2.24	1.08	1.42	–0.01	–2.41	25.3
His163	5.21	1.39	2.75	0.00	–0.01	–1.32	–0.04	198.6

Effect of Elimination of Covalent Bond between Cys145 and N3

Figure compares the MOE-based optimized structures around the N3 ligand with and without the covalent bond to Cys145. As illustrated in this figure, most of the structural changes before and after the Michael addition reaction took place on the ligand, while the protein was almost unchanged except for Cys145, where the S–C bond distance of Cys145 was 1.83 Å with the bond but 3.22 Å without it. The IFIE and PIEDA results for Cys145 are compiled in Table S1 for nonbonded Cys145, where the fragmentation points are illustrated in Figure S8. Interestingly, the interaction energies between the nonbonded Cys145 and Fragments 4 and 5 are not large before the Michael reaction. The remaining Fragments 1–3 show almost no interaction with the protein.

Figure 11

Pharmacophore structures around the N3 ligand with and without covalent bond to Cys145. Green and light blue frames correspond to the ligands with and without covalent bond, respectively.

For Fragments 4 and 5, the interaction energies are listed for important residues in Tables S2 and S3, respectively. It is notable that the conformational changes of the ligand lead to destabilization of the hydrogen bonding between Fragment 4 and the main chain > C=O of His164 (ΔE of +5.86 kcal/mol, relative to the corresponding entry of Met165 in Table ) and also stabilization (ΔE of −11.63 kcal/mol, similarly) of that with Gln166. About Fragment 5, there was a sizable conformational change for the terminal Phe ring, which destabilized the interaction with Thr25 (+2.00 kcal/mol, relative to the value in Table ) and that with Leu27 (+3.39 kcal/mol) mainly due to the loss of DI. For Fragments 1–3, no table difference was found with and without the covalent bond.

Summary of Interaction Analyses

In this section, a series of interaction analyses have been made in fragment-by-fragment fashion. Here, we would make several notes as follows. First, all five fragments of the N3 ligand have the hydrogen bonding with neighbored residues in the pharmacophore. In particular, Fragment 4 interacts strongly with His163, His164, and Glu166, and this fragment plays a leading role in stabilizing the ligand. Besides these three residues, His41, Gln189, and Thr190 have hydrogen bonds with fragments of the ligand. Second, the dispersion-type interactions assist the ligand binding (e.g., Phe140, Asn142, and Glu166 interacting with Fragment 4). Third, the phenyl ring of Fragment 5 and isoxazole ring of Fragment 1 do not have prominent interaction with surrounding residues. This fact suggests that if some chemical modifications are made on this ring to further optimize the inhibitor, interactions with Thr25, Thr26, and Asn142 for Fragment 5, and Ala191 and Pro168 for Fragment 1 may be induced. Fourth, the role of the peptide-like side chains in the ligand is illuminated. For example, the Leu-like side chain of Fragment 3 has a CH/π interaction with His41. His41 has a CH /π bond to Fragment 3 and also a hydrogen bond to fragment 5, and thus plays an important role in ligand binding. Supplemental modifications on the groups such as Ala-like side chain (currently not interacting with the neighboring residues) in Fragment 1 may also be effective. These findings obtained by the FMO calculation were in general agreement with the consideration given by X-ray crystallography.[4] Our discussion with quantitative interaction energy values could provide additional quantitative insights. In addition to confirming and reinforcing the experimental observations, the FMO calculation has newly found that His41 plays an important role in N3 binding.

Conclusion

In the present work, a series of FMO-based interaction analyses[10,12] with ABINIT-MP[12] have been made on the complex between the COVID-19 main protease and the peptide-like inhibitor N3, whose fundamental structure was obtained from PDB ID: 6LU7.[4] The object ligand was divided into five fragments, and their interactions were investigated by using the IFIE[12,25] and PIEDA[34] with and without the PB hydration.[21,22] It was found that the hydrogen bonds with surrounding residues in the pharmacophore were important throughout five fragments. The fourth fragment of the ligand was identified as the most important moiety in interactions with His163, His164 and Glu166. It was also noted that dispersion interactions provided additional stabilization of the ligand, such as the CH/π interaction with His41. A possibility of further optimization of ligand binding was suggested at the fifth fragment of the ligand for the interactions with Thr25, Thr26 and Asn142 in addition to His41. The influence of the covalent bond between the ligand and Cys145 was also investigated. Though some differences were observed for several residues, essential discussion on interaction energies were kept unchanged. It would be noted technically that these biochemical findings based on the IFIE and PIEDA results were carefully derived by avoiding the issue of assignment shift due to Cα fragmentation;[27] the function of peptide bond fragmentation[40] has been recently available in ABINIT-MP, but a cross-reference testing should be required before routine usages.

Virtual screening and high-throughput screening based on the crystal structure of the complex of Mpro and N3 have recently been performed.[4] In ref (4), the authors have discussed why cinanserin is a potential inhibitor for SARS-CoV-2. The present FMO analysis supports this idea because the significance of attractive interactions between the inhibitor and surrounding residues such as His41, Met49, Asn142, Met165 and Glu166 were identified quantitatively. Compounds that mimic many of these interactions are then expected to have higher inhibitory activity, indicating a direction of HIT-to-lead process. As a whole, our quantum-chemical analysis brings additional quantitative information to the recent comprehensive study[4] and can provide useful insights into the rational drug discovery based on the N3 pharmacophore. Thus, we hope that the present work by the FMO calculations for the SARS-CoV-2 Mpro inhibitor could contribute to the worldwide efforts toward the development of effective drugs against the COVID-19.[41]

From a stand-point of open science promotion against COVID-19, we have uploaded the data set of IFIE and PIEDA on the present system to the FMO database (FMODB).[42,43]