ligand - gated ion channels ( lgics ) are multimeric membrane proteins that play an important role in fast synaptic transmission and in the modulation of cellular activity .
they are involved in many neurological disorders , such as alzheimer s disease , anxiety , epilepsy , and learning and attention deficit , and are target sites for many drugs , including anesthetics and muscle relaxants.(1 ) lgics exist in multiple conformations , of which the closed ( inactive ) and open ( active ) are the best - studied .
the equilibrium between these two conformations is affected by the binding of an agonist , which stimulates the conformational change ( gating ) between the closed and the open states .
once the channel is open , ions can flow through the channel , altering the postsynaptic neuron membrane potential .
the 5-hydroxytryptamine 3 receptor ( 5-ht3r ) is a member of the cys - loop lgic superfamily , which includes nicotinic acetylcholine ( nach ) , -aminobutyric acid ( gabaa , gabac ) , and glycine receptors.(2 ) these receptors are composed of five subunits arranged in a pentagonal array around a central pore .
each subunit has a large extracellular n - terminal domain , which contains the ligand binding site , a transmembrane domain , and an intracellular domain .
the transmembrane domain is composed of four -helices : m1 , m2 , m3 and m4 ; the m2 helix of each subunit lines the interior of the channel pore .
mutagenesis experiments , in combination with various computational techniques , including homology modeling , ligand - docking , classical molecular dynamics , and first principle calculations , have supplied important information on the ligand - binding mechanism ( see , e.g. , refs ( 37 ) ) , but a clear atomistic picture of the successive gating mechanism is still missing .
data show that regions located at the interface of the extracellular and transmembrane domains are critical in linking binding with channel opening ; these include the cys - loop , the 12 loop , and the m2m3 loop .
a model has been proposed for 5-ht3r gating whereby transcis isomerization of a specific proline ( pro308 , alternatively labeled as pro8 * ) at the apex of the m2m3 loop would induce a conformational change in m2 , resulting in a displacement of the hydrophobic girdle that occludes the pore and opening the channel.(15 ) to examine the alleged switching role of proline in 5-ht3r and to provide atomistic insights to complement the mutagenesis experiments , we have performed a computational study of a proline dipeptide ( n - acetylproline methylamide ) , the simplest model to study proline isomerization within a peptide , shown in figure 1 .
we also investigated a range of dipeptides of the proline analogues that were studied experimentally.(15 ) proline dipeptide .
a number of previous studies of proline dipeptides have focused on single - point , total energy calculations by first principles methods in vacuo , for which the entropic contribution to the free energy was calculated within the harmonic approximation ; when included , the solvent effects were evaluated implicitly . here ,
using a classical force field scheme and the metadynamics method , we provide the free energy maps of the selected dipeptides and an evaluation of the free energy differences between the cis and trans isomers that takes into account the complexity of the basins of attraction without assumptions on their shape . at variance with previous free energy studies
, we include the water solvent explicitly so as to quantify the effects of both intra- and intermolecular hydrogen bonds .
our study is relevant not only for the ligand - gated ion channels that have inspired it , but also for other biomolecules in which the transcis isomerization of the prolyl peptide bond plays an important role.(24 ) it provides a reference database for experiments , where the properties of proline analogues have been coherently calculated under the same conditions .
moreover , it can be easily extended to other proline analogues , pressure and temperature conditions , and solvents .
we created the initial structures of proline dipeptide and its analogues using the program leap , contained in the amber simulation package,(25 ) by terminating the proline analogues with the two groups acetamide ( ace ) and n - methyl ( nme ) .
all the simulations were performed with the classical molecular dynamics amber(25 ) and orac(26 ) packages using the amber2003 force field.(27 ) we parametrized the nonstandard proline analogues by assigning to them partial charges from ab initio calculations ; for consistency and for a better agreement with ab initio calculations in vacuo , we used ab initio partial charges also for pro .
specifically , we calculated the esp atomic partial charges(28 ) of pro and its analogues at a density functional theory ( dft ) level with the cpmd package,(29 ) using the pbe gradient - corrected xc - functional,(30 ) norm - conserving martinstroulliers pseudopotentials,(31 ) a kinetic energy cutoff for the wave function expansion in plane waves of 70 ry , orthorhombic cells of size 15  14  14 a , and the hockney method for neglecting the electrostatic interactions between neighboring cells.(32 ) we determined the partial charges assigned to the force field for the trans isomer ; however , we also calculated the dftesp charges of proline dipeptide for the cis isomer and verified that the differences were small ( mainly confined to the redistribution of charge between n and h in the nh group , whose total charge remained however very similar )
. hence , the trans isomer charges are a good representation for all the conformers . for the simulations in aqueous solution , we surrounded the dipeptides with 431 tip3p water molecules(33 ) in a periodically repeated orthorhombic simulation cell .
this was initially equilibrated at room temperature and atmospheric pressure using a time step of 0.5 fs for the molecular dynamic and evaluating the electrostatic contributions with particle mesh ewald.(34 ) the volume of the cell was kept fixed at the equilibrated value during the metadynamics , and the temperature was controlled using the nos thermostat,(35 ) whereas in gas phase , we used a stochastic velocity rescaling.(36 ) to assess the reliability of the amber2003 force field with the dftesp charges for the description of the proline dipeptide transcis isomerization , we compared selected potential energy profiles in vacuo as a function of the torsional angles  ( defined by the atoms ch3o1cc in figure 1 ) and  ( defined by n1ccn ) , obtained with amber2003 and dft calculations ( figure 3 ) .
, we used the same protocol as for the esp charges ; each point corresponded to a geometry optimization at fixed  and . the amber2003 force field qualitatively reproduces the location of the minima , maxima , and saddle points obtained with dft and , most importantly for our study , the relative stability of the trans and cis isomers .
quantitatively , the differences in energy between the cis and trans absolute minima ( ect ) were 2.4 kcal / mol with amber and 3.6 kcal / mol with dft ; the barriers were at   90 and   0 15.3 kcal / mol with amber and 21.0 kcal / mol with dft and , similarly , at   + 90 and   0 15.7 kcal / mol with amber and 20.6 kcal / mol with dft .
the locations of our trans and cis minima and the corresponding etc are in fairly good agreement with previously published dft ( b3lyp ) and hartreefock results .
the larger error of the force field in the energy barrier with respect to the cistrans energy difference is likely to be due to the relevance of electronic structure effects for the transition states in which the partial double bond character of the n1c1 bond is broken .
although we were aware that the force field might not be able to capture subtle electronic structure effects , it was , indeed , able to reproduce the relevant features and allowed us to map the free energy surfaces ( fes ) in explicit solvent , a task that would be extremely computationally expensive with ab initio methods .
potential energy profiles at selected values of  and  ( in deg ) obtained with dft ( red ) and the amber2003 force field with dftesp partial charges ( black ) for proline dipeptide in vacuo .
the zero of the energy scale is set to the absolute minimum . to evaluate the free energy surfaces
, we used the metadynamics method , which is an efficient way to sample the free energy of complex polyatomic systems in the space of a few course - grained quantities ( collective variables , cvs).(37 ) metadynamics works by adding a history - dependent potential , built as a sum of gaussians , whose role is to discourage the system from revisiting regions in the cv space already visited , thus accelerating rare events involving barrier crossing .
when metadynamics is used with properly chosen collective variables , it provides an unbiased estimate of the free energy surface .
however , the efficiency of the metadynamics method strongly depends on which collective variables are chosen as reaction coordinates for a particular mechanism . in the case of pro transcis isomerization ,
the most intuitive choice of cv is the imide torsional angle  ( cn1c1ch3 in figure 1 ) , which is equal to 180 for the trans and to 0 for the cis isomer .
however , during the isomerization , this angle is coupled to the out - of - plane deformation ( pyramidalization ) of the imide nitrogen , described by the improper dihedral angle  ( c1cn1c ) : this should also be specified in addition to  to unambiguously identify the saddle points .
alternatively , fischer et al.(40 ) proposed the use of the improper dihedral angle  ( ch3o1cc ) , which at the same time takes into account both the transcis isomerization and the nitrogen pyramidalization . similarly to  ,  = 0 for the cis isomer and  = 180 for the trans isomer ; at variance from  ,  is around  90 at the saddle points , independently on the nitrogen pyramidalization , thus effectively removing the need to specify . in our metadynamics calculations , we used as cvs  and also the torsional angle  ( n1ccn ) , which describes the orientation of the c - terminal amide .
this is important because both experimental(41 ) and theoretical(40 ) studies have shown that the proline transcis isomerization is affected by the interaction between the c - terminal group and either the lone pair of the imide nitrogen ( n1 ) or the imide carbonyl oxygen ( o1 ) .
these interactions may , in fact , lower the torsional barrier by decreasing the double bond character of the bond between n1 and c1 in a process of autocatalysis .
we built the history - dependent potential by summing up gaussians of width 9 and height w = 0.07 kcal / mol , with a stride of  = 400 fs .
we then estimated the free energy from the negative of the history - dependent potential at the end of the simulation
. converged free energies required about 40  000 gaussians each , for a simulated time of 16 ns . to accelerate the free energy calculations
, we adopted the multiple walker scheme , which allows an efficient parallelization of the method.(42 ) as discussed in refs ( 38 ) and ( 39 ) , the statistical error on the estimated free energy depends on the rate , w/. our parameters were relatively conservative and allowed us to obtain accurate evaluations of the free energy profiles .
the statistical errors were lower than 0.2 kcal / mol for the simulations in vacuo and lower than 0.3 kcal / mol for the simulations in aqueous solution .
we further decreased the error by taking the average of the history - dependent potential over the last 4 ns of simulation.(43 )
in figure 4 , we show the free energy maps as a function of  and  for proline dipeptide in vacuo and aqueous solution ; the zero in each graph is the absolute minimum .
we obtained the free energy difference in vacuo between the cis and trans isomers by integrating over the relevant regions so as to account for the different shapes of the basin of attractions .
our calculated value ( 2.4  0.2 kcal / mol ) is in good agreement with previously published dft b3lyp/6 - 31g+d(d ) and b3lyp/6 - 311++g(d , p ) calculations , in which vibrational frequencies were calculated for the stationary points in the gas phase and then used to compute free energies . in vacuo , the trans global minimum at   180 and   60 corresponds to a configuration with an intramolecular hydrogen bond between the nh group and o1 .
the cis isomer has two minima ; no intramolecular hydrogen bonds were possible in this configuration ; hence , the transcis isomerization disrupts the intramolecular hydrogen bond within the dipeptide .
the presence of intramolecular hydrogen bonds was highly affected by the environment . in vacuo ,
90 <  < 90 ) , with  ranging from 0 to 90 ; in water , its occurrence was reduced to 16% .
this was due to the competition between the formation of the intramolecular hydrogen bond and that of intermolecular hydrogen bonds with the water molecules .
specifically , o and o1 each formed , on average , 2.1 and 1.9 hydrogen bonds with water .
n1 did not form any hydrogen bonds with water , whereas n formed only 0.8 hydrogen bonds ( 0.1 as acceptor and 0.7 as donor ) .
these averages were evaluated for any value of  ; restricting the averages to the trans isomer gave almost identical values ( 2.0 and 1.9 hydrogen bonds between water molecules and o and o1 , respectively ) .
the influence of the solvent is clearly reflected in the fes , where the unique trans minimum basin in vacuo is substituted by a pair of trans minima in water , at   145 and 25 , as shown in figure 4 .
hence , in water , the trans basin is more similar to the cis one , where no intramolecular hydrogen bonds can be formed . because of the importance of the intermolecular hydrogen bonds , it is essential that explicit water molecules were included in the simulations .
top : free energy maps at room temperature for proline dipeptide in vacuo ( left ) and water ( right ) as a function of the dihedral angles  and . bottom : the corresponding free energy profile integrated over . the zero of the energy scale is set at the absolute minimum , and the spacing between lines is 0.2 kcal / mol .
we calculated a free energy difference between the cis and trans isomers in water equal to 1.0 kcal / mol ; with an error of 0.3 kcal / mol , this is in good agreement with experiments , which estimates 0.6 kcal / mol .
the reduction of the cistrans free energy difference in aqueous with respect to that in vacuo is in line with dft / b3lyp and implicit solvent calculations .
the isomerization barrier , obtained from the fes profiles integrated over  ( figure 4 ) increased from 13.8 kcal / mol in vacuo to 15.6 kcal / mol in aqueous solution , thus confirming the mechanism of autocatalysis determined by the presence of the intramolecular hydrogen bond .
although the environment in the ion channel is much more complex , we expect that simulations of the transcis proline switch in water represent a realistic environment .
in fact , the m2m3 loop , where pro308 is situated , has been shown to be solvent - accessible in a variety of cys - loop receptors .
hence , we believe that the studies in aqueous solution are a good representation of the in vivo situation .
we therefore investigated the proline analogues ( figure 2 ) , which have been experimentally studied in the 5-ht3r , in aqueous solution .
these analogues exhibit various ring features and preference for either of the two isomers.(47 ) the corresponding free energy maps are shown in figure 5 , where , as in figure 4 , the zero of the scale is aligned with the absolute minimum of each fes .
the cistrans free energy differences with respect to pro and the lowest isomerization barriers are shown in table 1 .
all proline analogues favored the trans isomers , except dmp and tbp , which favored the cis .
free energy maps at room temperature of the proline dipeptide analogues in aqueous solution as a function of the dihedral angles  ( horizontal , in deg ) and  ( vertical , in deg ) .
the isomerisation barrier is the smallest between the two calculated for clockwise and counterclockwise transcis isomerisation . in pip and aze , the size of the ring was , respectively , increased and decreased as compared to pro ; the relative stability of the cis and trans configurations was very similar in pip and pro , whereas the cis stability was slightly increased in aze .
the steric interaction between the two methyl groups added to c and the ace methyl group disfavored the trans isomer . in tbp , the steric interaction between the bulky tert - butyl group attached to c and the ace group also disfavored the trans isomer ; however , the distance between the three methyl groups in the tert - butyl and the ace methyl group was larger than for dmp , resulting in a weaker interaction and a reduced destabilization of the trans .
both the isomerization barriers of dmp and tbp were substantially lower than that of pro .
the relative transcis stabilities of the fluoroprolines were also similar to that of pro , with t-4f - pro slightly stabilizing the cis more than both pro and c-4f - pro ; in aqueous solution , the fluorine atoms in both residues formed hydrogen bonds with water molecules , which neutralized the main effects of their presence on the ring .
specifically , the fluorine atom acted as acceptor for , on average , 0.9 hydrogen bonds with water in the t-4f - pro simulation and 0.7 in the c-4f - pro simulation .
2-me - pro favored the trans isomer , due to the methyl group attached to c , which sterically interacted with the ace methyl group , disfavoring the cis configuration .
the addition of a methyl group to c , as in the 3-me - pro , did not significantly affect the cis preference with respect to pro because it is farther away from the isomerization bond .
in addition , 2,4-ch2-pro favored the trans isomer due to the steric interaction between the ace methyl group and the methylene group attached to c and c ; the presence of the methylene group substantially decreased the isomerization barrier with respect to pro .
evaluation of the hydrogen bonds formed during the metadynamics for all the proline analogues revealed strong similarity to pro , with averages for the number of hydrogen bonds between water molecules and o1 , o , nh , and n1 ranging from 1.9 to 2.2 , from 1.6 to 1.8 , from 0.5 to 0.9 , and from 0.0 to 0.1 , respectively . the percentages of occurrence of the intramolecular hydrogen bond were also not dissimilar from that of pro , ranging from 12% for pip to 20% for c-4f - pro .
this indicates that the effect of the solvent for the proline analogues was as crucial as it was for pro .
a transcis isomerization of a proline peptide bond is an attractive mechanism to change protein conformation and has been shown to be important in the function of several proteins .
the data we present here support previous studies that show that such a mechanism can open the pore of the 5-ht3 receptor.(15 ) the comparison between the relative free energies of the cis and trans isomers of a range of proline analogue dipeptides and those inferred from the relative ec50s of the 5-ht3 receptor activation allowed us to isolate the isomerization of the proline peptide bond from all the other concurrent events that occur when 5-ht binds to this receptor and lead to channel gating .
the strong correlation we found supports the proposal that transcis isomerization is , indeed , a realistic switching mechanism . in the previous published experimental paper,(15 )
the cistrans free energy differences correlated with the experimental ec50s were obtained from the literature for a range of different peptides at different experimental conditions , and the comparison was done only for a subset of the investigated proline analogues . here
, we evaluated the relative cistrans free energy difference for all the experimentally investigated proline analogues , including those resulting in nonfunctional channels , using the same model system under the same conditions .
hence , we considered our data more accurate and consistent . the theoretical protocol that we have developed can be easily extended to other residues and conditions so as to provide a reference database for future experimental measurements done in the same spirit of those by lummis et al.(15 ) in other proteins or involving other prolines .
our simulations also revealed specific interactions within the dipeptides and with the solvent that explain the free energy profiles , thus providing additional information on the atomistic details of the isomerization mechanisms .
in particular , the isomerization from trans to cis disrupted the hydrogen bond between the residues immediately preceding and following pro308 ( specifically , between the nh group of nme and the o1 atom of ace ) .
such a hydrogen bond , which was present in the trans isomer both in vacuo and in solution in different percentages , was , in fact , absent in the cis .
this absence is likely to be mirrored in the m2m3 loop of 5-ht3r . through preliminary calculations on extended peptides mimicking the m2m3 loop
, we verified that the hydrogen bond between the threonine and leucine that sandwich pro308 was , indeed , broken during the transcis isomerization ; the occurrence of another hydrogen bond between this threonine and another threonine in the m2m3 loop also decreased substantially during the isomerization : modifications in the interaction network might facilitate the movement of the m2 helix leading to the opening of the channel .
indeed , the transcis proline isomerization has the potential to bring many changes in the interactions with surrounding residues in the 5-ht3 receptor , as can be inferred in figure 7 . here
, pro308 is shown in its trans and cis forms , together with nearby residues in homology models for the 5-ht3 receptor . pro308 and
the two adjacent residues thr307 and leu309 , highlighted with sticks , are in the m2m3 loop ; glu213 is in the f loop ; arg245 is at the boundary between the extracellular and the transmembrane domain ; and tyr250 is at the start of the m1 helix .
these homology models were generated using modeller 6v2(51 ) as previously described,(3 ) with the extracellular domain based on the crystal structure of the acetylcholine binding protein ( achbp ) at 2.7  resolution ( pdb i d : 1i9b(52 ) ) and the transmembrane domain on the structure of the nach receptor at 4  resolution ( pdb i d : 2bg9(53 ) ) .
the cistrans free energy differences relative to proline , when plotted against those obtained from the ec50 values ( figure 6 ) , cluster mainly around two regions for the proline analogues that behave similarly to proline ( pip , aze , t-4f - pro and c-4f - pro ) and those with a clear preference for the cis isomer ( tbp and dmp ) .
they reveal a correlation coefficient of 0.98 , similar to the experimental data previously reported(15 ) and indicating a tight relationship between the population of cis and trans isomers and the receptor response .
the slope of the best fit , however , is different : from our theoretical data , it is m = 0.66 , whereas the previous value was m = 1.(15 ) the latter suggests that the equilibrium between closed and open channels is purely due to the equilibrium between trans and cis pro308 , but the new data suggest the relationship may be more complex and possibly nonlinear . a linear relationship in fact
would also predict that all analogues , including 2-me - pro and 2,4-ch2-pro would be functional , with ec50s significantly less than the maximum concentration that has been tested .
correlation between cistrans free energy differences from the metadynamics simulations and from experimental data ( ec50s ) on the activation of 5-ht3r;(15 ) the data for the mutants corresponding to functional channels are indicated as circles ; those corresponding to nonfunctional channels ( hence with no available ec50 ) , with triangles .
homology models of the 5-ht3 receptors where pro308 , in the trans ( left ) and cis ( right ) configurations , and selected surrounding residues are shown explicitly .
pro308 and the adjacent residues thr207 and leu209 are in the m2m3 loop , glu213 is in the f loop , arg245 is at the boundary between the extracellular and the transmembrane domains , and tyr250 is at the start of the m1 helix . this does not , however , mean that the original hypothesis is incorrect .
for example , the observed differences could be due to the inaccuracy in the empirical force field , which can not capture subtle electronic effects that can alter quantitatively , rather than qualitatively , the isomer population .
the channel environment , which was absent in our calculations , is also likely to play an important role in triggering and catalyzing the isomerization so as to reconcile the observed gating time scale with that inferred from the isomerization barrier .
in particular , there could be steric hindrance effects that might prevent isomerization paths allowed in water , such as in the case of 3-me - pro , which does not behave in the dipeptide much differently from proline , but in the experiments produced a nonfunctional channel . in summary , using classical molecular dynamics combined with state - of - the - art metadynamics , we have characterized the free energy surfaces of a series of proline dipeptide analogues as the simplest model of transcis proline switch and deduced information about the relative stability of their trans and cis conformers . although a number of issues , which would benefit from simulations of more complex systems , remain open , the excellent correlation between our calculated free energy difference and those deduced from the experimental ec50s support the idea of a switching role for pro308 for the gating of 5-ht3r .
 transcis isomerization of a proline peptide bond is a potential mechanism to open the channel of the 5-ht3 receptor . here 
 , we have used the metadynamics method to theoretically explore such a mechanism . 
 we have determined the free energy surfaces in aqueous solution of a series of dipeptides of proline analogues and evaluated the free energy difference between the cis and trans isomers . 
 these theoretical results were then compared with data from mutagenesis experiments , in which the response of the 5-ht3 receptor was measured when the proline at the apex of the m2-m3 transmembrane domain loop was mutated . 
 the strong correlation between the experimental and the theoretical data supports the existence of a transcis proline switch for opening the 5-ht3 receptor ion channel . 
Introduction
Methods
Results
Discussion and Conclusions