Insights into the coordination chemistry of antineoplastic doxorubicin with 3d-transition metal ions Zn2+, Cu2+, and VO2+: a study using well-calibrated thermodynamic cycles and chemical interaction quantum chemistry models

Before delving into the analysis of each family of doxorubicin-metal complexes investigated, it is important to highlight that the computed entropy contribution to the solution phase Gibbs free energy of reaction displayed values close to zero for all the studied chemical processes (computed values were lower than 0.1 kcal·mol−1·K−1). This can be corroborated by comparing the results in Tables 2, 4, and 6 with those in Tables S2, S3, and S4 of the Supplementary Material, respectively. Nonetheless, it should be noted that the contribution of the total entropic term \(\left( \right)\) to the Gibbs free energy of reaction in the solution phase ranges between −6 and −18 kcal·mol−1, and it certainly contributes to the exergonicity of the studied reactions.

Cu(II)-doxorubicin complexes

Among the three families of metal-complexes investigated in this work, those comprised by Cu(II) ions have been extensively studied. Both, the 1:1 and 1:2 stoichiometries are known to be stable in aqueous solution, but there is some controversy about which one is predominant at near-neutral pH conditions. Dabrowiak and Greenaway observed the 1:2 Cu(II):DOX complex over a range of pH values from 7 to 8.2 [34]; conversely, Garnier-Suillerot et al. observed that the 1:1 species is the one predominant at a neutral pH condition, although none of these studies provided information about the corresponding thermodynamic equilibrium constants [95, 96]. In any case, the square planar coordination geometry was considered as the most prone to be observed in solution phase. It is known that, in general, Cu(II) complexes may display different coordination numbers, from four to six [97,98,99]. For instance, by UV–vis titration experiments, Salmon et al.[100] found that the hexahydrate \(\left[ }\left( }_ }} \right)_6 } \right]^\) complex is the most abundant in aqueous solution, while Tabouli et al.[101] reported the \(\left[ }\left( }_ } \right)_5 } \right]^\) (square pyramid) and \(\left[ }\left( }_ } \right)_6 } \right]^\) (octahedral) species, as the most stable coordination geometries for cooper-ammonia complexes, using a quantum chemistry theoretical protocol. We thus consider necessary to explore all the coordination numbers in our computations, for both the free Cu(II) (from tetrahydrate to hexahydrate complexes, m = 4, 5, 6) and the \(\left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^\) species, to find the chemical process which provides the lowest Gibbs free energy of complexation. An equivalent exploration will be performed for the remainder of the metal complexes under investigation.

Table 1 presents the results of the calibration study (Eq. (1)) for copper-doxorubicin metal-complexes. In addition to the electronic energy variations, we also report the relative reaction energies (in parentheses) referenced to the lowest energy complexation process, to facilitate comparison across different density functional approximations (DFAs). Remarkably, all the functionals here tested predicted the same trend as our DLPNO-CC standard; octahedral Cu(II)-DOX complexes (l = 4) are the most stable in gas phase, followed by the pentahedral (l = 3) and the tetrahedral (l = 2). For the octahedral and tetrahedral complexes, the HDOX molecule is preferable attached to the metallic center from the beta side, while for the pentahedral, the alpha orientation is preferred. Formation of complexes was generally enhanced when cooper tetrahydrate (square planar, (m = 4)) was used as starting material. It is worth noting that in some cases, our octahedral geometries evolved towards the penta-coordinated species after the geometry optimization process, due to a water molecule shifting towards the first solvation sphere. This was the case of the \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) and the \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) species, particularly when their geometries were optimized with the PBE and PBE0 DFAs, respectively. The PBE0 functional displayed the closest reaction energies to our DLPNO-CC standard, with an average deviation (AD) of 14.21 kcal·mol−1 followed by ωB97X-D (AD = 21.23 kcal·mol−1), M06 (AD = 24.41 kcal·mol−1) and PBE (AD = 43.93 kcal·mol−1), while the M06 DFA displayed the best (closest to our standard) relative energy differences (quantities in parentheses in Table 1) with an average deviation of 9.71 kcal·mol−1, followed by ωB97X-D (12.92 kcal·mol−1), PBE0 (18.83 kcal·mol−1) and PBE (20.40 kcal·mol−1). Therefore, results from the PBE0 and M06 functionals will be prioritized in the analysis of the complexation process of copper-doxorubicin species in solution phase. Nonetheless, it is important to acknowledge that even with functionals that exhibit low averaged deviations, there may still be cases where the predicted relative stability trend deviates from our CC standard. For instance, the M06 functional predicts a higher stability for the alpha pentahedral complex compared to its beta counterpart. Additionally, in our calculations using the PBE0 DFA, we encountered difficulties in obtaining a minimum energy structure for the \(\left[ }\left( }} \right)\left( }_ }} \right)_l } \right]^+\) complex, which is the second most stable species predicted by our CC protocol. Thus, for the analysis of chemical processes in the solution phase, we always prioritized the use of at least two different density functional approximations (DFAs), which enriched the discussion and provided more robust insights. Figure 4a, b, c–f show the BCPs and the NCII plots for the \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) (octahedral) and \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]^+\) (square planar) coordination geometries, respectively. The most favored species displays two OH-O interactions (green dots, Fig. 4a) against the only one observed for the square planar complex (Fig. 4c). The 3D-NCII surfaces calculated for both conformers (Fig. 4c, d), conjointly with the 2D-NCII maps (Fig. 4e, f) indicates the presence of more attractive strong non-covalent interactions (blue surfaces) for the \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) species and are mainly located at the coordination sphere for both compounds. This suggest that the metal–ligands attractive interactions may be stronger for the octahedral than for the tetrahedral complex.

Table 1 Computed energy variation corresponding to the reaction between \(\left[ }\left( }_ }} \right)_m } \right]^\) and HDOX ligand, forming \(\left[ }\left( }} \right)\left( }_ }} \right)_l } \right]^+\) complexes, in gas phase, using Eq. (11) relative values respect the most stable complexation process are reported in parenthesesFig. 4figure 4

a, b BCPs, c, d 3D-NCII plots and e, f 2D-NCII maps obtained for the \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]^+\) and \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) species, respectively, at the ωB97X-D/DGDZVP level of theory. Only BCPs between H2O molecules and DOX− are plotted, while BCPs between the metallic center and DOX− or H2O molecules were replaced by turquoise solid lines. “s” is the reduced gradient of the electron density

In Table 2 we report the Gibbs free energy of reaction in solution phase for the Cu(II)-DOX complexes. One may observe that both 1:1 and 1:2 stoichiometric species are prone to (co)exist in aqueous solution at the specified pH conditions (all the tested reactions, analyzed through the set of DFAs under consideration, were exergonic). Furthermore, these results indicate that the representative species at both stoichiometries may have a square planar geometry (l = 2 for n = 1 and l = 0 for n = 2). The formation of the former (1:1 stoichiometry) is particularly favored when the \(}^-\) species is attached to the metal from the beta side, \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]^+\), and also, if the free Cu(II) ion is found in its hexahydrate (m = 6) coordination geometry prior to the complexation process. It is worthy to recall that precisely the Cu(II) hexahydrate is the most abundant in aqueous solution. For the 1:2 complexes, results obtained with the PBE and PBE0 functionals indicate that the beta orientation is the one preferred by the HDOX molecule, \(\left[ }\left( }} \right)_2 } \right]\), while results from the M06 and the ωB97X-D DFAs rather suggest that is the alpha orientation the one attached to the metallic center, \(\left[ }\left( }} \right)_2 } \right]\). In any case, the Gibbs free energy difference between both orientations is less than 4 kcal·mol−1 and conjointly, these results suggest that for the 1:2 stoichiometry situation, both are prone to coexist (the reaction energy difference between the process with m = 5 and m = 6, using the M06 functional, is less than 0.3 kcal·mol−1). Results obtained with the PBE0 DFA, suggest that the 1:1 stoichiometry is slightly more stable than the 1:2 (by barely 4 kcal·mol−1). However, it is important to note that experimental validation is necessary to confirm this prediction, for instance, it is crucial to consider the potential for a higher overestimation of the loss of reaction entropy in the formation of the 1:2 complex. Figure 5a shows the geometry and BCP of the \(\left[ }\left( }} \right)_2 } \right]\) species. The complex is stabilized by five intermolecular BCPs, with one located on the left and four on the right-hand side. The 3D-NCII analysis shown in Fig. 5b supports this observation and shows that, both weak (green surface) and/or strong (blue surfaces) attractive non-covalent interactions are abundant in these molecular moieties. Interestingly, despite the size of the doxorubicin ligands, the 2D-NCII plot (Fig. 5c) indicate that repulsive interactions are the less abundant for this complex.

Table 2 Solution phase Gibbs free energy of complexation between \(\left[ }\left( }_ }} \right)_m } \right]^\) and n HDOX ligands, forming \(\left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^\) complexesFig. 5figure 5

a BCPs, b 3D-NCII and c 2D-NCII plots obtained for the \(\left[ }\left( }} \right)_2 } \right]\) species at the ωB97X-D/ DGDZVP level of theory. Only BCPs between H2O molecules and DOX− are plotted, while BCPs between the metallic center and DOX− or H2O molecules were replaced by turquoise solid lines. “s” is the reduced gradient of the electron density

Zn(II)-doxorubicin complexes

Zn(II) species are commonly found in their tetra, penta and/or hexa coordination modes in aqueous solution, for instance, the most abundant coordination geometries for the zinc-aquo complexes in solution phase are the tetrahedron and the octahedron, while the pentahydrate is practically absent [102,103,104,105,106,107,108,109,110,111,112,113]. Jabłońska-Trypuć et al. [37] reported that under near-neutral pH conditions, the Zn(II)-DOX complexes are predominantly found in a 1:2 stoichiometry, based on their partial molar ratio measurements. No further reports were found about the formation of these species. Therefore, like the analysis of Cu(II)-DOX complexes, we will explore all coordination patterns for both hydrated zinc and zinc-doxorubicin species, to propose plausible complexation processes in aqueous solution.

Our theoretical gas-phase benchmarking for zinc complexes (Table 3) reveals that the PBE0 functional provides the closest results to our DLPNO-CCSD standard (AD = 6.24 kcal·mol−1), followed by ωB97X-D (AD = 13.09 kcal·mol−1), M06 (AD = 14.65 kcal·mol−1) and PBE (AD = 21.61 kcal·mol−1); while, the M06 functional displays the best relative energy differences with an average deviation of 3.90 kcal·mol−1, followed by ωB97X-D (AD = 8.64 kcal·mol−1), PBE0 (AD = 12.59 kcal·mol−1) and PBE (AD = 16.76 kcal·mol−1). Thus, as previously considered for Cu(II)-DOX species, for zinc-doxorubicin metal complexes results obtained from PBE0 and M06 functionals will be prioritized in our analysis of results computed in solution phase. Octahedral Zn(II)-DOX geometries were predicted to be the most stable 1:1 stoichiometric complex under gas phase conditions, particularly if the \(}^-\) species is attached to the metal from the beta side. Also note that the spontaneity of these reactions is enhanced when the tetracoordinated Zn(II)-aquo species is considered as reactant (m = 4). In Fig. 6a, b we show the BCPs of the \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) and the \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]^+\) complexes, respectively, for comparison purposes. The complex with the higher stability displays two OH–O bonds between the H2O molecules that coordinates the metal center and the O atoms of the DOX− species. Conversely, the less stable tetrahedral complex only exhibited one of these OH-O bonds. Likewise, both the 3D-NCII plots (Fig. 6c, d) and the 2D-NCII maps (Fig. 6e, f) reveal that the attractive non-covalent interactions (blue and green surfaces) are substantially less abundant in the less stable complex, with particular emphasis on the region nearby the coordination sphere.

Table 3 Computed energy variation corresponding to the reaction between \(\left[ }\left( }_ }} \right)_m } \right]^\) and HDOX ligand, forming the \(\left[ }\left( }} \right)\left( }_ }} \right)_l } \right]^+\) complex, in gas phase, using Eq. (11)Fig. 6figure 6

a, b BCPs, c, d 3D-NCII plots and e, f 2D-NCII maps obtained for the \(\left[ }\left( }} \right)\left( }_ }} \right)_4 } \right]^+\) and \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]^+\) species, respectively, at the ωB97X-D/ DGDZVP level of theory. Only BCPs between H2O molecules and DOX− are plotted, while BCPs between the metallic center and DOX− or H2O molecules were replaced by turquoise solid lines. “s” is the reduced gradient of the electron density

The Gibbs free energy of complexation of the Zn(II)-DOX species in solution phase (Table 4) generally exhibit positive or slightly negative values for both stoichiometries, indicating that the formation of zinc complexes in aqueous solution is less favorable than that of copper complexes. The tested functionals consistently suggest that the 1:1 stoichiometric species may preferentially adopt a pentahedral coordination geometry (l = 3), and, according to results obtained with the ωB97XD-D, M06 and PBE0 functionals, the corresponding formation process is almost equally favored if the Zn(II)-aquo complex is in its penta- or hexahydrate forms (m = 5, 6; note that the latter is reported as the most stable in solution phase). All tested functionals indicate that this process is almost insensitive to whether the DOX− molecule is attached to the metal center from the alpha or beta side, and thus, both binding patterns are prone to be observed. Based on results obtained from the M06 functional, the formation of the octahedral complex is more likely, especially when the zinc pentahydrate (m = 5) is considered as reactant; nonetheless, one must recall that this Zn(II)-aquo complex is practically absent in aqueous solution, and thus, we rather preferred to consider the pentahedral zinc-doxorubicin complex as the representative 1:1 stoichiometric species, which can be formed from the hexahydrate species (which is abundant in aqueous solution), and has a complexation energy difference less than 1.0 kcal·mol−1 respect the octahedron.

Table 4 Solution phase Gibbs free energy of complexation between \(\left[ }\left( }_ }} \right)_m } \right]^\) and n HDOX ligands, forming \(\left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^\) complexes

Our findings suggest that the formation of the 1:2 stoichiometric Zn(II)-DOX complexes are not too favored; the PBE0 and ωB97X-D functionals predict that none of the complexation processes analyzed for these species will occur (all are endergonic), while M06 and PBE indicate that the octahedral (l = 2) or the tetrahedral (l = 0) complexes may exist in solution, respectively, with the drug attached to the metal from the alpha side, in any case. Based on our benchmarking study (Table 3), where we concluded that results obtained with the M06 functionals are of high priority, we selected the \(\left[ }\left( }} \right)_2 \left( }_ }} \right)_2 } \right]\) complex as the representative Zn(II)-DOX 1:2 species. This complex can be formed from three different zinc hydrates (m = 4, m = 5, m = 6). However, it should be noted that the pentahydrate is less abundant in aqueous solution, and more importantly, that the formation from the tetrahydrate is the most exergonic process. Now that we have identified the representative species for both stoichiometries, we can analyze their relative stability using the results obtained from the M06 functional. Both, the \(\left[ }\left( }} \right)\left( }_ }} \right)_3 } \right]^+\) and \(\left[ }\left( }} \right)_2 \left( }_ }} \right)_2 } \right]\) complexes are expected to be formed from the zinc tetrahydrate, which is abundant in solution phase. The relative stability of these species is thus directly determined by the difference in their molecular Gibbs free energies. Our calculations show that such energy difference is about 3 kcal·mol−1, favoring the 1:1 stoichiometric species and is too small to draw definitive conclusions about the predominance of one species over the other, for instance, recall that a higher overestimation of the loss of entropy is expected for the formation of the 1:2 complex. Figure 7a, b depict the geometry of the \(\left[ }\left( }} \right)\left( }_ }} \right)_3 } \right]^+\) and the \(\left[ }\left( }} \right)_2 \left( }_ }} \right)_2 } \right]\) species conjointly with the BCPs surrounding their corresponding coordination centers. The former display two BCP flanking the coordination sphere, while the second displays seven BCPs. These findings suggest that intramolecular interactions are more important for the stability of the 1:2 stoichiometry species, than in their 1:1 counterpart.

Fig. 7figure 7

BCPs for the a \(\left[ }\left( }} \right)\left( }_ }} \right)_3 } \right]^+\) and b \(\left[ }\left( }} \right)_2 \left( }_ }} \right)_2 } \right]\) species, obtained at the ωB97X-D/ DGDZVP level of theory. Only BCPs between H2O molecules and DOX− are plotted, while BCPs between the metallic center and DOX− or H2O molecules were replaced by turquoise solid lines

VO 2+doxorubicin complexes

It is documented that the VO2+ ion is found rather in acidic aqueous solutions at pH values up to 6.0 [114]. Therefore, the formation of the oxovanadium-doxorubicin complexes must occur at pH values below 6.0, however, the protonated (cationic) form of doxorubicin predominates under such conditions. Consequently, both doxorubicin and VO2+ species will carry positive electric charges, and their complexation may be hindered by electrostatic repulsion forces unless doxorubicin gets first deprotonated to its neutral form (for instance, due to the interaction with the metallic center). The whole complexation process can be thus described by the following Hess law,

$$\begin \left[ }\left( }_ }} \right)_5 } \right]^ + n} \to \left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^ + n}^+ \\ n}_2 }^+ \to n} + n}^+ \\ \_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\\ \left[ }\left( }_ }} \right)_5 } \right]^ + n}_2 }^+ \to \left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^ + 2n}^+ \\ \end$$

(12)

The Gibbs free variation of such process is given by,

$$\begin \Delta_r G_^* = \Delta_r G_^* + nRT \cdot }K_a \\ = \Delta_r G_^* + n\left[ } \cdot }^ } \right] \\ \end$$

(13)

where \(\Delta_r G_^*\) is the solution phase Gibbs free energy of reaction computed through our strategy based on TCs. In Eq. (13) we used the experimental pKa value of 7.84 corresponding to the deprotonation of the amine moiety of \(}_2 }^+\), while the value of the thermodynamic temperature T was set equal to 298.15 K. It should be noted that the addition of HDOX species to the VO(IV) center is more restricted than for Cu(II) and Zn(II), due to the presence of the oxo ligand. The coordination vacancies must be accurately classified to unambiguously identify the generated coordination geometries, which were established to be equatorial or axial, as depicted in Fig. 8. A HDOX molecule can be linked to the VO(IV) center through two equatorial vacancies or one equatorial and one axial. To indicate the former case, a superscript E will be included in our notation, whereas the superscript A will be used for the latter.

Fig. 8figure 8

Molecular structure of the \(\left[ }\left( }_ }} \right)_5 } \right]^\) species displaying the axial (ax) and equatorial (eq) coordination sites

Table 5 presents the results of our benchmarking study for the VO(IV)-doxorubicin complexes. All the theoretical methods here tested predicted that the octahedral complex \(\left[ }\left( }} \right)\left( }_ }} \right)_3 } \right]^+\) is the most stable 1:1 stoichiometric species in gas phase. DOX− was observed to prefer the beta orientation, occupying two of the equatorial coordination vacancies of the vanadium center. The PBE0 functional displayed the closest reaction energies to our DLPNO-CC standard, with an average deviation (AD) of 6.88 kcal·mol−1 followed by ωB97X-D (AD = 14.81 kcal·mol−1), M06 (AD = 15.62 kcal·mol−1) and PBE (AD = 81.40 kcal·mol−1). The ωB97X-D DFA in turn displayed the best (closest to our standard) relative energy differences (quantities in parentheses in Table 1) with an average deviation of 2.49 kcal·mol−1, followed by PBE0 (4.53 kcal·mol−1), M06 (5.28 kcal·mol−1) and PBE (5.98 kcal·mol−1). Figure 9a–c show the BCP, 3D-NCII and 2D-NCII plots for the \(\left[ }\left( }} \right)\left( }_ }} \right)_3 } \right]^+\) species, respectively displaying two attractive hydrogen bonding interactions (indicated by blue surfaces) that stabilize the complex. The most pronounced attractive non-covalent interactions are observed in the vicinity of the coordination sphere, primarily resulting from the interaction between solvent molecules and the metal center. This finding highlights the importance of adequately filling coordination vacancies to determine the most stable oxovanadium complex in solution phase. By ensuring the appropriate coordination environment, the complex can engage in favorable interactions with the surrounding solvent molecules, enhancing its stability and overall behavior.

Table 5 Computed energy variation corresponding to the reaction between \(\left[ }\left( }_ }} \right)_m } \right]^\) and HDOX species, forming \(\left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^\) complexes, in gas phase, using Eq. (11)Fig. 9figure 9

a BCPs, b 3D-NCII and c 2D-NCII plots obtained for the \(\left[ }\left( }} \right)\left( }_ }} \right)_3 } \right]^+\) species at the ωB97X-D/ DGDZVP level of theory. Only BCPs between H2O molecules and DOX− are plotted, while BCPs between the metallic center and DOX− or H2O molecules were replaced by turquoise solid lines. “s” is the reduced gradient of the electron density.

The Gibbs free energies of complexation of VO(IV) with doxorubicin in acidic aqueous solution were computed using Eq. (12) and are reported in Table 6. Except for the PBE functional, all the methods predicted that the \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]\) complex, which has a square pyramid coordination geometry, is the only predominant 1:1 stoichiometric species in solution. As in the gas phase, the oxygen atoms of DOX occupy equatorial vacancies of the oxovanadium center, but in solution phase, the alpha orientation was preferred by the drug. For the 1:2 stoichiometry case, we observed that the formation of the square pyramidal \(\left[ }\left( }} \right)_2 } \right]\) species was highly exergonic only when the drug binds to the oxovanadium center through the beta side. Nonetheless, our results strongly suggest that the 1:1 complex predominates over the 1:2 complex, as indicated by a significant Gibbs free energy difference of 53.78 kcal·mol−1 (computed using the PBE0 functional). This substantial energy difference is unlikely to be counteracted by errors arising from the overestimation of the loss of entropy. Figure 10a, b show the coordination geometries of the \(\left[ }\left( }} \right)\left( }_ }} \right)_2 } \right]\) and \(\left[ }\left( }} \right)_2 } \right]\) species. Both structures show only one BCP, which may be indicating that the intramolecular interactions may not be too relevant for the stability of VO(IV)-doxorubicin complexes in solution phase.

Table 6 Solution phase Gibbs free energy of complexation between \(\left[ }\left( }_ }} \right)_m } \right]^\) and n HDOX ligands, for \(\left[ }\left( }} \right)_n \left( }_ }} \right)_l } \right]^\) complexesFig. 10figure 10

BCPs obtained at the ωB97X-D/aug-cc-pVDZ level. (a) \(\left[ }\left( }} \right)(}_2 })_2 } \right]\) and (b) \(\left[ }\left( }} \right)_2 } \right]\). Only BCPs between H2O molecules and DOX− are plotted, while BCPs between the metallic center and DOX− or H2O molecules were replaced by turquoise solid lines

Conclusions

In this study, we used a hybrid thermodynamic cycle approach to analyze the complexation equilibrium between the antineoplastic doxorubicin and three transition metal ions, Cu2+, Zn2+, and VO2+, in aqueous solution. To improve accuracy, we performed a preliminary benchmarking study to minimize errors resulting from the approximated electronic structure of the involved species. Additionally, we used explicit solvation with solvent clusters to enhance errors cancellation between the solvation patterns of the participating solutes. Our analysis revealed the representative 1:1 and 1:2 stoichiometric species, as well as the most probable complexation process for each metal–complex family. For instance, we found that copper-doxorubicin complexes have a squared planar coordination geometry at both stoichiometries and are formed from Cu(II) in its hexahydrate form; doxorubicin is linked to the metal center from oxygens located at the beta side. The 1:1 complex was found to be slightly more stable than the 1:2 complex but the validity of this result strongly depends on the overestimation of the loss of entropy going from gas to solution phase. For zinc, we found that the most favored complexes were the penta- and tetrahedral species for the 1:1 and 1:2 stoichiometries, respectively. These complexes are formed from the tetra and hexahydrates of free zinc and chelated by oxygens from the alpha side of doxorubicin. We did not find evidence of one stoichiometry predominating over another. In the case of oxovanadium-doxorubicin complexes, both 1:1 and 1:2 stoichiometric species exhibited a squared pyramidal geometry and were formed from the free octahedral VO2+ species in acidic aqueous solution, with doxorubicin exhibiting the alpha and beta orientations, respectively. Our results suggest that the 1:1 stoichiometric oxovanadium complex is predominant, while the 1:2 species may be scarce in aqueous solution.

Interestingly, our thermodynamic approach predicted entropy contributions close to zero to the Gibbs free energy of reaction in all the studied processes. This consistent pattern across the different complexes would be attributed to the thermodynamic model employed in estimating solution phase quantities. Our computational strategy is carefully designed to minimize errors and facilitate the cancellation of potential inaccuracies. However, it is important to acknowledge that the cancellation of minority Gibbs free energy components can also impact the computed reaction entropies. The accurate determination of the reaction entropy in the solution phase remains a challenging task and can be one of the major sources of error in strategies that rely on continuum solvation models and thermodynamic cycles for calculating (contributions to) solvation Gibbs free energies of reactions in polar solvents. Unfortunately, the only way to quantify this error is by comparing the predicted values with experimental measurements.

The coordination numbers inferred in our thermodynamic study were corroborated by the analysis of corresponding BCPs, which also served as indicators of relevant inter- and intramolecular interactions involved in the stability of the studied complexes. NCI indexes allowed us to build plausible discussions regarding the stability of one complex over another, in terms of the non-covalent interactions exhibited by each compound. Together, they constituted valuable tools to rationalize the thermodynamic stability of the set of metal-doxorubicin complexes investigated here.

Among the four DFAs tested in our study, PBE0 consistently produced Gibbs free energies of reactions that were the closest to our DLPNO/aug-cc-pVDZ standard for all three families of metal complexes investigated, followed by ωB97X-D. Relative energies were reasonably well predicted with the M06 functional, while ωB97X-D was the second-best performing functional. Likewise, in solution phase, these three functionals consistently provided the same trends among the sets of metal complexes studied, which also can be attributed to the cancellation of errors propitiated by our hybrid thermodynamic approach. Based on these results, we conclude that PBE0, ωB97X-D and M06 are suitable for studying the complexation process between the three metal ions considered in this study and biomolecules containing the anthracycline moiety, in both, gas and solution phase, particularly using the hybrid thermodynamic approach described herein. Conversely, our findings suggest that the PBE functional may not provide reliable results in further studies and must be avoided.

Finally, our study demonstrates that with the aid of a proper thermodynamic theoretical protocol and adequate theoretical tools, valuable information about the complexation between metal ions and bioactive molecules can be inferred, even in the absence of experimental data. Our results also suggest that reasonable predictions can be made, which can subsequently be compared with experimental determinations to gain a deeper understanding of the coordination chemistry of drugs bound to transition metal ions. Overall, our study highlights the potential of computational methods based on thermodynamic cycles in the field of bioinorganic chemistry and provides a useful framework for future studies in this area.

留言 (0)

沒有登入
gif