Abstract

Indocyanine green is a great near-infrared fluorescence with good luminescent properties and important medical applications. In this paper, the theoretical spectrum and orbital model of its molecular level are established. The two most probable conformations were studied, and their energies, vibrational spectra, UV-Vis absorption spectra, frontier molecular orbitals (HOMO and LUMO), and energy gaps were obtained by density functional theory (DFT) calculations, respectively. This provides a theoretical and design basis for the development of novel dyes similar to indocyanine green dyes and a reference case for improved application methods and synthetic predesign of novel fluorescent dyes.

1. Introduction

With the increasing incidence of end-stage liver disease and liver cancer, the liver function evaluation is of great significance for clinical diagnosis and treatment [1]. Indocyanine green (ICG) is a commonly used diagnostic drug, mainly used as a safe and reliable approach to examine the liver function and effective blood flow. Obana et al. reported that 13 cases of adverse reactions were found in 3774 ICG angiography [2]. Indocyanine green can be ingested and excreted by the liver cells in a prototypical manner and can be used to directly evaluate the liver metabolism and thus liver function [3]. For patients preparing for hepatocellular carcinoma resection, the evaluation and prediction of residual liver tissue function with indocyanine green is of great clinical significance [4]. In the human body, light less than 700 nm is absorbed by tissues, light greater than 900 nm is blocked by water and fat, and infrared light with wavelengths between the two can effectively penetrate tissues, thereby obtaining corresponding biological information [5]. In fluorescence imaging, red light mainly has the advantages of low tissue absorption, low scattering, low autofluorescence background, and high signal ratio and can objectively reflect the metabolic state of living organisms [6, 7]. The optical properties of indocyanine green determine its good biological application prospects. In clinical application, it can truly reflect the metabolic status of the liver [8], predict the compensatory ability of postoperative residual liver tissue [9], and determine the number and boundary of liver tumors [10]. A more in-depth study of the luminescent properties of indocyanine green at the atomic level can make use of its advantages more efficiently, assist in the development of more efficient signal-receiving devices, reduce the dosage, and even explore new uses. At present, the luminescence mechanism at the molecular level is still unclear. It is hoped that by establishing a corresponding quantum mechanical calculation model, the luminescence mechanism can be accurately analyzed providing a theoretical basis for subsequent research and application.

Density functional theory [11] (DFT) is a quantum chemical calculation method developed in recent decades. It has the advantages of small computation, high relative precision, and being able to deal with large molecular systems. Among them, in most cases, the field changes of excited states are time-related, so people add time-related parameters into the density functional theory equation when dealing with the time-dependent system and derive the time-dependent density functional theory (TD-DFT). With the development of electronic structure theory and computer computing power, theoretical calculation can predict excitation energy and the molecular energy gap more and more accurately. The calculation of molecular properties and spectra of small organic systems has been in good agreement with the experimental results, which provides sufficient conditions for us to study the relevant information and properties of interested fluorescent molecules [1216].

Fluorescent dye molecules absorb light from the ground state to the excited state; on the one hand, the electronic structure and properties of the excited state are greatly different from the ground state molecules, and the electronic structure of the excited state, the mutual transformation of different electronic structures, and attenuation process of the excited state directly determine the photophysical process and the final application performance. On the other hand, relative to the ground state, the excited state is a transient state, from which the experiment on the collection of information is limited by the instruments and equipment directly, but in the design and synthesis of fluorescent dyes, this information is very important, and quantum chemistry calculations can, in theory, provide information on molecular excited states and provide the relevant data for the researchers. Therefore, it is very important to study the molecular properties and luminescence mechanism of fluorescent dyes.

In this paper, density functional theory (DFT) and time-dependent density functional theory (TD-DFT) were used to study the structure information, spectral properties, and luminescence mechanism of indocyanine green (ICG) dye. The precise structure and spectral information of ICG were obtained, the reliable luminescence model was established, and the clear luminescence mechanism was clarified. It provides valuable experience and examples for its future clinical application and theoretical design of new luciferin.

2. Results and Discussion

2.1. Determination of the Optimal Ground State Configuration of IGC

Two possible gas phase configurations of IGC molecules were optimized by DFT calculation, as shown in Figure 1. At the same time, structural optimization and frequency analysis were carried out to determine the minimum value near the potential energy surface. For conformation A, the sodium ion binds to the sulfonate ion at site 1, while the sulfonate at site 2 tends to be close to the positive ammonium ion in the spatial structure. In the case of conformation B, the sodium ion binds to two sulfonate ions to form a “bridge conformation.” The energy data of the two conformations are listed in Table 1. By comparing the zero-point energy, enthalpy, and free energy, it can be seen that the energy of conformation B is about 50 kcal·mol−1 lower than conformation A. This suggests that the possibility of IGC molecules in conformation B is much greater than in conformation A.

Then, after analyzing the energy, we compared the bond lengths of the two conformations as shown in Table 2. For the two sulfonate groups in the A conformation, the length of the S-O bond is about 1.48–1.51 Å, and the length of the S-O bond in the B conformation is 1.49–1.51 Å, and there is no significant difference between them. For the C-S bond, the length is 1.81 Å and 1.82 Å, respectively, but for the spatial structure, the sodium ion in the A conformation only coordinates with one of the sulfonate groups, while the sodium ion in the B conformation coordinates with both sulfonate groups. Then, the bond lengths of the conjugated chains between the two aromatic rings were compared, which was found that for C-C bonds, both single bonds (B1 and B3) and double bonds (B2, B4, and B5), the bond lengths showed an average distribution of about 1.40. The same is true for the C-N bond (B6 and B7), both of which are about 1.5 Å. Although the structures of the two conformations are quite different (one is a straight chain molecule and the other is a ring molecule), there is no significant difference in main bond lengths between the π system and the sulfonic group. Therefore, we can speculate that the energy difference may come from the deformation energy.

2.2. Noncovalent Interaction (NCI) Analysis of IGC Molecule

Many vibrational modes (infrared spectra) are greatly affected by the presence of hydrogen bonds [1719]. Therefore, a systematic noncovalent interaction analysis of the two conformations should be performed before the theoretical vibration spectra are analyzed. Figure 2 shows that in the conformation of A, the main noncovalent interactions are distributed between the conjugate chains and within the aromatic ring, and their types are van der Waals interactions and repulsive interactions, respectively. The van der Waals action is more distributed in B than in A, and there is repulsion in the aromatic ring. Although the distribution of noncovalent interactions between the two conformations is slightly different, the types of noncovalent interactions are van der Waals interactions and repulsive interactions, and there is no obvious hydrogen bond interaction. Therefore, in the following study, we did not consider the effect of hydrogen bonds on the vibration spectrum.

2.3. Establishment of Theoretical Model of Vibrational Absorption Spectrum

The vibrational absorption spectrum, also known as the infrared spectrum, is an important means to identify organic molecules, and its peak value depends on the structure and conformation of molecules. It is important to establish the vibrational absorption spectrum of the ICG molecule to determine its fine structure. After fully optimizing two conformations of ICG by DFT method, the vibration form of the spectral peak of the infrared spectrum was pointed out as marked in Figure 3.

By comparing the theoretical infrared spectra of the two conformations of IGC, it is not difficult to find that the two conformations obtained after optimization have weak vibrations near 682 cm−1 to 744 cm−1. However, there is a significant difference between the peaks near 948.8 cm−1. Conformation A has a “double peak,” while conformation B has a “single peak” with a narrow top and wide bottom. Both of them have a strong absorption peak at 1118.9 cm−1, and there is a combined low peak on the left side of A instead of B. In addition, there is no significant difference between the vibration spectra of the two conformations, indicating that there are two distant weak peaks at 1243.1 cm−1, 1335.1 cm−1 for A, and 1247.1 cm−1, 1331.1 cm−1 for B, and a shoulder peak at 1391.1 cm−1 and 1450.1 cm−1, then a strong narrow peak and a weak wide peak at 1450.1 cm−1 and 1562.0 cm−1. Therefore, it can be seen from the above that for ICG molecule, different conformations (although the structures of the two conformations are quite different) do not have a great influence on the main vibration peak of the corresponding vibration spectrum, no matter the position of the peak or the intensity of the peak.

2.4. Establishment of Theoretical Model of UV-Vis Absorption Spectrum

When the light wave acts on the material, the electrons inside the material will absorb the photon energy and produce the electron energy level transition from the ground state to the excited state. In the process of the electron energy level transition, there are vibrational and rotational energy level transitions. Therefore, the UV-Vis absorption spectrum is “band spectrum.” The electron energy level transition requires a large amount of energy, generally falling in the ultraviolet-visible region (200∼800 nm). UV-Vis spectrum can not only identify the species of compounds, but also reflect their photophysical properties. Therefore, it is very necessary and important to establish UV-Vis spectral model for organic molecules, especially for organic dyes.

Based on the IGC molecular structure obtained above, we calculated the excited state of the IGC molecule at a single point using TD-DFT method and obtained the prediction results of the UV-Vis spectra of the two conformations of the IGC molecule, as shown in Figure 4. Although the structures of the two conformations are very different, they all have similar peak types for the ultraviolet spectrum, with a strong narrow peak at about 200 nm, a low peak at about 290 nm, and a wide strong peak at about 570 nm. On the whole, conformation B has a weak redshift (about 5 nm) relative to conformation A, and the strength of conformation B is significantly greater than that of conformation A at the main peak (about 570 nm). In addition, some subtle differences also appear in the UV-Vis spectra of the two conformations. For example, the peak of conformation A at about 290 nm is obviously narrower and stronger, while the peak of conformation B is shorter and wider. Moreover, conformation A also has a low wide peak at about 360 nm, but this does not appear in the spectrum of conformation B.

2.5. Distribution of Frontier Molecular Orbital (HOMO and LUMO) and Energy Gap

The frontier molecular orbital theory was proposed in the 1950s. It was based on the idea that HOMO has the most energetic and least bound electrons, and so the liveliest and most volatile electrons. LUMO has the lowest energy among all unoccupied orbitals and is the easiest to accept electrons. Therefore, these two orbitals determine the gain and loss of molecules’ electrons and their transfer ability, determine the spatial orientation of intermolecular reactions and other important chemical properties, and are also important data for studying the photophysical properties of dyes. After discussing the theoretical UV-Vis spectra of the two conformations of IGC molecules, we also plotted the distribution of the frontier molecular orbitals and the HOMO/LUMO energy gaps of the two conformations of IGC molecules, in order to explore the information of ICG molecules’ orbitals and excited states. The first is the distribution of HOMO and LUMO. For conformation B (with an annular conformation), the distribution of HOMO and LUMO on the whole molecule is relatively uniform. A few of the highest occupied molecular orbitals (HOMO) are distributed in the aromatic ring region, mainly focusing on the double bond between the two aromatic rings. For the lowest unoccupied molecular orbital (LUMO) of the B conformation, the distribution is also concentrated on the double bond between the two aromatic rings, but unlike the HOMO, the LUMO is more concentrated on the aromatic ring. (See Figure 5)

For the chain conformation A, its frontier molecular orbital distribution is very different from that of the ring conformation B, and its LUMO distribution is similar to that of conformation B, almost all of which are distributed on the double bonds in the two aromatic rings. But the HOMO distribution of conformation A is very different from that of conformation B, all of which are distributed on the sulfonic acid group that is not bound to the sodium ion. Although the HOMO and LUMO distributions of conformation A and conformation B are quite different, the energy gap of them is 4.14 eV and 4.19 eV, respectively, there is not a lot of difference. In addition, we also compared the energy gap between HOMO and LUMO of ICG molecules with N-ethyl-N-nitrosourea (ENU, a potential anticancer agent) [20], cepharanthine (a natural alkaloid) [21], edotecarin (indolocarbazole anticancer agent) [22], gemcitabine (GEM, against a variety of solid tumors) [23], and bis-chloroethylnitrosourea (BCNU, anticancer drug) [24], the results are shown in Table 3.

3. Models and Method

All the calculations in this paper were carried out in Gaussian 09 D01 software [25]. B3LYP [26, 27] functional was used as the method, and 6–311 g (d, p) basis set was used as the basis group to fully optimize the ICG molecule in the ground state. Frequency analysis was used to determine the minimum point of the optimized conformation and obtain its energy at 298 K. The 3D structure of molecules is drawn by CYL-View [28] and VMD [29]. The obtained molecular structure is directly used for single point time-dependent density functional (TD-DFT) calculation at the level of CAM-B3LYP [30] /6–311++G (d, p). It should be noted that B3LYP is widely used in the structural optimization of small molecules [13, 15, 16, 3133], and for the basis set, we chose 6-311G (d, p), which is higher than the common standard. CAM-B3LYP is widely used to study the excited state of small organic molecules [3439], based on the structure optimization base set, we add the dispersion function to ensure more accurate excited state information. The solvent effect is represented by the implicit solvent model using the SMD-H2O model (developed by Truhlar’s group) [40]. The IGC model is derived from the PubChem database. In order to better fit the actual situation, we adopted its sodium salt molecule as the research model. NCI analyses [41, 42] were performed using Multiwfn (version 3.8) software [43].

4. Conclusions

In summary, the geometric structure, relative stability, vibration spectrum, UV-Vis spectrum, and frontier orbital information of the two conformations of ICG dye molecules are systematically and deeply studied by using density functional theory (DFT). The results show that the conformation B with annular conformation is significantly more stable than the conformation A with chain conformation from the point of view of energy, and the difference of zero energy, enthalpy, and free energy between them is 63.94 kcal·mol−1, 64.66 kcal·mol−1, and 59.47 kcal·mol−1, respectively, which can fully explain that the conformation B is more stable than the conformation A. Then, we analyzed the theoretical vibration spectra of the two conformations and found that although the spatial structure of the conformations was very different, the peak position of the spectra did not change significantly, only slightly different in some places. With this conclusion, we used TD-DFT theory to calculate the single point excited state of the above two conformations and plotted their theoretical UV-Vis spectra. However, the very different structures of the two conformations did not bring about significant differences in UV-Vis spectra, except that the peak of A conformation at 290 nm is sharper than that of B conformation. Conformation A has an extra wide peak at 360 nm and the spectrum of conformation B has a weak redshift (less than 5 nm) compared with A. In addition, the distribution of the frontier molecular orbital of two different conformations of IGC molecules was plotted to further investigate the orbital information and the source of the molecular photophysical properties. We learned that although the obvious differences in the spatial structure of the two conformations do not make a significant difference in the spectrum or distribution of LUMO, there is a significant difference in the distribution of HOMO, the distribution of B conformation on both sides of the aromatic ring and the double bond between the aromatic rings, and the distribution of A conformation on the sulfonic acid group without the sodium ion. However, unlike the frontier orbital distribution, the energy gap of the two conformations is not much different, only 0.05 eV.

Therefore, we believe that although the different conformations of the ICG molecule have great differences in spatial structure and energy, they do not have great differences in spectrum, LUMO distribution or energy gap. When we study the theoretical spectrum and orbital information of dye molecules, we may be able to reduce the attention to different conformations. Moreover, the theoretical study of common medical fluorescent dyes can not only satisfy our curiosity about luminescence mechanism and influencing factors but also provide an effective idea and solution for the adjustment of absorption bands of fluorescent dyes and the development of new fluorescent dyes. Based on this work, more theoretical and experimental studies on small molecular fluorescent dyes for medical diagnosis are being carried out in our laboratory.

Data Availability

The data used in this study are available as supporting information.

Conflicts of Interest

The authors declare that there are no conflicts of interest.

Acknowledgments

The authors gratefully acknowledge this research funded by National Natural Science Foundation of China, grant number 81960323.

Supplementary Materials

The following supporting information can be downloaded at https://www./xxx/s1, ball and stick model diagram, main excited state information and Cartesian coordinates for the two conformations of ICG molecule. (Supplementary Materials)