Transcript Slide 1
Monte Carlo simulation of ionic systems in the presence of dielectric boundaries: application to the selectivity of calcium channels Dezső Boda1,2, Mónika Valiskó2, Bob Eisenberg1, Doug Henderson3, Wolfgang Nonner4, Dirk Gillespie1 1 Rush University Medical Center, Chicago, USA; 2 Pannon University, Veszprém, Hungary; 3 Brigham Young University, Provo, USA; 4 University of Miami School of Medicine, Miami, USA Induced Charge Computation method Introduction Low resolution modeling of inhomogeneous ionic systems usually involves subsystems of different polarizabilities (for example, solution, membrane, protein, electrode, etc.). An electrostatically correct coarse grain modeling of such systems requires representing these subsystems with continuum dielectrics of different dielectric constants, while ions are treated explicitly (as charged hard spheres). The inhomogeneous dielectrics can be described by a position dependent dielectric coefficient e(r). We want to find the solution of Poisson’s equation: e (r) (r) 4(r) Where (r) is the source charge density (the point ions in the center of hard spheres) and f(r) is the electrostatic potential. The boundary condition (the tangential electric field and the normal electric displacement are continuous at the boundary) can be expressed in terms of charges induced on the dielectric boundary h(r). Computer simulation of such systems is problematic because the dielectric boundary forces are not pair-wise additive. Monte Carlo simulations require the calculation of the total electrostatic energy that includes the ioninduced charge interaction in addition to the ion-ion interaction (computed from pair-wise additive Coulombpotentials). Model of the calcium channel e (s) (s s) n(s) e (s) h(s) (s s) e (s) s s 3 ds e (s) E(s) n(s) B where e(s) is the jump in the dielectric coefficient at point s of the dielectric boundary (in the direction of the n(s) normal vector), e(s) is the average of the dielectric coefficients on the two sides, and E(s) is the electric field of the ions at the boundary. After discretization of the surface, this equation is transformed into a matrix equation Ah=c, where matrix A is defined by the expression in square brackets, h is the induced charge, and c is the right hand side that depends on the position of ions. The matrix depends only on the geometry of the dielectrics, so it need be filled and inverted only once at the beginning of the simulation. Moving ions in the simulation changes c, from which h can be computed by a matrix-vector multiplication in every simulation step. An incoming ion induces repulsive charge on the dielectric boundary (cations induce positive charge, anions induce negative charge) if the ion is in the higher dielectric regime (in the solution). In this work, we use the Induced Charge Computation (ICC) method to calculate the induced charges, and we give a biological example: the selectivity of a model calcium (Ca) channel. The induced charge can be calculated from this integral equation (Boda et al. Phys. Rev. E 69, 046702 (2004).): Selectivity of calcium channels Monte Carlo simulations Calcium channels are large membrane proteins that preferentially conduct Ca2+ ions into the cytoplasm (where Ca2+ ions actuate molecular switches) even if Ca2+ is present in micromolar quantity while monovalent ions (Na+ or K+) are present in 100-150 mM concentration. The temperature, the volume of the simulation cell, the number of Na+ ions, and the chemical potential of CaCl2 are the fixed parameters of our mixed thermodynamic ensemble. Selectivity of calcium channels varies in a wide range: from millimolar (the ryanodine receptor Ca channel) to micromolar (the L-type Ca channel). The different Ca2+ affinities of these channels are closely related to their different physiological functions. A: Cross-section of the cylindrical simulation cell. Two “bath” regions are connected via a molecular “pore”, a central cylinder flanked by vestibules. The “protein” forming the pore is a doughnut-shaped body. The protein is embedded in a “membrane” that separates the two baths. The figure is drawn to scale (units in Angstrom). B: Three-dimensional view of the protein surface. The surface grid outlines the “tiles” of the discretized surface used in the electrostatics computations. We study the physical mechanism by which genetic information can control Ca2+ affinity of Ca channels. We show that two physical parameters (the radius R of the selectivity filter and the dielectric coefficient e of the protein), at a first degree, can regulate Ca2+ affinity. Both of these factors (the shape of protein and the local polarizability of the channel wall) are determined by protein structure, e. g., by the genetic code. A common structural motif of Ca channels is four glutamate (Glu) or aspartate (Asp) amino acids in their selectivity filter. We model the negatively charged carboxyl groups of these amino acids by 8 half charged oxygen ions that are confined to the selectivity filter but otherwise free to move within it (see the figure above). Occupancy curves e (r s) Selectivity curves The occupancy of the pore is characterized by the average number of Na+ or Ca2+ ions in the cylindrical section (H=10 A) of the pore. Ca2+ affinity is characterized by the [CaCl2] where the number of Ca2+ and Na+ ions is equal in this region (crosspoint). Charge in the pore As e is decreased, more cations are attracted into the pore by the negative surface charge induced by the oxygens. Lower e favors better local electroneutrality in the pore. As e is decreased, more cations are absorbed in the pore. Ca2+ accumulates in the center, while Na+ also dwells at the pore entrances (R = 4.5 A). As e is decreased, the number of Ca2+ in the pore increases more steeply than the number of Na+. As e is decreased, the crosspoint is shifted towards lower [CaCl2] values, namely, the Ca2+ affinity is improved. As R is decreased, ionic density is increased in the pore. In the case of e=10, the system is trying to maintain charge neutrality, so Na+ is exchanged for Ca2+ as R is decreased. Effect of pore radius Effect of protein dielectric coefficient Concentration profiles The bath contains 100 mM NaCl, while CaCl2 is added to the bath gradually. Low Ca2+ concentrations are simulated in a grand canonical ensemble (insertion/deletion of neutral CaCl2 groups) that is coupled to a biased particle exchange between the channel and the bath to accelerate the convergence of the average number of Na+ and Ca2+ ions in the pore. As R is decreased, the number of Na+ decreases for both e=80 and 10. In the case of e=10, the number of Ca2+ increases with decreasing R. As R is decreased, the crosspoint is shifted towards lower [CaCl2] values, but this effect is stronger for e=10. Summary of the results Ca2+ concentration As R is decreased, increases in the center of the pore. The concentration of Na+ increases at the entrances of the pore. The crosspoint concentration plotted as a function of e and R. Ca2+ affinity is strongly improved by decreasing e and R simultaneously. As e is decreased, the separation of charge on the surface of the protein increases. Total induced charge on the protein surface is zero (Gauss’ law).