IWPLS'09 Potential of Mean Force of Ion Permeation through alpha7 nAChR Ion Channel Jens Krüger and Gregor Fels* Department of Chemistry, University of Paderborn, Warburger Strasse 100, 33100 Paderborn, Germany fels@uni-paderborn.de Associate Editor: Sandra Gesing and Jano van Hemert ABSTRACT opened, enabling the influx of Na+ and Ca2+ ions and the efflux of Many neuronal diseases such as Alzheimer’s dementia are related K+ ions (Fels et al. 1982). This causes a change in electrical poten- to a loss of inter-neuron communication. The nicotinic acetylcholine tial on the postsynaptic membrane that in turn enables the further receptor (nAChR) plays a crucial role in this process and is severely propagation the neuronal action potential (Buckingham et al. affected upon disease progression. Successful therapy approaches 2009). Individuals affected by Alzheimer’s dementia lack a suffi- rely on modulation of response signals, initiated by the flux of ions cient number of functional postsynaptic nAChR. Hence, the influx through the receptor integrated ion channel at the post-synaptic of Na+ is not sufficient to depolarize the membrane and to ignite a membrane. We here present the comparison of two methods for new action potential (Buckingham et al. 2009). Accordingly, un- calculating the potential of mean force of nAChR mediated ion per- derstanding the gating mechanism of nAChR ion channels should meation, in terms of accuracy and performance. be helpful in improving symptomatic Alzheimer’s therapy. We have evaluated and compared two methods for the calculation of the potential of mean force (PMF) of ion permeation through the 1 INTRODUCTION nAChR, Umbrella Sampling (US) and Steered Molecular Dynam- The nicotinic acetylcholine receptor (nAChR) is a cation-selective ics (SMD) pulling, in order to understand the changes in electrical channel made up from five homologous subunits symmetrically potentials connected to channel gating. The PMF corresponds to arranged around a central pore, thus forming a structure with five- the barrier a permeating ion has to overcome. It is crucial for the fold pseudo-symmetry (Lester et al. 2004). The nAChR at the neu- understanding of selectivity and conductivity of the nAChR. The romuscular endplate is a heteropentamer with α2βγδ subunits, required simulation steps as well as the quality of the resulting PMFs have been evaluated and are presented in a comparative while the neuronal receptor type used in this study is an α7 homo- manner. Special emphasis lies in the performance on two different pentamer. Cryo-electron microscopy studies of the Torpedo mar- high performance computing clusters (HPC) and the influence on morata nAChR have revealed three domains: an extracellular do- required CPU hours and real simulation time. main hosting the neurotransmitter binding sites; a transmembrane domain forming a channel across the cell wall and an intracellular domain providing binding sites for cytoskeletal proteins (Miyaza- 2 METHODS wa et al. 2003; Unwin 2005). The extracellular part of the protein A homology model for the transmembrane part of human α7 nAChR was is homologous to that of the acetylcholine binding protein found in constructed on the basis of an electron microscopy structure of torpedo the snail Lymnaea stagnalis (Brejc et al. 2001). The transmem- marmorata (PDB code 1OED) (Unwin 2005). The details of homology brane domains of each subunit consist of four membrane spanning modeling are described elsewhere (Kelly 1999; Wallace and Roberts 2004; helices named M1 to M4. The M2 helix is the pore lining part Gomaa et al. 2007). In brief, after aligning the sequence (SwissProt code while M4 faces the lipid environment. P36544) (Peng et al. 1994) with the template, backbone positions are as- The nAChR plays a crucial role in neuronal signaling and therefore signed for identical residues. Then loops are generated and candidates are chosen according to a Boltzman-weighted criterion. Then side chain data is holds a key function for learning and memory processes. It is se- assembled from an extensive rotamer library. The best intermediate model verely affected by neuronal diseases such as Alzheimer’s dementia is chosen based on electrostatic solvation energy and/or a packing score (Gotti and Clementi 2004). When an action potential arrives at the and finalized by energy minimization. The procedure is included in MOE- presynaptic membrane of a neuronal synapse, a biochemical cas- 2008.10. (Molecular Operating Environment Chemical Computing Group, cade leads to the release of the neurotransmitter acetylcholine. Inc., 1010 Sherbrooke St. W, Suite 910 Montreal, Quebec, Canada H3A Upon binding of two acetylcholine molecules to the extracellular 2R7) part of nAChR at the postsynaptic membrane, the integrated ion channel, located over 20 Å away from the ligand binding site, is All MD simulation were carried out using GROMACS-4.0.4 (van der Spoel et al. 2005; Hess et al. 2008) with the Gromos96 (ffG45a3) force field (Schuler et al. 2001). The temperature of the peptide, lipid and the solvent To whom correspondence should be addressed. * were separately coupled to a v-rescale thermostat with a coupling time of 2009 1 J. Krüger and G. Fels 0.1 ps. Semi-isotropic pressure coupling was applied with a coupling time Plots and pictures were made with Xmgrace-5.1.22, VMD-1.8.7 and MOE- of 1.0 ps and a compressibility of 4.5 · 10-5 bar-1 for the xy-plane as well as 2008.10. for the z-direction. Long range electrostatics were calculated using the particle-mesh Ewald (PME) algorithm with grid dimensions of 0.12 nm and interpolation order 4. Lennard-Jones and short-range Coulomb interactions were cut off at 1.4 and 0.8 nm respectively. The topology for the lipid bilayer (POPC (16:0-18:1 Diester PC, 1- Palmitoyl-2-Oleoyl-sn-Glycero-3-Phosphocholine) are the same used in earlier studies (Krüger and Fischer 2008) and were originally created on the basis of the parameters of Chandrasekhar et. al. (Chandrasekhar et al. 2003). Each pore model was embedded into a hydrated POPC bilayer system, by removing overlapping lipid and water molecules. After minimization the systems were equilibrated for 3x1 ns while initial position constraints were stepwise reduced (protein, backbone, C-alpha). In order to reduce stress induced by lateral pressure fluctuations each pore model was simulated during the equilibration phase using surface-tension pressure coupling with a tension of 37.5 mN/m (Krüger and Fischer 2009). The complete systems for US and pulling SMD consisted of the protein, 250 POPC, 17250 SPC- waters, 537 Na+ and 512 Cl-. US simulations follow the procedure described by Hub and de Groot (Hub and de Groot 2008). Only pressure coupling in the xy-plane was enabled while keeping the z-direction fixed. The starting configurations for the US simulations are based on three 10 ns equilibrium simulations. Three distinct channel conformations were prepared as described above and a full set of Fig. 1. Close view on the sodium ion at the gorge portion of the pore. US simulations carried out for each of them. The channels were divided While the hydration shell is slightly deformed in this confined space into 0.25 Å wide sections along the central pore axis. The ion was placed Thr244 and Ser249 show a stabilizing interaction. The hydrogen bond subsequently into the center of each section, removing overlapping waters, network is shown by dotted blue lines. followed by a thorough energy minimization. A harmonic restraint of 4000 kJ/mol was applied on the ion position along the pore axis. The sub- sequent simulations were carried out for 300 ps each. The PMFs were constructed with the WHAM procedure of Hub and de 3 RESULTS Groot (Hub and de Groot 2008) implemented into g_wham included in The sodium ion permeation through the nAChR was studied with GROMCAS. The first 50 ps were omitted and a cyclic correction for the two different methods for the construction of PMFs. The PMFs periodic system was applied, using an alpha of 1.75. The statistical error represent the energy barrier a passing ion has to overcome. This was estimated with a bootstrap analysis (N = 42). During pulling SMD semi-isotropic pressure coupling was used applying 1 barrier is directly correlated to conductivities accessible by expe- bar in z-direction while keeping the xy-plane fixed. The ion was placed at riments. Therefore PMFs are invaluable tools to access scientific the top or bottom of the simulation system on the central pore axis. The problems such as ion selectivity, channel open-/closing, mutations virtual spring attached to the ion had a constant of 100 kJ/mol/nm2 and was inside the pore or allosteric modulation. moved at 0.00375 nm/ns. Two times twelve independent simulations were carried out, pulling the ion through the pore from either side. Basically a PMF describes the probability to encounter the ion at a The construction of PMFs follows the method of Anishkin and Sukharev certain position along the reaction coordinate, compared to the (Anishkin and Sukharev 2004) considering the pulling as irreversible work bulk phase (equation 1). against the opposing friction forces and assuming a constant friction coeffi- cient. The friction coefficient was fitted for each independent simulation till PMF ( z ) = − k BT ln( P ( z ) / P0 ) (1) both ends equal a zero potential, which reflects the boundary condition of a free ion in bulk water. For US, which can be applied to a broad variety of other problems, The error estimate for the SMD PMFs was calculated using block averag- the probability is derived from umbrella histograms of the re- ing considering them to be correlated fluctuating quantity. The sets are strained atom (Kumar et al. 1992; Kumar et al. 1995; Beckstein divided in a number of blocks and averages are calculated for each block. and Sansom 2006; Hub and de Groot 2008). The error for the total average is calculated from the variance between averages of the m blocks Bi as follows: stderr2 = ∑(Bi-)2/(m·(m-1)). The The construction of PMFs for pulling SMD follows the Langevin complete derivation is given in the literature (Hess 2002). equation (equation 2) considering the forces exerted on the ion along the reaction coordinate as irreversible work against the fric- The simulations were run on a DELL Studio XPS (8 cores, i7 920) and on tion force (Gullingsrud et al. 1999; Anishkin and Sukharev 2004). facilities of the Paderborn Center for Parallel Computing PC2 Nz (http://wwwcs.uni-paderborn.de/pc2/, Arminius (400 cores, Xeon 3.2 GHz PMF ( z ) = ∑ ( ΔFn − γ Δzn / Δtn ) Δzn (2) EM64T, Infiniband) and Bisgrid (64 cores, dual-core Opteron 2.8 GHz, n =1 Infiniband)). 2 PMF of Ion Permeation through alpha7 nAChR Both approaches usually assume that the end points lie within the The conductive properties and time dependent characteristics of an bulk phase of the solvent. As the ion does not experience any po- ion channel are determined by various features of the pore, as well tential at these points, they have to be equal and zero (equation 3). as the influence of the extra- and intra-cellular parts. Beside gener- al aspects like the length of the pore and its diameter, the orienta- PMF ( z start ) = PMF ( zend ) = 0 (3) tion and polarity of side chains pointing into the lumen of the pore are highly relevant. The energy barrier of ion permeation not only has an enthalpic part but usually has also a large entropic part (Por- tella et al. 2008). As the hydrated ion enters the mouth of the pore the move ability of waters between the protein and the ion is more 70 limited. Consequently degrees of freedom are lost such as rotation- 60 al and translational degrees for the waters as well as the protein 50 side chains. The ion experiences repulsion and stabilization de- pending on the specific topology of the channel. Hydrophobicity 40 and hydrophilicity determine how well a passing ion is stabilized 30 within the narrow part of the pore (Krüger and Fischer 2009). It can be stated that first hydration shell of the sodium ion inside the 20 nAChR pore is never removed. At the gorge portion of the pore 10 deformations can be observed (Figure 1), which are directly stabi- 0 lized by serine and threonine residues (Thr244 and Ser249), which are found in rings along the pore (Bertrand et al. 1993). The so- -10 phisticated counter play between hydrophilic (e.g. Ser and Thr) -20 and hydrophobic (e.g. Leu and Ile) residues in the gorge area of the channel determines its conductivity and selectivity. 70 60 3.1 Umbrella Sampling 50 US simulations yield robust and well reproducible PMFs as long as 40 a high sampling density can be achieved. This can be a challenging G(z) (kJ mol ) -1 endeavor if the probability to encounter the ion at the gorge portion 30 of a pore is low. Related to this problem is the artificial (partial) 20 loss of the hydration shell of the ion during system preparation, when it is placed in hydrophobic areas of the pore. 10 Another commonly observed effect for pore US is due to small 0 sampling errors, that the ends of the PMF are not equal unlike being proposed in equation 3.. This effect can be overcome with -10 cyclic correction as described by Hub and de Groot (Hub and de -20 Groot 2008). 70 In this study the nAChR channel was prepared in three slightly 60 different ways, yielding three identical systems in terms of se- quence, secondary structure, lipid environment and solvent con- 50 tent. The difference lies within the small conformational diversity. 40 For each of them at least 300 restrained US were carried out. The single 300 ps simulations took 7.8 hours each using 4 cores. 30 The height of the PMFs was determined to be 51.4, 63.1 and 20 49.5 kJ/mol. The standard deviation for US PMFs determined by bootstrap analysis is in most cases around 2.5 kJ/mol, but in some 10 cases can reach as high as 7.6 kJ/mol (Figure 2). This is largely 0 influenced by the convergence of each simulation. Especially for the gorge portion of the channel larger deviations have to be ex- -10 pected. -20 An earlier study using US on restrained nAChR-M2 helices em- -4 -2 0 2 4 bedded into a bilayer-mimetic slab made from CH4 molecules z (nm) reports a barrier of 10.5 kbT for the sodium ion permeation (Beck- stein and Sansom 2006). These are just 25.5 kJ/mol or factor 2 less Fig. 2. PMF for the permeation of Na+ through the nAChR derived than the results presented in this study. The deviation has to be by US. The three curves are based on the same structural model, but contributed to the simulation conditions or more likely to the anal- each starting configuration was prepared in a slightly different way. ysis. The block data created by the modified WHAM analysis tool The permeation barrier lies for all three independent simulations be- from Alan Grossfield (http://membrane.urmc.rochester.edu/wham) tween 50 and 60 kJ/mol. The gray areas correspond to the standard used in that study can easily be misrepresented in Xmgrace by a deviation derived by bootstrap analysis. factor of 2. The misinterpretation has to be presumed comparing to 3 J. Krüger and G. Fels Fig. 3. Cut through the transmembrane part of the nAChR (orange, purple helices). Three independent permeation pathways of Na+ (blue) are indicated by thin gray lines, which are derived from different SMD pullings. the PMF heights determined in this study using two different me- estimate than having an unphysical flattening of the potential with thods. uncertain error. Ongoing research aims at overcoming this prob- lem. In terms of computational cost the two pulling analysis do not bear any differences. 3.2 Steered Molecular Dynamics - Pulling According to equation 2 different pulling velocities should not When an ion is pulled through a pore like the transmembrane part affect the final PMF. In order to evaluate the dependence of the of the nAChR it follows the natural path the ion would take in a friction coefficient on different pulling velocities a sodium ion was living organism. As shown in Figure 3 this path may vary signifi- pulled through the nAChR pore with 4 different velocities, each cantly for each individual permeation event. Hence different dis- with 6 repetitions. As shown in Figure 4 the absolute height of the tances were covered in different times, yielding slightly different PMFs is not affected within the margin of error. The increasing velocities. Therefore the friction coefficient has to be recalculated minimum around -2 nm has to be contributed to the assumption of for each simulation separately in order to fulfill equation 3. a constant friction coefficient. The higher the velocity of the ion The shape of pulling PMFs is slightly deformed as a constant fric- the less this assumption holds true. As stated above the deforma- tion coefficient was used. It was suggested by Gullingrud et. al. to tion of the PMF is related to this limitation. For the highest pulling use a velocity autocorrelation approach in order to apply a mod- velocities of 0.03 nm/ns a flattening of the energy barrier can be ulated friction coefficient along the reaction coordinate (Gul- observed. This effect is best described as ‘rupture’. The ion devel- lingsrud et al. 1999). As Patargias et. al. have demonstrated on ops such a high speed that it cannot develop all atomic interactions restraint Vpu pores, this approach negates the deformation (Patar- with surrounding residues. This can be compared to an experimen- gias et al. 2009). Unfortunately the ion experiences a complete tal condition under strong electric field called ‘electroporation’ different microenvironment at the mouth of the pore than it does showing analogy to some extend (Böckmann et al. 2008). within the gorge portion. Therefore the autocorrelation approach The PMF for pulling from the N- to the C-terminus is shown in the used in that study easily overestimates the friction within the pore upper panel of Figure 5, while the PMF for the opposite direction yielding a too low potential. For this study it was found to be more is shown below. Mimicking the influx of sodium into the cell acceptable to have an unphysical deformation with accurate error 4 PMF of Ion Permeation through alpha7 nAChR development the effective performance in ns/day differs only 30 % 70 on 8 cores. The sharing of the same system bus by multiple cores 60 on one node is presumably the main bottle neck hindering better 50 performance on more nodes. This effect is overlaid by the influ- 40 ence of the inter node communication via Infiniband thus leading 30 to the quasi parallel shift observed for large number of cores in G(z) (kJ/mol) 20 Figure 6. Arminus and Bisgird show a linear speedup on up to 32 10 and 40 cores, respectively. Using more cores leads to a fall-off of 0 performance. It has to be stated that the simulation system used in this studies was not optimized for highly parallel computing. As- -10 pects like the domain decomposition and grid spacing for PME -20 electrostatics have a considerable impact, especially when using -30 larger numbers of cores (Hess et al. 2008). Additional adjustments -40 would enable further scaling on more cores. -4 -2 0 z (nm) 2 4 To achieve an optimal investment of computational time a high Fig. 4. Dependency on pulling velocity. The darkest curve corres- performance per core was anticipated for US simulating. For SMD ponds to a pulling velocity of 0.00375 nm/ns, the velocity of the it was considered more important to minimize the simulated time spring. The speed was doubled with each grayscale. Each individual (wall clock). The relative short US simulations (300 ps) can be curve is based on 6 repetitions with a standard error (omitted for clari- finished using just 4 cores within 7.8 hours, while the longer SMD ty) similar to Figure 5. It can be stated that the minimum around - simulations (2100 ps) need 18 hours on 32 cores, both measured on 2 nm, which is considered to be an artificial deformation, increases Arminius. Considering this discrepancy and comparing to Figure 6 with higher speeds. For the highest speed of 0.03 nm/ns (lightest gray) (lower panel) the relative efficiency (ns per day and core) is the an effect best described as ‘rupture’ can be observed. highest on 2 cores for both clusters. The level remains high on up yields a maximum barrier height of 52.7±11.1 kJ/mol, while the 70 efflux has 57.3±6.9 kJ/mol. Assuming that the permeability corre- lates with exp(-∆Gmax/kBT) the influx is easier by a factor of 6 60 (Hub and de Groot 2008). This does not take the estimated error 50 into account which is of the same order of magnitude. This finding 40 is in very good agreement with experimental findings showing that G(z) (kJ/mol) 30 the nAChR is an rather unselective cation channel (Mishina et al. 1986). The evaluation of potassium, calcium and chlorine ion per- 20 meation barriers is subject of ongoing research. 10 The single 2100 ps pulling simulations took 18 hours each using 0 32 cores. In order to achieve a reasonable error estimate 12 repeti- tions for the influx and 12 for the efflux were prepared. -10 -20 3.3 Performance -30 The two methods US and SMD are carried out on identical simula- 70 tion systems, consisting of exactly the same number of atoms with 60 exactly the same atomic interactions. Although the principals for 50 PMF creation are fundamentally different, in terms of the simula- tion it is the same to constraint the ion to a certain position along 40 the reaction coordinate or the pull it along the same. G(z) (kJ/mol) 30 As stated above 300 separate simulations of 300 ps length each are 20 required for the construction of one US PMF. With three repeti- 10 tions a total of 9360 CPU hours are required. To construct a com- parable SMD pulling PMF twelve 2100 ps simulations with 0 6490 CPU hours are needed. On the first glance the SMD method -10 is one third less expensive than the US method. This conclusion -20 may only be partially true, strongly depending on the architecture of available compute resources. As shown in Figure 6 the simula- -30 -4 -2 0 2 4 tion setup used in this study scales superlinear on up to 8 cores. z (nm) Both HPC clusters are equipped with low latency infiniband inter- Fig. 5. PMFs generated by SMD pulling simulations. The upper panel connects, while Arminius has two single-core CPUs (Intel XeonDP shows the curve for the Na+ influx and the lower panel for the efflux. 3.2 GHz) on each node, Bisgrid has four dual-core CPUs on each The barrier height differs by 4.2 kJ/mol or in terms of permeability a node (AMD Opteron 8220, 2.8 GHz). Although the type of CPU factor of 6. used in both clusters differs by three generations of technological 5 J. Krüger and G. Fels to 8 cores and then drops sharply. The benefits of low latency infi- 4 CONCLUSION niband interconnects are not used to their full extend. This illu- US and SMD based PMFs have been investigated to help under- strates clearly that such kind of simulations with GROMACS standing cation flux of the nAChR ion channel. Both procedures would be best carried out using distributed computing on multi- yield comparable energy barriers for the permeation of Na+ core machines. through the nAChR. The errors estimated for both methods do not 12 differ considerably. From a technical point of view US tends to be more robust, while SMD enables the differentiation of in- and 10 efflux of ions. Both methods differ by 30 % with respect to the 8 computational cost. Due to the higher number of short simulations required for US this computational procedure is more recommend- ns/day 6 ed for low latency distributed computing than SMD. 4 2 ACKNOWLEDGEMENTS We thank the PC2 University of Paderborn for providing computer 0 time. JK gratefully acknowledges a research scholarship from the 0 4 8 12 16 20 24 28 32 36 40 44 48 Alexander von Humboldt-Foundation. Support of the e-Science cores Institute Edinburgh is acknowledged. 50 45 REFERENCES 40 Anishkin, A. and S. Sukharev (2004). Water dynamics and dewetting transitions in the 35 small mechanosensitive channel MscS. Biophys J 86(5): 2883-95. 30 Beckstein, O. and M. S. Sansom (2006). A hydrophobic gate in an ion channel: the speedup 25 closed state of the nicotinic acetylcholine receptor. Phys Biol 3(2): 147-59. 20 Bertrand, D., et al. (1993). Mutations at 2 Distinct Sites within the Channel Domain 15 M2 Alter Calcium Permeability of Neuronal Alpha-7 Nicotinic Receptor. Pro- 10 ceedings of the National Academy of Sciences of the United States of America 5 90(15): 6971-6975. 0 Böckmann, R. A., et al. (2008). Kinetics, statistics, and energetics of lipid membrane 0 4 8 12 16 20 24 28 32 36 40 44 48 electroporation studied by molecular dynamics simulations. Biophysical Journal 95(4): 1837-1850. cores Brejc, K., et al. (2001). Crystal structure of an ACh-binding protein reveals the ligand- binding domain of nicotinic receptors. Nature 411(6835): 269-276. 0,35 Buckingham, S. D., et al. (2009). Nicotinic Acetylcholine Receptor Signalling: Roles in Alzheimer's Disease and Amyloid Neuroprotection. Pharmacological Reviews 0,30 61(1): 39-61. 0,25 Chandrasekhar, I., et al. (2003). A consistent potential energy parameter set for lipids: ns/day and core dipalmitoylphosphatidylcholine as a benchmark of the GROMOS96 45A3 force 0,20 field. Eur Biophys J 32(1): 67-77. 0,15 Fels, G., et al. (1982). Equilibrium Binding of Acetylcholine to the Membrane-Bound Acetylcholine-Receptor. European Journal of Biochemistry 127(1): 31-38. 0,10 Gomaa, M. S., et al. (2007). Homology model of 1 alpha,25-dihydroxyvitamin D-3 0,05 24-hydroxylase cytochrome P450 24A1 (CYP24A1): Active site architecture and ligand binding. Journal of Steroid Biochemistry and Molecular Biology 104(1-2): 0,00 53-60. 0 4 8 12 16 20 24 28 32 36 40 44 48 Gotti, C. and F. Clementi (2004). Neuronal nicotinic receptors: from structure to cores pathology. Progress in Neurobiology 74(6): 363-396. Gullingsrud, J. R., et al. (1999). Reconstructing potentials of mean force through time Fig. 6. Performance (upper panel), speedup (middle panel) and rela- series analysis of steered molecular dynamics simulations. Journal of Computa- tive performance (lower panel) of the transmembrane nAChR system tional Physics 151(1): 190-211. presented in this study. The molecular setup consisting of 71994 atoms Hess, B. (2002). Determining the shear viscosity of model liquids from molecular has been simulated with GROMACS-4.0.4 on two clusters with low dynamics simulations. Journal of Chemical Physics 116(1): 209-217. latency Infiniband interconnects (blue – Bisgrid: 4 dual-core AMD Hess, B., et al. (2008). GROMACS 4: Algorithms for highly efficient, load-balanced, Opteron 8220, 2.8 GHz per node; red - Arminius: 2 single-core Intel and scalable molecular simulation. Journal of Chemical Theory and Computation XeonDP, 3.2 GHz per node). On up to 8 cores superlinear scaling 4(3): 435-447. could be observed for both systems. Up to 40 cores the speedup in- Hub, J. S. and B. L. de Groot (2008). Mechanism of selectivity in aquaporins and creases nearly linear on Arminius, while the performance breaks in at aquaglyceroporins. Proc Natl Acad Sci U S A 105(4): 1198-203. 32 cores on Bisgrid. In order to achieve an optimal investment of Kelly, K. (1999). "3D Bioinformatics and Comparative Protein Modeling in MOE." computational time a high performance per core was used for US Journal of the Chemical Computing Group, from simulating on 4 cores. For SMD it was considered more important to http://www.chemcomp.com/journal/bio1999.htm. minimize the simulated time (wall clock) and 32 cores were used. 6 PMF of Ion Permeation through alpha7 nAChR Krüger, J. and W. B. Fischer (2008). Exploring the conformational space of Vpu from HIV-1: a versatile adaptable protein. J Comput Chem 29(14): 2416-24. Krüger, J. and W. B. Fischer (2009). Structural implications of mutations assessed by molecular dynamics: Vpu(1-32) from HIV-1. Eur Biophys J. Kumar, S., et al. (1992). The Weighted Histogram Analysis Method for Free-Energy Calculations on Biomolecules.1. The Method. Journal of Computational Chemi- stry 13(8): 1011-1021. Kumar, S., et al. (1995). Multidimensional Free-Energy Calculations Using the Weighted Histogram Analysis Method. Journal of Computational Chemistry 16(11): 1339-1350. Lester, H. A., et al. (2004). Cys-loop receptors: new twists and turns. Trends in Neu- rosciences 27(6): 329-336. Mishina, M., et al. (1986). Molecular Distinction between Fetal and Adult Forms of Muscle Acetylcholine-Receptor. Nature 321(6068): 406-411. Miyazawa, A., et al. (2003). Structure and gating mechanism of the acetylcholine receptor pore. Nature 423(6943): 949-955. Patargias, G., et al. (2009). Reconstructing Potentials of Mean Force from Short Steered Molecular Dynamics Simulations of Vpu from HIV-1. J Biomol Struct Dyn 27(1): 1-12. Peng, X., et al. (1994). Human Alpha-7-Acetylcholine Receptor - Cloning of the Alpha-7-Subunit from the Sh-Sy5Y Cell-Line and Determination of Pharmacolog- ical Properties of Native Receptors and Functional Alpha-7-Homomers Expressed in Xenopus-Oocytes. Molecular Pharmacology 45(3): 546-554. Portella, G., et al. (2008). Not only enthalpy: Large entropy contribution to ion per- meation barriers in single-file channels. Biophysical Journal 95(5): 2275-2282. Schuler, L. D., et al. (2001). Molecular dynamics simulation of n-dodecyl phosphate aggregate structures. Eur Biophys J 30(5): 330-43. Unwin, N. (2005). Refined structure of the nicotinic acetylcholine receptor at 4A resolution. J Mol Biol 346(4): 967-89. van der Spoel, D., et al. (2005). GROMACS: fast, flexible, and free. J Comput Chem 26(16): 1701-18. Wallace, I. S. and D. M. Roberts (2004). Homology modeling of representative sub- families of Arabidopsis major intrinsic proteins. Classification based on the aro- matic/arginine selectivity filter. Plant Physiology 135(2): 1059-1068. 7