Identification of lead compounds through structure-based virtual screening

Structure-based virtual screening enables the prediction of optimal interactions between ligands and a macromolecular target for complex formation. The ligands are subsequently sorted according to their binding free energy for the target. This requires the three-dimensional structure of the target, with the compounds obtained from a database and categorised according to their affinity. In the present study, we downloaded a subset of natural compounds (n = 97,999) from the ZINC database for virtual screening against PRRSV Nsp4. We subsequently identified the top 10 compounds sorted according to their minimum binding free energy (range:− 10.0 to − 9.2 kcal/mol) for further analyses (Table 1).

Table 1 Binding free energies of the top 10 screened compounds along with the amino acid residues involved in interactions. The amino acid residues shown in bold are involved in hydrogen-bonding interactionsAnalysis and visualisation of the screened compounds

The active site residues of Nsp4 include His39, Asp64, and Ser118, as well as His133 and Ser136 that are reportedly essential for protein activity. The compound showing the optimal binding free energy (− 10.0 kcal/mol; ZINC38167083) demonstrated ligand interactions such as His39-mediated van der Waals interactions and Asp64-mediated Pi-anion interaction. Fig. 1 shows protein–ligand H-bond interactions involving the active site residues His39, Ser118, and Ser136, besides, Asp64 and His133 is interacted with van der Waals and pi-anion interaction with ZINC08877407. The interacting amino acid residues of top 10 compounds are shown in Table 1.

Fig. 1

2D representation of the binding interactions of top five screened natural compounds with Nsp4 depicted key amino acid residues contributed in protein-ligand interactions. A ZINC38167083, (B) ZINC16919178, (C) ZINC08792350, (D) ZINC01510656, and (E) ZINC08877407

Assessment of drug likeness through physicochemical property analysis

Physicochemical property analysis is one of the fundamental tasks in any drug discovery program. The top 10 screened compounds were then subjected to analysis of their physicochemical properties according to 7 principal descriptors (MW, logP value, status of HBDs and HBAs, 2D PSA, P, and VWSA). According to a previous study, a good drug should have an MW < 500 Da, an HBD < 5, and an HBA < 10. The MW, logP, HBD, and HBA of the selected compounds met the Lipinski’s rule. Additionally, PSA, P, and VWSA results displayed drug-like behaviour. (Table 2).

Table 2 Physicochemical properties of the top 10 screened compoundsMD simulation analysis

The structure of Nsp4 and top five screened compound-complex with Nsp4 was employed for 100 ns MD simulation study for predicting the dynamic changes during protein-ligand interaction and their nature of stability. The present study included various parameters i.e. RMSD, RMSF, Rg, SASA, H-bond, PCA, Gibbs free energy landscape, and binding free energy calculation.

Structural deviation analysis through RMSD

The RMSD value describes the dynamic behaviour among native structures to a new pose. After a 70 ns of simulation to obtain a stable trajectory, the RMSD values were 0.35, 0.25, 0.29, 0.23, 0.38, and 0.39 nm for Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407, respectively. These data suggest that Nsp4-ZINC08792350 and Nsp4-ZINC38167083 are highly stable complexes relative to the others. Because each Nsp4–compound complex demonstrated stability after the 70 ns simulation, we performed further evaluations on each for last 30 ns trajectory (Fig. 2A).

Fig. 2

Stability analysis (A) RMSD values for the Nsp4–compound complexes. Flexibility analysis (B) RMSF values for the Nsp4–compound complexes over the final 30 ns of the simulations. Compactness (C) Rg, and Solvent accessible surface area analysis (D) SASA values for the final 30 ns of the simulations. Black, red, green, blue, orange, and violet colours represent Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407, respectively. E Changes in the number of hydrogen bonds in each respective complex according to data from the final 30 ns of the simulations. Red, green, blue, orange, and violet colour represent Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407 respectively

Flexibility analysis through RMSF

Evaluation of the RMSF values used to assess structural rigidity revealed values of 0.08, 0.11, 0.12, 0.10, 0.11, and 0.11 nm for Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407, respectively (Fig. 2B). Higher RMSF values were due to ligand binding, causing alterations in protein geometry. Minimal fluctuations were observed in Nsp4-ZINC08792350 and Nsp4-ZINC38167083 complex compared with that in other complexes.

Radius of gyration (Rg) analysis

Assessment of complex compactness according to Rg calculation revealed values of 1.50, 1.27, 1.43, 1.46, 1.39, and 1.44 nm for Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407, respectively (Fig. 2C). The results indicate that the Nsp4-ZINC38167083 complex showed a more compact structure than the other complexes.

Solvent accessible surface area (SASA) analysis

To identify changes in the solvent-accessible regions of the complexes, we determined SASA values over the course of the final 30 ns of the simulation. Our study revealed values of 95.88, 98.33, 98.98, 97.13, 96.97, and 100.92 nm2 for Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407 (Fig. 2D), revealing relatively minimal changes after binding by each of the compounds.

Interaction analysis through hydrogen bonding

Hydrogen bonding is the most important bond for stabilizing protein–ligand interactions. The average number of hydrogen bonds for the complexes Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407 over the final 30 ns of the simulations was 0–1 and that for Nsp4-ZINC38167083 and Nsp4-ZINC16919178 was 0–2 and 0–3, respectively (Fig. 2E). Hence, these compounds interacted with Nsp4 and provided a stable complex during protein–ligand interactions.

Principal component analysis (PCA)

In PCA, the sum of the eigenvalues suggests the overall flexibility of a structure under different conditions. Therefore, the first 5 of 50 eigenvectors used to calculate eigenvalues from the final 30 ns of the simulation were used to determine the percentage change in structural movement. The results revealed that these five eigenvectors accounted for 42.85, 63.97, 63.27, 59.14, 64.83, and 71.05% of the motions for Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407 respectively (Fig. 3A), suggesting increased movement after the binding of each ligand. Moreover, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, and Nsp4-ZINC01510656 showed less overall motion relative to Nsp4-ZINC08877407. Additionally, generation of a 2D plot for assessing protein dynamics after ligand binding suggested the overall stability (lowcorrelated motions) of Nsp4, Nsp4-ZINC38167083, and Nsp4-ZINC08792350(Fig. 3B), indicating these compounds as possible leads for further evaluation as inhibitors.

Fig. 3

Principal component analysis (A) Eigenvalues derived from the final 30 ns of each simulation and used for PCA depicted Eigenvalues vs. first fifty eigenvector. B First two eigenvectors depicted Nsp4 motion in space for all the systems. Black, red, green, blue, orange, and violet colours represent Nsp4, Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407 respectively

Gibbs free energy landscape

We then calculated the Gibbs free energy landscape using the first two principal components (PC1 and PC2) in order to visualize the results. Fig. 4 shows the colour-coded plots generated for Nsp4 along with each complex. The lowest free energy values (≤9.08 kJ/mol) were observed for Nsp4-ZINC38167083, suggesting that this complex demonstrated overall thermodynamic stability. The other complexes (Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407) had values of to 11.4 kJ/mol, implying that these complexes have numerous high-energy minima.

Fig. 4

The color-coded illustration of the Gibbs free energy landscape plotted using PC1 and PC2. The lower energy systems are represented by the deeper blue color on the contour map. A Nsp4, (B) Nsp4-ZINC38167083, (C) Nsp4-ZINC16919178, (D) Nsp4-ZINC08792350, (E) Nsp4-ZINC01510656, and (F) Nsp4-ZINC08877407

Binding free energy

We then evaluate the binding free energy associated with each ligand through MM-PBSA using the final 10 ns of the simulation, for calculation of van der Waals and electrostatic interactions, Polar solvation, and SASA. The calculated binding free energy for Nsp4-ZINC38167083, Nsp4-ZINC16919178, Nsp4-ZINC08792350, Nsp4-ZINC01510656, and Nsp4-ZINC08877407 was − 124.54, − 128.44, − 159.33, − 122.50, and − 78.19 kJ mol− 1 respectively (Table 3).

Table 3 Average binding free energies of Nsp4 complexes in kJ mol − 1

The investigation of residual binding energy is a key method for identifying residues important to ligand binding. Fig. 5 shows that amino acid residues at positions 5 to 142 contributed significantly to binding of ZINC38167083, ZINC16919178, ZINC08792350, and ZINC01510656, which are the catalytic residues in the active site. Fewer contacts were observed in relation to ZINC08877407 binding, suggesting that ZINC38167083, ZINC16919178, ZINC08792350, and ZINC01510656 represent potential Nsp4 inhibitors.

Fig. 5

Plot depicting the amino acid residues of Nsp4 contributing to the binding with natural compounds. Red, green, blue, orange, and violet colours represent ZINC38167083, ZINC16919178, ZINC08792350, ZINC01510656, and ZINC08877407, respectively