Matta, Boyd. The quantum theory of atoms in molecules

282 10 Interpretation of Experimental Electron Densities by Combination of the QTAMC and DFT

Nguyen, S.J. Su, T.L. Windus, M. Dupuis, J.A. Montgomery, ‘‘General atomic and molecular electronic structure system’’. J. Comput. Chem., 1993, 14, 1347–1363.

78http://www.chemistry.mcmaster.ca/ aimpac

79V.G. Tsirelson, Yu.A. Abramov, V.E. Zavodnik, A.I. Stash, E.L. Belokoneva,

J.Stahn, U. Pietsch, D. Feil, ‘‘Critical points in a crystal and procrystal’’.

Struct. Chem., 1998, 9, 249–254.

80R.F.W. Bader, J.A. Platts, ‘‘Characterization of an F-center in an alkali halide cluster’’. J. Chem. Phys., 1997, 107, 8545–8553.

81L.J. Farrugia, P.R. Mallinson, B. Stewart, ‘‘Experimental charge density in the transition metal complex

Mn2(CO)10 a comparative study’’.

Acta Cryst., 2003, B59, 234–247.

82E.A. Zhurova, V.G. Tsirelson, ‘‘Electron density and energy density view on the atomic interactions in SrTiO3’’.

Acta Cryst., 2002, B58, 567–575.

83E.A. Zhurova, V.G. Tsirelson, A.I. Stash, M.V. Yakovlev, A.A. Pinkerton, ‘‘Electronic energy distributions in energetic materials NTO and the biguanidinium dinitramides’’. J. Phys. Chem., 2004, B108, 20173–20179.

84C.H. Hodges, ‘‘Quantum corrections to the Thomas–Fermi approximation and the Kirzhnits method’’. Can.

J.Phys., 1973, 51, 1428–1437.

85H. Lee, C. Lee, R.G. Parr, ‘‘Cojoint gradient correction to the Hartree– Fock kineticand exchange-energy density functionals’’. Phys. Rev., 1991, A44, 768–771.

86C.F. Weizsaecker, ‘‘Zur Theorie der Kernmassen’’. Z. Phys., 1935, 96, 431–444.

87T.A. Wesolowski, F. Tran, ‘‘Gradientfree and gradient-dependent approximations in the total energy bifunctional for weakly overlapping electron density’’. J. Chem. Phys., 2003, 118, 2072–2080.

88E.J. Baerends, O.V. Gritsenko, ‘‘A quantum chemical view of densityfunctional theory’’. J. Phys. Chem., 1997, A101, 5383–5403.

89A.D. Becke, ‘‘Density functional exchange-energy approximation with correct asymptotic behavior’’. Phys. Rev., 1988, A38, 3098–3100.

90J.P. Perdew, K. Burke, M. Ernzerhof, ‘‘Generelized Gradient Approximation Made Simple’’. Phys. Rev. Lett., 1996, 77, 3865–3868.

91J.P. Perdew, Y. Wang, ‘‘Accurate and simple analytic representation of the electron-gas correlation energy’’. Phys. Rev., 1992, B45, 13244–13249.

92V. Tsirelson, A. Stash, ‘‘On functions and quantities derived from the experimental electron density’’. Acta Cryst., 2004, A60, 418–426.

93E. Espinosa, M. Souhassou, H. Lachekar, C. Lecomte, ‘‘Topological analysis of the electron density in hydrogen bonds’’. Acta Cryst., 1999, B55, 563–572.

94P.R.T. Schipper, O.V. Gritsenko, E.J. Baerends, ‘‘Kohn–Sham potential and exchange and correlation energy densities from oneand twoelectron density matrices for Li2, N2, and F2’’. Phys. Rev., 1998, A57, 1729– 1742.

95V. Tsirelson, A. Stash, ‘‘Determination of the electron localization function from electron density’’.

Chem. Phys. Lett., 2002, 351, 142–148.

96V. Tsirelson, A. Stash, ‘‘Analyzing experimental electron density with localized-orbital locato’’. Acta Cryst., 2002, B58, 780–785.

97A. Stash, V. Tsirelson, ‘‘WinXPRO a program for calculating crystal and molecular properties using

multipole parameters of the electron density’’. J. Appl. Cryst., 2002, 35, 371–373.

98A.I. Stash, V.G. Tsirelson, ‘‘Modern opportunities of calculating physical properties of crystals using experimental electron density’’. Cryst. Report., 2005, 50, 209–216.

99P. Macchi, P. Coppens, ‘‘Relativistic analytical wavefunctions and scattering factors for neutral atoms beyond Kr and for all chemically important ions up to I ’’. Acta Cryst., 2001, A57, 656–662.

100S.S. Iyengar, M. Ernzerhof, S.N. Maximo , G.E. Scuseria, ‘‘Challenge of creating accurate and e ective kinetic-energy functionals’’. Phys. Rev., 2001, A63, 052508.

101C. Lee, Zh. Zhou, ‘‘Exchange-energy density functional reparametrization

References 283

of Becke’s formula and derivation of second-order gradient correction.’’ Phys. Rev., 1991, A44, 1536–1539.

102D.P. O’Neill, P.M.W. Gill, ‘‘Benchmark correlation energies for small molecules’’. Mol. Phys., 2005, 103, 763–766.

285

11

Topological Analysis of Proteins as Derived from Medium and High-resolution Electron Density: Applications to Electrostatic Properties

Laurence Leherte, Benoıˆt Guillot, Daniel P. Vercauteren,

Virginie Pichon-Pesme, Christian Jelsch, Ange´lique Lagoutte,

and Claude Lecomte

11.1 Introduction

The details observable in protein crystal structures depend on experimental diffraction resolution. In X-ray di raction (XRD), the experimental resolution (A˚ ) is defined by d ¼ l=ð2 sin yÞ where y is the Bragg angle and l the wavelength. At the usual resolution for macromolecular crystallography 1:6 < d < 3 A˚ , the number of Fourier components Fðh; k; lÞ of the electron density (ED) as obtained from:

rexpðx; y; zÞ ¼ V 1 XXXjFðh; k; lÞj exp ijðh; k; lÞexp ½2piðhx þ ky þ lzÞ&

h k l

ð1Þ

is not enough for estimation of all the atomic properties. Crystallographic refinement must generally be completed by using a standard stereochemistry for the macromolecule.

Another possibility, instead of using stereochemistry, is to use topological tools. Indeed, topological analysis of the ED can be applied at di erent resolution levels. For the lower-resolution structures the density maxima often do not correspond to atomic positions but should rather be associated with atomic groups. For example, at@3-A˚ resolution, we and others have shown that for protein crystallography there is a correspondence between ED peaks and the backbone or side-chain atoms of the amino acid residues [1–3]. For example, Guo et al. [2] showed that 3-A˚ resolution ED distributions could be reconstructed from structure factors generated by residue centers. Topological analysis of a medium resolution ED will, therefore, mostly deal with the ð3; 3Þ critical points (CPs), which are the maxima of the ED distribution, and their correspondence with chemical groups to give a geometrically meaningful representation of the protein. Such represen-

286 11 Topological Analysis of Proteins as Derived from Medium and High-resolution Electron Density

tations built on ð3; 3Þ CP properties have, for example, been used in protein– protein and protein–DNA docking applications [4].

At higher resolution, d < 2:5 A˚ , one can observe atomic details and an isotropic atomic Debye Waller factor Bj ¼ 8p2huj2i, where huj2i is the mean square atomic displacement, can be estimated. The dynamic structure factor can be written:

Nat
FdynðHÞ ¼ Xj 1	fj ðHÞ exp 2piðH:rjÞ exp ðBj sin2 y=lÞ	ð2Þ
¼

where Nat stands for the number of atoms in a unit cell, and fj ðHÞ stands for the atomic scattering factors, i.e. the Fourier transform of the free and neutral atom ED (independent atom model, IAM).

When d < @1:4 A˚ , the isotropic temperature factor may be replaced by an anisotropic factor for non hydrogen atoms:

Nat
FdynðHÞ ¼ Xj 1	fj ðHÞ exp 2piðH:rjÞ exp½ 2p2ðhihkU ijai :ak Þ&	ð3Þ
¼

where U ij are the tensor elements of the anisotropic atomic displacement and ai stands for the ith reciprocal dimension of the crystal unit cell.

The accuracy of the resulting ð3; 3Þ CPs located at the atomic positions is good enough to validate deviations from the standard geometry and evidence the ED of some hydrogen atoms, especially those with a small B factor. At a resolution value estimated to be d < 0:9 A˚ , di raction data, when accurate enough, contain information at a subatomic scale. Thus, information on valence ED distribution may be obtained when the anisotropic thermal displacement values are small enough, corresponding to an equivalent isotropic B factor lower than about 5 A˚ 2 [5].

The ð3; 3Þ CPs correspond to nuclei and therefore give an atomic description of the protein structure; hydrogen atoms clearly show up. Deviations from the spherical free atom ED, i.e. deviations from the IAM, appear in the dynamic experimental deformation ED maps as ED peaks located on the chemical bonds:

X	ð4Þ
DrdynðrÞ ¼ V 1 ðjFobsj jFcaljÞ expðijcalÞ expð 2piH:rÞ
H

At subatomic resolution, analysis of the total ED thus enables topological description of the atomic properties: atomic basins and their properties may be defined (charge, volume, dipole moment) as transferable pieces to evaluate chemical bonding and ligand or inhibitor–protein interactions [23].

In this review, we will describe topological approaches based on both medium and high-resolution ED representations and try to link them in a more general

11.2 Methodology and Technical Details 287

way by using topological results to calculate electrostatic properties of a highresolution protein model by using the transferability principle.

The chapter will be divided in four parts. The first part will focus on the technical details of the approaches used for topological analysis of high and mediumresolution ED distributions. The second part will summarize results related to the topological analysis approach based on our multipolar ED database fragment description (see below). To test both methods on a real example we will focus, in the last two sections, on results obtained by both approaches applied to the human aldose reductase (hAR) structure [6]. The analysis will focus on a selected subset of the protein active site involved in the binding of the NADPþ cofactor adenine moiety. The third part will be dedicated to results of topological analysis of ED at medium resolution and the fourth part will describe applications to the modeling of protein electrostatic properties.

11.2

Methodology and Technical Details

11.2.1

Ultra-high X-ray Resolution Approach

In contrast with the IAM model, in which all atoms are assumed to be neutral, spherical, and independent, the crystal static valence ED is modeled by a sum of multipolar pseudo atom density rvalat ðrÞ located at atomic centers [7–11], while the atomic core density rcoreðrÞ remains unchanged:

l max	X
Xl 0		ð5Þ
rvalat ðrÞ ¼ k3PvalrvalðkrÞ þ	k03Rnlðk0rÞ PlmHylmGðy; jÞ
¼	mH

where the term rval represents the spherically averaged free atom HF valence density. The second term of the summation describes the nonspherical part, in which the radial functions used are of Slater type: RnlðrÞ ¼ Nr nl expð xrÞ.

The functions ylm are spherical harmonics in real form; the ðy; jÞ coordinates are expressed in an atom-centered local axis system which facilitates application of chemical similarity. The refinable terms are the k and k0 coe cients, which describe the expansion–contraction of the perturbed valence ED [11], and the population terms Pval and Plm.

To apply this formalism to ultra high resolution protein di raction data we have proposed a multipolar data library [12, 13]: high resolution XRD data have been collected in the Laboratoire de Cristallographie et de Mode´lisation des Mate´riaux Mine´raux et Biologiques (LCM3B), Nancy, France, for a large group of mono or polypeptides to precisely determine the ED distribution of all natural amino acids (neutral or charged). Multipole refinement of the related structure factors enabled the building of an experimental database of multipolar properties for each

288 11 Topological Analysis of Proteins as Derived from Medium and High-resolution Electron Density

Fig. 11.1 Static deformation ED in the peptide-bond plane obtained with the multipolar ED database. Contour interval 0.05 e A˚ 3.

Full lines: positive values, dashed lines: negative and zero values.

protein-type atom in a given chemical environment [12, 13], from which a specific atom-type nomenclature has been developed. As an example, Fig. 11.1 gives the static electron deformation density:

Nat
DrstaticðrÞ ¼ Xj 1	½rjatmultðr rjÞ rjatIAMðr rjÞ&	ð6Þ
¼

calculated from the experimental database in the protein main chain HNaCbO peptide plane.

The valence ED distribution of covalent interactions and nonbonding electron pairs is clearly apparent. The database values were shown to be transferable to the protein amino acids and enable calculation of aspherical atomic scattering factors to be used for protein refinement. The validity of these aspherical scattering factors was checked more than 10 years ago [12] and has been confirmed since then by several studies. They have been successfully used to refine ultra- high-resolution protein structures – a scorpion toxin [14], crambin [15], and hAR [16]. As shown in these papers, the use of aspherical scattering factors improves all least-square statistical indices and consequently, leads to more physically meaningful bond distances and thermal anisotropic displacement data. Aspherical features, for example nonbonded density on CbO oxygen atoms, are also taken into account [17]. This finding, which was not surprising because it had already been reported by us and others from porphyrin [18] and naphthalene-type compound [19] charge density studies, led us to the development of the MoPro package of crystallographic programs [10]. In MoPro, for any chemical type of atom belonging to a protein type molecule the corresponding aspherical scattering factor is automatically assigned and these scattering factors can then be further used in the refinement process. Such generalized refinement, which is not

11.2 Methodology and Technical Details 289

necessarily a charge density refinement, because the aspherical scattering factors can be fixed, should be generalized to all small molecule refinements in the future. Furthermore, as shown below, this atom transferability enables estimation of electrostatic properties at least as accurate as when conventional force fields are used [16, 20].

Another approach, based on high-level theoretical calculations, was later proposed by P. Coppens’ group [21]. The objective of their theoretical database is different from ours because the also transferable multipolar data will be used to estimate properties such as protein–protein or protein–ligand interaction energies based on theoretical ED reconstruction whereas our approach is based on experimental ultra-high-resolution XRD data. Comparison of our experimental database and Coppens’ theoretical database is in progress and will be published soon. A similar project has been recently proposed by Lu¨ger and coworkers [22].

To summarize, we have at our disposal an experimental database which contains data for a set of transferable multipolar pseudo atoms that can be used for protein refinement and estimation of electrostatic properties. These transferable pseudo atoms have their own topology and may therefore be regarded as experimental topological atoms as defined by Bader [23]. In this work, topological analysis of the ED of those fragments has been performed according to the QTAIM theory. The corresponding data, i.e. CP, rb, ‘2rb and li calculated using the NewProp software [24], are discussed in Section 11.3.

11.2.2

Medium-resolution Approach

11.2.2.1 Promolecular Electron Density Distribution Calculated from Structure Factors

Several approaches are available for modeling protein ED maps at di erent levels of resolution. One consists in calculation of structure factors FðHÞ from the atomic coordinates r of the system and the IAM atomic scattering factors f ðHÞ and then in the application of a Fourier transformation (FT) algorithm to generate an ED distribution rðrÞ (Eq. 1). This procedure can be completed using any crystallographic refinement program package, for example XTAL [25].

The set of reciprocal space vectors H that is actually considered in the FT calculation is always finite and determines the crystallographic resolution d of the map.

As already mentioned, in a high-resolution map the density maxima actually correspond to atomic positions whereas in a lower-resolution map they are, instead, associated with groups of atoms. For a protein structure, at a resolution of approximately 3 A˚ , there is good correspondence between the ED peaks and the backbone or the side-chains of the amino acid residues [1–3]. Guo et al. [2], for example, showed that 3-A˚ ED distributions could be reconstructed from the socalled ‘‘globbic’’ structure factors associated with residue centers.

Following QTAIM, one can locate ð3; 3Þ CPs, i.e. ED maxima, using the program ORCRIT [26, 27]. In this sense, statistical studies that were conducted on a set of 140 highly idealized protein structures (resolution ¼ 2.85 A˚ , no hydrogen

290 11 Topological Analysis of Proteins as Derived from Medium and High-resolution Electron Density

atoms, XTAL default overall isotropic thermal variable hu2i ¼ 0:035 A˚ 2, overall atom occupancy ¼ 1) led to macromolecular models built on CP properties that could be used in protein–protein and protein–DNA docking applications [4].

11.2.2.2 Promolecular Electron Density Distribution Calculated from Atoms

Another approach used to model medium-resolution ED distributions consists in smoothing of an analytical representation of the ED. A simple way to model a molecular ED distribution function is to consider that it is the summation over individual atomic contributions ra.

The use of molecular properties such as the ED reconstructed in the framework of a promolecular representation has several essential advantages. First, as for the first method described above computation times are substantially reduced compared with, for instance, ab initio quantum-mechanical calculations. This is especially appealing when studying macromolecular systems [28].

In their works related to the promolecular atom shell approximation (PASA), Amat et al. [29] used atomic Gaussian ED functions that were fitted on 6-311G atomic basis set results (coe cients and exponents can be downloaded from http://iqc.udg.es/cat/similarity/ASA/funcset.html). A molecular or promolecular ED distribution is thus a sum over atomic Gaussian functions wherein expansion coe cients are positive to preserve the statistical meaning of the density function in the fitted structure:

raðr RaÞ ¼ Za X3 wa; i½ð2va; i=pÞ3=4e va; ijr Raj2 &2

i¼1

where wa; i and za; i are the fitted terms.

In our approach to generate smoothed 3D ED functions, pressed as the solution of the di usion equation according to sented by [30]:

ð7Þ

rM is directly exthe formalism pre-

3
ra; tðr RaÞ ¼ Xi 1	sa; i where sa; i ¼ aa; ie ba; ijr Ra j2	ð8Þ
¼

with:

aa; i ¼ Zawa; ið2va; i=pÞ3=2ð1 þ 8va; itÞ 3=2 and ba; i ¼				2va; i				ð9Þ
	ð		1	þ	8v	t	Þ
						a; i

where t is the ED smoothing degree. It has previously been demonstrated that t is actually equivalent to the overall isotropic temperature factor B, i.e., more precisely [31]:

2t ¼ B ¼ 8p2hu2i

ð10Þ

The peaks in a smoothed ED map can obviously be located with ORCRIT, in the same way as for maps generated with the software XTAL. Because the analytical expression of rðrÞ is known, however, peaks can also be located analytically,

11.2 Methodology and Technical Details 291

without calculation of an ED grid [32]. To follow the pattern of local maxima in a molecular ED distribution, as a function of the degree of smoothing, an algorithm initially described by Leung et al. [33] was implemented. The di erent steps of that merging/clustering algorithm, described elsewhere [32], consist in following the trajectory of the ED maxima, rpeakðtÞ, in a progressively smoothed ED distribution function:

D
rpeakðt þ DtÞ ¼ rpeakðtÞ þ rpeakðtÞ	‘rpeakðt þ DtÞ	ð11Þ

where ‘rpeak is the density variation and D is a predetermined constant value. The trajectory search is stopped when ‘rpeakðtÞ is lower than or equal to a limit value gradlim.

In all examples treated so far, the settings were always D ¼ 2:80 10 5 A˚ 2, the number of iterations 2000, and gradlim was set equal to 1:275 10 4 e A˚ 4. The results obtained using that algorithm are the location of the local maxima (peaks), their density, eigenvectors, and Laplacian values, and the atomic content of all fragments, at each value of t between 0 and tmax (Fig. 11.2).

11.2.3

A Test System – Human Aldose Reductase

The two topological approaches described above are applied in this section to the hAR structure, solved at a subatomic resolution of 0.66 A˚ as previously described by Howard et al. [6]. It consists of 316 amino acid residues bonded with the cofactor NADPþ and the inhibitor IDD594. There are 617 water sites and two citrate ions. The whole complex crystallizes in space group P21 with cell data a ¼ 49.43, b ¼ 66.79, c ¼ 47.40 A˚ , b ¼ 92:40 . The full multipolar refinement of hAR is being performed and will be reported elsewhere [34]. In cases of chronic hyperglycemia, aldose reductase is known to reduce part of the excess of glucose to sorbitol [35], accumulation of which in cells leads to long-term diabetes diseases, for example cataracts or nephropathies. hAR inhibition is, therefore, a well researched pharmacological target.

To test the di erent approaches of ED topological description and to apply these methods in the framework of analysis of electrostatic properties, a small structure subset was selected within the hAR model. It consists of the first amino acid layer surrounding the adenine moiety of the NADPþ cofactor (Fig. 11.3). This subset, referred to below as the ‘‘adenine binding site’’, consists of ten amino acids (160 atoms) of several types:

four nonpolar residues (LEU212; PRO215; LEU228; ALA245);

two polar residues (SER263 and ASN272);

two formally þ1 electron (e) positively charged amino acids (LYS262 and ARG268); and

two formally 1 e negatively charged amino acids (ASP216 and GLU271).