Identification of novel STAT3 inhibitors for liver fibrosis, using pharmacophore-based virtual screening, molecular docking, and biomolecular dynamics simulations (2024)

Abstract

The signal transducer and activator of transcription 3 (STAT3) plays a fundamental role in the growth and regulation of cellular life. Activation and over-expression of STAT3 have been implicated in many cancers including solid blood tumors and other diseases such as liver fibrosis and rheumatoid arthritis. Therefore, STAT3 inhibitors are be coming a growing and interesting area of pharmacological research. Consequently, the aim of this study is to design novel inhibitors of STAT3-SH3 computationally for the reduction of liver fibrosis. Herein, we performed Pharmacophore-based virtual screening of databases including more than 19,481 commercially available compounds and in-house compounds. The hits obtained from virtual screening were further docked with the STAT3 receptor. The hits were further ranked on the basis of docking score and binding interaction with the active site of STAT3. ADMET properties of the screened compounds were calculated and filtered based on drug-likeness criteria. Finally, the top five drug-like hit compounds were selected and subjected to molecular dynamic simulation. The stability of each drug-like hit in complex with STAT3 was determined by computing their RMSD, RMSF, Rg, and DCCM analyses. Among all the compounds Sa32 revealed a good docking score, interactions, and stability during the entire simulation procedure. As compared to the Reference compound, the drug-like hit compound Sa32 showed good docking scores, interaction, stability, and binding energy. Therefore, we identified Sa32 as the best small molecule potent inhibitor for STAT3 that will be helpful in the future for the treatment of liver fibrosis.

Subject terms: Cancer, Computational biology and bioinformatics, Drug discovery

Introduction

Many studies have shown that signal transducers and activators of transcription 3 (STAT3) a crucial member of the STAT family of cytoplasmic proteins that transmit signals to the nucleus (in response to cytokines and growth factors) there by regulating growth, development, and survival of the cell1. STAT3 was identified as an acute-phase response factor that is linked with cancer and other diseases. In healthy cells, the STAT3 remains inactive, while it activates when the cell is experiencing unfavorable circumstances including cancer and other disease conditions. It is a transcription factor and is encoded in humans as the STAT3 gene. STAT3 consists of seven members (STAT 1, 2, 3, 4, 5A, 5B, and 6) encoded by different chromosomal clustered genes2. In addition to this Jiang-Jiang Qinet al suggested that STAT also consists of seven domains such as coiled-coil, SH2, NH2-terminus, DNA-binding, linker, and C-terminal Transactivation3. Among the STAT3 domains, the SH2 domain is the most significant. At Try705, C-terminal Transactivation activates STAT3 by generating serine and tyrosine development at position 727. The SH2 domains consist of 100 amino acidsthat control the tyrosine phosphorylation of STAT3 receptors and homo-dimerization or hetero-dimerization of STAT3 monomer bindings. Both are vital for gene expression as well as for DNA binding. While, SH2-, DBD, and N-terminal all have many functions, and are therefore considered therapeutic targets for many neurotic and pathological conditions4.

STAT3 is essential for cell development, growth, and regulation. Therefore, there are potential ways to inhibit STAT3 protein directly. On the basis of their structure and function the N-terminal, DNA binding domain and SH2 domain can be selectively targeted in order to design the STAT3 inhibitors. The SH2 components of the STAT3 are considered to play the most fundamental role in STAT3 activation. Among the different domains, the SH2 domain is preferentially targeted, because SH2 enables STAT3 to recruit it to any activated receptor. Furthermore, the inhibition of SH2 blocks STAT3 dimerization, thereby inhibiting nuclear translocation and STAT3-dependent gene regulation5.

TGFβ1, a STAT3 target gene, promotes liver fibrosis by increasing the expression of TGF-β1 in a liver fibrosis mouse model and a rapid early STAT3 activation was observed in HSCs and fibroblasts, but not in normal hepatocytes. This suggests that STAT3 could be a target for liver disorder. Moreover, STAT3 inhibition might reduce angiogenesis and further fibrogenesis thus, indicate STAT3 as a potential therapeutic target in liver disease6,7.

To overcome STAT3-induced diseases, several classes of direct low molecular weight inhibitors of STAT3 have been discovered that prevent STAT3-DNA or STAT3 dimerization binding. These inhibitors include oligonucleotides, peptides and peptidomimetics. However, it is challenging to develop highly potent and selective inhibitors of STAT3 against the liver fibrosis. In an attempt, computational screening precisely identified novel STAT3 SH2 inhibitors that block dimerization of STAT3, consisting of various chemical compounds such as STA-219, S3I-201 and STX-011911 which relied on docking of chemical libraries1. For chemical induced liver injury, some STAT3 inhibitors designed as chemical probes which are more affected in therapeutic sense. STX-0119 is reported as a more potent in controlling the improvement of CCl4 induced liver fibrosis via decrease the activated HSCs. Some other inhibitor like HJC0123was proven to overcome liver fibrosis by inhibiting the phosphorylation and block the process of STAT3 activation8.

Sorafenib is multi-kinase inhibitor and is considered as a potential inhibitor for liver fibrosis. Sorafenib has the potential to inhibit STAT3 phosphorylation9. It is also reported that Sorafenib has inhibited HSC activation which is also helpful in anti-fibrotic activity, growth and collagen accumulation in vitro.

Sorafenib plays a vital role in inhibiting STAT3 by directly binding to the SH2 domain that is SHP-1 is consisting of 2 Src homology domains, phosphorylated tyrosine binding, C-terminal tail and PTP catalytic domain. SHP-1 activity is linked closely with its structural changeability. SHP-1 is having auto inhibited ability but when the N-terminal SH2 domain on surface into the catalytic domain, blocking the active site of the catalytic pocket. The WPD loop containing the remaining active site residue Asp 421 was affected with this residue and the catalytic activity of SHP-1. The negative regulator of STAT3 is SHP-1. The anti-HCC effect of Sorafenib is due to the formation of important salts in the N-SH2 domain, which becomes the open catalytic PTP domain and liberates SHP-1. Sorafenib and its derivatives SC-1, SC-40 and SC-43 showed similar SHP-1 reactivation and STAT3 signal inhibition in HCC cells10.

In this study, Pharmacophore model of Sorafenib was generated using MOE software. The validation of Pharmacophore model was done by Gunner Henry method. Validated model was further used for virtual screening of different databases like in-house, Zinc, Phytochemicals and South African natural compounds database. In order to evaluate drug-likeness, screened hits were filtered with the help of SwissADME (http://www.swissadme.ch/) online server. After that hits were docked against STAT3-SH2. From docking study, top 05 complexes were selected for MD simulation. The selection criterion of the top 05 complexes for MD simulation was based on the S score and interactions with receptor. Finally, MD simulation analysis such as RMSF, RMSD, Rg, and DCCM shows that the new predicted inhibitors of STAT3 like Sa32, ZINC47009207, 5459840, and SANC00347 are more stable than control.

Materials and methods

Pharmacophore modeling and virtual screening

A Pharmacophore model was developed in the MOE (Molecular Operating Environment) software11 for the virtual screening of synthetic and commercially available compounds databases. The model was then validated and the validated Pharmacophore model as 3D query was utilized for the virtual screening. To identify and retrieve the potent anti-STAT3 compounds, different databases including in-house, ZINC, Phytochemicals, and South African databases were screened12. The screened hits best fitted on the applied Pharmacophore were further filtered by ADMET properties to short-list compounds having drug-likeness. In order to determine the pharmacokinetic properties of the retrieved compounds, the Lipinski’s rule of five (RoF) was applied. Determination of drug-likeness property is particularly used to demonstrate new lead hits by screening of compound databases13.

Molecular docking of hit compounds and STAT3

The possible interaction of protein-ligands complex is explored by molecular docking study. Molecular docking helps topredict intermolecular framework formed among small molecule and protein, suggesting interaction modes leading to protein inhibition. In this study, all the drug-like retrieved hits were docked into the active site of STAT3. For docking purpose, the docking protocol implemented in MOE 2016 software was used14. For docking of STAT3, the structure of its SH2 domain (PDBID:1BG1) was taken from protein data bank https://www.rcsb.org/15 This SH2 domain (comprised of 499–688 amino acid sequence) of STAT3 is phosphorylated16. STAT3 receptor and its reference i.e., standard inhibitor (Sorafenib) was initially Protonated and subjected to energy minimization with default parameters (gradient: 0.05, Force Field: MMFF94x) of MOE 2016.08 software17. The energy converged system was subjected to docking. The result of docking was analyzed on the basis of protein/compounds interactions and docking scores. In order to validate the stability and interaction of drug-like hits in the binding site of STAT3-SH2 domain, the complexes were subjected to molecular dynamics simulations.

Molecular dynamics simulations

MD simulation of the best five docking score compounds and the standard compound (Sorafenib) were carried out to unveil the stability of the protein–ligand complexes. Molecular dynamic simulations and analysis for each complex was performed in AMBER 2018 using the force field (ff14SB). During the MD protocol, each system was neutralized by utilizing the LEAP module counter ions such as Na+ and Cl. TIP3P water model was filled in octahedral box with 10Å parameter to solvateeach system14. For long-range electrostatic interactions a distance of 10Å was used. Additionally, the Particle Mesh Ewald (PME) algorithm was used to deal with long-range electrostatic interactions. By using the SHAKE option and the Particle Mesh Ewald (PME) method, covalent bonds involving hydrogen atoms were optimized in each complex system16. The temperature was controlled by Langevin dynamics. Finally, 100ns MD simulations (two replicas) for each complex were performed on the GPU using CUDA version of PMEMD. For the analyses of MD trajectory, Amber 20's CPPTRAJ module was used.

ADME predictions

Evaluation of absorption, distribution, metabolism, and excretion (ADME) are the main criteria to recommend a molecule as a drug-like candidate. Hence, computational approaches are utilized to predict the ADME properties of chemical molecules and advocate their drug-likeness. In this study, we used SwissADME http://www.swissadme.ch online server19,20 and assessed the ADME properties such as GIT absorption, bioavailability, and solubility profile of the selected compounds21.

Binding free energy calculation

The molecular Mechanic/Poisson-Boltzmann Surface Area (MM-PBSA) approach22 is an efficient and reliable method to model molecular recognition, such as for protein–ligand binding interactions. To estimate the strength of binding of STAT3 in complex with either the reference compound or any of the candidate hits, the binding free energy was computed by employing the MMPBSA.py script. The MMPBSA was computed for the last 500 frames (10ns) of each simulated complex. The below equation was used to calculate the free energy18

Gbind=ΔGcomplex-[ΔG receptor+ΔG ligand]

The equations describe G bind shows total binding free energies, ΔG complex usefor complex free energies, and the remaining terms stand for the corresponding free energies of the receptor protein and ligand.

Results

Pharmacophore model generation

Ligand-based Pharmacophore model was generated. Sorafenib was selected as a reference for Pharmacophore model generation. Pharmacophorewas generated from a total of seven (7) features. Thus, the resulted Pharmacophore was comprised of two hydrogen bond acceptor (HBA), two hydrogen bond donor (HBD), one aromatic (Aro) and two hydrophobic (HY) features (Fig.1). Finally, in-house, ZINC, phytochemicals and South African compound libraries were screened based on the generated Pharmacophore. As a result, many compounds were found that interacted more strongly with the targeted compound than the reference compound.

Figure 1.

Identification of novel STAT3 inhibitors for liver fibrosis, using pharmacophore-based virtual screening, molecular docking, and biomolecular dynamics simulations (1)

Open in a new tab

Pharmacophore model validation

Internal database was used in this study as test set forPharmacophore model validation. The test set contained 852 inactive compounds for STAT3 and the rest of 48 molecules were known inhibitors. The validated model has the ability to differentiate between the active and inactive molecules. This model was used to investigateand carry out virtual screening. Many important parameters were calculated like total hits (Ht), active hits (Ha), % yield of actives, % ratio of actives, enrichment factor (E), and goodness-of-hit score (GH) which are placed in (Table 1).

Ha3A+Ht/4HtA1-Ht-Ha/D-A

Table 1.

Pharmacophore model evaluation based on Gunner-Henry (GH) scoring method.

S. noParametersPh4 model
1Number of compounds in database (D)900
2Number of active in database (A)48
3Total Hits (Ht)50
4Active Hits (Ha)40
5% Yield of actives[(Ha/Ht) × 100]80
6% Ratio of actives [(Ha/A) × 100]83
7Enrichment factor (E) [(Ha × D)/(Ht × A)]15
8False negatives [A − Ha]8
9False positives [Ht − Ha]10
10GH score0.79

Open in a new tab

GH Score of 0.7–0.8 indicates a validated model23. The model was used for further virtual screening purpose.

The higher the E value, the better the model's ability to identify active compounds. With an E-value of 15, the model identified 40 active hits out of 900 compounds screened, indicating that the model was able to differentiate active frominactive molecules. A Gunner Henry’s score of 0.7–0.8 preferred as a good model. So the goodness score for Pharmacophore model was 0.79. The validation of our generated Pharmacophore suggested that the resultant Pharmacophore was good enough to be used for virtual screening of different databases against the identification of STAT3 candidate hit compounds.

Virtual screening

The validated model was utilized to screen in-house database (1700 molecules), zinc database (12,000 compounds), phytochemicals (5000 compounds) and South African database (1000 compounds). The Pharmacophore model retrieved340 compounds as the best fitted compounds on the Pharmacophore features against the STAT3 candidate hits molecules.Identification of novel STAT3 inhibitors for liver fibrosis, using pharmacophore-based virtual screening, molecular docking, and biomolecular dynamics simulations (2)

Molecular docking

The ADMET evaluation of the Pharmacophore-retrieved hits filtered a total of 278 hits as drug-like molecules which were further docked in the active site ofSTAT3 which is comprised of Arg609, Arg595, Lys591, Arg609, Ser611, Glu612, Ser613, Gln633, Ile634, Gln635, Ser636, Val637, Gln638, Met586, Gly587, Phe588, Ilu589, and Ser590 residues. The top 25 ligands of STAT3 were selected according to the S-score and interaction type.

Table 2 suggests the S-Score and interacting residues of the top 25 hits with STAT3. Figure2 represents the interactions of best hits and the reference compound with STAT3 receptor.

Table 2.

The docking scores of top 25 ligands and their interactions with SH2 domain of STAT3.

S. noCompound IDDocking scoreInteracting residuesBond typeDistance
1Ref−7.2334GLU 638H-donor3.07
SER 613H-acceptor3.18
GLN 635H-acceptor3.53
25459840−11.2812ARG 518H-acceptor2.85
SER 521H-acceptor3.17
ARG 518H-acceptor3.04
ARG 593H-acceptor2.84
ARG 518Ionic2.85
ARG 518ionic3.98
3SANC00347−10.2663VAL 637H-donor2.89
LYS 591H-acceptor3.15
SER 613pi-H4.22
4ZINC47009207−9.2654GLU 523H-donor2.92
ARG 518H-acceptor3.15
SER 521H-acceptor3.17
ARG 593H-acceptor3.29
5sa32−9.1425GLU 523H-donor3.71
ARG 593H-donor4.38
ARG 593H-acceptor3.01
MET 586Pi-H3.1
6SANC00352−10.1882LYS 591H-acceptor3.48
SER 636H-acceptor3.23
GLU 638pi-H4.29
7SANC00512−9.9356GLU 594H-donor2.80
GLN 635H-acceptor3.21
LYS 626H-acceptor3.21
8442659−9.4941GLU 523H-donor3.58
ARG 593H-acceptor2.97
TYR 674Ionic2.92
910885340−9.7567GLU 612H-donor3.12
LYS 591H-acceptor2.86
SER 613H-acceptor2.99
10ZINC44978511−9.1258MET 586H-donor4.17
GLU 530H-donor3.07
ARG 593H-acceptor3.30
11ZINC23543539−9.0425GLU 523H-donor2.83
GLU 530H-donor3.28
ARG 593H-acceptor2.93
12SANC00856−9.4064GLU 594H-donor2.79
LYS 591H-acceptor3.05
SER 613H-acceptor3.03
13SANC01083−8.2785ILE 634H-donor3.01
GLU 594H-donor2.91
SER 636H-acceptor3.12
146-KA6-TP−8.1504MET 586H-donor3.81
GLU 523H-donor2.84
GLU 523pi-H3.73
15ZINC48237730−8.3923MET 586H-donor3.78
MET 586H-donor3.75
GLU 523pi-H4.36
1616736655−8.1621GLU 523H-donor3.36
ARG 593H-acceptor ionic3.17
GLU 523H-acceptor ionic3.81
179957384−7.9068ARG 595H-acceptor2.90
LYS 591Ionic2.86
SER 636pi-H4.32
18ZINC20540819−7.8776MET 586H-donor3.99
GLU 523pi-H3.76
MET 586pi-H4.60
197-KA7-TP−7.8152MET 586H-donor3.79
GLU 523H-donor2.83
GLU 523pi-H3.72
208.3-Formylcoumarin−7.7606GLU 594H-donor2.77
GLU 594H-donor2.89
GLU 638H-acceptor3.21
21BB-III-18−7.7005GLU 523H-donor2.98
ARG 593H-acceptor3.10
TYR 674pi-H3.92
22Sa5−8.6211GLU 594H-donor2.89
ILE 634Pi-H4.78
SER 636Pi-H4.38
23Sa10−7.7093MET 586H-donor3.51
GLU 523H-donor3.64
TYR 674Pi-H3.51
24Sa12−9.7339MET 586H-donor3.97
TYR 674H-pi4.26
GLU 523Pi-H3.70
25Sa15−8.06696TYR 674H-pi4.23
MET 586Pi-H4.15
PRO 675Pi-H3.94

Open in a new tab

Figure 2.

Open in a new tab

Among the 25 docked compounds, 5459840, SANC00347, ZINC47009207 and Sa32 were the most potent with docking score  −11.2812, −10.2663, −9.2654 and −9.1425, respectively. Most of the compounds showed the most accurate binding and interaction with the receptors compared to reference compound. The docking index for Sorafenib was −7.2334.

ADMET properties determination

Table 3 shows the drug likeness of the molecules 5459840; ZINC47009207, Sa32 and SANC00347 identified by online ADME/Tox tool of Swiss-ADME using the SMILES nomenclature. To show the ADME/Tox profile, we decided on maximum fundamental ADME/Tox functionality from online server. Table 3 displays the pharmacokinetic profiles of the selected compounds. Drug-likeness of a compound is strongly stimulated with the aid of 2 properties of water solubility and hydrophobicity (predicted using the ALI log S, log S, and ESOL methods)24. Figure3 shows the drug-like properties (molecular weight etc.) obtained from ADME/Tox tool. Orally active drugs must follow Lipinski rules five.

Table 3.

Drug likeness properties of the nominated compounds.

S. NoCompounds IDMWTPSAHB donorHB accepNo. Rot bondLogPLogS(ESOL)LogS(Ali)Violation
1Ref46492.353794−5.11−5.710
25459840478144.114751.94−2.77−3.050
3SANC0034743494.452703.46−5.49−5.960
4ZINC470092073101023481.06−2.16−2.600
5sa32298742543.02−3.72−4.230
6SANC003524361093713.51−5.34−5.950
7SANC0051266415147103.67−6.00−7.200
844265937276.360763.10−4.40−4.790
91088534048115741041.64−4.14−4.840
10ZINC4497851129682.072321.19−2.28−1.710
11ZINC235435393111362690.89−2.00−3.180
12SANC008563541003622.62−4.27−4.800
13SANC010834361884106−0.22−2.54−3.680
146-KA6-TP32183.812552.47−3.86−4.380
15ZINC4823773030883.160731.89−2.94−2.680
161673665543653.351351.19−3.93−3.040
17995738438366.460732.57−4.06−3.790
18ZINC20540819258851551.22−2.22−2.210
197-KA7-TP30954.351343.40−4.54−4.810
208.3-Formylcouma3741324741.96−3.81−4.730
21BB-III-18369872433.21−5.07−5.720
22Sa5325682463.193.27−4.230
23Sa10284742442.542.60−3.520
24Sa12344922762.402.45−3.700
25Sa1531483.642652.442.50−3.590

Open in a new tab

Figure 3.

Open in a new tab

Molecular dynamic simulation

To investigate the structural stability of the docked complexes, 100ns molecular dynamics simulations were executed to determine the stability of the top lead compounds in complex with receptor. MD simulations were performed by using Amber program software. Molecular dynamics (MD) simulations of the top five (05) hits with good interactions, binding energies, best docking scores, and have satisfied ADMET properties compared to the reference molecule was carried out. While, stability of the complexes ZINC47009207 (green) compound, Sa32 complex (Blue), and SANC00347 (sky-blue) were evaluate based on root mean square deviation (RMSD), root mean square fluctuation (RMSF), radius of gyration (Rg) and DCCM. Moreover, to recognize the structural variability in protein over a time of 100ns, the simulated complex with existence of lead compounds at particular trajectory frames that is frame 1, 20, 40, 60, 80, and 100 were analyzed and are shown in Figs.4, 5, 6 and 7.

Figure 4.

Open in a new tab

Figure 5.

Identification of novel STAT3 inhibitors for liver fibrosis, using pharmacophore-based virtual screening, molecular docking, and biomolecular dynamics simulations (6)

Open in a new tab

Figure 6.

Identification of novel STAT3 inhibitors for liver fibrosis, using pharmacophore-based virtual screening, molecular docking, and biomolecular dynamics simulations (7)

Open in a new tab

Figure 7.

Open in a new tab

RMSD analysis

The predicted complexes stability was further checkedby MD simulations. After 100ns MD simulations, the RMSD of the backbone of STAT3 and the ligands at 300K was plotted against time (ns). This was shown in Fig.4, the RMSDs of Sa32 and ZINC47009207 were found to be comparatively firm and stable at around 2.7 and 3Å, respectively. Sa32 wasstable from the beginning in replica 1 then there were some fluctuations from 40 to 60ns. The Sa32 in replica2 revealed minor deviations from 5–10 to 70–85ns but overall the complex RMSD was stable during the 100ns MD run.ZINC47009207 in replica 1has been stable since beginning with some unstable behavior from the 60ns to 70ns as compared to Ref complex. The ZINC47009207in replica2 revealed an unstable behavior during 21–70ns but then the system gain stability and remained stable till the end of MD run. The average RMSD value of the reference compound was found as 3.1Å in replica 1. The RMSD of Ref compound in replica2 was stable during the first 25ns and then the RMSD increases from 26ns to 50ns after that the system gain stability and remained stable till the 100ns. The compound545984 in replica 1around 20ns show fluctuation then get stable and afterward from 50ns to 65 unstable and from 85 to 90ns show fluctuation. The RMSD value of the 545984 was found as 2.8Å during simulation. The RMSD of compound 545984 in replica2 indicate stable during the first 50ns after that a slight increase in RMSD was observed and fluctuations were observed during 51–82ns after that, the system gain stability and remained stable. The compound SANC00347 shows stability from start to 20ns then gradually shows fluctuation from 21 to 80nsin replica 1. Then from 85 to 100ns they show stability. The RMSD value of SANC00347 was found as 4.5Å during MD simulation. The RMSD of SANC00347 in replica2 indicates a similar pattern to that of the replica1. The RMSD was stable during the first 20ns but major deviations were observed from 20 to70ns and then get stability and remained stable during the last 30ns. The RMSD plots of replica 2 are shown in Fig. S4 (Supplementary). Except the compound SANC00347 the RMSD of all other compounds indicated high stability during the entire 100ns MD simulation.

RMSF analysis

The RMSF is a significant factor that produces data regarding the structural flexibility of each and every residue i-e remains in the system. In replica 1 the RMSF values of the Ref, 5459840, ZINC47009207, Sa32 and SANC00347 complexes were calculated (Fig.5). In RMSF we found that fluctuation was found in the region including 18–21(THR516, LYS517, ARG518, GLY519) 39–41(VAL537, ASN538, TYR539) 43–51(SER540, GLY541, CYS542, GLN543, ILE544, THR545, TRP546, ALA547, LY548, PHE549, CYS550, LYS551) 122–127 (THR620, PHE621, THR622, TRP623, VAL624, GLU625) 160–164(LYS658, ILE659, MET660, ASP661, ALA662). The fluctuating residues were not present in the binding site and the active site residues showed stable behavior which indicates strong binding of compounds with the receptor. The RMSF value of SANC00347 was very high up to 12Å while the RMSF value of all other complexes was less than 4Å. RMSD analysis was in agreement with the RMSF analysis which indicated that among all the complexes SANC00347 was unstable as compared to the control and other complexes. In replica 2 a similar pattern of RMSF was observed as the replica 1. In all the systems the RMSF was found less than 4Å except the SANC00347. The RMSF plot of SANC00347 was observed as 13Å. The residues such as 25–40 and 125–135 in the complexes 5459840 and the SANC00347 revealed high fluctuations during the MD simulation. The residues of Sa32 and ZINC47009207 revealed minor fluctuations during MD simulation. The RMSF plots of all the systems in replica 2 are shown in Fig. S5 (Supplementary).

Radius of gyration (Rg)

STAT3 protein compactness with their complexes was evaluated by utilizing the radius of gyration (Rg). The analysis of Rg confirmed that all complexes of STAT3 i.e., Sa32 and ZINC47009207 while Fig.6 showed a divergent pattern of compactness. Extensively, the Rg of the Sa32 and ZINC47009207 complexes show stability throughout MD simulation and demonstrate a vigorous compact conformation. The Rg value of SA32 (replica 1) was found as 16.9–17.3Å. The Rg value of ZINC47009207 (replica 1) was found as 17.2–17.4Å. The Rg value of SANC00347 (replica 1) was 17.2–17.5Å. The Rg value of 5459840 (replica 1), and Ref (replica 1) were found to be 17.1–17.3Å and 17.2–17.41Å respectively. Whereas, other 5459840 and SANC00347 complexes were indicated to be less compact over time as compared to Ref. The dynamic behavior of protein–ligand complexes reveal that binding alters stability and remaining flexibility, thereby inducing therapeutic ability. In replica 2 the Rg of the SANC00347 was observed to be 17.1–17.5Å. The Rg of the complex 5459840 was observed to be 17.0–17.4Å. The Rg value of ZINC47009207, Sa32 and Ref was found to be 16.9–17.1, 16.6–17.1Å and 16.4–17.1Å respectively. The Rg plots of replica 2 of all the complexes are shown in Fig. S6 (Supplementary). As compared to replica 2 all the complexes were more compact and revealed high stability in replica 1.

DCCM analysis

The aim of the DCCM analysis is to study the changes in the binding pocket conformation of STAT3 protein when interacting with various ligands. This was accomplished through a post Molecular Docking analysis of the backbone of Ca during MD simulations lasting 100ns, to observe any fluctuations and correlated motions. To determine binding affinity of selected or targeted ligands in the active site of the correlation and anti-correlation motions of STAT3. The dynamics matrix of cross-correlation plots of 100ns MD simulations was plotted for every complex. The analysis of DCCM was to evaluate correlation between residuals and this was performed to elucidate the motion correlation of the residuals. Movements of amino acids are in parallel direction, indicating strongly correlated movements called positive correlation. On the other side, if the amino acids movements are anti-parallel then this shows movement of anti-correlation. Negative & positive correlations were calculated by DCCM analysis and represented as the cyan to dark blue region and the red to green region as shown in Fig.7A–E

Free energy binding

The binding energy of ligand-STAT3 complexes was calculated using MMPBSA methods. However, van der Waals and electrostatic interactions and non-polar solvation energy negatively contribute to the total interaction energy while only polar solvation energy positively contributes to total free binding energy. In terms of negative contribution, van der Waals interaction gives much larger contribution than electrostatic interactions for all the cases. The non-polar free energy contributes relatively less as compared to the total binding energy. This indicates that non-polar solvation energy, van der Waals and electrostatic interaction together contribute to the STAT3 polyphenol complex stability25. MMPBSA was calculated for the last 10ns from the total of 100ns MD trajectory. The binding free energies ΔG Bind of complexes Ref, 5459840, ZINC47009207, Sa32 and SANC00347 calculated during the 100ns simulation were found to be −9.8027kcal/mole, −11.5849kcal/mole, −10.6547kcal/mole and −24.0920kcal/mole respectively. The details of MMPBSA calculation of the complexes are summarized in Table 4. In Table 4 shown the calculated values of the molecular mechanic’s VDW, EFL, EPB, ENPOLAR and EDIPER under ΔPBSA between the protein and ligand during the simulation were evaluated from the MM-PBSA calculation. The van der Waals VDW interactions of ligands −29.1030kcal/mole (Reference), −29.6634kcal/mole (5459840), −43.4243kcal/mole (ZINC47009207), −34.0802kcal/mole (Sa32) and −35.0421kcal/mole (SANC00347), all the complex systems showed the negative values of VDWso all complexes illustrated sufficient hydrophobic interactions. The results of the free energy calculation showed that the complex having best binding energy, also best RMSD and best binding affinity toward receptor. As compared to Ref compound binding energy of SANC00347 was not good enough.

Table 4.

MMPBSA free binding energy of selected compounds.

S. NoComplexVDWEELEPBENPOLAREDIPERΔPBSA
1Reference−29.1030−1.246015.3688−17.921234.7609−9.8027
25459840−29.663439.2518−11.9422−18.772332.7110−11.5849
3ZINC47009207−43.4243−2.519416.7161−23.748443.3213−10.6547
4Sa32−34.0802−140.0321132.3842−24.822042.4582−24.0920
5SANC00347−35.0421−3.373115.7086−21.110337.0142−7.8595

Open in a new tab

Discussion

Liver fibrosis is an initial stage of cirrhosis and results from abundant chronic liver disorders and a wide process that presents with numerous depositions of extracellular matrix proteins that change the liver architecture and functionality resulting in cirrhosis and failed organ26. Liver disorders are responsible for approximately 2million fatalities each year worldwide, with 1million of these deaths attributable to cirrhosis severity. Liver fibrosis is now at the top of list among ten causes of death globally27. Therefore, some STAT3 inhibitors exist such as Sorafenib inhibits STAT3 phosphorylation in human tumors, HCC, and liver fibrosis etc. However, STAT3 played a beneficial role against liver fibrosis due to the proliferative function of STAT3.

STX-0119 is a STAT3 inhibitor of dimerization. STX-0119 mediates the development of CCl4 and thioacetamide- induced liver fibrosis by diminishing the activation of HSCs28.

There are some STAT3 SH2 domain inhibitors that are Shikonin a natural naphthoquinone derivative that has been recognized as STAT3 inhibitory potency. Napabucasin has been also identified as a STAT3 inhibitor which inhibited STAT3 by binding to its SH2 domain. Nifuroxazide and ODZ10117 have a greater inhibition on STAT3 activation via directly blocking the SH2 domain of STAT3 and blocking the transcription activity29.

This study has shown few STAT3 inhibitors that help to overcome STAT3 activity as well as anticancer effects. According to Tomohiro, Von Manstein, and Groner there are also some disadvantages, such as dose-limiting side effects, multidrug resistance, and extreme side effects. The result is a timetable for the development of essentially new STAT3 inhibitors with improved viability30.

The overall analysis of novel STAT3 inhibitors from different Databases such as In-house (Sa32), Zinc(ZINC47009207), Phytochemicals(5459840) and South African (SANC00347) were carried out. The docking interaction of compounds, RMSD, RMSF, DCCM, Rg, and MMPBSA based analysis are discussed below.

From previous study of docking results shows that N4 (inhibitor) which are STAT3 inhibitor interacted with STAT3-SH2 at multiple amino acid residues some are following LYS591, SER636, VAL637, AND GLU63831 STAT3 inhibitor i-e S31-201with compounds which improved bioavailability and anti-cancer properties. Their mechanism of action depends upon interaction with STAT3-SH2 and phosphor tyrosine 705 site. The paining with active dimmer from tyrosine peptide is represented together with inhibitor S31-201 from H-bind and residues including LYS591, SER611, SER613 and ARG60932.

In current study Sorafenib dock as a reference drug of STAT3 and Sorafenib interact with STAT3-SH2 with some common residues same as previous study that is shown in Fig.2 they also show structures similarity and binding interaction. And South African compounds SANC00347 show some similarity with S31-201 inhibitor with STAT3-SH2 interactions as shown in Fig.2 S31-201 and SANC00347 show interaction on same active site residue.

RMSD provides detailed information about protein and their structure predictions. In previous studies RMSD of STAT3 inhibitor (−)-Epigallocatechingallate-complex, Saikosaponin D complex, Kaempferol-3-O-rutinoside complexes and Ginsenoside RK1 and Picroside. On the other hand, the bound state of the protein can stabilize between 50 and 100ns. For the complex the RMSD was found to exceed 1.5nm, with multiple variations for Ginsenoside RK1 and a decrease in stability for Picroside I within 75ns. Both compounds showed lower stability and higher RMSD due to multiple amplitudes. A number of STAT3 complexes are observed by evaluating the active behavior of backbone residues and atomic position variations24. In the current study RMSD analysis of novel STAT3 inhibitors such as Sa32, ZINC47009207, 545984 and SANC00347 complexes show stability but also fluctuate at some positions as compared to Reference. But Sa32 and ZINC47009207 best RMSD and binding affinity shows in Fig.4.

The previous RMSF result provides valuable information about the dynamic behavior of the residues. A reduced flexibility of the complexes will occur due to a small change in atomic positional Hicks24.

In the current study, RMSF is the fundamental factor that provides data on the structural compliance of each residue in the system. The values of the following that is Ref, 5459840, ZINC47009207, Sa32 and SANC00347 complexes were calculated Fig.5 in RMSF we found that fluctuation was found in the region including 18–21, 39–41, 43–51, 122–127 and 160–164 respectively. Some hikes present which are not active site residues.

The current STAT3 inhibitors which are analyzed by Rg analysis and as an important residue amplitudes and backbone deviation, more comprehensive check, general compression in all complexes were needed. The analysis of Rg shows that all complexes of STAT3 i.e., Sa32 and ZINC47009207, 5459840 and SANC00347 showed a different pattern of compared to Ref, as shown in Fig.6. The dynamic behavior of protein–ligand complexes reveal that binding alters stability and residual flexibility, there by inducing therapeutic clouds.

To investigate the structural changes in the binding pockets of STAT3 protein, to determine the interacting effects of selected compounds in the active site of STAT3 on positively correlated and negatively correlated movements. Hence, for each complex system, the dynamic cross-plots of the correlation matrix are plotted for 100ns MD simulation, which analyzed the regions of negative and positive correlation, such as cyan to dark blue region and red to green, as shown in Fig.7A–E.

MMPBSA Binding Free Energy (Kcal/mol) is calculated for nominated complexes for estimated protein–ligand complexesof 100ns MD trajectory. The overall results of complexes were compared and on the basis of MMPBSA analysisSa32 hadthe best binding energy which is shown in Table 4. Sa32 is also best in all analyses like RMSD, and Rg indicating the strong binding affinity of the Sa32 compound toward the receptor.

Conclusion

In the current study, novel STAT3 inhibitors were identified with a comprehensive screening method. The Ph4-based virtual screening, docking, and MD simulations were carried out. Additionally, 25 complexes were selected based onbest docking scores which were further screened for the drug-likeness. ADMET drug similarity was utilized to filter the properties of pharmacokinetic of the designed ligand–protein complexes. Lastly, 5 molecules were selected having excellent properties as well as stable binding modes. In addition, Sa32 has a valuable and novel STAT3 inhibitor, which will further serve as a reference for the development of new effective STAT3 inhibitors. Sa32 compound is overall showing the best response toward STAT3 inhibition which is confirmed through RMSD, RMSF, RG and DCCM analysis. It is further recommended to perform Invitro study on the final predicted hits to evaluate the inhibitory potential of these compounds.

Supplementary Information

Supplementary Information. (3MB, docx)

Acknowledgements

The authors would like to thank Guangdong Basic and Applied Basic Research Foundation (2021B1515140026). Further, MAA would like to extend his appreciation to Prince Sattam bin Abdulaziz University for providing the necessary technical and computational tools under project number (PSAU/2023/R/1444).

Author contributions

Conceptualization, A.W. and X.U. methodology, H.R., A.W. and X.U. software, M.A.H., A.Z., M.A.A. and A.M. validation, B.S.A., A.M. and J.H. formal analysis, A.Z., J.H., M.A.H., H.A.A. and X.U. investigation, J.H., M.A.H., M.A. and B.S.A. resources, A.Z., M.A.A. and B.S.A. data curation, J.H., M.A.H., B.S.A., M.A.A. and A.Z., writing—original draft preparation, H.R. writing—review and editing, J.H., A.W., H.A.A. and X.U. visualization, M.A.H., A.M., H.A.A. and A.Z. supervision, A.W. and X.U. project administration, A.W. and Z.U. funding acquisition, X.U. All authors have read and agreed to the published version of the manuscript.

Data availability

All the data and its links are available in the manuscript.

Competing interests

The authors declare no competing interests.

Footnotes

Publisher's note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Contributor Information

Abdul Wadood, Email: awadood@awkum.edu.pk.

Xiaoyun Huang, Email: hxydg21@163.com.

Supplementary Information

The online version contains supplementary material available at 10.1038/s41598-023-46193-x.

References

  • 1.Poli G, et al. Identification of a new STAT3 dimerization inhibitor through a pharmacophore-based virtual screening approach. J. Enzyme Inhib. Med. Chem. 2016;31(6):1011–1017. doi: 10.3109/14756366.2015.1079184. [DOI] [PubMed] [Google Scholar]
  • 2.Tolomeo M, Cascio A. The multifaced role of STAT3 in cancer and its implication for anticancer therapy. Int. J. Mol. Sci. 2021;22(2):603. doi: 10.3390/ijms22020603. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 3.Qin J-J, et al. STAT3 as a potential therapeutic target in triple negative breast cancer: A systematic review. J. Exp. Clin. Cancer Res. 2019;38:1–16. doi: 10.1186/s13046-019-1206-z. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 4.Lakshmanan K, et al. Discovery of potential inhibitors for stat3: Ligand based 3D pharmacophore, virtual screening, molecular docking, dynamic studies and in vitro evaluation. J. Biomol. Struct. Dyn. 2022;40(21):11320–11338. doi: 10.1080/07391102.2021.1957717. [DOI] [PubMed] [Google Scholar]
  • 5.Wang X, et al. STAT3 inhibition, a novel approach to enhancing targeted therapy in human cancers. Int. J. Oncol. 2012;41(4):1181–1191. doi: 10.3892/ijo.2012.1568. [DOI] [PubMed] [Google Scholar]
  • 6.Wang Z, Long J, Zhang H. The STAT3 inhibitor S3I–201 suppresses fibrogenesis and angiogenesis in liver fibrosis. Lab. Invest. 2018;98(12):1600–1613. doi: 10.1038/s41374-018-0127-3. [DOI] [PubMed] [Google Scholar]
  • 7.Tang M, et al. Therapeutic targeting of STAT3 with small interference RNAs and antisense oligonucleotides embedded exosomes in liver fibrosis. FASEB J. 2021;35(5):e21557. doi: 10.1096/fj.202002777RR. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 8.Alkreathy HM, Esmat A. Lycorine ameliorates thioacetamide-induced hepatic fibrosis in rats: Emphasis on antioxidant, anti-inflammatory, and STAT3 inhibition effects. Pharmaceuticals. 2022;15(3):369. doi: 10.3390/ph15030369. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 9.Deng Y-R, et al. STAT3-mediated attenuation of CCl4-induced mouse liver fibrosis by the protein kinase inhibitor sorafenib. J. Autoimmun. 2013;46:25–34. doi: 10.1016/j.jaut.2013.07.008. [DOI] [PubMed] [Google Scholar]
  • 10.Hung M-H, et al. Downregulation of signal transducer and activator of transcription 3 by sorafenib: A novel mechanism for hepatocellular carcinoma therapy. World J. Gastroenterol. WJG. 2014;20(41):15269. doi: 10.3748/wjg.v20.i41.15269. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 11.Wadood A, et al. Machine learning-based virtual screening for STAT3 anticancer drug target. Curr. Pharm. Des. 2022;28(36):3023–3032. doi: 10.2174/1381612828666220728120523. [DOI] [PubMed] [Google Scholar]
  • 12.Irwin JJ, et al. ZINC: A free tool to discover chemistry for biology. J. Chem. Inf. Model. 2012;52(7):1757–1768. doi: 10.1021/ci3001277. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 13.Muegge I. Selection criteria for drug-like compounds. Med. Res. Rev. 2003;23(3):302–321. doi: 10.1002/med.10041. [DOI] [PubMed] [Google Scholar]
  • 14.Wadood A, et al. In silico drug designing for ala438 deleted ribosomal protein S1 (RpsA) on the basis of the active compound Zrl 15. ACS Omega. 2021;7(1):397–408. doi: 10.1021/acsomega.1c04764. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 15.Becker S, Groner B, Müller CW. Three-dimensional structure of the Stat3β homodimer bound to DNA. Nature. 1998;394(6689):145–151. doi: 10.1038/28101. [DOI] [PubMed] [Google Scholar]
  • 16.Park IH, Li C. Characterization of molecular recognition of STAT3 SH2 domain inhibitors through molecular simulation. J. Mol. Recogn. 2011;24(2):254–265. doi: 10.1002/jmr.1047. [DOI] [PubMed] [Google Scholar]
  • 17.Alotaibi BS, et al. New drug target identification in Vibrio vulnificus by subtractive genome analysis and their inhibitors through molecular docking and molecular dynamics simulations. Heliyon. 2023;9:e17650. doi: 10.1016/j.heliyon.2023.e17650. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 18.Ajmal A, et al. Computer-assisted drug repurposing for thymidylate kinase drug target in monkeypox virus. Front. Cell. Infect. Microbiol. 2023;13:618. doi: 10.3389/fcimb.2023.1159389. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 19.Daina A, Michielin O, Zoete V. SwissADME: A free web tool to evaluate pharmacokinetics, drug-likeness and medicinal chemistry friendliness of small molecules. Sci. Rep. 2017;7(1):42717. doi: 10.1038/srep42717. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 20.Ebrahimi KS, et al. In silico investigation on the inhibitory effect of fungal secondary metabolites on RNA dependent RNA polymerase of SARS-CoV-II: A docking and molecular dynamic simulation study. Comput. Biol. Med. 2021;135:104613. doi: 10.1016/j.compbiomed.2021.104613. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 21.Opo FADM, et al. Structure based pharmacophore modeling, virtual screening, molecular docking and ADMET approaches for identification of natural anti-cancer agents targeting XIAP protein. Sci. Rep. 2021;11(1):4049. doi: 10.1038/s41598-021-83626-x. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 22.Tabti K, et al. Identification of a potential thiazole inhibitor against biofilms by 3D QSAR, molecular docking, DFT analysis, MM-PBSA binding energy calculations, and molecular dynamics simulation. Phys. Chem. Res. 2023;11(2):369–389. [Google Scholar]
  • 23.Zhou Y, Di B, Niu M-M. Structure-based pharmacophore design and virtual screening for novel tubulin inhibitors with potential anticancer activity. Molecules. 2019;24(17):3181. doi: 10.3390/molecules24173181. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 24.Manoharan S, et al. Screening of potent STAT3-SH2 domain inhibitors from JAK/STAT compound library through molecular dynamics simulation. Mol. Divers. 2022;27:1297–1308. doi: 10.1007/s11030-022-10490-w. [DOI] [PubMed] [Google Scholar]
  • 25.Verma S, et al. Hydrophobic interactions are a key to MDM2 inhibition by polyphenols as revealed by molecular dynamics simulations and MM/PBSA free energy calculations. PLoS ONE. 2016;11(2):e0149014. doi: 10.1371/journal.pone.0149014. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 26.Erdem Koc G, Gokcimen A, Sahin F. The effect of boric acid and sodium pentaborate pentahydrate-treated foreskin derived mesenchymal stem cells on liver fibrosis. Biol. Trace Elem. Res. 2023;201:4834–4849. doi: 10.1007/s12011-023-03565-8. [DOI] [PubMed] [Google Scholar]
  • 27.Zhang L, Chan C. Isolation and enrichment of rat mesenchymal stem cells (MSCs) and separation of single-colony derived MSCs. JoVE. 2010;37:e1852. doi: 10.3791/1852. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 28.Zhao J, Qi Y-F, Yu Y-R. STAT3: A key regulator in liver fibrosis. Ann. Hepatol. 2021;21:100224. doi: 10.1016/j.aohep.2020.06.010. [DOI] [PubMed] [Google Scholar]
  • 29.Dong J, et al. Recent update on development of small-molecule STAT3 inhibitors for cancer therapy: From phosphorylation inhibition to protein degradation. J. Med. Chem. 2021;64(13):8884–8915. doi: 10.1021/acs.jmedchem.1c00629. [DOI] [PubMed] [Google Scholar]
  • 30.Chiba T. STAT3 inhibitors for cancer therapy-the rationale and remained problems. EC cancer. 2016;1(S1):S1–S8. [Google Scholar]
  • 31.Chen H, et al. Selectively targeting STAT3 using a small molecule inhibitor is a potential therapeutic strategy for pancreatic cancer. Clin. Cancer Res. 2023;29(4):815–830. doi: 10.1158/1078-0432.CCR-22-0997. [DOI] [PubMed] [Google Scholar]
  • 32.Fagard R, et al. STAT3 inhibitors for cancer therapy: Have all roads been explored? Jak-Stat. 2013;2(1):e22882. doi: 10.4161/jkst.22882. [DOI] [PMC free article] [PubMed] [Google Scholar]

Associated Data

This section collects any data citations, data availability statements, or supplementary materials included in this article.

Supplementary Materials

Supplementary Information. (3MB, docx)

Data Availability Statement

All the data and its links are available in the manuscript.

Identification of novel STAT3 inhibitors for liver fibrosis, using pharmacophore-based virtual screening, molecular docking, and biomolecular dynamics simulations (2024)
Top Articles
Latest Posts
Recommended Articles
Article information

Author: Tish Haag

Last Updated:

Views: 6264

Rating: 4.7 / 5 (47 voted)

Reviews: 94% of readers found this page helpful

Author information

Name: Tish Haag

Birthday: 1999-11-18

Address: 30256 Tara Expressway, Kutchburgh, VT 92892-0078

Phone: +4215847628708

Job: Internal Consulting Engineer

Hobby: Roller skating, Roller skating, Kayaking, Flying, Graffiti, Ghost hunting, scrapbook

Introduction: My name is Tish Haag, I am a excited, delightful, curious, beautiful, agreeable, enchanting, fancy person who loves writing and wants to share my knowledge and understanding with you.