In silico interaction analysis of selected natural compounds with bacteriophage-encoded hyaluronate lyase from Streptococcus pyogenes

TL;DR

This study uses computer modeling to find natural compounds that can inhibit a key virulence factor in Streptococcus pyogenes, offering a potential alternative to antibiotics.

Contribution

The paper introduces a novel in silico approach to identify natural HylP2 inhibitors with favorable pharmacokinetic properties.

Findings

Violacein, sulfangolid C, chlorotonil A, xiamycin, and kulkenon showed strong binding to HylP2.

MD simulations confirmed stable interactions of violacein, xiamycin, and kulkenon with key catalytic residues.

The compounds comply with Lipinski’s Rule of Five and show favorable ADMET profiles.

Abstract

The rising antibiotic resistance of Streptococcus pyogenes necessitates alternative anti-virulence strategies. Bacteriophage-encoded hyaluronate lyase (HylP2), a key virulence factor that promotes bacterial dissemination by degrading host extracellular matrix components, represents an attractive therapeutic target.

In this study, an integrated in silico approach was employed to identify potential HylP2 inhibitors from a library of 118 bioactive natural compounds. Following protocol validation through redocking of ascorbic acid (RMSD = 1.897 Å), virtual screening, ADMET prediction, molecular dynamics (MD) simulations, and per-residue energy decomposition analyses were performed.

Violacein (−7.7 kcal/mol), sulfangolid C (−7.427 kcal/mol), chlorotonil A (−7.4 kcal/mol), xiamycin (−7.3 kcal/mol), and kulkenon (−7.1 kcal/ mol) were identified as the most potent binders. ADMET analysis confirmed that these leads possess favorable pharmacokinetic properties and compliance with Lipinski’s Rule of Five. Subsequent 100-ns molecular dynamics (MD) simulations and per-residue energy decomposition revealed that violacein, xiamycin, and kulkenon formed stable, compact complexes by “trapping” catalytic residues Arg279 and Tyr264.

These findings suggest that these natural product scaffolds are promising anti-virulence leads that may limit S. pyogenes tissue invasion while minimizing selective pressure for resistance development.

Figures

The 3D structure of bacteriophage-encoded hyaluronate lyase derived from Streptococcus pyogenes (viewed using Biovia Discovery Studio Visualizer).

3-D (Left) and 2-D (Right) representations of the binding pose of ascorbic acid (A), bikaverin (B), violacein (C), xiamycin (D), chlorotinil (E), phenoxan (F), and kulkenon (G) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. 3-D (Left) and 2-D (Right) representations of the binding pose of ascorbic acid (A), bikaverin (B), violacein (C), xiamycin (D), chlorotinil (E), phenoxan (F), and kulkenon (G) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Radius of gyration (Rg) profiles from molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. Radius of gyration (Rg) profiles from molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Root-mean-square deviation (RMSD) profiles obtained from the molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. Root-mean-square deviation (RMSD) profiles obtained from the molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Root-mean-square fluctuation (RMSF) profiles obtained from MD simulations of violacein (A), bikaverin (B), chlorotonil-A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase derived from S. pyogenes.

3D free-energy landscapes (FEL) profiles obtained from MD simulations of violacein (A), bikaverin (B), chlorotonil-A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. 3D free-energy landscapes (FEL) profiles obtained from MD simulations of violacein (A), bikaverin (B), chlorotonil-A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Per-residue energy decomposition profiles from MD simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. Per-residue energy decomposition profiles from MD simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Tables

Table 1. Binding affinity of compounds derived from bacteria and cyanobacteria against bacteriophage-encoded hyaluronate lyase from Streptococcus pyogenes.

Compounds	Binding affinity (kcal/mol)	Compounds	Binding energy (kcal/mol)
Coformycin	−4.9	Abietane	−6.1
Capolactin B	−4.8	Ambigol A	−5.5
Formycin	−4.9	Ambigol B	−5.4
Resveratrol	−5.5	Malyngolide	−4.3
Phenalamide A2	−5.7	Lyngbyoic acid	−4.2
Scutellarein	−5.9	pitinoic acid A	−3.9
Siastatin B	−4.9	Hapalindole	−6.4
Xiamycin	−7.3	anaephene B	−4.5
Violacein	−7.7	anaephene A	−4.2
Sarkomycin	−4.5	Anaephene C	−5
Anthranoside C	−6.0	Cylindrofridin A	−4
Apigenin	−6.0	Capolactin A	−4.5
Baicalein	−6.2	Triton-x100#	−4.7
Pinocembrin	−6.1	Ascorbic acid#	−3.6
Rosmarinic acid	−5.9	Sodium dodecyl sulfate#	−4.0

Table 2. Binding affinity of compounds from Myxobacteria against bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Compounds	Binding affinity (kcal/mol)	Compounds	Binding energy (kcal/mol)
Epothilon D	−6.4	Ajudazol	−5
Kulkenon	−7.1	Althiomycin	−5
Myxochelin A	−4.7	Angiolactone	−5.5
Noricumazole C	−5.3	Aurachin E	−5.4
Phenoxan	−5.9	Carolacton	−5.3
Ratjadon	−5.2	Chlorotonil-A	−7.4
Thiangazole	−6.3	Corallorazine	−5
Sulfangolid C	−7.4	Cystobactamid	−6
Stipiamide	−4.7	Disciformycin A	−5.4
Soraphen F	−6.0	Enhygrolide B	−5.6
Labindole A	−4.4	Hyapyrone B	−6.7
Labindole B	−4.7	Hyalachelin	−4.8
3-chloro-9H-carbazole	−4.9	Hyaladione	−3.8
4-hydroxymethyl-quinoline	−6.4	Indiacen A	−5.2
Letermovir	−6.3	Indothiazinone	−4.8
Spirangien B	−6.0	Methyl indole-3-carboxylate	−4.4
Noricumazole A	−5.4	Melithiazol A	−4.4
Nannozinone B	−5.1	Nannozinone A	−5.5

Table 3. Binding affinity of fungi-derived compounds against bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Compounds	Binding affinity (kcal/mol)	Compounds	Binding energy (kcal/mol)
Aureonitol	−4.0	Cordycepin	−4.5
Ganoderic acid	−5.2	Curvularin	−5.3
Velutin	−4.5	Citreorosein	−5.1
Deoxyfunicone	−5.5	Griseofulvin	−4.6
Fuscinarin	−5.8	Brefeldin A	−5
Hinnuliquinone	−6.5	Bikaverin	−5.8
Griseoxanthone C	−6.1	Xantocillin	−5.5
Kaempferol	−5.9	Aphidicolin	−5.2
10-methoxydihydrofuscin	−5.6	Cajanol	−4.7
Quercetin	−5.2	Lysergic acid	−5.3
Cytosporin A	−5.4	Ergothioneine	−3.9
Mellein	−5.2	Kojic acid	−3.5
6-Hydroxymellein	−4.5	Glutathione	−3.9
Sabinene	−3.9	Camptothecin	−6.2
cyclo(L-Phe-L-Pro)	−3.9	Tyrosol	−3.5
6-Pentyl-2H-pyran-2-one (6-Amyl-2-pyrone)	−3.9	Xanthone	−5.1
Emericellin	−5.6	2,4-Diacetylphloroglucinol	−4.2
cyclo(L-Pro-L-Val)	−4.1	Ergosterol	−6.8
Patulin	−3.7	Ergosterol peroxide	−6.3
Emodin	−5.3

Table 4. Binding affinities of antibacterial compounds against bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Compounds	Binding affinity (kcal/mol)
2-(2-chlorophenyl)-6-methoxy-4H-chromen-4-one	−5.9
1,3-Dicyclohexylurea	−4.5
(2Z)-2-(4-ethoxybenzylidene)-6-hydroxy-1-benzofuran-3(2H)-one	−4.9
6-[(3-methoxyphenyl)methoxy]-2-[(3,4,5-trimethoxyphenyl)methylidene]-1-benzofuran-3-one	−5.8
6-[2-(4-methoxyphenyl)-2-oxoethoxy]-2-[(2,4,5-trimethoxyphenyl)methylidene]-1-benzofuran-3-one	−5.5
7-hydroxy-3-(2-methoxyphenoxy)-2-methylchromen-4-one	−5.0
(2E)-3-(2,4-dimethoxyphenyl)-1-(2,5-dimethoxyphenyl)prop-2-en-1-one	−4.7
[2-[(2,3-Dimethoxyphenyl)methylidene]-3-oxo-1-benzofuran-6-yl] 2-methylpropanoate	−5.3
3-(2-methoxyphenyl)-4-oxo-4H-chromen-7-yl dimethylcarbamate	−5.7
[5-Hydroxy-3-(4-methoxyphenyl)-4-oxochromen-7-yl] 2,2-dimethylpropanoate	−5.9
5-hydroxy-4-oxo-2-phenyl-4H-chromen-7-yl dimethylcarbamate	−5.9
Methyl (1S,3S,4R)-3-hydroxy-4,7,7-trimethylbicyclo[2.2.1]heptane-1-carboxylate	−4.2
Methyl 2-{[(2Z)-2-[(3-methylphenyl)methylidene]-3-oxo-1-benzofuran-6-yl]oxy}acetate	−5.4

Full text

Introduction

1

Streptococcus pyogenes is a clinically significant human pathogen responsible for a wide spectrum of diseases, ranging from mild skin and throat infections to severe, life-threatening conditions such as sepsis, streptococcal toxic shock syndrome, scarlet fever, and necrotizing fasciitis (1). The Gram-positive, β-hemolytic bacterium, also known as group A streptococci; GAS, commonly colonizes the throat, vaginal mucosa, and anorectal region (2). The global prevalence of S. pyogenes infections varies depending on the clinical manifestations and geographic factors. Infection involves multiple factors, including intricate interactions between host immune defenses and bacterial virulence determinants that collectively affect disease pathophysiology (3).

Bacteriophages are the most abundant viruses in the human microbiome and play a crucial role in maintaining or disrupting microbial ecological balance. They influence human health both directly through bacterial predation, and indirectly, by modulating host immune responses and metabolic processes (4). Interactions between bacteriophages and their bacterial hosts are central to the evolution and dissemination of virulence traits. In S. pyogenes, numerous toxins and virulence factors are bacteriophage-encoded, and phage-derived sequences constitute an integral component of the GAS genome (5–7). Among these factors, bacteriophage-encoded hyaluronidases (hyaluronate lyases) have been identified in S. pyogenes and related species such as S. equi (8). Notably, S. pyogenes strains SF370.1 and 10,403 harbor the phage-encoded hyaluronate lyases HylP1 and HylP2, respectively, which have been implicated in tissue invasion and bacterial dissemination (9, 10). Hyaluronate lyases such as HylP2 facilitate bacterial invasion by degrading components of the host extracellular matrix. These enzymes primarily target hyaluronic acid (HA) and chondroitin sulfates, while exhibiting limited activity toward dermatan sulfates (11). Degradation of these glycosaminoglycans compromises host connective tissue integrity, facilitating bacterial spread and enhancing the dissemination of other phage-encoded toxins. The resulting breakdown products may also serve as nutrient sources that support bacterial proliferation at the sites of infection (12, 13). HylP2 has been biochemically characterized, and ascorbic acid has been identified as an inhibitor of its enzymatic activity (11). Additionally, phage-encoded hyaluronate lyases interact with the human transmembrane glycoprotein CD44, further contributing to host–pathogen interactions (12, 13).

.Medical management of S. pyogenes infections primarily relies on antibiotic therapy; however, the emergence of strains with reduced penicillin sensitivity underscores the urgent need for alternative strategies (14–16). Conventional antibiotics exert high selective pressure, driving resistance (17). In contrast, targeting HylP2 represents an anti-virulence strategy. By inhibiting the degradation of the host extracellular matrix, these agents may restrict bacterial dissemination and “disarm” the pathogen without exerting the bactericidal pressure that leads to resistance (18, 19). Consequently, HylP2 is a highly attractive molecular target for drug discovery.

Although efforts are underway to combat S. pyogenes infections using anti-virulence strategies, challenges remain in identifying effective and selective therapeutic agents. Numerous bioactive substances from diverse sources have been investigated for the treatment of infectious diseases. Microbial metabolites provide a valuable starting point for drug development because microorganisms can produce various natural compounds with considerable health benefits; these compounds currently constitute approximately half of all approved therapeutic agents (20). Natural products derived from myxobacteria, cyanobacteria, and fungi are the promising sources of anti-streptococcal agents owing to their extensive metabolic diversity and evolutionarily optimized chemical defense mechanisms. Myxobacteria have yielded hundreds of structurally diverse secondary metabolites, many with antibacterial activities, making them powerful reservoirs for novel antibiotics (21, 22). Cyanobacterial compounds and fungal secondary metabolites also exhibit broad-spectrum antimicrobial activities (23–25). Thus, these microorganisms represent valuable, underexplored resources for discovering compounds active against streptococci.

To identify potential inhibitors from these diverse natural sources, structure-based drug discovery has emerged as a critical approach. In silico workflows allow for the rapid characterization of drug–target complexes, providing mechanistic insights into ligand-binding behavior (26, 27). In this study, we employ an integrated pipeline utilizing molecular docking to predict binding orientations, followed by ADMET filtering to prioritize compounds with favorable pharmacokinetic profiles (28, 29). To overcome the limitations of static docking, molecular dynamics (MD) simulations are further employed to evaluate the stability and conformational flexibility of the complexes under physiological conditions (30). This integrated approach enhances the likelihood of identifying safe, stable, and biologically active lead compounds against S. pyogenes.

Materials and methods

2

Retrieval of hyaluronate lyase and preparation for molecular docking

2.1

The three-dimensional (3D) structure of hyaluronate lyase (HylP2; PDB: 2DP5, resolution: 3.55 Å) was obtained from the RCSB protein data bank1 in PDB format (31). This enzyme forms a single chain of 332 amino acids, with Gln261, Tyr264, and Arg279 as the active site residues. The enzyme was visualized using UCSF Chimera version 1.16 after removal of crystallographic water molecules (32). Thereafter, it was imported into the PyRx virtual screening tool,2 saved as an AutoDock macromolecule, and subsequently converted to PDBQT format (33, 34).

Collection of ligand 3D structures and preparation for molecular docking

2.2

In total, 118 compounds were retrieved from PubChem3 in SDF format (35). The molecules were selected based on documented antimicrobial or anti-streptococcal activities and their biosynthetic origins, which are known sources of bioactive secondary metabolites. The dataset included 27 compounds from bacteria and cyanobacteria, 36 from myxobacteria, 39 from fungi, and 13 established antibacterial agents (21, 36, 37). These organisms were targeted because they are well-recognized producers of structurally diverse natural products with reported activities against Gram-positive pathogens, including Streptococcus spp. ascorbic acid, triton X-100, and sodium dodecyl sulfate, which are the known inhibitors of hyaluronate lyase (11, 38), were included as controls to validate the docking protocol.

All compounds were imported into the PyRx virtual screening tool (See footnote 2) after successful download. The molecular structures were energy-minimized using the Universal Force Field (UFF) and subsequently converted to the PDBQT format for docking. AutoDock Vina settings were configured to ensure reproducibility. The exhaustiveness parameter was set to 8 to provide a balanced trade-off between computational cost and comprehensive exploration of the binding site (33, 34).

Molecular docking, interaction analysis, and visualization

2.3

Molecular docking was performed using the Vina wizard when the enzyme and ligands were successfully imported into the PyRx virtual screening tool (33). The grid dimensions were set in the Vina search space (Center X: 27.5716, Y: 18.1412, Z: 31.4922; Dimensions [Å] X: 33.5308, Y: 23.1604, Z: 23.5902) to contain the active site. This configuration allowed ligands to interact freely with the active-site residues to adopt a favorable binding conformation.

Redocking validation was performed using ascorbic acid, a co-crystallized inhibitor, to assess protocol reliability. The co-crystal ligand was extracted from the protein structure, energy-minimized using the UFF, and converted to the PDBQT format in PyRx. Docking was performed in AutoDock Vina (via PyRx v0.8) with an exhaustiveness value of 8, using a grid box sized to fully encompass the catalytic cleft. The top-scoring redocked pose yielded a root-mean-square deviation (RMSD) of 1.897 Å, which is within the acceptable threshold (<2.0 Å) for successful pose reproduction, and was centered at X = 15.085187, Y = 16.508875, Z = 32.766250. This confirmed that the docking workflow accurately replicated the native binding mode and was suitable for subsequent virtual screening analyses. The binding energies (kcal/mol) generated from the interaction between the enzyme and the ligands were recorded in a Microsoft Excel spreadsheet. The Biovia Discovery Studio Visualizer4 was used to visualize interactions between the hyaluronate lyase and the ligands after docking.

Computational prediction of ADMET parameters

2.4

The compounds were further subjected to in silico ADMET evaluation using ADMETlab2.0.5 These compounds were individually uploaded into the database in either SDF or simplified molecular-input line-entry system format. Their pharmacokinetic properties and drug-likeness were subsequently evaluated in accordance with Lipinski’s Rule of Five. Physicochemical properties and toxicity parameters were also assessed.

System preparation and force field selection for MD simulation

2.5

MD simulations were performed for 100 ns to evaluate the structural stability and binding persistence of six lead candidates: violacein, bikaverin, chlorotonil A, kulkenon, phenoxan, and xiamycin.

The 3D structure of HylP2 (PDB: 2DP5) was obtained, and preparation steps were performed, including the elimination of non-critical water molecules, heteroatoms, and metals (39). To ensure structural correctness, missing loops were modeled, and the protein was protonated to represent a physiological pH of 7.0. The energy minimization of the ligand (violacein, bikaverin, chlorotonil A, kulkenon, phenoxan, and xiamycin) geometries was achieved using the UFF (40). Ligand topology files and partial charges were produced using allied force field builders such as the Automated Topology Builder. To ensure consistent force field treatment within the complex, Amber ff14SB or Amber ff19SB force fields were applied to the protein parts, whereas the General Amber Force Field was used for the ligands (41, 42).

Solvation, ionization, and equilibration

2.6

For each complex, the orthorhombic or triclinic simulation box centered on the protein–ligand complex had a minimum distance of 10 Å (1.0 nm) from the box wall to the protein to prevent self-interactions during periodic boundary conditions. The systems were solvated using a standard three-point water model (TIP3P or SPC) to simulate a biological aqueous environment. Counterions (Na^+^ or Cl^−^) supplemented with 0.15 M NaCl were added in the case of uneven systems (43, 44). The systems underwent energy minimization for 5,000–10,000 steps to remove steric clashes, either by the steepest descent or conjugate gradient methods. Equilibration was performed in two stages: the NVT ensemble (constant number of particles, volume, and temperature), which allowed the system to reach the desired temperature of 300 K or 310.15 K; and the NPT ensemble (constant number of particles, pressure, and temperature), which stabilized the system’s density and pressure at 1 bar (45).

Production of MD simulations and trajectory analysis

2.7

The final MD simulations were performed for 100 ns using the NPT ensemble with a time step of 2 fs (46). All bonds involving hydrogen atoms, were constrained using the SHAKE or LINCS algorithms (47). The particle mesh Ewald method was used to handle electrostatic interactions with a standard 10-Å cutoff for non-bonded interactions. The coordinates and energy data were sampled every 10–100 ps to generate a complete dataset for trajectory analysis. Stability and conformational behavior were quantified using RMSD and root-mean-square fluctuation (RMSF), respectively, to evaluate structural stability and residue flexibility. The radius of gyration (Rg) was calculated to evaluate protein compactness. The 3D free-energy landscapes (FELs) and per-residue energy decomposition were used to locate the lowest-energy conformations and identify residues with the highest contribution to binding free energy (48).

Results

3

Molecular docking of the compounds (ligands) against the crystal structure of S. pyogenes-derived hyaluronate lyase (Figure 1) was performed. The binding energies are presented in Tables 1–4.

The 3D structure of bacteriophage-encoded hyaluronate lyase derived from Streptococcus pyogenes (viewed using Biovia Discovery Studio Visualizer).

Binding energy is a function of the binding affinity of the ligand toward the enzyme. A lower value indicates a stronger affinity. The molecular docking of the compounds derived from bacteria and cyanobacteria against the S. pyogenes-derived hyaluronate lyase showed binding affinities ranging from −3.9 to −7.7 kcal/mol. The highest binding affinities were observed for violacein (−7.7 kcal/mol) and xiamycin (−7.3 kcal/mol). Baicalein, pinocembrin, apigenin, and anthranoside C showed binding affinities of −6.2, −6.1, and −6.0 kcal/mol, respectively. The positive controls (ascorbic acid, −3.6 kcal/mol; Triton X-100, −4.7 kcal/mol; sodium dodecyl sulfate, −4.0 kcal/mol) exhibited weaker binding, indicating that the predicted affinity of several natural compounds surpasses that of conventional inhibitors (Table 1). The compounds derived from myxobacteria revealed binding affinities ranging from −3.8 to −7.4 kcal/mol. Sulfangolid C, kulkenon, chlorotonil A, enhygrolide B, epothilon D, 4-hydroxymethyl-quinoline, and thiangazole showed binding affinities between −6.3 and −7.4 kcal/mol. Myxochelin A, noricumazole C, phenoxan, aurachin E, soraphen F, 3-chloro-9H-carbazole, melithiazol A, and carolacton showed binding affinities between −4.4 and −6.0 kcal/mol (Table 2). The binding affinities of fungi-derived compounds ranged from −3.5 to −6.8 kcal/mol. The highest binding affinities were observed for ergosterol (−6.8 kcal/mol), hinnuliquinone (−6.5 kcal/mol), ergosterol peroxide (−6.3 kcal/mol), and griseoxanthone C (−6.1 kcal/mol). In contrast, sabinene, cyclo(L-Phe-L-Pro), 6-Pentyl-2H-pyran-2-one, ergothioneine, and glutathione exhibited similar binding affinities of −3.9 kcal/mol, whereas brefeldin A, bikaverin, xantocillin, aphidicolin, lysergic acid, xanthone, kaempferol, 10-methoxydihydrofuscin, quercetin, and cytosporin A had binding affinities between −5.0 and −5.9 kcal/mol (Table 3). The binding affinities of the antibacterial compounds ranged from −4.2 to −5.9 kcal/mol. Within this group, the highest binding affinity (−5.9 kcal/mol) was observed for 5-hydroxy-4-oxo-2-phenyl-4H-chromen-7-yl dimethylcarbamate and 2-(2-chlorophenyl)-6-methoxy-4H-chromen-4-one (Table 4).

Visualization of compound–enzyme interactions revealed 0, 1, 2, 3, 4, 5, and 6 hydrogen bonds (H-bonds) in 26, 37, 24, 13, 8, 6, and 2 compounds, respectively. The most notable amino acid residues were Val185, Ser211, Met213, Leu215, Glu219, Thr224, Leu225, Lys226, Ile227, His229, Asn231, Asp239, Ala243, Ala244, Leu245, Leu247, Leu249, Gln261, Gly262, Ile263, Tyr264, Ile265, Leu275, Leu276, Arg279, Asn280, Leu281, Ser282, Phe286, Val288, Phe294, Ala296, Lys297, Glu298, Thr299, Ser300, and Gln301 (Supplementary Table S1). Among these, Gln261, Tyr264, and Arg279 constituted the active-site residues and were involved in interactions with kulkenon, noricumazole C, phenoxan, soraphen F, labindole B, aureonitol, velutin, ascorbic acid, kojic acid, emodin, citreorosein, griseofulvin, brefeldin A, bikaverin, xantocillin, lysergic acid, nannozinone B, anaephene B, xanthone, ajudazol, althiomycin, aurachin E, carolacton, cystobactamid, disciformycin A, and melithiazol A; comparable bond lengths were observed across these complexes. Violacein, xiamycin, sulfangolid C, and chlorotonil A formed between 0 and 4 H-bonds but did not interact directly with the active-site residues (Supplementary Table S1). The docking results indicated that the top-scoring ligands interacted with key catalytic residues, including Arg279 and Tyr264. These interactions suggest a possible active-site binding mode, which should be validated through enzyme kinetics studies.

The detailed interactions of compounds forming at least two H-bonds with the active-site residues are presented in Figure 2. Noricumazole C formed six H-bonds with the enzyme, including conventional H-bonds with Thr224, Leu225, Ile227, and the active-site residue Arg279 (bond length: 2.16 Å). A similar pattern was observed for ascorbic acid, which formed six H-bonds—three of these bonds involved Arg279 (bond lengths: 2.02, 2.70, and 3.01 Å). Disciformycin A formed five H-bonds with Gln261, Gly262, and Tyr264, with bond lengths of 2.62, 2.49, and 2.10–2.98 Å, respectively. Similarly, brefeldin A formed two H-bonds with the active-site residues Gln261 (2.42 Å) and Tyr264 (2.06 Å). In addition to hydrogen bonding, these compounds exhibited hydrophobic and electrostatic interactions with the enzyme (Figure 2).

3-D (Left) and 2-D (Right) representations of the binding pose of ascorbic acid (A), bikaverin (B), violacein (C), xiamycin (D), chlorotinil (E), phenoxan (F), and kulkenon (G) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. 3-D (Left) and 2-D (Right) representations of the binding pose of ascorbic acid (A), bikaverin (B), violacein (C), xiamycin (D), chlorotinil (E), phenoxan (F), and kulkenon (G) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

The in silico ADMET prediction analysis of the 26 compounds (including ascorbic acid) that interacted with at least one active-site residue is shown in Supplementary Table S2. Most compounds satisfied Lipinski’s rule and were not P-glycoprotein or CYP3A4 substrates. They were also predicted not to penetrate the blood–brain barrier. These compounds exhibited sufficient water solubility (−4.0–0.5 log mol/L) and lipophilicity (logD: 1–3), with moderate (5–10) clearance (Supplementary Table S2). Only eight compounds had suitable (0–3) n-octanol/water distribution coefficients. All compounds exhibited excellent (0.04–20 L/kg) volume distributions, with most having medium (0.3–0.7) to excellent (0–0.3) half-lives. Although a few compounds were predicted to cause eye irritation, many compounds were predicted to lack respiratory toxicity, human hepatotoxicity, or drug-induced liver injury and are non-carcinogenic. Only 10 compounds showed excellent (<90%) plasma protein-binding capacity.

MD simulations were conducted for 100 ns to gain insight into the dynamic stability and conformational behavior of the protein–ligand complexes. Rg analysis revealed initial expansion followed by structural compaction of the protein–ligand complexes during the simulation. The violacein- and bikaverin-bound complexes exhibited Rg peaks of approximately 2.8 nm, followed by a downward trend with noticeable fluctuations during the first 20 ns. Subsequently, they stabilized at approximately 1.8 nm and maintained this compact conformation throughout the 100-ns simulation. Chlorotonil A showed a high initial Rg peak of approximately 3.2 nm, which gradually decreased to approximately 1.8 nm over during the simulation, indicating progressive structural stabilization. Similarly, the kulkenon complex displayed an elevated Rg peak near 2.5 nm during the first 20 ns, before converging toward a stable, compact state. Phenoxan- and xiamycin-bound complexes followed comparable trends, beginning at around 2.4 nm and reaching transient peaks near 2.8 nm within the first 20 ns. These complexes subsequently maintained stable Rg values of approximately 1.8 nm for the remainder of the simulation period (Figure 3).

Radius of gyration (Rg) profiles from molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. Radius of gyration (Rg) profiles from molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

RMSD analysis revealed distinct equilibration and stability patterns among the protein–ligand complexes during the 100-ns simulation. The violacein-bound complex exhibited an RMSD peak of approximately 2.0 nm at around 20 ns, followed by pronounced fluctuations up to approximately 70 ns. Thereafter, the RMSD showed a downward trend and stabilized at approximately 1.5 nm toward the end of the simulation, indicating the attainment of a stable conformational state. The bikaverin complex displayed an early RMSD peak of approximately 1.5 nm within the initial phase of the simulation, followed by increased fluctuations and an upward peak exceeding 2.5 nm at approximately 20 ns. The complex then maintained a relatively stable RMSD signal around 2.5 nm throughout the remainder of the 100-ns period.

A similar overall stabilization pattern was observed for chlorotonil A, kulkenon, phenoxan, and xiamycin complexes, which converged to stable RMSD values of approximately 1.7–1.8 nm over during the simulation. Specifically, the chlorotonil A complex exhibited a high initial RMSD peak at approximately 3.2 nm, followed by a gradual decrease to approximately 1.8 nm at the end of the simulation. Kulkenon showed a comparable trend, with a pronounced RMSD peak observed near 2.5 nm during the first 20 ns before stabilization. Phenoxan- and xiamycin-bound complexes followed similar trajectories, beginning at approximately 2.4 nm, reaching transient peaks near 2.8 nm within the first 20 ns, and subsequently stabilizing around 1.8 nm for the remainder of the simulation (Figure 4).

Root-mean-square deviation (RMSD) profiles obtained from the molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. Root-mean-square deviation (RMSD) profiles obtained from the molecular dynamics simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

RMSF analysis revealed that key residues of the protein, including those located within the catalytic region, exhibited consistently low fluctuations throughout the simulation period (Figure 5).

Root-mean-square fluctuation (RMSF) profiles obtained from MD simulations of violacein (A), bikaverin (B), chlorotonil-A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase derived from S. pyogenes.

3D FEL analysis was used to evaluate the stability and dynamics of protein–ligand complexes during the simulation. FEL revealed well-defined low-energy basins corresponding to thermodynamically favorable stable conformational states. These minima represent the dominant conformations of the complexes, which exhibited limited conformational space and reduced flexibility owing to ligand binding. The lack of significant high-energy barriers suggested smooth transitions and overall stability during the simulation (Figure 6).

3D free-energy landscapes (FEL) profiles obtained from MD simulations of violacein (A), bikaverin (B), chlorotonil-A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. 3D free-energy landscapes (FEL) profiles obtained from MD simulations of violacein (A), bikaverin (B), chlorotonil-A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

Per-residue energy decomposition analysis was performed to identify key amino acid residues that contribute to ligand binding. Residues Phe339, Ser353, and Gln354 showed favorable binding contributions in the violacein complex, each exhibiting interaction energies of approximately −2.0 kcal/mol. The bikaverin complex displayed notable contributions from Ile331, Asn338, and Lys338, with interaction energies of approximately −2.0 kcal/mol, indicating stable residue–ligand interactions. The chlorotonil A complex showed a dominant contribution from Phe347, with a comparatively higher interaction energy of approximately −4.0 kcal/mol, suggesting strong residue involvement in ligand stabilization (Figure 7).

Per-residue energy decomposition profiles from MD simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes. Per-residue energy decomposition profiles from MD simulations of violacein (A), bikaverin (B), chlorotonil A (C), kulkenon (D), phenoxan (E), and xiamycin (F) in complex with bacteriophage-encoded hyaluronate lyase from S. pyogenes.

The kulkenon complex exhibited minimal residue-specific contributions, with interaction energies close to 0 kcal/mol, indicating weaker or more distributed interactions across the binding site. Phe339 was a major contributor to the phenoxan complex, with an interaction energy of approximately −4.0 kcal/mol, whereas Phe347 showed a moderate contribution (−2.0 kcal/mol). Similarly, residues Phe339, Ile331, and Phe347 contributed favorably to binding in the xiamycin-bound complex, each displaying interaction energies of approximately −1.0 kcal/mol. The ligand atom LIG387 consistently exhibited highly negative energy contributions across all complexes, highlighting its critical role in the stabilization of protein–ligand interactions (Figure 7).

Discussion

4

In this study, a few compounds had slightly lower binding affinities than previously reported inhibitors (ascorbic acid, Triton X-100, and sodium dodecyl sulfate). Most compounds interacted with S. pyogenes-derived hyaluronate lyase with higher binding affinities. The top-performing compounds (violacein, −7.7 kcal/mol; sulfangolid C, −7.4 kcal/mol; chlorotonil A, −7.4 kcal/mol; xiamycin, −7.3 kcal/mol; kulkenon, −7.1 kcal/mol) exhibited binding affinities broadly comparable to natural-product docking results reported in other antibacterial target studies. For instance, Chao et al. (49). evaluated natural compounds against three Streptococcus-associated targets (GlmU, PPAT, and RpoD). They reported binding affinities of −9.263 to −8.622 kcal/mol for GlmU, −8.038 to −7.25 kcal/mol for PPAT, and −7.154 to −6.521 kcal/mol for RpoD. The binding energies of the lead compounds identified in the present study align most closely with the affinity range reported for PPAT (−8.038 to −7.25 kcal/mol) and overlap with the upper range of affinities observed for RpoD (−7.154 to −6.521 kcal/mol). Although our compounds did not reach the extremely strong affinities obtained for GlmU (−9.263 to −8.622 kcal/mol), they demonstrated clinically relevant binding strengths, particularly considering that hyaluronate lyase is an anti-virulence rather than a bactericidal target.

Anti-virulence inhibitors are often reported to exhibit moderate binding energies. This ensures stable engagement with the catalytic pocket without exerting bactericidal pressure, which accelerates resistance development. Hassan et al. (50) identified a top hit of −7.2 kcal/mol [Sodium (1,5-dihydroxy-2-oxopyrrolidin-3-yl)-hydroxy-dioxidophosphanium] against S. pneumoniae alpha-enolase. Violacein, sulfangolid C, and chlorotonil A, the strongest compounds in this study, showed affinities that slightly exceeded or matched the highest affinities reported by Hassan et al. (50). This indicates that the natural molecules identified in this study interact with hyaluronate lyase at a binding strength comparable to that of validated in silico hits against other streptococcal targets.

Rehman et al. (51) have reported extremely high binding affinities for natural compounds docked to two essential S. pyogenes proteins. The interactions of DnaA–Sophorastilbene A and TCR protein–Aloin B had binding affinities of −21.31 kcal/mol and −18.02 kcal/mol, respectively. These strong binding scores are associated with essential replication and transcription proteins that possess highly conserved deep nucleotide-binding pockets. These proteins often display high docking scores owing to large hydrophobic cavities. In contrast, hyaluronate lyase is a secreted protein with a surface-exposed catalytic cleft designed for processing long glycosaminoglycan chains rather than for tight ligand encapsulation (52, 53). A similar pattern was observed by Rivera-Pérez et al. (54), who recorded markedly strong docking scores for inhibitors of the phosphotransferase system of S. mutans. Of these, ZINC15958489 and ZINC15968797 exhibited the highest binding affinities at −13.0 kcal/mol and −12.9 kcal/mol, respectively. These significantly stronger affinities correspond to molecules interacting with deeply recessed binding pockets in a membrane-associated transferase, which are considerably different from the topology of the hyaluronate lyase cleft.

Hydrogen bonding provides clear insights into the mechanisms whereby a molecule binds to an active site. Hydrogen bonding is instrumental in complex molecule inhibition because it ensures the stability of structures and functions (55). In this study, 26 compounds formed at least one H-bond with an active-site residue. Several top-ranking compounds interacted directly within the defined active-site pocket of hyaluronate lyase, including key catalytic and substrate-recognition residues (Gln261, Tyr264, and Arg279). Noricumazole C, disciformycin A, althiomycin, melithiazole A, griseofulvin, citreorosein, and bikaverin formed multiple H-bonds and hydrophobic contacts with these residues (particularly Arg279). These results suggest a stable and well-anchored binding mode. Their interaction patterns closely resembled those observed for ascorbic acid (11), supporting the relevance of their predicted orientations.

Engagement of active-site residues is especially important for hyaluronate lyase because inhibition requires obstruction of the catalytic groove responsible for glycosaminoglycan cleavage. Therefore, the ability of these compounds to interact with residues located within or immediately adjacent to the catalytic site provides a mechanistic rationale for their predicted inhibitory potential rather than relying solely on docking scores (55). In addition to hydrogen bonding, the stability and specificity of the ligand binding to hyaluronate lyase are also significantly affected by hydrophobic interactions, π–π stacking, and Van der Waals forces. Many of the top compounds, including violacein, xiamycin, sulfangolid C, chlorotonil A, and kulkenon, possess aromatic and nonpolar moieties that interact with the hydrophobic pockets within the catalytic cleft (11, 55). Aromatic rings in violacein and xiamycin facilitate π–π stacking interactions with aromatic residues, further stabilizing the ligand orientation and contributing to proper alignment with the catalytic residues Arg279 and Tyr264 (9, 12, 13).

Van der Waals forces fine-tune the ligand–enzyme complex by enabling close-range, weakly attractive interactions that complement H-bonds and hydrophobic contacts (56). These non-covalent interactions act synergistically to secure the ligand within the substrate-binding cleft and potentially enhance competitive inhibition by mimicking the interactions of natural substrates (11). Therefore, evaluating both H-bonds and these additional interactions elucidates the molecular determinants underlying the efficacy of hyaluronate lyase inhibitors. The predicted interactions of the top hits with the catalytic residues Arg279 and Tyr264 suggested that these compounds occupy the active site and potentially interfere with substrate binding. As molecular docking remains a predictive approach, biochemical assays such as enzyme inhibition kinetics are required to validate the inhibitory mechanisms and binding modes of the identified compounds.

To evaluate the biological activity of these compounds and predict whether they have favorable or harmful effects when utilized in pharmaceutical applications, ADMETLab2.0 was used to assess drug-like properties. Most of the compounds that engaged in H-bond interactions with at least one of the active-site residues passed Lipinski’s Rule of Five. Lipinski’s Rule of Five, which was used to determine whether the active compounds were orally active (57). The compounds were predicted to have good solubility, although only 10 had a plasma protein-binding capacity <90%, which may affect bioavailability. Early evaluation of drug solubility is crucial for drug development because low solubility hinders effective and complete oral absorption (58). Plasma protein-binding capacity can directly impact oral bioavailability because the free concentration of the drug is at risk when it binds to serum proteins (59). Overall, the compounds possessed good toxicological profiles, as they did not induce respiratory toxicity, human hepatotoxicity, or drug-induced liver injury, and were predicted to be non-carcinogenic.

Drug-induced liver damage has emerged as the most prevalent safety issue associated with drug withdrawal over the past 50 years and remains a primary reason medication are removed from the market to ensure patient safety (60). Minimal or no blood–brain barrier penetration may be necessary for medications with peripheral targets to prevent adverse effects on the central nervous system (61). The late and expensive discontinuation of drug development initiatives is frequently caused by adverse liver consequences in clinical trials (62, 63).

A smaller subset of compounds emerged as the most promising leads when the molecular docking results were integrated with ADMET and drug-likeness assessments. Violacein, xiamycin, kulkenon, chlorotonil A, and bikaverin demonstrated high binding affinities and favorable pharmacokinetic properties, including full compliance with Lipinski’s rule. Therefore, these molecules represent the lead candidates, combining strong predicted target engagement with acceptable drug-likeness and safety profiles. Phenoxan also exhibited good ADMET characteristics, although its binding energy was slightly lower than that of the top-performing compounds.

MD simulations added a temporal dimension to the docking results and permitted assessment of complex stability under simulated physiological conditions. The reliability of violacein, bikaverin, chlorotonil A, kulkenon, phenoxan, and xiamycin was dynamically validated by monitoring their interactions with the protein–ligand complex during a 100-ns MD trajectory.

The RMSD plots typically showed initial variation followed by convergence. This stabilization indicated that the ligands were firmly attached to the catalytic cleft and did not dissociate. Violacein and xiamycin showed synchronized trajectories after equilibration, indicating that their predicted binding orientations were reliable inhibitory poses. Rg plots for violacein and kulkenon showed a downward trend followed by a plateau, indicating that the protein adopted a denser form upon ligand binding. This structural compaction may isolate active-site residues from the surroundings, hindering natural substrate access to the production machinery. Local flexibility analysis via RMSF showed that important residues in the catalytic region exhibited low fluctuations throughout the simulation. The decrease in flexibility indicates that the presence of these natural compounds effectively “locks” the catalytic region, possibly preventing the conformational changes necessary for enzymatic activity.

The 3D FEL provided a strong thermodynamic argument for these interactions. The complexes were located at considerably deep and distinct energy minima. Per-residue energy decomposition indicated that the ligand molecules (LIG387) were primarily responsible for the largest negative total energy contribution to system stabilization. Residues such as Phe339 and Phe347 are pivotal anchors at the energy level of several complexes.

MD simulations identified violacein, xiamycin, and kulkenon as the most favorable compounds because of their highly stable and consistent convergence throughout the 100-ns trajectory. Rg profiling substantiated the maintenance of the compact structural state of the protein. FEL confirmed thermodynamic favorability; These compounds had markedly deep energy minima, further supported by their significant negative total energy contributions. These results indicate that these compounds were the most promising candidates for further experimental validation.

This study had some limitations: First, it was based entirely on in silico computational approaches. Therefore, the findings should be interpreted as predictive rather than conclusive. Molecular docking provides a theoretical estimate of ligand–protein binding affinity, but these scores do not directly correspond to actual enzymatic inhibition. Docking also assumes a largely rigid protein structure and may not fully capture the conformational flexibility of hyaluronate lyase or the dynamic nature of the enzymes long surface-exposed catalytic cleft. This was complemented by MD simulations. Additionally, ADMET predictions are based on quantitative structure–activity relationship models that carry inherent uncertainty and may not fully represent the actual pharmacokinetic behavior in biological systems. Consequently, the identified lead compounds require experimental validation using biochemical assays, enzyme inhibition studies, and cytotoxicity profiling to confirm their biological activity, potency, and safety. Further in vitro and in vivo studies are required before these candidates can be considered as viable therapeutic inhibitors.

Conclusion

5

The increasing prevalence of antibiotic resistance poses a major global health challenge, highlighting the need for alternatives to conventional bactericidal agents. Targeting virulence factors, such as the hyaluronate lyase derived from S. pyogenes, offers a promising anti-virulence approach that may reduce pathogenicity without exerting selective pressure for resistance. In this study, structure-based molecular docking, ADMET analysis, and 100-ns MD simulations were employed to evaluate the binding potential of 118 natural compounds.

Virtual screening identified violacein, chlorotonil A, sulfangolid C, xiamycin, and kulkenon as high-affinity candidates with favorable ADMET properties. Mechanistic analysis revealed that these top hits consistently interacted with the key catalytic residues Arg279 and Tyr264. MD simulations further validated these findings. Violacein, xiamycin, and kulkenon demonstrated exceptional structural stability and convergent trajectories, effectively “locking” the catalytic region in a compact non-functional state. The direct engagement with the catalytic pocket and stable energetic minima identified through the FEL suggest that these compounds act as competitive inhibitors, thereby preventing the enzyme from processing host glycosaminoglycans. This mechanism aligns with those of anti-virulence therapy strategies and has the potential to reduce host tissue degradation and bacterial dissemination without imposing selective pressure on the bacteria.

Although these findings offer a robust foundation for the development of hyaluronate lyase-targeted inhibitors, they remain predictive and are derived from computational methods. Therefore, in vitro enzyme assays, cytotoxicity screening, and in vivo infection models are required to validate these molecules as true inhibitors and to establish their therapeutic potential. Overall, this study identified promising natural product scaffolds and provided dynamic mechanistic insights into their modes of action, supporting the continued exploration of anti-virulence strategies against S. pyogenes.