A thermodynamic description for physiological transmembrane transport

Work required for transmembrane molecular fluxes

Consider a system consisting of a biological membrane surrounded by two aqueous compartments (e.g. extracellular and intracellular). Assume, to start with, that the compartments contain molecules of a single type s (e.g. Na⁺, K⁺, glucose), possibly in different concentrations. Let ΔG_s be the energy required for the transport of the molecules across the membrane in a specific direction (e.g. inside to outside). To write an expression for ΔG_s it is necessary to take the direction of motion of the s-molecules into account. To do so, label the extracellular and intracellular compartments as 0 and 1, respectively, and let c_s and d_s represent the source and the destination compartments for the transport of the s-molecules. The pair (c_s, d_s)=(0,1) represents inward transport and the pair (c_s, d_s)=(1,0) represents outward transport. The work required to transport n_s molecules of type s from compartment c_s to compartment d_s can then be written as

(Aidley, 1998; Blaustein et al., 2004; De Weer et al., 1988) where q, z_s, [s]₀, and [s]₁ represent the elementary charge, the valence, the extracellular, and the intracellular concentrations for the molecules of type s, respectively. Two particular cases are worth noticing. First, if s is an ion, then z_s ≠ 0 and Equation (1) becomes

where v_s is the Nernst potential for the s-molecules¹ (Nernst, 1888). Second, if the s-molecules are not charged, then z_s = 0 and the work required to move the s-type molecules from c_s to d_s simplifies to

If ΔG_s < 0, then the molecules can be transported passively (e.g. electrodiffusion), decreasing the electrochemical gradient for s across the membrane. In contrast, if ΔG_s > 0, the transmembrane transport of s from c_s to d_s is not thermodynamically favorable, which means the transport from c_s to d_s requires energy that is not available in the electrochemical gradient for s (active transport). As a consequence, active transport of s would increase the driving force for the motion of s across the membrane.

Joint transmembrane transport of different types of molecules

To find an expression for ΔG that describes a more general transport mechanism, assume that transport takes place as single events in which molecules of m different types move in parallel, or possibly sequentially (e.g. first Na⁺, then K⁺), across the membrane. Let S be a set that represents the types of molecules that are jointly transported in a single event. For instance, for Na⁺-H⁺ exchangers, S = {N a, H}, with m = 2. The energy required to transport these molecules is the sum of the energies required to transport each of the molecules in S. In other words,

As before, transport is thermodynamically favorable when ΔG_s ≤ 0. If not, extra energy is required. To distinguish between these two cases, define the total energy of the transport mechanism as

where δ_Ext = 1 if ΔG_s > 0, and 0 otherwise. In particular, for ATP-driven transport, the extra energy supplied by hydrolysis of ATP (De Weer et al., 1988; Tanford, 1981) is

where v_ATP ≈ −450 mV (Endresen et al., 2000), but could vary depending on the amounts of ATP, ADP, and P_i (De Weer et al., 1988). Similar expressions could be derived for active transport driven by light, or other sources of energy.

Flux due to transmembrane transport

The formulation in Equation (10) can be combined with Equation (1) to derive a generalized expression for flux and model different known mechanisms of physiological transmembrane transport, possibly combining the transport of different molecules simultaneously (e.g. Na-H exchange). In this case, the forward direction of the transport would be described by the combined forward transport of each of the different molecules under consideration. For instance, the source and target compartments for Na⁺ and Ca²⁺ are different in Na-Ca exchangers. The stoichiometry for the transport mediated by Na-Ca exchangers in the forward direction involves three Na⁺ molecules moving inward (along their electrochemical gradient) in exchange for one Ca²⁺ molecule moving outward (against their electrochemical gradient) (Mullins, 1979; Venetucci et al., 2007).

Let α and β be the flux rates in the forward and backward directions, in units of molecules per ms per µm⁻². These rates depend, a priori, on the energy required for the transport of the molecules in S. The net flux rate associated to the net transmembrane transport, can then be written as

How do α and β depend on ΔG? The steady state relationship between the energy ΔG and the the forward and backward flow rates, hereby represented by α and β, can be written as

where k is Boltzmann’s constant, and T the absolute temperature.

Assuming that α and β are continuous functions, the rates α and β can be rewritten as

where r (molecules per ms per µm⁻²) may depend on temperature (Sen & Widdas, 1962), the transmembrane potential (Starace et al., 1997), the concentrations inside and outside of the membrane (Yue et al., 1990), and other factors (Novák & Tyson, 2008). Note that the functional form of the rates in Equations (9) are similar to those by Butler (1924); Erdey-Grúz & Volmer (1930). Also, notice that the steady state relationship between α and β in Equation (8) can be obtained from Equations (9), for any r and any b. However, it should be the case that r and b vary in specific ranges depending on the physico-chemical characteristics of the pore through which molecules cross the membrane, and in general, on the transport mechanism. As mentioned earlier, the rate r should be larger for electrodiffusive transport in comparison to the slower transport rates for pumps and other carrier proteins. If the parameter b ∈ [0, 1], then bΔG and (b−1)ΔG have opposite signs and can be thought of as the energies required to the transport of the molecules in S in the forward and backward directions, respectively, with b biasing the transport in the forward direction when close to 1, and in the backward direction when close to 0 (Figure 1).

Figure 1. Fluxes biased in the target→source (backward, b=0.1, black dashed line), source→target (forward, b=0.9, red dash-dot line), or showing no rectification (b=0.5, blue solid line).

See Table 1 for examples.

The flux can then be written explicitly combining Equation (7) and Equations (9) to obtain,

Taking the above observations into account, it is possible to combine Equation (4) and Equation (5), to write an expression similar to Equation (8) for the steady state balance between the forward and backward transport of all the molecules in S.

Flux and current

Substitution of the formulas for ΔG from Equation (4) Equation (5) into Equation (10), the flux rate resulting from simultaneously transporting molecules in S across the membrane can be written explicitly

where v_T = kT/q and

represents the net number of charges moved across the membrane.

If the transport is electrogenic, then the product qη (in Coulombs) represents the net charge moved across the membrane, relative to the extracellular compartment. Non electrogenic transport yields η = 0, which means the flow does not depend on the transmembrane potential, and

If only ions are involved in the transport, the flux simplifies to

where

The quantity v_o/η can be thought of as a reversal potential. If η < 0, then positive charge is transported inward, or negative charge is transported outward. In contrast, η > 0 means that positive charge is transported outward or negative charge transported inward. For instance, inward electrodiffusion of single Na⁺ ions gives an η = −1, which can be thought of as loosing one positive charge from the extracellular compartment in each transport event (see Table 1). In particular, for electrodiffusive (passive) transport of ions of type l, v_o reduces to n_lz_l (c_l − d_l) v_l. A list with examples of energies and total charge movements for different transport mechanisms can be found in Table 1.

The first, more complex, form of the flux in Equation (13) could be useful when working with models for which changes in the concentrations of different molecules are relevant.

Transmembrane current. The flux that results in electrogenic transport (Equation (13) and Equation (14)) can be converted to current density after multiplication by qη. In short form,

with qr in amperes/m² or equivalent units.

Substitution of Equation (13) or Equation (14) into Equation (16) yields a general formula for the current generated by transmembrane ionic flux (Figure 2), that uses the same functional form for channels (protein or lipid) and pumps. Recall that Equation (16) can also be written explicitly in terms of the transmembrane concentrations of one or more of the ions involved using Equation (13). It is possible to derive expressions for r that take into account biophysical variables like temperature and the shape and length of the pore through which the molecules cross (Endresen et al., 2000; Pickard, 1969).

Figure 2. Fluxes for K-electrodiffusion for b_K ∈ {0.1, 0.5, 0.9} and their cubic approximations.

Inward rectification occurs for b_K < 1/2 and outward rectification for b_K > 1/2 and qr_K N_K = 1.

Special cases and examples

A number of nontrivial and important properties of transmembrane ionic currents, including rectification, are also described by Equation (16). Also, different models for current already in the literature can be obtained by making approximations or setting particular cases from Equation (16). Examples include electrodiffusive currents that result from integration of the Nernst-Planck equation along the length of membrane pore (Jacquez & Schultz, 1974; Johnston et al., 1995; Pickard, 1969). Of particular interest, conductance-based currents are linear approximations of the formulation (16), around the reversal potential for the current.

Lower order approximations to the generic formulation and conductance based models. Conductance-based currents (Hodgkin & Huxley, 1952) are linear approximations of the generic current from Equation (16), around the reversal potential v_o/η. To see this, use Taylor’s theorem (Courant & John, 2012; Spivak, 2018) to rewrite the generic current from Equation (16) as a series around v_o

Truncation of the series to first order gives

where g = η²qr/v_T is in units of nS/µm², which has the functional form of the conductance-based current used in the Hodgkin & Huxley (1952) model. For instance, the linear approximation for the current through an open sodium channels around v_Na in Equation (18) gives g_Na = qr_Na/v_T, and v_o = η_Nav_Na, with η_Na = −1, so that i_Na ≈ g_Na(v − v_Na).

Notice that third order approximations to Equation (14) can also capture rectification. In contrast, first order approximations (conductance-based models) cannot capture rectification.

Rectification results from asymmetric bidirectional flow. The flux of molecules across the membrane can be biased in either the outward or the inward direction when mediated by proteins. This was first called ”anomalous rectification” by Katz (1949), who noticed that K⁺ flows through muscle membranes more easily in the inward, than in the outward direction (Adrian, 1969; Armstrong & Binstock, 1965). It was later found that some K⁺ channels display the bias in the opposite direction (Woodbury, 1971). The former type of K⁺ current rectification is called inward, and the latter outward.

Rectification is a bias in either of the two directions of transport, which may result from changes in the structure of the proteins or pores through which the molecules cross the membrane (Hollmann et al., 1991; Riedelsberger et al., 2015). The type of rectification (inward or outward) depends on what molecules are being transported and on the structure of the proteins mediating the transport. Rectification is therefore not only a property of ions, as shown by molecules like glucose, which may cross membranes via GLUT transporters bidirectionally, but asymmetrically, even when the glucose concentration is balanced across the membrane (Lowe & Walmsley, 1986).

Rectification can be described by Equation (13) by setting b to values different from 1/2, and becomes more pronounced as b is closer to either 0 or 1. These values represent biases in the transport toward the source, or the target compartment, respectively. As a consequence, rectification yields an asymmetry in the graph of α − β as a function of ΔG (Figure 2). For electrogenic transport, rectification can be thought of as an asymmetric relationship between current flow and voltage, with respect to the reversal potential v_o. The particular case b = 1/2 (non-rectifying) yields a functional form for current similar to that proposed by Pickard (1969), and later reproduced by (Endresen et al., 2000), namely

From here on, subscripts will be used to represent different transport mechanisms. For instance, the current for a Na-Ca pump will be written as i_NaCa.

Electrodiffusion of K⁺ through channels (η = 1 and v₀ = v_K), is outward for v > v_K, and inward for v < v_K. The K⁺ current through the open pore is therefore

Current flow through inward rectifier channels (Riedelsberger et al., 2015) can be fit to values of b_K < 1/2. For instance,

describes a current with limited flow of K⁺ in the outward direction, similar to the currents described originally by Katz (1949). Analogously, b_K > 1=2 limits the inward flow. For example, the current

describes outward rectification (Riedelsberger et al., 2015).

Based on the work of Riedelsberger et al. (2015) on K⁺ channels, inward (outward) rectification arises when the S4 segment in K⁺ channels is located in the inner (outer) portion of the membrane. These two generic configurations can be thought of in terms of ranges for the parameter b_K, namely, b_K < 1/2 for inward, and b_K > 1/2 for outward rectification.

In general, ion channels are typically formed by different subunits, that may combine in different ways, resulting in structural changes that may restrict the flow of ions through them, causing rectification. For instance, non-NMDA glutamatergic receptors that can be activated by kainic acid and α-amino-3-hydroxy-5-methyl-4-isoxazole propionic acid (AMPA) conduct Na⁺, K⁺, and Ca²⁺, with different permeabilities depending on the subunits that form the receptor (Hollmann et al., 1991). The reason is that the specific combination of GluR subunits forming the receptor restrict ionic flow in different ways. In particular, the currents recorded in oocytes injected with combinations of GluR1 and GluR3 cRNA have different steady state amplitudes and show different levels of rectification (Figure 3).

Figure 3. Currents recorded from oocytes injected with GluR3 cRNA (blue), or a combination of GluR1 and GluR3 cRNA (orange), in Ca²⁺ ringer solution after activation by AMPA.

The curves were fit with (v_o, b, rq)=(-30,0.45,21), for GluR3 and (v_o, b, rq)=(-35,0.35,20), for GluR1+GluR3. The data was digitized from Figure 3B in the article by Hollmann et al. (1991), using the ginput function from the python module matplotlib (Hunter, 2007) - see Supplementary File 2.

Primary active transport. The Na-K ATPase is a primary active transporter that uses the energy from the hydrolysis of one molecule of ATP for the uphill transport of Na⁺ and K⁺ (De Weer et al., 1988). The kinetics of the Na-K ATPase can be assumed to translocate 3 Na⁺ ions outward and 2 K⁺ ions inward (η_NaK = 1) with a reversal potential ν_NaK = v_ATP + 3ν_Na − 2v_K (see Table 1) in a single transport event (Chapman, 1973; Garrahan & Glynn, 1967; Gadsby et al., 1985; Post & Jolly, 1957). Importantly, the transport kinetics of the Na-K ATPase and by extension, the current, reverse for potentials smaller than v_NaK (De Weer et al., 1988).

The current-voltage relationships recorded from Na-K ATPases in guinea pig ventricular cells are shaped as hyperbolic sines (Gadsby et al., 1985). Those currents would be fit with b_NaK ≈1/2, yielding currents of the form

The voltage-dependence of the Na-K ATPase currents is reported to show a plateau as v increases past the reversal potential for the current, in response to steroids like strophandin (Nakao & Gadsby, 1989). In such cases, the Na-K ATPase current can be assumed to be inwardly rectifying and fit with values of b_NaK ≈ 0, so that,

or alternatively,

The rectification for the Na-K pump ATPase has also been reported to occur in small neurons of the dorsal root ganglion in rats (Hamada et al., 2003). The alternative expression (25) also explains qualitatively different behaviors of the Na-K current as a function of the transmembrane concentrations of Na⁺ and K⁺. For instance, if either [Na]₁ or [K]₀ increase and v > v_NaK, then the amplitude of i_NaK would increase at a smaller rate of change in comparison to when v < v_NaK, which grows exponentially in size. This is also in line with reports of non significant changes in the transport by Na-K ATPases in response to elevated intracellular Na⁺ during heart failure (Despa et al., 2002), in which the transmembrane potential is likely to be depolarized.

Secondary active transport. An example of a pump that mediates secondary active transport is the Na-Ca exchanger, which takes 3 Na⁺ ions from the extracellular compartment in exchange for one intracellular Ca²⁺ ion in forward mode (Pitts, 1979; Reeves & Hale, 1984). The reversal potential for the current is v_NaCa = 2v_Ca − 3v_Na, with η_NaCa = 1. Assuming b_NaCa = 1/2, the Na-Ca current is

The driving force v − v_NaCa could reverse in sign with large enough increases in the intracellular concentration of Ca²⁺, or in the membrane potential. As a result, the current could have a dual contribution to the change in transmembrane potential, as predicted by some theoretical models of cardiac pacemaker activity (Rasmusson et al., 1990a; Rasmusson et al., 1990b).

Electrodiffusive transport. Consider transmembrane electrodiffusive transport of a single ionic type x, with z_x and v_x representing the valence and the Nernst potential for x-ions, respectively. In this case, the reversal potential satisfies

and the generic expression (16) can be rewritten as

In the absence of rectification (b_x = 0.5),

For calcium channels,

See Jacquez & Schultz (1974); Pickard (1969) Table 1 for other examples.

The applicability of the general formulations described above is illustrated next with models of cardiac and neuronal membrane potential.

Cardiac pacemaking in the sinoatrial node

The pacemaking dynamics of cells in the rabbit sinoanatrial node (Figure 4) can be modeled using low dimensional dynamical systems based on the assumption that v changes as a function of K⁺, Ca²⁺, and some Na⁺ transmembrane transport, (Herrera-Valdez & Lega, 2011; Herrera-Valdez, 2014). Transmembrane currents are assumed to be mediated by a combination of channel-mediated electrodiffusion and pumping mechanisms. Explicitly, Ca²⁺ transport is mediated by L-type Ca_v13 channels (Mangoni et al., 2003) and Na⁺-Ca²⁺ exchangers (Sanders et al., 2006). K⁺ transport is mediated by delayed-rectifier voltage-activated channels (Shibasaki, 1987), and Na⁺-K⁺ ATPases (Herrera-Valdez & Lega, 2011; Herrera-Valdez, 2014). In this model, the activation for the L-type Ca²⁺ channels is fast, and assumed to be at steady state (Herrera-Valdez & Lega, 2011). The proportion of activated K⁺ channels and the proportion of inactivated Ca²⁺ channels are both represented by a variable w (Av-Ron et al., 1991; Herrera-Valdez & Lega, 2011). The activation phase of currents recorded in voltage-clamp experiments often displays sigmoidal time courses (Covarrubias et al., 1991; Hodgkin & Huxley, 1952; Tsunoda & Salkoff, 1995). Therefore, the activation dynamics represented by w are described by solutions to equations of the form

Figure 4. Central sinoatrial node pacemaking dynamics using the system (39)-(34).

A. Transmembrane potential and the reversal potential v_NaCa as a function of time. B,C. Dynamics of large currents and small currents, respectively.

where F_w and R_w represent the voltage-dependent steady state and rate (1/ms) for the opening of Kd channels (Willms et al., 1999). The steady state for the activation of voltage-dependent channels is described by the function

which has a graph with increasing sigmoidal shape as a function of v. The parameters g_u and v_u control the steepness and the half-activation potential, for u ∈ {m, w}. The activation rate for K⁺ channels is a voltage-dependent function of the form

where b_w represents a bias in the conformational change for activation. The function R_w has the shape of a hyperbolic cosine when b_w is 1/2. The transmembrane (normalized) currents are given by

where c represents the intracellular Ca²⁺ concentration, and φ_x is a difference of exponential functions as defined above, with x ∈ {NaK, NaCa, KD, CaL}. The temporal evolution for v is then described by

During pacemaking the concentrations of Na⁺ and K⁺ across the membrane are assumed to change negligibly, but the Ca²⁺ concentration does change at least 10-fold (Herrera-Valdez & Lega, 2011; Rasmusson et al., 1990a; Rasmusson et al., 1990b). Therefore, the system includes an equation for the dynamics for c in which c converges to a steady state c_∞ in the absence of Ca²⁺ fluxes, and increases proportionally to the total transport of Ca²⁺ ions via L-type channels and Na⁺-Ca²⁺ exchangers (Figure 4). Explicitly,

where k_c (µM/mV) represents the impact of the transmembrane Ca²⁺ fluxes on the free intracellular Ca²⁺ concentration. The minus sign in front of k_c accounts for the fact that the sign of the J_CaL is negative. The sign in front of J_NaCa is because the forward flux of Ca²⁺ mediated by the Na-Ca exchanger is opposite to that of electrodiffusive Ca²⁺.

Figure 5. Dynamics of the calcium current and double activation.

A. Behavior of the inactivating L-type Ca²⁺ current with respect to the transmembrane potential (blue line) and a non-inactivating current (gray line), and (B) with respect to the time-dependent change in v.

The solutions of Equation (32)–Equation (40) with parameters as in Table 2 reproduce important features of the membrane dynamics observed in the rabbit’s central sinoatrial node, including the period (ca. 400 ms), amplitude (ca. 70 mV), and maximum ∂_t v (<10 V/s) of the action potentials (Zhang et al., 2000).

The solutions of the system show a number of interesting features related to ionic fluxes. First, the Na-Ca current reverses when v = v_NaCa (Figure 4A, blue line). During the initial depolarization and until the maximum downstroke rate, approximately, v_NaCa < v, which means J_NaCa > 0, so that Ca²⁺ extrusion by the Na-Ca exchanger occurs only for a brief period of time during the downstroke and also after each action potential (Figure 4C, blue line). Second, as previously reported in different studies involving spiking dynamics, the time course of the Ca²⁺ current shows a partial inactivation with a double peak (Figure 4B, blue line) around a local minimum (Carter & Bean, 2009; Rasmusson et al., 1990a; Rasmusson et al., 1990b), and in agreement with data from voltage-clamp experiments (Mangoni et al., 2006). A number of models have made attempts to reproduce the double activation by making extra assumptions about gating (Rasmusson et al., 1990a; Rasmusson et al., 1990b). For instance, some models include a second activation variable, or the multiple terms in the steady state gating, or in the time constant for activation or inactivation. However, the explanation for the double peak can be much simpler. The calcium current J_CaL is a negative-valued, non monotonic function for v < v_Ca, which can be thought of as a product of a amplitude term that includes gating and the function φ_CaL. The normalized current J_CaL has a local minimum (maximum current amplitude) around -10 mV (Figure 4B, blue line and Figure 5A, blue line), after which the current decreases, reaching a local maximum as the total current passes through zero, at the peak of the action potential around 10 mV, (Figure 4B, where ∂_t v = 0). The first peak for the Ca²⁺ current occurs when v reaches the maximum depolarization rate (Figure 5B). As v increases (e.g. upstroke of the action potential). The second peak for the current occurs as the membrane potential decreases, and passes again through the region where the maximal current occurs (local minimum for J_Ca). The two local minima for J_CaL represent peaks in the Ca²⁺ current that have different amplitudes due the difference in the time course of v during the upstroke and the downstroke of the action potential (Figure 5A, blue line, and Figure 5B, where ∂_t v = 0). It is important to remark that the dual role played by w is not the cause of the double activation. This is illustrated by analyzing the behavior of a non-inactivating J_CaL without the inactivation component, (Figure 5A, gray line). The double activation can also be observed in models in which the activation of K⁺ channels and the inactivation of Ca²⁺ or Na⁺ channels are represented by different variables (Rasmusson et al., 1990a) and in dynamic voltage clamp experiments on neurons in which there are transient and persistent sodium channels (Carter & Bean, 2009).

The double peak in the Ca²⁺ current is reflected in the intracellular Ca²⁺ concentration (Figure 6, gray line), and by extension, on the Nernst potential for Ca²⁺ (Figure 6, blue line), which display two increasing phases and two decreasing phases, respectively. The first and faster phase in both cases occur during the initial activation of the L-type channels. The second phase occurs during the downstroke, as second peak of the Ca²⁺ current occurs. As a consequence, the reversal potential for the Na-Ca exchanger, v_NaCa = 3v_Na − 2v_Ca (Figure 6, orange line) also has two phases, this time increasing. Increasing the intracellular Ca²⁺ (Figure 6, gray line) concentration decreases the Nernst potential for Ca²⁺, and viceversa. By extension, the reversal potential for the Na-Ca exchanger, v_NaCa = 3v_Na − 2v_Ca becomes larger when c increases. Ca²⁺ enters the cell in exchange for Na⁺ that moves out when v > v_NaCa, during most of the increasing phase and the initial depolarization phase of the action potential (blue lines in Figure 4A and C, and Figure 6).

Figure 6. Calcium dynamics during pacemaking.

Time courses of the intracellular calcium concentration (gray, left axis), the Nernst potential for Ca²⁺ (orange, right axis), and the reversal potential for the Na-Ca exchanger (blue, right axis).