Multimodal characterization of a pyrazolo[3,4-b]pyridine aldehyde: quantum chemical insights, NLO activity and VEGFR-2 docking/molecular dynamics

Abstract

The aldehyde-functionalized pyrazolo[3,4-b]pyridine derivative CJ129 was examined in the solid state for the first time using single-crystal X-ray diffraction. The analysis showed that the molecule adopts a monoclinic P2₁/c arrangement, in which the packing is supported by π–π stacking interactions of 3.683 Å between the molecules. To gain a broader picture of its behavior, we combined several quantum-chemical approaches—DFT using M06-2X/def2-TZVPP and TD-DFT with CAM-B3LYP/6-311++G(d,p)—to explore its electronic properties, chemical reactivity, and nonlinear optical (NLO) response. CJ129 features a moderate HOMO–LUMO separation (6.17 eV), a notably high electrophilicity index (ω = 70.31 kcal mol⁻¹), and clear evidence of intramolecular charge transfer, all of which contribute to its strong NLO behavior (β = 26.81 × 10⁻³⁰ esu; γ = 118.34 × 10⁻³⁶ esu), which exceeds the values reported for chalcone-based systems. To complement the electronic analysis, we performed 500-nanosecond molecular dynamics simulations of the complex with VEGFR-2. Throughout the trajectory, CJ129 remained stably positioned in the active site, mainly through hydrophobic interactions and recurring hydrogen bonds involving Asp1046 and Asn923. Taken together, these results provide the first integrated crystallographic and computational description of CJ129 and suggest its potential as an NLO-active molecule and promising scaffold for targeting VEGFR-2 in antiangiogenic research.

---

1. Introduction

Angiogenesis is an essential process for the formation, remodeling, and maintenance of the vascular system, contributing to physiological homeostasis, including embryonic development, wound healing, and tissue regeneration, and to clinically significant pathologies, such as cancer, diabetic retinopathy, chronic inflammation, and metastasis. In this context, the vascular tyrosine kinase receptor VEGFR-2 (Vascular Endothelial Growth Factor Receptor-2) is recognized as the principal mediator of VEGF-induced pro-angiogenic signaling, playing a decisive role in endothelial cell proliferation, migration, and survival. Its activation arises from the binding of VEGF-A, -C, or -D to the receptor's IgD2 and IgD3 domains, promoting dimerization and autophosphorylation at key tyrosine residues (Tyr801, Tyr951, Tyr1054, Tyr1059, Tyr1175, and Tyr1214), which in turn triggers cascades such as PLCγ–PKC–MAPK, SRC–PI3K–Akt, and p38–MAPK. These pathways converge to drive new vessel formation, increase vascular permeability, and remodel microenvironment functions, over which VEGFR-2 exerts dominant control over other members of the receptor family. These mechanisms explain why VEGFR-2 has become a priority therapeutic target in oncology, supporting the development of tyrosine kinase inhibitors such as sunitinib, sorafenib, and axitinib. However, the toxicity, acquired resistance, and limited selectivity of current inhibitors have motivated the search for new compounds that are structurally more efficient and possess tunable electronic profiles.^1,2

The need to expand the repertoire of inhibitors has driven the systematic exploration of nitrogen-containing heterocycles as bioactive scaffolds capable of occupying the VEGFR-2 catalytic pocket with a high affinity. Recent efforts in medicinal chemistry have demonstrated that pyrazoles, triazoles, indoles, quinazolinones, and benzimidazoles can exhibit significant inhibition, with IC₅₀ values in the nano- to micromolar range, supported by π–π interactions and hydrogen bonds within the hinge domain.² For example, Marques et al. reviewed a broad set of heterocyclic VEGFR-2 inhibitors, highlighting pyrazole derivatives with IC₅₀ values of 0.008–0.05 µM and indole–triazole hybrids capable of outperforming sorafenib in vitro. Likewise, various synthetic strategies have enabled the design of molecules with optimized hydrophobic “tails”, flexible urea/thiourea linkers, and rigid aromatic cores that favor selective receptor engagements. These findings underscore the importance of electronic modulation of the central ring, which, together with the nature of the heterocycle, is a key factor in the performance of type I or type II VEGFR-2 inhibitors.³

Among these scaffolds, pyrazolopyridines have emerged as particularly relevant family of compounds due to their ability to act as purine bioisosteres, efficiently occupying the hinge region of various kinases.^4–6 Gheidari et al. conducted an extensive computational study of 37 pyrazolo[3,4-b]pyridine derivatives as potential TRKA inhibitors, reporting docking energies between −12.6 and −14.1 kcal mol⁻¹, with specific hydrogen-bonding interactions involving Glu546, Met620, and Lys627, along with substantial complex stability in 100-ns molecular dynamics simulations. These findings validate the high affinity of this scaffold for catalytic sites enriched in H-bond donor/acceptor residues and demonstrate that its rigid and aromatic nature enables consistent interactions with conserved kinase structures.⁷

In parallel, experimental studies have highlighted the bioactive potential of pyrazolopyridines against VEGFR-2 and other cancer-relevant targets in preclinical models. In a structural review, Donaire-Arias et al. described 1H-pyrazolo[3,4-b]pyridine derivatives with direct inhibitory activity on VEGFR-2, including one compound (“Compound 28”) with an IC₅₀ of 1.46 µM, along with CDK1 inhibition, highlighting the multitarget versatility of this scaffold. At the structural level, these compounds exhibit clear SAR patterns: substituents at C3 and C5 enhance interactions within the hydrophobic region, whereas halogenated groups strengthen binding through σ-hole-type interactions with residues inside the catalytic pocket. The study also reported antiproliferative activity in tumor cell lines sensitive to angiogenesis, consolidating pyrazolopyridines as strong candidates for the development of novel anti-angiogenic agents.⁸

On the other hand, comprehensive mechanistic study, Michailidou et al. developed a series of pyrazolopyridine derivatives as direct inhibitors of VEGFR-2–mediated angiogenesis, demonstrating that several compounds were able to inhibit proliferation, migration, and capillary tube formation in human endothelial cells (HUVECs) under both basal and VEGF-A–stimulated conditions. In particular, derivatives 17a, 18b, and 18c exhibited significant inhibition of VEGFR-2 phosphorylation at Tyr1175, a key site for downstream activation of the ERK1/2 and AKT pathways, without affecting total receptor expression. In murine Lewis lung carcinoma models, these compounds significantly reduced tumor microvessel density and tumor growth, with no evidence of systemic toxicity, providing the first in vivo validation of the pyrazolopyridine scaffold as a functional inhibitor of VEGFR-2-dependent tumor angiogenesis.⁴

Considering these precedents, the present study aims to evaluate the potential of 6-chloro-3-methyl-1-phenyl-1H-pyrazolo[3,4-b]pyridine-5-carbaldehyde as a VEGFR-2 modulator, integrating a multidisciplinary approach that includes (i) detailed electronic analysis through density functional theory (DFT) and evaluation of nonlinear optical (NLO) properties; (ii) structural characterization based on quantum descriptors and comparisons with crystallographic data; (iii) molecular docking and molecular dynamics studies to explore its interaction with VEGFR-2 at the atomic level; and (iv) comparison of its computational properties with examples reported in the literature. This approach not only broadens the understanding of the capabilities of the pyrazolopyridine core in kinase inhibition but also proposes new perspectives for the development of selective VEGFR-2 inhibitors with potential applications in targeted antiangiogenic therapies.

2. Methodology

2.1. Synthesis

The synthetic procedure was carried out following a methodology previously developed and optimized by our research group.⁹ The synthesis was performed in two consecutive stages. First, a mixture of 5-amino-3-methyl-1-phenylpyrazole 1 (1 g, 5.8 mmol) in an AcOH/Ac₂O mixture (10 mL) was refluxed for 3 h. The reaction mixture was cooled, and the resulting solid was filtered off and recrystallized from ethanol to obtain product 2. Acetamide 2 (0.9902 g, 4.6 mmol) was dissolved in DMF (5.0 mL), and the solution was cooled to 0 °C. POCl₃ (2.58 mL, 6 mmol) was added slowly under vigorous stirring at 0 °C. The reaction mixture was then allowed to reach room temperature and refluxed at 80 °C for 8 h. After completion of the reaction, the mixture was cooled and neutralized with a saturated aqueous sodium acetate (NaOAc) solution. The resulting product was extracted with ethyl acetate (EtOAc, 20 mL), washed with water (100 mL), and dried over anhydrous Na₂SO₄. The solvent was evaporated under reduced pressure, and the crude material was purified by column chromatography (30% EtOAc/hexane) to obtain pure CJ129 with a yield of 40% (Scheme 1).

Scheme 1 Synthesis of acetamide and final product (CJ129).

2.2. Computational analyses

2.2.1. Electronic structure calculations and molecular properties. All molecular geometries were optimized using the ORCA 6.0.1 software package¹⁰ and visualized using ChimeraX, incorporating the SEQCROW extension for enhanced molecular analysis.¹¹ The optimized structures were confirmed as true minima on the potential energy surface through harmonic vibrational analysis, with no imaginary frequencies.¹²

Calculations were performed using the M06-2X hybrid meta-GGA functional in combination with the def2-TZVPP basis set, which offers a reliable description of non-covalent interactions, thermodynamic properties, and dispersion effects.¹³ The ionization potential (IP) and electron affinity (EA) were determined exclusively using the adiabatic approach, which accounts for structural relaxation upon electron removal or addition.^12,14,15 The computed vibrational frequencies were subsequently scaled by a factor of 0.945 to account for the systematic anharmonicity and basis set effects inherent to the harmonic approximation.¹⁶

Global reactivity descriptors, including chemical hardness, electronegativity, and electrophilicity index, were derived from the adiabatic IP and EA values according to the definitions proposed by Pearson.¹⁷ Local reactivity was assessed by computing the condensed Fukui functions for both nucleophilic (f⁺) and electrophilic (f⁻) attacks using the finite difference method applied to Hirshfeld atomic charges.^18,19

Spin density distributions were also calculated for the radical anion and radical cation species, providing insights into the localization and delocalization of the unpaired electron and its implications for radical stability.¹⁹

2.2.2. Electronic and nonlinear optical properties. The electronic excitation properties were investigated using time-dependent density functional theory (TD-DFT) at the CAM-B3LYP/6-311++G(d,p) level of theory. The solvent effects were modeled using the polarizable continuum model (PCM) with methanol as the solvent. The first 50 excited states were computed, focusing on the UV-vis absorption region.¹²

The nonlinear optical (NLO) response parameters were evaluated using Gaussian 16, employing the same CAM-B3LYP functional and basis set. The finite field approach was used to calculate the average dipole moment, total linear polarizability, and first- and second-order hyperpolarizabilities using well-established methodologies reported in the literature.^20,21 The (hyper)polarizabilities were calculated using the CPHF/CPKS response formalism implemented in Gaussian (polar = DCSHG, Cubic; CPHF = RdFreq) within the finite field framework. Convergence for the applied field strength was verified through the successful convergence of the linear response equations under strict criteria, as well as by the numerical stability of the polarizability, first-order, and second-order hyperpolarizability tensors (α, β, and γ) obtained at the selected field (ω = 0.0428 a.u.).

2.2.3. Molecular dynamics. Molecular dynamics simulations were performed to evaluate the conformational flexibility, binding stability, and interaction dynamics of CJ129 in the VEGFR-2 active site. To ensure reproducibility, three independent 500 ns trajectories were generated using GROMACS 2023.3.²² Post-simulation analysis employed a Python-based pipeline that integrated the protein–ligand interaction profiler (PLIP) for the systematic characterization of the binding interactions.²³

Molecular system construction: the crystallographic coordinates of VEGFR-2 (PDB ID: 4ASE) were retrieved and preprocessed by eliminating nonstructural water molecules and bound ligands while preserving the crystallographic water molecules relevant to catalysis.²⁴ Protonation states compatible with pH 7.4 were assigned using pdb2gmx software. The protein was parameterized using the CHARMM36m force field. Compound CJ129's topology and parameters were generated using the CGenFF server and subsequently adapted to the GROMACS format using cgenconv. The initial binding pose was obtained using molecular docking calculations performed with SMINA to ensure that the ligand occupied the orthosteric binding pocket.^25,26 The assembled complex was immersed in a dodecahedral solvation box that extended 1.2 nm beyond the protein surface. The box was filled with TIP3P explicit water molecules and ionized with sodium and chloride ions to achieve a physiological ionic strength of 0.15 M NaCl.

Thermodynamic equilibration and production phase: the solvated system underwent energy minimization via the steepest descent algorithm until convergence (maximum force <1000 kJ mol⁻¹ nm⁻¹). Equilibration proceeded through sequential canonical (NVT, 100 ps) and isothermal-isobaric (NPT, 500 ps) ensembles at 310 K and 1 bar. Velocity rescaling was employed for temperature coupling, and the Parrinello–Rahman barostat was used for the pressure control. The production dynamics were extended for 500 ns under NPT conditions (310 K and 1 bar) with a 2 fs integration step. Covalent bonds involving hydrogen were constrained using the LINCS algorithm. Electrostatic interactions beyond 1.0 nm were computed using the Particle Mesh Ewald method. Lennard-Jones and short-range electrostatic terms employed a cutoff radius of 1.0 nm. Coordinate snapshots were recorded at 10 ps intervals for downstream analysis.²²

2.3. Instrumentation used for crystal measurement

Single crystals of CJ129 were grown by slow evaporation from a hexane/ethyl acetate mixture at room temperature. A well-formed crystal was mounted on top of the micromount in a random orientation. Intensity diffraction data were collected at 296 K using a Bruker D8 VENTURE diffractometer equipped with a bidimensional CMOS Photon III detector, employing IµS 3.0 microfocus Mo-Kα radiation (λ = 0.71073 Å). The diffraction frames were integrated using the APEX-4 package and corrected for absorption using SADABS.²⁷ The structure of CJ129 was solved via the Intrinsic Phasing method using ShelXT²⁸ and refined by full-matrix least-squares based on F² with ShelXL.²⁹ All solutions and refinements of the datasets were performed using Olex2 software 1.5.³⁰ Hydrogen atoms were calculated and refined with isotropic thermal parameters using riding coordinates, whereas the rest of the atoms were refined with anisotropic parameters. All CH and CH₂ groups were fixed at Uiso values of 1.2×, and the CH₃ groups were assigned Uiso values of 1.5×. A summary of the crystal data, collection parameters, and refinement details is provided in Table 5, with additional crystallographic information available in the CIF files. ORTEP views were generated using the Olex2 software.³⁰

3. Results and discussion

3.1. DFT analyses

Computational analyses indicated that CJ129 had a HOMO energy of −176.33 kcal mol⁻¹ and an LUMO energy of −34.07 kcal mol⁻¹, resulting in a HOMO–LUMO energy gap (ΔE_H–L) of 142.26 kcal mol⁻¹ (≈6.17 eV). This value reflects a moderate orbital separation, suggesting good thermodynamic stability of the ground state (Table 1).

Table 1 The electronic properties of CJ129, including the HOMO (H) and LUMO (L) energies, HOMO–LUMO energy gap (ΔE_H–L), ionization potential (IP), electron affinity (EA), chemical hardness (η), and electrophilicity index (ω), are expressed in kcal mol⁻¹

Properties	Values
H	−176.33
L	−34.07
ΔEH–L	142.26
IP	182.17
EA	26.84
η	77.67
ω	70.31

---

The ionization potential (IP = 182.17 kcal mol⁻¹) and electron affinity (EA = 26.84 kcal mol⁻¹) result in a chemical hardness of η = 77.67 kcal mol⁻¹, indicating an electronically resilient system. Additionally, the global electrophilicity index (ω = 70.31 kcal mol⁻¹) is notably high, suggesting that the molecule can efficiently behave as an electron-accepting species in intermolecular interactions or under external-field perturbation.

Analysis of the electron density surfaces associated with the frontier molecular orbitals revealed a favorable spatial separation between the HOMO and LUMO. The HOMO is primarily distributed over the pyrazole ring and phenyl moiety, whereas the LUMO is mainly localized on the pyridine fragment and carbonyl group (Fig. 1). This arrangement delineates a clear intramolecular push–pull charge-transfer pathway, which is a hallmark of molecular systems with nonlinear optical (NLO) properties.²⁰

Fig. 1 HOMO and LUMO frontier molecular orbitals (in kcal mol⁻¹).

Such electronic redistribution enhances the ability of the molecule to undergo field-induced polarization, thereby facilitating stronger NLO responses. In practice, this donor–acceptor architecture promotes efficient charge migration upon excitation, directly contributing to the observed increases in the polarizability and hyperpolarizability parameters. This behavior positions the compound as a promising candidate for photonic and optoelectronic applications that require robust charge-transfer dynamics.

The nonlinear optical (NLO) responses of the series were first benchmarked against urea, a classical reference material in organic NLO chemistry (Table 2). As expected, urea exhibited the lowest values in the set, with a first hyperpolarizability β of only 0.41 × 10⁻³⁰ esu and a modest second hyperpolarizability γ of 2.94 × 10⁻³⁶ esu. In contrast, all pyrazolo[3,4-b]pyridine-based systems and chalcone derivatives displayed markedly enhanced NLO parameters, confirming the decisive role of extended π-conjugation and intramolecular charge transfer (ICT) in boosting the molecular response.

Table 2 Nonlinear properties of the studied compound CJ129 at the CAM-B3LYP/6-311++G(d,p) level

Compound	µ (Debye)	α (×10^−24 esu)	β (×10^−30 esu)	γ (×10^−36 esu)
Urea^31	3.88	4.823	0.41	2.94
Ref. 1 (ref. 31)	5.01	50.83	27.46	137.39
Ref. 2 (ref. 32)	2.44	45.89	11.82	3.06
Ref. 3 (ref. 33)	3.01	31.34	8.344	67.87
Ref. 4 (ref. 34)	3.54	27.62	22.41	44.18
This work	2.13	41.25	26.81	118.34

---

Among the pyrazolo[3,4-b]pyridine analogues, ref. 1 and ref. 2 constitute direct structural comparators for the title compound.^31,32 Ref. 1, bearing a highly conjugated benzodioxole–chalcone fragment attached to the pyrazolopyridine core, displayed the largest overall NLO response in the series (β = 27.46 × 10⁻³⁰ esu, γ = 137.39 × 10⁻³⁶ esu). This can be rationalized by the presence of an extended donor–π–acceptor framework, in which the benzodioxole ring acts as an efficient electron donor and the pyrazolopyridine–carbonyl segment behaves as an acceptor, enabling a long-range ICT pathway. Ref. 2 also preserves the pyrazolopyridine scaffold but incorporates a benzofuran unit bearing amino, methoxy, and hydroxy substituents, which modulate the electron density and partially disrupt the symmetry of the system. Consistent with this more complex and less strict linear conjugation, ref. 2 showed intermediate NLO values (β = 11.82 × 10⁻³⁰ esu, γ = 3.06 × 10⁻³⁶ esu), which were clearly higher than those of urea but significantly lower than those of ref. 1.

In this context, CJ129 emerges as a particularly interesting case. Despite having a more compact π-system than that of ref. 1 and lacking strongly electron-donating substituents such as dialkylamino or benzodioxole groups, it exhibits a first hyperpolarizability β = 26.81 × 10⁻³⁰ esu and a second hyperpolarizability γ = 118.34 × 10⁻³⁶ esu, values that are essentially comparable to those of ref. 1 and clearly superior to those of ref. 2. This result indicates that the combination of the pyrazolopyridine core, formyl group at C-5, and chloro-substituted heteroaromatic framework provides an efficient, structurally compact push–pull arrangement. The phenyl–pyrazole fragment behaves as the main electron-rich domain, whereas the pyridine–carbonyl region constitutes the electron-deficient counterpart, as supported by frontier orbital analysis and Fukui indices. The relatively moderate dipole moment of the title compound (2.13 D) compared to that of CJ129 (5.01 D) suggests that the enhancement of β and γ does not rely solely on a large permanent dipole but rather on an efficient field-induced redistribution of charge along the conjugated skeleton.

The comparison with the chalcone-type references (ref. 3 and ref. 4) is particularly revealing, since chalcones are widely recognized as robust NLO building blocks.^33,34 Ref. 3 and 4, both benzodioxole-based chalcones with different π-bridge architectures, present β values of 8.34 × 10⁻³⁰ esu and 22.41 × 10⁻³⁰ esu, and γ values of 67.87 × 10⁻³⁶ esu and 44.18 × 10⁻³⁶ esu, respectively. These numbers confirm their good NLO performance, which is driven by a classical donor–π–acceptor pattern across a flexible conjugated bridge. Notably, the title pyrazolopyridine derivative outperformed both chalcones in terms of the first and second hyperpolarizabilities: its β value exceeded that of compound ref. 4 and was more than three times that of compound ref. 3, while its γ value was almost twice that of compound ref. 3 and nearly three times that of compound 4. This indicates that the rigid, fused heteroaromatic framework of pyrazolopyridine, combined with a strategically placed formyl group and a chloro substituent, can rival—and even surpass—the performance of more extended but less rigid chalcone architectures.

From a molecular design perspective, these trends highlight several key structure–property relationships. First, the pyrazolo[3,4-b]pyridine core behaved as an efficient NLO platform, capable of delivering high β and γ values even when the π-system was relatively compact, provided that an asymmetric distribution of donor and acceptor regions was preserved. Second, the introduction of a formyl group at C-5 in the title compound appears to be crucial for maximizing the ICT toward the carbonyl–pyridine region, which in turn amplifies the hyperpolarizability without necessarily increasing the dipole moment. Third, the comparison between ref. 1/ref. 2 and the title compound suggests that excessive functionalization or overly complex π-extensions do not automatically translate into superior NLO responses; instead, an optimal balance between the conjugation length, electronic asymmetry, and molecular rigidity is required.

Overall, the quantitative comparison with structurally related pyrazolopyridine derivatives (ref. 1 and 2) and high-performing chalcones (ref. 3 and 4) demonstrates that 6-chloro-3-methyl-1-phenyl-1H-pyrazolo[3,4-b]pyridine-5-carbaldehyde is a competitive and, in several aspects, superior candidate for NLO applications. Its high β and γ values, combined with a relatively simple and synthetically accessible scaffold, make this molecule an attractive starting point for the rational design of new pyrazolopyridine-based materials for second- and third-order nonlinear optical devices applications. Although experimental NLO measurements were not performed in this study, the reliability of the predicted hyperpolarizability values was supported by several factors. The calculations were performed using a long-range corrected functional (CAM-B3LYP) combined with the finite field approach, which has been extensively validated for charge-transfer-driven organic NLO systems. Moreover, the relative magnitudes of the calculated β and γ values are consistent with those reported for structurally related pyrazolopyridine and chalcone derivatives studied experimentally, reinforcing the predictive value of the present computational analysis.

On the other hand, the local reactivity of the system was examined using the dual Fukui functions (f⁺ and f⁻), which identify the preferred regions for nucleophilic and electrophilic attack, respectively (Table 3 and Fig. 2). As shown in Fig. 2, the blue isosurfaces corresponding to the nucleophilic Fukui function (f⁺) were predominantly located around the formyl fragment, with atoms C16, C18, and O19 being the most reactive electrophilic centers (f⁺ = 0.13, 0.09, and 0.10, respectively). This marked accumulation of electrophilic character aligns with the electron-accepting function of the carbonyl-pyridine region in the push–pull architecture of the molecule.

Table 3 Fukui index and spin densities ρ_spin of CJ129

Atom number	ρspin (RA)	ρspin (RC)	(f^+)	(f^−)
1	0.00	0.00	0.02	0.04
2	0.00	0.19	0.03	0.09
3	0.00	0.00	0.02	0.04
4	0.00	0.10	0.01	0.05
5	0.00	0.10	0.00	0.04
6	0.00	0.09	0.01	0.05
7	0.00	0.23	0.03	0.06
8	0.07	0.00	0.03	0.02
9	0.03	0.04	0.02	0.03
10	0.06	0.12	0.02	0.07
11	0.06	−0.02	0.05	0.03
12	0.00	0.01	0.01	0.02
13	0.06	0.05	0.05	0.04
14	−0.01	−0.01	0.01	0.02
15	0.09	0.05	0.05	0.04
16	0.31	0.00	0.13	0.02
17	0.00	0.00	0.08	0.05
18	0.13	0.00	0.09	0.02
19	0.14	0.02	0.10	0.05

---

Fig. 2 Isosurfaces of the Fukui dual descriptors of CJ129. The yellow regions correspond to the electrophilic Fukui function (f⁻), and the blue areas correspond to the nucleophilic Fukui function (f⁺).

Conversely, the yellow isosurfaces in Fig. 2 represent regions with high values of the Fukui electrophilic function (f⁻), identifying C2 (0.09), N7 (0.06), and C10 (0.07) as the main nucleophilic sites capable of donating electron density. This clear spatial segregation between donor- and acceptor-rich regions reflects a strong intrinsic electronic asymmetry, which defines the potential intramolecular charge transfer (ICT) pathways of the molecule. These pathways are crucial for nonlinear optical (NLO) responses, as they enhance the ability of the system to undergo field-induced polarization.

Spin density analyses of the radical cation and anion further support this interpretation. As illustrated in Fig. 3, the α- and β-spin isosurfaces show distinct localization patterns depending on the species charge. In the radical cation (Fig. 3a), the unpaired electron density is mainly concentrated on the pyrazole ring, indicating that this domain functions as the primary electron reservoir under oxidative conditions. In contrast, the radical anion (Fig. 3b) exhibited spin accumulation in the carbonyl and pyridine regions, revealing that the negative charge was preferentially stabilized by the electron-deficient portion of the molecule.

Fig. 3 Spin density isosurfaces of the radical cation (a) and radical anion (b) of CJ129. The blue and purple regions represent the spatial distributions of the α- and β-spin densities, respectively, illustrating the localization of the unpaired electron after single-electron reduction.

This asymmetric spin-density redistribution between the oxidized and reduced states confirmed the presence of a strong electronic gradient along the molecular framework. Such behavior is highly desirable in NLO materials, as the ability to stabilize charge-separated electronic states is directly associated with enhanced hyperpolarizability and efficient ICT under external electric-field stimulation. These characteristics also contribute to improved performance in optoelectronic processes, including harmonic generation, optical switching and photonic signal modulation.

Taken together, the results depicted in Fig. 2 and 3 demonstrate that the architecture of 6-chloro-3-methyl-1-phenyl-1H-pyrazolo[3,4-b]pyridine-5-carbaldehyde promotes (i) a well-defined separation between electron-donating and electron-accepting sites, (ii) robust ICT pathways, and (iii) efficient stabilization of the charged and excited states. This combination of electronic asymmetry, directional charge-transfer capability, and spin density delocalization is a hallmark of high-performance organic NLO materials, further supporting the potential of this compound as a competitive candidate for advanced optoelectronic applications.

In contrast, the spectroscopic simulation in methanol predicted three main electronic transitions at S₃, S₄, and S₇, corresponding to 273, 249, and 222 nm, respectively (Fig. 4 and Table 4). These transitions are dominated by HOMO → LUMO (95%), HOMO → LUMO + 1 (91%), and HOMO − 1 → LUMO (58%) contributions, with significant oscillator strengths (f = 0.42–0.80), confirming their allowed nature. These computational results closely match the experimental UV-vis spectrum, which exhibits intense absorption bands in the 250–280 nm region, characteristic of the π → π* transitions. The strong agreement between theory and experiment reinforces the reliability of the electronic model used and supports the structural assignment of this compound. Furthermore, the intensity and positioning of the absorptive bands highlight the molecule's capability to respond efficiently to optical stimuli, which is an essential feature for materials designed for optical switching, photonic sensing, or other light-driven functional applications.

Fig. 4 (a) Simulated UV-vis spectra and (b) experimental UV-vis spectrum of the CJ129 compound.

Table 4 Main electronic transitions calculated in methanol

Excited state	Excitation vertical (nm)	Ecalc^a (eV)	fos^b	Major contribution (%)	λabs^c (nm)
a Calculated transition energy.b Oscillator strength.c Experimental absorption maxima. H: HOMO, L: LUMO.
S3	273	4.54	0.42	H + 0 → L + 0 (95%)	269
S4	249	4.97	0.80	H + 0 → L + 1 (91%)	—
S7	222	5.58	0.20	H − 1 → L + 0 (58%)	223

---

Finally, the experimental IR spectrum was compared with the DFT-simulated spectrum calculated at the M06-2X/def2-TZVPP level (Fig. 5). The vibrational frequencies were uniformly scaled by a factor of 0.945, as recommended for this functional and basis set, to correct for the systematic overestimation inherent in DFT calculations.¹⁶

Fig. 5 Comparison of the experimental (black) and DFT-calculated (blue) IR spectra of the studied compound.

A strong band was observed at approximately 1700 cm⁻¹, consistent with the C=O stretching of the aldehyde group. Likewise, the signals at 2913 and 2969 cm⁻¹ correspond to the symmetric and asymmetric C–H sp³ stretching modes, respectively, in line with the presence of the methyl group attached to the pyrazole ring. Above 3000 cm⁻¹, the characteristic bands of C–H sp² stretching associated with the aromatic framework are clearly distinguished. Overall, excellent qualitative agreement was observed between the experimental and calculated spectra of the samples. The small discrepancies in band positions and intensities between the experiment and theory can be explained by the different environments studied experimentally and computationally, as well as by the level of quantum-chemical treatment employed.

Taken together, the agreement between the experimental and calculated spectra supports the correct structural assignment of the compound and reinforces the reliability of the employed computational model, providing confidence in both the vibrational description and the resulting structural interpretation.

3.2. Crystal structure description

A suitable single crystal of CJ129 was analyzed at room temperature, and its molecular structure is depicted in Fig. 6. The organic compound crystallized in the non-centrosymmetric monoclinic space group P2₁/c with one molecule per asymmetric unit (Table 5). The crystal structure revealed a predominant π⋯π stacking interaction (3.683 Å, slipping angle 22.56(6)°) that propagated along the crystallographic c-axis (Fig. 6). These interactions form one-dimensional chains with an ABAB molecular arrangement (Fig. 7). All the atoms in the asymmetric unit occupied the general Wyckoff position 4e, and no atoms were found to occupy special positions.

Fig. 6 (a) View of the asymmetric unit of CJ129 with atom numbering scheme. The atomic displacement ellipsoids are drawn at the 50% probability level. (b) Perspective view of the CJ129 crystal packing highlighting the π⋯π interactions propagating along the c-axis (yellow dashed lines).

Table 5 Crystal data, data collection, and structure refinement for the CJ129 compound

Compound	CJ129
Empirical formula	C14H10ClN3O
Formula weight (g mol^−1)	271.70
Temperature (K)	296
Crystal system	Monoclinic
Space group	P21/c
a (Å)	9.2405(2)
b (Å)	20.5457(5)
c (Å)	6.9132(2)
α (°)	90
β (°)	109.941(1)
γ (°)	90
Volume (Å^3)	1233.80(5)
Z	4
ρcalc g (cm^3)	1.463
µ (mm^−1)	0.30
F (000)	560.0
Crystal size (mm^3)	0.12 × 0.058 × 0.015
Radiation	Mo Ka (λ = 0.71073)
2θ Range for data collection (°)	2.34 to 26.37
Index ranges	−11 ≤ h ≤ 11/−25 ≤ k ≤ 25/−8 ≤ l ≤ 8
Reflections collected	23 955
Independent reflections	2156 [Rint = 0.0445]
Data/restraints/parameters	2520/0/173
Goodness-of-fit on F^2	1.099
Final R indexes [all data]	R^1 = 0.0403/wR^2 = 0.1182
Wyckoff position (all atoms)	4e
Site symmetry (all atoms)	1
Largest diff. peak/hole/e Å^−3	0.25/−0.26

---

Fig. 7 Chains of CJ129 projected on the ab-plane, showing its ABAB order along the c-axis.

3.3. Molecular dynamic simulations

The 500-ns molecular dynamics simulation revealed the remarkable structural stability of the VEGFR-2/compound CJ129 complex (Fig. 8A and B). The root-mean-square deviation (RMSD) of the protein backbone stabilized at 2.90 ± 0.21 nm (29.0 ± 2.1 Å) after removing the first 10 ns of equilibration (the first 500 frames). This stabilization indicates that the complex rapidly reached a conformationally stable state without exhibiting significant structural transitions during the production stage of the simulation. Meanwhile, the ligand RMSD showed minimal fluctuations, averaging 0.27 ± 0.06 nm (2.7 ± 0.6 Å). These low-amplitude deviations confirmed that CJ129 remained firmly anchored within the binding pocket without undergoing meaningful reorientation. This demonstrates a highly stable and thermodynamically favorable binding mode.

Fig. 8 Structural stability and interaction profile of the VEGFR-2–ligand complex. (A) Average profile of noncovalent interactions per residue between the ligand and key residues of the catalytic domain. The bars indicate the mean number of contacts by interaction type: hydrophobic (orange), hydrogen bonding (blue), π–π stacking (green), and cation–π (red). (B) Time evolution of the root-mean-square deviation (RMSD) for the protein backbone (pink) and ligand heavy atoms (green) over a 500-ns molecular dynamics simulation. (C) Root-mean-square fluctuation (RMSF) per residue plot (D) Three-dimensional representation of the final conformation of the VEGFR-2 (PDB ID: 4ASE)–ligand complex after 500 ns simulation.

The residue-wise RMSF profile revealed the dynamic characteristics of the kinase proteins (Fig. 8C), displaying a rigid catalytic core and flexible peripheral regions. The analysis yielded an average RMSF of 1.07 ± 0.78 nm (10.7 ± 7.8 Å), with most residues fluctuating within the 0.4–1.2 nm (4–12 Å) range. The binding-site residues, including Asp1046, Ala1050, and Asn923, exhibited fluctuations below 2 nm (20 Å). This indicates moderate mobility and suggests that the interaction with CJ129 does not substantially alter the local architecture of the catalytic pocket. This localized stability contributes directly to the overall conformational robustness of the complex. A set of 18 residues contributed to the stabilization of the complex. The predominant stabilizing mechanism was hydrophobic contact. Ile888, which formed approximately 460 contacts, emerged as the primary hydrophobic node. This created a favorable apolar environment that accommodated the ligand's aromatic core within the catalytic pocket. These interactions play a central role in molecular recognition and ligand affinities. In contrast, hydrogen bonding interactions presented an average of 18 ± 35 contacts per residue, exhibiting variable occupancy patterns. Asp1046 and Asn923 formed highly persistent hydrogen bonds (>85%), indicating their fundamental role in stabilizing the active site. However, Glu885 and Ser884 exhibited more dynamic interactions with an occupancy of ∼40%, likely associated with fine modulation of local flexibility and conformational adjustment during ligand-binding.

Although less frequent, cation–π interactions display notable specificity. A persistent cation–π contact involving Lys868 was observed, which may be associated with ligand selectivity. The residues critical for stabilizing the complex, Asp1046, Ala1050, Ile888, Glu885, Lys868, and Asn923, are located in proximity to the catalytic site of VEGFR-2 and the functional core formed by Cys919, Lys868, and Glu885 (Fig. 8D). The interaction network in this region suggests that the ligand stabilizes the inactive kinase conformation, preventing the catalytically competent configuration required for substrate phosphorylation.

Importantly, the strong nonlinear optical (NLO) response of CJ129 is not an isolated physicochemical feature but is intrinsically connected to the same electronic properties that govern its interaction with biological targets. The pronounced intramolecular charge transfer (ICT), evidenced by the spatial separation of the HOMO and LUMO, high electrophilicity index, and enhanced first- and second-order hyperpolarizabilities, reflects a highly polarizable electronic framework capable of efficient charge redistribution under external perturbations. In a biological context, this electronic flexibility is directly relevant to protein–ligand recognition, as it facilitates adaptive electrostatic complementarity within the VEGFR-2 binding pocket.

Specifically, the same donor–acceptor architecture that amplifies the NLO response enhances the ability of CJ129 to stabilize noncovalent interactions, such as hydrogen bonding, cation–π contacts, and hydrophobic packing, within the kinase active site. The carbonyl–pyridine fragment, identified as the primary electron-accepting region by frontier orbital and Fukui analyses, participates in persistent hydrogen bonding with key catalytic residues (Asp1046 and Asn923), whereas the electron-rich pyrazole–phenyl domain favors π–π and hydrophobic interactions. This dual electronic behavior mirrors the push–pull mechanism responsible for the NLO enhancement.

From a broader perspective, these results suggest that the NLO descriptors reported here can be interpreted as quantitative proxies for electronic softness, polarizability, and charge-transfer efficiency, which are increasingly recognized as critical for kinase inhibition. Therefore, the strong NLO performance of CJ129 is not only relevant for optoelectronic applications but also provides mechanistic insights into its favorable antiangiogenic profile as a VEGFR-2 modulator. This dual functionality highlights the value of integrating NLO analysis into the rational design and pre-screening of bioactive heterocycles that target angiogenesis-related pathways.

4. Conclusions

The present study demonstrates that 6-chloro-3-methyl-1-phenyl-1H-pyrazolo[3,4-b]pyridine-5-carbaldehyde (CJ129) is a paradigmatic example of how a compact, electronically well-balanced, and structurally rigid heterocyclic scaffold simultaneously integrates advanced physicochemical properties and biologically relevant behavior. Crystallographic characterization confirmed a highly ordered architecture with efficient π–π interactions, while DFT analyses revealed a clear donor–π–acceptor electronic segregation, responsible for its high global electrophilicity, radical stability, and pronounced nonlinear optical response. In this regard, the hyperpolarizability values obtained position CJ129 as a competitive NLO material when compared to more extended reference systems, providing key insight that NLO efficiency does not depend exclusively on conjugation length, but rather on proper intramolecular electronic organization.

In contrast, molecular dynamics simulations suggested that CJ129 exhibited a robust and consistent affinity for the catalytic domain of VEGFR-2, with a stable binding orientation throughout the 500 ns trajectory. The interaction pattern is dominated by hydrophobic contacts and reinforced by highly persistent hydrogen bonds. This pattern supports a stable and energetically favorable association that does not disrupt the structural integrity of the catalytic pocket. This cooperative network of interactions suggests an inhibition mechanism that stabilizes an inactive kinase conformation. Although these findings are based on computational modeling and require experimental validation, the overall dynamic and energetic profiles position CJ129 as a promising modulator of VEGFR-2 activity and a promising structural scaffold for developing novel antiangiogenic agents.

Conflicts of interest

There are no conflicts to declare.

Data availability

CCDC 2479116 (CJ129) contains the supplementary crystallographic data for this paper.³⁵

References

1. X. Wang, A. M. Bove, G. Simone and B. Ma, Front. Cell Dev. Biol., 2020, 8, 599281.
2. P. J. Chaudhari, A. R. Nemade and A. A. Shirkhedkar, RSC Adv., 2024, 14, 33384–33417.
3. C. S. Marques, P. Brandão and A. J. Burke, Molecules, 2024, 29, 5341.
4. M. Michailidou, V. Giannouli, V. Kotsikoris, O. Papadodima, G. Kontogianni, I. K. Kostakis, N. Lougiakis, A. Chatziioannou, F. N. Kolisis, P. Marakos, N. Pouli and H. Loutrari, Eur. J. Med. Chem., 2016, 121, 143–157.
5. W. M. Eldehna, M. H. Abdulla, M. S. Nafie, A. E. Elsawi, S. Ayman, M. I. Shahin, N. S. Alhassan, A. M. Zubaidi, H. A. Ghabbour, M. Elaasser, A. A. Al-Karmalawy and H. A. Abdel-Aziz, Bioorg. Chem., 2024, 153, 107778.
6. S. M. Alhamaky, N. A. Khalil, A. K. A. Bass, N. Osama and M. S. A. Hassan, Drug Dev. Res., 2025, 86, e70056.
7. D. Gheidari, M. Mehrdad, K. jani and A. Feizkhah, Sci. Rep., 2025, 15, 27367.
8. A. Donaire-Arias, A. M. Montagut, R. P. de la Bellacasa, R. Estrada-Tejedor, J. Teixidó and J. I. Borrell, Molecules, 2022, 27, 2237.
9. E. Polo-Cuadrado, K. Ferrer, J. Sánchez-Márquez, A. Charris-Molina, Y. A. Rodríguez-Núñez, L. Espinoza-Catalán and M. Gutiérrez, Heliyon, 2024, 10, e32573,  DOI:10.1016/j.heliyon.2024.e32573.
10. F. Neese, Wiley Interdiscip. Rev.:Comput. Mol. Sci., 2012, 2, 73–78.
11. E. C. Meng, T. D. Goddard, E. F. Pettersen, G. S. Couch, Z. J. Pearson, J. H. Morris and T. E. Ferrin, Protein Sci., 2023, 32, e4792.
12. E. Polo-Cuadrado, E. Osorio, K. Acosta-Quiroga, P. A. Camargo-Ayala, I. Brito, J. Rodriguez, J. B. Alderete, O. Forero-Doria, E. F. Blanco-Acuña and M. Gutiérrez, RSC Adv., 2024, 14, 10199–10208.
13. R. G. Parr, R. A. Donnelly, M. Levy and W. E. Palke, J. Chem. Phys., 1978, 68, 3801–3807.
14. K. Acosta-Quiroga, E. Rocha-Valderrama, M. Zúñiga-Bustos, R. Mera-Adasme, G. Cabrera-Barjas, C. Olea-Azar and M. Moncada-Basualto, Antioxidants, 2025, 14, 479.
15. E. Polo-Cuadrado, K. Acosta-Quiroga, C. Rojas-Peña, Y. A. Rodriguez-Nuñez, E. F. Blanco-Acuña, J. J. Lopez, I. Brito, J. Cisterna, J. B. Alderete and M. Gutiérrez, RSC Adv., 2025, 15, 8609–8621.
16. I. M. Alecu, J. Zheng, Y. Zhao and D. G. Truhlar, J. Chem. Theory Comput., 2010, 6, 2872–2887.
17. R. G. Pearson, J. Chem. Sci., 2005, 117, 369–377.
18. P. Fuentealba, E. Florez and W. Tiznado, J. Chem. Theory Comput., 2010, 6, 1470–1478.
19. E. Osorio, E. G. Pérez, C. Areche, L. M. Ruiz, B. K. Cassels, E. Flórez and W. Tiznado, J. Mol. Model., 2013, 19, 2165–2172.
20. K. Acosta-Quiroga, E. Polo-Cuadrado, M. J. Percino, E. Blanco-Acuña, D. Villaman, E. Pérez-Gutiérrez, M. E. Patiño, J. B. Alderete, R. Gonzalez, J. Saavedra-Olavarría, E. G. Pérez, C. Olea-Azar, M. Moncada-Basualto and C. Rojas-Peña, J. Mater. Chem. C, 2025, 13, 9083–9098.
21. C. Rojas-Peña, E. Polo-Cuadrado, D. Villaman, A. Liempi, M. Zúñiga-Bustos, F. Martínez-Gómez, J. Saavedra-Olavarría, E. G. Pérez, J. D. Maya, C. Olea-Azar and M. Moncada-Basualto, Spectrochim. Acta, Part A, 2026, 347, 126984.
22. M. Abraham, A. Alekseenko, C. Bergh, C. Blau, E. Briand, M. Doijade, S. Fleischmann, V. Gapsys, G. Garg, S. Gorelov, G. Gouaillardet, A. Gray, M. E. Irrgang, F. Jalalypour, J. Jordan, C. Junghans, P. Kanduri, S. Keller, C. Kutzner, J. A. Lemkul, M. Lundborg, P. Merz, V. Miletić, D. Morozov, S. Páll, R. Schulz, M. Shirts, A. Shvetsov, B. Soproni, D. van der Spoel, P. Turner, C. Uphoff, A. Villa, S. Wingbermühle, A. Zhmurov, P. Bauer, B. Hess and E. Lindahl, GROMACS 2023.3 Manual (2023.3), Zenodo, 2023,  DOI:10.5281/ZENODO.10017699.
23. P. Schake, S. N. Bolz, K. Linnemann and M. Schroeder, Nucleic Acids Res., 2025, 53, W463–W465.
24. M. McTigue, B. W. Murray, J. H. Chen, Y. L. Deng, J. Solowiej and R. S. Kania, Proc. Natl. Acad. Sci. U. S. A., 2012, 109, 18281–18289.
25. O. Trott and A. J. Olson, J. Comput. Chem., 2010, 31, 455–461.
26. J. Eberhardt, D. Santos-Martins, A. F. Tillack and S. Forli, J. Chem. Inf. Model., 2021, 61, 3891–3898.
27. G. Sheldrick, SADABS, Program for Empirical Absorption Correction of Area Detector Data, 1996.
28. G. M. Sheldrick, Acta Crystallogr., Sect. A:Found. Adv., 2014, 71, 3–8.
29. G. M. Sheldrick, Acta Crystallogr., Sect. C:Struct. Chem., 2015, 71, 3–8.
30. O. V. Dolomanov, L. J. Bourhis, R. J. Gildea, J. A. K. Howard and H. Puschmann, J. Appl. Crystallogr., 2009, 42, 339–341.
31. E. Polo-Cuadrado, K. Ferrer, E. Osorio, I. Brito, J. Cisterna, L. Espinoza, J. B. Alderete and M. Gutiérrez, RSC Adv., 2023, 13, 5197–5207.
32. S. A. Halim and M. A. Ibrahim, RSC Adv., 2022, 12, 13135–13153.
33. S. Shettigar, P. Poornesh, G. Umesh, B. K. Sarojini, B. Narayana and K. Prakash Kamath, Opt. Laser Technol., 2010, 42, 1162–1166.
34. Y. Sheena Mary, K. Raju, C. Y. Panicker, A. A. Al-Saadi, T. Thiemann and C. Van Alsenoy, Spectrochim. Acta, Part A, 2014, 128, 638–646.
35. CCDC 2479116: Experimental Crystal Structure Determination, 2026,  DOI:10.5517/ccdc.csd.cc2p6qh4.

---