2.1 Mathematical description of transport processes

A theoretical/mathematical description of membrane transport requires to quantify the transport processes involved. For this purpose, the biophysical knowledge on transport (Hille, Reference Hille2001) needs to be translated into the universal language of mathematics. Here, only basic knowledge on the membrane transport is often sufficient to obtain rather simple but very powerful equations. The art of mathematical modelling lies in the elimination of redundancies, which means that many physiological parameters can be combined into a few, non-redundant model parameters. As this manuscript will illustrate, such steps simplify both the equations and the computational handling of the model. However, in order to draw conclusions from the model calculations, the model output must be translated from mathematical language to biology. The challenge in understanding and the art of physiological interpretation lies particularly in the fact that independent physiological parameters in biology can be redundant in mathematics. Having this in mind, we may now tackle the question:

Steady state corresponds to the state in which the system does not change any longer. Indeed, homeostatic conditions are steady state conditions because the system exhibits a certain stability and maintains its parameters constant. The homeostatic properties of a system can be analysed in a two-step process. First, the steady state is analysed and then, in the second step, the system is destabilised in order to investigate its resilience. In this study, we begin with the first step and will present the second in a follow-up paper.

Figure 2. Modelling of fluxes through a transporter in a steady state. (a) The flux curve of a transporter J_X (V) (grey line) is zero at $V={E}_X^{ss}$ ( ${J}_X\left({E}_X^{ss}\right)=0$ , blue point). (b) In a steady state, flux and voltage are constant, and the curve is represented by the red point $\left({V}_{ss}|{J}_X^{ss}\right)$ . The abscissa of this point is V_ss , while the ordinate, ${J}_X^{ss}$ , can be expressed as m ⋅ ( ${V}_{ss}-{E}_X^{ss}$ ) using triangulation (m = tan(α)). (c) Changes in the activity are mirrored by changes in the slope m and the angle α. An increase in the transporter activity (dashed grey line), for example, by a higher expression or activation by phosphorylation, results in a larger ${J}_X^{ss}$ -value (purple point), which is represented by a larger angle (α1) and hence a larger slope (m1). A decrease in transporter activity (dotted grey line), on the other hand, results in a smaller ${J}_X^{ss}$ -value (pink point), smaller angle (α2) and smaller slope (m2). Thus, the slope m contains the regulatory features of the transporter.

In the beginning, we consider passive or secondary active transporters, that is, membrane proteins that do not couple ATP hydrolysis to transport activity (pumps will be considered later). The flux through a passive transporter (J_X ) depends on the electrochemical gradient ( $\Delta \mu$ ) of the transported solute across the membrane (Kedem & Katchalsky, Reference Kedem and Katchalsky1963). For a permeating substrate X, this gradient is given by:

(1a) $$\begin{align}{\Delta \mu}_X&={\mu}_X^{\mathrm{in}}-{\mu}_X^{\mathrm{out}}= RT\cdot \ln \left({\left[X\right]}_{\mathrm{in}}\right)+{z}_XF{\psi}_{\mathrm{in}}- RT\cdot \ln \left({\left[X\right]}_{\mathrm{out}}\right)\nonumber\\&\quad-{z}_XF{\psi}_{\mathrm{out}}=F\cdot \left({z}_XV-\frac{RT}{F}\cdot \ln \left(\frac{{\left[X\right]}_{\mathrm{out}}}{{\left[X\right]}_{\mathrm{in}}}\right)\right),\end{align}$$

with the gas constant R, the Faraday constant F, the absolute Temperature T, the internal and external electrical potentials ${\psi}_{\mathrm{in}}$ and ${\psi}_{\mathrm{out}}$ ( ${\psi}_{\mathrm{in}}-{\psi}_{\mathrm{out}}=V$ , the membrane voltage), the valence z_X and the concentrations [X]_in and [X]_out, of the substrate X, respectively (Carpaneto et al., Reference Carpaneto, Geiger, Bamberg, Sauer, Fromm and Hedrich2005; Gerber et al., Reference Gerber, Fröhlich, Lichtenberg-Fraté, Shabala, Shabala and Klipp2016; Nour-Eldin et al., Reference Nour-Eldin, Andersen, Burow, Madsen, Jørgensen, Olsen, Dreyer, Hedrich, Geiger and Halkier2012). The net flux of X through the transporter (J_X ) is a function of the gradient Δμ_X , that is, the concentrations and, if ${z}_X\ne 0$ , also the membrane voltage. If the gradient is zero, $\Delta \mu =0$ , there is no net flux, that is, J_X = 0. Thus, if J_X depends on the membrane voltage, it is zero at $V={E}_X$ , with ${E}_X=\frac{RT}{z_XF}\cdot \ln \left(\frac{{\left[X\right]}_{\mathrm{out}}}{{\left[X\right]}_{\mathrm{in}}}\right)$ , also known as the Nernst (or zero-flux) potential for X. Otherwise, it is zero if ${\left[X\right]}_{\mathrm{out}}={\left[X\right]}_{\mathrm{in}}$ .

In the case of a secondary active co-transporter that transports n_X molecules/ions of type X (valence z_X ) along with n_Y molecules/ions of type Y (valence z_Y ), the electrochemical gradient is given by:

(1b) $$\begin{align}&{\Delta \mu}_{XYs}={n}_x\cdot \Delta {\mu}_X+{n}_Y\cdot \Delta {\mu}_Y\\&\quad=F\cdot \left[\left({n}_X{z}_X+{n}_Y{z}_Y\right)\cdot V-{n}_X\frac{RT}{F}\cdot \ln \left(\frac{{\left[X\right]}_{\mathrm{out}}}{{\left[X\right]}_{\mathrm{in}}}\kern-1.2pt\right)-{n}_Y\frac{RT}{F}\cdot \ln \left(\frac{{\left[Y\right]}_{\mathrm{out}}}{{\left[Y\right]}_{\mathrm{in}}}\kern-1.2pt\right)\kern-1.2pt\right],\nonumber\end{align}$$

for a symporter, and by

(1c) $$\begin{align}&{\Delta \mu}_{XYa}={n}_x\cdot \Delta {\mu}_X-{n}_Y\cdot \Delta {\mu}_Y\\&\quad=F\cdot \left[\left({n}_X{z}_X-{n}_Y{z}_Y\right)\cdot V-{n}_X\frac{RT}{F}\cdot \ln \left(\frac{{\left[X\right]}_{\mathrm{out}}}{{\left[X\right]}_{\mathrm{in}}}\right)+{n}_Y\frac{RT}{F}\cdot \ln \left(\frac{{\left[Y\right]}_{\mathrm{out}}}{{\left[Y\right]}_{\mathrm{in}}}\right)\right],\nonumber\end{align}$$

for an antiporter. In the case of symporters, the flux J_XY is zero at $V={E}_{XYs}$ , with ${E}_{XYs}=\frac{n_X}{n_X{z}_X+{n}_Y{z}_Y}\cdot {E}_X+\frac{n_Y}{n_X{z}_X+{n}_Y{z}_Y}\cdot {E}_Y$ . In the case of antiporters, it is zero at $V={E}_{XYa}$ , with ${E}_{XYa}=\frac{n_X}{n_X{z}_X-{n}_Y{z}_Y}\cdot {E}_X-\frac{n_Y}{n_X{z}_X-{n}_Y{z}_Y}\cdot {E}_Y$ , if ${n}_X{z}_X-{n}_Y{z}_Y\ne 0$ ; otherwise it is zero if ${n}_X\cdot {E}_X={n}_Y\cdot {E}_Y$ .

As we will show now, this rudimentary information on Δμ = 0 and J = 0 can be used to mathematically describe the flux in the steady state even without further knowledge of the function J_X (V). The function J_X (V) describes a curve in a 2-dimensional space (Figure 2a). We know from the above consideration that J_X (E_X ) = 0 (Figure 2, blue point). This point is strategic because at E_X the direction of the flux changes, for example, if J_X < 0 at V < E_X , then J_X > 0 at V > E_X (Figure 2, grey regions). In the steady state, both the membrane voltage and the flux through the transporter are constant. Thus, in this case, we need to describe just one point of the J _X(V) curve (Figure 2b, red point). For such a description two coordinates are sufficient: One is the membrane voltage in steady state, V_ss , and the other is the flux in steady state, ${J}_X^{ss}$ . Using triangulation, ${J}_X^{ss}$ can be replaced by (Figure 2b, orange triangle)

(2) $$\begin{align}{J}_X^{ss}=m\cdot \left({V}_{ss}-{E}_X^{ss}\right),\end{align}$$

with slope m of magnitude $m=\tan \left(\alpha \right)$ . Flux values, J_X (V) have the unit s⁻¹, while membrane voltage and E_X have the unit V. Therefore, the slope m has the unit V⁻¹⋅s⁻¹. We can consider it as a free parameter (ranging from zero to infinity), which allows us to reach any point on the vertical green line in Figure 2b. In this way, any J_X (V) curve is represented by ${E}_X^{ss}$ , V_ss and the parameter m, even without knowing the exact shape of J_X (V). m combines several physiologically important parameters because a change in the activity of the transporter protein is reflected by a change of the slope m (i.e. the angle α). An increase in activity, for instance by enhanced expression of the transporter or its activation by posttranslational modifications, results in higher J-values (Figure 2c, dashed grey line, purple point). This is reflected by a larger angle (α1) and a larger slope (m1). Conversely, a reduction in activity (Figure 2c, dotted grey line, pink point) results in smaller J-values, a smaller angle (α2) and a smaller slope (m2). Thus, the parameter m gathers all regulatory features of the respective transporter. It can take on any value in the interval [0,∞). It should be noted, though, that m is generally not linear to biological parameters. The rather simple mathematical description of equation (2) therefore covers the entire range of possibilities without any constraints of approximations.

To illustrate this approach in practical terms, we develop it in the following, step by step, for (i) K⁺ channels, (ii) H⁺/K⁺ symporters and (iii) H⁺/K⁺ antiporters. Finally, we will also introduce a suitable model for the (iv) H⁺-ATPase that primarily energises all transport processes.

(i) K⁺ channels

The electrochemical gradient (compare equation (1a)) that determines the flux across K⁺ channels is given in the case of highly selective K⁺ channels by:

(3) $$\begin{align}{\Delta \mu}_{\mathrm{KC}}={\mu}_K^{\mathrm{in}}-{\mu}_K^{\mathrm{out}}= RT\cdot \ln \left(\frac{{\left[{K}^{+}\right]}_{\mathrm{in}}}{{\left[{K}^{+}\right]}_{\mathrm{out}}}\right)+F\cdot V,\end{align}$$

with the internal and external K⁺ concentrations [K⁺ ]_in and [K⁺ ]_out, respectively. This gradient is zero at the zero-flux potential for K⁺ channels, V = E_K , which is the Nernst potential for potassium

(4) $$\begin{align}{E}_K=\frac{RT}{F}\cdot \ln \left(\frac{{\left[{K}^{+}\right]}_{\mathrm{out}}}{{\left[{K}^{+}\right]}_{\mathrm{in}}}\right).\end{align}$$

Thus, at $V={E}_K$ the function of the flux, ${J}_{K,\mathrm{KC}}(V)$ (unit s⁻¹), has the value ${J}_{K,\mathrm{KC}}\left({E}_K\right)=0$ . Applying the same considerations that led to equation (2), the flux in steady state through the K⁺ channel can be described by:

(5) $$\begin{align}{J}_{K,\mathrm{KC}}^{ss}={J}_{K,\mathrm{KC}}\left({V}_{ss}\right)={m}_{\mathrm{KC}}\cdot \left({V}_{ss}-{E}_K^{ss}\right).\end{align}$$

The parameter m _KC (unit V⁻¹⋅s⁻¹) determines the activity of the channels and comprises all potential regulatory features, such as functional protein expression, post-translational modifications and/or voltage dependence. However, it should be noted again that m _KC is in most cases not linear to these biological parameters. The net K⁺ flux through the K⁺ channels causes a transmembrane electric current (unit A) given by:

(6) $$\begin{align}{I}_{\mathrm{KC}}={e}_0\cdot {J}_{K,\mathrm{KC}},\end{align}$$

with the elementary (positive) charge e₀ (= 1.6 × 10⁻¹⁹ A⋅s).

(ii) H⁺/K⁺ symporters

The electrochemical gradient (compare equation (1b)) that determines the flux across a symporter with n_S H⁺/1 K⁺ stoichiometry is given by:

(7) $$\begin{align}&{\Delta \mu}_{\mathrm{HKs}}={n}_s\cdot \left({\mu}_H^{\mathrm{in}}-{\mu}_H^{\mathrm{out}}\right)+\left({\mu}_K^{\mathrm{in}}-{\mu}_K^{\mathrm{out}}\right)\\&\quad= RT\cdot \left[{n}_s\cdot \ln \left(\frac{{\left[{H}^{+}\right]}_{\mathrm{in}}}{{\left[{H}^{+}\right]}_{\mathrm{out}}}\right)+\ln \left(\frac{{\left[{K}^{+}\right]}_{\mathrm{in}}}{{\left[{K}^{+}\right]}_{\mathrm{out}}}\right)\right]+F\cdot \left({n}_s+1\right)\cdot V,\nonumber\end{align}$$

with the internal and external proton concentrations [H⁺ ]_in and [H⁺ ]_out, respectively. This gradient is zero at the zero-flux potential V = E _HKs with

(8) $$\begin{align}&\left({n}_s+1\right)\cdot {E}_{\mathrm{HKs}}\\&\quad=\frac{RT}{F}\cdot \left[{n}_s\cdot \ln \left(\frac{{\left[{H}^{+}\right]}_{\mathrm{out}}}{{\left[{H}^{+}\right]}_{\mathrm{in}}}\right)+\ln \left(\frac{{\left[{K}^{+}\right]}_{\mathrm{out}}}{{\left[{K}^{+}\right]}_{\mathrm{in}}}\right)\right]={n}_s\cdot {E}_H+{E}_K,\nonumber\end{align}$$

and the Nernst potentials for potassium (E_K ) and protons (E_H ).

Hence, the net fluxes for K⁺ and protons, respectively, at V = E _HKs are J _K,HKs(E _HKs) =0 and J _H,HKs(E _HKs) = n_s ⋅ J _K,HKs(E _HKs) = 0. Applying the same considerations that led to equation (2), the fluxes in steady state can be described by:

(9) $$\begin{align}{J}_{H,\mathrm{HKs}}^{ss}={J}_{H,\mathrm{HKs}}\left({V}_{ss}\right)={m}_{\mathrm{HKs}}\cdot {n}_s\cdot \left(\left({n}_s+1\right)\cdot {V}_{ss}-{n}_s\cdot {E}_H-{E}_K^{ss}\right),\end{align}$$

(10) $$\begin{align}{J}_{K,\mathrm{HKs}}^{ss}={J}_{K,\mathrm{HKs}}\left({V}_{ss}\right)={m}_{\mathrm{HKs}}\cdot \left(\left({n}_s+1\right)\cdot {V}_{ss}-{n}_s\cdot {E}_H-{E}_K^{ss}\right).\end{align}$$

Also here, the parameter m _HKs gathers all regulatory features of the symporter. The transmembrane electric current (unit A) that is caused by the net H⁺ and K⁺ fluxes (unit s⁻¹) through the H⁺/K⁺ symporters is given by:

(11) $$\begin{align}{I}_{\mathrm{HKs}}={e}_0\cdot \left({J}_{H,\mathrm{HKs}}+{J}_{K,\mathrm{HKs}}\right).\end{align}$$

(iii) H⁺/K⁺ antiporters

The electrochemical gradient (compare equation (1c)) that determines the flux across an antiporter with n_a H⁺/1 K⁺ stoichiometry is given by:

(12) $$\begin{align}&{\Delta \mu}_{\mathrm{HKa}}={n}_a\cdot \left({\mu}_H^{\mathrm{in}}-{\mu}_H^{\mathrm{out}}\right)-\left({\mu}_K^{\mathrm{in}}-{\mu}_K^{\mathrm{out}}\right)\\&\quad= RT\cdot \left[{n}_a\cdot \ln \left(\frac{{\left[{H}^{+}\right]}_{\mathrm{in}}}{{\left[{H}^{+}\right]}_{\mathrm{out}}}\right)-\ln \left(\frac{{\left[{K}^{+}\right]}_{\mathrm{in}}}{{\left[{K}^{+}\right]}_{\mathrm{out}}}\right)\right]+F\cdot \left({n}_a-1\right)\cdot V.\nonumber\end{align}$$

This gradient is zero at E _HKa with

(13) $$\begin{align}&\left({n}_a-1\right)\cdot {E}_{\mathrm{HKa}}\\&\quad=\frac{RT}{F}\cdot \left[{n}_a\cdot \ln \left(\frac{{\left[{H}^{+}\right]}_{\mathrm{out}}}{{\left[{H}^{+}\right]}_{\mathrm{in}}}\right)-\ln \left(\frac{{\left[{K}^{+}\right]}_{\mathrm{out}}}{{\left[{K}^{+}\right]}_{\mathrm{in}}}\right)\right]={n}_a\cdot {E}_H-{E}_K.\nonumber\end{align}$$

If we apply again the same considerations that led to equation (2), the net fluxes for K⁺ and H⁺ in steady state can be expressed as:

(14) $$\begin{align}{J}_{H,\mathrm{HKa}}^{ss}={m}_{\mathrm{HKa}}\cdot {n}_a\cdot \left(\left({n}_a-1\right)\cdot {V}_{ss}-{n}_a\cdot {E}_H+{E}_K^{ss}\right).\end{align}$$

(15) $$\begin{align}{J}_{K,\mathrm{HKa}}^{ss}=-{m}_{\mathrm{HKa}}\cdot \left(\left({n}_a-1\right)\cdot {V}_{ss}-{n}_a\cdot {E}_H+{E}_K^{ss}\right).\end{align}$$

Here as well, the parameter m _HKa combines all regulatory characteristics of the antiporter. The transmembrane electric current (unit A) that is caused by the net H⁺ and K⁺ fluxes (unit s⁻¹) through the H⁺/K⁺ antiporters is given by:

(16) $$\begin{align}{I}_{\mathrm{HKa}}={e}_0\cdot \left({J}_{H,\mathrm{HKa}}+{J}_{K,\mathrm{HKa}}\right).\end{align}$$

It is zero, if n_a = 1. In this case, the antiport would be electroneutral and the flux would not depend on the membrane voltage, but would only be proportional to the difference E_H − E_K between the Nernst potentials of H⁺ and K⁺.

(iv) H⁺-ATPase

Finally, we consider the proton pump as an active membrane transporter. The derivation of a mathematical description of the pump current of H⁺-ATPases was presented in detail in the Supplementary Material by Reyer et al. (Reference Reyer, Häßler, Scherzer, Huang, Pedersen, Al-Rascheid, Bamberg, Palmgren, Dreyer, Nagel, Hedrich and Becker2020). The mechanistic model of a pump cycle (Dreyer, Reference Dreyer2017; Reyer et al., Reference Reyer, Häßler, Scherzer, Huang, Pedersen, Al-Rascheid, Bamberg, Palmgren, Dreyer, Nagel, Hedrich and Becker2020) results in a sigmoidal curve for pump currents (unit A) that describes very well the experimental findings (Lohse & Hedrich, Reference Lohse and Hedrich1992):

(17) $$\begin{align}{I}_{\mathrm{HATPase}}={I}_{\mathrm{Hmax}}\cdot \frac{1-{e}^{-\left(\frac{F}{RT}\cdot \left(V-{V}_{0,\mathrm{pump}}\right)\right)}}{1+{e}^{-\left(\frac{F}{RT}\cdot V+B\right)}+{e}^{-\left(D\cdot \frac{F}{RT}\cdot V+C\right)}}={I}_{\mathrm{Hmax}}\cdot {i}_p(V),\end{align}$$

with ${i}_p(V)=\frac{1-{e}^{-\left(\frac{F}{RT}\cdot \left(V-{V}_{0,\mathrm{pump}}\right)\right)}}{1+{e}^{-\left(\frac{F}{RT}\cdot V+B\right)}+{e}^{-\left(D\cdot \frac{F}{RT}\cdot V+C\right)}}$ . The proton flux mediated by the H⁺-ATPase (J _HATPase; unit s⁻¹) can be determined by dividing the current by the elementary (positive) charge e₀ :

(18) $$\begin{align}{J}_{\mathrm{HATPase}}=\frac{I_{\mathrm{HATPase}}}{e_0}&=\frac{I_{\mathrm{Hmax}}}{e_0}\cdot \frac{1-{e}^{-\left(\frac{F}{RT}\cdot \left(V-{V}_{0,\mathrm{pump}}\right)\right)}}{1+{e}^{-\left(\frac{F}{RT}\cdot V+B\right)}+{e}^{-\left(D\cdot \frac{F}{RT}\cdot V+C\right)}}\nonumber\\&=\frac{I_{\mathrm{Hmax}}}{e_0}\cdot {i}_p(V).\end{align}$$

Here, the numerical values of B, C and D are of minor importance for the shape of the curve. Suitable values are B = 5.4, C = 2.5, D = 0.1 (Dreyer, Reference Dreyer2021a). Most important are the parameters I _Hmax and V _0,pump. I _Hmax determines the maximum pump current and is proportional to the number of active proton pumps in the membrane, while V _0,pump is the voltage at which the pump current is zero [i_p (V _0,pump) = 0]. It depends in a complex manner on the cytosolic and apoplastic proton concentrations and the cytosolic ATP, ADP and Pi concentrations (Rienmüller et al., Reference Rienmüller, Dreyer, Schönknecht, Schulz, Schumacher, Nagy, Martinoia, Marten and Hedrich2012). For suitable simulations, in which the energy status of the cell does not change in the considered time interval, V _0,pump can be set and kept constant (e.g. V _0,pump = −200 mV). This value determines the most negative value that the membrane voltage can attain.

2.2 Mathematical description of changes in concentrations and voltage

The proton and K⁺ fluxes change the concentrations of H⁺ and K⁺ on both sides of the membrane. The combined net K⁺ efflux (J_K = J _K,KC + J _K,KHs + J _K,KHa; unit s⁻¹) increases [K⁺ ]_out and decreases [K⁺ ]_in, while the combined H⁺ efflux (J_H = J _HATPase + J _H,KHs + J _H,KHa; unit s⁻¹) increases [H⁺ ]_out and decreases [H⁺ ]_in. Additionally, the buffer capacities of the internal and external compartments can partially mitigate the changes in proton concentration, described by the general buffer reaction ${H}^{+}+{Buf}^{-}\rightleftarrows HBuf$ with an anionic base Buf^- and its conjugate acid, HBuf. All these changes are represented in mathematical terms by ordinary differential equations. The changes in K⁺ concentrations are governed by:

(19) $$\begin{align}\frac{d{\left[{K}^{+}\right]}_{\mathrm{in}}}{dt}=\frac{-{J}_K}{{\mathrm{Vol}}_{\mathrm{in}}\cdot {N}_A}=\frac{-1}{{\mathrm{Vol}}_{\mathrm{in}}\cdot {N}_A}\cdot \left[{J}_{K,\mathrm{KC}}+{J}_{K,\mathrm{HKs}}+{J}_{K,\mathrm{HKa}}\right].\end{align}$$

(20) $$\begin{align}&\frac{d{\left[{K}^{+}\right]}_{\mathrm{out}}}{dt}=\frac{J_K}{{\mathrm{Vol}}_{\mathrm{out}}\cdot {N}_A}=\frac{1}{{\mathrm{Vol}}_{\mathrm{out}}\cdot {N}_A}\cdot \left[{J}_{K,\mathrm{KC}}+{J}_{K,\mathrm{HKs}}+{J}_{K,\mathrm{HKa}}\right].\end{align}$$

Here N_A is the Avogadro constant and Vol _in and Vol _out are the volumes of the compartments on the respective side of the membrane. Changes in proton concentrations are determined by:

(21) $$\begin{align}&\!\!\!\!\!\!\frac{d{\left[{H}^{+}\right]}_{\mathrm{in}}}{dt}=\frac{-1}{{\mathrm{Vol}}_{\mathrm{in}}\cdot {N}_A}\cdot \big\{\left[{J}_{\mathrm{HATPase}}+{J}_{H,\mathrm{HKs}}+{J}_{H,\mathrm{HKa}}\right]\quad\\&\quad-{k}_{v,\mathrm{in}}\cdot \left({\left[{H}^{+}\right]}_{\mathrm{in}}\cdot {\left[{Buf}^{-}\right]}_{\mathrm{in}}-{K}_{s,\mathrm{in}}\cdot {\left[ HBuf\right]}_{\mathrm{in}}\right)\big\}.\nonumber\end{align}$$

(22) $$\begin{align}&\!\!\!\!\!\!\frac{d{\left[{H}^{+}\right]}_{\mathrm{out}}}{dt}=\frac{1}{{\mathrm{Vol}}_{\mathrm{out}}\cdot {N}_A}\cdot \big\{\left[{J}_{\mathrm{HATPase}}+{J}_{H,\mathrm{HKs}}+{J}_{H,\mathrm{HKa}}\right]\quad\\&\quad-{k}_{v,\mathrm{out}}\cdot \left({\left[{H}^{+}\right]}_{\mathrm{out}}\cdot {\left[{Buf}^{-}\right]}_{\mathrm{out}}-{K}_{s,\mathrm{out}}\cdot {\left[ HBuf\right]}_{\mathrm{out}}\right)\big\}.\nonumber\end{align}$$

Here, ${K}_{s,\mathrm{in}}={10}^{-{pK}_{s,\mathrm{in}}}$ and ${K}_{s,\mathrm{out}}={10}^{-{pK}_{s,\mathrm{out}}}$ indicate the dissociation constants of the buffer reactions that buffer pH _in and pH _out to pK _s,in and pK _s,out, respectively. The buffer concentrations, in turn, change according to:

(23) $$\begin{align}&\frac{d{\left[{Buf}^{-}\right]}_{\mathrm{in}}}{dt}=-\frac{d{\left[ HBuf\right]}_{\mathrm{in}}}{dt}\\&\quad=\frac{-1}{{\mathrm{Vol}}_{\mathrm{in}}\cdot {N}_A}\cdot {k}_{v,\mathrm{in}}\cdot \left({\left[{H}^{+}\right]}_{\mathrm{in}}\cdot {\left[{Buf}^{-}\right]}_{\mathrm{in}}-{K}_{s,\mathrm{in}}\cdot {\left[ HBuf\right]}_{\mathrm{in}}\right).\nonumber\end{align}$$

(24) $$\begin{align}&\frac{d{\left[{Buf}^{-}\right]}_{\mathrm{out}}}{dt}=-\frac{d{\left[ HBuf\right]}_{\mathrm{out}}}{dt}\\&\quad=\frac{-1}{Vol_{\mathrm{out}}\cdot {N}_A}\cdot {k}_{v,\mathrm{out}}\cdot \left({\left[{H}^{+}\right]}_{\mathrm{out}}\cdot {\left[{Buf}^{-}\right]}_{\mathrm{out}}-{K}_{s,\mathrm{out}}\cdot {\left[ HBuf\right]}_{\mathrm{out}}\right).\nonumber\end{align}$$

Compared to the transport reactions, the buffering reactions take place very quickly, which is reflected in large k _v,in and k _v,out values. Because of this difference in time scale, the parameters k _v,in and k _v,out do not need to be known in detail. In simulations, it is sufficient to set them 2–3 magnitudes larger than the other parameters (e.g. k _v,[in,out] /(I _max/e₀ ) = 100 μM⁻²) to obtain almost instantaneous buffer reactions in comparison to the time-scale of transport.

The net fluxes usually cause a net charge transport across the membrane. This non-zero current provokes a change in the membrane voltage. The membrane is actually comparable to a dielectric in an electric capacitor that separates the charges in the aqueous interior from the aqueous exterior. Charge transport from one side to the other changes the membrane voltage according to:

(25) $$\begin{align}\frac{dV}{dt}=\frac{-1}{C}\cdot \left[{I}_{\mathrm{HATPase}}+{I}_{\mathrm{KC}}+{I}_{\mathrm{HKs}}+{I}_{\mathrm{HKa}}\right],\end{align}$$

with the membrane capacitance C (unit F).

2.3 Modelling the steady state with fixed transporters activities

At this point, we have translated Figure 1 into the language of mathematics. For each transporter, we expressed the fluxes and currents in equations (5), (6), (9)–(11) and (14)–(18) with a few parameters. The dynamic properties of the system are determined by equations (19)–(25). In the following, we will solve the differential equations for homeostatic conditions. If the system is considered in a steady state, neither the membrane voltage nor the concentrations change with time, which means that the left sides of equations (19)–(25) are zero. In this condition, the equations (19)–(25) simplify to two non-redundant equations taking also into account that I_X = e₀ ⋅J_X :

(26) $$\begin{align}0={J}_{K,\mathrm{KC}}^{ss}+{J}_{K,\mathrm{HKs}}^{ss}+{J}_{K,\mathrm{HKa}}^{ss}.\end{align}$$

(27) $$\begin{align}0={J}_{\mathrm{HATPase}}^{ss}+{J}_{H,\mathrm{HKs}}^{ss}+{J}_{H,\mathrm{HKa}}^{ss}.\end{align}$$

All the other equations are either zero or linear combinations of the equations (26) and (27). In the next step, the ${J}_X^{ss}$ were replaced by the expressions deduced in equations (5), (9), (10), (14), (15) and (18) yielding:

(28) $$\begin{align}0&={m}_{\mathrm{KC}}\cdot \left({V}_{ss}-{E}_K^{ss}\right)+{m}_{\mathrm{HKs}}\cdot \left(\left({n}_s+1\right)\cdot {V}_{ss}-{n}_s\cdot {E}_H-{E}_K^{ss}\right)\\&\quad-{m}_{\mathrm{HKa}}\cdot \left(\left({n}_a-1\right)\cdot {V}_{ss}-{n}_a\cdot {E}_H+{E}_K^{ss}\right).\nonumber\end{align}$$

(29) $$\begin{align}0&=\frac{I_{\mathrm{Hmax}}}{e_0}\cdot {i}_p\left({V}_{ss}\right)+{n}_s\cdot {m}_{\mathrm{HKs}}\cdot \left(\left({n}_s+1\right)\cdot {V}_{ss}-{n}_s\cdot {E}_H-{E}_K^{ss}\right)\\&\quad+{n}_a\cdot {m}_{\mathrm{HKa}}\cdot \left(\left({n}_a-1\right)\cdot {V}_{ss}-{n}_a\cdot {E}_H+{E}_K^{ss}\right).\nonumber\end{align}$$

Although the mathematical parameters m_X can depend on relevant biological parameters such as the K⁺ and H⁺ concentrations or V (in a linear and non-linear way), they are constant in the steady state and are therefore fixed numbers with units. Then, we multiplied both equations with the factor ${e}_0/{I}_{\mathrm{Hmax}}$ and define ${{g}_X:= {m}_X\cdot {e}_0/{I}_{\mathrm{Hmax}}}$ (unit V⁻¹). This operation eliminated a redundant parameter and reduced the set of four parameters (I _Hmax, m _KC, m _HKs, m _HKa) to a set of three (g _KC, g _KHs, g _KHa) representing the relative activity of the transporter proteins (relative to the maximal activity of the H⁺-ATPase).

(30) $$\begin{align}0&={g}_{\mathrm{KC}}\cdot \left({V}_{ss}-{E}_K^{ss}\right)+{g}_{\mathrm{HKs}}\cdot \left(\left({n}_s+1\right)\cdot {V}_{ss}-{n}_s\cdot {E}_H-{E}_K^{ss}\right)\\&\quad-{g}_{\mathrm{HKa}}\cdot \left(\left({n}_a-1\right)\cdot {V}_{ss}-{n}_a\cdot {E}_H+{E}_K^{ss}\right).\nonumber\end{align}$$

(31) $$\begin{align}0&={i}_p\left({V}_{ss}\right)+{n}_s\cdot {g}_{\mathrm{HKs}}\cdot \left(\left({n}_s+1\right)\cdot {V}_{ss}-{n}_s\cdot {E}_H-{E}_K^{ss}\right)\\&\quad+{n}_a\cdot {g}_{\mathrm{HKa}}\cdot \left(\left({n}_a-1\right)\cdot {V}_{ss}-{n}_a\cdot {E}_H+{E}_K^{ss}\right).\nonumber\end{align}$$

Additionally, the equations (30) and (31) could be arranged in vector and matrix form:

(32) $$\begin{align}M\bullet \left(\begin{array}{c}{V}_{ss}\\ {}{E}_K^{ss}\end{array}\right)=\overrightarrow{A},\end{align}$$

with:

(33) $$\begin{align}M=\left(\begin{array}{@{}cc@{}}{g}_{\mathrm{KC}}+\left({n}_s+1\right)\cdot {g}_{\mathrm{HKs}}-\left({n}_a-1\right)\cdot {g}_{\mathrm{HKa}}& -{g}_{\mathrm{KC}}-{g}_{\mathrm{HKs}}-{g}_{\mathrm{HKa}}\\ {}{n}_s\cdot \left({n}_s+1\right)\cdot {g}_{\mathrm{HKs}}+{n}_a\cdot \left({n}_a-1\right)\cdot {g}_{\mathrm{HKa}}& -{n}_s\cdot {g}_{\mathrm{HKs}}+{n}_a\cdot {g}_{\mathrm{HKa}}\end{array}\right),\end{align}$$

and:

(34) $$\begin{align}\overrightarrow{A}=\left(\begin{array}{c}{n}_s\cdot {g}_{\mathrm{HKs}}\cdot {E}_H-{n}_a\cdot {g}_{\mathrm{HKa}}\cdot {E}_H\\ {}{n}_s^2\cdot {g}_{\mathrm{HKs}}\cdot {E}_H+{n}_a^2\cdot {g}_{\mathrm{HKa}}\cdot {E}_H-{i}_p\left({V}_{ss}\right)\end{array}\right).\end{align}$$

This might be considered as just a formalism and not needed for two equations. However, the conversion should illustrate that we have opened the powerful toolbox of linear algebra. More complex homeostats with more than two equations can be handled in the same way and the tools from linear algebra for solving this type of equation can be applied in the same manner. Equation (32) has exactly one solution for the $\left({V}_{ss}|{E}_K^{ss}\right)$ pair that depends on the activities of the transporter proteins (g_X ) and the proton gradient (E_H ):

(35) $$\begin{align}{V}_{ss}\kern1.2pt{-}\kern1.2pt{E}_H+{i}_p\left({V}_{ss}\right)\cdot \frac{g_{\mathrm{KC}}+{g}_{\mathrm{HKs}}+{g}_{\mathrm{HKa}}}{{\left({n}_s+{n}_a\right)}^2\cdot {g}_{\mathrm{HKs}}\cdot {g}_{\mathrm{HKa}}+{g}_{\mathrm{KC}}\cdot \left({n}_s^2\cdot {g}_{\mathrm{HKs}}+{n}_a^2\cdot {g}_{\mathrm{HKa}}\right)}\kern-1pt=0.\end{align}$$

(36) $$\begin{align}{E}_K^{ss}={E}_H-{i}_p\left({V}_{ss}\right)\cdot \frac{g_{\mathrm{KC}}+\left({n}_s+1\right)\cdot {g}_{\mathrm{HKs}}-\left({n}_a-1\right)\cdot {g}_{\mathrm{HKa}}}{{\left({n}_s+{n}_a\right)}^2\cdot {g}_{\mathrm{HKs}}\cdot {g}_{\mathrm{HKa}}+{g}_{\mathrm{KC}}\cdot \left({n}_s^2\cdot {g}_{\mathrm{HKs}}+{n}_a^2\cdot {g}_{\mathrm{HKa}}\right)}\end{align}$$ $$\kern0.96em ={V}_{ss}-{i}_p\left({V}_{ss}\right)\cdot \frac{n_s\cdot {g}_{\mathrm{HKs}}-{n}_a\cdot {g}_{\mathrm{HKa}}}{{\left({n}_s+{n}_a\right)}^2\cdot {g}_{\mathrm{HKs}}\cdot {g}_{\mathrm{HKa}}+{g}_{\mathrm{KC}}\cdot \left({n}_s^2\cdot {g}_{\mathrm{HKs}}+{n}_a^2\cdot {g}_{\mathrm{HKa}}\right)}.$$

Please note

1. (1) Equations (35) and (36) are not defined for cases in which two of the three parameters g _KC, g _HKs and g _HKa are zero. In cases, in which only one g_X is different from zero, the equations (30) and (31) result in:

   g _KC = g _HKs = 0 (only active H⁺/K⁺ antiporters): ${V}_{ss}={V}_{0,\mathrm{pump}}$ , ${E}_K^{ss}={n}_a\cdot {E}_H-\left({n}_a-1\right)\cdot {V}_{0,\mathrm{pump}}$

   g _KC = g _HKa = 0 (only active H⁺/K⁺ symporters): ${V}_{ss}={V}_{0,\mathrm{pump}}$ , ${E}_K^{ss}=\left({n}_s+1\right)\cdot {V}_{0,\mathrm{pump}}-{n}_s\cdot {E}_H$

   g _HKs = g _HKa = 0 (only active K⁺ channels): ${V}_{ss}={V}_{0,\mathrm{pump}}$ , ${E}_K^{ss}={V}_{0,\mathrm{pump}}$

2. (2) Equation (35); of type ${V}_{ss}-A+{i}_p\left({V}_{ss}\right)\cdot B=0$ , is an implicit solution. It cannot be explicitly solved analytically for V_ss because i_p (V_ss ) is a function of V_ss . Nevertheless, mathematical tools also allow us to bypass this difficulty. The left-hand side of equation (35) before the equals sign is a strictly monotonically increasing function of V_ss , which means that the first derivative after V_SS is > 0. Or in other words, with increasing V_ss , the left-hand side of the equation also increases, starting from negative values at negative V_ss and reaching positive values at positive V_ss . The corresponding curve therefore cuts the zero axis only once. Thus, for given values of $A={E}_H$ and $B=f\left({g}_{\mathrm{KC}},{g}_{\mathrm{HKs}},{g}_{\mathrm{HKa}}\right)$ there is exactly one V_ss that obeys this equation. The system is free from potential bifurcation. If E_H , g _KC, g _HKs and g _HKa are known, the value of V_ss can be determined numerically by root-finding algorithms, such as the Newton–Raphson method.