Assessment of the Validity of Some Minnesota Density Functionals for the Prediction of the Chemical Reactivity of the SYBR Green I and Ethidium Bromide Nucleic Acid Stains

Abstract

This research work has assessed many Minnesota density functionals to find their molecular structure and electronic properties possessed by SYBR green I (SYBRGI) and ethidium bromide (EtBr) nucleic acid stains. In the determination of the global descriptors that come up from conceptual density functional theory (CDFT), the processes include: Self-Consistent Field Energy Differences (∆SCF) and higher occupied molecular orbital (HOMO) and lower unoccupied molecular orbital (LUMO) frontier orbitals energies. Regarding the deduced outcomes for the conceptual DFT indices, many of the descriptors have been adjusted to achieve the “Koopmans in DFT (KID)” process. It has also been shown that the only density functionals that confirm this approximation are the range-separated hybrids (RSH).

1. Introduction

The chemical reactivity theory [also known as the conceptual density functional theory (CDFT)] is a vital technique that is used to predict, evaluate, and interpret the results from chemical processes [1, 2, 3, 4].

Research done by Parr and his associates [1] reveals that several theories and models have been discovered after the evaluation of the molecular system with the use of DFT. Almost all the discovered theories are helpful in research because they enable scholars to achieve quantitative forecasts of a chemical reactivity system. In addition to this, the theories can further be quantified and are generally termed as conceptual DFT descriptors.

To obtain the quantitative figures of conceptual descriptors, it may be necessary to analyze the Kohn-Sham theory by calculating the energy system, molecular density, and the orbital energies associated with frontier orbitals [i.e., higher occupied molecular orbital (HOMO) and lower unoccupied molecular orbital (LUMO)] [5, 6, 7, 8, 9, 10].

For research on a molecular system, the first activity before proceeding is selection of the model chemistry. Model chemistry is the collection of basic set, density function, together with an implied solvent model that is known to be consistent for the problem under research. Several studies provide insights on the way to choose the model chemistry. A researcher may also decide to preview past studies when choosing the model.

Even though the fundamentals of DFT reveal that universal functional density is present and that computations using this function can be used to obtain all the features of the system, it is always necessary in practical cases that one refers to the estimated density functionals that have been established for the past 3 decades. For the approximate functionals, almost all of them are perfectly fit to be used in estimating some features, while some can be used for estimating other features. In separate scenarios, you can encounter density functionals that are perfectly fit for estimating the features of a given molecular system and a functional group. It is also important to assess separate density functionals for a separate functional group which can be added to the molecular system under research.

When researching on chemical reactivity (which is a process that entails the transfer of electrons), a person performs computations for both ground and open systems, i.e., cation and anion. It is not easy to obtain consistent outcomes using these computations (when diffuse functions should be a part of the basis set) [5, 6, 7, 8, 9, 10]. This necessitates adoption of a more consistent technique that provides all the data that a person will require directly from the outcomes of the computations at ground state in the molecular system under research. In addition to this, a person may also want to find the deionization ability together with the electron affinity of any system being researched without having to calculate the radical cation and anion. This can be determined by the Koopmans’ theorem [7, 8, 9, 10] that relies on Hartree-Fock Theory, which states that the energy of the HOMO (i.e., I = −εH) can be used to estimate the ionization potential. Alternatively, the electron affinity can be estimated using the minus the energy of the LUMO (i.e., A = −εL).

The legitimacy of the Koopmans’ theorem is yet again a contentious issue because of the existing difference between the fundamental band gap and the HOMO and the LUMO gaps. This can be termed as derivative discontinuity. It has again been discovered that an exact physical description may be assigned to Kohn-Sham HOMO using “the Kohn-Sham analogue of Koopmans’ theorem in Hartree-Fock theory” (this theory explains that in the exact theory, the KS HOMO is opposite and same as the ionization potential) [11, 12, 13, 14]. The effects brought about by the difference between the fundamental band gap and the HOMO and the LUMO gaps have ensured that no Koopmans’ theorem creates a direct relationship between the LUMO energy and the electron affinity. To eliminate these effects, a suggestion has been made by scholars to conceive that the ionization potential of the N + 1 electron system (anion) is almost equal to the electron affinity of N electron system [15]. Regarding the range-separated hybrids (RSH) density functionals [16, 17, 18], e.g., that the repulsive coulomb potential has to be separated in the long-range (LR) and short-range (SR) terms, e.g., via r⁻¹ = r⁻¹ erf (γr) + r⁻¹ erfc (γr), with γ representing the range-separation parameter, it was highlighted by Kronik et al. [15] that the legitimacy of the Koopmans’ theorem may be approximately approved with a well-informed choice of this final parameter.

Use of the γ-tuning technique can be useful in upgrading the features projected by the density functionals. This is due to better utilization of the of Koopmans’ theorem, leading to a superior understanding of frontier orbitals energies with I and A. An example can be cited from Lima et al. [19], who just demonstrated a better explanation of the optical features of carotenoids through tuning of some LR density functionals that are linked.

This therefore illustrates that conformity of any given density functional can be determined by assessing how it has adopted the “Koopmans’ in DFT” (KID) process, which guides its behavior to be almost equal to the ideal density functional. This is essential for any precise computation of the conceptual DFT descriptors that help in forecasting and analysis of chemical reactivity in molecular systems. Still, the γ-tuning technique for range-separated hybrids density functional is system dependent. This implies that separate density functionals are to be used in the computations of the descriptors for separate molecular features. We are then going to concentrate on part of the density functionals that have displayed the required precision in physics and chemistry [20].

The main aim of this study is to do a comparative research relating to the performance of the just identified Minnesota family of density functionals for the account of the chemical reactivity of two nucleic acids intercalating stains, SYBR green I (SYBRGI) [21] and ethidium bromide (EtBr) [22]. The molecular structures of the two are shown in Figure 1.

Figure 1.

Molecular structures of (a) ethidium bromide (EtBr) and (b) SYBR green I (SYBRGI).

2. Theoretical background

Within the context of DFT [2, 23], the chemical potential μ, which estimates the escaping tendency of the electron from the equilibrium point is stated as follows:

$$\mu ={\left(\frac{\partial E}{\partial N}\right)}_{v\left(r\right)}$$E1

where χ represents the electronegativity.

Chemical hardness is represented by η, which is the opposition to charge transfer:

$$\eta ={\left(\frac{{\partial }^{2}E}{{\partial N}^2}\right)}_{v\left(r\right)}$$E2

Employing a finite difference approximation and the Koopmans’s theorem [7, 8, 9, 10], we can write the above equation as:

$$\mu =-\frac{1}{2}\left(I+A\right)\approx \frac{1}{2}\left({\epsilon }_L+{\epsilon }_H\right)=-{\chi }_K$$E3

$$\eta =\left(I-A\right)\approx \left({\epsilon }_L-{\epsilon }_H\right)={\eta }_K$$E4

where ε_H is the HOMO energy and ε_L is the LUMO energy.

An expression for the electrophilicity index ω is as below:

$$\omega =\frac{{\mu }^2}{2\eta }=\frac{{\left(I+A\right)}^2}{4\left(I-A\right)}\approx \frac{{\left({\epsilon }_L+{\epsilon }_H\right)}^2}{4\left({\epsilon }_L-{\epsilon }_H\right)}={\omega }_K$$E5

Expressions for electrodonating ω⁻ and electroaccepting ω⁺ powers are as below [24]:

$${\omega }^{-}=\frac{{\left(3I+A\right)}^2}{4\left(I-A\right)}\approx \frac{{\left(3{\epsilon }_H+{\epsilon }_L\right)}^2}{16{\eta }_K}={\omega }_K^{-}$$E6

and

$${\omega }^{+}=\frac{{\left(I+3A\right)}^2}{4\left(I-A\right)}\approx \frac{{\left({\epsilon }_H+3{\epsilon }_L\right)}^2}{16{\eta }_K}={\omega }_K^{+}$$E7

To obtain a comparison for ω⁺ and − ω⁻, the explanation below for net electrophilicity has been suggested [25]:

$$∆{\omega }^{\pm }={\omega }^{+}-\left(-{\omega }^{-}\right)={\omega }^{+}+{\omega }^{-}\approx {\omega }_K^{+}-\left(-{\omega }_K^{-}\right)={\omega }_K^{+}+{\omega }_K^{-}=∆{\omega }_K^{\pm }$$E8

3. Settings and computational details

In this research project, each of the computations was done using Gaussian 09 programs [26] and the density functional methods as compelled in the computational package. The gradient method was used to obtain the equilibrium geometries of molecules in this research. Additionally, vibration frequencies and the force parameters were estimated through computation of analytical frequencies on still areas after optimization to check whether they were the actual minima. Def2SVP was used in this research project as the basic set for optimization of geometry and frequencies. Computation of the electronic features was achieved using Def2TZVP basic set [27, 28].

We chose a group of Minnesota density functionals, which give consistent outcomes when computing the molecular structures and systems during the research activities. The group below were selected: M11, which falls under RSH meta-generalized gradient approximation (GGA) [29]; M11L, which falls under double-range local meta-GGA [30]; MN12L, which falls under nonseparable meta-nonseparable gradient approximation (NGA) [31]; MN12SX, which falls under nonseparable hybrid nonseparable meta-NGA [32]; N12, which falls under nonseparable gradient estimation [33]; N12SX, which falls under RSH nonseparable gradient estimation [32]; SOGGA11, which falls under generalized gradient approximation (GGA) density functional [34]; and SOGGA11X, which falls under generalized gradient approximation density functional [35]. GGA can be explained as a functional whereby the functional relies on the both the up down spin densities and the minimized gradient. Nonseparable gradient approximation (NGA) can be explained as a functional that relies on up down spin densities and minimized gradient and falls under nonseparable. In the current research, each of the computations was done where the solvent used was water and by doing the calculations in conformity to the Solvation Model based on Density (SMD) solvation standard [36].

4. Results and discussion

Firstly, the molecular structures of SYBRGI and EtBr were first optimized by MOL structures, then by finding the most stable conformers through Avogadro program [37, 38]. This was done through random sampling with molecular mechanic techniques. After the optimization, the resulting conformers were then re-optimized with MN12SX, M11L, M11, N12, SOGGA11X, and SOGGA11 density functionals. In addition, Def2SVP basic set and SMD solvation model were used, water being used as a solvent.

The HOMO and LUMO energies (in eV), the ionization potential I and electron affinity A (in eV), electronegativity χ, chemical hardness η, global electrophilicity ω, electrodonating ω⁻, and electroaccepting ω⁺ powers, and the net electrophilicity ∆ω^± of the EtBr and SYBRGI molecules calculated with the same density functionals and solvation model are presented in Tables 1 and 2, respectively. The upper part of the tables shows the results derived assuming the validity of the Kid procedure (hence the subscript K), and the lower part of the tables shows the results derived from the calculated vertical I and A through a ∆SCF technique.

Property	M11	M11L	MN12L	MN12SX	N12	N12SX	SOGGA11	SOGGA11X
HOMO	−7.535	−4.956	−4.684	−5.187	−4.164	−4.951	−4.151	−6.108
LUMO	−0.535	−3.028	−2.613	−2.576	−2.674	−2.511	−3.121	−1.808
χK	4.035	3.992	3.649	3.882	3.419	3.731	3.636	3.958
ηK	7.000	1.928	2.071	2.611	1.490	2.440	1.030	4.300
ωK	1.163	4.133	3.213	2.886	3.924	2.852	6.420	1.822
ω^−K	4.781	10.383	8.380	7.876	9.651	7.723	14.722	5.891
ω^+K	0.746	6.391	4.732	3.994	6.232	3.992	11.086	1.933
∆ω^±K	5.527	16.774	13.112	11.869	15.883	11.714	25.808	7.824
I	5.585	5.183	4.883	5.216	4.520	4.970	4.819	5.385
A	2.711	2.782	2.408	2.624	2.344	2.579	2.788	2.663
χ	4.148	3.983	3.646	3.920	3.432	3.775	3.804	4.024
η	2.874	2.401	2.475	2.592	2.177	2.392	2.031	2.721
ω	2.994	3.303	2.685	2.965	2.706	2.978	3.562	2.975
ω^−	8.241	8.748	7.348	8.051	7.263	7.993	9.152	8.132
ω^+	4.093	4.786	3.702	4.131	3.831	4.219	5.348	4.108
∆ω^±	12.335	13.514	11.050	12.182	11.094	12.212	14.500	12.241

Table 1.

HOMO and LUMO energies (in eV), the ionization potential I and electron affinity A (in eV), electronegativity χ, chemical hardness η, global electrophilicity ω, electrodonating ω⁻ and electroaccepting ω⁺ powers, and the net electrophilicity ∆ω^± of the EtBr molecule.

Property	M11	M11L	MN12L	MN12SX	N12	N12SX	SOGGA11	SOGGA11X
HOMO	−7.593	−5.132	−4.924	−5.325	−4.587	−5.160	−4.779	−6.099
LUMO	−0.558	−2.933	−2.569	−2.545	−2.618	−2.460	−2.959	−1.768
χK	4.075	4.033	3.747	3.935	3.603	3.810	3.869	3.933
ηK	7.034	2.199	2.355	2.780	1.970	2.700	1.820	4.331
ωK	1.181	3.698	2.980	2.785	3.295	2.688	4.112	1.786
ω^−K	4.839	9.550	7.981	7.711	8.514	7.449	10.273	5.809
ω^+K	0.763	5.517	4.234	3.776	4.911	3.639	6.404	1.876
∆ω^±K	5.602	15.067	12.214	11.486	13.425	11.088	16.677	7.685
I	5.407	5.302	5.062	5.263	4.767	5.074	4.980	5.252
A	2.643	2.747	2.400	2.589	2.358	2.521	2.690	2.564
χ	4.025	4.024	3.731	3.926	3.563	3.797	3.835	3.908
η	2.764	2.555	2.661	2.674	2.409	2.553	2.290	2.688
ω	2.931	3.169	2.615	2.882	2.635	2.824	3.210	2.841
ω^−	8.048	8.509	7.262	7.894	7.202	7.706	8.481	7.803
ω^+	4.023	4.485	3.531	3.968	3.639	3.909	4.646	3.895
∆ω^±	12.071	12.994	10.972	11.862	10.841	11.615	13.128	11.699

Table 2.

HOMO and LUMO energies (in eV), the ionization potential I and electron affinity A (in eV), electronegativity χ, chemical hardness η, global electrophilicity ω, electrodonating ω⁻ and electroaccepting ω⁺ powers, and the net electrophilicity ∆ω^± of the SYBRGI molecule.

For examining the outcomes to determine if the KID process is fulfilled, and the drive from past works [15, 19], we have come up with descriptors having the ability to compare the outcomes from HOMO and LUMO computations with those attained using vertical I and A and a ΔSCF technique. It should again be known that we have no plans to form a gap fitting by reducing the descriptor. We plan to determine if the density functionals employed in this research contain the fixed range parameter γ that helps in effective execution of the KID process. It is somehow astonishing that our research at present lacks the parameter γ. We also included a minus of the energy of the LUMO of the neutral system instead of using A as minus of HOMO of the electron system [15, 19].

The initial three descriptors are associated with the basic accomplishment of “Koopmans in DFT” estimation by associating ε_H with −I, ε_L with −A, and their responses in explaining the HOMO-LUMO gap:

$$J_I=\left|{\epsilon }_H+E_{gs}\left(N-1\right)-E_{gs}\left(N\right)\right|$$E9

$$J_A=\left|{\epsilon }_L+E_{gs}\left(N1\right)-E_{gs}\left(N+1\right)\right|$$E10

$$J_{HL}=\sqrt{J_I^2+J_A^2}$$E11

Four separate descriptors will then be used to examine how the density functionals under research will help in forecasting the electronegativity χ, the chemical hardness η, the global electrophilicity ω, and the collection of conceptual DFT descriptors through deliberation of the energies of the HOMO and LUMO or the vertical I and A:

$$J_{\chi }=\left|\chi -{\chi }_K\right|$$E12

$$J_{\eta }=\left|\eta -{\eta }_K\right|$$E13

$$J_{\omega }=\left|\omega -{\omega }_K\right|$$E14

$$J_{D1}=\sqrt{J_{\chi }^2+J_{\eta }^2+J_{\omega }^2}$$E15

D1 represents the initial collection of conceptual DFT descriptors.

Finally, we came up with four extra descriptors to determine the success of the density functionals under research in forecasting of electrodonating power ω−, the electroaccepting power ω⁺, the net electrophilicity ∆ω^±, together with the four descriptors combined and considering the HOMO and LUMO energies or the vertical I and A:

$$J_{{\omega }^{-}}=\left|{\omega }^{-}-{\omega }_K^{-}\right|$$E16

$$J_{{\omega }^{+}}=\left|{\omega }^{+}-{\omega }_K^{+}\right|$$E17

$$J_{∆\omega \pm }=\left|∆\omega \pm -∆\omega {\pm }_K\right|$$E18

$$J_{D2}=\sqrt{J_{\omega -}^2+J_{\omega +}^2+J_{∆\omega \pm }^2}$$E19

D2 represents the second collection of conceptual DFT descriptors.

The results of the calculations of J_I, J_A, J_HL, J_χ, J_η, J_ω, J_D1, J_ω −, J_ω +, J_∆ω±, and J_D2 for the EtBr and SYBRGI are displayed in Tables 3 and 4, respectively.

Descriptor	M11	M11L	MN12L	MN12SX	N12	N12SX	SOGGA11	SOGGA11X
J1	1.950	0.228	0.199	0.028	0.356	0.020	0.668	0.724
JA	2.176	0.246	0.205	0.048	0.331	0.068	0.333	0.855
JHL	2.922	0.335	0.285	0.055	0.486	0.071	0.747	1.120
Jχ	0.113	0.009	0.003	0.038	0.013	0.044	0.168	0.066
Jη	4.127	0.473	0.403	0.019	0.687	0.048	1.001	1.579
Jω	1.831	0.830	0.528	0.079	1.219	0.126	2.858	1.154
JD1	4.516	0.955	0.665	0.090	1.399	0.142	3.033	1.957
Jω−	3.460	1.635	1.033	0.175	2.388	0.271	5.570	2.242
Jω+	3.347	1.626	1.030	0.137	2.401	0.227	5.738	2.176
J∆ω+	6.808	3.260	2.062	0.313	4.788	0.498	11.307	4.417
JD2	8.338	3.993	2.526	0.384	5.864	0.611	13.849	5.410

Table 3.

Descriptors J_I, J_A, J_HL, J_χ, J_η, J_ω, J_D1, J_ω −, J_ω+, J_∆ω±, and J_D2 for the ethidium bromide (EtBr) molecule calculated from the results of Table 1.

Descriptor	M11	M11L	MN12L	MN12SX	N12	N12SX	SOGGA11	SOGGA11X
J1	2.185	0.170	0.137	0.062	0.180	0.086	0.201	0.847
JA	2.085	0.187	0.169	0.044	0.259	0.061	0.269	0.796
JHL	3.020	0.252	0.218	0.076	0.316	0.105	0.336	1.162
Jχ	0.050	0.009	0.016	0.009	0.040	0.012	0.034	0.025
Jη	4.270	0.356	0.306	0.106	0.439	0.147	0.470	1.643
Jω	1.751	0.529	0.365	0.097	0.660	0.136	0.902	1.055
JD1	4.615	0.638	0.477	0.144	0.794	0.201	1.018	1.952
Jω−	3.209	1.041	0.719	0.183	1.312	0.257	1.792	1.994
Jω+	3.259	1.032	0.703	0.192	1.272	0.270	1.758	2.019
J∆ω+	6.469	2.073	1.422	0.376	2.584	0.527	3.550	4.013
JD2	7.923	2.539	1.742	0.460	3.165	0.646	4.348	4.915

Table 4.

Descriptors J_I, J_A, J_HL, J_χ, J_η, J_ω, J_D1, J_ω−, J_ω +, J_∆ω± and J_D2 for the SYBR green I (SYBRGI) molecule calculated from the results of Table 2.

As shown in Tables 1 and 2, and the outcomes from Tables 3 and 4, the precision provided by the KID process is outstanding for the MN12SX, which falls under RSH meta-NGA, and N12SX, which falls under RSH NGA density functionals. In reality, values for J_I, J_A and J_HL is not zero. However, the values found can satisfactorily be likened to the past studies of Lima et al. [19], whereby the minima were found by selecting a parameter that imposes such a trend.

The outcomes are necessary because they reveal that we should not depend on J_I, J_A, and J_HL alone, i.e., if we depend on outcomes from J_χ, alone, almost all the values will near zero. For the remaining descriptors, only MN12SX and N12SX reveal such trends. This shows that outcomes for J_χ can be due to elimination of errors.

Authentication of the KID process is not done correctly by the GGA (SOGGA11) and hybrid-GGA (SOGGA11X). Local density functionals like M11L, MN12L, and N12 are also inappropriate.

It is vital to know that even though the RSH hybrid NGA and RSH meta-NGA density functionalities are necessary when computing the conceptual DFT descriptors, it is a different case for RSH GGA (M11) density functional. According to Tables 1 and 2, this functional doesnot provide enough explanation concerning LUMO energy, and this can be due to an inaccurate figure of γ in the functional. A fine tuning of γ can handle the issue.

5. Conclusions

Weighing on the outcomes from this research work, DFT-based reactivity descriptors like electronegativity, chemical hardness, global electrophilicity, electrodonating, and electroaccepting powers, and net electrophilicity can be used to forecast EtBr’s chemical reactivity.

It has also been illustrated that the KID process can effectively be implemented by the RSH meta-NGA (MN12SX) and the RSH NGA (N12SX) density functionalities. They can then be used in place of the tuned density functionals using a gap-fitting process, and we believe that such a trend can be helpful when analyzing the chemical reactivity of bigger molecular systems.

Acknowledgments

This work has been partially supported by CIMAV, SC, and Consejo Nacional de Ciencia y Tecnología (CONACYT, Mexico) through Grant 219566/2014 for Basic Science Research and Grant 265217/2016 for a Foreign Sabbatical Leave. Daniel Glossman-Mitnik conducted this work while a Sabbatical Fellow at the University of the Balearic Islands from which support is gratefully acknowledged. This work was also funded by the Ministerio de Economía y Competitividad (MINECO) and the European Fund for Regional Development (FEDER) (CTQ2014-55835-R).