Computational investigation of covalent organic pyrgos[n]cages as potential dual modulators of the VEGFR-2 D2–D3 domain and DNA structural interfaces

Abstract

Angiogenesis, driven by vascular endothelial growth factor (VEGF) and its receptor VEGFR-2 domain D2–D3 (VEGFR-2 D23) signaling, plays a central role in tumor vascularization and progression. In this multiscale computational study, we investigate pyrgos[n]cages (PC[n], n = 1–4), purely covalent organic cages with electronically adaptive π-stacked systems, capable of engaging the extracellular VEGFR-2 D23 ligand-binding interface (wild type and mutant) and DNA (duplex and fork) structures. Molecular dynamics simulations and electronic structure analyses indicate that the electronically soft and polarizable PC[2–4] frameworks preferentially localize within the VEGFR-2 D23 interdomain cleft, stabilized by hydrogen bonding, π⋯π stacking, and π⋯H interactions. This multivalent interaction pattern suggests potential steric interference with VEGF binding while preserving overall receptor structural stability, including mutant VEGFR-2 D23 variants. Complementary DNA-binding analyses reveal structurally stable groove-associated binding stabilized by directional hydrogen bonds, phosphate backbone interactions, π⋯π stacking, and van der Waals contacts. Among the studied systems, PC[4] exhibits the most favorable binding enthalpy estimates for VEGFR-2 D23 (−48.3 ± 3.5 kcal mol⁻¹), mutant M-VEGFR-2 D23 (−24.3 ± 5.8 kcal mol⁻¹), and duplex DNA (−44.8 ± 3.8 kcal mol⁻¹). At the same time, PC[2] shows favorable interaction with fork DNA (−46.3 ± 4.5 kcal mol⁻¹). These findings provide mechanistic insight into the organic cage PC[n] in engaging biologically relevant protein-DNA interfaces and highlight its potential to modulate the VEGFR-2 D23 domain and DNA binding.

---

1. Introduction

Angiogenesis is crucial for the growth of new blood vessels from existing ones¹ and is essential for development, tissue growth, healing, and reproduction.^2,3 While they are vital for normal body functions and embryonic development, their disruptions lead to various health issues.⁴ Abnormal angiogenesis allows tumours to secure a steady supply of oxygen and nutrients, which promotes their growth and can lead to the spread of cancer through the lymphatic system.⁵ Understanding the molecular factors that control blood vessel growth is, therefore, a key focus in the fields of vascular biology and translational medicine.⁶ The major contributing factor is the vascular endothelial growth factor VEGF-A and its associated receptor VEGFR-2.^7,8 Pioneering domain deletion studies established that VEGFR-2 requires the cooperative action of domains 2 and 3 (D23) for high-affinity VEGF-A recognition. Additionally, VEGFR-2 D23 binds VEGF-C/D/E, which have direct and indirect effects on angiogenesis.^9,10 Traditional small molecules are often ineffective against such shallow protein–protein interaction (PPI) interfaces, demanding the need for biologics, including monoclonal antibodies, engineered proteins, and nucleic acid aptamers, that can mimic or disrupt extended ligand-binding epitopes.^11–13 Consequently, designing molecular scaffolds that structurally resemble VEGF-like supramolecular architectures offers a promising strategy to interfere with VEGF/VEGFR-2 D23 interactions, thereby blocking downstream signalling and selectively restricting tumour-associated angiogenesis.

Supramolecular scaffolds, including macrocycles, calixarenes, cucurbiturils (CB), and cyclodextrin cages, have long been explored for biomedical applications, including antimicrobial agents, drug delivery vehicles, and modulators of PPIs.^14–17 These systems excel in molecular encapsulation and host–guest recognition. Li et al. demonstrated the effectiveness of CB[7] as an inhibitor of MPTP/MPP⁺-induced neurodegeneration in zebrafish, suggesting a promising strategy for preventing toxin-induced Parkinson's disease and expanding the potential biological and biomedical applications of CB[7].¹⁸ Dockerill et al. developed a supramolecular drug strategy in which peptide nucleic acid (PNA) hybridized drug fragments enable potent yet fully reversible activity.¹⁹ Their bivalent thrombin inhibitors showed strong anticoagulant effects that could be rapidly neutralized by a complementary PNA antidote. However, their reliance on flexible linkages or polymeric backbones often imparts conformational plasticity and heterogeneous charge distributions.^20,21 Craig et al. demonstrated the binding of tetra cationic metallo-supramolecular pillarplexes to both 4-way (4WJ) and 3-way (3WJ) DNA junctions.²² Similarly, Boer et al. investigated the self-assembly of functionalizable two-component 3D DNA arrays facilitated by Fe-based supramolecular cylinders that induce the formation of 3WJ branch points.²³ Hence, these findings provide a new roadmap for targeting higher-order junction structures using a metallo-supramolecular approach. Macrocyclic drugs have gained clinical approval, with a prominent presence among antiviral agents, most notably those targeting the hepatitis C virus NS3/4A protease, a challenging enzyme featuring a shallow, and solvent-exposed groove that typically resists conventional small-molecule inhibitors.²⁴ Beyond infectious disease, macrocycles have proven invaluable for modulating PPIs across large, flat, and dynamic protein surfaces, where rigidity and shape complementarity confer distinct advantages.²⁵ Very recently, Zhang et al. pioneered a new class of cationic multidecker covalent organic cages, termed pyrgos[n]cages (PC[n]s, n = 1–4).²⁶ These tower-like supramolecules exhibit exceptional rigidity and a uniform distribution of positive charges, unlike traditional flexible scaffolds. The PC[n]s exerted a synergistic antibacterial effect by collapsing the bacterial membrane and interacting with vital cellular components. Remarkably, their bactericidal efficiency increased with the decker number, and they proved highly effective against both extracellular and intracellular drug-resistant pathogens across diverse bacterial infection models. This landmark study demonstrated the potential of rigid, preorganized, and charge-defined supramolecules as a new frontier in biomedical research. Interestingly, the molecular architectures of PC[n]s are similar to those of pillarplexes,²⁷ which are known to modulate DNA recombination and repair upon binding to the 4WJ and 3WJ DNA, and are potential anticancer agents.²²

Here, we computationally investigate a unique class of metal-free supramolecular covalent organic cages, PC[n], to evaluate their interaction behavior and binding stability at the VEGFR-2 interface. Specifically, we assess their association with the extracellular VEGFR-2 D23 domain in both wild-type and mutant forms, as well as their interaction with duplex and forked DNA structures. Three-way junction (3WJ) DNA motifs arise transiently during replication and homologous recombination, and their selective engagement may provide mechanistic insight into cellular regulation and potential therapeutic intervention. We perform high-throughput in silico analysis to determine the relative binding energetics and structural stability of PC[n] macrocycles at protein–nucleic acid interfaces. Molecular dynamics (MD)²⁸ were conducted to examine structural persistence and intermolecular interactions with VEGFR-2 D23 and DNA systems. Complementary density functional theory (DFT)²⁹ calculations were used to probe the electronic structure and intrinsic interaction preferences of PC[n] toward individual nucleobases (adenine, guanine, thymine, and cytosine). Additional analyses, including functional group symmetry-adapted perturbation theory (F-SAPT),^30,31 quantum theory of atoms in molecules analysis (QTAIM),^32,33 noncovalent interaction analysis (NCI),³⁴ and energy decomposition analysis (EDA) and natural bond orbital (NBO)^35,36 calculations, were employed to elucidate the electronic origins, energetic contributions, and charge-transfer characteristics governing these associations. Together, the multiscale computational results provide mechanistic insight into how PC[n] frameworks engage protein–DNA interfaces through cooperative noncovalent interactions. These findings establish a computational foundation for supramolecular interface recognition and support future experimental evaluation of purely organic macrocyclic architectures for angiogenesis-related and nucleic acid-targeted applications.

2. Results and discussion

2.1. Structure and reactivity of PC[n] (n = 1–4)

PC[n] (n = 1–4) represents a π-stacked cationic cage (Fig. 1) featuring imidazolium-based vertical walls capped by phenyl rings. Notably, in PC[4], the terminal wall is distinguished by a benzimidazolium moiety that completes the cage architecture. The axial cage height H1 of PC[2] (Fig. 1a) is 10.093 Å, and that of PC[3] (Fig. 1b) and PC[4] (Fig. 1c) are 15.666 and 20.981 Å, respectively. PC[2] was used as a standard model for computational benchmarking due to its simple stacked structure with an inter-ring spacing (H2) of 5.073 Å and an average lateral window (H3) of 8.197 Å. The top H₂C⋯CH₂ distance across the PC[1–4] structures is ∼5.8 Å. Various methods, B3LYP, PBE0, B3PW91, BLYP, B97, PBE, BP86, BPBE, M05, M06, M052X, M062X, TPSS, CAM-B3LYP, and LC-wPBE, including three-body dispersion D3 corrected functionals and def2TZVPP basis set, are used for the geometry optimization followed by frequency analysis (Section S1, SI). The three diagnostic metrics H1 (opt/cry: 10.096/10.093 Å), H2 (opt/cry: 5.048/5.073 Å), average H3 (opt/cry: 8.183/8.197 Å), and the top H₂C⋯CH₂ distance (opt/cry: 5.822/5.814 Å) of the optimized PC[2] structure (Fig. 1(d)) align well with its experimental structure, claiming the benchmarking of the B3PW91-D3/def2TZVPP level of theory. The comparison of detailed metrics defining the structural features, derived using various levels of theory, suggests that M05-2X, M06-2X, CAM-B3LYP, M06, and PBE0 also lie in the optimal range, as shown in Fig. 1(e). The benchmarking dataset reveals that BLYP (H1: 10.309 Å, H2: 5.155 Å) results in the largest overestimates, while LC-wPBE (H1: 10.031 Å, H2: 5.016 Å) underestimates the most. For H3, all functionals show similar deviations from the experimental structure. The structural RMSD analysis between the experimental and optimized geometry (Fig. 1(e)) shows that all functionals exhibit a similar deviation of ∼0.317 Å. Specifically, M06 yields the lowest RMSD of 0.311 Å, while BLYP yields the highest at 0.326 Å, and B3PW91 remains relatively low at 0.315 Å. In the monomer PC[1] (Fig. 1(f)), the metric H1 = H2 = 4.891 Å is lower than the H2 of PC[2] (opt/cry: 5.048/5.073 Å), PC[3] (5.198 Å), and PC[4] (5.165 Å). This is due to the flattening of each cage through the dispersion-driven π–π stacking, leading to an increase in the interplanar gap with the peripheral region.

Fig. 1 Geometric comparison of PC[n] structures and DFT benchmark against crystal structures, highlighting the three diagnostic metrics H1–H3 (axial cage height H1 (green arrow), inter-ring spacing H2 (blue arrow), and lateral window H3 (black arrow)) and the H₂C⋯CH₂ distance (brown arrow) for (a) PC[2], (b) PC[3], and (c) PC[4]. (d) Gas-phase optimized PC[2] with H1–H3 measured identically for direct comparison with (a). (e) Method benchmarking for PC[2]: H1, H2, and average H3 predicted by a panel of DFT functionals (rows) plotted against the experimental target, along with the corresponding RMSD of the optimized crystal structure. (f) Gas-phase optimized PC[1] (side/front views) displayed with the same geometric references.

Frontier molecular orbital (FMO) analysis reveals that the π-stacked architecture of PC[n] exhibits tunable electronic softness and reactivity governed by π-delocalization and orbital distribution. The energy gap (ΔE_gap) between the highest occupied and the lowest unoccupied molecular orbitals (HOMO–LUMO) systematically decreases from PC[1] (115.836 kcal mol⁻¹, kinetically stable) to PC[4] (8.465 kcal mol⁻¹, highly reactive), reflecting enhanced polarizability and electronic flexibility with increasing decker number (Fig. S1). This reduction originates from progressive π-electron delocalization that narrows orbital separation, rendering PC[4] the most chemically soft and reactive species. The FMO lobes in PC[2–4] are localized at the terminal imidazole and benzimidazole region, indicating site-selective reactivity. Higher orbital transitions (HOMO−5 → LUMO+5 and HOMO−7 → LUMO+7) display larger gaps, corresponding to energetically less relevant excitations; however, in PC[4], the slight difference between HOMO–LUMO (8.465 kcal mol⁻¹) and HOMO−5 → LUMO+5 (42.551 kcal mol⁻¹) denotes dense orbital packing and strong energetic delocalization. Collectively, these results identify PC[4] as the most electronically adaptive and chemically responsive macrocycle, followed by PC[3] and PC[2], ideally suited for dynamic interactions within biomolecular environments.

2.2. Structural analysis of VEGF-A/C/E@VEGFR-2 D23

VEGFR-2 D23 preserves a standard binding blueprint that adapts to each of VEGF-A, VEGF-C, and VEGF-E (Fig. S2), coupling conserved recognition to ligand-tailored activation (see Section S2, computational methodology, SI). In the VEGF-A (PDB ID: 3V2A)³⁷ bound complex (Table S1), a hydrophobic core formed by Val217/219 and Ile256 establishes the primary anchoring interface, while flexible glycine (Gly220/255) and polar residues (Asn274, Ser311) reinforce the interaction through complementary hydrogen bonding and conformational adaptability. The VEGF-C (PDB ID: 2X1X)³⁸ bound complex retains this hydrophobic framework through Ala195, Gly196, and Met197, further stabilized by aromatic stacking from Tyr165/194, and additional polar contributions from Asn253 and Gln132, refining ligand selectivity (Table S2). In the VEGF-E bound (PDB ID: 3V6B)³⁷ complex, similar residues (Ser193, Gly196, Met197, Ala195) dominate the buried interface, supported by Tyr165/194 and Asn253, which mediate π-hydrophobic and hydrogen bonding interactions for structural stabilization (Table S3). As demonstrated by Leppänen et al., the VEGF-D monomeric structure exhibits a high degree of structural conservation with the RMSD of 1.300 Å and 1.100 Å across 96 C_α atoms of VEGF-C and 91 C_α atoms of VEGF-A, respectively.³⁹ Interestingly, VEGF-C and VEGF-D homology domains share approximately 60% sequence identity, and VEGF-D conserves nearly all residues implicated in VEGFR-2 D23 binding activation. Therefore, any ligand capable of inhibiting VEGF-A/C association with VEGFR-2 D23 is also expected to impede VEGF-D binding effectively. Collectively, these interaction patterns reveal a modular binding mechanism in which the D2 domain provides hydrophobic anchoring. In contrast, the D3 domain contributes electrostatic complementarity and aromatic stabilization, together enabling high-affinity ligand recognition, receptor dimerization, and the initiation of angiogenic signalling.

2.3. MD simulation of PC[n] with the VEGFR-2 D23

Microsecond-scale all-atom MD simulations (see Section S3, computational methodology, SI) suggest that PC[n] macrocycles (n = 1–4) may associate with the canonical VEGF docking cleft located between VEGFR-2 domains D2 and D3. The ligands appear to localize near an electrostatic interaction region involving Glu284 on D3, while forming loop-resolved contacts across loops L1–L5 (Fig. 2a). These interactions may involve a combination of hydrogen-bonding, polar contacts, and hydrophobic interactions that contribute to interfacial stabilization. Across all complexes, the D2–D3 bending angle remains within the experimentally observed range (≈122–149°),³⁸ suggesting that ligand engagement does not significantly perturb the native interdomain mechanics, while potentially contributing to steric and electrostatic changes at the interface that could influence VEGF accessibility. PC[1] does not exhibit a persistent association with VEGFR-2 D23, which may be attributed to its smaller macrocyclic scaffold and limited geometric adaptability, as suggested by FMO analysis, reducing its ability to span the shallow D2–D3 interface. In contrast, PC[2] remains localized at the D2–D3 junction for approximately 100–800 ns (Fig. 2b) and, following a transient increase in bending angle (136.14° at 900 ns), adopts a more stabilized configuration near loop L5 by 1000 ns, with a mean angle of 131.85 ± 6.47° (defined using the Val219–Gly220–Tyr221 linker). A cooperative interaction network may form involving residues from L2 (Tyr190, Ser193, Tyr194), L3 (Val218), and L5 (Leu252, Gly282, Ser283, Glu284), including short polar contacts with Glu140 (OE1, 2.74 Å), Ser193 (OG, 2.87 Å), and Glu284 (OE2, 2.78 Å), as well as carbonyl-proximal interactions with Gly282 (O, 3.71 Å) and Ser283 (O, 3.91 Å) (Fig. 2c; Table S4). These interactions fall within typical hydrogen-bonding and van der Waals interaction distances (≤4 Å) and may contribute to favorable interaction energetics through shape complementarity and reduced solvent exposure. PC[3] adopts a vertically tilted orientation (Fig. S3(a)) primarily associated with Glu284 (OE2, 3.01 Å), with intermittent contacts involving Tyr194, Ser193, and Val218 (Fig. S3(b), and Table S5). This configuration is accompanied by greater bending angle variability (128.40 ± 7.80°), particularly during the 100–200 ns and 700–800 ns intervals, suggesting comparatively weaker interfacial retention relative to PC[2]. In contrast, PC[4] adopts a more extended orientation across the cleft (Fig. S4a), transiently occupying the D2–D3 seam during the early simulation period (≤600 ns), while maintaining an average bending angle of 131.66 ± 5.48°. This configuration is associated with hydrophobic contacts involving Tyr165 (CB, 3.49 Å), Tyr194 (CB, 3.65 Å), and Leu252 (CD2, 3.15 Å), as well as polar interactions with Gly282 (O, 3.17 Å) and Glu284 (OE2, 2.97 Å) (Fig. S4b; Table S6), suggesting favorable structural complementarity at the interdomain interface.

Fig. 2 MD simulation analysis of the VEGFR-2 D2–D3 domain complexed with the PC[2]. (a) Structural architecture of VEGFR-2 D2–D3 showing β-strands and loop regions (L1–L5). (b) Time-resolved conformational transitions of VEGFR-2@PC[2] over a 1000 ns MD trajectory, with the inner plot depicting the bending angle fluctuation. (c) Binding interface of the PC[2] with VEGFR-2, highlighting key interacting residues and hydrogen bond distances (Å). The top panel shows views of the complex with 180° rotation for spatial clarity. (d) RMSD plots of the VEGFR-2 D23@PC[2] complex, VEGFR-2 D23, and PC[2] illustrating structural stability and equilibration throughout the simulation.

The VEGF-A–VEGFR-2 D23 interaction has been previously experimentally characterized using alanine-scanning mutagenesis, which identified key residues involved in ligand recognition. Based on these findings, the VEGFR-2 Leu252Ala/Asn253Ala double mutant, located within the D2–D3 interdomain binding cleft, was selected to computationally examine the structural and electrostatic contributions of these residues to ligand recognition. The highly reactive PC[4] macrocycle was used as a representative probe to evaluate how these mutations may influence supramolecular binding behavior and associated interaction patterns at the interface. This analysis suggests a structure–function relationship between receptor site perturbation and macrocyclic adaptability. PC[4] may re-anchor deeper within the D3 domain while maintaining the native interdomain geometry (128.96 ± 6.78°) (Fig. S5a). It potentially forms a short hydrogen bond with Asp257 (OD1, 2.75 Å) and a secondary interaction with Ser310 (O, 3.12 Å), along with additional longer-range contacts (Fig. S5b and Table S7) involving Gly255 (O, 3.61 Å), Phe258 (O, 3.66 Å), and hydrophobic interactions with Leu272 (CD1, 3.35 Å), with proximity to Ile256 (C, 3.70 Å). The residue Glu284 appears to serve as a key electrostatic interaction site contributing to a stable polar interaction network; under mutation, this interaction may shift toward adjacent D3 residues such as Asp257 and Ser310. Surrounding this region, the ligands may engage the amphiphilic interface through polar interactions involving Glu140, Ser193, Gly282, Ser283, Glu284, and Asp257, as well as hydrophobic contacts with Tyr165, Tyr194, Leu252, Ile256, and Leu272. These interactions may contribute to favorable desolvation and improved shape complementarity. Throughout the simulations, the native D2–D3 bending angle remains preserved, while ligand positioning at the interface may contribute to structural and electrostatic changes that could influence VEGF docking. These observations provide a computational framework for understanding macrocycle-mediated modulation of shallow PPI interfaces.

RMSD and MM/GBSA (see Section S4, computational methodology, SI) analyses collectively suggest that an interplay between receptor flexibility and ligand electronic adaptability may govern VEGFR-2 D23@PC[n] recognition. Across all trajectories, conformational dynamics appear to be dominated by intrinsic D2–D3 hinge bending. At the same time, PC[n] ligands remain consistently localized near an electrostatic interaction region centered on Glu284, with ligand RMSD values ranging from 0.51 to 1.66 Å. This localization is supported by amphiphilic interface engagement and hydrophobic contacts involving residues such as Tyr165, Tyr194, and Leu252. PC[2] maintains a relatively stable association (RMSD: 4.97 Å, Fig. 2(d)), potentially facilitated by hydrogen bonding with Glu284, Glu140, and Ser193, along with interactions involving the L2–L5 loop regions, and is associated with a favorable binding enthalpy (ΔH_total = −47.7 ± 6.2 kcal mol⁻¹). In contrast, PC[3] appears to be stabilized primarily through interactions with Glu284, permitting greater interdomain flexibility (RMSD: 7.08 Å, Fig. S3(c)) while maintaining interfacial association, consistent with a slightly less favorable enthalpy (ΔH_total = −43.9 ± 6.5 kcal mol⁻¹). The π-extended PC[4] exhibits the most extensive interface engagement, combining hydrophobic packing with directional hydrogen bonding involving residues such as Glu284 and Gly282, which may contribute to localized structural stabilization despite overall receptor flexibility (RMSD: 9.18 Å, Fig. S4(c)). Energetically, MM/GBSA decomposition analysis supports these interaction trends. PC[4] shows the most favorable binding enthalpy (ΔH_total = −48.3 ± 3.5 kcal mol⁻¹, Table S8) arising from a balance between van der Waals interactions (ΔE_vdW = −32.8 ± 3.5 kcal mol⁻¹) and electrostatic contributions (ΔE_elec = 654.5 ± 43.4 kcal mol⁻¹), partially offset by polar solvation effects (ΔE_polar = −616.7 ± 42.7 kcal mol⁻¹). Mutation of Leu252/Asn253 to Ala alters local packing and hydrogen-bonding interactions (RMSD: 8.10 Å, Fig. S5(c)). It is associated with a reduced calculated binding enthalpy to −24.3 ± 5.8 kcal mol⁻¹, suggesting the importance of these residues in stabilizing ligand association. Collectively, these computational observations suggest that PC[4], due to its extended π-system and conformational adaptability, may achieve enhanced interaction complementarity with VEGFR-2 D23 through favorable dispersion and electrostatic interactions. This mode of interface engagement could influence the structural and electrostatic environments of the VEGF-binding region. It may provide valuable insights for the computational design of supramolecular modulators targeting shallow PPI interfaces.

2.4. MD simulations of PC[n] with the DNA

DNA represents a critical therapeutic target because it encodes, organizes, and regulates the cell's genetic information. Small molecules that interact with DNA may influence genome processing, as demonstrated by clinically used duplex-binding agents such as cisplatin and doxorubicin. Importantly, the genome does not exist solely in the canonical B-form duplex. During essential cellular processes, including replication, repair, and transcription, DNA may transiently adopt alternative conformations, such as replication forks and related intermediates. These non-duplex structures are increasingly recognized as important regulatory and structurally dynamic regions involved in genome maintenance. Accordingly, the rational design of molecules capable of selectively recognizing duplex and fork-like DNA structures may provide opportunities to modulate DNA-associated processes. It could contribute to the development of targeted anticancer strategies.

Across the MD trajectories, the PC[n] (n = 2–4) display a unified amphiphilic strategy for DNA (nucleobases are numbered as per Scheme S1) recognition, characterized by persistent groove association with minimal structural distortion. This interaction mode appears to be supported by a combination of directional hydrogen bonding, electrostatic interactions with phosphate groups, and dispersion-driven π-surface complementarity. PC[2], the smallest macrocycle, preferentially localizes within the minor groove, adopting a compact orientation centered on DC4–DG5–DA6 and cross-strand DG25 (Fig. 3(a)). Multiple donor–acceptor interactions are observed (Table S9), including contacts with DC4 C2(=O)⋯H–C (2.95 Å), DG5 sugar(O)⋯H–C (2.76 Å), DG5 sugar(H)⋯N–C (2.42 Å), and DA6 phosphate(=O)⋯H–C (2.45–2.50 Å) and sugar (H)⋯N–C (2.68 Å). Additional electrostatic interactions involving the DG25 phosphate group (2.51/2.67 Å) suggest a multivalent interaction network that may contribute to favorable interaction energetics and reduced solvent exposure. PC[3] (Fig. S6(a)) preferentially associates with the major groove, forming an extended hydrogen bonding network with phosphate oxygens (Table S10) of DC4, DG5, DG17, and DA18 (2.36–3.04 Å). These interactions appear to provide electrostatic stabilization while preserving overall helical geometry, with no substantial perturbation to twist or rise parameters. The complex exhibits relatively limited translational displacement and localized conformational fluctuations near the binding region, consistent with multivalent interaction stabilization. The π-expanded PC[4] (Fig. S7(a)) demonstrates more extensive major groove association (Table S11), forming a dense network of hydrogen bonds involving DG5 phosphate oxygen (DG5 phosphate (=O)⋯H–C 2.47 Å; DT20 C4(=O)⋯H–C 2.18 Å; DA6 N7⋯H–C 2.93 Å). These interactions are complemented by π–π stacking (5.2 Å) and C–H⋯π contacts (2.18–3.04 Å) involving DC4–DG5–DA6–DA7 and cross-strand DG17–DT20. The macrocyclic framework appears to undergo conformational adjustment to improve interaction complementarity and surface burial, which may contribute to localized structural stabilization and reduced groove flexibility. In the case of forked 3WJ DNA, PC[2] (Fig. S8(a)) associates near the junction region, bridging two arms while interacting with the third through π⋯π and CH⋯π contacts involving DA5, DT10, and DT16 (3.7–5.8 Å). Additional hydrogen-bonding interactions are observed with the sugar and phosphate groups of DA11 and the C2(=O) group of DT10, with interaction distances of 2.47, 2.99, and 3.06 Å, respectively (Table S12). This interaction pattern may contribute to localized stabilization of the junction architecture while preserving native base-pairing interactions, suggesting favorable structural and electrostatic complementarity between the macrocycle and the DNA junction environment.

Fig. 3 MD simulation analysis of the DNA@PC[2] complex. (a) Structural representation of the ligand bound within the DNA double helix, highlighting key hydrogen-bonding interactions with nucleobases. The distances are given in Å. The inset shows the overall orientation of DNA@PC[2] binding, highlighting the groove-specific interaction mode. (b) RMSD profiles of the complex, DNA, and ligand over a 500 ns MD trajectory.

RMSD and MM/GBSA analyses collectively suggest conformational stability and enthalpically favorable interaction patterns in the PC[n]–DNA assemblies, highlighting their structural compatibility and electrostatic complementarity with the DNA framework. Across all systems, the PC[n] macrocycles exhibit low ligand RMSDs (≤1 Å), indicating persistent groove engagement and favorable geometric alignment with the DNA structure. PC[2] localizes within the minor groove (complex: 5.48 Å, DNA: 2.18 Å, and ligand: 0.65 Å, Fig. 3(b)), where its association appears to be supported by hydrogen bonding and phosphate-mediated electrostatic interactions. Similarly, PC[3] (complex: 2.73 Å, DNA: 2.03 Å, and ligand 0.58 Å, Fig. S6(b)) adopts a strand-bridging configuration that may contribute to duplex stabilization through multivalent contacts and localized conformational fluctuations. The π-extended PC[4] exhibits somewhat greater ligand flexibility (complex: 3.13 Å, DNA: 1.01 Å, and ligand: 2.53 Å, Fig. S7(b)), which may reflect adaptive major-groove association facilitated by cooperative hydrogen bonding and dispersion-mediated π interactions. In contrast, PC[2] associated with forked 3WJ DNA (Fig. S8(b)) shows elevated RMSD values for the complex (5.08 Å) and DNA (4.86 Å), likely reflecting intrinsic flexibility of the junction arms. Despite this, the ligand itself remains conformationally stable (0.52 Å), suggesting persistent localization near the branchpoint region. Energetically, MM/GBSA analysis indicates that PC[2–4] binding is predominantly enthalpy-favored, with calculated total enthalpies (ΔH_total) of −26.0 ± 5.3, −17.0 ± 4.6, and −44.8 ± 3.8 kcal mol⁻¹ for PC[2], PC[3], and PC[4], respectively (Table S13). These interaction energies arise primarily from favorable van der Waals contributions (−19.9 to −62.1 kcal mol⁻¹) counterbalanced by electrostatic-solvation compensation. Notably, PC[2] binding to forked 3WJ DNA (ΔH_total = −46.3 ± 4.5 kcal mol⁻¹) exhibits comparable calculated interaction energetics to PC[4], which may be associated with cooperative hydrogen bonding, π-mediated stacking interactions, and electrostatic contacts with phosphate groups at the junction. Together, these computational results suggest a trend in interaction energetics, with PC[4] showing the most favorable calculated binding enthalpy, followed by PC[2] and PC[3]. This trend may reflect differences in macrocycle size, electronic structure, and conformational adaptability, which influence interaction complementarity with the DNA surface. Collectively, these observations indicate that DNA recognition by PC[n] macrocycles may be governed by electrostatic complementarity, dispersion interactions, and geometric compatibility with the DNA grooves and junction regions. These findings provide a computational basis for understanding macrocycle–DNA interaction mechanisms and may inform the rational design of supramolecular systems targeting diverse DNA structural motifs.

2.5. Odd–even effect in the energy analysis

From the simulation, even-membered PC[n] frameworks exhibit more favorable binding enthalpies and reduced ligand RMSD compared to odd-membered analogues, indicating enhanced structural stability and interaction persistence at biomolecular interfaces. MD trajectory analysis suggests that even-membered macrocycles adopt more symmetric, conformationally stable geometries, as evidenced by ligand RMSD, enabling improved spatial alignment of functional interaction sites and greater buried surface area upon binding. This geometric preorganization appears to facilitate more persistent intermolecular contacts, particularly in DNA–PC[n] systems, where the odd–even effect is most pronounced. At the same time, a similar, though less distinct, trend was observed in VEGFR-2 D23 interactions. FMO analysis further reveals that even-membered systems exhibit more delocalized and symmetrically distributed HOMO–LUMO orbitals, consistent with enhanced electronic polarizability and interaction adaptability. Collectively, these structural and electronic features suggest that geometric symmetry and electronic delocalization contribute to the improved interaction stability of even-membered PC[n] macrocycles.

2.6. Binding energy analysis

The isolated PC[1]–nucleobase complexes were examined to establish the intrinsic electronic interaction preferences between the PC[n] cage scaffold and nucleobase functional groups in the absence of structural and environmental complexity. These quantum-chemical models enable direct quantification of directional hydrogen bonding, lone-pair donation, polarization, and dispersion interactions, which are fundamental contributors to supramolecular recognition. Although duplex DNA introduces additional stabilizing factors such as base stacking, phosphate backbone electrostatics, and conformational constraints, the isolated nucleobase calculations provide a well-defined electronic reference framework for interpreting macrocycle–DNA interactions. Notably, the base-dependent interaction hierarchy and preferred hydrogen-bonding geometries identified in the DFT calculations are consistent with the persistent nucleobase-groove contacts observed in the explicit-solvent molecular dynamics simulations of full DNA duplex and junction structures. This agreement supports the role of intrinsic electronic complementarity as a contributing factor underlying the structurally resolved binding modes. Accordingly, the isolated-base DFT analysis provides mechanistic insight into the electronic origin of recognition. At the same time, the atomistic MD simulations capture the structural, dynamical, and environmental determinants of binding in the complete DNA system.

Each nucleobase was positioned along the horizontal and vertical axes of the macrocycle, revealing a consistent preference for corner-stacking configurations rather than direct placement over the phenyl caps or within the imidazolium cavity (Fig. S9). The resulting local minima (Fig. S10) show that nucleobases anchor to the terminal periphery through directional C–H⋯N/O hydrogen bonds, forming discrete, site-specific complexes. Adenine interacts through N3⋯H60 (2.53 Å), N3⋯H65 (2.16 Å), and N7⋯H70 (2.08 Å) (Fig. S9), yielding a binding energy (ΔE_b) of −38.796 kcal mol⁻¹ (Table 1), while guanine engages via N3⋯H28 (2.07 Å), N7⋯H43 (2.18 Å), and N7⋯H48 (2.52 Å), producing ΔE_b of −30.645 kcal mol⁻¹. Thymine exhibits the weakest ΔE_b (−23.632 kcal mol⁻¹), stabilized primarily by 4C=O⋯H hydrogen bonds (2.07–2.50 Å). In contrast, cytosine exhibits the strongest ΔE_b (−43.326 kcal mol⁻¹) through a combination of 2C=O⋯H47 (2.37 Å), 2C=O⋯H74 (2.04 Å), and N1⋯H72 (2.01 Å) contacts, highlighting its dual donor–acceptor capacity and superior electrostatic complementarity. The thermodynamic parameters (ΔH and ΔG of Table 1) corroborate that binding is strongly enthalpy-driven, with ΔG values rendered less negative due to expected entropic penalties from loss of translational and rotational freedom (ΔS ≈ −0.030 to −0.040 kcal mol⁻¹ K⁻¹). At the B3PW91-D3/def2-TZVPP level, the enthalpies of association (ΔH of −36.488, −22.161, −41.017, and −29.272 kcal mol⁻¹ for A, T, C, and G, respectively) and corresponding ΔG (−26.100, −11.530, −30.247, and −18.742 kcal mol⁻¹) reveal cytosine's superior binding stability, consistent with its bifunctional hydrogen bonding pattern. To evaluate methodological consistency, the interaction energetics were further examined using multiple dispersion-corrected density functionals, including B3PW91-D3, BLYP-D3, PBE0-D3, M06-2X, TPSS-D3, and CAM-B3LYP-D3, along with anharmonic corrections. All computational approaches reproduced a consistent interaction energy trend (C > A > G ≫ T), suggesting that the observed recognition pattern is robust across different exchange–correlation functionals (Table S14). Collectively, these results provide computational evidence that PC[1] exhibits base-dependent interaction preferences, primarily governed by hydrogen bonding and electrostatic complementarity, supporting its role as a modular recognition unit within higher-order PC[n] supramolecular assemblies.

Table 1 Binding energetics of PC[1] macrocycle with individual DNA nucleobases computed using the B3PW91-D3/def2-TZVPP level of theory. Reported values include binding energy (ΔE_b), zero-point energy (ZPE)-corrected binding energy (E_b-ZPE), basis set superposition error (BSSE)-corrected energy (E_b-BSSE), Gibbs free energy change (ΔG), and enthalpy change (ΔH). All energies are in kcal mol⁻¹

Systems	ΔEb	ΔEb-ZPE	ΔEb-BSSE	ΔG	ΔH
PC[1]@adenine	−38.796	−37.301	−39.510	−26.100	−36.488
PC[1]@thymine	−23.632	−22.623	−25.650	−11.530	−22.161
PC[1]@cytosine	−43.326	−41.695	−44.930	−30.247	−41.017
PC[1]@guanine	−30.645	−29.976	−31.870	−18.742	−29.272

---

2.7. Qualitative and quantitative interaction analysis

To elucidate the nature of the interaction, NCI analysis was performed based on the wavefunction obtained at the B3PW91-D3/def2TZVPP level of theory (see Section S5, computational methodology, SI). The NCI plots reveal that the nucleobases are strongly influenced by the hydrogen bonding, as visualized through the dark blue spattered region (Fig. 4). The green surfaces claim the strong, attractive interactions arising from electrostatic contributions. The ring strain of the aromatic core regions localizes the red-spattered regions. To reinforce the nature of these interactions, the interaction strengths were estimated from the electron density, ρ(r), at the BCP and its corresponding Laplacian, ∇²ρ(r), via QTAIM analysis. The quantitative and qualitative analyses reveal that the C–H⋯N/C–H⋯O and H⋯H hydrogen bonds between the nucleobase and the PC[1] form a triangular bond path (Fig. S11). Furthermore, Espinosa–Molins–Lecomte⁴⁰ (EML) a semi-empirical equation was used to determine the hydrogen bond strength. Adenine binds with PC[1] with a total hydrogen bond strength of −14.375 kcal mol⁻¹, attributed to N⋯H hydrogen bonds like those of guanine, −14.346 kcal mol⁻¹. Thymine has a weaker hydrogen bond strength of −13.382 kcal mol⁻¹ predominantly from O⋯H hydrogen bonds, while cytosine binds with a maximum of −17.166 kcal mol⁻¹ via O⋯H and N⋯H hydrogen bonds. The strength of hydrogen bonding is also linearly correlated to the total interaction energy (Fig. S12). Across all systems, ρ(r) falls within 0.009–0.028 a.u., a regime characteristic of weak but significant closed-shell interactions, including hydrogen bonding and H⋯H contacts (Table S15). The uniformly positive Laplacian, ∇²ρ(r) > 0, further corroborates the electrostatic and noncovalent nature of these interactions, distinguishing them from shared-shell (covalent) interactions. The negative potential energy density, V(r), observed for all contacts indicates net stabilizing interactions, consistent with hydrogen bonding or van der Waals attraction. In contrast, the corresponding positive kinetic energy densities, G(r), reflect the dominance of electron density depletion at these critical points, a feature typical of closed-shell interactions. The total energy density K(r) remains slightly negative or near zero, suggesting interactions that lie on the boundary between pure hydrogen bonds and dispersion-dominated contacts, commonly observed in supramolecular host-nucleobase assemblies. Bond ellipticity (ε) values provide further insight into the anisotropy of the electron density distribution. Higher ε values (>1), particularly in H⋯H contacts, indicate significant anisotropy and reduced directional stability, reflecting the flexible nature of these interactions. Conversely, low ε values (<0.1) for H⋯N and H⋯O contacts indicate more directionally constrained and structurally stable interactions, consistent with classical hydrogen-bond behavior. Collectively, these topological descriptors clearly establish that PC[1] engages nucleobases through a cooperative network of weak hydrogen bonds and H⋯H contacts. Although individually weak, these interactions cumulatively contribute to the overall stability of the PC[1]–DNA binding interface, aligning with the supramolecular architecture's ability to recognize and anchor nucleobases.

Fig. 4 NCI analysis of the complex measured at an isosurface value of 0.5 a.u., illustrating the spatial distribution and strength of weak intermolecular interactions. The green regions represent van der Waals interactions, blue regions indicate strong attractive interactions such as hydrogen bonding, and red regions denote steric repulsion.

2.8. Energy decomposition and natural bond orbital analyses

EDA and NBO analyses together elucidate the electronic origins of PC[1]–nucleobase recognition, revealing that host–guest stabilization arises from a synergy of electrostatics, polarization, and lone-pair-driven charge delocalization (see Section S6 Computational methodology in the SI). FSAPT0-based EDA shows that electrostatic energy (E_elec) dominates the attractive potential, supplemented by moderate, induction energy (E_ind) and dispersion energy (E_disp) components (Fig. 5a), with total interaction energies following the hierarchy PC[1]–C (−46.340 kcal mol⁻¹) > A (−36.840 kcal mol⁻¹) > G (−28.750 kcal mol⁻¹) > T (−23.660 kcal mol⁻¹), consistent with hydrogen bond propensities and QTAIM evidence of charge accumulation at N⋯H/O⋯H contacts. Cytosine exhibits the strongest electrostatic stabilization due to its localized lone pairs, while guanine shows enhanced dispersion (18.31%, Table S16) from π–π complementarity. Exchange repulsion (E_exc) scales with electron density localization, peaking in cytosine (22.509 kcal mol⁻¹). Across all complexes, the nucleobases exhibit pronounced polarization within the electrostatic field of PC[1]. The density difference maps (Fig. 5b) visualize this effect, showing electron depletion at the nucleobase O/N lone pairs (red) and electron accumulation around the donor H atoms (blue), evidencing charge displacement along the H⋯O/H⋯N hydrogen bonds. Charge-transfer analysis confirms predominant electron flow from PC[1] to the nucleobases (−11.279 kcal mol⁻¹, except thymine −8.787 kcal mol⁻¹), with minor back-donation, maximized for cytosine (−3.731 kcal mol⁻¹). Thymine accepts the least electron density because its strongly resonance-stabilized C=O groups (at C2 and C4) and C5-methyl substitution limit the availability of its lone pairs for charge accommodation. In contrast, adenine, guanine, and cytosine exhibit more accessible and less stabilized heteroatom lone pairs, allowing them to participate in more feasible donor–acceptor recognition. This structural flexibility makes the other nucleobases inherently better electron acceptors. NBO analysis confirms that lone-pair (LP → σ*) delocalization governs donor–acceptor strength, with cytosine (E(2) = 12.310 kcal mol⁻¹, Table S17) and adenine (8.130 kcal mol⁻¹) exhibiting the strongest interactions via pyrimidine and imidazole nitrogen donors, respectively. These results demonstrate that, rather than purely Coulombic attraction, lone-pair-driven polarization and charge transfer dictate the selective, robust binding of PC[n] macrocycles, underpinning their dual recognition of nucleobases and VEGFR-2 D23 through electronically adaptive, π-polarized interfaces.

Fig. 5 EDA of PC[1]–nucleobase interactions. (a) Decomposition of the total interaction energy (E_total) into electrostatic (E_elec), exchange (E_exc), induction (E_ind), and dispersion (E_disp) components for adenine, thymine, cytosine, and guanine complexes, illustrating the relative contributions of each term to the overall stabilization. (b) Electronic polarization density plots of PC[1] (fragment A) interacting with each nucleobase (fragment B), visualized at a contour level of 0.001 a.u. The red and blue lobes represent the redistribution of polarization density between the fragments. Red arrows indicate electronic polarization from A to B, while green arrows denote back-polarization from B to A, with interaction energies (in kcal mol⁻¹) annotated for each direction.

3. Comparison with experimental structures

The experimental crystal structures of the VEGFR-2 D2–D3 domain in complex with VEGF ligands, including VEGF-A (PDB ID: 3V2A), VEGF-C (PDB ID: 2X1X), and VEGF-E (PDB ID: 3V6B), provide structural insight into ligand recognition based on buried accessible surface area (ASA) analysis. For the VEGF-A@VEGFR-2 D23 complex, residues such as Val217, Gly220, Gly255, Ile256, Phe288, and Ser311 exhibit high buried ASA values (>90%), indicating their importance in ligand engagement (Table S1). Similarly, in the VEGF-C@VEGFR-2 D23 complex, residues Ala195, Gly196, Met197, and Tyr165 show significant buried ASA contributions (Table S2), while Tyr194, Ala195, Gly196, Leu252, and Asn253 are prominently involved in VEGF-E binding (Table S3). In comparison, the computationally simulated VEGFR-2 D23@PC[2–4] systems reveal that several of these residues, including Tyr165, Tyr194, Leu252, Gly182, Ser283, and Glu284, overlap with experimentally identified hotspot regions involved in VEGF-C and VEGF-E recognition (Tables S4–S7), suggesting that PC[2–4] macrocycles may associate with structurally relevant regions of the VEGFR-2 binding interface. Although the PC[2–4] binding sites do not directly coincide with the primary VEGF-A hotspot residues, their persistent localization near the D2–D3 interdomain cleft and interactions with adjacent structural elements may influence the local steric and electrostatic environment of the receptor interface.

Given the cationic nature and membrane-associating properties of PC[n], their potential to interact with negatively charged intracellular targets, such as DNA, has been experimentally demonstrated by Zhang et al.²⁶ Consistent with these observations, our computational analyses indicate that PC[n] macrocycles exhibit favorable electrostatic complementarity with DNA, driven primarily by interactions with the negatively charged phosphate backbone and nucleobase heteroatoms. Experimental gel retardation and fluorophotometric competition assays have reported enhanced DNA-binding propensity with increasing cage size, showing comparatively stronger binding. In agreement with this experimental evidence, our simulations and energy decomposition analyses suggest that higher-order PC[n] architectures exhibit increased interaction surface area, enhanced π-electronic delocalization, and greater multivalent interaction capability, which may contribute to stronger and more persistent association with DNA. These computational findings are consistent with experimentally observed size-dependent binding trends and support the role of multivalent electrostatic and dispersion interactions in governing PC[n]–DNA recognition.

4. Conclusions

This study computationally investigates covalent organic pyrgos[n]cages, PC[n] (n = 1–4), as a class of electronically adaptive supramolecular macrocycles capable of interacting with both protein and DNA interfaces. Using an integrated multiscale simulation framework, the results provide mechanistic insights into their supramolecular recognition behavior and potential relevance for biomolecular modulation. The progressively π-stacked PC[n] architectures exhibit increased electronic delocalization and molecular softness as the decker number increases from 1 to 4. A distinct odd–even trend is observed, with even-numbered deckers showing more favorable calculated interaction energetics compared to odd-numbered analogues. The systematic narrowing of the HOMO–LUMO gap from 115.836 kcal mol⁻¹ in PC[1] to 8.465 kcal mol⁻¹ in PC[4] suggests an enhanced capacity for charge redistribution and electronic adaptability in response to local electrostatic environments. At the VEGFR-2 D2–D3 interface, PC[2–4] exhibits persistent interfacial association mediated by hydrogen bonding, electrostatic interactions, and π-mediated contacts involving residues such as Tyr165, Tyr194, Ser193, Gly255, Ile256, and Glu284. These interactions appear to localize the macrocycles near the canonical ligand-binding region while preserving the overall structural stability of both wild-type and mutant receptor conformations. In DNA systems, PC[n] macrocycles preferentially associate within groove regions, forming interaction networks supported by hydrogen bonding, electrostatic contacts with phosphate groups, and dispersion-driven surface complementarity, while maintaining overall helical geometry. F-SAPT and NBO analyses suggest that nucleobase lone-pair donation, particularly from nitrogen atoms, contributes significantly to interaction strength and directionality, consistent with polarization-mediated stabilization. QTAIM analysis further supports the presence of well-defined hydrogen bonding interactions through characteristic electron density distributions at bond critical points. DFT and energy decomposition analyses indicate an enthalpy-favored interaction hierarchy (C > A ≈ G ≫ T), dominated primarily by electrostatic and induction contributions, with additional stabilization arising from dispersion interactions. Collectively, these computational findings suggest that PC[n] macrocycles exhibit electronically cooperative and structurally adaptable interaction behavior, enabling multivalent recognition at both protein and DNA interfaces. Their intrinsic electronic flexibility and tunable supramolecular architecture provide a computational design framework for developing macrocyclic systems with enhanced interface complementarity and interaction specificity, with potential relevance for targeting protein–protein and protein–DNA interaction regions.

Author contributions

Conceptualization, B. M.; methodology, B. M.; investigation, B. M., and P. M.; writing – original draft, B. M.; writing – review & editing, B. M. and P. M.; funding acquisition, P. M.; resources, P. M.; supervision, P. M.

Conflicts of interest

There are no conflicts to declare.

Data availability

All data generated or analyzed during this study are included within the published article.

Supplementary information (SI): computational methodologies, Fig. S1–S12 and Tables S1–S17, and the coordinates of the optimized geometries of the PC[1], PC[2], and PC[1]–nucleobases are included. See DOI: https://doi.org/10.1039/d5cp04983k.

Acknowledgements

Shiv Nadar Institution of Eminence Deemed to be University is thanked for providing fellowship to Biswajit Mohanty, computational resources through MAGUS-02, and overall research support. Dr Kshatresh Dutta Dubey is thanked for their inspiring discussions.

Notes and references

1. Z.-L. Liu, H.-H. Chen, L.-L. Zheng, L.-P. Sun and L. Shi, Signal Transduction Targeted Ther., 2023, 8, 198.
2. L. A. DiPietro, J. Leucocyte Biol., 2016, 100, 979–984.
3. S. B. Fox, G. Gasparini and A. L. Harris, Lancet Oncol., 2001, 2, 278–289.
4. R. Kerbel and J. Folkman, Nat. Rev. Cancer, 2002, 2, 727–739.
5. P. Carmeliet and R. K. Jain, Nature, 2000, 407, 249–257.
6. M. Potente, H. Gerhardt and P. Carmeliet, Cell, 2011, 146, 873–887.
7. C. Lee, M.-J. Kim, A. Kumar, H.-W. Lee, Y. Yang and Y. Kim, Signal Transduction Targeted Ther., 2025, 10, 170.
8. M. Shibuya, Genes Cancer, 2011, 2, 1097–1105.
9. M. L. George, M. G. Tutton, F. Janssen, A. Arnaout, A. M. Abulafi, S. A. Eccles and R. I. Swift, Neoplasia, 2001, 3, 420–427.
10. M. Shibuya, Cancer Sci., 2003, 94, 751–756.
11. J. Zhou and J. Rossi, Nat. Rev. Drug Discovery, 2017, 16, 181–202.
12. X. Wang, D. Ni, Y. Liu and S. Lu, Front. Chem., 2021, 9, 682675.
13. L. Wang, W.-Q. Liu, S. Broussy, B. Han and H. Fang, Front. Pharmacol., 2024, 14, 1307860.
14. L. Mabonga and A. P. Kappo, Biophys. Rev., 2019, 11, 559–581.
15. P. B. Crowley, Acc. Chem. Res., 2022, 55, 2019–2032.
16. L. Armstrong, S. L. Chang, N. Clements, Z. Hirani, L. B. Kimberly, K. Odoi-Adams, P. Suating, H. F. Taylor, S. A. Trauth and A. R. Urbach, Chem. Soc. Rev., 2024, 53, 11519–11556.
17. C. Muankaew and T. Loftsson, Basic Clin. Pharmacol. Toxicol., 2018, 122, 46–55.
18. S. Li, H. Chen, X. Yang, D. Bardelang, I. W. Wyman, J. Wan, S. M. Y. Lee and R. Wang, ACS Med. Chem. Lett., 2015, 6, 1174–1178.
19. M. Dockerill, D. J. Ford, S. Angerani, I. Alwis, L. J. Dowman, J. Ripoll-Rozada, R. E. Smythe, J. S. Liu, P. J. B. Pereira and S. P. Jackson, Nat. Biotechnol., 2025, 43, 186–193.
20. T. Haino, Polym. J., 2013, 45, 363–383.
21. K. Velmurugan, M. Mohan, B. Li, K. Wang, M. Zuo and X.-Y. Hu, Mater. Adv., 2020, 1, 2646–2662.
22. J. S. Craig, L. Melidis, H. D. Williams, S. J. Dettmer, A. A. Heidecker, P. J. Altmann, S. Guan, C. Campbell, D. F. Browning, R. K. O. Sigel, S. Johannsen, R. T. Egan, B. Aikman, A. Casini, A. Pöthig and M. J. Hannon, J. Am. Chem. Soc., 2023, 145, 13570–13580.
23. D. R. Boer, J. M. C. A. Kerckhoffs, Y. Parajo, M. Pascu, I. Usón, P. Lincoln, M. J. Hannon and M. Coll, Angew. Chem., Int. Ed., 2010, 49, 2336–2339.
24. X. Yang, Z. Jiao, C. Fan, W. Ye, S. Lv, X. Wang, C. Wang, K. Zhang, X. Ke and W. Zhou, J. Med. Chem., 2025, 68, 21035–21071.
25. P. G. Dougherty, Z. Qian and D. Pei, Biochem. J., 2017, 474, 1109–1125.
26. Y. Zhang, M. Luo, X. Shi, A. Li, W. Zhou, Y. Yin, H. Wang, W.-L. Wong, X. Feng and Q. He, Sci. Adv., 2024, 10, eadp4872.
27. P. J. Altmann and A. Pöthig, J. Am. Chem. Soc., 2016, 138, 13171–13174.
28. S. A. Hollingsworth and R. O. Dror, Neuron, 2018, 99, 1129–1143.
29. W. Kohn, A. D. Becke and R. G. Parr, J. Phys. Chem., 1996, 100, 12974–12980.
30. M. Hapka, M. Przybytek and K. Pernal, J. Chem. Theory Comput., 2021, 17, 5538–5555.
31. R. M. Parrish, T. M. Parker and C. D. Sherrill, J. Chem. Theory Comput., 2014, 10, 4417–4431.
32. R. F. W. Bader, Chem. Rev., 1991, 91, 893–928.
33. C. F. Matta and A. A. Arabi, Future Med. Chem., 2011, 3, 969–994.
34. S. Tretiakov, A. Nigam and R. Pollice, Chem. Rev., 2025, 125, 5776–5829.
35. J. P. Foster and F. Weinhold, J. Am. Chem. Soc., 1980, 102, 7211–7218.
36. A. E. Reed, L. A. Curtiss and F. Weinhold, Chem. Rev., 1988, 88, 899–926.
37. M. S. Brozzo, S. Bjelić, K. Kisko, T. Schleier, V.-M. Leppänen, K. Alitalo, F. K. Winkler and K. Ballmer-Hofer, Blood, 2012, 119, 1781–1788.
38. V.-M. Leppänen, A. E. Prota, M. Jeltsch, A. Anisimov, N. Kalkkinen, T. Strandin, H. Lankinen, A. Goldman, K. Ballmer-Hofer and K. Alitalo, Proc. Natl. Acad. Sci. U. S. A., 2010, 107, 2425–2430.
39. V.-M. Leppänen, M. Jeltsch, A. Anisimov, D. Tvorogov, K. Aho, N. Kalkkinen, P. Toivanen, S. Ylä-Herttuala, K. Ballmer-Hofer and K. Alitalo, Blood, 2011, 117, 1507–1515.
40. E. Espinosa, E. Molins and C. Lecomte, Chem. Phys. Lett., 1998, 285, 170–173.

---