Chemical Senses Advance Access originally published online on October 21, 2005
Chemical Senses 2005 30(9):781-792; doi:10.1093/chemse/bji070
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Structural Activation Pathways from Dynamic Olfactory Receptor–Odorant Interactions
1 Department of Molecular and Cell Biology, University of Connecticut, Storrs, CT, USA, 2 Massachusetts Eye and Ear Infirmary, Harvard Medical School, Boston, MA, USA and 3 Department of Neurobiology and Yale Center for Medical Informatics, Yale University School of Medicine, New Haven, CT, USA
Correspondence to be sent to: Chiquito Crasto, Department of Neurobiology and Yale Center for Medical Informatics, Yale University School of Medicine, New Haven, CT, USA. e-mail: chiquito.crasto{at}yale.edu
 |
Abstract |
|---|
 Top Abstract IntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
We have simulated an odor ligand's dynamic behavior in the binding region of an olfactory receptor (OR). Our short timescale computational studies (up to 200 ps) have helped identify unprecedented postdocking ligand behavior of ligands. From in vacuo molecular dynamics simulations of interactions between models of rat OR I7 and 10 aldehyde ligands, we have identified a dissociative pathway along which the ligand exits and enters the OR-binding pocket—a transit event. The ligand's transit through the receptor's binding region may mark the beginning of a signal transduction cascade leading to odor recognition. We have graphically traced the rotameric changes in key OR amino acid side chains during the transit. Our results have helped substantiate or refute previously held notions of amino acid contribution to ligand stability in the binding pocket. Our observations of ligand activity when compared to those of experimental (electroolfactogram response) OR-activation studies provide a view to predicting the stability of ligands in the binding pocket as a precursor to OR activation by the ligand.
Key words: ligand binding, molecular dynamics, olfactory receptors
|
Introduction |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
An olfactory receptor (OR) is a G protein–coupled receptor (GPCR)—a superfamily of proteins identified by seven transmembrane (TM) helical domains connected by intracellular and extracellular loops. GPCRs respond to extracellular stimuli with a wide range of intracellular and intercellular functions (Ji et al., 1998
; Muller, 2000). The interactions of the odor ligands with ORs in the olfactory epithelium result in a signal transduction cascade, a key step in olfactory sensory perception (Shepherd, 1995).
After the publication of early drafts of the human genome (Lander et al., 2001; Venter et al., 2001), four groups have independently identified human ORs (Glusman et al., 2001; Zozulya et al., 2001; Niimura and Nei, 2003; Malnic et al., 2004). ORs from the mouse genome were also identified (Young et al., 2002; Zhang and Firestein, 2002; Godfrey et al., 2004). More recently, work has also been done in identifying dog ORs (Quignon et al., 2003; Olender et al., 2004b; Parker et al., 2004). The Olfactory Receptor Database (ORDB) (http://senselab.med.yale.edu/senselab/ordb) and its companion website OdorDB (http://senselab.med.yale.edu/senselab/odordb) list OR gene sequences for more than 35 organisms, in addition to a list of functionally characterized ORs (by odors).
Experimental odorant binding to ORs
Electroolfactograms (EOGs) and calcium imaging techniques are typically used to functionally characterize odor–OR interactions. Malnic et al. (1999) used calcium imaging studies to study the activation of a set of 14 mouse olfactory neurons to a combinatorial panel of odors of functional groups: alcohols, carboxylic acids, bromo- and dicarboxylic acids of carbon chain lengths varying between five and nine, and at varying concentrations. Araneda et al. (2000) also used EOGs to identify the responses of rat OR I7 to a panel of 90 odors. Their best results were obtained for eight-carbon atom chain aldehydes. Other functional activation studies have shown that phenylic aldehydes and ketones activate mouse OR M71 (Bozza et al., 2002) and the activation of mouse OR OR912-93 follows binding by 2-heptanone and 3-heptanone (Rouquier et al., 1999). Luu et al. (2004) have shown that key receptor residues are responsible for interacting with amino acid odor ligands in zebra fish and goldfish. Recent reports have pointed to ORs that are broadly tuned versus those whose responses are restricted to a few odors (Sanz et al., 2005).
Katada and coworkers have carried out comprehensive functional analyses of mouse receptors activated by lyral, eugenol, vanillin, and ethyl vanillin. Their experimental, modeling, and docking studies have shown that ligand dimerization contributes to inhibition or competitive antagonism (Katada and Touhara, 2004; Katada et al., 2005). Additionally, Shirokova et al. (2005) have recently shown that the type of G protein in the transduction cascade plays an important role in determining whether an odor will be an agonist or an antagonist.
We use a model of rat I7 (ORL11 in ORDB; entry M64386 in GenBank) OR in our studies. Rat I7 had been previously functionally characterized as responsive to octanal (eight-carbon atom chain aldehyde group) (Krautwurst et al., 1998; Zhao et al., 1998). Of all the odor molecules tested, rat I7 is known to be activated only by straight-chain aldehydes of length ranging between seven and nine carbon atoms. Mouse I7, which has very high sequence identity with its rat homolog but has a key residue mutation (V206I) that complements the shortage of one methyl group between octanal and heptanal, was earlier shown to be activated only by heptanal (Krautwurst et al., 1998). However, Bozza and coworkers (2002) later showed that both heptanal and octanal activated rat and mouse I7, experimentally. Hall et al. (2004), from theoretical binding-energy calculations following computational docking of these ligands in the receptor, confirmed these findings.
Most studies of functionally characterized ORs show that the binding region interactions are governed in part by one or more electrostatic interactions between the charged functional group and a residue on the OR TM domains. Additionally, van der Waals interactions also support the ligand's presence in the OR-binding region. Uchida et al. (2000) agree with the view that functional groups and additional nonbonded interactions are two determinants in activation. With both mouse and rat I7, the lysine on residue 164 on TM4 likely forms a Schiff base with the carbonyl oxygen of the aldehyde (Singer, 2000). Correlated mutational analyses showed that Lys164 and Asp204 were critical to ligand binding (Singer and Shepherd, 1994). In OR912-93, a histidine on TM4 finds itself in a position similar to K164 on TM4 and is likely the primary locus of electrostatic interactions. A similar interaction is found in rhodopsin where a lysine (TM7) is the primary covalent interaction with its sole activating ligand, retinal (Nakayama and Khorana, 1991). The role of a counter ion, that is, the role of Asp204, is played by a glutamate residue on TM3 (Sakmar et al., 1989) of rhodopsin.
The need for modeling studies
The first model of rat OR I7 was created by Singer (2000). We used this model in our studies. A description of the creation of this model appears in Materials and Methods. Singer (2000) also studied the docking and short timescale dynamics of aldehydes of varying chain lengths in the I7-binding pocket. Lys164 and Asp204, and to a lesser extent F262, were identified as key determinants in binding (and perhaps activating) I7. Results suggested that the carbonyl oxygen of octanal tended toward the formation of a Schiff base with the Lys164 nitrogen atom. Singer (2000) also identified competitive binding between the ligand oxygen and aspartate (D204) oxygen atoms toward the K164 nitrogen.
Except for bovine rhodopsin, no GPCR has been structurally characterized (Palczewski et al., 2000). The most recent rhodopsin structure has a resolution of 2.2 Å (Okada et al., 2004). The lack of a crystal structure makes OR–odor dynamic modeling challenging. Attempts to computationally elucidate OR structures thus far have involved homology modeling (Man et al., 2004), which uses bovine rhodopsin as a structural template, ab initio techniques (Vaidehi et al., 2002), or the creation of idealized canonical helical domains from diffraction-derived rhodopsin electron densities. Previous, related studies have described stable odor docking configurations (Floriano et al., 2000, 2004; Vaidehi et al., 2002). ORs are activated following binding with odor molecules. The binding pocket of ORs has been established from previous experimental studies supported by computational studies—the latter provide a view of the odor–OR interaction at the molecular level (Singer, 2000; Floriano et al., 2004).
OR–odor interactions are complicated: many odors activate one receptor, and vice versa (Floriano et al., 2000; Hall et al., 2004). The discriminatory nature of ORs toward odor molecules (being able to identify several thousand odors with far fewer functional receptors) further complicates the identification of a mechanism for olfactory perception.
In order for OR activation to occur, therefore, as a first step, an odor ligand interacts with the OR-binding pocket, activates the OR, initiating the signal transduction cascade, and then departs from the receptor-binding region, inactivating the OR. There is currently no study that identifies the activation–inactivation of ORs at a molecular level. Little is known of interactions of the ligand with ORs at longer timescales. Even less is known of the ligand's transit through the OR-binding region during OR activation. These studies discussed specific interactions of ligands with key binding-pocket amino acid residues, which mediate ligand–receptor binding and selectivity. Clearly then, longer timescale simulations of ORs and odors are needed. Increasingly available information necessitates efforts to understand the mechanism of olfaction at the molecular level.
Araneda et al. (2000) tested the activities of 90 ligands with rat I7 via EOG recordings after expressing OR I7 from olfactory sensory neurons. In further studying a narrower list of 10 aldehydes, they reported results that for OR activation molecular length, functional group and level of unsaturation in the carbon chain of the ligands was critical. We use models of 10 aldehyde ligands (chemical names of these compounds are identified in Materials and Methods) experimentally studied by Araneda and coworkers (2000) in our studies. Our simulations allow visualization of the ligand's behavior in the I7's binding pocket. Using these, we can compare how ligand stability (based on size and chemical character) might dictate OR activation.
Through our observations from molecular dynamics simulations of 10 odor ligands with the most studied and one of the first-identified OR, rat I7 (Buck and Axel, 1991), we identified a transit—dissociative—pathway within the OR's binding region and key amino acid residues that facilitated this transit. Additionally, our observations reveal certain trends in binding that correlate well with experimental results of EOG studies of I7 OR activation. Our computationally derived results will serve as predictors of OR–odor interactions.
|
Materials and methods |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
We carried out computational (molecular dynamics) OR–ligand simulations for longer durations than had been carried out previously to address the following key issues related to ligand binding: (1) substantiate or repudiate previously held notions of binding in terms of specific amino acid contributions, (2) recognize if short-term interactions and static docking configurations are sustained over longer time periods, (3) identify a pathway by which the ligand enters and leaves the OR-binding region—a transit pathway, (4) identify how (and which) amino acid residues in the binding pocket facilitate this transit event, and (5) predict, if possible, whether ligand stability in the OR-binding pocket can be compared to experimental OR activation leading to olfactory perception.
OR model used in this study
The OR model used in this study is one of the rat ORs, I7 (Singer, 2000), one of the first ORs cloned and identified (Buck and Axel, 1991). Singer (2000) created this OR model using a 7.5-Å resolved electron diffraction structure (Schertler, 1998) of rhodopsin. We summarize here the steps that were taken in creating this model. A detailed description of the model creation, docking, and preliminary dynamics results, which also include a rationale for the model-creation methodology (experimental support for modeling strategies) can be found in Singer (2000). Singer used hydrophobicity moments to identify the helical TM regions from rat I7's primary structure. Canonical helices were then created, and the energy of each helix was individually minimized. Each side chain on the helices was cast in its lowest energy conformation by accessing a rotamer library. TM domains were assembled by positioning these helices by using the electron diffraction–derived densities for rhodopsin while setting a dielectric constant (equal to 1) to mimic the surrounding membrane. The helices were then rotated to preserve the hydrophobicity construct of a GPCR assembly; that is, hydrophobicity moment was pointed directly away from the TM assembly. The entire structure thus created was then energy minimized and refined using molecular dynamics.
Validation of I7 model to match the high-resolution structure of rhodopsin
Singer (2000) validated the above methodology by constructing a bacteriorhodopsin model and comparing it favorably to bacteriorhodopsin's crystal structure [2BRD in the Protein Data Bank (PDB)] (Goto and Iwamoto, 1997) with acceptable root mean squared deviations (RMSDs).
We have revalidated the structure of our model, comparing it to the currently available highest resolution structure of rhodopsin, a 2.2-Å resolved X-ray diffraction structure by Okada et al. (2004) [1U19 in the PDB]. Since our I7 model lacks loops, the rhodopsin structure–PDB file was edited to contain only the TM regions. The homology modeling software Modeller8v2 (Marti-Renom et al., 2000) was used to perform the superimposition of the backbone atoms for each of the seven TM helices. The computed positional RMSD was 2.2 Å. The rhodopsin X-ray crystallographically derived structure also shows kinks, which are artifacts of its primary structure. The idealized helices in our model ignore these kinks, which would include rhodopsin structure–specific features and influence the structure. Figure 6 describes the two superimposed structures.
Docking
We built 10 ligand models whose binding results had been studied experimentally (Araneda et al., 2000). Each ligand was docked in the binding pocket of the OR. This was followed by dynamic simulation of the system for up to 200 ps. The ligands: n-octanal, trans-2-octenal, 2-octynal, trans,trans-2,4-octadienal, tetrahydrocitral (3,7-dimethyloctanal), 7-methyloctanal, 3-methyl-4-(4-methyl-cyclohexyl)-propanal, citronellal (3,7-dimethyl-6-octenal), citral (3,7-dimethyl-2,6-octadienal), and 2,5,7-trimethyl-2,6-octadienal were constructed using the InsightII suite of software (http://www.accelrys.com/insight/). Ligand conformational energies were minimized using the Discover module in InsightII. We added hydrogen atoms to our I7 OR model to create a system of pH 7.0. Atomic charges were assigned using consistent valence force field (CVFF).
We used DOCK (Gschwend and Kuntz, 1996; Gschwend et al., 1996) to identify the ideal binding configurations of the ligands in the binding pocket of the I7 OR. Using every atom in the OR model, the DMS (dot molecular surface) module (Richards, 1977) was used to calculate a solvent-accessible molecular surface area for the I7 model, and DOCK's SPHGEN (SPHere GENerator) module identified cavity site points in the receptor. Only spheres associated with I7 atoms that represented biologically relevant docking regions were retained for final docking. We discarded spheres that were inaccessible to the extracellular surface of the OR; these spheres were structurally "below" the TM3 and TM4 crossover plane in the I7 model. This strategy does not spatially impact docking; it also saves time and computational resources. GRID, another module in DOCK, was used to generate force fields and interaction parameters to compute intermolecular binding. DOCK used spheres that were retained to compute spatial restraints based on van der Waals interactions. Flexible_Ligand a module in DOCK allowed the modification of torsion angles in the ligand.
Ligand positioning in DOCK is spatial. It is based primarily on van der Waals interactions and energy calculations from assigned force fields. DOCK scores docking conformations on an increasing docking energy scale. Additionally, DOCK allows users to recognize polar regions in molecules by upscaling the electrostatic component of the interatomic interactions. The highest scoring (lowest energy) docked configuration was selected for molecular dynamics simulation studies. For every ligand tested, the lowest energy-docking configuration tethered the Lys164 basic nitrogen to the aldehyde carbonyl oxygen in support of Singer's (2000) contention of the Lys164–octanal interaction based on correlated mutational analysis.
Molecular dynamics simulations
Each docked system was introduced into the InsightII Discover module. Ligand–OR system energies were minimized using a conjugate gradient and CVFF. In order to preserve the integrity of the helical regions, the positional parameters of alpha-carbon atoms in the helical domains were restrained. While introduction of a simulated bilayer (Vaidehi et al., 2002) results in a more realistic representation of the membrane OR protein system, we believe that it influences neither the nature of docking (restricted to the binding pocket) nor the OR–ligand binding energies—the latter are used to establish levels of OR activation (Floriano et al., 2004). Allowing dynamic motion for every atom of the bilayer in the simulation would add time and consume considerable computational resources, which in our work—10 simulations of I7-ligand systems—is forbidding. Fixing the positional parameters of the alpha-carbon atoms and the terminal amino acids of each helix preserves the integrity of the helical domains and the binding pocket, yet provides unrestricted movement for the amino acid side chains. Introducing a bilayer would preclude the necessity to constrain the alpha-carbon atoms of the TM domains. Each ligand–OR simulation was performed for durations between 100–200 ps at 300 K with 100 fs of equilibration time.
During simulation, trajectories were stored every 10 fs. Relevant distances (between amino acids and ligands) required to identify rotameric behavior of amino acids were traced and graphed as the simulation proceeded. Graphical data were acquired in InsightII and postprocessed using Microsoft Excel and Ploticus (http://ploticus.sourceforge.net). The graphs in Figures 3–5 and the supplemental information available online (http://senselab.med.yale.edu/senselab/OR-model) depict the behavior of the ligand with respect to key amino acid residues in the binding pocket over the duration of the simulation. The graphical traces are of distance or angle changes in amino acids with ligand motion during the simulation. The blue lines in the graphs represent positional changes. The orange characteristics overlaying the blue lines are a measure of the moving averages every picosecond. Moving averages smooth out random peaks that result from high-frequency conformational fluctuations. They also help to identify characteristic behavioral trends.
The positional parameters at every time step (in PDB format) were rendered using the visualization program Pymol (http://pymol.sourceforge.net). These files were used to create movies of the simulations. The total energy of the system, resolved into van der Waals and Coulombic energies calculated at every time step of the simulation were also graphed. These energy graphs for the simulations of I7 with octanal (experimentally known to activate I7 strongly) and I7 with 2,5,7-trimethyl-2,6-octadienal (a nonactivator) are described in Figures 7 and 8. The graphing methods used are the same as for traces of positional changes: dark blue lines trace the energies at for every time step; orange lines show moving averages. In Figures 7 and 8, the distance graph between the ligand carbonyl oxygen and Lys164 nitrogen is reintroduced to provide a perspective on how the system energy changes as the simulation proceeds—especially during the transit peaks.
For each of the 10 ligands, the supplemental information includes links to the PDB-formatted positional parameters of the ligand docked in the receptor-binding pocket (energy minimized position), figures depicting the ligand in its initial docked position, the ligand in its maximum (exit) displacement from the energy-minimized docked (EMD) position, graphs of relevant rotameric changes in key amino acids as the simulation proceeds, and graphs of energy calculations during the simulation and a link to a movie (MPEG) of the ligand during the exit event.
|
Results |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
Transit event
Our simulation results for 10 ligands (from a docked position) with the rat OR I7 show that, within the limits of our model and simulation parameters, the ligand is associated with a "transit event" along an exit (and possible entry) pathway. Movies of most of the simulations show that the ligand backs away from its primary electrostatic tether with Lys164 toward an opening on the extracellular side of the receptor. We define the transit event as follows. The ligand backs away by at least 1 Å (for most ligands, more than 2.0 Å) from its EMD position and its electrostatic interaction with the Lys164 (TM4).
We quantified an exit as the ligand carbonyl oxygen atom's greatest RMSD from its starting energy-minimized position. The transit event correlates with structural changes in the binding pocket. The interaction between Asp204 and Lys164 increases, decreasing their interresidue distance. The dihedral angle between the planes of the two carboxylic groups of the aspartate residue also changes to maximize its interaction with the Lys164 amine. The ligand distance from Val285 in the "rear" of the binding pocket decreases commensurately.
Figure 1a (EMD position) and b (maximum displacement from EMD position) show the initial and transit event configurations for octanal in the binding pocket of I7. The top–down perspective shows an opening, despite the close packing of receptor atoms. Figure 1b shows the ligand (in its exit position) nearer to this opening, indicating that the exit or dissociative path leads to this opening. This is line with the exit pathway from the electrostatic carbonyl–Lys164 interaction. Figure 2 shows the same transit pathway (for octanal) from a side-view perspective with both the EMD position (blue) and the position of maximum distance from Lys164 (orange). The arrows in Figure 2 show the potential exit path from the receptor.
|
|
Figure 3 traces the ligand (carbonyl oxygen) and Lys164 nitrogen atom distance over the course of a simulation for all 10 ligands. The peaks which depict the transit event are broad. The periodicity of the peaks, however, does not follow any particular trend. This is despite doubling the simulation time for some ligands. Sharper peaks indicate high-frequency conformational changes, which the ligand and residues constantly undergo during the simulation. The amplitudes of exit peaks do not necessarily depict the maximum ligand displacement from K164.
|
Response of amino acids in the OR-binding pocket during simulations
Graphs in Figures 4 and 5 depict the position and angle changes in key residues of I7's binding pocket during simulations with two ligands: n-octanal and 2,5,7-trimethyl-2,6-octadienal. Experimental results (Araneda et al., 2000) have shown that the former is a strong activator for I7, while the latter does not activate the receptor at any concentration. Similar graphs for simulations involving the other eight ligands are available online with the supplemental information.
|
|
In Figures 4 and 5, the first graph traces the ligand–Lys164 distance. Peaks in the graph indicate (though not always) transit events marked by increasing distance between the ligand carbonyl group and the Lys164. The next graph traces the distance between the terminal carbon atom (atoms if the terminal is branched) and V285, opposite to K164 in the binding pocket. To test Singer's (2000)
assertion of competitive binding, the dihedral angle of the aspartic acid would change to engender an electrostatic interaction with the cationic amine of K164. The third graph traces this angle change. Graphs that trace the RMSD displacements of the carbonyl oxygen atom from the minimized position of the ligand in the binding pocket prior to the simulation run follow. The positions of the phenylalanines F205 (TM5) and F262 (TM6) have given rise to conjectures that their phenyl rings act as "lids" over the binding pocket, preventing the premature exit of the ligand. Figures 1a,b and 2 show that the phenylalanine contributions to the exit are indirect at best. F262 was one of the sites identified by Singer and Shepherd (1994) as having a marginal contribution to ligand binding. The last two graphs trace the motion of the phenyl rings of F205 and F262. Ligand behavior during simulations
From the graphs for every simulation, a composite of ligand behavior begins to emerge. The trends from these interactions, with regard to ligand stability in the binding pocket and OR activation, can be compared to experimental results (Araneda et al., 2000). The figures and movies representing the transit events and the graphs of the simulation are at http://senselab.med.yale.edu/senselab/OR-Model. What follows is a description of binding trends during the simulations based on the structures of each ligand. The graphs of the interactions of amino acid residues with the ligand provide a view to the stability of the ligand in the binding pocket.
N-octanal, trans-2-octanenal, and 2-octynal are dissimilar only in the degree of unsaturation. Two transit peaks are observed for each ligand. The experimental results indicate that the activation strengths decrease with increase in unsaturation. The graphs for octynal indicate that the first peak is not a transit peak because the Lys164–ligand distance increase is not accompanied by a decrease in Lys164–Asp204 distance. The movie shows that at this peak, the ligand undergoes a conformational adjustment. This increases the distance between the ligand and Lys164. The second peak, however, meets the criteria for a transit event. Ligand oxygen RMSD displacements are consistent with a transit event for all three ligands. The phenylalanines, F205 and F262, show strong displacement, but their motions are correlated with each other than with the exit event. The movies indicate that F205 and F262 are not in the exit path. Their displacements, while lagging behind transit events, indicate that these amino acid residues move away from the exit path to make the exit sterically feasible.
Similar behavior is observed in F205 and F262 during the simulation of trans,trans-2,4-octadienal. EOG responses indicate that this ligand is the strongest activator for I7 (Araneda et al., 2000). Two peaks are observed during the simulation. The true transit peak is observed at 
70 ps; the first is not a true transit peak—as previously defined. When methyl groups are added to octadienal as in the case of citral (3,7-dimethyl-2,6-octadienal), the bulky groups at the carbon chain termini and the restricted rotation due to the double bonds prevents the ligand from stabilizing in the binding pocket. The movie and the graphs indicate that the ligand tries to exit almost immediately and spends most of the simulation time displaced from its electrostatic tether with Lys164. Toward the end of the simulation, the ligand finds a more stable conformation and enters the binding pocket. But an attempt to exit does not follow the exit path observed for the ligands that have been experimentally observed to activate the I7 OR. The movie shows that the ligand burrows deeper into the binding pocket. Experimental results indicated that citral minimally activates I7 and only when it is used in significantly higher concentrations than other ligands.
On addition of one more methyl group, the bulkiness and unsaturation prevents the ligand 2,5,7-trimethyl-2,6-octadienal from fitting in the binding pocket. Experimentally, this ligand does not activate I7 at even very high concentrations. During docking, this ligand was forced into the binding pocket at higher docking energies than the other ligands. The movie indicates that the ligand puckers and is occasionally positioned vertically in the binding pocket. Similar puckering and a vertical configuration are observed in 3-methyl-4-(4-methyl-cyclohexyl)-propanal from steric repulsions due to a six-member ring in the carbon chain. This ligand also activates I7 poorly. The movie for this simulation indicates that none of the peaks in the graphs coincide with exit events.
For ligands with a methyl substituent on the terminal carbon atom such as citral, distance traces for this methyl carbon (designated C7-M) are also included in the graphs. Tetrahydrocitral (3,7-dimethyloctanal) is one of the stronger activators of I7 (Araneda et al., 2000). At docking, this ligand, which has two methyl groups and can freely rotate along the axis of the carbon chain due to an absence of unsaturation, docked more vertically. Two transit events (peaks) are observed. The exit events are relatively not as long lived as those for other ligands because as soon as the 3,7-dimethyloctanal backs off from its EMD position, the presence of an additional terminal methyl substituent is repelled causing the ligand to quickly regain its EMD position. For 7-methyloctanal and citronellal (3,7-dimethyl-6-octenal), which are both good activators, one definitive transit peak is observed. Because of the free rotation in the former, the terminal C8 carbon atom behaves differently from C7-M. For dimethyloctanal, after the first peak, the other smaller, random peaks during the rest of the simulation do not depict exit events. During the citronellal simulation, F205 shows a positional change at the end of the transit event and remains in that position till the end of the simulation. The movie of this simulation shows that at 100 ps, the ligand abruptly changes position and is thus stabilized for 10 ps. No correlated positional changes are seen in other amino acids monitored during the simulation.
The movies of the simulations indicate that the phenylalanine rings while not in the exit path do perform a lidlike motion. These positional changes in the phenyl rings are to facilitate the movement of the ligand by minimizing steric repulsion from the aromatic electron cloud.
Energy calculations during simulation
The total energy of the OR–odor system was traced during each of the 10 simulations. Energies were traced in order to ensure that the exit events were not associated due to system destabilization due to increased energy. The electrostatic interaction between Lys164 and the octanal carbonyl group has been implicated as a primary interaction. Therefore, the Coulombic and van der Waals components of the total energy were graphed separately. The exit peaks in the Lys164–odor distance graph are also plotted. Figures 7 and 8 show the energy traces for octanal and 2,5,7-trimethyl-2,6-octadienal. In all but one case, the energy graphs show increasingly lower energy as the molecular dynamics simulations involve energy stabilization. Similar energy graphs for the other ligand are available on the OR-Model web page.
The characteristics of the total energy graphs are strongly weighted toward Coulombic interactions, while the van der Waals interactions are relatively stabilized. The Coulombic interactions show very small (and mostly insignificant) peaks that coincide with the exit events. But no significant instances of destabilization of the system can be found to coincide with the exit events. The flatness of the graphs after 100 ps indicates that our OR–odor system is stabilized. Only the energy graphs for citronellal are associated with a significant drop in energy at the beginning of a well-defined (in graphs and movies) exit peak. The energy during the exit peak however remains steady. This energy does not change for the duration of the simulation.
Simulations for octanal, octanal, and octynal were carried out for 100 ps. As mentioned previously, each shows two exit peaks. The energy at the simulation temperature, however, shows that 100-ps simulation time is insufficient for complete energy stabilization of our OR–odor system. The remaining seven simulations were carried out for twice that time.
|
Discussion |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
Our work has been carried out using the first published computational model of rat I7 (Singer, 2000
). Barring the availability of an OR crystal structure, the only structurally comparable aspects of all GPCRs are helical TM domains. Our I7 model replicates the TMs—helical densities are adequately resolved even from low-resolution electron densities and include amino acid side chains in energetically stable conformations. As described in Materials and Methods, we have validated our model by comparing it to the highest resolution structure of rhodopsin currently available (Okada et al., 2004). The RMSD between the positions of backbone atoms of rhodopsin and our model is 2.2 Å. Figure 6 shows that the only discrepancy between the two structures is in TM1. Our observations (from movies and graphs of the simulations) show that TM1 is not involved in the odorant binding.
|
Limitations and rationale
In this work, simulations have been carried out on a model without interhelical loops. The absence of loops does not interfere with the observed ligand–receptor dynamics and transit events within the binding pocket. This is because each simulation proceeds from an initial docked position, which is within the binding pocket. It is possible that the dynamic shifting of loops during a simulation would influence the behavior of binding-pocket residues due to a push–pull effect on the TM domains. The dynamic behavior (time lag to transit events) of the F262 and F205 is possibly an artifact of the lack of loops in our model. On the other hand, since phenylalanine movements are not correlated with the transit peaks, it is possible that the dynamics of these residues are indicative of possible active and inactive states of the OR. Our modeling system and simulation conditions are inadequate to identify this.
Our simulations were also carried out in vacuum. It would benefit studies of this nature to perform dynamic simulations of ligands in ORs in an aqueous medium. Increasing the parameters and degrees of freedom by including the loops and the bilayer and modeling several thousand molecules of water in the system, however, are prohibited by limited available computational facilities. The transit pathway identified in our study is, however, not directly affected by either phenylalanine.
In a recent publication, Man et al. (2004) have used a human OR model created using homology modeling. This model has the advantages of using a high-resolution, X-ray–determined structural model of rhodopsin. A study of the results of simulating OR–ligand interactions using a homology-derived model would be interesting to identify whether similar transit pathways are observed. Such a model, based on structural template matching (the sequence homology between I7 and rhodopsin is low) may, however, introduce rhodopsin structure–specific biases into the model. These include differences in lengths of loops (Otaki and Firestein, 2001) and kinks (Olender et al., 2004a; Yohannan et al., 2004) in TM helical domains. Olender and coworkers (2004a) have studied the kinks in homology-derived OR structures to address some of these issues. Variations in sequence similarities across ORs and other GPCRs also introduce rotameric bias into homology-modeled structures of ORs (Visiers et al., 2002). The MEMBSTRUK (Floriano et al., 2000; Vaidehi et al., 2002) protocol devised at the California Institute of Technology uses the low-resolution electron densities of rhodopsin followed by ab initio techniques to model the amino acid residue side chains. Their results have often been used to successfully mimic experimental results (Floriano et al., 2004).
Comparing dynamic behavior and experimental OR activation
In discussing receptor activation trends, one should dissociate from the notion of direct comparisons of experimental results with our observations. EOGs reflect the recorded activity of a cascade of events in vitro or in vivo. They measure very different responses (Watt and Storm, 2001). It should be recognized that experimental (wet) and computational modeling and simulation (dry) results can be used in conjunction to better understand the mechanism leading to olfactory perception. Araneda et al. (2000) graded the activation of I7 by each ligand graphically. The activation peaks indicate, however, that the relative activations of the more strongly activating ligands—octanal, octenal, octadienal, tetrahydrocitral, methyloctanal, and citronellal—are statistically indistinguishable (figure 4a in Araneda et al., 2000). Relative activations then cannot be quantified.
Our results (graphs in Figures 3–5
and on the OR-Model website) show that theoretical intrabinding pocket interactions are predictive of strong experimental binding and activation of a ligand in the OR (Araneda et al., 2000). These interactions govern the possibility of, and isolate the pathway by which, the ligand entering the pocket, stabilizing in it, and then exiting. Ligand stability (and possibly, activation) is restricted by steric factors such as type of substituent, branching, as well as degree of unsaturation in the carbon atom chain. Our studies show that citral and 2,5,7-trimethyl-2,6-octadienal fail to stabilize in the binding pocket of the I7 model or show a definite exit pathway. Experimentally, they fail to activate rat I7.
Energy profiles
The energy profiles in Figures 7 (for n-octanal) and 8 (for 2,5,7-trimethyl-2,6-octadienal) indicate that the system stabilizes as the molecular dynamics simulation proceeds without any instability that causes the ligand to escape the OR-binding pocket. The stability of the energy profile indicates that the system is stable at the simulation temperature. During our short runs, we have identified a definitive dissociative pathway through the ligand's binding pocket. This dissociation is nonrandom. Our energy profiles of the simulations show that the exit events take place even after the system's energy has stabilized. We have resolved energies into Coulombic and van der Waals to identify that the delinking of the Lys164–ligand electrostatic interaction is not an artifact of increased system energy. We have shown this consistently for 10 ligands of varying substituents and levels of unsaturation.
|
|
Since, in most cases, the energy profile levels off after 120 ps (most of the exit peaks having occurred before that), one might surmise, that under the given simulation conditions for out OR–odor system longer simulation times would not add to the conclusions drawn. The energy characteristics of octanal, octenal, and octynal show a decline in energies but no leveling off. Clearly, the simulation times of 100 ps are inadequate. The isolation of the exit pathway that is repeated for two exit events for each of these ligands and the other seven is the primary conclusion of this work.
|
Conclusions |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
We carried out simulation studies for 10 eight-carbon atom chain aldehydes in the rat I7 OR. I7 was the first OR identified. It is also the most studied, both experimentally and computationally. The experimental precedence and the impetus for our work, was set by the results of Araneda and coworkers (2000)
. Our results have shown that for a ligand to activate a receptor, it should be dynamically stable in the receptor-binding region. Steric factors play an important role in such stabilizations. Our simulation runs begin with a ligand docked in its minimum energy conformation in the binding pocket. Our results attest to an exit pathway from the receptor. In most of our simulations, even if the ligand is not docked in a position for a facile exit, when conformational changes allow it (after as long as 100 ps) to be in the exit pathway, an exit occurs. While, with increasing energy (temperature) within the docked system during a simulation, an exit is inevitable, the pathway for exit and the influence of specific amino acid residues from the binding pocket during the exit is noteworthy. Our simulations clearly show that such a path in and out of the binding pocket exists.
One of the comparable trends between the experiments of Araneda et al. (2000) and our work is from observations of ligands that failed to activate the OR. Simulation trajectories for these (experimentally nonactivating) ligands in our model of I7 show the lack of a clearly defined transit path and transit event, indicating that the ligand will minimally activate the receptor. Studies such as this are likely predictors of OR binding by ligands and consequently, OR activation.
|
Acknowledgements |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
This work was generously supported by grant 2 P01 DC 004732-05 from the National Institute for Deafness and Communicative Disorders, National Institutes of Health. P.L. is supported as a part-time undergraduate researcher at the Department of Neurobiology, Yale University School of Medicine. The authors would like to thank Professor Gordon Shepherd at the Department of Neurobiology, Yale University School of Medicine, for his critical evaluation of this manuscript.
|
References |
|---|
|
TopAbstractIntroductionMaterials and methodsResultsDiscussionConclusionsAcknowledgementsReferences |
|---|
Araneda, R.C., Kini, A.D. and Firestein, S. (2000) The molecular receptive range of an odorant receptor. Nat. Neurosci., 3, 1248–1255.[CrossRef][ISI][Medline]
Bozza, T., Feinstein, P., Zheng, C. and Mombaerts, P. (2002) Odorant receptor expression defines functional units in the mouse olfactory system. J. Neurosci., 22, 3033–3043.
Buck, L. and Axel, R. (1991) A novel multigene family may encode odorant receptors: a molecular basis for odor recognition. Cell, 65, 175–187.[CrossRef][ISI][Medline]
Floriano, W.B., Vaidehi, N. and Goddard, W.A. 3rd (2004) Making sense of olfaction through predictions of the 3-D structure and function of olfactory receptors. Chem. Senses, 29, 269–290.
Floriano, W.B., Vaidehi, N., Goddard, W.A. 3rd, Singer, M.S. and Shepherd, G.M. (2000) Molecular mechanisms underlying differential odor responses of a mouse olfactory receptor. Proc. Natl Acad. Sci. USA, 97, 10712–10716.
Glusman, G., Yanai, I., Rubin, I. and Lancet, D. (2001) The complete human olfactory subgenome. Genome Res., 11, 685–702.
Godfrey, P.A., Malnic, B. and Buck, L.B. (2004) The mouse olfactory receptor gene family. Proc. Natl Acad. Sci. USA, 101, 2156–2161.
Goto, K. and Iwamoto, M. (1997) Evidence of alpha-helix slidings during bacteriorhodopsin photocycle-energetics coupling. Tohoku J. Exp. Med., 182, 15–33.[CrossRef][ISI][Medline]
Gschwend, D.A., Good, A.C. and Kuntz, I.D. (1996) Molecular docking towards drug discovery. J. Mol. Recognit., 9, 175–186.[CrossRef][ISI][Medline]
Hall, S.E., Floriano, W.B., Vaidehi, N. and Goddard, W.A. 3rd (2004) Predicted 3-D structures for mouse I7 and rat I7 olfactory receptors and comparison of predicted odor recognition profiles with experiment. Chem. Senses, 29, 595–616.
Ji, T.H., Grossmann, M. and Ji, I. (1998) G protein-coupled receptors. I. Diversity of receptor-ligand interactions. J. Biol. Chem., 273, 17299–17302.
Katada, S., Hirokawa, T., Oka, Y., Suwa, M. and Touhara, K. (2005) Structural basis for a broad but selective ligand spectrum of a mouse olfactory receptor: mapping the odorant-binding site. J. Neurosci., 25, 1806–1815.
Katada, S. and Touhara, K. (2004) A molecular basis for odorant recognition: olfactory receptor pharmacology. Nippon Yakurigaku Zasshi, 124, 201–209.[Medline]
Krautwurst, D., Yau, K.W. and Reed, R.R. (1998) Identification of ligands for olfactory receptors by functional expression of a receptor library. Cell, 95, 917–926.[CrossRef][ISI][Medline]
Lander, E.S., Linton, L.M., Birren, B. et al. (2001) Initial sequencing and analysis of the human genome. Nature, 409, 860–921.[CrossRef][Medline]
Luu, P., Acher, F., Bertrand, H.O., Fan, J. and Ngai, J. (2004) Molecular determinants of ligand selectivity in a vertebrate odorant receptor. J. Neurosci., 24, 10128–10137.
Malnic, B., Godfrey, P.A. and Buck, L.B. (2004) The human olfactory receptor gene family. Proc. Natl Acad. Sci. USA, 101, 2584–2589.
Malnic, B., Hirono, J., Sato, T. and Buck, L.B. (1999) Combinatorial receptor codes for odors. Cell, 96, 713–723.[CrossRef][ISI][Medline]
Man, O., Gilad, Y. and Lancet, D. (2004) Prediction of the odorant binding site of olfactory receptor proteins by human-mouse comparisons. Protein Sci., 13, 240–254.
Marti-Renom, M.A., Stuart, A.C., Fiser, A., Sanchez, R., Melo, F. and Sali, A. (2000) Comparative protein structure modeling of genes and genomes. Annu. Rev. Biophys. Biomol. Struct., 29, 291–325.[CrossRef][ISI][Medline]
Muller, G. (2000) Towards 3D structures of G protein-coupled receptors: a multidisciplinary approach. Curr. Med. Chem, 7, 861–888.[ISI][Medline]
Nakayama, T.A. and Khorana, H.G. (1991) Mapping of the amino acids in membrane-embedded helices that interact with the retinal chromophore in bovine rhodopsin. J. Biol. Chem., 266, 4269–4275.
Niimura, Y. and Nei, M. (2003) Evolution of olfactory receptor genes in the human genome. Proc. Natl Acad. Sci. USA, 100, 12235–12240.
Okada, T., Sugihara, M., Bondar, A.N., Elstner, M., Entel, P. and Buss, V. (2004) The retinal conformation and its environment in rhodopsin in light of a new 2.2 A crystal structure. J. Mol. Biol., 342, 571–583.[CrossRef][ISI][Medline]
Olender, T., Feldmesser, E., Atarot, T., Eisenstein, M. and Lancet, D. (2004a) The olfactory receptor universe—from whole genome analysis to structure and evolution. Genet. Mol. Res., 3, 545–553.[Medline]
Olender, T., Fuchs, T., Linhart, C., Shamir, R., Adams, M., Kalush, F., Khen, M. and Lancet, D. (2004b) The canine olfactory subgenome. Genomics, 83, 361–72.[CrossRef][ISI][Medline]
Otaki, J.M. and Firestein, S. (2001) Length analyses of mammalian G-protein-coupled receptors. J. Theor. Biol., 211, 77–100.[CrossRef][ISI][Medline]
Palczewski, K., Kumasaka, T., Hori, T., Behnke, C.A., Motoshima, H., Fox, B.A., Le Trong, I., Teller, D.C., Okada, T., Stenkamp, R.E., Yamamoto, M. and Miyano, M. (2000) Crystal structure of rhodopsin: a G protein-coupled receptor. Science, 289, 739–745.
Parker, H.G., Kim, L.V., Sutter, N.B., Carlson, S., Lorentzen, T.D., Malek, T.B., Johnson, G.S., DeFrance, H.B., Ostrander, E.A. and Kruglyak, L. (2004) Genetic structure of the purebred domestic dog. Science, 304, 1160–1164.
Quignon, P., Kirkness, E., Cadieu, E., Touleimat, N., Guyon, R., Renier, C., Hitte, C., Andre, C., Fraser, C. and Galibert, F. (2003) Comparison of the canine and human olfactory receptor gene repertoires. Genome Biol., 4, R80.[CrossRef][Medline]
Richards, F.M. (1977) Areas, volumes, packing and protein structure. Annu. Rev. Biophys. Bioeng., 6, 151–176.[CrossRef][ISI][Medline]
Rouquier, S., Stubbs, L., Gaillard-Sanchez, I. and Giorgi, D. (1999) Sequence and chromosomal localization of the mouse ortholog of the human olfactory receptor gene 912-93. Mamm. Genome, 10, 1172–1174.[CrossRef][ISI][Medline]
Sakmar, T.P., Franke, R.R. and Khorana, H.G. (1989) Glutamic acid-113 serves as the retinylidene Schiff base counter ion in bovine rhodopsin. Proc. Natl Acad. Sci. USA, 86, 8309–8313.
Sanz, G., Schlegel, C., Pernollet, J.C. and Briand, L. (2005) Comparison of odorant specificity of two human olfactory receptors from different phylogenetic classes and evidence for antagonism. Chem. Senses, 30, 69–80.
Schertler, G.F. (1998) Structure of rhodopsin. Eye, 12(Pt 3b), 504–510.
Shepherd, G.M. (1995) The cognitive neurosciences. In Gazzaniga, M.S. (ed.), The Cognitive Neurosciences: A Handbook for the Field. MIT Press, Cambridge, MA, pp. 105–102.
Shirokova, E., Schmiedeberg, K., Bedner, P., Niessen, H., Willecke, K., Raguse, J.D., Meyerhof, W. and Krautwurst, D. (2005) Identification of specific ligands for orphan olfactory receptors. G protein-dependent agonism and antagonism of odorants. J. Biol. Chem., 280, 11807–11815.
Singer, M.S. (2000) Analysis of the molecular basis for octanal interactions in the expressed rat 17 olfactory receptor. Chem. Senses, 25, 155–165.
Singer, M.S. and Shepherd, G.M. (1994) Molecular modeling of ligand-receptor interactions in the OR5 olfactory receptor. Neuroreport, 5, 1297–1300.[ISI][Medline]
Uchida, N., Takahashi, Y.K., Tanifuji, M. and Mori, K. (2000) Odor maps in the mammalian olfactory bulb: domain organization and odorant structural features. Nat. Neurosci., 3, 1035–1043.[CrossRef][ISI][Medline]
Vaidehi, N., Floriano, W.B., Trabanino, R., Hall, S.E., Freddolino, P., Choi, E.J., Zamanakos, G. and Goddard, W.A. 3rd (2002) Prediction of structure and function of G protein-coupled receptors. Proc. Natl Acad. Sci. USA, 99, 12622–12627.
Venter, J.C., Adams, M.D., Myers, E.W. et al. (2001) The sequence of the human genome. Science, 291, 1304–1351.
Visiers, I., Ballesteros, J.A. and Weinstein, H. (2002) Three-dimensional representations of G protein-coupled receptor structures and mechanisms. Methods Enzymol., 343, 329–371.[ISI][Medline]
Yohannan, S., Faham, S., Yang, D., Whitelegge, J.P. and Bowie, J.U. (2004) The evolution of transmembrane helix kinks and the structural diversity of G protein-coupled receptors. Proc. Natl Acad. Sci. USA, 101, 959–963.
Young, J.M., Friedman, C., Williams, E.M., Ross, J.A., Tonnes-Priddy, L. and Trask, B.J. (2002) Different evolutionary processes shaped the mouse and human olfactory receptor gene families. Hum. Mol. Genet., 11, 535–546.
Zhang, X. and Firestein, S. (2002) The olfactory receptor gene superfamily of the mouse. Nat. Neurosci., 5, 124–133.[ISI][Medline]
Zhao, H., Ivic, L., Otaki, J.M., Hashimoto, M., Mikoshiba, K. and Firestein, S. (1998) Functional expression of a mammalian odorant receptor. Science, 279, 237–242.
Zozulya, S., Echeverri, F. and Nguyen, T. (2001) The human olfactory receptor repertoire. Genome Biol., 2, 18.
Accepted September 30, 2005
CiteULike 
Connotea 
Del.icio.us What's this?
This article has been cited by other articles:
|
|
 |
|
  N. Benbernou, S. Tacher, S. Robin, M. Rakotomanga, F. Senger, and F. Galibert Functional Analysis of a Subset of Canine Olfactory Receptor Genes J. Hered., August 3, 2007; (2007) esm054v2. [Abstract] [Full Text] |