Calculating the Aromaticity of Heterocycles

2 downloads 17 Views 3MB Size Report
nucleus-independent chemical shift, abbreviated as NICS (1996JA6317), has found great .... in the following sequence: thymine < cytosine < guanine < adenine .... Aromaticity of guanine and cytosine tautomers strongly depends on the.

Calculating the Aromaticity of Heterocycles Halina Szatylowicz1, *, Olga A. Stasyuk1, Tadeusz M. Krygowski2 1

Faculty of Chemistry, Warsaw University of Technology, Warsaw, Poland Department of Chemistry, University of Warsaw, Warsaw, Poland *Corresponding author: E-mail: [email protected]


Contents 1. Introduction 2. Outline of the Most Frequently Used Aromaticity Descriptors 2.1 Energetic Aromaticity Descriptors 2.2 Geometry-Based Descriptors 2.3 Magnetic Measures 2.4 Electronic Indices 3. Differences between Carbo- and Heterocyclic p-Electron Systems 4. Aromaticity of Nucleobases 5. Cyclic Hetero-p-Electron Systems as Test Probes for Statistical Multidimensionality of Aromaticity Characteristics 6. Conclusions Acknowledgments References

302 304 304 305 306 307 308 311 316 323 323 324

Abstract Heterocycles play a very important role in chemistry, biochemistry, pharmacy, and related fields; therefore, the ability to describe their p-electron delocalization is essential. The description and the application of the most effective and popular aromaticity indices based on energetic, structural, magnetic as well as electronic properties are presented herein for heterocycles. A multidimensional character of aromaticity has again been confirmed for five- and six-membered heterocyclic rings and for their fused systems indicating that utilization of one index only cannot be sufficient for the estimation of the aromatic character of a given system. In turn, the usefulness of the Kohonen neuron networks for the classification of the p-electron cyclic systems into three sets, namely aromatic, nonaromatic, and antiaromatic, has been shown. Although for all these groups some linear correlations between various aromaticity descriptors may exist, for each separated subgroup these relationships fail. Advances in Heterocyclic Chemistry, Volume 120 ISSN 0065-2725

© 2016 Elsevier Inc. All rights reserved.




Halina Szatylowicz et al.

Keywords: Aromatic character; Aromaticity indices; Azine; Azole; Heterocycles; Nucleobases; Pyridine; Pyrrole; Tautomers; Thiophene


Aromaticity index Aromatic stabilization energy Delocalization index Deoxyribonucleic acid Extra cyclic resonance energy Energy decomposition analysis Electron localization function Aromatic fluctuation index Harmonic oscillator model of aromaticity Harmonic oscillator model of electron delocalization Isomerization stabilization energy Multicenter bond index Magnetic resonance energy Natural bond orbital Nucleus-independent chemical shift Nuclear magnetic resonance Principal component analysis para-Delocalization index p-Electron-donor-acceptor descriptor Quantum theory of atoms in molecules Ring current Resonance energy s-Electron-donor-acceptor descriptor Self-organizing map Topological resonance energy

1. INTRODUCTION Aromaticity is a very important phenomenon but somehow mysterious and very often understood intuitively. This is one of the most frequently used terms in organic chemistry and related fieldsdevery day over 100 papers appear with “aromatic/aromaticity” words in a title, an abstract or keywords (ISI Web of Science, retrieved in February 2016). On the other hand, a definition of aromaticity is enumerative in nature (2000T1783), i.e., an aromatic compound has to be considered as fulfilling some particular physicochemical properties. It must be stressed that from the very beginning there have been two viewpoints for interpretation of the term “aromaticity”: either as a structural similarity to benzene (1865BSF98) or as a similarity of chemical properties of a system to the properties of benzene (1866ACP327); for a review, see Balaban et al. (2005CRV3436). Initially, the phrase “aromatic

Calculating Aromaticity of Heterocycles


compound” was related to its resistance in chemical reactions or, in other words, to a demonstration of a relatively high stability. For aromatic systems, this observation has allowed for their distinction from olefinic analogs. However, this fundamental property observed in aromatic systems has been for a long time only a qualitative feature. Only in the first half of twentieth century did the first quantitative criterion of aromaticity appear: Pauling and Sherman (1933JCP606) introduced the concept of resonance energy (RE), indicating the stability of the aromatic system (benzene) in comparison to its olefinic analog. Since that time many various quantitative measures of aromaticity based on structural (2001CRV1385, 2014CRV6383), magnetic (2005CRV3842, 2015CSR6597), energetic (2005CRV3773), and electronic (2005CRV3812, 2005CRV3911, 2015CSR6434) properties have been presented. Subsequently, the aromaticity concept built on carbocycles has been applied successfully to heterocyclic organic compounds, i.e., to a family of aromatic compounds which contain one or more heteroatoms as a part of their cyclic p-electron system. However, some of the measures developed for aromatic hydrocarbons, particularly those based on geometrical parameters, have been modified to be used for heterocycles (1993JCI70, 2001T5715, 2010SY1485). Indeed, a heteroatom embedded into a ring certainly affects its s- and p-electron system. Such incorporation effects can be evaluated by s-electron-donor-acceptor descriptor [sEDA(II)] and p-Electron-donoracceptor descriptor [pEDA(II)] characteristics, which have been developed especially for unsaturated five- and six-membered cycles with heteroatoms (2012JOC2608), or their later modification for fused heterocycles, sEDA(III) and pEDA(III) (2015JPO290). Many important biomolecules, including building blocks of DNA/ RNA and aromatic amino acids, contain heteroaromatic rings. Their characteristic structures contribute to their biological individuality and play an important role in the processes occurring with their participation. For example, it has been found that tautomeric preferences of systems which can exist in different tautomeric forms can often be explained by their pelectron delocalization. However, functional groups can significantly modify aromaticity of heterocycles, unlike p-electron system of benzene which is almost resistant to substituent effects (2004JOC6634). The concept of aromaticity is constantly evolving with the development of new compounds which can be considered as aromatic, and also with the discovery of new types of aromaticity, such as s-, p- (1988JSTT93),


Halina Szatylowicz et al.

d-aromaticity (2007AGE4277), as well as metalloaromaticity (2001CCR957, 2013WCMS105) and others. Despite the large number of methods applicable to measure aromaticity, any particular one is not sufficient to characterize an aromatic system unequivocally. The main thing is that they do not always lead to the same results within a series of compounds. Therefore, a complete description is usually provided by taking into account properties such as stabilization energy, equalization of bond lengths, specific magnetic behavior, etc. The most commonly used descriptors of heteroaromaticity are described in Section 2.

2. OUTLINE OF THE MOST FREQUENTLY USED AROMATICITY DESCRIPTORS Some difficulties in the quantitative description of aromaticity of heterocyclic systems result from the disparities between carbon atoms and heteroatoms involved in the p-electron system, occurring mostly due to the differences in the electronegativity, size, and valency. This is particularly important for energy- and geometry-based characteristics, because both of them use reference systems and hence encounter some obstacles in their proper choice.

2.1 Energetic Aromaticity Descriptors Consider first the concept of aromatic stabilization energy (ASE) (2005CRV3773), which is a modern representation of the Pauling RE. The ASE is based on the energy of isodesmic or homodesmotic reactions (Scheme 1), and the idea can be presented by a simple Eqn (1). The systems with high positive ASE values are aromatic, whereas these with high negative ASE values are antiaromatic. A key problem is a proper choice of the reference systems which are not defined in a general way. ASE ¼ Eðaromatic systemÞ  Eðnonaromatic reference analogueÞ















Scheme 1 Isodesmic reactions for (a) five- and (b) six-membered heterocycles.

Calculating Aromaticity of Heterocycles


Depending on the reference reaction, nonaromatic molecules can significantly differ from each other, hence the energy difference between the conjugated and the reference system can contain not only stabilization effects due to the p-electron delocalization but also strain, hyperconjugation, and differences in the types of bonds and hybridization (2002OL2873). It has been shown that ASE values obtained by utilization of different schemes of isodesmic reactions may give completely uncorrelated results (2003T1657); therefore, the application of ASE values requires substantial caution. To avoid extraneous discordances, Schleyer and P€ uhlhofer (2002OL2873) have recommended use of the “isomerization stabilization energy (ISE),” which is based on the differences between total energies computed for only two species, a methyl derivative of the aromatic system and its nonaromatic exocyclic methylene isomer. Recently, the results for a set of heterobenzenes (2010JSTT47, 2014CTC20) and metallobenzenes (2014CEJ14885) have been described. In turn, Fernandez and Frenking (2007FD403) have proposed the use of the difference in p-components of the energy (DEp), calculated by energy decomposition analysis (EDA); for a review, see Bickelhaupt and Baerends (2000RCC1). The application of this method for metallobenzenes in the case of a cyclic molecule and its acyclic analog with the same number of diene conjugations has been described (2007CEJ5873).

2.2 Geometry-Based Descriptors The geometry-based aromaticity indices (AIs) are calculated by a direct use of either the bond lengths (R) or the Gordy’s bond orders defined as a function of R2 (1947JCP305) and then applied in the Bird index (I) (1985T1409). The Bird index depends on ring size and can be calculated for five- and six-membered rings (I5 and I6, respectively), as well as for the fused entities. Up to now, this structural index has been employed for a great number of hetero-p-electron systems with different heteroatoms (1996T10255, 2001T5715, 2011PCCP20536). Bond lengths were used by Julg and Françoise (1967TCA249) for the first time and subsequently in the harmonic oscillator model of aromaticity (HOMA) index (1972TL3839, 1993JCI70) or its further modifications, harmonic oscillator model of electron delocalization (HOMED) (2010SY1485) and harmonic oscillator model of heterocycle electron delocalization (HOMHED) (2012STC375) with special parametrization for heterocyclic systems; for a review, see Krygowski et al. (2014CRV6383). The main


Halina Szatylowicz et al.

idea of the geometry-based AIs is a numerical description of bond lengths alternation following the general Eqn (2): X 2 AI ¼ 1  a Rref  Ri (2) where a is a constant taking into account number of bonds in a system and a normalization coefficient to give AI ¼ 1.0 for a fully aromatic system and AI ¼ 0 for the nonaromatic ones. AI < 0 qualifies a molecule as an antiaromatic system. Ri denotes bond lengths taken into account for a given pelectron system, whereas Rref means bond length taken as a reference. The choice of the reference bond length is the main characteristic of various modifications of AI. The commonly used HOMA index has a drawback in the case of “unusual” heterocycles, because it must be additionally parametrized taking into account experimental bond lengths, which are not always available. Nevertheless, Zborowski et al. (2012STC595) have successfully used computationally obtained bond lengths for parametrization of the carboneboron bond. In the case of the HOMED index, bond lengths for the reference molecules and compounds under consideration are calculated by quantum-chemical methods at the same level of theory that may cancel out computational errors in the procedure of the HOMED estimation (2010SY1485). However, if a given AI is applied for a series of structurally similar systems, then the dependence on the choice of the reference systems becomes less important. It is noteworthy that geometry-based AIs, unlike other descriptors of aromaticity, may be applied to any p-electron fragment of any system, hence being an indispensable tool for studying aromaticity or more generally the p-electron delocalization of a particular structural fragment of the complex system. Detailed information on geometry-based AIs is presented in a review by Krygowski et al. (2014CRV6383).

2.3 Magnetic Measures Magnetism-based AIs result from the property of electrons to form a ring current when a molecule is exposed to a perpendicularly directed magnetic field (2000PNMCS1, 2001CRV1349). Since the character of the current depends in some way on the level of the delocalization of p-electrons in the system, experimentally and theoretically accessible magnetic properties, such as nuclear magnetic resonance (NMR) chemical shifts (1961JCS859, 2001CRV1301), magnetic susceptibility (1969JA1991, 1974CRV663), and others (2001CRV1349, 2005CRV3889), may serve for the estimation

Calculating Aromaticity of Heterocycles


of aromaticity of molecules. In the last two decades Schleyer’s concept of the nucleus-independent chemical shift, abbreviated as NICS (1996JA6317), has found great popularity and numerous applications. NICS is defined as the negative value of an isotropic magnetic shielding and can be computed at any point of space. Negative NICS values indicate the presence of a diatropic ring current, and therefore the aromatic character of the system. In contrast, rings with positive NICS values are regarded as antiaromatic. NICS(0) and NICS(1), calculated at the geometric center of the ring and 1 Å above, respectively, are the most widely used aromaticity descriptors (2005CRV3842). Different variations of NICS have been later suggested, like dissected NICS (1997JA12669), NICSp and NICSs, or its out-of-plane component NICS(1)zz and its p contribution NICS(1)pzz (2004PCCP273). The latter indices have been found to be more accurate in the description of the p-electron delocalization because they contain “pure” contribution of p-electrons. A further step in the development of NICS has been the NICS scan method (2006JOC883), also applicable to heterocyclic aromatic systems. It has been shown that the results obtained from the NICS curves are more informative and consistent with other AIs than a single isotropic NICS values. In 2014, a new NICS-based methodology (NICS-XY-scan) for the identification of local and global induced ring currents in polycyclic systems was introduced, and it has found an application for the description of heteroacenes (2014JOC11644, 2016JA1792). Details on the NICS evolution and applications can be found in a comprehensive review by Chen et al. (2005CRV3842). Recent developments in the field of magnetic criteria of aromaticity have been presented by Gershoni-Poranne and Stanger (2015CSR6597). Note particularly that NICS indices depend on the ring size of the system in question and cannot be used for a quantitative comparison of the aromaticity degree for molecules with different number of atoms in the ring (2006JOM4359).

2.4 Electronic Indices Slightly earlier than the appearance of NICS, attempts to construct numerical characteristics of the electron delocalization based on Bader’s “quantum theory of atoms in molecules” (QTAIM) (1990MI) or “electron localization function” (ELF) (1990JCP5397) had begun. The most important for aromaticity description are electron sharing indices. For example, the delocalization index (DI) provides a value d(A,B) which is the number of electrons delocalized or shared between atoms A and B. The next useful aromaticity


Halina Szatylowicz et al.

descriptor is the para-delocalization index (PDI) (2003CEJ400), which is an average value of all DIs of para-related atoms in the ring, but it can be applied only to six-membered rings. For planar systems, PDI values can be separated into s and p contributions. Benzene was found to be the most aromatic ring with PDI and PDIp values 0.101 and 0.093, respectively, whereas s contributions (obtained as a difference between PDI and PDIp) amount to ca. 0.008 electrons. The multicenter delocalization index (Iring) (1990STC423, 2000 PCCP3381), and its other versiondmulticenter bond index (MCI) (2005JPO706), can also be applied to heterocyclic compounds with a different size of the ring to quantify their aromaticity. For planar species, MCI can be additionally split into s- and p contributions (MCIs and MCIp, respectively). When computed in an aromatic ring, more positive value of the MCI means more aromatic character of the ring. Another popular electronic measure of aromaticity is the fluctuation index (FLU) (2005JCP014109). It measures electron fluctuation differences with respect to an aromatic reference system by comparing the contiguous electron delocalization indices (1999JPCA304) along cyclic structure. The FLU index is close to 0 in aromatic species and departs from this value in nonaromatic ones. The drawback in this case is similar to that for geometry-based indices, namely the dependence on the reference systems. The application of the electron delocalization indices as an aromaticity measure has been recently amply reviewed by Sola et al. (2010SY1156) and Feixas et al. (2015CSR6434).

3. DIFFERENCES BETWEEN CARBO- AND HETEROCYCLIC p-ELECTRON SYSTEMS Undoubtedly, the most common heteroaromatics are nitrogen-containing five- and six-membered rings, namely azoles and azines. Apart from the individual molecules, they are often present as a part of fused ring systems. Further, they may contain several nitrogen atoms, which can be categorized in one of two types: (1) azine-like, i.e., containing a lone electron pair lying in the plane of a molecule carrying its basicity and one 2pz electron which is a part of the whole p-electron system resembling the p-electron structure of benzene (Figure 1(a)) and (2) azole-like, where, besides the aforementioned type of nitrogen atom/atoms, an acidic NH unit exists (Figure 1(b)). The investigations of a sequential replacement of CH units in benzene by nitrogen atoms leading to hexazine have shown that aromaticity of the

Calculating Aromaticity of Heterocycles


Figure 1 Schematic comparison between orbitals involved in the p-electron system of (a) six- and (b) five-membered rings for carbo- and heterocycles.

entire set of azines estimated by NICS and extra cyclic resonance energies (ECREs) “is nearly the same and resembles that of benzene” (2010OL4824). The application of HOMA and HOMED indices fully support the former conclusion (2010SY1485, 2014THC129). In a series of sixmembered monoheterocycles, it has been found that different AIs based on structural, magnetic, energetic, and electronic properties indicate significant aromaticity of all heterocycles under consideration (2010JSTT47, 2011PCCP20536, 2014CTC20). Despite the important inconsistency between the different indices, benzene and pyridine have been found to be the most aromatic species. Recently, Stojanovic and Baranac-Stojanovic (2016JOC197) have studied the effect of the same replacement (CH unit by the nitrogen atom) on the aromatic character of three isomeric azaborines (Scheme 2). Based on the HOMA, NICS(0)pzz, PDI, and ECRE AIs, they have stated that the replacement slightly affects aromaticity of the system. Other well-known aromatic heterocyclic compounds, like pyrrole, furan, and thiophene, can be regarded as a butadiene unit and a heteroatom bearing lone pairs built-in a five-membered ring. Hence their electronic structure is similar to the cyclopentadienyl anion (Figure 1(b)). Nevertheless, the cyclopentadienyl anion has the highest aromaticity compared to the heteroatom analogs. According to the results reported by Najmidin et al. (2013JMM3529), neither pyrrole, furan, and thiophene nor their aza-derivatives are stabilized by the replacing of the CH unit in the cyclopentadienyl anion with one or more nitrogen atoms.


Halina Szatylowicz et al.

Scheme 2 Effect of CH substitution with N atom on aromaticity of azaborines. Reprinted with permission from Stojanovic and Baranac-Stojanovic (2016JOC197). Copyright 2016 American Chemical Society.

The difference in p-electron systems of benzene and six-membered heterocycles can be determined through the examination of substituent effects. For example, substituents in pyridine (2013STC725) affect its aromaticity much more strongly than in benzene for which aromaticity weakly depends on the substituent (2004JOC6634). The same can also be stated for N-substituted azoles. The less aromatic pyrrole ring is much more sensitive to substituent effects than the more aromatic benzene ring (2007STC797). Similarly, the N-substitution of imidazole and pyrazole dramatically influences the aromaticity of their rings (2007STC965, 2011JPCA8571). A common feature is the large decrease of aromaticity in the case of substituents with electron-withdrawing properties, whereas the effect of electron-donating substituents is less pronounced (2015MI1). Alonso et al. (2011PCCP20564) have successfully applied the neural network methodology in the study of the substituent effect on aromaticity for a set of pyrimidine derivatives with a potential pushepull character. The interplay between aromaticity, planarity, steric effect, and charge transfer properties of all substituted pyrimidine derivatives has been also discussed (2011PCCP20564). To study the difference between p-electron systems of hetero- and carbocycles, the pEDA(II) descriptor (2012JOC2608), which measures heteroatom incorporation effects on the p-electron system, is very helpful. Using natural bond orbital (NBO) analysis, it is possible to estimate the changes in population of s- and p-electron in a given system (2009JPO769), that is, the amount of electrons shifted to or withdrawn from the p valence orbitals of the core molecule by the heteroatom/heteroatomic group. The obtained

Calculating Aromaticity of Heterocycles


correlations show that the greater the p-electron withdrawing effect of the incorporated heteroatom the less aromatic is the system, and vice versa, and the greater the donating effect the system becomes more aromatic. The broad investigations of aromaticity in heterocyclic analogs of benzene have pointed out a lack of correlations between different aromaticity characteristics based on structural aspects, the electron delocalization, magnetic and energetic properties (2010JSTT47, 2011PCCP20536). In contrast, a comparative study by Najmidin et al. (2013JMM3529) of aromaticity in pyrrole, furan, thiophene, and their aza-derivatives has shown that among seven indices applied, only the NICS(1) index does not correlate with the others. In the connection with such contradictory information, a detailed study on interrelation between AIs or, in other words, on statistical multidimensionality of the aromaticity characteristics is given in the following section.

4. AROMATICITY OF NUCLEOBASES It is well known that aromatic stacking interaction is one of the most important forces stabilizing the double-stranded DNA structure. Therefore, aromaticity of individual nucleobases has attracted attention for many years. Since the early 2000s different AIs, based on structural, magnetic and electronic properties, have been applied to the most known canonical forms of DNA bases present in the B-DNA structure. It has been found that the aromatic character of nucleobases differ from each other. According to the structural index HOMA and the magnetic index NICS, aromaticity increases in the following sequence: thymine < cytosine < guanine < adenine (2003JOC8607, 2005JSTT29). The same order has been also obtained using the electronic AIs FLU and PDI (2006JPCA12249). The canonical 9H form of adenine is the most aromatic species among all DNA bases. In the WatsoneCrick pair, aromaticity of the whole adenine molecule as well as of its six- and five-membered rings does not change in any significant way (2003JOC8607). In turn, in DNA bases containing at least one carbonyl group the formation of WatsoneCrick pairs leads to an aromaticity increase, because this group participates directly in Hbonding which weakens its p-electrons localization capability. Hydrogen bonds do not have a large effect on the aromaticity of the five-membered rings, whereas the p-electron structure of the six-membered rings with exocyclic functional groups undergoes perceptible changes. Scheme 3 shows the effect of the functional groups on the aromaticity of the most stable tautomers of DNA bases and their parent molecules


Halina Szatylowicz et al.

Scheme 3 Effects of fusion and substituents on aromaticity of nucleobases; red numbersdHOMA index. HOMA, harmonic oscillator model of aromaticity.

(2015MI). Obviously, the aromatic character of imidazole and pyrimidine rings hardly changes due to their fusion into one purine molecule if the most stable 9H tautomer is considered. The effect of the amino group on aromaticity of both rings in adenine is negligible. In contrast, introduction of the double-bonded carbonyl group induces an outflow of the electron density from the ring and leads to a dramatic decrease of aromaticity in guanine and cytosine. This effect is more pronounced in thymine where two carbonyl groups are present. If the molecule can exist in several tautomeric forms, aromaticity is one of the most important factors that can dictate the tautomeric preference (2005CRV3561). However, the first results for selected tautomers of DNA bases has shown that for such complex systems there is no simple relation between the p-electron delocalization and the stability of tautomers (2003JOC8607). Recently, Raczy nska et al. (2013CTC35) have studied the role of the p-electron delocalization in the tautomerism of purine (2010JPO828, 2013JMM3947), adenine (2013CTC35, 2014JMM2234) and, earlier uracil (2009STT103). All possible major, minor and rare tautomers in neutral, oxidized, and reduced forms have been investigated. Aromaticity of the neutral NH tautomers is the greatest, whereas the CH tautomers are less delocalized. The one-electron oxidation/reduction only slightly affects the p-electron delocalization of the NH forms, whereas aromaticity of the CH forms varies significantly upon electron transfer. This indicates

Calculating Aromaticity of Heterocycles


Figure 2 Correlations between the global HOMED indices and the relative energies (DE in kcal/mol) estimated for neutral isomers of adenine and its building blocks in the gas phase (a) and in water solution (b). HOMED, harmonic oscillator model of electron delocalization; DFT, density functional theory; PCM, polarizable continuum model. Reprinted ska et al. (2015RSCA36587) with permission from The Royal with permission from Raczyn Society of Chemistry.

that the p-electron system of CH tautomers is more sensitive to electron transfer than in the case of NH forms. It is worth mentioning that the results for adenine and its parent compounds (purine, imidazole, 4-aminopyridine) support the statement that aromaticity of the system determines its tautomeric preference (2015RSCA36587). The relative energies of prototropic tautomers correlate well with the global HOMED indices (Figure 2). The preferred tautomers are always highly aromatic, with HOMED values close to 1. Similar relationships have been found in the gas phase and in aqueous (polarizable continuum model) solution. The correlations are better for imidazole and purine because of the absence of exocyclic functional groups, which bring additional intramolecular interactions with the p-electron system. Interestingly, the presented dependences are steeper for imidazole and 4-aminopyrimidine than for systems with fused rings, i.e., adenine and purine. However, the thymine tautomers, with lack of correlation between HOMA index and relative energies of the tautomers, are an instructive


Halina Szatylowicz et al.

Figure 3 Dependence of HOMA on relative energy, Erel, for thymine tautomers. Group (a) blue (black in print versions) circles, (b) red (light gray in print versions) diamonds, and (c) crosses. HOMA, harmonic oscillator model of aromaticity. Reproduced by permission of The Royal Society of Chemistry from Stasyuk et al. (2014OBC6476).

example of an exception to the rule that the most stable tautomer is the most aromatic species (2014OBC6476). When six-membered heterocyclic systems are substituted, then interactions between substituents and ring nitrogen atoms (N or NH) lead to unusual energetic effects (Figure 3). The stability of the tautomers is determined by additional intramolecular interactions (attractive or repulsive) between neighboring bond dipoles, usually with participation of the functional groups (Figure 4). Higher numbers of attractive interactions in the tautomer leads to its higher stability. This case can also be extended to other nucleobases to describe the stability of their particular tautomers. In the 9H and 7H tautomers of purine and adenine (the NH unit is located in the five-membered ring), both rings fulfill the H€ uckel rule and, hence, their aromatic character is greater than in the 1H and 3H tautomers, as supported by HOMA and NICS values (2012JOC4035, 2014OBC456).

Calculating Aromaticity of Heterocycles


Figure 4 Structures of the most stable thymine tautomers and natural bond orbital atomic charges.

Aromaticity of guanine and cytosine tautomers strongly depends on the presence of C]O or C]NH groups. The tautomers with hydroxyl and amino groups are significantly more aromatic than their keto and imino analogs (2014CCA335, 2016STC111, 2016STC133). The main conclusion drawn from the study of the relationships between tautomeric equilibria and the p-electron delocalization of nucleobases is that their most stable tautomers are not always the most aromatic ones. This difference obviously originates from the effects caused by functional groups attached to the aromatic rings. Only exocyclic groups attached via a double bond (carbonyl and imino) induce a significant decrease in the p-electron delocalization, whereas tautomers without these groups are usually highly aromatic. Thus, it can be stated that aromaticity of the nucleobases depends not on the number of the functional groups but rather on their electronic character. Aromaticity of particular rings in several most stable tautomers of DNA bases has been analyzed not only by the structural HOMA index but also by the magnetic NICS(0) and NICS(1) indices (2014CCA335, 2014OBC456, 2014OBC6476, 2016STC111). Although correlations between different aromaticity descriptors for particular nucleobases are not good, they are quite well correlated if all DNA bases are taken into account (Figure 5). This is consistent with the view that aromaticity is a multidimensional phenomenon; therefore, various aromaticity criteria for a restricted group of compounds may not correlate with each other (2002JOC1333). The effect of intermolecular interactions on the changes in aromaticity of nucleobases has been also investigated (2014CCA335, 2014OBC456, 2014OBC6476, 2016STC111). Table 1 demonstrates the ranges in which the HOMA index for the whole cyclic system (HOMAtot) of the most stable tautomers of DNA bases can change due to the formation of various types of complexes. Formation of neutral H-bonds, investigated by interactions with HF, does not much influence aromaticity of nucleobases. However, participation in charge-assisted H-bonds (with F) and complexation with metal


Halina Szatylowicz et al.

Figure 5 Correlations between two aromaticity indices calculated for the five- (empty signs) and six-membered rings of studied tautomers of DNA bases: (a) NICS(0) vs HOMA, cc(6) ¼ 0.889; (b) NICS(1) vs HOMA, cc(6) ¼ 0.931. HOMA, harmonic oscillator model of aromaticity; NICS, nucleus-independent chemical shift. Table 1 HOMAtot values for the most stable tautomers of DNA bases and their changes upon intermolecular interactions, M ¼ Li, Na, K* Mþ HF F Tautomer











9H-adenine 9H-guanine 1H, 3H-thymine 1H-cytosine

0.926 0.761 0.490

0.922 0.748 0.527

0.933 0.810 0.574

0.011 0.062 0.047

0.825 0.759 0.369

0.947 0.788 0.726

0.122 0.029 0.357

0.861 0.694 0.450

0.986 0.874 0.727

0.125 0.180 0.277











* Data taken from Stasyuk (2015MI).

cations (Liþ, Naþ, Kþ) usually leads to a relatively large magnitude of aromaticity variation. Thus, the strongest intermolecular interactions usually lead to the greatest aromaticity changes. Effects of intermolecular interactions are very similar for purine and pyrimidine bases and are more pronounced in the six-membered rings due to the functional groups attached to these rings.

5. CYCLIC HETERO-p-ELECTRON SYSTEMS AS TEST PROBES FOR STATISTICAL MULTIDIMENSIONALITY OF AROMATICITY CHARACTERISTICS The first widely planned study of heteroaromaticity originated from the Katritzky group (1989JA7). At the beginning, for this purpose they

Calculating Aromaticity of Heterocycles


selected molecules (benzene and 15 five- and six-membered hetero-p-electron systems) for which several AIs based on geometrical, energetic, and magnetic characteristics were available. Using a principal component analysis (PCA) they demonstrated that description of 83% of the total variance of the AIs needed three orthogonal vectors PC1, PC2, and PC3, accounting for 47%, 22%, and 13%, respectively. The first two of them correlated well with geometry- and magnetism-based indices and were named as the “classical” and “magnetic” characteristics of aromaticity. The need to use three orthogonal PC vectors indicates unambiguously a statistical multidimensionality of aromaticity described by these indices. The study has been later extended on wider group of heteroaromatics, including the benzofused systems (1990JPR853, 1990JPR870, 1990JPR885). The above mentioned works have contributed to the beginning of a vivid disputation concerning an ability for describing the p-electron delocalization in a ring by one descriptor of aromaticity only, which would imply that the AIs are mutually correlated. Herein it should be stressed that heterocyclic systems have been mostly selected as objects of the investigation. At the beginning, the results obtained by Jug and Koester (1991JPO163) for heteroaromatic systems, as well as by Krygowski et al. (1995JCI203) for the polycyclic benzenoid hydrocarbons, have agreed with those reported previously by Katritzky et al. (1989JA7, 1990JPR853, 1990JPR870, 1990JPR885). Schleyer et al. (1995AGE337) have presented a different view. They have argued that the geometry-based Julg index (A), ASE, and magnetic susceptibility exaltation (L) are highly correlated (with R2 w 0.99) for fivemembered ring systems (C4H4X with X ¼ CHþ, SiHþ, BH, AlH, CH2, PH, SiH, O, S, NH, CH). Thus, the applied indices yielded quantitatively the same order of aromaticity (antiaromaticity) for the studied systems. However, the Julg index has been calculated by applying only three CC bonds of the molecule containing five bonds (two CX bonds). A further important consequence has been a suggestion (1996JA6317) that magnetic criteria are the most useful for describing aromaticity as “the diamagnetic susceptibility exaltation (L) is uniquely associated with aromaticity” and the introduction of a new magnetism-based indexdNICS. A plot of NICS against ASE for 10 C4H4X ring systems, with cc ¼ 0.966, was shown as a support of this new aromaticity descriptor. Katritzky et al. (1998JOC5228), inspired by the paper by Schleyer et al. (1995AGE337), have extended the set of 11 C4H4X ring systems by other 9. The obtained results show no correlation between L and ASE (R2 ¼ 0.034).


Halina Szatylowicz et al.

Figure 6 Dependence between L, NICS, NICS(1), and HOMA vs ASE for all 105 structures: (a) exaltation of magnetic susceptibility vs ASE; (b) NICS vs ASE; (c) NICS(1) vs ASE; (d) HOMA vs ASE. ASE, aromatic stabilization energy; HOMA, harmonic oscillator model of aromaticity; NICS, nucleus-independent chemical shift; L, magnetic susceptibility exaltation. Reprinted with permission from Cyranski et al. (2002JOC1333). Copyright 2002 American Chemical Society.

Undoubtedly, the result of correlation between AIs clearly depends on the choice of systems taken into consideration. Finally, the study of this problem has been undertaken in a joint work by the authors of the above controversy (2002JOC1333). AIs such as ASE, L, NICS and NICS(1) as well as HOMA were applied to a set of 75 polyhetero-fivemembered p-electron systems (including aza- and phospha-derivatives of furan, thiophene, pyrrole, and phosphole) and a set of 30 mono-heterofive-membered rings. The obtained results (Figure 6) show the character of the mutual relationships, which is confirmed by the correlation coefficients


Calculating Aromaticity of Heterocycles

Table 2 Correlation coefficients for correlation between L, NICS, NICS(1), HOMA, and ASE for all compounds* ASE L NICS NICS(1)


0.828 0.941 0.922 0.900

(102) (102) (102) (39)

0.891 (102) 0.881 (102) 0.827 (39)

0.975 (102) 0.845 (41)

0.856 (41)

* The sample size is given in parentheses; data taken from Cyra nski et al. (2002JOC1333).

(Table 2). Therefore, for all the compounds in question, some kind of a general trend in correlations between AIs exists. Hence, aromaticity may be regarded in some way as the one-dimensional phenomenon. What is more, the whole set of compounds can be divided into three major groups: aromatic, nonaromatic, and antiaromatic. When one narrow set of compounds is considered, then the results are much worse, as shown by the data for aromatic compounds only (with ASE >5 kcal/mol) in Table 3. Ramsden (2010T401) has detected statistically significant relationships between aromaticity characterized by ASE, HOMA, or NICS(1) indices and the number and positions of nitrogen atoms in the ring of aza-derivatives of pyrrole, furan, and thiophene. It has been demonstrated that azasubstitution at positions 2 and/or 5 increases aromaticity measured by all three indices. However, for aza-substitution at positions 3 and/or 4, two trends are observed. The NICS(1) index suggests an increase but ASE and HOMA indices indicate a decrease of aromaticity of the ring (Figure 7) with the growing number of nitrogen atoms. Najmidin et al. (2013JMM3529) have examined aromaticity for the same aza-derivatives using topological resonance energy (TRE), magnetic resonance energy (MRE), ring current (RC), and ring current diamagnetic susceptibility (cG) methods. An analysis of the results and the ones obtained by others [ASE, HOMA and NICS(1)] has been performed, but Table 3 Correlation coefficients for correlation between L, NICS, NICS(1), HOMA, and ASE for all compounds with ASE >5 kcal/mol* ASE L NICS NICS(1)


0.062 0.611 0.405 0.733

(66) (66) (66) (27)

0.383 (66) 0.207 (66) 0.700 (27)

0.761 (66) 0.809 (27)

0.516 (27)

* The sample size is given in parentheses; data taken from Cyranski et al. (2002JOC1333).


Halina Szatylowicz et al.

Figure 7 The variation of (a) NICS(1) and (b) ASE with the number and position of ring nitrogen atoms in pyrrole and its aza-derivatives. ASE, aromatic stabilization energy; NICS, nucleus-independent chemical shift. Reprinted with permission from Ramsden (2010T401). Copyright 2009 Elsevier Ltd.

excellent correlations are only found between TRE, MRE, RC, and cG (Table 4); that is, between energetic and magnetic criteria of aromaticity based on graph-theoretical approach (2006JA2873). Furthermore, Aihara (2008BCJ241), based on this approach, has postulated the one-dimensional character of aromaticity although these indices are physically dependent. Additionally, no correlation between NICS(1) and TRE has been observed. Among a set of 10 indicators of aromaticity (structural, magnetic, and electronic), the indices based on electron delocalization have been found to be the most precise for describing aromaticity in a series of 15 tested molecular systems (2008JCC1543). The results obtained show that an expected order of aromaticity is better reproduced for the five-membered heterocyclic species (C4H4X, with X ¼ CH, NH, O, CH2, BH, and CHþ) than for the six-membered ones (differing by the number and positions of nitrogen Table 4 Correlation matrix for aromaticity indices estimated for pyrrole, furan, thiophene, and their aza-derivatives* ASE HOMA NICS(1) TRE MRE RC


0.754 0.517 0.794 0.786 0.787 0.787

0.434 0.871 0.871 0.871 0.868

0.666 0.668 0.668 0.673

1.000 1.000 0.999

* The data for cc calculations taken from Najmidin et al. (2013JMM3529).

1.000 1.000



Calculating Aromaticity of Heterocycles

Table 5 Correlation coefficients between FLU1/2, MCI, INB, Iring, HOMA, ASE, NICS(0), NICS(1), and NICS(1)zz for five-membered heteroaromatic species* FLU1/2 MCI INB Iring ING HOMA ASE NICS(0) NICS(1) MCI INB Iring ING HOMA ASEx NICS(0) NICS(1) NICS(1)zz

0.137 0.083 0.883 0.230 0.963 0.802 0.041 0.840 0.631 0.713 0.095 0.948 0.818 0.934 0.958 0.059 0.931 0.975 0.824 0.783 0.868 0.114 0.850 0.920 0.685 0.812 0.770 0.116 0.828 0.899 0.653 0.816 0.749 0.102 0.843 0.906 0.673 0.822 0.761

0.968 0.953 0.998 0.959 0.999


* Values of aromaticity indices taken from Feixas et al. (2014THC129). x Data ASE taken from Cyra nski et al. (2002JOC1333).

atoms). Moreover, these findings are supported by mutual correlation between aromaticity descriptors (Tables 5 and 6). An application of different aromaticity descriptors (e.g., ASE, rRCP, I6, and NICS(1)zz) to the eight heterocyclic analogs of benzene (C5H5X with X ¼ SiH, GeH, N, P, As, Oþ, Sþ, and Seþ) leads to the conclusion of a lack of mutual correlations between them (2011PCCP20536). Ebrahimi et al. (2010JSTT47) have come to the same results for aromaticity of groups IIIA to VIA heterobenzenes (24 systems), although, within these groups, the correlation between aromaticity parameters is satisfactory or even good. The first successful classification of compounds according to aromaticity criteria has been provided by Alonso and Herradon (2007CEJ3913, 2010JCC917), who used for this purpose the Kohonen neural networks, a type of artificial neural networks applied also in chemistry (for detail see (2003MI) and (2005JCIM264)). The application of self-organizing maps (SOMs) (2001MI) shows a good correlation between aromaticity of a compound and its placement in a particular neuron. Since the position of the Table 6 Correlation coefficients between PDI, FLU1/2, HOMA, NICS(0), NICS(1), and RE for six-membered nitrogen derivatives of pyridine* PDI FLU1/2 HOMA NICS(0) NICS(1)


0.038 0.479 0.347 0.844 0.361

0.185 0.942 0.141 0.677

0.338 0.856 0.384

* Values of aromaticity indices taken from Feixas et al. (2014THC129). x Data RE taken from Cyranski (2005CRV3773).

0.413 0.519



Halina Szatylowicz et al.

Figure 8 Sammon map obtained for the training data set, showing the relative distances between the input variables [ASE, L, NICS(1)zz, and HOMA] in the original space. The color scale indicates the Euclidean distances between the weight vector of each neuron and the neuron activated by benzene. ASE, aromatic stabilization energy; HOMA, harmonic oscillator model of aromaticity; NICS, nucleus-independent chemical shift; L, magnetic susceptibility exaltation. Reprinted with permission from Alonso and Herradon (2010JCC917). Copyright 2009 Wiley Periodicals, Inc.

compound in the map depends on its aromatic character so a new quantitative scale of aromaticity, based on the Euclidean distance between neurons in an SOM, has been established. At the beginning, as input data for the training of the network, ASE, magnetic susceptibility exaltation (L), NICS, and HOMA indices for five-membered heterocycles (105 molecules) and the cyclopentadienyl cation were used (2007CEJ3913). Then, the trained Kohonen network was applied to a few selected five-membered rings and sixmembered rings, i.e., the compounds which were not used for the training of the neural network. This methodology was then expanded to a full set of organic compounds ranging from highly aromatic to highly antiaromatic rings (2010JCC917). Hence, this method can be used to quantify the p-electron delocalization for molecules of different sizes (from five- to eightmembered rings), both carbocyclic and heterocyclic, including neutral, cationic, or anionic species. Sammon map (1969IEEETC11644) in Figure 8 illustrates the capability of this method. Following the conclusion of Cyranski et al. (2002JOC1333), the computed correlations between various aromaticity descriptors for subgroups (aromatic, nonaromatic, and antiaromatic) are much worse than for the whole set of data (see Supporting Information by Alonso and

Calculating Aromaticity of Heterocycles


Herradon (2010JCC917)). The application of PCA for all data leads to one principal component explaining 90.7% of the whole variance. When only aromatic compounds are considered, four principal components are necessary to describe 90.5% of the whole variance.

6. CONCLUSIONS Many important chemical compounds such as biomolecules, medicines, explosives, and others contain p-electron fragments with one or more heteroatoms. Various properties of such molecules hinge on the pelectron delocalization, which in turn is strongly dependent on the presence of heteroatoms, their numbers, and character. The level of the p-electron delocalization in cyclic systems is usually estimated by AIs, which, however, not always are mutually correlated. For more homogeneous systems and large variability of the indices, higher mutual correlations may be found. However, for a smaller set of compounds (for example, classified as aromatic, nonaromatic, or antiaromatic), linear regressions between aromaticity parameters characterize rather small slopes, which are strongly biased by the random errors of the data, and hence the mutual correlations are worse. Here we have described the most popular aromaticity descriptors applicable to heteroaromatic systems and have considered possible correlations between them for five- and six-membered heterocycles. The relationship of aromaticity with stability of tautomers, substituent effects and intermolecular interactions in the case of nucleobases have been also discussed. Since many descriptors of aromaticity are obtained by quantum chemistry methods, their values depend on the computational level. Moreover, energetic AIs based on different schemes of isodesmic or homodesmotic reactions may give completely unrelated results. In the case of geometrybased indices, the final results are directly connected with the choice of reference bond lengths. Magnetism-based indices depend on molecular dimension (size of the ring) and hence may be applied only to molecules of the same size. Thus, for a good description of the aromatic character of any molecule and a comparison of the aromaticity degree in a set of compounds, the best way is to combine several aromaticity descriptors based on different physicochemical properties of molecules.

ACKNOWLEDGMENTS We thank Warsaw University of Technology and University of Warsaw for financial support.


Halina Szatylowicz et al.

REFERENCES 1865BSF98 1866ACP327 1933JCP606 1947JCP305 1961JCS859 1967TCA249 1969IEEETC11644 1969JA1991 1972TL3839 1974CRV663 1985T1409 1988JSTT93 1989JA7 1990JCP5397 1990MI 1990JPR853 1990JPR870 1990JPR885 1990STC423 1991JPO163 1993JCI70 1995AGE337 1995JCI203 1996JA6317 1996T10255 1997JA12669 1998JOC5228 1999JPCA304 2000PCCP3381 2000PNMCS1 2000RCC1

A. Kekulé, Bull. Soc. Chim. Fr., 3, 98e110 (1865). E. Erlenmayer, Ann. d. Chem. u. Pharm., 137, 327 (1866). L. Pauling and J. Sherman, J. Chem. Phys., 1, 606e617 (1933). W. Gordy, J. Chem. Phys., 15, 305e310 (1947). J.A. Elvidge and L.M. Jackman, J. Chem. Soc., 859e866 (1961). A. Julg and P. Françoise, Theor. Chim. Acta, 7, 249e259 (1967). J.W. Sammon, IEEE Trans. Comput., C-18, 401e409 (1969). H.J. Dauben, J.D. Wilson, and J.L. Laity, J. Am. Chem. Soc., 91, 1991e 1998 (1969). J. Kruszewski and T.M. Krygowski, Tetrahedron Lett., 13, 3839e 3842 (1972). W.H. Flygare, Chem. Rev., 74, 653e687 (1974). C.W. Bird, Tetrahedron, 41, 1409e1414 (1985). V.I. Minkin, M.N. Glukhovtsev, and B.Ya Simkin, J. Mol. Struct.: THEOCHEM, 181, 93e110 (1988). A.R. Katritzky, P. Barczy nski, G. Mussumura, D. Pisano, and M. Szafran, J. Am. Chem. Soc., 111, 7e15 (1989). A.D. Becke and K.E. Edgecombe, J. Chem. Phys., 92, 5397e 5403 (1990). R.F.W. Bader, Atoms in Molecules. A Quantum Theory, Clarendon: Oxford (1990). A.R. Katritzky, V. Feygelman, G. Musumarra, P. Barczynski, and M. Szafran, J. Prakt. Chem., 332, 853e869 (1990). A.R. Katritzky, V. Feygelman, G. Musumarra, P. Barczynski, and M. Szafran, J. Prakt. Chem., 332, 870e884 (1990). A.R. Katritzky and P. Barczynski, J. Prakt. Chem., 332, 885e 897 (1990). M. Giambiagi, M.S. de Giambiagi, and K.C. Mundim, Struct. Chem., 1, 423e427 (1990). K. Jug and A.M. Koester, J. Phys. Org. Chem., 4, 163e169 (1991). T.M. Krygowski, J. Chem. Inf. Comput. Sci., 33, 70e78 (1993). P.v.R. Schleyer, P.K. Freeman, H. Jiao, and B. Goldfuss, Angew. Chem. Int. Ed. Engl., 34, 337e340 (1995). T.M. Krygowski, A. Ciesielski, C.W. Bird, and A. Kotschy, J. Chem. Inf. Comput. Sci., 35, 203e210 (1995). P.v. R. Schleyer, C. Maerker, H. Dransfeld, H. Jiao, and N.J.R.v. E. Hommes, J. Am. Chem. Soc., 118, 6317e6318 (1996). T.M. Krygowski and M.K. Cyra nski, Tetrahedron, 52, 10255e 10264 (1996). P.v. R. Schleyer, H. Jiao, N.J.R.v. E. Hommes, V.G. Malkin, and O.L. Malkina, J. Am. Chem. Soc., 119, 12669e12670 (1997). A.R. Katritzky, M. Karelson, S. Sild, T.M. Krygowski, and K. Jug, J. Org. Chem., 63, 5228e5231 (1998). X. Fradera, M.A. Austen, and R.F.W. Bader, J. Phys. Chem. A, 103, 304e314 (1999). M. Giambiagi, M.S. de Giambiagi, C.D. dos Santos Silva, and A.P. de Figuetredo, Phys. Chem. Chem. Phys., 2, 3381e3392 (2000). P. Lazzeretti, Prog. Nucl. Mag. Res. Sp., 36, 1e88 (2000). F.M. Bickelhaupt and E.J. Baerends, Rev. Comput. Chem., 15, 1e86 (2000).

Calculating Aromaticity of Heterocycles

2000T1783 2001CCR957 2001CRV1301 2001CRV1349 2001CRV1385 2001MI 2001T5715 2002JOC1333 2002OL2873 2003CEJ400 2003JOC8607 2003MI 2003T1657 2004PCCP273 2004JOC6634 2005CRV3436 2005CRV3561 2005CRV3773 2005CRV3812 2005CRV3842 2005CRV3889 2005CRV3911 2005JCIM264 2005JCP014109 2005JPO706 2005JSTT29 2006JA2873 2006JOC883 2006JOM4359 2006JPCA12249


T.M. Krygowski, M.K. Cyra nski, Z. Czarnocki, G. Haefelinger, and A.R. Katritzky, Tetrahedron, 56, 1783e1796 (2000). H. Masui, Coord. Chem. Rev., 219, 957e992 (2001). R.H. Mitchell, Chem. Rev., 101, 1301e1315 (2001). J.A.N.F. Gomes and R.B. Mallion, Chem. Rev., 101, 1349e 1383 (2001). T.M. Krygowski and M.K. Cyra nski, Chem. Rev., 101, 1385e 1419 (2001). T. Kohonen, Self-Organizing Maps, 3rd ed., Springer: Berlin (2001). S.I. Kotelevskii and O.V. Prezhdo, Tetrahedron, 57, 5715e 5729 (2001). M.K. Cyra nski, T.M. Krygowski, A.R. Katritzky, and P.v. R. Schleyer, J. Org. Chem., 67, 1333e1338 (2002). P.v. R. Schleyer and F. P€ uhlhofer, Org. Lett., 4, 2873e2876 (2002). J. Poater, X. Fradera, M. Duran, and M. Sola, Chem., Eur. J., 9, 400e 406 (2003). M.K. Cyra nski, M. Gilski, M. Jask olski, and T.M. Krygowski, J. Org. Chem., 68, 8607e8613 (2003). J. Zupan and J. Gasteiger, Neural Networks for Chemists. An Introduction, Wiley-VCH: Weinheim (2003). M.K. Cyra nski, P.V.R. Schleyer, T.M. Krygowski, H. Jiao, and G. Hohlneicher, Tetrahedron, 59, 1657e1665 (2003). C. Corminboeuf, T. Heine, G. Seifert, P.v. R. Schleyer, and J. Weber, Phys. Chem. Chem. Phys., 6, 273e276 (2004). T.M. Krygowski, K. Ejsmont, B.T. Stepien, M.K. Cyra nski, J. Poater, and M. Sola, J. Org. Chem., 69, 6634e6640 (2004). A.T. Balaban, P.v. R. Schleyer, and H.S. Rzepa, Chem. Rev., 105, 3436e3447 (2005). E.D. Raczy nska, W. Kosi nska, B. Osmialowski, and R. Gawinecki, Chem. Rev., 105, 3561e3612 (2005). M.K. Cyra nski, Chem. Rev., 105, 3773e3811 (2005). G. Merino, A. Vela, and T. Heine, Chem. Rev., 105, 3812e 3841 (2005). Z.F. Chen, C.S. Wannere, C. Corminboeuf, R. Puchta, and P.v. R. Schleyer, Chem. Rev., 105, 3842e3888 (2005). T. Heine, C. Corminboeuf, and G. Seifert, Chem. Rev., 105, 3889e 3910 (2005). J. Poater, M. Duran, M. Sola, and B. Silvi, Chem. Rev., 105, 3911e 3947 (2005). J.J. Panek, A. Jezierska, and M. Vracko, J. Chem. Inf. Mod., 45, 264e 272 (2005). E. Matito, M. Duran, and M. Sola, J. Chem. Phys., 122, 014109 (2005). P. Bultinck, R. Ponec, and S. Van Damme, J. Phys. Org. Chem., 18, 706e718 (2005). P. Cysewski, J. Mol. Struct.: THEOCHEM, 714, 29e34 (2005). J.-I. Aihara, J. Am. Chem. Soc., 128, 2873e2879 (2006). A. Stanger, J. Org. Chem., 71, 883e893 (2006). J.O.C. Jiménez-Halla, E. Matito, J. Robles, and M. Sola, J. Organomet. Chem., 691, 4359e4366 (2006). O. Huertas, J. Poater, M. Fuentes-Cabrera, M. Orozco, M. Sola, and F.J. Luque, J. Phys. Chem. A, 110, 12249e12258 (2006).


2007AGE4277 2007CEJ3913 2007CEJ5873 2007FD403 2007STC797 2007STC965 2008BCJ241 2008JCC1543 2009JPO769 2009STT103 2010JCC917 2010JPO828 2010JSTT47 2010OL4824 2010SY1156 2010SY1485 2010T401 2011JPCA8571 2011PCCP20536 2011PCCP20564 2012JOC2608 2012JOC4035 2012STC375 2012STC595 2013CTC35 2013JMM3947 2013JMM3529 2013STC725 2013WCMS105

Halina Szatylowicz et al.

H.-J. Zhai, B.B. Averkiev, D.Y. Zubarev, L.-S. Wang, and A.I. Boldyrev, Angew. Chem. Int. Ed., 46, 4277e4280 (2007). M. Alonso and B. Herradon, Chem., Eur. J., 13, 3913e3923 (2007). I. Fernandez and G. Frenking, Chem., Eur. J., 13, 5873e5884 (2007). I. Fernandez and G. Frenking, Faraday Discuss, 135, 403e421 (2007). K.K. Zborowski, I. Alkorta, and J. Elguero, Struct. Chem., 18, 797e 805 (2007). F. Blanco, I. Alkorta, K. Zborowski, and J. Elguero, Struct. Chem., 18, 965e975 (2007). J.-I. Aihara, Bull. Chem. Soc. Jpn., 81, 241e247 (2008). F. Feixas, E. Matito, J. Poater, and M. Sola, J. Comput. Chem., 29, 1543e1554 (2008). W.P. Ozimi nski and J.C. Dobrowolski, J. Phys. Org. Chem., 22, 769e 778 (2009). E.D. Raczy nska, K. Zientara, K. Kolczy nska, and T. Ste˛ pniewski, J. Mol. Struct.: THEOCHEM, 894, 103e111 (2009). M. Alonso and B. Herradon, J. Comput. Chem., 31, 917e928 (2010). E.D. Raczy nska and B. Kami nska, J. Phys. Org. Chem., 23, 828e 835 (2010). A.A. Ebrahimi, R. Ghiasi, and C. Foroutan-Nejad, J. Mol. Struct. THEOCHEM, 941, 47e52 (2010). Y. Wang, J.I. Wu, Q.-S. Li, and P.v.R. Schleyer, Org. Lett., 12, 4824e 4827 (2010). M. Sola, F. Feixas, J.O.C. Jiménez-Halla, E. Matito, and J. Poater, Symmetry, 2, 1156e1179 (2010). E.D. Raczy nska, M. Hallman, K. Kolczy nska, and T. Ste˛ pniewski, Symmetry, 2, 1485e1509 (2010). C.A. Ramsden, Tetrahedron, 66, 2695e2699 (2010). C. Curutchet, J. Poater, M. Sola, and J. Elguero, J. Phys. Chem. A, 115, 8571e8577 (2011). I.V. Omelchenko, O.V. Shishkin, L. Gorb, J. Leszczynski, S. Fias, and P. Bulting, Phys. Chem. Chem. Phys., 13, 20536e20548 (2011). M. Alonso, C. Miranda, N. Martin, and B. Herradon, Phys. Chem. Chem. Phys., 13, 20564e20574 (2011). A. Mazurek and J.C. Dobrowolski, J. Org. Chem., 77, 2608e 2618 (2012). O.A. Stasyuk, H. Szatylowicz, and T.M. Krygowski, J. Org. Chem., 77, 4035e4045 (2012). C.P. Frizzo and M.A.P. Martins, Struct. Chem., 23, 375e380 (2012). K.K. Zborowski, I. Alkorta, J. Elguero, and L.M. Proniewicz, Struct. Chem., 23, 595e600 (2012). E.D. Raczy nska, K. Kolczy nska, T.M. Ste˛ pniewski, and B. Kami nska, Comput. Theor. Chem., 1022, 35e44 (2013). E.D. Raczy nska and B. Kami nska, J. Mol. Model., 19, 3947e 3960 (2013). K. Najmidin, A. Kerim, P. Abdirishit, H. Kalam, and T. Tawar, J. Mol. Model., 19, 3529e3535 (2013). I.V. Omelchenko, O.V. Shishkin, L. Gorb, F.C. Hill, and J. Leszczynski, Struct. Chem., 24, 725e733 (2013). F. Feixas, E. Matito, J. Poater, and M. Sola, WIREs Comput. Mol. Sci., 3, 105e122 (2013).

Calculating Aromaticity of Heterocycles

2014CCA335 2014CEJ14885 2014CRV6383 2014CTC20 2014JMM2234 2014JOC11644 2014OBC456 2014OBC6476 2014THC129

2015CSR6434 2015CSR6597 2015JPO290 2015MI 2015MI1 2015RSCA36587 2016JA1792 2016JOC197 2016STC111 2016STC133


O.A. Stasyuk, H. Szatylowicz, and T.M. Krygowski, Croat. Chem. Acta, 87, 335e342 (2014). R. Lin, K.H. Lee, K.C. Poon, H.H.Y. Sung, I.D. Williams, Z. Lin, and G. Jia, Chem., Eur. J., 20, 14885e14899 (2014). T.M. Krygowski, H. Szatylowicz, O.A. Stasyuk, J. Dominikowska, and M. Palusiak, Chem. Rev., 114, 6383e6422 (2014). B. Hou, P. Yi, Z. Wang, S. Zhang, J. Zhao, R.L. Mancera, Y. Cheng, and Z. Zuo, Comput. Theor. Chem., 1046, 20e24 (2014). E.D. Raczy nska and M. Makowski, J. Mol. Model., 20, 2234 (2014). M. Schaffroth, R. Gershoni-Poranne, A. Stanger, and U.H.F. Bunz, J. Org. Chem., 79, 11644e11650 (2014). O.A. Stasyuk, H. Szatylowicz, and T.M. Krygowski, Org. Biomol. Chem., 12, 456e466 (2014). O.A. Stasyuk, H. Szatylowicz, and T.M. Krygowski, Org. Biomol. Chem., 12, 6476e6483 (2014). F. Feixas, J. Poater, E. Matito, and M. Sola, Aromaticity of Organic and Inorganic Heterocycles, In F. De Proft and P. Geerlings, editors: Structure, Bonding and Reactivity of Heterocyclic Compounds, Topics in Heterocycl. Chem., Vol. 38 (2014), pp 129e160. F. Feixas, E. Matito, J. Poater, and M. Sola, Chem. Soc. Rev., 44, 6434e6451 (2015). R. Gershoni-Poranne and A. Stanger, Chem. Soc. Rev., 44, 6597e 6615 (2015). A. Mazurek and J.C. Dobrowolski, J. Phys. Org. Chem., 28, 290e 297 (2015). O.A. Stasyuk, Effect of Intermolecular Interactions on the p-electron Structure and Tautomerism of Nucleobases (Ph.D. Thesis), Warsaw University of Technology: Warsaw (2015). H. Szatylowicz, O.A. Stasyuk, and T.M. Krygowski, Advances in Heterocyclic Chemistry, In E.F.V. Scriven and C.A. Ramsden, editors, Academic Press, Elsevier, Vol. 116 (2015), pp 137e192. E.D. Raczy nska, M. Makowski, M. Hallmann, and B. Kami nska, RSC Adv., 5, 36587e36604 (2015). A.H. Endres, M. Schaffroth, F. Paulus, H. Reiss, H. Wadepohl, F. Rominger, R. Kramer, and U.H.F. Bunz, J. Am. Chem. Soc., 138, 1792e1795 (2016). M. Stojanovic and M. Baranac-Stojanovic, J. Org. Chem., 81, 197e 205 (2016). O.A. Stasyuk, H. Szatylowicz, and T.M. Krygowski, Struct. Chem., 27, 111e118 (2016). E.D. Raczy nska, M. Sapu1a, K. Zientara-Rytter, K. Kolczy nska, T.M. Ste˛ pniewski, and M. Hallmann, Struct. Chem., 27, 133e 143 (2016).