GATING CURRENT MODELS COMPUTED WITH CONSISTENT INTERACTIONS

1. BIRS Workshop: Ion Transport: Electrodiffusion, Electrohydrodynamics and Homogenization, May 29-June 3, 2016 GATING CURRENT MODELS COMPUTED WITH CONSISTENT INTERACTIONS Allen Tzyy-Leng Horng1 (洪子倫), Robert S. Eisenberg2, Chun Liu3, Francisco Bezanilla4 1Feng Chia Univ, Taichung, Taiwan, tlhorng123@gmail.com http://newton.math.fcu.edu.tw/~tlhorng 2 Rush Univ., Chicago, IL, USA, 3Penn State Univ., State College, PA, USA, 4Univ. Of Chicago, Chicago, IL, USA

2. ABSTRACTGating currents of the voltage sensor involve back-and forth movements of positively charged arginines through the hydrophobic plug of the gating pore. Transient movements of the permanent charge of the arginines induce structural changes and polarization charge nearby. The moving permanent charge induces current flow everywhere. Everything interacts with everything else in this structural model so everything must interact with everything else in the mathematics, as everything does in the structure. Energy variational methods EnVarA are used to compute gating currents in which all movements of charge and mass satisfy conservation laws of current and mass. Conservation laws are partial differential equations in space and time. Ordinary differential equations cannot capture such interactions with one set of parameters. Indeed, they may inadvertently violate conservation of current. Conservation of current is particularly important since small violations (<0.01%) quickly (microseconds) produce forces that destroy molecules. Our model reproduces signature properties of gating current: (1) equality of ON and OFF charge (2) saturating voltage dependence and (3) many (but not all) details of the shape of charge movement as a function of voltage, time, and solution composition. The model computes gating current flowing in the baths produced by arginines moving in the voltage sensor. The movement of arginines induces current flow everywhere producing ‘capacitive’ pile ups at the ends of the channel. Such pile-ups at charged interfaces are well studied in measurements and theories of physical chemistry but they are not typically included in models of gating current or ion channels. The pile-ups of charge change local electric fields, and they store charge in series with the charge storage of the arginines of the voltage sensor. Implications are being investigated.

3. INTRODUCTION Voltage Sensors Much of biology depends on the voltage across cell membranes. The voltage across the membrane must be sensed before it can be used by proteins. Permanent charges move in the strong electric fields within membranes, so carriers of gating and sensing charge were proposed as voltage sensors even before membrane proteins were known to span lipid membranes [1]. Movement of permanent charges of the voltage sensor is gating current and movement is the voltage sensing mechanism. Measurements Measurement of gating currents (inside the hydrophobic pore) by voltage clamp (with current subtracted from its counterpart of mutated neutral side chain) is possible through conservation law: [1] Hodgkin, A.L. and Huxley, A.F. (1952) “A quantitative description of membrane current and its application to conduction and excitation in nerve”. J. Physiol. 117:500-544.

4. INTRODUCTION, Cont. Structure and sensors Knowledge of membrane protein structure has allowed us to identify and look at the atoms that make up the voltage sensor. Protein structures do not include the membrane potentials and macroscopic concentrations that power gating currents, therefore simulations are needed. Atomic level simulations Molecular (really atomic) dynamics do not provide an easy extension from the atomic time scale 2×10-16 sec to the biological time scale of gating currents that reaches 50×10-3sec. Calculations of gating currents from simulations must average the trajectories (lasting 50×10-3sec sampled every 2×10-16 sec) of ~106atoms all of which interact through the electric field to conserve charge and current, while conserving mass. It is difficult to enforce continuity of current flow in simulations of atomic dynamics because simulations compute only local behavior while continuity of current is global, involving current flow far from the atoms that control the local behavior.

5. Our Modeling Approach for this Inverse Problem (a much more difficult problem due to its dynamics) A hybrid approach is needed, starting with the essential knowledge of structure, but computing only those parts of the structure used by biology to sense voltage. In close packed (‘condensed’) systems like the voltage sensor, or ionic solutions, ‘everything interacts with everything else’ because electric fields are long ranged as well as strong. In ionic solutions, ion channels, even in enzyme active sites, steric interactions are also of great importance that prevent the overfilling of space. Closely packed charged systems are best handled mathematically by variational methods. Variationalmethods guarantee that all variables satisfy all equations (and boundary conditions) at all times and under all conditions. We have then used the energy variational approach developed by Chun Liu [2,3], to derive a consistent 1D reduced model of gating charge movement, based on the basic features of the structure of crystallized channels and voltage sensors. [2] Eisenberg, B., et al. (2010). "Energy Variational Analysis EnVarA of Ions in Water and Channels: Field Theory for Primitive Models of Complex Ionic Fluids." Journal of Chemical Physics 133: 104104  [3] Horng, T.-L., et al. (2012). "PNP Equations with Steric Effects: A Model of Ion Flow through Channels." The Journal of Physical Chemistry B 116(37): 11422-11441.

6. From [8] Jensen, M. O., et al. (2012). “Mechanism of voltage gating in potassium channels.” Science, 336(6078):229-233.

7. One of the 4 subunit VSDs: 10 nonpolar residues from S1, S2 and S3 constructing the hydrophobic pore for arginines(charged) to go through during S4 inward/outward movement (driven by voltage). extracellular (Phe) intracellular From [7] Lacroix, J. J., et al. (2014). “Moving gating charges through the gating pore in a Kv channel voltage sensor.” PNAS, 111:1950-1959.

8. extracellular extracellular depolarization activating ion channel conductance deactivating ion channel conductance polarization intracellular intracellular

9. intracellular extracellular Figure 1. Geometric configuration of the reduced mechanical model including the attachments of arginines to the S4 segment.

10. hydrophobic pore, channel extracellular intracellular vesituble, antechamber, hydrated lumen vesituble, antechamber, hydrated lumen Figure 2: Bead and spring model for argininesand S4. In current continuum model, Ci,i=1, 2, 3, 4 will be represented by concentration distributions not particle, while S4 will be treated as a particle here.

11. Mathematical Description Figure 2. Axisymmetric geometric configuration for mathematical model. Note that pore diameter is 0.3nm (arginine diameter). LR=1.5nm, L=0.7nm. Vestibule radius is R=1nm.

12. The reduced 1D dimensionless PNP-steric equations are expressed as below. The first one is Poisson equation: The second equation is the transport equation based on conservation law: Na and Cl in vesitbules are EDL capacitors. Quick response to voltage change in a very short time compared with slow arginines moving. Quasi-steady assumption for Na and Cl here. steric effect in the form of cross diffusion [3] spring connected to S4 energy barrier in pore

13. Vi , i=1, 2, 3 and 4 being the trap potential for ci representing a spring connecting ci to S4. where K is the spring constant, zi is the fixed anchoring position of spring for cion S4 relative to ZS4(t), which is the center-of-mass z position of S4. Here we set z1=-0.6nm, z2=-0.2nm, z3=0.2nm, z4=0.6nm from the structure showing the arginine anchoring interval on S4 is 0.4nm. ZS4(t) follows the motion of equation based on spring-mass system: overdamped dragging by arginine 0 where mS4, bS4 and KS4 are mass, damping coefficient and restraining spring constant for S4. ZS4,0 is the natural position of ZS4(t) when the net force at the RHS of (6) vanishes.zi,CMis the center of mass for ci ,

14. The free energy barrier Vb , universal to all arginines, in (4) is only non-zero in hydrophobic pore (zone 2), which generally needs to be obtained by potential mean force (PMF) from the structure of protein. Here it is dominated by hydration energy due to hydrophobic pore consisting of surrounding nonpolar residues [4]. Without exact shape of Vb known, here we just assume [4] Zhu, F. and Hummer, G. (2012). “Drying transition in the hydrophobic gate of the GLIC channel blocks ion conduction.” Biophys. J., 103:219-227.

15. Boundary and interface conditions for electric potential ϕ are intracellular continuity of electric potential and displacement (flux) extracellular Boundary and interface conditions for arginine are continuity of concentration and flux Boundary conditions for Na and Cl are Dirichlet and no-flux BCs for Na and Cl Key parameters: Di=50, i=1,2,3,4, K=3, KS4=3, bS4=1.5. best fit with experiment so far with these 3 sensitive parameters.

16. Numerical Method The framework is method of lines (MOL). High-order multi-block pseudospectralmethod is used here to discretize (1) and (2) in space. The resultant semi-discrete system is then a set of coupled ODEs in time, chiefly from (2), and algebraic equations, chiefly from (1) and boundary/interface conditions (9-11). This system is further integrated in time by an ODAE solver (ODE15S in MATLAB) together with the initial condition. High-order pseudospectral method provides good accuracy with economic resolutions. ODE15S is a variable-step-variable-order (VSVO) solver, which is highly efficient in time integration. With these two highly efficient techniques, we can do numerous numerical experiments economically to search parameters that fit experiments best.

17. Results: QV curve (biased) sigmoidal ZS4,0=3.441nm V1/2=-48mV equilibrium state Experimental data: [5] Bezanilla, F., et al. (1994). “Gating of K+ channels : II. The components of gating currents and a model of channel.” Biophy. J., 66:1011-1021.

18. Results: Case 1, V (mV) changed from -90 to -8 and back to -90 holding potential command voltage pulse duration in 140 10 On and off gating currents. I measured right at the middle of channel. Arginine #4 goes first, back last. All arginines move back and forth forming individual currents. area under total current: number of charges moved

19. Results: Case 1, V (mV) changed from -90 to -8 and back to -90

20. Results: Case 1, V (mV) changed from -90 to -8 and back to -90 watch for slight bulge in zone 1 and 3

21. Results: Case 1, V (mV) changed from -90 to -8 and back to -90 S4 movement contributes significantly to the total movement of arginines. Cooperativity of S4.

22. Results: Case 1, V (mV) changed from -90 to -8 and back to -90 Here arginines moved quickly. Equilibrium reached quickly after change of voltage. Very low arginine concentrations inside channel all the time due to energy barrier.

23. Results: Case 1, V (mV) changed from -90 to -8 and back to -90 The trajectory cartoon looks reasonable but somewhat unreal, because arginines jump through hydrophobic pore quickly instead of passing through it slowly.

24. Results: Case 1, V (mV) changed from -90 to -8 and back to -90 The stochastic trajectory of each arginine is calculated based on its probability of existence at zones 1 and 3 to determine using Zi,CMat zone 1 or 3 respectively.

25. Results: Case 2, V (mV) changed from -90 to -50 and back to -90 close to V1/2 fast movement of arginine due to full cooperativity of S4. slow movement due to poor cooperativity of S4 since -50 mV close to V1/2.

26. Results: Case 2, V (mV) changed from -90 to -50 and back to -90

27. Results: Case 2, V (mV) changed from -90 to -50 and back to -90 Will increase to 2 if voltage pulse duration is long enough according to QV curve.

28. Results: Case 2, V (mV) changed from -90 to -50 and back to -90 Still in relaxation, not reaching equilibrium yet.

29. Results: Case 2, V (mV) changed from -90 to -50 and back to -90

30. Results: Family of gating currents for a range of command voltages, V (mV) changed from -90 to -62, …,-8 and back to -90 Area saturates as voltage increases. Longest time span happens (largest τ2) at V=-44 mV.

31. Results: Family of gating currents for a range of command voltages, V (mV) changed from -90 to -62, …,-8 and back to -90 A replot of previous figure.

32. Results: Family of gating currents for a range of command voltages, V (mV) changed from -90 to -62, …,-8 and back to -90. Comparison with experimental data from [5]. Measured right at the middle of channel. Spike is caused by EDL capacitor current at vestibulessubject to sudden voltage change. Not cleanly subtracted in voltage clamp experiment. Pure and non-pure exponential decay.

33. Results: Kinetics, τ2 curve, comparison with experimental data from [5]. Less steepness in the right branch of simulation results compared with experiment is revealed in QV curve as well. Gating on current fit with after the peak.

34. Results: Effect of pulse duration, comparison with experimental data from [5]. Shape saturates as duration increases. Shape saturates as duration increases.

35. Discussions and conclusions • The present 1D continuum model of the voltage sensor is an attempt at capturing the essential structural details that are necessary to reproduce the basic features of experimentally recorded gating currents. After finding appropriate parameters, we found that the general kinetic and steady-state properties are well represented by the simulations. This indicates that this continuum approach, which takes in account all interactions, and satisfies conservation of current, is a good model of voltage sensors. (Pancho is convinced and satisfied by this continuum approach!) • The current result has a qualitatively good agreement with the experiment. Quantitative agreement with experiment is not expected, since (1) exact PMF is not known, (2) the present model is a single-channel model while experimental recordings were from ensemble average of many ion channels (and VSDs).

36. Discussions and conclusions, cont. • The pore hydration energy (not arginine solvation energy) dominates the energy barrier Vb(z) in the present model [4]. Here, Vb(z) is independent of time. However, once the first arginine enters the hydrophobic pore by carrying some water with it cooperatively. This partial wetting of pore will lower Vband enable the next arginine to enter the pore with less difficulty. This will explain the cooperativitybetween arginineswhen they jump through the pore. This remains to be fixed in the future modeling. • It has been reported that a very strong electric field might affect the hydration equilibrium of the hydrophobic pore and would lower its hydration energy barrier [6]. This may help explain the quick saturation in the upper right branch of QV curve (and smaller τ2). This also remains to be fixed in the future modeling. [6] Vaitheeswaran, S., et al. (2004). “Electric field and temperature effects on water in the narrow nonpolar pores of carbon nanotubes.” J. Chem. Phys., 121(16):7955-7965.

37. Discussions and conclusions, cont. • The next step is to model the details of interactions of the moving arginines with the wall of the channel. There is plenty of detailed information on the amino acid side chains in the channel and how each one of them have important effects in the kinetics and steady-state properties of gating charge movement [7]. The experimentally studied side chains reveal steric as well as dielectric interactions with the arginines that the present 1D model over-simplifies by just lumping in Vb. This can only be possibly fixed by a 3D model in the future. • Further work must address the mechanism of coupling between the voltage sensor movements and the ion conduction pore. It seems likely that the classical mechanical models will need to be extended to include coupling through the electrical field. It is possible that the voltage sensor modifies the stability of conduction current (consolidation at activating and deactivating states).

38. References [1] Hodgkin, A.L. and Huxley, A.F. (1952) “A quantitative description of membrane current and its application to conduction and excitation in nerve”. J. Physiol. 117:500-544. [2] Eisenberg, B., et al. (2010). "Energy Variational Analysis EnVarA of Ions in Water and Channels: Field Theory for Primitive Models of Complex Ionic Fluids." Journal of Chemical Physics 133: 104104  [3] Horng, T.-L., et al. (2012). "PNP Equations with Steric Effects: A Model of Ion Flow through Channels." The Journal of Physical Chemistry B 116(37): 11422-11441. [4] Zhu, F. and Hummer, G. (2012). “Drying transition in the hydrophobic gate of the GLIC channel blocks ion conduction.” Biophys. J., 103:219-227. [5] Bezanilla, F., et al. (1994). “Gating of K+ channels : II. The components of gating currents and a model of channel.” Biophy. J., 66:1011-1021. [6] Vaitheeswaran, S., et al. (2004). “Electric field and temperature effects on water in the narrow nonpolar pores of carbon nanotubes.” J. Chem. Phys., 121(16):7955-7965. [7] Lacroix, J. J., et al. (2014). “Moving gating charges through the gating pore in a Kv channel voltage sensor.” PNAS, 111:1950-1959. [8] Jensen, M. O., et al. (2012). “Mechanism of voltage gating in potassium channels.” Science, 336(6078):229-233.

39. Besides thanking the long-term collaborators in US, thanks for the support of our group members: Tai-Chia Lin (NTU), Ren-Hsiang Chen (THU), Chun-HaoTeng (NCHU), Son-Reng Hwang (NTUST), Chiun-Jiang Lee (NHCNEU), Former post-doc (Yi-Ping Lo) and many graduate students. About PNP/PB, besides the present gating current project, I have done Selectivity of Ca channel by PNP-Steric model. (with Eisenberg, Liu and Lin) Nanofiltration in desalination of sea water with nanopore membrane: over-lapping EDL modeled by Navier-Stokes equations + PNP-steric + dielectric exclusion (supported by ITRI). Supercapicitor made of carbon with numerous pores of nanosize: PNP-steric model + dielectric exclusion. (with HsisengTengat NCKU) and is now working on (4) 3D simulation of KcsA (supported by MOST). (5) Exploration of gA by its 3D simulations (with ShixinXuat UC Riverside). Look forwards to have collaborations with people of this workshop in the future.