Computational analysis of dynamic allostery and control in the SARS-CoV-2 main protease

Abstract

The COVID-19 pandemic caused by the novel coronavirus SARS-CoV-2 has no publicly available vaccine or antiviral drugs at the time of writing. An attractive coronavirus drug target is the main protease (M^pro, also known as 3CL^pro) because of its vital role in the viral cycle. A significant body of work has been focused on finding inhibitors which bind and block the active site of the main protease, but little has been done to address potential non-competitive inhibition, targeting regions other than the active site, partly because the fundamental biophysics of such allosteric control is still poorly understood. In this work, we construct an elastic network model (ENM) of the SARS-CoV-2 M^pro homodimer protein and analyse its dynamics and thermodynamics. We found a rich and heterogeneous dynamical structure, including allosterically correlated motions between the homodimeric protease's active sites. Exhaustive 1-point and 2-point mutation scans of the ENM and their effect on fluctuation free energies confirm previously experimentally identified bioactive residues, but also suggest several new candidate regions that are distant from the active site, yet control the protease function. Our results suggest new dynamically driven control regions as possible candidates for non-competitive inhibiting binding sites in the protease, which may assist the development of current fragment-based binding screens. The results also provide new insights into the active biophysical research field of protein fluctuation allostery and its underpinning dynamical structure.

Keywords: SARS-CoV-2; allostery; elastic network model; protein dynamics.

Figure 1.

The crystal structure of the SARS-CoV-2 M^pro with N3 inhibitor. Chain A of the dimer is shown in blue, while chain B, in orange. Domains are labelled with boxed Roman numerals (I, II and III). Amino acid residues of the catalytic dyad are indicated as yellow spheres for H41 and magenta spheres for C145. Asterisks mark residues from chain B (orange). Chain termini are shown as spheres and labelled N and C for chain A (blue) and N* and C* for chain B (orange). The N3 inhibitor is shown as green sticks. Experimentally identified dynamically allosteric residues N214 and SAL284-286 are shown as labelled black spheres. Identified control residue candidates for dynamic allostery in the SARS-CoV-2 M^pro that are distant from the catalytically active residues H41 and C145 are coloured in black on both homodimeric chains. These residues can be found in table 1.

Figure 2.

Constructing ENM of SARS-CoV-2 M^pro step-by-step. (a) SARS-CoV-2 M^pro secondary structure cartoon. (b) Elastic model of M^pro generated with PyANM package in PyMOL. C_α atoms are shown in blue; while node-connecting springs (black) are shown only for one chain for comparison. (c) The first real vibrational mode eigenvectors (yellow) visualization. For clarity, displacement vectors are scaled five times.

Figure 3.

The cross correlation of the motion of 6lu7 ENM. (a) The cross-correlation maps calculated for the first real 25 modes (bottom-right region of the plots) and spacing between residues (C_α nodes) and ligand nodes (top-left region of the plots) for apo, holo1 and holo2 forms. The first colour scale shows the extent of cross correlation, with a cross correlation of 1 (red) indicating perfectly correlated motion, −1 (blue) showing perfectly anti-correlated motion and 0 (white) no correlation. The second colour scale (black to white) depicts the Euclidean distance between two ENM nodes in the Cartesian space in 0–16 Å range. The secondary structure of M^pro is indicated along the horizontal residue axes, with cyan waves indicating alpha helices, and magenta triangles indicating beta sheets. Coloured lines on the vertical axes point to the boxed regions on the plots mentioned in the main text. The green ticks on the axis indicate the location of the biologically active residues (table 1). (b) A real-space representation of the correlations in 6lu7 ENM with respect to residue 145 on chain A (green sphere). The dashed squares indicate the anti-correlation enhanced in apo–holo1 and reduced in holo1–holo2 transitions. The cross correlation matrix was calculated using only the C_α atoms for the protein and all heavy atoms for the ligand (N3 inhibitor).

Figure 4.

Mutation scan maps for thermodynamic control of M^pro calculated from the ENM over the first real 25 fluctuation modes. (a) A map for the fluctuation free energy change. The map plots the relative change in free energy to the wild-type ((G_mut − G_wt)/|G_wt|) due to the dimensionless change in the spring constant (k_R/k) for the mutated residue with the amino acid number shown. White corresponds to values of free energy predicted by the wild-type ENM. Red corresponds to an increase in (G_mut − G_wt)/|G_wt| (decreased value of G_mut comparing to G_wt), whereas blue corresponds to a decrease in (G_mut − G_wt)/|G_wt| (increased value of G_mut comparing to G_wt). (b) The map for the vibrational free energy change plotted in real space onto the wild-type M^pro homodimer structure at k_R/k = 4.00. (c) A map for the global control space of allostery in M^pro. The map plots the change in cooperativity coefficient (K₂/K₁) due to the dimensionless change in the spring constant (k_R/k) for the mutated residue with the amino acid number shown. White corresponds to values of K₂/K₁ predicted by the wild-type ENM. Red corresponds to an increase in K₂/K₁ (stronger negative cooperativity), whereas blue corresponds to a decrease in K₂/K₁ (weaker negative cooperativity or positive cooperativity). (d) The global map plotted in real space onto the wild-type M^pro homodimer structure at k_R/k = 0.25.

Figure 5.

Two-point mutational maps for 6lu7 ENM with all possible pairwise combinations of residue mutations with equal spring constant change k_R/k equal to 0.25 and 4.00 over the first real 25 fluctuation modes. (a,b) 2-point mutational maps for 6lu7 ENM with all possible pairwise combinations of residue mutations with equal spring constant change (a) k_R/k = 4.00 (spring stiffening) and (b) k_R/k = 0.25 (spring relaxation). (c,d) Maps for the two-dimensional global control space of allostery in M^pro for (c) k_R/k=4.00 and (d) k_R/k = 0.25 with colour bar centre at K₂/K₁ = 1.032. Black solid lines separate two homodimer chains, while dashed lines represent 1-point mutational scan results for the given spring constant change.