1. Introduction
Hantaviruses are a group of rodent-borne viruses that have caused severe, often fatal outbreaks globally for decades, yet remain without any specific antiviral treatment [1]. Depending on the region and the specific viral strain, infection produces one of two devastating syndromes: hemorrhagic fever with renal syndrome (HFRS)—first recognized among American soldiers near the Hantaan River during the Korean War, and later linked to the Finnish Puumala virus—or hantavirus (cardio)pulmonary syndrome (HPS/HCPS), first identified after the 1993 Four Corners outbreak of Sin Nombre virus (SNV) in the United States [2] [3]. HPS, in particular, is considered alarming due to its high case fatality rate (CFR), with the rate approaching 29% - 40%: among the highest of any known pulmonary viral disease [4] [5]. Andes virus (ANDV), the primary virological focus of this study—and not in the least due to its relevance with the recent outbreak—causes HPS across South America, and is of notable concern amongst the New World hantaviruses due to it being the only hantavirus with documented person-to-person transmission, giving it a higher epidemic potential than many related strains such as SNV [1] [6]. Despite the pathophysiological severity of HPS, contemporary treatment remains entirely supportive; no drug has been marketed which specifically targets its etiological viruses or their mechanisms of entry [7].
In this regard, pathogenic New World hantaviruses like ANDV and SNV hijack a receptor called αvβ3 integrin, which normally assists in anchoring vascular endothelial cells to the surrounding extracellular matrix through arginine-glycine-aspartate (RGD)-motif interactions with structural proteins such as vitronectin and fibronectin [8] [9]. Furthermore, αvβ3 appears to play a critical role in facilitating vascular endothelial growth factor receptor-2 (VEGFR-2)-mediated angiogenesis, abetting cancer metastasis (i.e. prostate cancer), as well as mediating osteoclast-bone adsorption for proper bone resorption [10]-[12]. The ANDV/SNV viral glycoprotein (Gn/Gc) binds the aforementioned receptor and effectively occupies it, preventing the endothelial cell from adhering properly to the surrounding ECM, and compromising vascular structural integrity [13]. This results in “leaky” blood vessels, which—when occurring in the alveolar capillaries of the lungs—produces pulmonary edema: a leading cause of death in HPS [8] [13]. Because this attachment step happens at the cell surface, prior to viral replication yielding progeny virions, it remains an attractive interest point for pharmacological inhibition—if the receptor itself can be competitively inhibited, the pathophysiological proceedings of HPS clinical infection can potentially be halted in the earliest stages [8]. In particular, while ECM proteins (i.e. vitronectin and fibronectin) have been shown to bind the RGD binding site, New World/HPS-etiological hantaviruses, such as New York 1 virus (NY-1V) and SNV—as well as Old World viruses such as Hantaan virus (HTNV)—bind to a distinct region of the β3 subunit of the α-β heterodimer known as the pleckstrin-semaphorin-integrin domain (PSI), independent of the RGD binding site [13]-[15] (Figure 1). While New World hantaviruses bind the aforementioned αvβ3 subtype—and thus foremost cause pulmonary symptomatology where capillary integrity is among the most delicate—Old World hantaviruses preferentially bind the αIIβ3 subtype found exclusively on platelets, thus correlating with thrombocytopenic, hemorrhagic manifestations [13] [15] [16].
Notably, Raymond et al. [15] report that the PSI domain of αvβ3 remains exposed at the apical surface when the integrin is in an inactive, ‘bent’ conformation, and that this calcium-driven conformation facilitates cellular infection by hantaviruses such as NY-1V and HNTV (Figure 2).
Serris et al. [14] further determined the crystal structure of the ANDV Gn head/Gc heterodimer, and used these structures together with a cryo-electron tomography map of the related Tula virus (TULV) to build a complete atomic model of how Gn/Gc dimers assemble into the tetrameric spikes that cover the surface of the virus. This work likewise furthered the theory that Gc is a fusion protein that needs to undergo a large conformational change to drive membrane fusion during entry, while Gn actively restrains this change, acting as a ‘safety-catch’ that keeps Gc in its pre-fusion form until the virus reaches the acidic environment of the host-cell endosome during endocytosis [14]. This structural detail, while concerning the viral aspect of the interaction between hantaviruses and host endothelial cells expressing αvβ3, remains useful to this study by implying which surface of Gn engages the αvβ3 receptor, giving a clearer indication of what receptor-binding molecules requires to structurally/functionally compete against.
In the above regard, while no studies have explicitly depicted a 3D crystallographic hantavirus Gn/β3 interaction from any strain on its own, the ample data provided by Serris et al. [14] and Guo et al. [21] have allowed for a demonstrative protein-protein docking to be conducted, portraying the putative interactions between ANDV Gn and the β3 subunit of αvβ3 via the PSI domain (Figure 3).
There remains, moreover, a strong precedent showing that a general strategy
Figure 1. Crystallographic ribbon structures of the αvβ3 heterodimer depicting highlighted αv subunit (green) (a), the RGD binding site in complex with cyclo-RGD (b), and the pleckstrin-semaphorin-integrin (PSI) domain, spanning amino acid residues 1 - 53 on the β3 subunit (brown), used by pathogenic hantavirus attachment and entry (c) [PDB IDs: 1L5G (1a, 1b) & 1U8C (1c)]. Crystallographic structures were adapted from Xiong et al. [17] [18].
of disrupting the above hantavirus/β3 interaction with small molecules that target the receptor itself yields favorable results in literature. For instance, Hall et al. [22] began this endeavor from a phage-display-derived cyclic nonapeptide, cyclo-[CPFVKTQLC], and used alanine scanning and stepwise deletion to truncate it to a minimally-effective cyclic pentapeptide, cyclo-[CPFVC], that retained most of the inhibitory activity against SNV, specifically. Hall et al. [23] furthered their findings, using the peptide as a template for ligand-based virtual screening (LBVS) against a large chemical library, refining a four-point pharmacophore model over two screening rounds. This eventually produced peptidomimetic small molecules active in the nanomolar (nM) range against ANDV, SNV, and HTNV infection
![]()
Figure 2. The three primary conformations of αvβ3: inactive (bent), hybrid, and active (open), with the PSI domain labeled on the β3 subunit. Reproduced with permission from Agarwal et al. [19] (CC-BY 4.0).
Figure 3. Ribbon diagrams of Zhang et al.’s [20] HawkDock-computed protein-protein docking of β3—Gn, depicting: (a) putative hydrophobic interactions between the PSI domain of β3 (red—Phe53 [auth. Phe1056] and Ile51 [auth. Ile1054]) and Phe440/Leu443 on ANDV Gn (blue ribbon), and (b) putative salt-bridge interaction between the PSI domain and β3 (red—Glu42 [auth. 1042]) and Lys2 on ANDV Gn (blue ribbon), adapted from Serris et al. [14] and Guo et al. [21]. PDB IDs: 1U8C (β3) and 9P31 (Gn).
in vitro through the β3-inhibition pathway, while remaining inactive against the non-pathogenic Prospect Hill virus (PHV) which uses a different integrin type for cellular entry (β1), corroborating the proposed mechanism of action [23]. Among these compounds was G319-0078, the most efficacious and one of the most potent of the hits identified, with an inhibition rate of 46.5% in αvβ3-SNV bioassay, and IC50 = 1.6 nM (±0.6 nM): 100-fold more potent than the original cyclo-[CPFVC] parent compound [23] (Figure 4).
The binding site of cyclo-[CPFVC] and G319-0078 on αvβ3, moreover, was discovered to be at a distinct portion of the β3 subunit previously known to bind the ReoPro (Abciximab, murine 7E3 IgG) chimeric monoclonal antibody,
Figure 4. A schematic representation of the structural formulae of cyclo-[CPFVC] and G319-0078 (a), and ribbon structure obtained via mcule.com’s AutoDock Vina of cyclo-[CPFVC] (purple) and G319-0078 (yellow) binding to the β3 subunit (brown Gaussian volume) of αvβ3 integrin (b), adapted from Hall et al. [22] [23]. Note the crucial interaction between Trp1129 (highlighted green), and the phenylalanine ring of cyclo-[CPFVC] and benzyl ring of G319-0078 [23].
originally developed to inhibit thrombotic platelet activity via αIIβ3 integrin (though being cross-reactive with αvβ3 at the identical binding site) during percutaneous coronary procedures and acute coronary syndrome [24] [25]. Notably, the binding site of Abciximab was found by Nešić et al. [26] to be adjacent to the binding site of the γ-chain (Res. 404-411) of fibrinogen on αIIβ3, with the binding site being corroborated by van den Kerhof et al. [27] the following year (Figure 5).
The rationale for this study, as such, was to directly edify the above precedent: rather than targeting the virus’ own glycoprotein (Gn/Gc), of which information
Figure 5. A 3D crystallographic ribbon diagram of ReoPro (Abciximab) antibody (pink) complexed onto αvβ3—with the αv subunit in green and β3 subunit in brown—and fibrinogen γ-chain binding (highlighted green/labeled red), adapted from Nešić et al. [26] (PDB ID: 2VDO).
in literature concerning a viable antiviral target site in isolation—while available in the form of fusion-inhibiting stem peptides, mimicking domain III of ANDV Gc—was scarce, the targeting the αvβ3 integrin receptor itself was deemed the more promising option in order to prevent hantavirus Gn/Gc from binding to it [28]-[30]. This, in part, is due to the different serotypes of hantaviruses, which vary in their epitope amino acid sequences by up to 44% (particularly comparing New and Old World viruses), with cross-reactive efficacy of direct Gn/Gc inhibitors being observed notably in binding the more-conserved Gc domains as compared to more variable receptor-binding Gn [31] [32]. The above sequence differences likewise vary in ANDV and SNV Gn/Gc by up to 22%, possibly explaining the different pathophysiologies and infectivities of the two HPS-causing New World strains [33] [34]. The β3 subunit, on the other hand, remains conserved across αvβ3 and αIIβ3 alike, thus offering potential further druggability extended to other, non-HPS-etiological, Old World hantavirus strains alike [35].
In the aforementioned regard, cyclophosphinamide structural analogues of G319-0078, the lead compound identified by Hall et al. [23], were designed in silico with the aim of improving its potency at the β3 subunit of αvβ3. Cyclophosphinamides were chosen, during the course of this study for optimized in silico lead compounds, due to the discovery that cyclization of the peptidomimetic amide tether scaffold afforded better docking scores than uncyclized derivatives, or parent compound G319-0078. Moreover, phosphinamides have long been known to be good bioisosteres of amides and sulfonamides in vitro and in vivo, showing optimized lipophilicity and improved solubility compared to the former [36] [37]. Phosphinamide-based hydroxamic acids have, moreover, been shown to effectively inhibit matrix metalloprotease-1 (MMP-1), a collagenase implicated in tumor growth and metastasis, arthritis, and multiple sclerosis (MS) [36]. Additionally, the cyclization of peptides and peptidomimetics has been found to render them more resistant to enzymatic hydrolysis and more membrane permeable than acyclic peptides, improving bioavailability and overall pharmacokinetic scores [38] [39].
2. Materials and Methods
2.1. Blind Docking and Ligand-Based Virtual Screening (LBVS)
A two-stage in silico screening strategy was adopted, combining an initial unbiased blind docking step with a subsequent ligand-based virtual screening (LBVS) round, in order to both independently corroborate the previously reported peptidomimetic binding site on αvβ3, and to generate a structurally diverse panel of bioisosteric analogues for comparative evaluation. Blind docking of the parent ligand, G319-0078, to the crystallographic structure of αvβ3 integrin (PDB ID: 1U8C) was first performed using the SwissDock web server [40] [41], which replaced the previous EADock DSS engine with a choice of Attracting Cavities and an integrated AutoDock Vina (ADVina) engine for pose scoring. As SwissDock performs its search across a large protein surface (search-box size up to 20 Å), this blind, hypothesis-free approach was deliberately chosen over a pre-defined local search so that the resulting binding cavity could be treated as an independent, computationally-derived confirmation of the orthosteric site, rather than an assumption carried over strictly from prior literature. Because SwissDock’s docking poses are scored using ADVina—the same engine subsequently used for the LBVS docking of analogues—both stages of the screening pipeline shared a common, directly-comparable scoring function throughout.
The docking cavity with the most favorable docking score returned by SwissDock was found to correspond to a binding pocket with coordinates (X: 27.502, Y: 22.528, Z: 49.248), matching the site previously localized experimentally and computationally by Hall et al. [23], thereby validating the use of this cavity as the fixed search space for all subsequent docking runs. Notably, whereas the original identification of this site by Hall et al. [23] relied on a qualitative, four-point pharmacophore model derived from LBVS, the present blind docking step provided a quantitative, physics-based binding-affinity estimate, offering independent energetic corroboration of the site beyond the pharmacophore-level agreement previously available.
Having established and fixed the coordinates of the orthosteric binding cavity, a focused LBVS round was carried out to evaluate the effect of systematic bioisosteric substitution on the benzyl and sulfonamide moieties of G319-0078. G319-0078 was re-docked via mcule.com’s “1-Click Docking” tool [42], which similarly runs ADVina under default exhaustiveness (i.e. 8) and pH 7.4, and additionally cross-validates the identity of each output conformer against the input ligand’s InChI string to exclude docking artifacts arising from unintended structural modification during pose generation. A total of 39 structural analogues (L2 - L40) were designed around the fixed 1,3-benzodioxole and alkyl-amide backbone of the parent compound: 27 monosubstituted analogues (L2 - L28) bearing halogen, alkyl, heteroaryl, and polar-group substitutions at the benzylic position, and 12 disubstituted analogues (L29 - L40) combining a benzylic heteroaryl substitution with an additional modification of the central sulfonamide linker (i.e. selenonylamide, phosphonamide, or extended fused-heteroaryl systems). Each analogue was docked individually against the fixed binding cavity via the same software, and the four top-scoring docking poses retained, with the best/largest (most negative) docking score recorded and collated against the docking score of the re-docked parent ligand, G319-0078, allowing for comparative analysis of each substitution to predicted binding affinity.
2.2. Computational Lead Optimization Strategy
Building on the initial LBVS round (2.1), a second round of in silico lead optimization was undertaken to refine the two top-performing hits, L34 and L30, into a focused set of cyclophosphinamide candidates. L34 and L30 were selected as optimization templates because they represented a disubstituted heteroaryl/sulfonamide scaffold and a phosphonamide-substituted scaffold, respectively, offering two structurally-distinct starting points from which to explore the cyclophosphinamide bioisostere rationale outlined in the Introduction. Modifications to these templates were designed systematically along four axes, each grounded in an established medicinal chemistry rationale.
First, a hydrogen-bond donor/acceptor (HBD/HBA)-orthogonal substituent—specifically a phthalimide group—was introduced at the amide N-alkyl position, benchmarked against nitrogen-free and non-orthogonal alternatives. This design choice was informed by the well-documented asymmetry between HBDs and HBAs in ligand-receptor recognition, whereby orthogonally positioned HBD/HBA pairs can engage independent, non-competing polar contacts on a target surface without incurring the de-solvation penalties typically associated with donor-heavy substituents [43]. Second, the heteroaryl versus aliphatic character of ring systems introduced at both the N-alkyl and benzylic positions were systematically varied, guided directly by the substitution trends observed during the initial LBVS round (2.1). Third, the peptidomimetic amide tether was cyclized via incorporation of an aryl phosphinamide functional group, converting the acyclic sulfonamide/amide backbone of the parent scaffold into a cyclic phosphinamide ring system. This modification was grounded in the established use of phosphonamides as amide/sulfonamide bioisosteres offering favorable lipophilicity and solubility profiles, combined with the well-documented pharmacokinetic advantages of macrocyclization (see Introduction) [36] [37]. Finally, ring expansion/contraction and fluorination/trifluoromethylation of heteroaryl substituents at the benzylic and N-alkyl positions were explored, reflecting the established capacity of fluorine substitution to modulate lipophilicity, metabolic stability, and binding conformation in drug design [44].
In addition, an alternative terminal functionalization at the benzylic position was explored by substituting the carbamate group present in the initial designs with a cyclic hydroxamic acid moiety, in order to assess the effect of introducing an additional terminal HBD at this position; hydroxamic acid functional groups have a well-established history as potent HBD/metal-chelating pharmacophores in medicinal chemistry, including within phosphinamide-based scaffolds [36]. Combining these four design axes yielded five candidate cyclophosphinamide leads, designated X47-A, X47-B, X47-C, X47-D, and X47-E (Figure 6), each incorporating a distinct combination of the above modifications while retaining the core phthalimide-bearing cyclophosphinamide scaffold. Each of the five candidates was subsequently docked against the same fixed αvβ3 orthosteric binding cavity used throughout this study, following the mcule.com “1-Click Docking” ADVina protocol described (2.1). For each candidate, the four top-scoring docking poses were retained per ligand, and both the single best-scoring pose and the average docking score across all four poses were recorded in order to provide a more robust estimate of predicted binding affinity.
2.3. Molecular Dynamics Simulation (MDsim)
Given that the platform used for this study, SiBioLead [45], could not reliably parameterize and simulate the complete αvβ3 heterodimer owing to the size (>1500 Da), multidomain architecture, and conformational flexibility of the full integrin ectodomain, molecular dynamics simulation (MDsim) was restricted to the isolated β3 subunit alone, which contains both the PSI domain and the orthosteric binding cavity characterized throughout Sections 2.1 and 2.2. This scope restriction reflects a common practical trade-off encountered in MD simulation of large, multidomain protein complexes, where computational cost and force-field parameterization reliability generally scale unfavorably with system size, motivating the use of a reduced, computationally-tractable subsystem centered
![]()
Figure 6. Skeletal structures of the five lead ligands developed: X47-A, X47-B, X47-C, X47-D, and X47-E. The variable benzylic groups are designated in red, and the phosphinamide heteroaryl moiety in blue; the constant phthalimide scaffold is depicted in black.
on the region of direct pharmacological interest [46]. Because only the isolated β3 subunit was simulated, all conclusions drawn from this stage of the study characterize the stability of the ligand-β3 subunit complex, independent of any potential allosteric or quaternary-complex effects that the αv subunit might exert within the intact heterodimer (see 2.6).
For each lead candidate (X47-A—E), a corresponding re-docked ligand-β3 complex was used as the starting structure for a single run of MDsim, performed using GROMACS software hosted on SiBioLead. The protein-ligand system was parameterized using the OPLS all-atom force field [47] for optimal small-molecule parameterization, and solvated within a dodecahedral simulation box using a three-point simple point charge (SPC) water model; the dodecahedral geometry was selected to minimize solvent volume relative to a conventional cubic box, while maintaining a sufficient buffer around the solute [48]. The system was subsequently brought to physiological ionic strength through the addition of 150 mM NaCl. Prior to production, energy minimization (EM) was carried out over a maximum of 5000 steps using the Steepest Descent integrator, followed by sequential NVT (constant volume/temperature) and NPT (constant pressure/temperature) equilibration phases of 1 nanosecond (ns) each in order to stabilize the system. The pre-production temperature/pressure was, moreover, set to 310 K and atmospheric pressure (1.0 bar), respectively, to approximate human physiological conditions in vivo. Production MD simulations were then carried out for 50 ns per system using the Leap-Frog integrator. Three result criteria were established thereof for comparative analysis: root mean square deviation (RMSD, measuring the average structural deviation of the ligand-protein complex over the MDsim timeframe), root mean square fluctuation (RMSF, recording average individual amino acid fluctuations), and solvent-accessible surface area (SASA, indicating the surface area of the protein in direct contact with the solvent).
2.4. Normal Mode Analysis (NMA) and Dynamic Cross-Correlation Matrix (DCCM)
To complement the explicit-solvent, time resolved assessment of protein-ligand stability provided by MDsim (2.3), normal mode analysis (NMA) and dynamic cross-correlation matrix (DCCM) computation was performed for each of the five lead ligand-β3 complexes using DynaMut [49]. DynaMut employs the Elastic Network Contact Model (ENCoM), a coarse-grained NMA approach which—in contrast to conventional distance-based elastic network models—additionally weights inter-residue interactions by atom-type-dependent potentials, thereby capturing the anisotropic contribution of side-chain packing to predicted global protein motion [49]. This distinguishes ENCoM-based NMA from a mere geometric approximation of flexibility, allowing ligand-dependent differences in the binding cavity to propagate into the predicted low-frequency mode ensemble of the full β3 subunit.
Each of the five re-docked ligand-β3 complexes (X47-A—E, as generated via mcule.com’s 1-Click Docking, 2.2) were submitted as a static input structure to DynaMut, from which two per-residue outputs were retained for comparative, Elastic Network Contact Model (ENCoM)-based analysis across the five leads: a) deformation energy, representing the local energetic cost of displacing a given residue from its equilibrium position, with higher values indicating greater local rigidity, and b) atomic fluctuation (Å), representing the predicted amplitude of residue motion under the same conditions. In parallel, DCCM was computed for each complex across all Cα-Cα amino acid residue pairs, generating a residue-by-residue matrix of correlation coefficients (−1 to +1), wherein positive values denote residue pairs predicted to move in-phase (correlated), negative values denote out-of-phase (anticorrelated) motion, and near-zero values denote dynamically independent motion. Given the established role of coordinated conformational transitions in β3-mediated hantavirus attachment—specifically, the calcium-driven bent-to-active transitions that modulate PSI domain (Gly1001—Phe1056) exposure at the integrin apex [15] [19]—DCCM analysis was used specifically to assess whether ligand binding at the orthosteric cavity was predicted to correlate with coordinated motion at this distal, hantavirus-binding domain, offering an indirect readout of potential allosteric communication thereof.
2.5. ADME-Tox Computation
Physicochemical, pharmacokinetic, and drug-likeness properties of the lead compounds were computed using SwissADME [50], a freely accessible web tool that derives a panel of descriptors directly from a molecule’s structure without requiring physical samples. In addition to core physicochemical descriptors such as molecular weight (MW), number of rotatable bonds (NRotB), topological polar surface area (TPSA), hydrogen bond donor (HBD) and acceptor (HBA) counts, and consensus octanol-water partition coefficient (CLogPo/w, averaged across five separate lipophilicity prediction methods), SwissADME was used to assess Lipinski’s Rule of Five (RO5)/Veber’s rule (ROV) compliance, P-glycoprotein (P-gp) substrate status, and molar refractivity (MRef) for each lead. Gastrointestinal absorption and blood-brain barrier permeation were additionally evaluated using SwissADME’s Brain Or IntestinaL EstimateD permeation (BOILED-Egg) model, a classification tool that plots each compound within a two-dimensional space defined by Wildman-Crippen lipophilicity (WLOGP) (y axis) and apparent polar surfaces (TPSA) (x axis): compounds falling within the model’s white region are predicted to be well absorbed by the gastrointestinal tract, while those falling within the yolk (yellow) region are predicted to permeate the blood-brain barrier. Because the cyclophosphinamide leads under evaluation are structurally large, polar, multi-ring peptidomimetics, RO5 compliance was not treated as a strict pass/fail druggability filter; rather, RO5/ROV violation counts were interpreted in the context of the well-established “beyond Rule of Five” (bRO5) chemical space, within which macrocyclic and peptidomimetic drugs routinely and successfully violate one or more classical RO5 thresholds, while retaining oral bioavailability through alternative mechanisms such as conformation-assisted intramolecular hydrogen bonding [51].
Physicochemical descriptors obtained from SwissADME were cross-verified using mcule.com’s “Property Checker” tool, and toxicological properties of the lead compounds were independently assessed using ProTox 3.0 [52]-[54], a webserver that combines molecular similarity analysis, pharmacophore-based methods, and machine learning classification models to predict toxicity across a broad panel of endpoints including: acute oral toxicity (PO LD50), organ-specific toxicity (hepatotoxicity, cardiotoxicity, nephrotoxicity, immunotoxicity, neurotoxicity, and respiratory toxicity), carcinogenicity, and relevant Tox21 stress-response/nuclear-receptor pathway activity for a given input structure.
Finally, phase I (oxido-reductive/hydrolytic) metabolic transformation of the lead compounds was predicted using BioTransformer 3.0 [55] [56], a metabolism-prediction platform that combines a knowledge-based system of literature-derived biotransformation rules with machine learning classification models (including a random forest-based module for predicting cytochrome P450 substrate selectivity) to generate metabolite structures for a given parent compound. For this study, the phase I prediction module of BioTransformer 3.0 was applied to each lead compound to generate hepatic metabolite structures, and their corresponding predicted biotransformation pathways.
2.6. Limitations
Several limitations of the present study were to be acknowledged. First—and most fundamentally—this study proceeds from the long-standing pathophysiological premise that αII/αvβ3 integrins serve as the principal hantavirus entry receptor via PSI-domain engagement: a premise that has been challenged by recent genetic depletion studies implicating protocadherin-1 (PCDH1), rather than β3 integrin, as the cardinalmost entry determinant for several New World (though not Old World) pathogenic hantavirus strains in human endothelial cells [57] [58]. Consequently, the in silico leads reported here should be interpreted as candidate β3-binding inhibitors validated at the structural/computational level only, and their translation to demonstrable antiviral efficacy would require experimental confirmation of receptor dependence in the relevant viral strain and cell lines in vitro prior to further conclusions.
Second, all binding-affinity estimates reported in this study derive from empirical, force field-based docking scoring functions (AutoDock Vina, via SwissDock and mcule.com’s 1-Click Docking), which—in our assessment—only partially capture the true, in vitro thermodynamics of protein-ligand binding. Docking scores reported throughout this study should therefore be interpreted as relative, comparative rankings among structurally-related analogues, rather than as quantitative predictors of true binding affinity in absolute terms.
Third, as detailed in Section 2.3, molecular dynamics simulation (MDsim) was restricted to the isolated β3 subunit rather than the complete αvβ3 heterodimer, owing to parameterization and computational constraints associated with simulating the full, multidomain integrin on the SiBioLead platform. As a result, the MDsim-derived stability assessments reported in this study characterize the ligand receptor complex specifically within the isolated ligand-proximal β3 subunit, and do not capture potential allosteric or quaternary-context effects that the αv subunit, or the full bent/hybrid/open conformational transitions of the intact heterodimer, might exert on binding-site dynamics in vitro/in vivo. Additionally, each system was simulated as a single 50 ns trajectory; while sufficient to assess short-timescale local stability of the docked pose, in our assessment this simulation length and lack of independent runs may not fully capture rarer conformational transitions or long-timescale binding-site dynamics.
Moreover, it should be noted as a target-liability limitation that—due to the inability of this study on its own to establish selectivity for αvβ3 or αIIβ3 over the other, potential cross-reactivity should be noted, with the added consequences of unforeseen off-site effects on platelet aggregation and hemodynamics remaining a notable topic of contention for such an in silico-only study.
3. Results and Discussion
3.1. Blind Docking and Ligand-Based Virtual Screening (LBVS)
The initial blind docking via SwissDock of G319-0078 to αvβ3 (PDB ID: 1U8C) in order to re-visualize and discover the coordinates of the peptidomimetic binding site yielded a distinct binding cavity at coordinates X: 27.502, Y: 22.528, Z: 49.248, with a docking score of −6.9 kcal/mol, matching the binding site depicted by Hall et al. [23] [38] [39]. Following re-docking of the G319-0078-αvβ3 complex with ADVina, L2-40 were designed rationally, and produced varying top-scoring docking scores (kcal/mol) relative to G319-0078 (Figure 7).
All of the initial hits were comparable to G319-0078 in terms of docking score, with no hit docking with a difference of more than 0.7 kcal/mol from the original ligand. As depicted in Figure 7, G319-0078 docked with a score of −7.3 kcal/mol, with 26 hits (L2, L3, L5-10, L12, L14, L16-22, L24-28, L33, L37, L39, and L40) showing comparable or weaker predicted binding affinity via a smaller (more positive) docking score. On the other hand, 13 hits (L4, L11, L13, L15, L23, L29-32, L34-36, and L38) showed a stronger predicted binding affinity to αvβ3 than the original ligand, depicted via larger (more negative) docking scores. Of these, L34—characterized by a cinnoline-7-sulfonamide and benzylic 6H-cyclohepta[c]pyridine—had the strongest predicted binding affinity to αvβ3 (−8.0 kcal/mol, −0.7 kcal/mol than G319-0078). The two analogues with a selenonylamide (L29) and phosphonamide (L30) substitution of the central sulfonamide in G319-0078, alongside a benzyl-to-indole substitution, likewise depicted larger (more negative) docking scores of −7.8 kcal/mol; a 0.5 kcal/mol increase in predicted binding affinity compared to G319-0078.
Overall, the disubstituted ligands (L29-40) were found to exert stronger predicted binding affinity to the orthosteric binding site of αvβ3 than the monosubstituted derivatives (L2-28). Remarkably, heterocyclic benzylic moiety substitutions did not show a clear pattern of preference for a specific heteroatom over another in the same isosteric class, however certain trends could be observed, such as meta- and para-substitutions yielding larger docking scores in fluoro- and chlorobenzyl, and the para-trifluoromethylbenzyl substituent increasing docking scores up to −7.8 kcal/mol (L15). On the other hand, analogous aliphatic substitutions resulted in a markedly decreased predicted binding affinity (i.e. L10, −6.2 kcal/mol) compared to G319-0078.
3.2. Computational Lead Optimization Strategy
Further in silico lead optimization was conducted based on bioisosteric modifications of L34 and L30, given the rationale described for the incorporation of
Figure 7. A chart depicting G319-0078, alongside L2-40—hit compounds docked at the orthosteric site. Note the monosubstituted (blue) bioisosteric modifications (L2-28) on the benzyl of G319-0078, whilst disubstituted ligands (L29-40) contain both benzylic (blue) and sulfonamidic (red) modifications. The 1,3-benzodioxole and alkyl amide backbone have remained unchanged.
cyclophosphinamides (see Introduction). Notably, it was discovered that incorporating orthogonal hydrogen bond donors and acceptors (HBD/HBA), such as in phthalimide, at an N-alkyl position increased predicted binding affinity of the ligand to αvβ3, whilst removing nitrogen-containing heterocyclic groups at the N-alkyl position had the converse effect. Moreover, the heterocyclic moiety had to be of a heteroaryl class, as aliphatic heterocycles likewise decreased predicted binding affinity, much as was the case during initial docking studies (3.1).
It was likewise discovered that cyclizing the alkyl tether scaffold via an aryl phosphinamide functional group increased predicted binding affinity by up to 2 kcal/mol at the designated binding site, with the top four docking scores of five in silico lead compounds (X47-A—E) depicted below (Table 1).
Table 1. A table of the top four docking scores (kcal/mol) of five in silico lead ligands (X47-A, X47-B, X47-C, X47-D, X47-E) with αvβ3 (PDB ID: 1U8C), computed via mcule.com’s AutoDock Vina software: 1ClickDocking [42].
Ligand |
Pose 1 |
Pose 2 |
Pose 3 |
Pose 4 |
Average Docking Score (kcal/mol) |
X47-A |
−8.0 |
−7.9 |
−7.8 |
−7.6 |
−7.83 |
X47-B |
−8.2 |
−8.1 |
−7.3 |
−7.2 |
−7.70 |
X47-C |
−8.8 |
−7.7 |
−7.6 |
−7.2 |
−7.83 |
X47-D |
−9.0 |
−8.6 |
−8.3 |
−7.7 |
−8.40 |
X47-E |
−8.5 |
−8.1 |
−8.0 |
−7.5 |
−8.03 |
The lead ligand with the largest top docking score (-9.0 kcal/mol) was X47-D, whilst the smallest top docking score was observed with ligand X47-A (−8.0 kcal/mol). Nevertheless, X47-B retained the smallest average docking score, at -7.70 kcal/mol, indicating the lowest average predicted binding affinity to αvβ3 at its top-four-scoring poses. Conversely, X47-D reported an average docking score of −8.40 kcal/mol, making it the ligand with the strongest affinity to the orthosteric binding site of G319-0078 on αvβ3. Notably, a ligand which binds to a receptor/enzyme with an ADVina docking score larger (more negative) than −6.0 kcal/mol is considered, according to Shityakov and Förster [59], to render a compound “active” in binding its target, thus implicating all of the five aforementioned leads in viable hypothetical druggability of αvβ3.
It was observed, furthermore, that the primary structure-activity relationship (SAR) interactions responsible for the favorable predicted binding affinity of the leads consisted of parallel displaced π-stacking between the phthalimide moiety and Trp1129 (i.e. X47-D), or the latter and the variable benzylic groups attached to the carbonyl-adjacent alpha carbon of the central ring (Figure 8). Additional SAR factors included cation-π interactions between the terminal ammonium of Lys1125 at the phosphinamidic heteroaryl moiety. Moreover, it was observed that transient H-bonding interactions of the HBA/HBD-orthogonal phthalimide moiety with correspondent, proximal HBAs/HBDs on either the backbone (i.e. amide of Lys1208) or side chain (i.e. OH of Ser1211) could play a role in stabilizing the ligand-receptor complex further, contributing to the larger docking scores of the phthalimide-bearing leads as compared to the initial hits. Configuring overall charge from cationic [+1] (X47-A—C) to neutral (X47-D/E), however, produced no specific correlation with increased predicted binding affinity of the latter compounds per se, nor remarkable SAR interactions.
![]()
Figure 8. A ribbon model of X47-D (a) binding and interacting with the β3 subunit of αvβ3, with overt cation-π interactions (Lys1125), putative π-π stacking (Trp1129), and H-bonding (Ser1211), and a subsequent model of X47-B; (b) engaging with the orthosteric binding site via isotypical side-chain interactions (with the exception of Lys1208, instead of Ser1211, H-bonding with phthalimide moiety).
Additionally, it was discovered that inserting a cyclic hydroxamic moiety at the benzylic position (as in X47C—E)—as opposed to a carbamate—notably increased predicted binding affinity (Table 1). This was speculated to be due to the increased HBD characteristic of hydroxamic acids interacting with hydrophilic side chains in the region/their backbone amides, which was confirmed via additional distance computation of the hydroxamic hydroxyl and the nearest HBA in the X47-C/αvβ3 complex (Figure 9).
The computed distances between the hydroxamic hydroxyl and the nearest HBAs (amidic oxygen of Lys1208 and Val1207) spanned 3.0 Å to 3.4 Å, placing the ligand within the feasible H-bonding range of <4 Å [60] (with the optimal range being 2.6 - 3.1 Å [61] [62]), corroborating our hypothesis for favoring terminal HBD functionalizations at the benzylic position over carbamates. In the above way, the modification carried by introducing cyclic hydroxamic acid moieties pertained to a ring-contraction from a six- to five-membered system, whilst still maintaining heteroaryl character of the ring, which likewise correlated with an improvement in predicted binding affinity of the lead ligands to αvβ3, speculated to be—in addition to decreased steric hindrance allowing further cavity penetration—due to added H-bonding with the proximal, flexible amine side chain of Lys1208 (Figure 9). This is supported by a study by Gallivan and Dougherty [63], who found that aromatic moieties with directly-attached/embedded HBD functional groups (i.e. tyrosine phenol, and tryptophan amine) markedly strengthen cation-π interactions to create an empirically observable difference compared to HBD-lacking aromatic systems (i.e. phenylalanine benzene, or X47-A/B with solely tertiary carbamates at the benzylic position) [64].
![]()
Figure 9. A ribbon model depicting the distances between the hydroxyl of the cyclic hydroxamic acid moiety of X47-C (yellow) and the amide oxygen on the backbones of Lys1208 (3.4 Å) and Val1207 (3.0 Å) of αvβ3 (cyan) post-docking, computed and visualized with Schrodinger’s PyMOL.
Nevertheless, in order to examine in situ stability computationally, further molecular dynamics simulation (MDsim) was required to understand receptor-ligand interactions over time, an expansive solvent/salt environment, and temperature mimicking that of the human physiological system itself.
3.3. Molecular Dynamics Simulation (MDsim)
Preliminary GROMACS MDsim runs of X47-A—E yielded variable—though at times consistent—results of RMSD and RMSF (Figure 10).
Over the duration of the simulation (50 ns), RMSD tended to yield an average of ~0.6 nm over the five lead ligands, with the apex deviation observed ~0.9 nm with ligand X47-C at the 32.5 ns mark (Figure 10(c)). Correspondingly, only X47-B showed a stable plateau (~0.65 nm) across the simulation duration, reached after 7 ns, while X47-C reported the second-most stabilizing RMSD plateau (~0.8 nm), reached after approximately 25 ns, respectively. On the other hand, X47-A
Figure 10. RMSD (left) and RMSF (right) graphs of lead ligands X47-A (a), X47-B (b), X47-C (c), X47-D (d), and X47-E (e), as computed via SiBioLead’s GROMACS MDsim software [45].
was observed with impermanent, ternary plateaus ranging from durations: 1 - 15 ns (~0.5 nm), 17 - 33 ns (~0.6 nm), and 40 - 50 ns (~0.7 nm). Similarly, X47-D and X47-E did not produce meaningful plateauing in the given simulation duration, indicating unstable in silico binding to the solvated target receptor subunit. Considering optimal RMSD values should plateau below 0.2 nm (<2 Å) [65] [66], none of the ligands show favorability in maintaining a stable protein conformation during MDsim. Nevertheless, the stable plateauing of X47-B and X47-C at the aforementioned deviation distances over the duration of MDsim indicated that the β3 subunit underwent a major—albeit nonetheless stabilizing—conformational change [67].
The RMSF values, moreover, collated in Figure 10 allowed for the inference that residues ~Leu1185 consistently showed among the largest fluctuations in all five leads, corresponding to a loop motif proximal to the orthosteric binding site (16.2 Å), with the largest value observed being 0.9 nm (X47-D). Concurrently, residues from to Ser1020 to Phe1056 consistently showed higher than average fluctuation orthogonal to the binding site, indicating instability at the distal, hantavirus-binding PSI domain for all but X47-C (<0.4 nm fluctuation, aside from the 0.7 nm fluctuation at ~Lys1125 binding site). The helix-loop domain at residues Cys1523-Asp1692, corresponding to helix-loop domain which is—in the heterodimer—stabilized by the neighboring αv subunit—was not considered due to the inherent thermodynamic instability of said domain in the isolated β3 subunit as such (Figure 11).
![]()
Figure 11. A ribbon-Gaussian volume diagram of the β3 subunit of αvβ3, depicting the high-RMSF Cys-1523-Asp1692 (red) domain, binding cavity (Met1180-Tyr1190, encircled violet and labeled), and PSI domain loop motif (Ser1020-Phe1056, blue) [PDB ID: 1U8C].
On average, RMSF values of the main protein structure were lowest with X47-B (~0.15 nm), whilst they were the highest with X47-C (~0.2 nm). Nevertheless, X47-A and X47-C had the lowest fluctuation variation among the leads between their peaks and troughs, with their average RMSF values of ~0.2 nm putting them within the range of stabilizing ligand-protein interaction in MDsim (<0.4 nm/4 Å [68]). Solvent-accessible surface area (SASA) computation during the MDsim runs, furthermore, indicated that—while all lead ligands demonstrated an decrease from the beginning value of 330 nm2—the greatest decrease was observed with X47-A and X47-D (330 to 290 nm2); the smallest decrease, conversely, was observed with X47-E (330 to 310 nm2) (Figure 12).
Figure 12. Solvent accessible surface area (SASA, nm2) of the five lead ligands (X47-A—E), computed via SiBioLead’s GROMACS MDsim software [45].
The most stable, linear decrease in SASA, however, was observed with X47-C and X47-D, whilst the other ligands subscribed to vast irregularities in SASA variation over the 50 ns MDsim duration. Thus, it was inferred that X47-C and X47-D were most effective at burying hydrophobic residues from the water solvent cage, and stabilizing the protein in a new conformation in silico.
3.4. Normal Mode Analysis (NMA) and Dynamic Cross-Correlation Matrix (DCCM)
Further normal mode analysis (NMA) and dynamic cross-correlation matrix (DCCM) computation via DynaMut [49] yielded varying results amongst the five lead ligands (Figure 13).
Figure 13. Normal mode analysis (NMA) of the β3 subunit of αvβ3 in complex with ligands: X47-A, X47-B, X47-C, X47-D, and X47-E, depicting deformation energies (a) and atomic fluctuations (b). Computed via DynaMut [49] and visualized via Google Sheets.
As depicted in Figure 13(a), deformation energies reached an apex of ~27.5 at ~Gln1440 (X47-C), corresponding to a loop segment connecting the αv-stabilized Cys1523-Asp1692 domain of the β3 subunit (Figure 11). Atomic fluctuations (Figure 13(b)), on the other hand, tended to correspond across all five leads at the aforementioned residue (maximum ~2.2 Å, X47-C), with a notable further fluctuation increase at the ~Met1180 proximal to the orthosteric binding site (maximum ~1.4 Å, X47-C), and from ~Ser1020-Phe1056 (PSI domain loop, maximum ~2.2, X47-A). While deformation energies, indicative of the energy expense to reposition a specific protein residue from its reference structure, remained highest at a set value for X47-C, the largest average overall deformation energy across residues was deduced for X47-D, indicating the least local flexibility and highest protein rigidity with the latter ligand-protein complex [49]. This correlates with the observed SASA trendline of the X47-D—β3 complex (Figure 12), with its most linear, stable decrease over the MDsim duration compared to the other leads. Moreover, atomic fluctuations likewise favored X47-D as the most stabilizing of ligands due to its lower fluctuation values compared to the others.
Additional dynamic cross-correlation matrix (DCCM) computation detailed the strongest positive correlation scores (darkest red) between residues Gly1001-Ala1050 (res. 1 - 50), corresponding to part of the hantavirus-binding PSI domain, most notably with ligand X47-A, while additional strong positive correlation was observed around residues Gln1440 (res. 440), most notably with X47-E (Figure 14). Further positive correlation is observed at the terminal residues corresponding to the Cys1523-Asp1692 domain (Figure 11) stabilized by the αv subunit in vivo, most notably with X47-C.
Interestingly, X47-C depicted the least negative correlation amongst residues (especially Lys1350-Asp1550) compared to the other four ligand-protein complexes, indicating low oppositional residue movement [49]. Nevertheless, it was deduced that the propensity for all five ligands throughout the NMA/DCCM simulation to invoke a positive correlation at the hantavirus-binding PSI domain may indicate conformational change that could play a hypothetical role in terms of viral inhibitory feasibility, with X47-A and X47-C implicated most therein in the aforementioned change.
3.5. ADME-Tox Computation
Further ADME-Tox computation conducted by SwissADME revealed that leads X47-A—E, while having variable Brain or IntestinaL EstimateD permeation (BOILED-Egg), Wildman-Crippen water partition coefficient (WLogP), and topological polar surface area (TPSA, Å2) values, invariably showed sub-optimal intestinal absorption and no blood-brain barrier (BBB) permeation (Figure 15). Moreover, all of the lead ligands were found to be P-glycoprotein (P-gp) substrates, implicating them in further reduced intestinal absorption, BBB permeation, and overall bioavailability [69].
In the above regard, X47-C, X47-D, and X47-E were found to have too large a
Figure 14. Dynamic cross correlation matrices (DCCM) of the β3 subunit of αvβ3 in complex with ligands: X47-A (a), X47-B (b), X47-C (c), X47-D (d), and X47-E (e), computed via DynaMut [49]. Red indicates positive correlation between residue movement, whilst blue indicates negative correlation thereof. Note that residue no. “1” corresponds to authentic residue no. “1001”, for analytical purposes.
Figure 15. A Brain or IntestinaL EstimateD permeation (BOILED-Egg) model of the five lead compounds (X47-A to E), with the white region corresponding to high intestinal absorption (HIA), yellow “yolk” region to the zone of blood brain barrier (BBB) permeability, and the blue color of the X47 compound bank representing positive P-glycoprotein (P-gp) substrate status (computed via SwissADME [50]). The x-axis consists of topological polar surface area (TPSA, Å2) while the y-axis consists of fragment-based Wildman-Crippen partition coefficient (WLogP) values.
TPSA (162.70 Å2, 182.67 Å2, and 169.78 Å2, respectively) to achieve effective intestinal absorption, or permeate the BBB, with X47-A and X47-B positioned at the borderline of HIA (TPSA: 140.50 Å2 and 138.81 Å2, respectively) (Table 2). Furthermore, all leads maintained unfavorable WLogP values for either HIA or BBB permeation, ranging from ~2.8 (X47-C) to ~4.5 (X47-E).
Additional factoring of Lipinski’s “Rule of Five” (RO5)—which states that effective bioavailability via administration per orum (PO) entails no more than one violation of: <5 consensus LogP (CLogPo/w), <500 Da molecular weight (MW), <10 HBAs, and <5 HBDs—as well as Veber’s rules (ROV, where good PO bioavailability coincides with: TPSA < 140 Å2, and rotatable bond number (NRotB) ≤ 10), however, provided variable data [70] [71]. Namely, X47-B was found to have only one RO5/ROV violation—that being its MW > 500 Da (672.67 Da)—making it the only lead displaying good predicted RO5/ROV bioavailability; other ligands contained two (X47-A) or three (X47-C—E) RO5/ROV violations, rendering them suboptimal for PO administration according to the aforementioned criteria. Nonetheless, when considering that molar refractivity (MRef) values for optimal PO bioavailability range between 40 to 130 cm3∙mol−1, none of the lead ligands score favorably due to their high polarizability and steric bulk [72].
Table 2. A table depicting the rotatable bond number (NRotB), total polar surface area (TPSA, Å2), number of hydrogen bond donors (HBD) and hydrogen bond acceptors (HBA), consensus octanol-water partition coefficient (CLogPo/w)1, P-glycoprotein (P-gp) substrate status, Lipinski’s “Rule of Five” (RO5) and Veber’s rule (ROV) violations, and molecular refractivity (MRef) of X47-A to E leads. These results were obtained from SwissADME [50], and verified by ProTox-3.0 [52] and mcule.com’s “Property Checker” software [42].
Compound |
NRotB |
TPSA (Å2) |
HBD |
HBA |
CLogPo/w |
P-gp substrate |
RO5/ROV Violations |
MRef (cm3∙mol−1) |
X47-A |
9 |
140.50 |
2 |
8 |
2.86 |
Yes |
2 (MW > 500, TPSA > 140) |
194.51 |
X47-B |
10 |
138.81 |
2 |
8 |
3.21 |
Yes |
1 (MW > 500) |
191.96 |
X47-C |
8 |
162.70 |
3 |
11 |
2.05 |
Yes |
3 (MW > 500, HBA > 10, TPSA > 140) |
181.22 |
X47-D |
8 |
182.67 |
4 |
13 |
2.86 |
Yes |
3 (MW > 500, HBA > 10, TPSA > 140) |
178.45 |
X47-E |
8 |
169.78 |
4 |
12 |
3.57 |
Yes |
3 (MW > 500, HBA > 10, TPSA > 140) |
180.66 |
With regard to toxicological analyses of the five lead compounds, X47-B, X47-C, and X47-D displayed markedly low predicted oral acute toxicity, measured by high median lethal dose (LD50) values of 4000 mg/kg (X47-B) and 6000 mg/kg (X47-C & X47-D). On the other hand, X47-A and X47-E displayed lower oral LD50 values of 420 mg/kg and 502 mg/kg, respectively, indicating greater oral acute toxicity. Nevertheless, all compounds had a toxicity class of >4 (>IV), indicating an absence of “toxic”/“fatal” classification PO, and a mere “harmful” status when administered orally [52]. Additionally, all compounds displayed predicted inactivity against a panel of nuclear receptors [aryl hydrocarbon receptor (AhR), androgen receptor (AR), androgen receptor ligand-binding domain (AR-LBD), aromatase, estrogen receptor alpha (ERα), and peroxisome proliferator-activated receptor gamma (PPARγ)], indicating a lack of predicted (xeno)androgenic/(xeno)estrogenic effects leading to reproductive issues, hormone-sensitive cancers (i.e. prostate and breast cancer), and metabolic syndrome and insulin resistance (i.e. via added PPARγ antagonism) [73] [74]. Moreover, none of the lead ligands displayed predicted cardiotoxicity, carcinogenicity, nephrotoxicity, hepatotoxicity, or neurotoxicity, although all leads displayed a high probability for immunotoxicity (>0.92) and respiratory toxicity (>0.78); of these, the highest probability for respiratory toxicity was found in X47-D (0.84), while the highest probability for immunotoxicity was found in X47-E (0.99) [52] [57] [58].
Further analyses of phase I metabolism by the cytochrome P450 (CYP) enzyme superfamily indicated a specific predicted preference for the CYP1A2 subtype shared by all five leads. Namely, the primary phase I biotransformations—as computed by BioTransformer 3.0 [53] [54]—by the aforementioned enzyme uniform across the lead ligands consisted of N-dealkylation of the phthalimide moiety, aromatic/benzylic hydroxylation, and terminal aliphatic hydroxylation (Figure 16).
Of notable concern to potential in vitro/in vivo toxicity is the formation of 1,3-dioxoisoindole-5-carbaldehyde consequent to phase I biotransformation (via N-dealkylation, Figure 16). Following additional ProTox 3.0 predicted toxicity
Figure 16. Three computed phase I (oxido-reductive/hydrolytic) metabolic pathways of the X47 class, here represented by X47-A, including: N-dealkylation of the phthalimide group, aromatic/benzylic hydroxylation, and terminal aliphatic hydroxylation at the carbamic alkyl chain. The pathways were computed with BioTransformer 3.0 [53] [54].
profiling, the 1,3-dioxoisoindole-5-carbaldehyde metabolite—while maintaining a predicted oral LD50 value of (Class V toxicity)—showed putative binding to flavin-containing monoamine oxidase A (MAO-A) responsible for the breakdown of indoleamine (i.e. serotonin) and catecholamine (i.e. dopamine and norepinephrine) neurotransmitters in the central nervous (CNS) and cardiovascular systems [75]. This, coupled with the predicted ability of 1,3-dioxoisoindole-5-carbaldehyde to readily cross the BBB (probability = 0.93), poses potential off site effects hypothesized to pose notable clinical risk of sympathomimetic toxidrome [52] [76].
Further computational analysis showed that the free pyridine at the aryl phosphinate moiety of X47-C undergoes a predicted pyridine N-oxidation into a pyridine N-oxide species (Figure 17). While studies have demonstrated that certain pyridine N-oxides, such as those which are chloropyridine-derived, may invoke a dose-dependent cytotoxicity/clastogenicity, 2,6-dimethyl pyridine (2,6-lutidine) N-oxides were found to have antimutagenic properties in vitro via Ames test [77] [78]. Hence, a definitive implication for X47-C’s phase I biotransformation into the corresponding N-oxide cannot be ascertained from in silico studies alone.
Figure 17. The additional computed phase I pathway of the “free pyridine” lead, X47-C, metabolized into pyridine N-oxide, computed with BioTransformer 3.0 [53] [54].
With the above in silico prediction methods considered (LBVS, docking, MDsim, NMA/DCCM, and ADME-Tox computation), it was concluded that the most promising lead drug candidate appeared to be X47-C out of the five presented ligands (X47-A, X47-B, X47-C, X47-D, and X47-E). This was deduced based across its second-largest docking score of −8.8 kcal/mol (0.2 units below -9.0 kcal/mol of top-scoring X47-D), an average docking score of −7.83 kcal/mol across the four top-scoring poses (rivaled by X47-A), its stable RMSD plateauing (~0.8 nm) within ~25 ns, coupled with its stable, linearly-decreasing SASA (330 - 305 nm2), and low-variation RMSF (~0.2 nm) compared to other leads (with the exception of X47-A). NMA and DCCM, moreover, supported high deformation energies, low atomic fluctuations, and dense positive dynamic cross-correlation, particularly at the hantavirus-binding PSI domain. Further ADME-Tox analysis likewise provided optimal predicted PO LD50 toxicological data for X47-C (6000 mg/kg) compared to the other in silico leads (with the exception of X47-D in terms of absolute LD50 values), with a further lower predicted risk of respiratory toxicity and immunotoxicity compared to X47-D and X47-E, respectively. This factor remains paramount due to the nature of hantavirus infection, particularly the respiratory infection of HPS, and the immune response that follows suit which frequently leads to enhanced vascular permeability—and worsened clinical outcomes—due to a “cytokine storm” [4].
While a valid concern about the above conclusion pertains to X47-C’s notable three RO5/ROV violations, rendering it of low PO bioavailability by classical pharmacokinetic standards, the introduction of concepts such as the “beyond Rule of Five” (bRO5) to account for the high bioavailability of natural products and peptides that traditionally violate RO5/ROV standards have provided additional paradigms to view pharmacokinetic effects of X47-C and its analogues [51]. Notably, transient intramolecular hydrogen bonding of the phosphinamide backbone with other HBA/HBD residues may reduce the overall polar surface area (TPSA) of the molecule enough to render optimal intestinal epithelium permeation, while macrocyclization of the phosphinamide—as demonstrated with X47-C and the other leads—reduces NRotB, much to the same effect on bioavailability in theory [51].
Although X47-C and its analogues have been rationally-designed keeping metabolic stability of phosphinamide peptidomimetics in mind, it still remains a future task to prove in vivo metabolic stability of the compounds in question, in addition to in vitro assays to evaluate cytotoxicity (i.e. MTT) and antiviral inhibitory activity against a large panel of hantaviruses (Old and New World alike) which bind to β3 integrins collectively (i.e. αv/αII). The fact that the above lead compounds bind to the conserved β3 subunit, as opposed to the variable α subunit, renders—in our opinion—further in vitro research justification for cross-hantavirus inhibitory activity all the more relevant for potential clinical candidacy, alongside screening for off-site effects involving cardiovascular and respiratory systems due to the nature of β3 integrin binding in question.
4. Conclusion
With the above analyses, we have concluded that cyclophosphinamide lead X47-C presented the most optimal in silico candidate for hantavirus inhibition at the β3 subunit of αvβ3 integrin. Building on the pioneering cyclopeptide and peptidomimetic research conducted by Hall et al. [22] [23], X47-C has consistently scored well across molecular docking (−8.8 kcal/mol top-scoring pose, −7.83 kcal/mol average across four poses), RMSD (plateau at ~0.8 nm, 25 ns), RMSF (average ~0.2 nm, least variation), steady and linear SASA decrease over 25 nm2 (330 - 305 nm2), favorable NMA and DCCM outcomes, and toxicological analyses. Further research directions should point to viable synthetic strategies, coupled with in vitro and in vivo testing to empirically validate the inhibitory potential and predicted mechanism of action (MOA) of X47-C, principally due to the wide cross-hantavirus inhibition theorized due to the conservation of the β3 subunit of New World hantavirus-binding (αvβ3) and Old World hantavirus-binding (αIIβ3) integrins.
Data Availability
The data supporting the findings of this study are available from the corresponding repository at DOI: https://doi.org/10.5281/zenodo.22276619.
Acknowledgements
The authors would like to acknowledge the facilitation of this study by the Future Generation of Scientists (FGS) team, without whose resources this project would not have been possible.
Author Contributions
Filip Latkovic: Writing—Original Draft, Writing—Review & Editing, Conceptualization, Methodology, Investigation, Visualization, Supervision, Marlowe Samuel: Writing—Original Draft, Investigation, Formal Analysis, Mary Amgad Adly Wasef Haroon: Investigation, Formal Analysis, Validation, Abdelrahman Shreef Soubi Ragab: Investigation, Formal Analysis, Zsofia Szegedi: Resources, Project Administration, Sviatoslav Gladyshev: Resources.
NOTES
1Across fragment-based and atomistic calculations: WLogP, XLogP3, iLogP, WLogP, and SILICOS-IT.