Ab Initio Fragment Molecular Orbital Study of Molecular Interactions in

J. Phys. Chem. B 2008, 112, 12081–12094

12081

Ab Initio Fragment Molecular Orbital Study of Molecular Interactions in Liganded Retinoid X Receptor: Specification of Residues Associated with Ligand Inducible Information Transmission Mika Ito,*,†,‡ Kaori Fukuzawa,§ Takeshi Ishikawa,| Yuji Mochizuki,‡,⊥ Tatsuya Nakano,‡,# and Shigenori Tanaka†,‡ Graduate School of Human DeVelopment and EnVironment, Kobe UniVersity, 3-11 Tsurukabuto, Nada-ku, Kobe 657-8501, Japan, CREST Project, Japan Science and Technology Agency, 4-1-8 Honcho, Kawaguchi, Saitama 332-0012, Japan, Mizuho Information and Research Institute, Inc., 2-3 Kanda Nishiki-cho, Chiyoda-ku, Tokyo 101-8443, Japan, CEID, Gifu UniVersity, Yanagido 1-1, Gifu 501-1194, Japan, Department of Chemistry, Faculty of Science, Rikkyo UniVersity, 3-34-1 Nishi-Ikebukuro, Toshima-ku, Tokyo 171-8501, Japan, and DiVision of Medicinal Safety Science, National Institute of Health Sciences, 1-18-1 Kamiyoga, Setagaya-ku, Tokyo 158-8501, Japan ReceiVed: April 18, 2008; ReVised Manuscript ReceiVed: July 8, 2008

The ab initio fragment molecular orbital calculations were performed for the R-subtype of the human retinoid X receptor (hRXRR) complex with its natural ligand 9-cis retinoic acid (9cRA) to quantitatively specify the key residues with important roles for the ligand inducible information transmission of RXR. In the RXR9cRA complex, the transactivation helix 12 (H12) adopts a canonical agonist conformation, which just corresponds to the transcriptional activation function 2 activating domain core (AF2C). Through the analyses of molecular interactions by the second-order Møller-Plesset perturbation (MP2) method, it was proved that Trp305 and Leu436 of the AF2C binding pocket would be important for the stabilization of the H12 canonical agonist conformation, and, at the same time, for the recognition of the 9cRA molecule. Besides, through the analyses of orbital interactions by the local MP2 (LMP2) method, it was found that Trp305 and Leu436 would recognize the 9cRA molecule especially at its C19 methyl group, which has been most notably targeted to modify for agonist and antagonist design. Moreover, on the basis of the relationships of molecular interactions, it was suggested that the interactions of Trp305 and Leu436 with AF2C residues would be significantly influenced by the interactions of Trp305 and Leu436 with 9cRA. Taken together, our findings quantitatively demonstrated that Trp305 and Leu436 would be the possible key residues for the information transmission in liganded RXR, accounting for their importance suggested by experiments. Altogether, these results substantiated that our approach is useful for the understanding of the detailed molecular mechanism underlying the transcriptional regulation of RXR and related nuclear receptors at the quantum mechanical level. 1. Introduction The retinoid X receptor (RXR) is a member of the nuclear receptor (NR) superfamily which regulates expression of many genes involved in various physiological actions of its ligands at the transcriptional level. RXR not only forms a homodimeric DNA complex, but also can form heterodimeric DNA complexes with various NRs.1 Because RXR thus has diverse important biological roles associated with human life and diseases, it has been one of the primary targets of drug discovery. As well as the functions of other NRs, the functions of RXR are induced by the binding of ligands. A natural ligand of RXR is one of carboxylic acid derivatives of vitamin A, 9-cis retinoic acid (9cRA), which controls morphogenesis, differentiation, and homeostasis during embryonal development and postnatal life. 9cRA is also an effective inhibitor of tumor cell growth, and * To whom correspondence should be addressed. E-mail: ito@ insilico.h.kobe-u.ac.jp. Tel: +81-78-803-7991. Fax: +81-78-803-7761. † Kobe University. ‡ Japan Science and Technology Agency. § Mizuho Information and Research Institute, Inc. | Gifu University. ⊥ Rikkyo University. # National Institute of Health Sciences.

this antitumor activity is useful in therapy and prevention of cancers such as human immunodeficiency virus (HIV) associated Kaposi’s sarcoma.2,3 To date, many experimental studies have been devoted to elucidate the transcriptional activation mechanism of NRs so as to efficiently exploit the functions of NRs.4,5 In an early experimental study6 by X-ray crystallography, a “mouse trap” mechanism was proposed for the transcriptional activation mechanisms of NRs. In this mechanism, conformational changes of NR ligand-binding domain (LBD) occur especially in its C-terminal helix 12 (H12) induced by the ligand binding, and NR LBD traps its ligand, just like a “mousetrap” does. H12 possesses the ligand-dependent transcriptional activation function 2 activating domain core (AF2C) which is thought to be a crucial region for the ligand-dependent transcriptional activation of NRs. On the basis of the early experimental study, H12 was proposed as a ligand-dependent “switch” of the transcriptional activation mechanism, which generates the surface for coactivator binding or recognition in its ligand-induced functional conformation.4 For the transcriptional activation of NRs, the interactions between the AF2C regions of NR LBDs and LXXLL motifs (L represents leucine; X represents any amino

10.1021/jp803369x CCC: $40.75  2008 American Chemical Society Published on Web 08/27/2008

12082 J. Phys. Chem. B, Vol. 112, No. 38, 2008 acid) of coactivators are known to be essential. It is now widely accepted that the transcriptional activity of RXR, as well as that of other NRs, is regulated by the three steps, that is, the binding of ligand, the repositioning of H12, and the exchange of the binding of transcriptional coregulators such as coactivators and corepressors.5 It is well known that the ligand-induced conformations of H12 are different by the types of ligands, including agonists, antagonists, partial agonists, and partial antagonists.4,5 Agonists induce the formation of H12 agonist conformations which allow the binding of coactivators and cause the transcriptional activation of NRs, whereas antagonists induce the formation of H12 antagonist conformations which inhibit the binding of coactivators and cause the transcriptional inactivation of NRs. Agonists and antagonists thus have the similar activities and the opposite activities, respectively, to the activities of endogenous natural ligands of NRs. Partial agonists and partial antagonists have lower activities than those of agonists and antagonists, respectively. Interestingly, a particular conformation, so-called “agonist conformation”, of H12, is essential to support the effective interactions between the AF2C region of NR LBD and the LXXLL motif of a coactivator and to subsequently induce the transcriptional activation of NR.4,5 It is considered that the information of ligands is translated into the particular conformations of H12 and recognized by the coregulators. However, it remains mysterious how the information of a ligand is transmitted to H12 without direct contacts between the ligand and the transactivation helix H12. We suppose that, between a ligand and H12, there may be some residues via which the ligand information can be transmitted to H12. To better understand the transcriptional regulation mechanism of NRs including RXR, the specification of such key residues for the ligand inducible information transmission is necessary. On the other hand, structural features of RXR LBD were well analyzed in an experimental study7 by X-ray crystallography. In the R-subtype of the human retinoid X receptor (hRXRR) complex with its natural ligand 9cRA, H12 adopts a canonical agonist conformation upon the 9cRA binding, and it just corresponds to the region of AF2C. As shown in Figure 1a, AF2C, which is composed of seven residues of H12,8 is packed into the AF2C binding pocket (AF2CBP), which is composed of ten residues of H3, H4, H5, H10, and H11, in the hRXRRcRA complex. The residues of AF2CBP are located near the residues of AF2C within 4.2 Å of the closest distances between their heavy atoms, whereas 9cRA is away from the residues of AF2C more than 5.7 Å of the distances of their closest heavy atoms. This experimental study7 suggested that H12 is held in place mainly by a set of tight hydrophobic interactions involving nonpolar residues (Leu451, Met452, Met454, and Leu455) of AF2C and residues (Cys269, Ala272, Leu276, Arg302, Trp305, Leu433, Leu436, Phe437, and Lys440) of AF2CBP, and additionally by the hydrogen bonds between the backbone amide groups of nonpolar residues (Phe450 and Leu451) of AF2C and the side chain carboxyl group of a acidic polar residue (Asp273) of AF2CBP. Structural features of ligands of RXR were also well analyzed in experimental studies.9,10 As shown in Figure 1b, the natural ligand 9cRA is an L-shaped large molecule composed of fortynine atoms including twenty-two heavy atoms, and it has two regions: the hydrophobic region and the hydrophilic region including a carboxyl group. One of the most characteristic moieties of the hydrophobic region is the C19 methyl group, at which 9cRA has a sharp turn. Because of the characteristics of the ligand, skeletal modifications of 9cRA have been notably

Ito et al.

Figure 1. (a) Ribbon display of the X-ray crystal structure (PDB code 1FBY) of the hRXRR LBD complex with 9cRA (purple). The positions of AF2C (red) in H12 and AF2CBP (yellow) are also displayed. (b) Chemical structure of 9cRA. The numbers of some carbon atoms are indicated.

performed for the hydrophobic region,9 especially for the C19 moiety,10,11 in order to synthesize agonists and antagonists. Some experimental studies10-12 show that ligands of RXR including agonists and antagonists have diverse hydrophobic regions especially at C19, though they have similar hydrophilic regions to the hydrophilic region of 9cRA including a carboxyl group. These experimental data suggest that the hydrophobic region of 9cRA, especially the C19 methyl group of 9cRA, may be more important for the recognition of the 9cRA molecule, than the hydrophilic region of 9cRA, that is, the carboxyl group. In addition, it is known that the ligand binding pocket (LBP) of RXR adopts an L-shaped pocket which can accommodate the L-shaped 9cRA, and thus, RXR LBP has preference of 9cRA over the other ligands.9 In some experimental studies7,13 by X-ray crystallography, nineteen residues surrounding 9cRA were regarded as the pertinent residues of LBP in the hRXRR-9cRA complex. These residues of LBP are located near 9cRA within 4.2 Å of the closest distances between their heavy atoms. One of the residues of LBP, Trp305, is located near the C19 methyl group of 9cRA, and the L-shaped pocket adopts a sharp turn at Trp305. An experimental study7 suggested that the antagonistic effect on hRXRR may result from the displacement of Trp305 and the subsequent disruption of the hydrophobic interactions that stabilize the agonist position of H12. The other one of the residues of LBP, Leu436, is also located near the C19 methyl

J. Phys. Chem. B, Vol. 112, No. 38, 2008 12083 group of 9cRA. An experimental study14 by point mutations reported that substitution of Val for Leu436 decreased the transcriptional activation of hRXRR. It should be noted that these two residues, Trp305 and Leu436, are not only the members of LBP but also the member of AF2CBP. Together with above experimental data, it is suggested that RXR, Trp305 and Leu436, which are in contract with H12 in its agonist conformation and whose conformations can be affected by the bound ligand, should be considered as target residues for the design of new NR modulators.10 However, it remains unclear whether Trp305 and Leu436 are the key residues for the ligand inducible information transmission. To precisely specify the key residues, detailed analyses of molecular interactions in the RXR-9cRA complex should be performed. By taking advantage of ab initio quantum mechanical (QM) calculations, we can perform detailed and accurate analyses of molecular interactions as well as molecular structures and properties including charge redistributions, which cannot be described by classical molecular mechanical (MM) methods. Moreover, by using the fragment molecular orbital (FMO) method,15-20 we can perform the ab initio QM calculations of biomacromolecules, though the conventional molecular orbital (MO) method is limited to the ab initio QM calculations of small molecules. In particular, when the fragmentations are performed according to the amino acid unit for the protein in the FMO calculations, we can estimate molecular interaction energies between residues using the interfragment interaction energies (IFIEs).21-25 Furthermore, by means of electron-correlation methods beyond the Hartree-Fock (HF) method such as the second-order Møller-Plesset perturbation (MP2) method26,27 in the FMO procedure, we can also appropriately describe the dispersion interactions, which are known to be important for interactions in biomacromolecules.24,28 Through previous theoretical studies,24,25,28,29 it is demonstrated that the FMO calculations using the MP2 method can provide valuable data, even though they were performed under gas-phase conditions. In addition to the conventional MP2 method, the local MP2 (LMP2) method30 was recently implemented in the FMO procedure of our ABINIT-MP program.31 By performing the fragment interaction analysis based on local MP2 (FILM),30 we can further describe the dispersion interactions in detail at the orbital level. In this work, in an effort to provide insight into the detailed molecular mechanism of the transcriptional regulation of RXR, we have attempted to clarify how the information of a ligand is transmitted to H12. To address this question, we have tried to quantitatively specify the key residues for the ligand inducible information transmission of RXR by the ab initio FMO method. The FMO calculations were performed for the hRXRR-9cRA complex, which have the H12 canonical agonist conformation, at the MP2 level under gas-phase conditions. In this investigation, the residues which have two important roles of the stabilization of the H12 canonical agonist conformation, and, at the same time, the recognition of the 9cRA ligand molecule were probed. Since dispersion interactions are considered to be essential for interactions between hydrophobic regions to stabilize the H12 canonical agonist conformation and to recognize the 9cRA molecule, the MP2 level of the FMO calculations is required in this investigation. To understand the entire transcriptional regulation mechanism of RXR, the residues which have roles of the formation of the H12 canonical agonist conformation and at the same time the recognition of 9cRA molecule should be probed by taking into account the ligandinduced conformational changes of H12. However, no halfway

X-ray structures of the RXR complexes between the ligand bound and unbound forms have been reported. Moreover, it is difficult to calculate numerous structures of the RXR complexes by the FMO method at present, even if the halfway X-ray structures were reported. Thus, in this work the molecular interactions in the ligand bound form of the RXR complex were analyzed. It is expected that the residues specified in this work may hold the key to understanding the transcriptional regulation mechanism of RXR, since they may affect the ligand-induced conformation of H12. To specify the key residues, we adopted one strategy that is a trial to find the key residues from the residues directly surrounding AF2C (the AF2CBP residues) instead of the residues directly surrounding 9cRA (the LBP residues). This strategy was planned based on the following reasons. Most of the residues in LBP have strong electrostatic interactions with 9cRA, because 9cRA is a polar ligand. On the other hand, the molecular recognition of the 9cRA molecule can be assumed to depend on dispersion or van der Waals interactions between the hydrophobic region of the 9cRA molecule and its surrounding residues, which are fairly weaker than electrostatic interactions between the hydrophilic region of the 9cRA molecule and its surrounding residues. Consequently, it is difficult for us to find the key residues for the molecular recognition of the ligand molecule from the residues directly surrounding 9cRA by the analyses of molecular interactions. Besides, the residues directly surrounding 9cRA may be capable of the information transmission from 9cRA to various directions, whereas the residues directly surrounding AF2C may be capable of the information transmission from 9cRA to AF2C that is the direction of interest. Therefore, we tried to find the key residues from the residues directly surrounding AF2C. 2. Theoretical Calculations The molecular interactions in a RXR-9cRA complex were examined in this study. The sample structures of the RXR-9cRA complex were prepared by the molecular dynamics (MD) structural relaxation and the molecular mechanical (MM) geometry optimization in aqueous solution. The ab initio FMO calculations were carried out for the sample structures under gas-phase conditions, because the ab initio FMO calculations including numerous water molecules in aqueous solution are difficult to perform at present. However, previous FMO calculations23-25,28,29 under gas-phase conditions have provided much valuable information on molecular interactions in biomacromolecules using X-ray crystal structures. Details of the MD and MM calculations to prepare the sample structures and the FMO calculations for the sample structures are described below. 2.1. Preparation of Molecular Structures. The initial atomic coordinate of a wild type hRXRR-9cRA complex was obtained from the Research Collaboratory for Structural Bioinformatics (RCSB) Protein Data Bank (PDB),32 PDB code 1FBY (Figure 1a).7 The entire hRXRR LBD consisting of 232 amino acid residues (residues 227-458) and a 9cRA molecule were employed for simulations. Missing hydrogen atoms, main chains, and side chains in the PDB file were complemented manually by using the molecular graphic software Molecular Operating Environment (MOE), Version 2006.08.33 The X-ray crystal structures of wild type hRXRR-9cRA complexes without missing main chains at loops have not been reported to date, though a similar X-ray crystal structure of a wild type hRXRR9cRA complex with its coactivator (SRC1; steroid receptor coactivator-1), PDB code 1FM9,34 has been reported. Therefore, missing main chains at a loop (residues 245-262) of the wild

12084 J. Phys. Chem. B, Vol. 112, No. 38, 2008 type hRXRR-9cRA complex, PDB code 1FBY, were complemented using the same region (residues 245-262) of the wild type hRXRR-9cRA-SRC1 complex, PDB code 1FM9. Hydrogen atoms were added to both N- and C-terminal residues of the peptide chains and to all of the dissociative and associative residues in their charged states. The total number of atoms in the RXR-9cRA complex is 3733 including hydrogen atoms. The charges of RXR and 9cRA in the RXR-9cRA complex are 0e and -1e, respectively, and the total charge of the RXR-9cRA complex is -1e. The AMBER99 force field35 was used to describe the RXR protein. To describe the 9cRA molecule, the atomic charges of 9cRA were obtained after geometry optimization and subsequent single-point calculation of electrostatic potential at the HartreeFock (HF) level with the 6-31G(d) basis set by using the Gaussian03 program,36 and were fitted by using the restrained electrostatic potential (RESP) procedure.37,38 These force field parameters of 9cRA were assigned on the basis of the atom types of the force field model developed by Cornell et al.39 By employing the initial atomic coordinate prepared above, the atomic coordinate of the RXR-9cRA complex was build by using the LEaP module of the Amber7 program40 on Intel Pentium 4 CPU 3.0 GHz (2 CPUs). After the coordination and energy minimization of hydrogen atoms under gas-phase conditions by using the AMBER99 force field for RXR and the force field for 9cRA prepared by the above procedure, the RXR-9cRA complex was solvated in a rectangular box of TIP3P water molecules41 with a minimum solute-wall distance of 10 Å and neutralized by adding one sodium counterion. The total number of atoms in the solvated system including the RXR9cRA complex, the water molecules, and the sodium counterion is 33482. The following MM and MD calculations were performed for the RXR-9cRA complex under the periodic boundary condition with the Sander module of the Amber8 program42 using MDGRAPE-3 system43-45 on Intel Xeon CPU 2.7 GHz (2 CPUs). Energy minimizations were performed by the steepest descent (SD) and the conjugate gradient (CG) methods. The solvated system was optimized prior to the MD simulation through two steps. At the first step, the RXR-9cRA complex was frozen, and the solvent water molecules and sodium counterion were optimized by 500 steps of the SD energy minimization followed by 500 steps of the CG energy minimization. At the second step, the entire solvated system was optimized by 1000 steps of the SD energy minimization followed by 1500 steps of the CG energy minimization. The optimized system was gradually heated from 0 to 300 K within 20 ps by a constant volume MD simulation. A harmonic restraint of 10 kcal/mol/Å2 was imposed on solutes while heating the system. Then, a production MD simulation was carried out for 10 ns under a periodic boundary condition in the NPT ensemble at constant pressure (1 atm) with isotropic position scaling and at 300 K with the Berendsen temperature coupling46 (using a time constant of 0.5 ps for heat bath coupling). No harmonic restraints were imposed during the production MD simulation. The SHAKE algorithm47,48 was applied to fix all covalent bonds containing a hydrogen atom. The particle mesh Ewald (PME) method49 was used to treat the long-range electrostatic interactions. A cutoff of 10 Å was applied to the noncovalent interactions and a time step of 2 fs was used. During the 10 ns production MD simulation, the coordinates of the simulated complex were saved every 1 ps. The root-mean-squared deviation (RMSD) measurements of the RXR-9cRA complex along the simulation time are available as Supporting Informa-

Ito et al. tion. The RMSD measurements calculated for all residues show that the structure of the RXR-9cRA complex was structurally relaxed or equilibrated after about 7 ns by the MD simulation. Additionally, the RMSD measurements calculated for all residues except for residues in loops, that is, residues in R-helices and β-sheets, show that the structures of these residues were structurally relaxed or equilibrated after about 3 ns by the MD simulation. Because of the elimination of the RMSD values for residues in flexible loops, the line of the latter RMSD measurements is flatter than that of the former RMSD measurements. The snapshots after 3 ns of the MD simulation were thus assumed to be good samples for the FMO calculations, because the residues of interest are located in R-helices. For the FMO calculations, five snapshots of equilibrated structures were taken from the last 2 ns (8.0, 8.5, 9.0, 9.5, and 10.0 ns) of the 10 ns production MD simulation. To remove irrelevant contacts between atoms and bonds in the RXR-9cRA complex, each snapshot in the MD trajectory was annealed and energy minimized prior to the FMO calculations by the following procedures. The system was gradually cooled from 300 to 0 K within 20 ps by a constant volume MD simulation. A harmonic restraint of 10 kcal/mol/Å2 was imposed on solutes while cooling the system. Then, the system was optimized by 2 ps of quenched dynamics (at a time step of 0.2 fs and at 0.1 K), followed by the CG energy minimization until the threshold energy decreased to less than 0.1 kcal/mol. 2.2. FMO Calculations. The ab initio FMO15-20 calculations were performed for the five sample structures of the RXR-9cRA complex prepared by the procedure mentioned above. The results of the FMO calculations for the five snapshots of the MD simulations at 10.0, 9.5, 9.0, 8.5, and 8.0 ns, respectively, are shown as samples 1, 2, 3, 4, and 5 in figures and tables of the following sections. All of the FMO calculations were carried out under gas-phase conditions, that is, the conditions without water, by using the ABINIT-MP program31 on AMD Opteron Processor 246 CPU 2.0 GHz (16 CPUs), and the visualizations of the results obtained by the FMO calculations were carried out using the BioStation Viewer.50 Since the FMO calculations were carried out under gas-phase conditions primarily due to computational limitations, long-range electrostatic interactions would be overestimated. However, short-range dispersion or van der Waals interactions which were especially focused on in this study would be adequately estimated by the FMO calculations even under gas-phase conditions.24,25,28,29 To save computational time without losing significant accuracy, the approximations of electrostatic potentials considered as the Mulliken orbital charge (esp-aoc) and the fractional point charge (esp-ptc) were applied to fragments whose separations of the closest contact atoms were more than 0.0 and 2.0 in units of van der Waals (vdW) radii, respectively.18 The Coulomb interaction approximation (dimeres) was also applied to fragments whose separation was more than 2.0 in vdW units.18 The fragmentation was performed according to the amino acid unit for the protein, and each amino acid residue of RXR and the 9cRA molecule were treated as a single fragment. The number of fragments in the RXR-9cRA complex is 233. The IFIEs in the RXR-9cRA complex were calculated at the residue level by the MP2 method26,27 with the 6-31G(d) basis set. The IFIEs of the N- and C-terminal residues of RXR with the other residues are not listed in the tables of the following sections, because these residues are not ends of the peptide chains in the actual system. In addition, the FMO calculations were also carried out by the LMP2 method30 with the 6-31G(d) basis set to examine the molecular interactions in detail at the orbital level, since the

J. Phys. Chem. B, Vol. 112, No. 38, 2008 12085 conventional MP2 method cannot describe the molecular interactions at the orbital level. The LMP2 method30,51,52 has been developed with the achievement of linear scaling, based on the local correlation approach which was originally proposed by Pulay.53,54 In the LMP2 method, the total correlation energy is described as the sum of the pair correlation energies based on localized molecular orbitals. By using the LMP2 method in the FMO procedure, we can obtain IFIE which is decomposed into pair correlation energies among localized orbitals. Moreover, by using the Pipek-Mezey localization scheme55 which can retain the σ-π separation, we can understand which orbital interacts with a orbital of interest. Therefore, the Pipek-Mezey localization scheme was adopted for the FILM calculations.30 Formulations of the energies obtained by the FMO calculations are described below. The total energy at the HF level (EHFtotal) in the FMO calculation is given as follows

EHFtotal )

∑ EHFIJ - (Nf - 2)∑ EHFI I>J

(1)

I

where and EHFI and EHFIJ are energies of fragment monomer I and dimer IJ, respectively. The IFIE at the HF level (∆EHFIJ) is defined as follows

∆EHFIJ ) (E′HFIJ - EHFI - E′HFJ) + Tr(∆PIJVIJ) ∆PIJ

(2)

VIJ

where is a difference density matrix, is an environmental electrostatic potential for fragment dimer IJ from other fragments, and E′HFI and E′HFIJ are energies of fragment monomer I and dimer IJ without environmental electrostatic potential. The many-body effects are considered through the environmental electrostatic potentials. By using ∆EHFIJ, EHFtotal is rewritten as follows

EHFtotal )

∑ ∆EHFIJ + ∑ EHFI I>J

(3)

I

In the FMO calculation at the MP2 level, the total energy (EMP2total) is given as follows

EMP2total )

∑ EMP2IJ - (Nf - 2)∑ EMP2I I>J

(4)

I

where and EMP2I and EMP2IJ are the MP2 energies of fragment monomer I and dimer IJ, respectively. The IFIE at the MP2 level (∆EMP2IJ) is defined as follows

∆EMP2IJ ) ∆EHFIJ + ∆EcorrIJ

(5)

where ∆EcorrIJ is the electron correlation energy by the MP2 calculation, which is obtained by the following equation

∆EcorrIJ ) EcorrIJ - EcorrI - EcorrJ Ecorr

(6)

Ecorr

In this equation, I and IJ are electron correlation energies of fragment monomer I and dimer IJ, respectively. In addition, the total energy at the LMP level (E˜LMP2total) is given as follows,

E˜LMP2total )

∑ E˜LMP2IJ - (Nf - 2)∑ E˜LMP2I I>J

(7)

I

where and E˜LMP2I and E˜LMP2IJ are the LMP2 energies of fragment monomer I and dimer IJ, respectively. The IFIE at the LMP2 level (∆E˜LMP2IJ) is defined as follows

∆E˜LMP2IJ ) ∆EHFIJ + ∆E˜corrIJ

(8)

where ∆E˜corrIJ is the electron correlation energy by the LMP2 calculation, which is obtained by the following equation

∆E˜corrIJ )

∑ εi′j′.

(9)

i′gj′

In this equation, εi′j′ are the pair correlation energies obtained from the dimer calculation by the LMP2 method. Because occupied orbitals in LMP2 are highly localized, these orbitals (i′ and j′) can be assigned to either fragment. Then, the summation of εi′j′ is limited within the span of only interfragment orbital pairs.30 On the basis of previous theoretical studies,30,51,52 it is assumed that due to the restriction of virtual spaces and the selection of correlated orbital pairs the correlation energies calculated by the LMP2 method would be estimated somewhat smaller than those calculated by the conventional MP2 method. In the previous studies,30,51,52 it was shown that the restriction of virtual spaces yields the reduction of the basis set superposition error (BSSE) from interaction energies. 3. Results and Discussion 3.1. Structural Features. At the beginning of this investigation, the structural features of the samples 1, 2, 3, 4, and 5 of the RXR-9cRA complex obtained from the snapshots at 10.0, 9.5, 9.0, 8.5, 8.0 ns, respectively, of the MD simulation were analyzed. Figure 2a shows a superposition of the five samples of RXR. As shown in Figure 2a, five samples of RXR have similar structures and retain their H12 canonical agonist conformations, which are similar to the conformation observed in the X-ray crystal structure (Figure 1a). In all of the five samples, as well as in the X-ray crystal structure, AF2C in H12 is surrounded by AF2CBP, and it is not in contact with 9cRA. 9cRA is away from the AF2C residues more than 5.9 Å of the distance of the closest heavy atoms in all of the five samples, while this distance measured for the five samples is almost the same as the distance (5.7 Å) measured for the X-ray crystal structure. On the other hand, the AF2CBP residues are close to the AF2C residues within 4.2 Å of the distance of the closest heavy atoms in all of the five samples, as well as in the X-ray crystal structure. The conformations of the AF2C residues (Phe450, Leu451, Met452, Glu453, Met454, Leu455, and Glu456) in H12 were also analyzed. Because the five samples have similar conformations of the AF2C residues, which are also similar to the conformations observed in the X-ray crystal structure, the conformation of the AF2C residues of only sample 1 of RXR9cRA complex is shown in Figure 2b. As shown in Figure 2b, side chains of two acidic polar residues (Glu453 and Glu456) of AF2C are located outside of AF2CBP. In our previous works,28,29 it has already been shown that interactions of these two acidic polar residues of AF2C significantly contribute to the coactivator binding. On the other hand, side chains of the other five nonpolar residues (Phe450, Leu451, Met452, Met454, and Leu455) of AF2C are located inside of AF2CBP as shown in Figure 2b. AF2CBP is composed of ten residues including six nonpolar residues (Ala272, Leu276, Trp305, Leu433, Leu436, and Phe437), one neutral polar residue (Cys269), one acidic polar residue (Asp273), and two basic polar residues (Arg302 and Lys440). Thus, it is supposed that interactions of AF2C nonpolar residues with AF2CBP residues may significantly contribute to the stabilization of the H12 canonical agonist conformation. In the following section, interaction energies for the stabilization of the H12 canonical agonist conformation of the RXR9cRA complex are discussed. 3.2. Interaction Energies of AF2C. Interaction energies, which are described by IFIEs, and electron correlation energies in the interaction energies, which correspond to dispersion

12086 J. Phys. Chem. B, Vol. 112, No. 38, 2008

Ito et al. TABLE 1: IFIEs (∆EMP2IJ) and Its Electron Correlation Energies (∆EcorrIJ)a of the Whole AF2C and the Whole AF2C Hydrophobic Region (AF2Cpho) with the Whole RXR or the Whole AF2CBP Calculated for the RXR-9cRA Complex sample 1

MD 10.0 ns

2

MD 9.5 ns

3

MD 9.0 ns

4

MD 8.5 ns

5

MD 8.0 ns averageb SDc

interaction

∆EMP2IJ

∆EcorrIJ

INT1 INT2 INT3 INT1 INT2 INT3 INT1 INT2 INT3 INT1 INT2 INT3 INT1 INT2 INT3 INT1 INT2 INT3 INT1 INT2 INT3

-322.18 -146.88 -98.97 -339.63 -170.13 -101.13 -247.35 -143.46 -91.31 -356.95 -165.16 -113.21 -338.92 -168.74 -121.22 -321.01 -158.87 -105.17 38.44 11.36 10.67

-75.16 -64.24 -36.35 -75.02 -65.91 -34.67 -57.85 -54.43 -28.70 -78.27 -65.11 -41.25 -76.78 -65.69 -45.03 -72.62 -63.08 -37.20 7.48 4.36 5.60

a Energies (kcal/mol) are calculated at the MP2/6-31G(d) level. INT1, 2, and 3 are AF2C-RXR (except for AF2C), AF2Cpho-RXR (except for AF2C), and AF2Cpho-AF2CBP interactions, respectively. Interaction energies of the N- and C-terminal residues of RXR are not included. b Average of samples 1-5. c Standard deviation of samples 1-5.

Figure 2. (a) Side view of the superposition of the five samples of the RXR LBD obtained from the MD simulation. The positions of AF2C (red) in H12 and AF2CBP (yellow) of the five samples and the position of 9cRA (purple) of sample 1 are also displayed. (b) Another side view of sample 1 of the RXR-9cRA (purple) complex. The positions of AF2C (red) in H12 and AF2CBP (yellow) are also displayed. All of residues in AF2C including two charged residues (Glu453 and Glu456) are indicated.

interaction energies, of whole AF2C in H12 were calculated for the five samples of the RXR-9cRA complex by the FMO method at the MP2/6-31G(d) level as shown in Table 1. The MP2 interaction energy (∆EMP2IJ) and MP2 electron correlation energy (∆EcorrIJ) are given by eq 5 and eq 6, respectively, as described above. INT1, INT2, and INT3 of Table 1 are the interaction between whole AF2C and RXR except for AF2C, that between whole AF2C hydrophobic region (AF2Cpho) and RXR except for AF2C, and that between whole AF2Cpho and whole AF2CBP, respectively. The AF2Cpho is composed of the five nonpolar residues of AF2C. Note that the calculated interaction energies (∆EMP2IJ) are perhaps too large to be considered literally. This is caused by the lack of effects of solvent and entropy and by BSSE. A previous theoretical study56

on binding affinities of various ligands for the estrogen receptor (ER) by MD simulations has been shown that, though the ligand binding energies (∆E) calculated without effects of solvent and entropy are larger than the ligand binding free energies (∆G) calculated with these effects or those derived from experiments by about 50 kcal/mol, the qualitative tendency of the former ligand binding affinities is consistent with that of the latter ligand binding affinities. Our previous work29 on the coactivator binding, which is mainly stabilized by the interactions between the acidic polar residues of AF2C and basic polar resides of a coactivator, has also shown that although the binding energies (∆E) are larger than the binding free energies (∆G) by about 1 order of magnitude, the qualitative tendency of the former binding affinities is the same as that of the latter binding affinities. In addition, previous theoretical studies25,57 on binding energies of several nucleic acid base pairs have shown that although the binding energies calculated at the MP2/6-31G level are overestimated by about 20% compared with the values including the BSSE corrections calculated at the MP2/cc-pVTZ level, the qualitative tendency of the former binding affinities is consistent with that of the latter binding affinities. Therefore, it is considered that the qualitative tendency of the interaction energies (∆EMP2IJ) would not change, even if the effects of solvent and entropy were included in this calculation and BSSE was removed from this calculation. Besides, in this work the 6-31G(d) basis set was used to reduce the BSSE. As shown in Table 1, the average magnitude of interaction energies (∆EMP2IJ) of the AF2Cpho-RXR interaction (INT2) is smaller than that of the AF2C-RXR interaction (INT1) by 162.1 kcal/mol, whereas the fluctuation of ∆EMP2IJ of INT1 is larger than that of INT2 by the standard deviation (SD) of 27.1 kcal/ mol. These differences are caused by the interactions of the two acidic polar residues (Glu453, Glu456) of AF2C. It is indicated

J. Phys. Chem. B, Vol. 112, No. 38, 2008 12087

Figure 3. Visualization of the IFIE (∆EMP2IJ) between whole AF2C hydrophobic region (AF2Cpho) and the other regions calculated for sample 1 of the RXR-9cRA complex at the MP2/6-31G(d) level. AF2Cpho is colored in yellow, and the attractive and repulsive interactions are colored in red and blue, respectively.

that even though the fluctuation of the interactions of the two acidic polar residues of AF2C are large, H12 is kept with its canonical agonist conformation in all of the five samples as shown in Figure 2a. Additionally, since side chains of these two acidic polar residues are located outside of AF2CBP as shown in Figure 2b, and these two acidic polar residues would be fully hydrated in aqueous solution, the interactions of the two acidic polar residues of AF2C would be substantially reduced by solvent effects in aqueous solution. Thus, it is suggested that the H12 canonical agonist conformation would be stabilized by the interaction of AF2Cpho with the other regions of RXR. Figure 3 shows the visualization of the interaction energy (∆EMP2IJ) between AF2Cpho and the other regions calculated for sample 1, while similar visualizations were obtained from the other four samples as well. As shown in Figure 3, AF2Cpho has significantly attractive interactions with the surrounding AF2CBP residues, which are colored in yellow in Figure 2b. As well as this visualization, ∆EMP2IJ values of the AF2CphoRXR interaction (INT2) and the AF2Cpho-AF2CBP interaction (INT3) show that AF2Cpho has significantly attractive interactions with AF2CBP residues. As shown in Table 1, the average magnitude of ∆EMP2IJ of INT3 holds about 66% of that of INT2. In INT3, the average electron correlation energy (∆EcorrIJ) holds about 35% of ∆EMP2IJ, indicating that the contribution of dispersion interactions to the stabilization of the H12 canonical agonist conformation is substantial. Besides, dispersion force is known to work effectively for short-range interactions of especially nonpolar residues, and therefore the dispersion interactions would work effectively even in aqueous solution. Together with above results, it is indicated that interaction between AF2Cpho and AF2CBP would significantly contribute to the stabilization of the H12 canonical agonist conformation. This AF2Cpho-AF2CBP interaction was analyzed in detail at the residue level, and it is discussed in the following section. 3.3. Interaction Energies of AF2C at the Residue Level. Average values of interaction energies (∆EMP2IJ) and the electron correlation energies (∆EcorrIJ) between residues of AF2Cpho and residues of AF2CBP were calculated at the MP2/6-31G(d) level as shown in Table 2. The five samples of the RXR-9cRA complex were used to obtain these average values. As mentioned above, AF2Cpho is composed of five nonpolar residues (Phe450,

Leu451, Met452, Met454, and Leu455), and AF2CBP is composed of ten residues (Cys269, Ala272, Asp273, Leu276, Arg302, Trp305, Leu433, Leu436, Phe437, and Lys440) including six nonpolar residues, one neutral polar residue, one acidic polar residue, and two basic polar residues. The numbering of these residues is not consecutive. Consequently, the residues for listed IFIEs are not located directly next to each other. The total 1 and 2 values are the total energy of interactions between all AF2Cpho residues and each AF2CBP residue and that between each AF2Cpho residue and all AF2CBP residues, respectively. The standard deviations (SDs) of ∆EMP2IJ and ∆EcorrIJ were also calculated. Some of the SDs of the residues are relatively high values of larger than 1 kcal/mol. It is indicated that the fluctuations of interactions between AF2C and AF2CBP are relatively large, because the AF2C is located outside of RXR and was strongly affected by the aqueous solution in the MD simulation. Table 2 shows that seven residues of AF2CBP have large attractive total 1 values of ∆EMP2IJ of lower than -2.0 kcal/ mol including three polar residues (Asp273, Arg302, and Lys440) and four nonpolar residues (Leu276, Trp305, Leu436, and Phe437). In the polar residues of AF2CBP, Asp273 has the largest (in magnitude) total 1 value of ∆EMP2IJ (-58.7 kcal/ mol), most of which stems from the interaction to Leu451 (-24.8 kcal/mol) of AF2Cpho. As shown in Figure 4, Asp273 is located near Leu451. The average distance between the O atom of the side chain carboxyl group of Asp273 and the H atom of the main chain amide group of Leu451 is 2.32 Å, indicating that Asp273 forms a hydrogen bond with Leu451. Thus, this strong Asp273-Leu451 interaction is due to the close conformation of Asp273 and Leu451 with the hydrogen bond. Asp273 also has large attractive ∆EMP2IJ values for the interactions with Phe450 (-12.3 kcal/mol) and Met452 (-13.0 kcal/ mol) of AF2Cpho, since Asp273 is located near these residues of AF2Cpho. In the nonpolar residues of AF2CBP, Trp305 has the largest (in magnitude) total 1 value of ∆EMP2IJ (-12.7 kcal/ mol), most of which stems from the interactions to Leu455 (-6.2 kcal/mol) and Met454 (-4.2 kcal/mol) of AF2Cpho. Additionally, it is noteworthy that, unlike the other residues of AF2Cpho, Leu455 of AF2Cpho is more stabilized by the interaction of the AF2CBP nonpolar residues than by that of the AF2CBP polar residues with the calculated ∆EMP2IJ value of -12.2 kcal/mol. Because Leu455 is located at the end of the AF2Cpho, which is very close to the C-terminal residue of RXR, the interaction of Leu455 with the nonpolar residues of AF2CBP may have a great influence on the direction of H12 canonical agonist conformation. Table 2 also shows that six residues of AF2CBP have large attractive total 1 values of ∆EcorrIJ of lower than -2.0 kcal/mol including two polar residues (Asp273 and Arg302) and four nonpolar residues (Leu276, Trp305, Leu436, and Phe437), and they are included in the above seven residues which have large attractive total 1 values of ∆EMP2IJ. On the other hand, the seventh residue (Lys440) has a smaller total 1 value (-1.8 kcal/ mol) of ∆EcorrIJ than the other six residues, because this residue is a basic polar residue and is located relatively away from the residues of AF2Cpho as compared to the other six residues. In the polar residues of AF2CBP, Asp273 has the largest (in magnitude) total 1 value of ∆EcorrIJ (-7.3 kcal/mol), most of which stems from the interaction to Leu451 (-5.4 kcal/mol) of AF2Cpho. This strong dispersion interaction is due to the close conformation of Asp273 and Leu451, as mentioned above. In the nonpolar residues of AF2CBP, Leu276 and Trp305 have the largest (in magnitude) total 1 value of ∆EcorrIJ (-6.7 kcal/

12088 J. Phys. Chem. B, Vol. 112, No. 38, 2008

Ito et al.

TABLE 2: Averagea IFIEs (∆EMP2IJ) and Its Electron Correlation Energies (∆EcorrIJ)b between the Residues of the AF2C Hydrophobic Region (AF2Cpho) and the Residues of AF2CBP Calculated for the RXR-9cRA Complex AF2CBP

AF2Cpho

position

polarityc

residue

Phe450

Leu451

Met452

Met454

Leu455

total 1d

SDe

∆EMP2

H3 H3 H3 H3 H4 H5 H10 H11 H11 H11

neutral non acidic non basic non non non non basic

∆EcorrIJ

H3 H3 H3 H3 H4 H5 H10 H11 H11 H11

neutral non acidic non basic non non non non basic

Cys269 Ala272 Asp273 Leu276 Arg302 Trp305 Leu433 Leu436 Phe437 Lys440 polarf nonpolarg total 2h Cys269 Ala272 Asp273 Leu276 Arg302 Trp305 Leu433 Leu436 Phe437 Lys440 polarf nonpolarg total 2h

0.26 0.03 -12.27 -1.46 -1.36 -0.22 -0.01 -0.11 0.08 -3.63 -16.99 -1.69 -18.68 0.00 0.00 -1.76 -2.18 -0.02 0.00 0.00 0.00 0.00 0.00 -1.78 -2.18 -3.96

-0.28 -0.49 -24.76 -1.93 -2.90 -1.12 -0.03 -0.51 -0.18 -2.06 -30.01 -4.26 -34.27 -1.33 -1.43 -5.40 -2.48 0.00 -1.62 0.00 -0.42 0.00 0.00 -6.73 -5.96 -12.69

-0.17 0.06 -13.03 1.49 -1.20 -0.96 0.23 -1.45 0.83 -3.73 -18.13 0.21 -17.92 -0.64 0.00 -0.13 -0.19 0.00 -0.20 0.00 -1.00 -0.11 -0.82 -1.58 -1.51 -3.09

-0.13 0.13 -5.07 -0.76 -10.41 -4.19 -0.25 0.02 -0.45 1.99 -13.62 -5.50 -19.11 0.00 0.00 0.00 -1.82 -2.81 -3.35 0.00 0.00 0.00 -0.03 -2.84 -5.17 -8.01

-0.04 0.04 -3.61 0.34 -1.89 -6.17 -0.99 -1.77 -3.66 2.57 -2.97 -12.22 -15.19 -0.01 0.00 0.00 -0.04 -0.14 -1.54 -1.23 -2.16 -3.38 -0.95 -1.10 -8.35 -9.45

-0.36 -0.23 -58.74 -2.33 -17.76 -12.65 -1.04 -3.83 -3.38 -4.86 -81.71 -23.45 -105.17 -1.98 -1.43 -7.29 -6.71 -2.97 -6.71 -1.23 -3.58 -3.49 -1.79 -14.03 -23.17 -37.20

0.24 0.16 8.16 0.47 3.07 2.90 0.15 1.07 0.99 2.27 7.29 3.91 10.67 0.44 0.14 1.65 0.63 0.84 2.55 0.30 1.02 1.10 1.31 2.47 3.41 5.60

IJ

a Average of samples 1, 2, 3, 4, and 5 obtained by the snapshot at 10.0, 9.5, 9.0, 8.5, and 8.0 ns of the MD simulation, respectively. Energies (kcal/mol) are calculated at the MP2/6-31G(d) level. c Nonpolar residues are indicated by “non”, and neutral, acidic, and basic polar residues are indicated by “neutral”, “acidic”, and “basic”, respectively. d Total energy of the interactions between all AF2Cpho residues and each AF2CBP residue. e Standard deviation of the total 1 values calculated for samples 1-5. f Total energy of the interactions between each AF2Cpho residue and polar AF2CBP residues. g Total energy for the interactions between each AF2Cpho residue and nonpolar AF2CBP residues. h Total energy for the interactions between each AF2Cpho residue and all AF2CBP residues. b

Figure 4. Positions of Asp273 and Leu451 of sample 1 of the RXR9cRA complex. The positions of AF2C (red) in H12 and AF2CBP (yellow) are also displayed. The average distances between Asp273 and Leu451 calculated for the five samples are indicated.

mol). Most of ∆EcorrIJ of Leu276 stems from the interactions to Leu451 (-2.5 kcal/mol) and Phe450 (-2.2 kcal/mol) of AF2Cpho. Most of ∆EcorrIJ of Trp305 stems from the interaction to Met454 (-3.4 kcal/mol) of AF2Cpho. In addition, Phe437 and Leu436 of AF2CBP have relatively large ∆EcorrIJ values (-3.4 and -2.2 kcal/mol, respectively) for the interactions with Leu455 of AF2Cpho. Note that, generally, dispersion energy is known to be sensitive to detailed conformations of molecules. Then, dispersion energy could have a great influence on detailed conformations of molecules. In addition, dispersion interactions between the AF2Cpho and AF2CBP residues are essential for hydrophobic interactions between the AF2Cpho and AF2CBP

residues, which are suggested to be a main force to hold the proper position of H12 in aqueous solution.7 Thus, it is suggested that these AF2CBP residues (Asp273, Leu276, Arg302, Trp305, Leu436, and Phe437) would significantly contribute to stabilization of the H12 canonical agonist conformation not only by their strong electrostatic interactions but also by their dispersion interactions with the AF2Cpho residues. From the above results about the AF2Cpho-AF2CBP interaction, it is suggested that the above seven AF2CBP residues (Asp273, Leu276, Arg302, Trp305, Leu436, Phe437, and Lys440) which have stronger interactions with the AF2Cpho residues would be essential for the stabilization of the H12 canonical agonist conformation. In the following section, the 9cRA-AF2CBP interaction is also discussed in detail. 3.4. Interaction Energies of 9cRA at the Residue Level. Average values of ∆EMP2IJ and ∆EcorrIJ between 9cRA and the residues of AF2CBP were calculated at the MP2/6-31G(d) level as shown in Table 3. The five samples of the RXR-9cRA complex were used to obtain these average values. The SDs of ∆EMP2IJ and ∆EcorrIJ were also calculated. All of the SDs of the residues are relatively low values of smaller than 1 kcal/mol. It is indicated that, unlike the fluctuations of interactions between AF2C and AF2CBP, the fluctuations of interactions between 9cRA and AF2CBP are relatively small, because the 9cRA is located inside of RXR and was not strongly affected by the aqueous solution in the MD simulation. 9cRA is an acidic molecule with one carboxyl group, and it is a large molecule composed of 49 atoms including 22 heavy atoms (Figure 1b). It is known that skeletal modifications of 9cRA have been notably performed for the hydrophobic region,9 especially for

J. Phys. Chem. B, Vol. 112, No. 38, 2008 12089 TABLE 3: Averagea IFIEs (∆EMP2IJ) and Its Electron Correlation Energies (∆EcorrIJ)b between the 9cRA Ligand and the Residues of AF2CBP Calculated for the RXR-9cRA Complex AF2CBP position H3 H3 H3 H3 H4 H5 H10 H11 H11 H11

polarityc neutral non acidic non basic non non non non basic

residue

∆EMP2IJ

SD 1d

∆EcorrIJ

SD 2e

Cys269 Ala272 Asp273 Leu276 Arg302 Trp305 Leu433 Leu436 Phe437 Lys440 polarf nonpolarg totalh

2.77 2.47 23.56 -3.79 -13.42 1.56 -2.76 -1.93 -0.07 -14.32 -1.42 -4.52 -5.94

0.24 0.63 0.66 0.22 0.17 0.67 0.97 0.42 0.07 0.48 0.68 1.42 2.09

-0.72 -2.70 -0.34 0.00 0.00 -0.77 -2.14 -1.68 -0.14 0.00 -1.07 -7.42 -8.49

0.11 0.16 0.04 0.00 0.00 0.30 0.37 0.45 0.02 0.00 0.14 0.81 0.93

a Average of samples 1, 2, 3, 4, and 5 obtained by the snapshot at 10.0, 9.5, 9.0, 8.5, and 8.0 ns of the MD simulation, respectively. Energies (kcal/mol) are calculated at the MP2/6-31G(d) level. c Nonpolar residues are indicated by “non”, and neutral, acidic, and basic polar residues are indicated by “neutral”, “acidic”, and “basic”, respectively. d Standard deviation of the IFIEs calculated for samples 1-5. e Standard deviation of the correlation energies calculated for samples 1-5. f Total energy for the interactions between 9cRA and polar AF2CBP residues. g Total energy for the interactions between 9cRA and nonpolar AF2CBP residues. h Total energy for the interactions between 9cRA and all AF2CBP residues. b

the C19 moiety,10,11 in order to synthesize agonists and antagonists. Therefore, the hydrophobic region of 9cRA, especially the C19 methyl group of 9cRA, can be regarded as the more important region for the recognition of the 9cRA molecule than the hydrophilic region of 9cRA, that is, the carboxyl group. Thus, analyses of dispersion interactions, which can be evaluated by ∆EcorrIJ, between 9cRA and the residues are necessary to specify the important residues for the recognition of the 9cRA molecule. As shown in Table 3, the acidic polar residue, Asp273, has large repulsive ∆EMP2IJ value (23.6 kcal/mol), and the basic polar residues, Arg302 and Lys440, have large attractive ∆EMP2IJ value (-13.4 and -14.3 kcal/mol, respectively). Because 9cRA is an acidic molecule, acidic and basic polar residues have large repulsive and attractive ∆EMP2IJ values, respectively, for their interactions with 9cRA, even though these residues are located apart from 9cRA. The ∆EcorrIJ values of almost zero calculated for Asp273 (-0.3 kcal/mol), Arg302 (-0.0 kcal/mol), and Lys440 (-0.0 kcal/mol) show that these acidic and basic polar residues are located apart from 9cRA, and that it is difficult for them to recognize details of the molecular structure of 9cRA. On the other hand, the neutral polar residue, Cys269, and the four nonpolar residues, Ala272, Trp305, Leu433, and Leu436, have larger attractive ∆EcorrIJ values than the other residues. The ∆EcorrIJ values of lower than -0.7 kcal/mol calculated for Cys269 (-0.7 kcal/mol), Ala272 (-2.7 kcal/mol), Trp305 (-0.8 kcal/mol), Leu433 (-2.1 kcal/mol), and Leu436 (-1.7 kcal/ mol) show that these neutral polar and nonpolar residues are located near 9cRA, and that they can recognize details of the molecular structure of 9cRA. Thus, it is suggested that the conformations of these five AF2CBP residues (Cys269, Ala272, Trp305, Leu433, and Leu436) could be affected by the conformation of 9cRA. Note that, as well as the other criteria used in this study, the limit (-0.7 kcal/mol) of correlation energy was used to only select the residues which show relatively strong interactions, though this limit may be small. As shown in the previous section, the seven AF2CBP residues (Asp273, Leu276, Arg302, Trp305, Leu436, Phe437, and Lys440) would be essential for the stabilization of the H12 canonical agonist conformation, and, as shown in this section, the above five AF2CBP residues (Cys269, Ala272, Trp305,

Leu433, and Leu436) would be essential for the recognition of the 9cRA ligand molecule. Only Trp305 and Leu436 are included in both the former seven AF2CBP resides and the latter five AF2CBP residues. In other words, the other three residues of the latter five AF2CBP residues (Cys269, Ala272, and Leu433) would be essential for the recognition of the 9cRA ligand molecule, but they would not be essential for the stabilization of the H12 canonical agonist conformation. Thus, through the analyses of interactions of AF2C and 9cRA with AF2CBP, it is suggested that Trp305 and Leu436 of AF2CBP are important residues to stabilize the H12 canonical agonist conformation and, at the same time, to recognize the 9cRA ligand molecules. However, it remains unclear which moieties of 9cRA are recognized by Trp305 and Leu436. This question is solved in the following section. 3.5. Interaction Energies of 9cRA at the Orbital Level. To understand which moieties of 9cRA are recognized by Trp305 and Leu436, the Trp305-9cRA and Leu436-9cRA interactions were further analyzed in detail at the orbital level. Table 4 shows the average interaction energies and correlation energies of the Trp305-9cRA and Leu436-9cRA interactions calculated at the LMP2/6-31G(d) level. The five samples are used to obtain this average values. The LMP2 interaction energy (∆E˜LMP2IJ) and the LMP2 electron correlation energy (∆E˜corrIJ) are given by eq 8 and eq 9, respectively, as described above. The ∆E˜corrIJ values were calculated for the interactions of all 9cRA moieties with Trp305 and Leu436, while the ∆E˜LMP2IJ values were calculated only for the interactions of the whole 9cRA with Trp305 and Leu436. Only selected moieties with larger electron correlation energies than the other moieties are listed in Table 4. For the interactions of the whole 9cRA with Trp305 and Leu436, ∆E˜corrIJ values of -0.6 and -1.4 kcal/mol, respectively, were calculated at the LMP2/6-31G(d) level as shown in Table 4, which are slightly smaller than their ∆EcorrIJ values (-0.8 and -1.7 kcal/mol, respectively) calculated at the MP2/6-31G(d) level. Consequently, the ∆ E˜LMP2IJ values of 1.7 and -1.7 kcal/ mol were calculated at the LMP2/6-31G(d) level, while the repulsive Trp305-9cRA interaction energy and the attractive Leu436-9cRA interaction energy are slightly larger and smaller than their ∆EMP2IJ values (1.6 and -1.9 kcal/mol, respectively)

12090 J. Phys. Chem. B, Vol. 112, No. 38, 2008

Ito et al.

TABLE 4: Averagea IFIEs (∆E˜LMP2IJ) and Its Electron Correlation Energies (∆E˜corrIJ)b of the Interactions between Trp305 or Leu436 and 9cRA Ligand Moieties Calculated for the RXR-9cRA Complex residue

moiety of ligandc

Trp305

whole 9cRA C19 methyl group C9-C10 double bond C10 methine group

Leu436

whole 9cRA C19 methyl group C17 methyl group C16 methyl group

∆E˜LMP2IJ 1.72

-1.65

SD 1d

∆E˜corrIJ

SD 2e

ratiof

0.63

-0.62 -0.27 -0.10 -0.05

0.26 0.15 0.03 0.02

0.43 0.16 0.09

-1.40 -0.35 -0.28 -0.15

0.36 0.15 0.11 0.07

0.25 0.20 0.11

0.31

a Average of samples 1, 2, 3, 4, and 5 obtained by the snapshot at 10.0, 9.5, 9.0, 8.5, 8.0 ns of the MD simulation, respectively. b Energies (kcal/mol) are calculated at the LMP2/6-31G(d) level. c Only selected moieties with larger electron correlation energies than the other moieties are shown. d Standard deviation of the IFIEs calculated for samples 1-5. e Standard deviation of the correlation energies calculated for samples 1-5. f Ratio of the average correlation energy of each ligand moiety to that of the whole ligand molecule.

Figure 5. Visualizations of the orbital pairs for the interactions between the moieties of 9cRA and (a1-3) Trp305 or (b1-3) Leu436 of AF2CBP calculated for sample 1 of the RXR-9cRA complex at the LMP2/6-31G(d) level. Only selected orbital pairs with the largest electron correlation energies for the interactions of the moieties of 9cRA are shown. The positive and negative orbital phases of each moiety of 9cRA are colored in red and blue, respectively, and those of each moiety of Trp305 or Leu436 are colored in green and yellow, respectively.

calculated at the MP2/6-31G(d) level. These values show that similar tendency of the interaction energies and the electron correlation energies can be obtained from the MP2 and LMP2 methods, and thus the validity of the LMP2 calculations of these energies for the Trp305-9cRA and Leu436-9cRA interactions was confirmed. Table 4 shows that, in the whole Trp305-9cRA interaction, the listed three hydrophobic moieties of 9cRA (the C19 methyl group, the C9-C10 double bond, and the C10 methine group) have stronger interactions with Trp305 than the other moieties of 9cRA. Table 4 also shows that, in the whole Leu436-9cRA interaction, the listed three hydrophobic moieties of 9cRA (the C19, C17, C16 methyl groups) have stronger interactions with Leu436 than the other moieties of 9cRA. Figure 5 shows the visualizations of the orbital pairs for the interactions between these moieties of 9cRA and Trp305 or Leu436 calculated for sample 1 of the RXR-9cRA complex, while similar visualizations were obtained from the five samples. Only selected orbital pairs with the largest εi′j′ values for the interactions of these moieties of 9cRA are shown. In the Trp305-9cRA interaction, the σ-orbital of the C19 methyl group of 9cRA interacts with the π-orbital of the indole ring of Trp305 (Figure 5a1), the

π-orbital of the C9-C10 double bond of 9cRA interacts with the σ-orbital of the indole ring of Trp305 (Figure 5a2), and the σ-orbital of the C10 methine group of 9cRA interacts with the σ-orbital of the indole ring of Trp305 (Figure 5a3). In the Leu436-9cRA interaction, the σ-orbital of the C19 methyl group of 9cRA interacts with the σ-orbital of the methyl group of Leu436 (Figure 5b1), the σ-orbital of the C17 methyl group of 9cRA interacts with the lone pair orbital of the main chain N atom of Leu436 (Figure 5b2), and the σ-orbital of the C16 methyl group of 9cRA interacts with the σ-orbital of the methylene group of Leu436 (Figure 5b3). These visualizations show that the Trp305-9cRA and Leu436-9cRA dispersion interactions are shared primarily by the orbital pairs formed by the hydrophobic moieties of 9cRA and the side chain hydrophobic moieties of Trp305 and Leu436, respectively. As shown in Table 4, Trp305 interacts with the C19 methyl group, the C9-C10 double bond, and the C10 methine group of 9cRA holding 43, 16, and 9% of the whole Trp305-9cRA interaction, respectively. Leu436 interacts with the C19, C17, and C16 methyl groups of 9cRA holding 25, 20, and 11% of the whole Leu436-9cRA interaction, respectively. In all of the five samples, Trp305 and Leu436 have larger correlation

J. Phys. Chem. B, Vol. 112, No. 38, 2008 12091

Figure 6. Positions of selected AF2CBP residues, (a) Leu276, (b) Trp305 and Leu436, of sample 1 of the RXR-9cRA complex. The positions of 9cRA and the AF2C (red) residues which are the closest to the selected AF2CBP (yellow) residues are also displayed. The average distances between closest heavy atoms of residues and 9cRA molecule calculated for the five samples are indicated.

energies for the interactions with these hydrophobic moieties of 9cRA than the other moieties of 9cRA, though the correlation energies of these interactions are very small. Thus, it is indicated that the 9cRA molecule is recognized by Trp305 and Leu436 at these hydrophobic moieties, especially at the C19 methyl group, of 9cRA. It is noteworthy that the C19 methyl group of 9cRA, which has been most notably targeted to modify for agonist and antagonist design,10,11 more strongly interacts with both of Trp305 and Leu436 than the other groups of 9cRA. 3.6. Relationships of Interactions of AF2C and 9cRA. From the above analyses, Trp305 and Leu436 of AF2CBP were suggested to be the possible key residues of the ligand inducible information transmission of RXR with the dual roles of stabilizing the H12 canonical agonist conformation and at the same time recognizing the 9cRA ligand molecule particularly at the C19 methyl group. In order to confirm whether the Trp305–AF2C and Leu436–AF2C interactions are influenced by the Trp305–9cRA and Leu436–9cRA interactions, respectively, the relationships between the Trp305–9cRA and Trp305–AF2C interactions and between the Leu436–9cRA and Leu436–AF2C interactions were analyzed using the interaction energies (∆EMP2IJ) calculated at the MP2/6-31G(d) level. Leu276 of AF2CBP was also analyzed in contrast to Trp305 and Leu436. For this analysis, Leu451, Met454, and Leu455 of AF2C, which have stronger interactions with Leu276, Trp305, and Leu436 of AF2CBP (see Table 2), were selected. The positional features of the selected residues of AF2CBP, Leu276, Leu436, and Trp305, are described below. Figure 6 shows the positions of the selected residues of sample 1 of the RXR-9cRA complex, while similar positions were observed in the five samples. The average distances between closest heavy atoms of residues and 9cRA were calculated using the five samples, and are also shown in Figure 6. As shown in Figure 6a, Leu276 is located near Leu451 of AF2C with the calculated average distance of closest heavy atoms of 3.7 Å, but is not located near 9cRA. The closest heavy atoms of Leu276 and 9cRA are the main chain N atom and the carboxyl O atom,

respectively, and the distance between these two atoms is 7.9 Å. Additionally, Leu276 is away from the C19 atom of the C19 methyl group of 9cRA, which is regarded as the one of the most important moieties of 9cRA for the 9cRA molecular recognition,10,11 with the longer average distance (9.9 Å between the C2 atom of Leu276 and the C19 atom of 9cRA). On the other hand, as shown in Figure 6b, Leu436 is located near Leu455 of AF2C and 9cRA with the calculated average distances of 4.1 and 4.6 Å, respectively. The C19 atom of 9cRA is in the closest position to side chain C atom of Leu436. Trp305 is located near Met454 of AF2C, Leu455 of AF2C, and 9cRA with the calculated average distances of 4.0, 4.6, and 4.8 Å, respectively. The C19 atom of 9cRA is in the closest position to side chain C1 atom of Trp305. Figure 6b shows that the possible key residues, Trp305 and Leu436, are located in the top and bottom sides of AF2C, respectively. The pair of Trp305 and Leu436 appears to interact with the AF2C residues by sandwiching AF2C. The energetic features of the selected residues of AF2CBP are described below. As shown in Table 2, Leu436 interacts with Leu455 of AF2C by the average interaction energy (∆EMP2IJ) of -1.8 kcal/mol which is about the same as the average Leu276-Leu451 interaction energy of -1.9 kcal/ mol, and both of the Leu436-Leu455 and Leu276-Leu451 interactions have the large correlation energies (∆EcorrIJ) of -2.2 and -2.5 kcal/mol, respectively. In addition, Leu436 interacts with 9cRA by the average interaction energy (∆EMP2IJ) of -1.9 kcal/mol which is similar to the average Leu276-9cRA interaction energy of -3.8 kcal/mol. The Leu436-9cRA interaction has the large correlation energy (∆EcorrIJ) of -1.7 kcal/mol, whereas the Leu276-9cRA interaction has virtually no correlation energy (∆EcorrIJ). This is because Leu436 is located near 9cRA, whereas Leu276 is away from 9cRA as shown in Figure 6. On the other hand, Trp305 interacts with Met454 and Leu455 of AF2C by the average interaction energy (∆EMP2IJ) of -4.2 and -6.2 kcal/ mol, respectively, and the Trp305-Met454 and Trp305-Leu455 interactions have the large correlation energies (∆EcorrIJ) of -3.4 and -1.5 kcal/mol, respectively. In addition, Trp305

12092 J. Phys. Chem. B, Vol. 112, No. 38, 2008

Ito et al.

Figure 7. Relationships between AF2CBP-9cRA (x-axis) and AF2CBP-AF2C (y-axis) interaction energies calculated at the MP2/6-31G(d) level. The five samples of the RXR-9cRA complex are used for the plot. The regression lines and the correlation coefficients (R) based on the least-squares method are also shown. (a) Leu276-9cRA (x-axis) and Leu276-Leu451 (y-axis) interaction energies are plotted. (b) Leu436-9cRA (x-axis) and Leu436-Leu455 (y-axis) interaction energies are plotted. (c) Trp305-9cRA (x-axis) and Trp305-Met454 (y-axis) interaction energies are plotted. (d) Trp305-9cRA (x-axis) and Trp305-Leu455 (y-axis) interaction energies are plotted.

interacts with 9cRA by the average interaction energy (∆EMP2IJ) of 1.6 kcal/mol with the correlation energy (∆EcorrIJ) of -0.8 kcal/mol. It should be noted that the Leu436-9cRA and Trp305-9cRA interactions have opposite interaction energies, that is, the Leu436-9cRA and Trp305-9cRA interaction energies are attractive and repulsive, respectively. The relationship of the interaction energies (∆EMP2IJ) of the selected residues of AF2CBP (Leu276, Leu436, and Trp305), which have positional and energetic features above, were analyzed. The ∆EMP2IJ values of the Leu276-Leu451, Leu436Leu455, Trp305-Met454, and Trp305-Leu455 interactions obtained from the five samples are plotted against the ∆EMP2IJ values of the Leu276-9cRA, Leu436-9cRA, Trp305-9cRA, and Trp305-9cRA interactions, respectively, in Figure 7. The regression lines and the correlation coefficients (R) based on the least-squares method are also shown in Figure 7. As shown in Figure 7a, the ∆EMP2IJ value of the Leu276-Leu451 interaction is hardly related to the ∆EMP2IJ value of Leu276-9cRA interaction; the correlation coefficient (R) is 0.12. On the other hand, as shown in Figure 7b, the ∆EMP2IJ value of the Leu436-Leu455 interaction is closely related to the ∆EMP2IJ value of Leu436-9cRA interaction; the correlation coefficient (R) is 0.93. Thus, it is indicated that the Leu436-Leu455 attractive interaction becomes stronger as the Leu436-9cRA attractive interaction becomes stronger. In addition, as shown in Figure 7c,d, the ∆EMP2IJ value of the Trp305-Met454 and Trp305-Leu455 interactions are also related to the ∆EMP2IJ value of Trp305-9cRA interaction; the correlation coefficients (R) are -0.74 and -0.62, respectively. Thus, it is also indicated that the Trp305-Met454 and Trp305-Leu455 attractive interaction becomes stronger as the Trp305-9cRA repulsive interaction becomes stronger. It is noteworthy that the absolute values of the correlation coefficients calculated for the interactions

of Trp305 (-0.74 and -0.62 for the Trp305-Met454 and Trp305-Leu455 interactions, respectively) and Leu436 (0.93 for the Leu436-Leu455 interaction) are larger than 0.5. This result suggests that the Trp305-AF2C and Leu436-AF2C interactions would be influenced by the Trp305-9cRA and Leu436-9cRA interactions, respectively. From the above results, it is suggested that the stability of the H12 canonical agonist conformation of RXR may be affected by the interactions of Trp305 and Leu436 with 9cRA. Therefore, it has been confirmed that the Trp305 and Leu436 would be the possible key residues for the ligand inducible information transmission of RXR. 4. Conclusions The ab initio FMO calculations were performed for the hRXRR-9cRA complex with the canonical agonist conformation of H12, which corresponds to AF2C, to quantitatively specify the key residues with important roles for the ligand inducible information transmission of RXR. First, the interaction between AF2C and AF2CBP was analyzed. It was indicated that the H12 canonical agonist conformation is stabilized mainly by the interaction between the AF2C hydrophobic region and AF2CBP. It was shown that seven (Asp273, Leu276, Arg302, Trp305, Leu436, Phe437, and Lys440) of AF2CBP residues have strongly attractive interactions (of lower than -2.0 kcal/mol) with the AF2C hydrophobic region and significantly contribute to stabilize the H12 canonical agonist conformation. Second, the interaction between 9cRA and AF2CBP was analyzed. It was shown that five (Cys269, Ala272, Trp305, Leu433, and Leu436) of AF2CBP residues have stronger dispersion interactions (of lower than -0.7 kcal/mol) with 9cRA,

J. Phys. Chem. B, Vol. 112, No. 38, 2008 12093 indicating that these AF2CBP residues are in contact with 9cRA and their conformations can be affected by the conformation of 9cRA. From the residue level analyses of interactions of AF2C and 9cRA with AF2CBP, it was proved that Trp305 and Leu436 of AF2CBP would be important residues to stabilize the H12 canonical agonist conformation and, at the same time, to recognize the 9cRA ligand molecules. By the orbital level analysis of the interactions of Trp305 and Leu436 with 9cRA using the LMP2 method, it was shown that Trp305 and Leu436 have stronger dispersion interactions with the C19 methyl group of 9cRA than with the other functional or substituted groups of 9cRA. Thus, it was found that Trp305 and Leu436 would recognize the 9cRA ligand molecules especially at its C19 methyl group, which has been most notably targeted to modify for agonist and antagonist design.10,11 In addition, the relationships between the Trp305-9cRA and Trp305-AF2C interactions and between the Leu436-9cRA and Leu436-AF2C interactions were analyzed. In this analysis, AF2C residues (Met454 and Leu455) which have strong interactions with Trp305 or Leu436 were selected. It was shown that the Trp305-Met454 and Trp305-Leu455 interaction energies are closely related to the Trp305-9cRA interaction energy with large absolute values of correlation coefficients (R ) -0.74 and -0.62) and that the Leu436-Leu455 interaction energy is closely related to the Leu436-9cRA interaction energy with large absolute value of correlation coefficient (R ) 0.93). Therefore, it was suggested that the interactions of Trp305 and Leu436 with AF2C residues (Met454 and Leu455) would be influenced by the interactions of Trp305 and Leu436 with 9cRA. Taken together, our findings quantitatively demonstrated that Trp305 and Leu436 of AF2CBP would be the possible key residues for the ligand inducible information transmission of RXR with the dual roles of the stabilization of the H12 canonical agonist conformation and, at the same time, the recognition of the 9cRA ligand molecule particularly at the C19 methyl group, providing a comprehensive picture consistent with experimental suggestions.7,10,14 Hence, it was suggested that the 9cRA ligand information might be transmitted to AF2C in the H12 canonical agonist conformation via Trp305 and Leu436 by their molecular interactions. In the present paper, the structures for protein systems were generated using classical MD simulations, while their electronic properties were analyzed on the basis of ab initio FMO method. It may thus be remarked that more consistent treatments58-60 concerning the molecular interactions would be desirable in future investigations. Besides, to obtain fully converged results on structural and energetic aspects, the FMO-MP2 calculations should be performed for a lot of samples, though the FMOMP2 calculations were performed for the only five samples in this study due to computational constraints. Additionally, to further understand the mechanism of the ligand inducible information transmission of RXR, molecular interactions with the progression of the ligand-induced conformational changes of H12 in the RXR complex should be analyzed by performing numerous FMO-MP2 calculations for various structures of the RXR complexes under aqueous-solution conditions, while such calculations could be performed on supercomputers like the Earth Simulator in the future.61 However, the present QM calculations based on the FMO method could provide us with valuable information about the molecular interactions in liganded RXR with its H12 canonical agonist conformation at the residue and orbital levels. It is expected that if our approach is performed for other NRs, the key residues for the information transmission

in liganded NRs can also be specified. Altogether, the knowledge obtained from this work would be helpful for our better understanding of the transcriptional regulation mechanism of RXR and related NRs at the QM level. Acknowledgment. We thank Dr. Hirofumi Watanabe for technical assistance. This work was supported by the “Core Research for Evolutional Science and Technology” (CREST) project of the Japan Science and Technology Agency (JST). Supporting Information Available: Root-mean-squared deviation (RMSD) measurements of the MD simulation. This information is available free of charge via the Internet at http:// pubs.acs.org. References and Notes (1) Mangelsdorf, D. J.; Thummel, C.; Beato, M.; Herrlich, P.; Schu¨tz, G.; Umesono, K.; Blumberg, B.; Kastner, P.; Mark, M.; Chambon, P.; Evans, R. M. Cell 1995, 83, 835–839. (2) Ross, S. A.; McCaffery, P. J.; Drager, U. C.; Luca, L. M. D. Physiol. ReV. 2000, 80, 1021–1054. (3) Miles, S. A.; Dezube, B. J.; Lee, J. Y.; Krown, S. E.; Fletcher, M. A.; Saville, M. W.; Kaplan, L.; Groopman, J.; Scadden, D. T.; Cooley, T.; Von Roenn, J.; Friedman-Kien, A. AIDS 2002, 16, 421–429. (4) Moras, D.; Gronemeyer, H. Curr. Opin. Cell. Biol. 1998, 10, 384– 391. (5) Nagy, L.; Schwabe, J. W. Trends Biochem. Sci. 2004, 29, 317– 324. (6) Renaud, J. P.; Rochel, N.; Ruff, M.; Vivat, V.; Chambon, P.; Gronemeyer, H.; Moras, D. Nature 1995, 378, 681–689. (7) Egea, P. F.; Mitschler, A.; Rochel, N.; Ruff, M.; Chambon, P.; Moras, D. EMBO J. 2000, 19, 2592–2601. (8) Chambon, P. FASEB J. 1996, 10, 940–954. (9) Lera, A. R.; Bourguet, W.; Altucci, L.; Gronemeyer, H. Nature ReV. Drug Disc. 2007, 6, 811–820. (10) Nahoum, V.; Pe´rez, E.; Germain, P.; Rodriguez-Barrios, F.; Manzo, F.; Kammerer, S.; Lemaire, G.; Hirsch, O.; Royer, C. A.; Gronemeyer, H.; Lera, A. R.; Bourguet, W. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 17323– 17328. (11) Otero, M. P.; Torrado, A.; Pazos, Y.; Sussman, F.; Lera, A. R. J. Org. Chem. 2002, 67, 5876–5882. (12) Bourguet, W.; Vivat, V.; Wurtz, J.-M.; Chambon, P.; Gronemeyer, H.; Moras, D. Mol. Cell 2000, 5, 289–298. (13) Lu, J.; Cistola, D. P.; Li, E. Biochemistry 2006, 45, 1629–1639. (14) Peet, D. J.; Doyle, D. F.; Corey, D. R.; Mangelsdorf, D. J. Chem. Biol. 1998, 5, 13–21. (15) Kitaura, K.; Sawai, T.; Asada, T.; Nakano, T.; Uebayasi, M. Chem. Phys. Lett. 1999, 312, 319–324. (16) Kitaura, K.; Ikeo, E.; Asada, T.; Nakano, T.; Uebayasi, M. Chem. Phys. Lett. 1999, 313, 701–706. (17) Nakano, T.; Kaminuma, T.; Sato, T.; Akiyama, Y.; Uebayasi, M.; Kitaura, K. Chem. Phys. Lett. 2000, 318, 614–618. (18) Nakano, T.; Kaminuma, T.; Sato, T.; Fukuzawa, K.; Akiyama, Y.; Uebayasi, M.; Kitaura, K. Chem. Phys. Lett. 2002, 351, 475–480. (19) Fedorov, D. G.; Kitaura, K. Theoretical Development of the Fragment Molecular Orbital (FMO) Method. In Modern Methods for Theoretical Physical Chemistry of Biopolymers; Starikov, E. B., Lewis, J. B., Tanaka, S., Eds.; Elsevier: Amsterdam, The Netherlands, 2006; pp. 3-38. (20) Nakano, T.; Mochizuki, Y.; Fukuzawa, K.; Amari, S.; Tanaka, S. Developments and Applications of ABINIT-MP Software Based on the Fragment Molecular Orbital Method. In Modern Methods for Theoretical Physical Chemistry of Biopolymers; Starikov, E. B., Lewis, J. B.; Tanaka, S., Eds.; Elsevier: Amsterdam, The Netherlands, 2006; pp. 39-52. (21) Amari, S.; Aizawa, M.; Zhang, J.; Fukuzawa, K.; Mochizuki, Y.; Iwasawa, Y.; Nakata, K.; Chuman, H.; Nakano, T. J. Chem. Inf. Model. 2006, 46, 221–230. (22) Nemoto, T.; Fedorov, D. G.; Uebayasi, M.; Kanazawa, K.; Kitaura, K.; Komeiji, Y. Comput. Biol. Chem. 2005, 29, 434–439. (23) Fukuzawa, K.; Kitaura, K.; Uebayasi, M.; Nakata, K.; Kaminuma, T.; Nakano, T. J. Comput. Chem. 2005, 26, 1–10. (24) (a) Fukuzawa, K.; Mochizuki, Y.; Tanaka, S.; Kitaura, K.; Nakano, T. J. Phys. Chem. B 2006, 110, 16102-16110; (b) ibid, 24276. (25) (a) Fukuzawa, K.; Komeiji, Y.; Mochizuki, Y.; Kato, A.; Nakano, T.; Tanaka, S. J. Comput. Chem. 2006, 27, 948–960. (b) Fukuzawa, K.; Komeiji, Y.; Mochizuki, Y.; Kato, A.; Nakano, T.; Tanaka, S. J. Comput. Chem. 2007, 28, 2237–2239. (26) Mochizuki, Y.; Nakano, T.; Koikegami, S.; Tanimori, S.; Abe, Y.; Nagashima, U.; Kitaura, K. Theor. Chem. Acc. 2004, 112, 442–452.

12094 J. Phys. Chem. B, Vol. 112, No. 38, 2008 (27) Mochizuki, Y.; Koikegami, S.; Nakano, T.; Amari, S.; Kitaura, K. Chem. Phys. Lett. 2004, 396, 473–479. (28) Ito, M.; Fukuzawa, K.; Mochizuki, Y.; Nakano, T.; Tanaka, S. J. Phys. Chem. B 2007, 111, 3525–3533. (29) Ito, M.; Fukuzawa, K.; Mochizuki, Y.; Nakano, T.; Tanaka, S. J. Phys. Chem. A 2008, 112, 1986–1998. (30) Ishikawa, T.; Mochizuki, Y.; Amari, S.; Nakano, T.; Tokiwa, H.; Tanaka, S.; Tanaka, K. Theor. Chem. Acc. 2007, 118, 937–945. (31) ABINIT-MP, Version 4.1 is available from the website of RSS21 project www.ciss.iis.u-tokyo.ac.jp/rss21/ (accessed Aug. 2008). (32) (a) Berman, H. M.; Westbrook, J.; Feng, Z.; Gilliland, G.; Bhat, T. N.; Weissig, H.; Shindyalov, I. N.; Bourne, P. E. Nucleic Acids Res. 2000, 28, 235–242. (b) The RCSB Protein Data Bank (http://www.rcsb. org/). (33) MOE, Version 2005.06; Chemical Computing Group: Montreal, Canada, 2005. (34) Gampe, R. T.; Montana, V. G.; Lambert, M. H.; Miller, A. B.; Bledsoe, R. K.; Milburn, M. V.; Kliewer, S. A.; Willson, T. M.; Xu, H. E. Mol. Cell 2000, 5, 545–555. (35) Wang, J.; Cieplak, P.; Kollman, P. A. J. Comput. Chem. 2000, 21, 1049–1074. (36) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A., Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision D.01; Gaussian, Inc.: Wallingford, CT, 2004. (37) Bayly, C. I.; Cieplak, P.; Cornell, W.; Kollman, P. A. J. Phys. Chem. 1993, 97, 10269–10280. (38) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Kollmann, P. A. J. Am. Chem. Soc. 1993, 115, 9620–9631. (39) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. J. Am. Chem. Soc. 1995, 117, 5179–5197. (40) Case, D. A.; Pearlman, D. A.; Caldwell, J. W.; Cheatham, T. E., III; Wang, J.; Ross, W. S.; Simmerling, C. L.; Darden, T. A.; Merz, K. M.; Stanton, R. V.; Cheng, A. L.; Vincent, J. J.; Crowley, M.; Tsui, V.; Gohlke, H.; Radmer, R. J.; Duan, Y.; Pitera, J.; Massova, I.; Seibel, G. L.; Singh, U. C.; Weiner, P. K.; Kollman, P. A. AMBER 7; University of California: San Francisco, CA, 2002.

Ito et al. (41) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926–935. (42) Case, D. A.; Darden, T. A.; Cheatham, T. E., III; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R.; Merz, K. M.; Wang, B.; Pearlman, D. A.; Crowley, M.; Brozell, S.; Tsui, V.; Gohlke, H.; Mongan, J.; Hornak, V.; Cui, G.; Beroza, P.; Schafmeister, C.; Caldwell, J. W.; Ross, W. S.; Kollman, P. A. AMBER 8; University of California: San Francisco, CA, 2004. (43) Taiji, M.; Narumi, T.; Ohno, Y. Petascale Special-purpose Computer for Molecular Dynamics Simulations. In Petascale Computing; CRC Press, 2008; pp 183-210. (44) Taiji, M.; Narumi, T.; Ohno, Y.; Futatsugi, N.; Suenaga, A.; Takada, N.; Konagaya, A. A Petaflops Special-Purpose Computer System for Molecular Dynamics Simulations. In Proc. Supercomputing 2003, Phoenix, AZ, November 15–21, 2003; on CD-ROM, 2003. (45) Narumi, T.; Ohno, Y.; Okimoto, N.; Koishi, T.; Suenaga, A.; Futatsugi, N.; Yanai, R.; Himeno, R.; Fujikawa, S.; Ikei, M.; Taiji, M. A 185 TFLOPS Simulation of Amyloid-forming Peptides from Yeast Prion Sup35 with the Special-purpose Computer System MDGRAPE-3. In Proc. Supercomputing 2006, Tampa, FL, Nov. 11–17, 2006; on CD-ROM, 2006. (46) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684–3690. (47) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327–341. (48) Miyamoto, S.; Kollman, P. A. J. Comput. Chem. 1992, 13, 952– 962. (49) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. J. Chem. Phys. 1995, 103, 8577–8593. (50) BioStation Viewer, Version 8.0 is available from the website of RSS21 project http://www.ciss.iis.u-tokyo.ac.jp/rss21/ (accessed Aug. 2008). (51) Schu¨tz, M.; Hetzer, G.; Werner, H.-J. J. Chem. Phys. 1999, 111, 5691–5705. (52) Hetzer, G.; Schu¨tz, M.; Stoll, H.; Werner, H.-J. J. Chem. Phys. 2000, 113, 9443–9455. (53) Pulay, P. Chem. Phys. Lett. 1983, 100, 151–154. (54) Pulay, P.; Saebø, S. Theor. Chim. Acta 1986, 357–368. (55) Pipek, J.; Mezey, P. G. J. Chem. Phys. 1989, 90, 4916–4926. (56) van Lipzig, M. M. H.; ter Laak, A. M.; Jongejan, A.; Vermeulen, N. P. E.; Wamelink, M.; Geerke, D.; Meerman, J. H. N. J. Med. Chem. 2004, 47, 1018–1030. (57) Sˇponer, J.; Jurecjeka, P.; Hobza, P. J. Am. Chem. Soc. 2004, 126, 10142–10151. (58) Xie, W.; Gao, J. J. Chem. Theory Comput. 2007, 3, 1890–1900. (59) Komeiji, Y.; Inadomi, Y.; Nakano, T. Comput. Biol. Chem. 2004, 28, 155–161. (60) Okiyama, Y.; Watanabe, H.; Fukuzawa, K.; Nakano, T.; Mochizuki, Y.; Ishikawa, T.; Tanaka, S.; Ebina, K. Chem. Phys. Lett. 2007, 449, 329– 335. (61) Mochizuki, Y.; Yamashita, K.; Murase, T.; Nakano, T.; Fukuzawa, K.; Takematsu, K.; Watanabe, H.; Tanaka, S. Chem. Phys. Lett. 2008, 457, 396–403.

JP803369X