key: cord-0696738-i9ar8to1
authors: Sk, Md Fulbabu; Jonniya, Nisha Amarnath; Roy, Rajarshi; Poddar, Sayan; Kar, Parimal
title: Computational Investigation of Structural Dynamics of SARS-CoV-2 Methyltransferase-Stimulatory Factor Heterodimer nsp16/nsp10 Bound to the Cofactor SAM
date: 2020-11-24
journal: Front Mol Biosci
DOI: 10.3389/fmolb.2020.590165
sha: d674e8a0650dd20d95df5909ca8bb9fa3c2da7e2
doc_id: 696738
cord_uid: i9ar8to1

Recently, a highly contagious novel coronavirus disease 2019 (COVID-19), caused by SARS-CoV-2, has emerged, posing a global threat to public health. Identifying a potential target and developing vaccines or antiviral drugs is an urgent demand in the absence of approved therapeutic agents. The 5′-capping mechanism of eukaryotic mRNA and some viruses such as coronaviruses (CoVs) are essential for maintaining the RNA stability and protein translation in the virus. SARS-CoV-2 encodes S-adenosyl-L-methionine (SAM) dependent methyltransferase (MTase) enzyme characterized by nsp16 (2′-O-MTase) for generating the capped structure. The present study highlights the binding mechanism of nsp16 and nsp10 to identify the role of nsp10 in MTase activity. Furthermore, we investigated the conformational dynamics and energetics behind the binding of SAM to nsp16 and nsp16/nsp10 heterodimer by employing molecular dynamics simulations in conjunction with the Molecular Mechanics Poisson-Boltzmann Surface Area (MM/PBSA) method. We observed from our simulations that the presence of nsp10 increases the favorable van der Waals and electrostatic interactions between SAM and nsp16. Thus, nsp10 acts as a stimulator for the strong binding of SAM to nsp16. The hydrophobic interactions were predominately identified for the nsp16-nsp10 interactions. Also, the stable hydrogen bonds between Ala83 (nsp16) and Tyr96 (nsp10), and between Gln87 (nsp16) and Leu45 (nsp10) play a vital role in the dimerization of nsp16 and nsp10. Besides, Computational Alanine Scanning (CAS) mutagenesis was performed, which revealed hotspot mutants, namely I40A, V104A, and R86A for the dimer association. Hence, the dimer interface of nsp16/nsp10 could also be a potential target in retarding the 2′-O-MTase activity in SARS-CoV-2. Overall, our study provides a comprehensive understanding of the dynamic and thermodynamic process of binding nsp16 and nsp10 that will contribute to the novel design of peptide inhibitors based on nsp16.

Coronaviruses (CoVs) are considered as an etiological agent for causing severe acute respiratory syndrome (SARS) in humans. In the past two decades, SARS-CoV and MERS-CoV (middle east respiratory syndrome coronavirus) were responsible for the epidemic in 2003 and 2012, respectively. Recently, in December 2019, the novel coronavirus pandemic, caused by SARS-CoV-2, first broke out in Wuhan city of China and has been spreading worldwide (Bogoch et al., 2020; Guan et al., 2020; Lin X. et al., 2020) . So far, ∼42 million people worldwide are infected by SARS-CoV-2, including ∼1.1 million deaths. In India, the COVID-19 tally has already crossed 7.7 million, including ∼0.1 million deaths. However, neither prophylactic vaccines nor any direct antiviral drugs are available to effectively treat the human and animal coronavirus disease (Lu, 2020; Pillaiyar et al., 2020; Sheahan et al., 2020) .

CoVs are single-stranded positive-sense RNA viruses (Eckerle et al., 2010; Fehr and Perlman, 2015) belonging to the family Coronaviridae and possess the largest genome (26.4-31.7 kb) (Woo et al., 2010) . They are classified into four genera, namely Alphacoronaviruses, Betacoronaviruses, Gammacoronaviruses, and Deltacoronaviruses (King et al., 2011) . CoVs can infect both humans and animals (Lun and Qu, 2004; Coleman and Frieman, 2014) , and can cause diseases like Hepatitis and Pneumonitis in mouse (Weiss and Leibowitz, 2011) and neurologic & respiratory diseases in humans (Arabi et al., 2015; Huang et al., 2020) . So far, seven distinctive strains of human coronaviruses (HCoV) have been disclosed that includes 229E and NL63 (Alphacoronaviruses), and OC43, HKU1, SARS-CoV (2002) , MERS-CoV (2012) , and SARS-CoV-2 (2019) (Betacoronaviruses) (Weiss and Navas-Martin, 2005; Zeng et al., 2018; Singh et al., 2020) . SARS-CoV-2 shares a strong correlation with the bat CoV RaTG13 having a 96.2% genome sequence identity, suggesting its possible evolution from bat CoVs Zhou P. et al., 2020) . SARS-CoV-2 displays a sequence similarity of 80% with SARS-CoV, whereas only 50% with MERS-CoV Wu et al., 2020; Zhou Y. et al., 2020) . Currently, most of the therapeutic options that are available for controlling COVID-19 is based on the previous knowledge and information gathered from SARS-CoV and MERS-CoV.

Most viruses or eukaryotic cellular mRNA possess the 5 ′end capping mechanism that plays a vital role in mRNA splicing, translation initiation, stability, and intracellular RNA transport (Furuichi and Shatkin, 2000) . The capping of the 5 ′ -end occurs through a sequential enzymatic process. This involves three enzymes, such as RNA guanylyltransferase (GTase), RNA triphosphatase (TPase), and RNA guanine-N7methyltransferase (N7-MTase). This generates a cap-0 structure (m7GpppN). It is further methylated at the 2 ′ -O position of mRNA by 2 ′ -O-methyltransferase (2 ′ -O-MTase) and generates the cap-1 (m7GpppNm) and cap-2 (m7GpppNmNm) structures. This mimicking of the eukaryotic mRNA capping mechanism ultimately helps the virus evade the host innate immune system (Furuichi and Shatkin, 2000; Wang et al., 2015) . Both the MTase uses S-adenosyl-L-methionine (SAM or AdoMet) as a donor FIGURE 1 | Crystal structure of SARS-CoV-2 nsp10-nsp16 complexed with cofactor S-Adenosyl Methionine (SAM) (PDB: 6W4H). The nsp10, interacting part of nsp10, interacting part of nsp16, nsp16, and the cofactor SAM in a complex is shown in color teal, gold, dark sky blue, light green, and yellow, respectively. The surface represents an electrostatic surface.

of methyl and gives a by-product, S-adenosyl-L-homocysteine (SAH or AdoHcy).

The genome of SARS-CoV-2 is comprised of 10 ORFs (Open Reading Frame). ORF1ab encodes the replicase polyprotein 1ab (PP1ab), which gets cleaved by the two viral proteases, PL pro (papain-like protease) and 3CL pro (3-C like protease) at the Nterminus and C-terminus, respectively, to form all the 16 nonstructural proteins (nsp1, nsp2, nsp3 by PL pro , and nsp4-nsp16 by 3CL pro ). The remaining ORFs are associated with encoding structural proteins, such as spike (S), envelope (E), membrane (M), and nucleocapsid (N) proteins, as well as other accessory proteins (Harcourt et al., 2004; Yang et al., 2005; Chen et al., 2020; Wu et al., 2020) . Previous structural and biochemical characterization studies on SARS-CoV (2002) showed that in the nsp10-nsp16/SAM complex, nsp10 enhances the methylase activity of nsp16 by increasing the stability of the SAM-binding pocket (Chen et al., 2011; Wang et al., 2015) .

The present study involves analyzing protein-protein interactions in the heterodimeric structure formed by nsp10 and nsp16 of SARS-CoV-2. The recently solved crystal structure of SARS-CoV-2 nsp16/nsp10 heterodimer bound to SAM (PDB: 6W4H) is used in the current study. The three-dimensional structure of the complex is shown in Figure 1 . Here, we have used the standard molecular dynamics (MD) simulations in conjunction with the molecular-mechanics Poisson-Boltzmann surface area (MM-PBSA) method to identify the critical residues involved in the complex formation. Further, we conducted MD simulations of nsp16/SAM, elucidating the importance of nsp10 in the binding of SAM to nsp16. The detailed structural analysis and inter-molecular interactions reveal the interface of nsp16/ns10, which provides a distinctive feature for coronaviruses, can be exploited as an attractive target in developing specific antiviral drugs for controlling COVID-19. MD simulations were conducted to elucidate the role of nsp10 in the interaction of SAM to nsp16. Also, the conformational changes in nsp10 and nsp16 due to their association were investigated.

We retrieved the experimental coordinate of the SARS-CoV-2 nsp16/nsp10 complex from Protein Data Bank (PDB), which was crystallized at 1.8 Å (PDB ID 6W4H) . The complex structure includes S-adenosyl-L-methionine (SAM). The protonation states were determined using PROPKA 3.1 (Olsson et al., 2011) , and the corresponding residues were modified accordingly. The following three systems were simulated for the current study: complex (nsp16/nsp10/SAM), nsp16 SAM (nsp16/SAM), and nsp10 apo (apo nsp10). The starting configurations of nsp16 SAM and nsp10 apo were constructed manually from nsp16/nsp10/SAM (PDB: 6W4H).

Missing hydrogens in the crystal structure were built using the Leap module of AmberTools19 (Case et al., 2018) , and a proper amount of Na + was added to neutralize the system. All the systems were solvated in a periodic octahedron TIP3P water box (Price and Brooks, 2004 ) with a 10 Å buffering distance from all directions. The proteins (nsp10 and nsp16) were assigned the Amber ff14SB force field parameters (Maier et al., 2015) while the cofactor SAM was modeled using the updated Generalized Amber Force Field (GAFF2) (Wang et al., 2004) and AM1-bcc (Jakalian et al., 2002) charges. The bond lengths having hydrogen atoms were kept fixed using the SHAKE algorithm (Kräutler et al., 2001) , which allowed 2 fs time-step for MD simulations. The Particle-Mesh Ewald (PME) summation (Darden et al., 1993) scheme was employed to compute the long-range interactions. The non-bonded cut-off was set to 10 Å.

Firstly, the solvated systems were subjected to an energy minimization using 500 steps of steepest descent, followed by another 500 steps of the conjugant gradient algorithm. During the minimization, the solute atoms were kept fixed with a restraint force of 2.0 kcal mol −1 Å −2 . The second stage of minimization was carried out without any restraint force. After the minimization, each system was gradually heated from 0 to 300 K in the NVT ensemble, and the solute atoms were restrained with force constant of 2.0 kcal mol −1 Å −2 . The temperature was maintained using a Langevin thermostat (Pastor et al., 1988; Loncharich et al., 1992) . Subsequently, each system was simulated for 50 ps at 300 K at a constant pressure of 1 atm using the Berendsen barostat (Berendsen et al., 1984) . In this stage also, the same restraint force was applied to the solute atoms. Before the production run, we equilibrated each system by conducting 1 ns MD simulation in the NPT ensemble without restraining the system. Finally, the production simulation was carried out for 1 µs using the pmemd.cuda module of AMBER18. Overall, 100,000 snapshots were generated and used for the analysis using the Cpptraj module (Roe and Cheatham, 2013) of AMBER18 (Case et al., 2018) .

To study the effect of salt concentration on the binding, we simulated the complex system at 300 K for 1 µs at two different salt concentrations (0.15 and 0.25 M). The desired salt concentration was achieved by adding an appropriate number of monovalent ions (Na + and Cl − ) obtained with the protocol suggested by Machado and Pantano (2020) . Similarly, we have also investigated the effect of temperature on the binding by conducting MD simulation of the complex at a relatively elevated temperature of 310 K.

All analyses, including root-mean-square deviations (RMSDs), root-mean-square fluctuations (RMSFs), radius of gyration (R g ), solvent accessible surface area (SASA), etc. were performed using the Cpptraj module (Roe and Cheatham, 2013) of AMBER18. A distance of ≤ 3.5 Å and an angle cut-off of ≥120 • were used for hydrogen bond calculations. The same criterion was used in earlier studies (Jeffrey, 1997; Hu et al., 2016; Shi et al., 2018; Sanachai et al., 2020) .

One of the widely used unsupervised data reduction schemes is principal component analysis (PCA) (Ichiye and Karplus, 1991) . First, a covariance matrix (C) is constructed from the atomic fluctuations of C α -atoms of each residue (Amadei et al., 1996) . The elements C ij of the covariance matrix C are defines as

where x i and x j denotes the instant coordinates of the ith or jth atoms, while < x i > and < x j > denote the mean of ith or jth atoms over the ensemble. The diagonalization of the covariance matrix yields orthogonal eigenvectors and the corresponding eigenvalues. The resulting eigenvectors are termed as principal components (PCs) and indicate the direction of the movement, while the corresponding eigenvalue describes the amplitude of motions. We adopted the same methodology for PCA, as discussed in our earlier studies (Jonniya et al., 2019; Sk et al., 2020c) . The free energy landscape (FEL) was generated based on the following Equation (Frauenfelder et al., 1991) :

where k B represents the Boltzmann constant, and T is the absolute temperature. N i is the population of the ith bin, and N m is the population of the most populated bin. The 2dimensional FEL was constructed using PC1 and PC2 as the reaction coordinate.

The residual network analysis approach is widely used to explore the viral fitness and resistance development of protein structure. The network analysis of protein structures (NAPS) (Chakrabarty et al., 2019) server (http://bioinf.iiit.ac.in/NAPS/) was used to identify key residue interactions in the residual network and the network-based hydrophobic contacts from the simulation trajectories. Here, we used the protein-protein complex option, followed by the C α network type and unweighted edge weight. An edge represented the distance between a pair of C α atoms within the lower and upper thresholds (default upper threshold = 7 Å; lower threshold = 0 Å). We considered the default parameters for the long-range interaction networks and minimum residue separation of 1. The hydrophobicity indices of 20 amino acids range from −4.5 to 4.5 (Kyte and Doolittle, 1982) . These values were colored in gradients of red to show the hydrophobic residues, while the gradients of blue were used to display the hydrophilic residues.

The interaction energy between nsp10 and nsp16, and between the cofactor SAM and nsp16 in its both monomer and the complex were computed using the molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) methodology (Kollman et al., 2000; Wang et al., 2006; Kar et al., 2007a Kar et al., ,b, 2011 Kar et al., , 2013 Hou et al., 2011) . MMPBSA.py script available in the AmberTools19 was used for the analysis. Details of the MM-PBSA protocol were provided in our previous studies (Kar and Knecht, 2012a,b,c,d; Kar et al., 2013; Jonniya et al., 2019; Roy et al., 2020; Sk et al., 2020a,d) , and the same protocol was adopted here. The binding free energy is estimated by using the following equations of the MM-PBSA scheme:

where E internal , G solv , and T S represent the total internal energy, desolvation free energy, and conformational entropy, respectively. Further, the internal energy is composed of E covalent (bond, dihedral, and angle), E elec (electrostatic) and E vdW (van der Waals) and the desolvation free energy is composed of polar ( G pol ) and non-polar solvation energies ( G np ). Here, G pol was estimated from the Poisson-Boltzmann (PB) equation. The dielectric constant of the solute and solvent was set to 1.0 and 80.0, respectively. G np was estimated from the following Equation (6),

Here, SASA represented the solvent-accessible surface area and was estimated using the LCPO (linear combination of pairwise overlap) algorithm (Weiser et al., 1999) with a probe radius of 1.4 Å. The surface tension coefficient, γ and offset (b) value were set to 0.00542 kcal.mol −1 .Å −2 and 0.92 kcal.mol −1 , respectively (Gohlke et al., 2003) . For the MM-PBSA calculation, 2000 snapshots selected from the last 300 ns trajectory with a frequency of 15 ps were employed. The configurational entropy was calculated using the normal mode analysis (NMA) method (Karplus and Kushick, 1981; Rempe and Jónsson, 1998; Xu et al., 2011) , and ∼300 snapshots were used in the calculation due to high computational cost. Each configuration was energy minimized in an implicit solvent (nmode_igb = 1) using a maximum of 50,000 steps, and a target root-mean-square gradient of 10 −4 kcal mol −1 Å −1 via mmpbsa_py_nabnmode (Hawkins et al., 1995 (Hawkins et al., , 1996 Kar and Knecht, 2012a,b,c; Kar et al., 2013) . NMA was applied to calculate the vibrational entropy (S vib ) using the following equation (Carlsson and Åqvist, 2006) ;

Where R denotes the universal gas constant, k is the Boltzmann constant, and T is the absolute temperature. The Planck constant and the vibrational frequency are denoted as h and ν, respectively. Further, the decomposition of the total binding free energy at the residue level was conducted by using the Molecular Mechanics Generalized Born Surface Area (MM-GBSA) scheme (Kar and Knecht, 2012a,b,c; Kar et al., 2013; as per the following equation.

The total contribution from each residue can also be defined as the sum of van der Waals ( E vdW ), electrostatic ( E ele ), polar ( G pol ), and non-polar ( G np ) terms. This method was proposed by Gohlke et al. (2003) . Finally, the Computational Alanine Scanning (CAS) was performed for some essential residues. This method yields the energy difference between the wild type and mutant (alanine) variants.

The basic principle involves in AS (alanine scanning) is the substitution of a residue with alanine that has an impact on the side-chain beyond C β and not in the main-chain. In vitro, AS has been proven an advantageous mutagenesis method in finding hotspot residues in protein-protein interfaces. CAS is an excellent alternative approach to the in vitro experimental alanine scanning (Massova and Kollman, 1999; Moreira et al., 2007) . Therefore, in this study, CAS was applied to illustrate further the importance of specific residues except alanine mentioned in the binding decomposition free energy of nsp16 of COVID-19.

To explore the mechanism underlying the dimerization of nsp16 and nsp10 as well as the preferential binding of the cofactor SAM to the nsp16/nsp10 heterodimer compared to nsp16 (nsp16 SAM ), a conformational free energy landscape (FEL) and binding free energy calculations were performed using the MD/MM-PBSA scheme.

To explore the thermodynamic stability of each system and to ensure the rationality of the sampling method, we monitored the structural and energetic properties during the entire 1 µs production simulation. The time evolution of the root-meansquare deviations (RMSDs) of the protein backbone atoms, which reflects the stability of the system, were calculated relative to the initial configuration and shown in Figure 2A . It is worth mentioning here that a stable RMSD does not always provide stable energy profiles; hence we have also verified the potential and total energy of each system from the respective MD trajectory (data not shown). The RMSD plots indicated that all the studied systems reached equilibrium after ∼100 ns. The average RMSD value varies between 2.8 and 4.0 Å for all systems (see Table 1 ). The highest and lowest deviations were obtained for nsp10 apo and complex, respectively. An intermediate RMSD value of 3.1 Å was obtained for nsp16 SAM . Overall, the nsp16/nsp10/SAM complex showed more stability throughout the simulations than nsp16 SAM or nsp10 apo . This suggests that the dimerization leads to the overall stabilization of the complex.

Moreover, we also calculated the temporal RMSDs of heavy atoms of the cofactor SAM (see Supplementary Figure 1A ) and the backbone atoms of residues within 5 Å around SAM in the binding pocket (see Supplementary Figure 1B ). It is evident from Supplementary Figure 1B that the RMSD values of SAM and the binding pocket for the complex were stable up to 850 ns of the simulation. After that, in the last 150 ns, a sudden increase in RMSD was observed. Overall, the average RMSD value was 0.7 Å and 0.8 Å for SAM and the binding pocket, respectively (see Supplementary Table 1 ). It indicates that SAM binds strongly onto the cofactor binding cavity of the heterodimer nsp16/nsp10. However, in the case of nsp16 SAM , the RMSD values of SAM and the binding cavity were relatively stable up to 100 ns. After that, the RMSD value increased by ∼3 times compared to the complex system (see Supplementary Table 1 and Supplementary Figure 1) . These high values of RMSDs suggest that SAM does not bind to nsp16 monomer (nsp16 SAM ), which is in agreement with the experimental study by Chen and coworkers for SARS-CoV (Chen et al., 2011) . Additionally, three loops, namely 71-79, 100-108, and 130-148, which comprises the SAM binding pocket as shown in the crystal structure, were also analyzed and found to be stable for both complex and nsp16 SAM (see Supplementary Figures 2A-C) . However, the cap-binding groove of nsp16 in SARS-CoV, as well as SARS-CoV-2, is mainly composed of two flexible loops, comprised of residues 26-38 and 130-148, while in the case of flavivirus, the NS5 MTase was relatively stable with α-helices (A1, A2, and half of αD) along the cap-binding groove (Egloff et al., 2002) . Among the two loops involved in the cap-binding, only the loop 26-38 exhibited differences in the complex and nsp16 SAM system, as revealed from our simulations. The 26-38 loop is located near the interface of nsp16 and nsp10. On the other hand, the 130-148 loop is found near the SAM binding pocket, away from the nsp16/nsp10 interface. As shown in Supplementary Table 1 , the average RMSD value of this loop was higher in nsp16 SAM (3.1 Å) compared to the complex (1.4 Å). The time evolution of RMSD of the 26-38 loop and the corresponding potential mean force (PMF) were displayed in Supplementary Figures 3A,B . From Supplementary Figure 3A , it was observed that in the case of nsp16 SAM , the RMSD value of the loop 26-38 increased during the initial 300 ns and remained stable thereafter. On the other hand, in the complex case, the RMSD of the loop was found to increase during the initial 300 ns and gradually decreased up to 600 ns of the simulation and finally reached equilibrium. We computed the potential of mean force (PMF), taking the loop's RMSD as a reaction coordinate and shown in Supplementary Figure 3B . The primary low energy structure of the loop was observed at ∼3.2 Å for nsp16 SAM . However, for the complex system, the global energy minimum was found at a relatively lower RMSD (1.4 Å). Besides, we detected two secondary minima (∼1.9 Å, and ∼2.8 Å) for the 26-38 loop in the complex system, and the energy barriers of the adjacent states were ∼1.1 and 1.2 kcal/mol, respectively. Overall, our results suggest that the dimerization resulted in the loop rearrangement near the interface.

Next, we measured the root-mean-square fluctuations (RMSFs) of Cα atoms for all systems (see Figure 2B ) to elucidate the effect of dimerization on the residual flexibility. It is evident from Figure 2B that all three systems displayed a similar trend in RMSF patterns. However, we found some dynamic fluctuations in different loop regions, including the N-and C-terminals. A relatively large fluctuation was observed for various loop regions located at residues 26-38, 74-80, 104-115, and 130-148 . The loop regions 26-38 and 104-115, which contribute to the interface of the heterodimer, showed higher fluctuations in complex compared to nsp16SAM, which suggests that the dimerization leads to the loop rearrangements. In contrast, the 74-80 loop region, located near SAM binding pocket, stabilized after the heterodimerization. However, the cap-binding groove of the 130-148 loop region exhibited no differences in both complex and apo forms, which agrees with the RMSD analysis (Supplementary Figure 2C) . Overall, it depicts the binding of nsp10 to nsp16 favors and stabilizes the binding pocket of SAM, leading to its high activity.

Finally, we also measured the solvent-accessible surface area (SASA) of all systems and reported in Table 1 . The average SASA values were found to be 18928.1, 13467.3, and 6663.8 Å 2 , for the complex, nsp16 SAM, and nsp10 apo , respectively. This suggested that the dimerization of nsp16 and nsp10 covered ∼1,203 Å 2 surface area in total, indicating a very stable interaction. Our simulation results also favor an experimental study that showed the value of the solvent-exposed surface area for the heterodimer complex of nsp16/nsp10 as 19710 Å 2 (Rosas-Lemus et al., 2020).

The RMSD profile of SAM and the binding pocket (Supplementary Figure 1) for the complex and nsp16 SAM systems suggested that SAM binds more strongly to the heterodimer nsp16/nsp10 than the monomer nsp16. The 2-dimensional free energy landscape (FEL) was generated for the complex and nsp16 SAM systems to explore the strong and weakly bound states of the SAM molecule. The RMSD of heavy atoms of SAM and the center of mass (CoM) distance between SAM and nsp16, which characterizes the displacement of SAM from the substrate-binding pocket of nsp16, were used as reaction coordinates for the construction of FEL. The value of a CoM distance < 14.5 Å denotes the strong binding of SAM toward nsp16. The two-dimensional FEL of the complex (nsp16/nsp10/SAM) and monomer (nsp16 SAM ) systems were shown in Figures 3A,B , and the corresponding positions of SAM in the binding cavity were depicted in Figures 3C,D. FEL of the complex and monomer systems suggested that the conformational state of SAM in two systems prefer disparate configurations. In the case of the heterodimeric complex system, the underlying FEL was characterized by a single global minimum (see Figure 3A) , which corresponds to the strongly bound state of SAM. On the other hand, in the case of nsp16 SAM , the global free energy minimum was obtained at a CoM distance of ∼17.5 Å, which corresponds to a weakly or unbound state of SAM. The secondary minimum was detected at a CoM distance of ∼13.5 Å, which corresponds to a strongly bound state of SAM to nsp16 SAM . Overall, our results suggested that although nsp16 monomer (nsp16 SAM ) favored two binding states, the weakly bound or unbound state is energetically more favorable than the strongly bound state. These results illustrated a stronger binding of SAM with the nsp16/nsp10 heterodimer than nsp16 alone (nsp16 SAM ), which agrees with the experimental results obtained for SARS-CoV (Chen et al., 2011) .

The principal component analysis (PCA) was carried out for both nsp16 and nsp10 when they are in complex and monomeric states. FELs of sub-units nsp16 (nsp16 complex ) and nsp10 (nsp10 complex ) in the complex were compared with the corresponding monomer simulations (nsp16 SAM and nsp10 apo ). Each eigenvector and eigenvalue was plotted in the decreasing order, as shown in Supplementary Figure 4 . For all cases, the first few eigenvectors describe the collective motion of local fluctuations. Comparing the four systems indicated that the first few PCs that describe the properties of movements were not the same. The first two eigenvectors encapsulated for 42, 67, 46, and 70% of overall movements in nsp16 complex (nsp16 in complex), nsp10 complex (nsp10 in complex), nsp16 SAM , and nsp10 apo , respectively. Similarly, the first ten eigenvectors accounted for 80-90% of the total motions in all four systems.

Next, we constructed the 2-dimensional FEL at 300 K for each system using the first two principal components (PC1 and PC2) as reaction coordinates and shown in Figures 4A-D. The sampling space of four systems was different, as evident from Figure 4 . The conformational sampling space of nsp16 complex (nsp16 in nsp16/nsp10/SAM), depicted in Figure 4A , was found to be restricted compared to the other three systems, which sampled more expansive conformational space. From Figure 4A , we observed a global minimum of 62.7% occupancy and a secondary minimum of 37.3% occupancy, which suggested that the nsp16 complex system is more stabilized in the presence of nsp10. In the monomeric form of nsp16 SAM, a wider conformational basin was sampled (see Figure 4C) , and the energy barriers between adjacent minima were ∼1.5 kcal/mol. On the other hand, nsp10 sampled more phase space than nsp16 in the complex as well in its monomer apo form, as shown in Figures 4B,D. However, nsp10 apo has three distinct conformations with an energy barrier of 4.0 kcal/mol suggesting structural disparity in both nsp10 complex and nsp10 apo as compared to nsp16.

The Network Analysis of Protein Structure (NAPS) server provided a visual examination of sub-network based on the physicochemical properties of the protein residues to find more details about the interaction between nsp16 and nsp10 (see Figure 5) . Herein, we explored the 3D interaction network of hydrophobic residues of nsp16 and nsp10. The results showed that the number of hydrophobic interaction networks between nsp16 and nsp10 was initially very high, and after 100 ns, the number of networks reduced to ∼2-3. The hydrophobic interaction networks, such as (V104, A71) and (P80, V42), were found as strong and stable contacts throughout the simulation time (see Supplementary Table 2) . We also calculated the total number of interactions and noticed that the number of interaction networks was more or less conserved throughout the first 500 ns (see Supplementary Figure 5 ). In general, it can be concluded from our analyses that nsp16 and nsp10 interact with a very high affinity.

The differences in binding affinity of the SAM molecule toward nsp16 alone (nsp16 SAM ) and nsp16/nsp10 heterodimer (complex) were calculated by utilizing the MM-PBSA scheme. The MM-PBSA scheme provides various components contributing to the total binding energy ( G bind ), such as van der Waals interactions ( E vdW ), electrostatic interactions ( E ele ), polar solvation energy ( G pol ), non-polar solvation free energy ( G np ), and configurational entropy (T S). All these components of the binding free energy of SAM to nsp16 SAM or nsp16/nsp10 (complex) were shown in Table 2 . The binding free energy of SAM to the heterodimer nsp16/nsp10 was estimated as −6.8 (± 0.9) kcal/mol, which agrees with the experimental result (−7.4 ± 0.1 kcal/mol) (Lin S. et al., 2020) . On the other hand, SAM was found to bind very weakly ( G bind = −0.6 ± 0.8 kcal/mol) to nsp16 monomer (nsp16 SAM ), as evident from Table 2 . This is consistent with Figures 3A,B . A similar mode of SAM binding was observed for SARS-CoV nsp16/nsp10 or nsp16 alone. For both cases, the intermolecular van der Waals interactions ( E vdW ), electrostatic interactions ( E ele ), and non-polar solvation free energy ( G np ) favored the binding of SAM. In contrast, polar solvation energy ( G pol ) and entropy (T S) disfavored the complexation. The electrostatic interaction energy was more favorable than the van der Waals interactions for both cases. In nsp16 SAM , although E ele was much more favorable (−76.8 kcal/mol) than E vdW (−28.6 kcal/mol). However, the disfavorable polar solvation energy ( G pol = 84.7 kcal/mol) overcompensated the intermolecular electrostatic interaction energy. Hence, the overall polar contribution ( E ele + G pol ) was found to be unfavorable (7.9 kcal/mol) to the SAM binding. In contrast, the overall non-polar contribution ( E vdW + G np ) was found to be −32.1 kcal/mol. It implies that hydrophobic interactions drive the binding. In the case of the complex system (nsp16/nsp10/SAM), the binding affinity of SAM to nsp16 increases, as evident from Table 2 . Here, the electrostatic interactions energy ( E ele = −181.1 kcal/mol) was more favorable than the van der Waals energy ( E vdW = −40.4 kcal/mol). The total polar contributions ( E ele + G pol ) value was found to be favorable (−2.1 kcal/mol) to the binding of SAM, which is in contrast to what had been observed for nsp16 SAM . Further, the overall polar energy was less favorable than the net non-polar contribution ( E vdW + G np = −44.5 kcal/mol). Hence, in the complex, the main force behind the binding of SAM to nsp16 is hydrophobic interactions. Overall, the binding free energy analysis showed that the binding affinity of the SAM molecule to nsp16 increases with the presence of nsp10. This implies that nsp10 acts as a stimulator to bind SAM to nsp16/nsp10 of SARS-CoV-2, which agrees with the experiment (Viswanathan et al., 2020) . Wang and coworkers obtained a similar result for SARS-CoV and MERS-CoV (Wang et al., 2015) .

Further, to explore the significant residues of nsp16 in the binding of SAM and to evaluate the differences in the complex (nsp16/nsp10) and monomer (nsp16 SAM ), the perresidue decomposition of the binding free energy was performed using the MM-GBSA approach. All the residues having > 1.5 kcal/mol of energetic contributions were considered important and shown in Figure 6 . As seen in Table 3 , the number of amino acids contributing to binding with SAM is high in the complex as compared to nsp16 SAM . Residues such as Leu100, Asp99, Cys115, Met131, Phe149, and Asp114 are common in cases of complex and nsp16 SAM . Apart from these residues, additional residues, such as Asn43, Tyr47, Gly71, Tyr132, and Ser74, also played a significant role in the binding of SAM to the complex heterodimer. These results are consistent with the higher binding affinity of SAM in the complex than the nsp16 monomer. All these residues were also considered as SAMengaging in an experimental study by Lin S. et al. (2020) . It is worth noting here that all these SAM-interacting residues are conserved both in SARS-CoV and MERS-CoV, suggesting a conserved binding mode of SAM in these viruses. The conserved SAM-binding mechanism to nsp16/nsp10 also opens the possibility of developing a pan-CoV inhibitor by targeting this SAM-binding pocket.

Next, we investigated the hydrogen bond (H-bond) interactions between SAM and nsp16 in the complex and apo forms (see Table 4 ). The % occupancy calculated from the respective MD trajectories reflected the stability of H-bonds. As seen in Table 4 , residues with > 50% H-bond occupancy in the MD simulations were recognized for the complex (nsp16/nsp10) compared to monomer nsp16 SAM . The residues, including Asp130, Tyr47, Gly71, and Asp99, showed more H-bond stability in the complex. Overall, a stable H-bond is formed between nitrogen and oxygen atoms of SAM with Asp130, Tyr47, Gly71, and Asp99 residues of nsp16, respectively. Hence, the identified significant residues like Asp130, Tyr47, Gly71, and Asp99 of nsp16 may aid in the development of SAM competitive inhibitors.

The changes in the distance of the atoms forming H-bond between SAM and nsp16 were also determined and shown in Supplementary Figure 6 . In the complex, the distance between FIGURE 5 | Sub-network representation of network 3D view of (A,B) initial structure and (C,D) final structure of the complex. Each residual C α atom is represented by a sphere where the red sphere indicates hydrophobic residues. Blue and green color spheres correspond to residues of nsp16 and nsp10, respectively. The yellow line edges represent the contact network of hydrophobic residues between nsp16 and nsp10. the oxygen atom of Asp130, Tyr47, and Gly71 of nsp16 and the nitrogen atom of SAM illustrated the average distance of 3.01, 2.89, and 3.01 Å, respectively, suggesting strong H-bond interactions. In the case of nsp16 SAM , the average distances for these atoms were 12.48, 14.64, and 11.06 Å, respectively, indicating a lack of formation of H-bonds in the nsp16 monomer. However, for Asp99, the distance between its oxygen atom (OD1 and OD2) and the oxygen atom (O3) of SAM illustrated an average distance of 2.9 Å and 3.2 Å for the complex and monomer, respectively. In contrast, the average distance between Asp99 (OD2) and SAM (O2) was higher for monomer (4.2 Å) than the complex (3.1Å). All these results were consistent with the occupancy analysis. Further, they emphasized that Asp130, Tyr47, and Gly71 of nsp16 were key residues in the binding of SAM and for the 2 ′ -O-MTase activity of nsp16.

To evaluate further the binding free energy of the heterodimer (nsp16/nsp10), G bind was estimated using the MD/MM-PBSA approach and reported in Table 5 . From Table 5 , the binding free energy ( G bind ) between nsp16 and nsp10 was found to be −47.4 kcal/mol. The various components of the total binding free energy indicate that the intermolecular van der Waals interactions ( E vdW ), electrostatic interactions ( E ele ), and nonpolar solvation free energy ( G np ) favors the binding of nsp10 and nsp16. While polar solvation free energy ( G pol ) disfavor the binding. The electrostatic interactions ( E ele ) energy is higher (−429.4 kcal/mol) than the van der Waal interactions ( E vdW ) (−90.4 kcal/mol) (see Supplementary Figure 7) . However, the disfavouring components of polar solvation energy ( G pol ) compensate for the E ele being a value of 481.7 kcal/mol. Hence, the total polar interactions ( E ele + G pol ) disfavor the binding between nsp10 and nsp16 with a value of 52.3 kcal/mol. In contrast, the total non-polar contributions ( E vdW + G np ) favor the binding between nsp16 and nsp10 (−100.0 kcal/mol). Therefore, the binding between nsp10 and nsp16 is mainly driven by hydrophobic interactions. Overall, the binding free energy analysis depicts, the binding affinity of the SAM in the complex nsp16/nsp10 ( G bind = −46.6 kcal/mol) is similar to that between nsp10 and nsp16 ( G bind = −47.4 kcal/mol) of SARS-CoV-2. These observations were in agreement with its similar homology virus, SARS-CoV (2002) . A similar binding affinity was seen between nsp10 and nsp16, as well as between SAM and the nsp16/nsp10 complex of SARS-CoV (Chen et al., 2011) .

To further investigate the critical residues involved in the binding between nsp10 and nsp16, the binding free energy was decomposed at the individual residue level using the MM-GBSA approach and shown in Figure 7 . Here, we considered only residues having the energy value ≥ 1.0 kcal/mol and listed in Supplementary Table 3 . The key residues from nsp10 include Leu45, Ala71, Val42, Met44, Tyr96, Gly69, Thr47, Arg78, Gly70, Gly94, and Pro59. Similarly, the critical residues from nsp10 include Ile40, Val104, Ala83, Val78, Met247, Val44, Gln87, Arg86, Lys76, Lys38, Val84, and Met41. It indicates that hydrophobic interactions significantly control the binding between nsp10 and nsp16. Our findings agree with a recent experimental study (Lin S. et al., 2020) where it was reported that a total of 31 

residues (Asn40, Val42, Lys43, Met44, Leu45, Cys46, Thr47, Pro59, Gly70, Ala71, Cys77, Arg78, Lys93, Gly94, and Tyr96 in nsp10 and Lys38, Gly39, Ile40, Met41, Val44, Lys76, Val78, Pro80, Ala83, Arg86, Gln87, Val104, Ser105, Asp106, Leu244, and Met247 in nsp16) mediates the intimate nsp16-nsp10 interaction. The per-residue decomposition of the binding free energy was further validated by other computational approaches, such as FTmap (Kozakov et al., 2015) . The final complex structure obtained from the MD trajectory was analyzed using FTmap, and Gln87, Arg86, Val104, Met247 of nsp16, and Leu45, Met44, Pro59, Arg78, Tyr96 of nsp10 were identified as the hotspot residues, which again agrees with the MM-GBSA analysis.

In the above finding of binding free energy analysis, we have seen that nsp10 acts as a stimulator in the binding of the SAM to nsp16. Hence, to further explore the binding interactions between nsp16 and nsp10, H-bond and hydrophobicity analysis was calculated. As seen in Table 6 , the H-bond occupancy reflects the stability of H-bond formation between protein-protein in the MD simulations. The strong H-bonds were formed between residues Ala83 (nsp16) and Tyr96 (nsp10) and between Gln87 (nsp16) and Leu45 (nsp10). Other residues, including Asp106 (nsp16), formed two H-bonds with Ala71 and Gly94 of nsp10, Lys38 (nsp16) to Lys43 (nsp10), and Ser105 (nsp16) to Lys93 (nsp10). Hydrophobic interactions also played an important role in protein-protein or protein-ligand interactions Roy et al., 2020; Sk et al., 2020d) . Different H-bonds and hydrophobic interactions from the stable structure of the nsp16-nsp10 obtained from the MD simulations were plotted via Ligplot (Wallace et al., 1995) and shown in Figure 8 . The percentage occupancy of the residual contacts with the cut-off 3.9 Å for the protein-protein is also listed in Supplementary Table 4 . Overall, these results highlighted the significant residues forming strong interactions between nsp16 and nsp10, which may aid in the development and design of inhibitors that could block these protein-protein interactions by inhibiting the 2 ′ -O-MTase activity.

The changes in the binding free energy were computed as in Equation (4) (see Method section) after replacing the residue of WT to alanine. The impact of a single mutation in the proteinprotein interaction is mainly reflected by the more positive value of the G bind . In our study, we have conducted CAS for nsp16 by considering those residues with the decomposition free energy > 1.5 kcal/mol, as given in Supplementary Table 3 . The binding free energy components calculated from the CAS mutagenesis for residues I40A, V104A, R86A, V78A, V44A, M247A, and Q87A of nsp16 are listed in Table 5 and compared with the WT. As seen in Supplementary Table 3 , although the G bind of WT is higher than mutants, different energy components of mutants follow the same trend as WT. The E ele is higher than E vdW, but disfavouring polar solvation energy G pol compensates E ele . As shown in Figure 9A , for all the systems, the binding free energy is mainly coming from the hydrophobic interactions, where the total non-polar energy ( E vdW + G np ) is higher than the total polar energy ( E ele + G pol ). Figure 9B reflects significant residues in the protein interface calculated from CAS mutagenesis. It shows that G bind for mutants, I40A, V104A, and R86A are comparatively high, suggesting that primarily these residues play a significant role in the heterodimer formation, which is in accordance with the decomposition of energy. The CAS mutagenesis results depict that in addition to hydrophobic residues I40 and V104, the hydrophilic residue R86 also plays a vital role in the binding of nsp10-nsp16.

Next, we investigated the effect of temperature and salt concentration in the binding of nsp16/nsp10 and nsp16/nsp10/SAM using the MM-PBSA method. This method presents itself to be a perfect model to study the effect of conformational dynamics of the same protein evolved in different environmental conditions to the fate of the ligand, ultimately decided by the binding affinity or dissociation.

Supplementary Tables 5, 6 displayed the binding free energy components for nsp16-nsp10 and nsp16/nsp10/SAM, respectively, at two different temperatures. As shown in Table 5 , the predicted binding free energies between nsp10 and nsp16 ( G bind ) was found to be −47.7 kcal/mol, and −52.2 kcal/mol at 300 and 310 K, respectively. It suggests that nsp16 binds more strongly to nsp10 at 310 K than 300 K. It can also be noted from Supplementary Table 5 that at the elevated temperature (310 K), the magnitude of total polar interactions ( E ele + G pol ) decreases (52.3 to 49.6 kcal/mol) and oppose less unfavorably to the proteinprotein association. At a higher temperature, the non-polar interactions ( E vdW + G np ) decreased from −100 to −102 kcal/mol and contributed more favorably to the association of nsp16-nsp10. Similarly, Supplementary Table 6 suggested that the estimated binding free energies (without entropy) of the cofactor SAM with nsp16/nsp10 ( G bind ) were −46.6 and −13.2 kcal/mol at 300 and 310 K, respectively. It was further evident from Supplementary Table 6 that after increasing the temperature, the binding affinity of SAM toward nsp16/nsp10 decreased. At room temperature (300 K) both total polar ( E ele + G pol = −2.1 kcal/mol) and non-polar ( E vdW + G np = −44.5 kcal/mol) interactions favored the binding of SAM with nsp16. After increasing the temperature to 310 K, the total polar interaction ( E ele + G pol ) became unfavorable to the binding of SAM to nsp16/nsp10. However, the total non-polar contribution ( E vdW + G np ) still favored the binding. As compared to room temperature, the strength of the non-polar ( E vdW + G np ) interactions decreases significantly (−44.5 to −23.5 kcal/mol). At 310 K, we have noticed that up to 400 ns, SAM binds strongly to nsp16. After that, SAM got detached from the binding cavity (see Supplementary Figure 8) .

The salt concentration also influenced the protein-protein interactions as well as the binding of SAM to nsp16/nsp10, which are evident from Supplementary Tables 5, 6. As shown in Supplementary Table 5, the predicted binding free energies of nsp10 with nsp16 ( G bind ) were found to be −44.8 and −37.7 kcal/mol, for the salt concentration of 0.15 and 0.25 M, respectively. After increasing the salt concentration, it was observed that the total polar ( E ele + G pol ) and non-polar interactions ( E vdW + G np ) increased from 54.1 to 55.3 kcal/mol and −98.9 to −93 kcal/mol, respectively. The rise in unfavorable interactions caused the destabilization of the protein-protein interactions. Similarly, the binding of SAM to the heterodimer was found to be affected by the increased salt concentration (see Supplementary Table 6 ). The estimated binding free energy of SAM ( G bind ) was −23.7 and −22.1 kcal/mol for the salt concentration of 0.15 and 0.25 M, respectively. These values are significantly higher than what was obtained for the neutral salt concentration (−46.6 kcal/mol). E vdW and E elec were found to contribute less favorably to the binding of SAM to nsp16/nsp10 in comparison to neutral simulations. A similar trend was observed when the temperature was increased from 300 to 310 K. The time evolution of CoM distance between nsp16 and SAM was displayed in Supplementary Figure 9 . It is evident from Supplementary Figure 9 that with the increased salt concentration, the CoM distance between nsp16 and SAM increased and remained stable around 17 Å and 18 Å for 0.15 and 0.25 M, respectively. This implies that the cofactor SAM moves away from the binding cavity with the increased salt concentration.

Herein, we have employed extensive MD simulations of 1 µs along with the molecular mechanics/Poisson-Boltzmann surface area (MM/PBSA) method to study the binding mechanism of the cofactor SAM to the nsp16/nsp10 heterodimer and sub-unit nsp16 alone (nsp16 SAM ) of SARS-CoV-2. Our MD/MMPBSA calculations suggested that SAM binds strongly to the nsp16/nsp10 heterodimer and fails to bind to the nsp16 monomer. It emphasized that nsp10 helps in the strong interaction between SAM and nsp16 to execute the 2 ′ -O-MTase activity. This finding agrees well with the experimental study, and a similar observation was made for other coronaviruses, including SARS-CoV and MERS-CoV. We have also investigated the interactions between nsp16 and nsp10. Our study revealed that the binding of nsp10 stabilizes the complex (nsp16/nsp10) structure. Further, our study showed that hydrophobic interactions are critical for the heterodimer association. Thus, apart from the active site of nsp16, the interface of nsp16/nsp10 can also be considered as a potential target site in the design of antiviral drugs such as peptide inhibitors. It includes Ile40, Val104, Ala83, Val78, Met247, Val44, Gln87, Arg86, Lys76, Lys38, Val84, and Met41 from nsp16, and Leu45, Ala71, Val42, Met44, Tyr96, Gly69, Thr47, Arg78, Gly70, Gly94, and Pro59 from nsp10. Besides, the stable hydrogen bond between Ala83 (nsp16) and Tyr96 (nsp10), and between Gln87 (nsp16) and Leu45 (nsp10) were important in the nsp16-nsp10 interface. The CAS study reveals that residues I40A, V104A, R86A, V78A, V44A, M247A, and Q87A of nsp16 were considered as hot spot residues for the association of nsp16-nsp10.

G bind for mutants, I40A, V104A, and R86A are comparatively high, suggesting the significance of these residues. Finally, we investigated the effect of temperature and salt concentration on the binding of SAM to the heterodimer. Our study revealed that the SAM binding was impaired due to an increase in temperature or salt concentration. Finally, our study provides a comprehensive understanding of the dynamic and thermodynamic process of binding nsp16 and nsp10 that will contribute to the rational design of inhibitors targeting nsp16/nsp10.

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

An efficient method for sampling the essential subspace of proteins

Severe neurologic syndrome associated with Middle East respiratory syndrome corona virus (MERS-CoV)

Molecular dynamics with coupling to an external bath

Potential for global spread of a novel coronavirus from China

Calculations of solute and solvent entropies from molecular dynamics simulations

NAPS update: network analysis of molecular dynamics data and protein-nucleic acid complexes

Biochemical and structural insights into the mechanisms of SARS coronavirus RNA ribose 2 ′ -O-methylation by nsp16/nsp10 protein complex

Prediction of the SARS-CoV-2 (2019-nCoV) 3C-like protease (3CL pro) structure: virtual screening reveals velpatasvir, ledipasvir, and other drug repurposing candidates

Coronaviruses: important emerging human pathogens

Particle mesh Ewald: An N· log (N) method for Ewald sums in large systems

Infidelity of SARS-CoV Nsp14-exonuclease mutant virus replication is revealed by complete genome sequencing

An RNA cap (nucleoside-2 ′ -O-)-methyltransferase in the flavivirus RNA polymerase NS5: crystal structure and functional characterization

Coronaviruses: an overview of their replication and pathogenesis

The energy landscapes and motions of proteins

Viral and cellular mRNA capping: past and prospects

Insights into protein-protein binding by binding free energy calculation and free energy decomposition for the Ras-Raf and Ras-RalGDS complexes

Clinical Characteristics of Coronavirus Disease 2019 in China

Identification of severe acute respiratory syndrome coronavirus replicase products and characterization of papain-like protease activity

Pairwise solute descreening of solute charges from a dielectric medium

Parametrized models of aqueous free energies of solvation based on pairwise descreening of solute atomic charges from a dielectric medium

Assessing the performance of the MM/PBSA and MM/GBSA methods. 1. The accuracy of binding free energy calculations based on molecular dynamics simulations

Computational studies of a mechanism for binding and drug resistance in the wild type and four mutations of HIV-1 protease with a GRL-0519 inhibitor

Clinical features of patients infected with 2019 novel coronavirus in Wuhan

Collective motions in proteins: a covariance analysis of atomic fluctuations in molecular dynamics and normal mode simulations

Fast, efficient generation of highquality atomic charges. AM1-BCC model: II. Parameterization and validation

An Introduction to Hydrogen Bonding

Investigating specificity of the antihypertensive inhibitor WNK463 against With-No-Lysine kinase family isoforms via multiscale simulations

Investigating phosphorylationinduced conformational changes in WNK1 kinase by molecular dynamics simulations

A comparative study of structural and conformational properties of WNK kinase isoforms bound to an inhibitor: insights from molecular dynamic simulations

Energetic basis for drug resistance of HIV-1 protease mutants against amprenavir

Energetics of mutation-induced changes in potency of lersivirine against HIV-1 reverse transcriptase

Origin of decrease in potency of darunavir and two related antiviral inhibitors against HIV-2 compared to HIV-1 protease

Mutation-induced loop opening and energetics for binding of tamiflu to influenza N8 neuraminidase

Importance of polar solvation for cross-reactivity of antibody and its variants with steroids

Importance of polar solvation and configurational entropy for design of antiretroviral drugs targeting HIV-1 protease

Dispersion terms and analysis of size-and charge dependence in an enhanced poisson-Boltzmann approach

Systematic study of the boundary composition in Poisson Boltzmann calculations

Method for estimating the configurational entropy of macromolecules

Virus Taxonomy: 9th Report of the International Committee on Taxonomy of Viruses

Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models

The FTMap family of web servers for determining and characterizing ligand-binding hot spots of proteins

A fast SHAKE algorithm to solve distance constraint equations for small molecules in molecular dynamics simulations

A simple method for displaying the hydropathic character of a protein

Crystal structure of SARS-CoV-2 nsp10/nsp16 2 ′ -O-methylase and its implication on antiviral drug design

Novel coronavirus pneumonia outbreak in 2019: computed tomographic findings in two cases

Langevin dynamics of peptides: the frictional dependence of isomerization rates of N-acetylalanyl-N ′ -methylamide

Drug treatment options for the 2019-new coronavirus (2019-nCoV)

Genomic characterisation and epidemiology of 2019 novel coronavirus: implications for virus origins and receptor binding

Animal-to-human SARS-associated coronavirus transmission?

Split the charge difference in two! a rule of thumb for adding proper amounts of ions in MD simulations

ff14SB: improving the accuracy of protein side chain and backbone parameters from ff99SB

Computational alanine scanning to probe protein-protein interactions: a novel approach to evaluate binding free energies

Angstrom Resolution Crystal Structure of NSP16 -NSP10 Complex from SARS-CoV-2

Computational alanine scanning mutagenesis-an improved methodological approach

PROPKA3: consistent treatment of internal and surface residues in empirical p K a predictions

An analysis of the accuracy of Langevin and molecular dynamics algorithms

Recent Discovery and Development of Inhibitors Targeting Coronaviruses

A modified TIP3P water potential for simulation with Ewald summation

A computational exercise illustrating molecular vibrations and normal modes

PTRAJ and CPPTRAJ: software for processing and analysis of molecular dynamics trajectory data

The crystal structure of nsp10-nsp16 heterodimer from SARS-CoV-2 in complex with S-adenosylmethionine

Investigating conformational dynamics of lewis Y oligosaccharides and elucidating blood group dependency of cholera using molecular dynamics

Insights into the binding recognition and susceptibility of tofacitinib toward janus kinases

Comparative therapeutic efficacy of remdesivir and combination lopinavir, ritonavir, and interferon beta against MERS-CoV

Insight into binding mechanisms of inhibitors MKP56, MKP73, MKP86, and MKP97 to HIV-1 protease by using molecular dynamics simulation

Plant-derived natural polyphenols as potential antiviral drugs against SARS-CoV-2 via RNA-dependent RNA polymerase (RdRp) inhibition: an in-silico analysis

Exploring the energetic basis of binding of currently used drugs against HIV-1 subtype CRF01_AE protease via molecular dynamics simulations

Computational investigation of structural dynamics of SARS-CoV-2 methyltransferasestimulatory factor heterodimer nsp16/nsp10 bound to the cofactor SAM

Elucidating biophysical basis of binding of inhibitors to SARS-CoV-2 main protease by using molecular dynamics simulations and free energy calculations

Exploring the potency of currently used drugs against HIV-1 protease of subtype D variant by using multiscale simulations

Structural basis of RNA cap modification by SARS-CoV-2

LIGPLOT: a program to generate schematic diagrams of protein-ligand interactions

Automatic atom type and bond type perception in molecular mechanical calculations

Development and testing of a general amber force field

Coronavirus nsp10/nsp16 methyltransferase can be targeted by nsp10-derived peptide in vitro and in vivo to reduce replication and pathogenesis

Fast, approximate algorithm for detection of solvent-inaccessible atoms

Coronavirus pathogenesis

Coronavirus pathogenesis and the emerging pathogen severe acute respiratory syndrome coronavirus

Coronavirus genomics and bioinformatics analysis. Viruses

Analysis of therapeutic targets for SARS-CoV-2 and discovery of potential drugs by computational methods

Fast and accurate computation schemes for evaluating vibrational entropy of proteins

Design of wide-spectrum inhibitors targeting coronavirus main proteases

Epidemiology and clinical characteristics of human coronaviruses OC43, 229E, NL63, and HKU1: a study of hospitalized children with acute respiratory tract infection in Guangzhou, China. Eur

A pneumonia outbreak associated with a new coronavirus of probable bat origin

Networkbased drug repurposing for novel coronavirus 2019-nCoV/SARS-CoV-2

This manuscript has been released as a pre-print at ChemRxiv (Sk et al., 2020b) .

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb. 2020.590165/full#supplementary-material Conflict of Interest: The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.Copyright © 2020 Sk, Jonniya, Roy, Poddar and Kar. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.