The Boron & De Weer Model of Intracellular pH Regulation

1Introduction¶

In 1976 Boron & De Weer published their landmark paper on “Intracellular $\mathrm{pH}$ transients in squid giant axons caused by $\mathrm{CO_2}$, $\mathrm{NH_3}$, and metabolic inhibitors” Boron & De Weer, 1976. The authors used a squid giant axon preparation and a mathematical model of $\mathrm{pH}$ buffering and the transport of protons, bicarbonate ($\mathrm{HCO_3^-}$) and $\mathrm{CO_2}$ to establish the experimental evidence for active regulation of intracellular $\mathrm{pH}$ ($\mathrm{pH_i}$) by a transporter in the plasma membrane that — at the expense of energy — either moves acid out of the cell, or base into the cell. Today, we refer to such a transporter generically as an acid-extrusion mechanism. For simplicity, Boron & De Weer modelled it as a proton pump, although the result would have been almost indistinguishable had they modelled it as the uptake of $\mathrm{HCO_3^-}$ or carbonate ($\mathrm{CO_3^{2-}}$). The paper reported on the consequences of adding and then removing extracellular $\mathrm{CO_2}/\mathrm{HCO_3^-}$, $\mathrm{NH_3/NH_4^+}$ (where $\mathrm{NH_4^+}$ is ammonium), or the metabolic inhibitors, cyanide, azide and dinitrophenol (DNP).

In the first experiment, following exposure of the cell to elevated $\mathrm{CO_2}$ and $\mathrm{HCO_3^-}$, $\mathrm{CO_2}$ rapidly enters the cell and intracellular $\mathrm{CO_2}$ equilibrates with the extracellular $\mathrm{CO_2}$, and generates intracellular $\mathrm{H^+}$ and $\mathrm{HCO_3^-}$ via the $\mathrm{CO_2}$ hydration reaction ($\mathrm{CO_2}+\mathrm{H_2O} \longrightarrow \mathrm{H^+}+\mathrm{HCO^-_3}$). The accumulating $\mathrm{H^+}$ results in a rapid fall of $\mathrm{pH_i}$ (see Figure 1A & Figure 1B). To the extent that the membrane is permeable to $\mathrm{HCO_3^-}$ as well as to $\mathrm{CO_2}$, $\mathrm{HCO_3^-}$ will initially enter the cell passively, along its electrochemical gradient. Soon, however, the accumulation of intracellular $\mathrm{HCO_3^-}$ reverses the $\mathrm{HCO_3^-}$ electrochemical gradient and would be expected to lead to the passive efflux of $\mathrm{HCO_3^-}$. This loss of cellular $\mathrm{HCO_3^-}$ would tend to acidify the cell because — to replenish the lost $\mathrm{HCO_3^-}$ — additional $\mathrm{CO_2}$ would enter the cell and form even more $\mathrm{H^+}$ and $\mathrm{HCO_3^-}$ (the passive $\mathrm{CO_2}/\mathrm{HCO_3^-}$ shuttle). Thus, the expectation was that prolonged exposure to $\mathrm{CO_2}$ would cause $\mathrm{pH_i}$ to fall rapidly (passive influx of $\mathrm{CO_2}$) and then to drift more slowly in the acidic direction (passive efflux of $\mathrm{HCO_3^-}$). In fact, Boron & De Weer observed an alkaline drift, leading to the postulate of active extrusion of $\mathrm{H^+}$ — or an equivalent process^[1] — at a rate that exceeds the passive shuttling by the $\mathrm{CO_2}/\mathrm{HCO_3^-}$ couple (see Figure 1A & Figure 1C).

Following removal of external $\mathrm{CO_2}$, intracellular $\mathrm{CO_2}$ diffuses out, while intracellular $\mathrm{HCO_3^-}$ combines with $\mathrm{H^+}$ to leave the cell as $\mathrm{CO_2}$. Thus, the entire intracellular $\mathrm{H^+}$ load associated with $\mathrm{CO_2}$ entry would be removed, returning $\mathrm{pH_i}$ to its value before the addition of $\mathrm{CO_2}/\mathrm{HCO_3^-}$ (or to a slightly lower value, to the extent that $\mathrm{HCO_3^-}$ had passively exited during the $\mathrm{CO_2}/\mathrm{HCO_3^-}$ exposure). In fact, Boron & De Weer observed that $\mathrm{pH_i}$ overshoots its resting value by an amount consistent with the net removal of $\mathrm{H^+}$ by the active, acid-extrusion mechanism during the $\mathrm{CO_2}/\mathrm{HCO_3^-}$ exposure (see Figure 1A & Figure 1C).

In the second experiment, following exposure of the cell to extracellular $\mathrm{NH_3}/\mathrm{NH_4^+}$ in the form of $\mathrm{NH_4Cl}$ (ammonium chloride), the intracellular environment rapidly becomes alkaline as $\mathrm{NH_3}$ enters and combines with $\mathrm{H^+}$ to form $\mathrm{NH_4^+}$ (equivalent to the hydration of $\mathrm{NH_3}$ to form $\mathrm{NH_4^+}$ and $\mathrm{OH^-}$). If this were the entire story, then $\mathrm{pH_i}$ would rise monotonically to a relatively alkaline value, and then the subsequent removal of $\mathrm{NH_4Cl}$ would cause $\mathrm{pH_i}$ to fall to precisely its initial value. In fact, Boron & De Weer observed that the exposure to $\mathrm{NH_4Cl}$ causes $\mathrm{pH_i}$ to rise rapidly and then fall slowly. Moreover, the subsequent removal of $\mathrm{NH_3}/\mathrm{NH_4^+}$ causes $\mathrm{pH_i}$ to undershoot its original value (see Figure 2A). Thus, Boron & De Weer postulated that, during the $\mathrm{NH_3}/\mathrm{NH_4^+}$ exposure, $\mathrm{NH_4^+}$ passively enters the cell down its electrochemical gradient. Early during the exposure, this $\mathrm{NH_4^+}$ influx would oppose the $\mathrm{NH_3}$ entry and slightly reduce the $\mathrm{pH_i}$ increase. Later during the $\mathrm{NH_4Cl}$ exposure, after intracellular $\mathrm{[NH_3]}$ rises to match extracellular $\mathrm{[NH_3]}$ ($\mathrm{[NH_3]_o}$), the continued passive influx of $\mathrm{NH_4^+}$ would generate intracellular $\mathrm{H^+}$ and $\mathrm{NH_3}$. The result would be a slow fall of $\mathrm{pH_i}$ and a rise in intracellular $\mathrm{[NH_3]}$, the latter leading to the passive exit of $\mathrm{NH_3}$ (the passive $\mathrm{NH_3}/\mathrm{NH_4^+}$ shuttle, see Figure 2A & Figure 2B). The experimental data are consistent with the proposed mathematical model.

Finally, exposure of the cells — in turn — to cyanide, DNP and azide resulted in intracellular acidosis, consistent with the accumulation of acid metabolites.

In the present paper, we re-formulate the models from Boron & De Weer (1976), henceforth referred to as ‘BDW’, and specify the simulation using the Physiome modelling standards CellML Cuellar et al., 2003 and SED-ML Bergmann et al., 2017 in order to ensure that the model reproduces the graphs in the original paper and that the model is fully curated.^[2] Note that this effort requires the specification of some parameters used in BDW’s simulations, but not described in the BDW paper. The curated and annotated model is made available in a form that users can run with OpenCOR^[3] to understand the model and to explore the effect of changes in parameter values.

2pH Buffering by Weak Acids and Weak Bases¶

We begin by reviewing a few rudimentary concepts of $\mathrm{pH}$ buffering by weak acids and bases Roos & Boron, 1981Bevensee & Boron, 2013Boron & Boulpaep, 2016 to provide the background for understanding the derivation and implementation of the BDW model.

Buffers.¶

According to Brönsted’s definition Brönsted, 1923, an acid is any substance that can donate a $\mathrm{H^+}$. Conversely, a base is any substance that can accept a $\mathrm{H^+}$. A buffer is any substance that can reversibly consume or produce $\mathrm{H^+}$, thereby minimising changes in $\mathrm{pH}$.

The dissociation of the uncharged weak acid ($\mathrm{HA}$) to the anionic weak base ($\mathrm{A^-}$) is described by the equilibrium reaction:

which is governed by the equilibrium constant^[4]

An example is the carbonic acid ($\mathrm{H_2CO_3}$) dissociation reaction,

The total weak acid concentration, $\mathrm{[TA]}$, is the sum of $\mathrm{[HA]}$ and $\mathrm{[A^-]}$. Note that $\mathrm{[TA]}$ is one of the two main unknowns in the BDW model for weak acids.

The dissociation of the cationic weak acid ($\mathrm{BH^+}$) to the uncharged weak base ($\mathrm{B}$) is described by the equilibrium reaction,

where the equilibrium constant is

An example is the $\mathrm{NH_4^+}$ dissociation reaction,

The total weak base concentration, $\mathrm{[TB]}$, is the sum of $\mathrm{[BH^+]}$ and $\mathrm{[B]}$. Note that $\mathrm{[TB]}$ is one of the two main unknowns in the BDW model for weak bases.

The $\mathrm{CO_2}/\mathrm{HCO_3^-}$ buffer pair.¶

The formation of $\mathrm{HCO_3^-}$ and $\mathrm{H^+}$ from $\mathrm{CO_2}$ by hydration is given by the equilibrium reaction

where the equilibrium constant is

Taking logarithms of both sides of (5), and recognising from Henry’s law that

where $s$ is the solubility coefficient for $\mathrm{CO_2}$ and $p_\mathrm{CO_2}$ is the partial pressure of $\mathrm{CO_2}$, we obtain the familiar Henderson-Hasselbalch equation

Here, $\mathrm{pH}=-\log[\mathrm{H^+}]$^[5] and $\mathrm{pK_{CO_2}}=-\log(K_\mathrm{CO_2})$.

In terms of the nomenclature above, one might regard $\mathrm{CO_2}$ as the weak acid $\mathrm{HA}$^[6], and $\mathrm{HCO_3^-}$ as its conjugate base $\mathrm{A^-}$.

Buffering power ($\beta$).¶

By definition, $\beta$ is the amount of strong base (e.g., $\mathrm{NaOH}$), or the negative of the amount of strong acid (e.g., $\mathrm{HCl}$), that one must add to $1~\mathrm{L}$ of solution to raise $\mathrm{pH}$ by one $\mathrm{pH}$ unit:

The units of $\beta$ are $\mathrm{mM}$. For additional details, refer to Roos & Boron (1981)Boron & Boulpaep (2016). Note that BDW defined $\beta$ as a negative number, as did Koppel and Spiro in their original definition of buffering Koppel, 1914Roos & Boron, 1980, rather than as a now-conventional positive quantity, as did Van Slyke in his later work Van Slyke, 1922. BDW’s definition, which they consistently applied, has no effect on the outcome of their simulations. In the present paper, we will follow the definition of Van Slyke — defining $\beta$ as a positive number — and make appropriate sign changes to the derived equations.

3The Boron & De Weer Model for the Permeation by an Uncharged Weak Acid and its Conjugate, Anionic Weak Base¶

The BDW mathematical model consists of two time-dependent ordinary differential equations (ODEs), one describing the time-course of the concentration of total intracellular buffer ($\mathrm{[TA]_i} = \mathrm{[HA]_i}+\mathrm{[A^-]_i}$) and the other the time-course of the intracellular free $\mathrm{H^+}$ concentration (which is related to $\mathrm{pH_i}$). BDW derived these two equations for the general cases in which any buffer pair $\mathrm{HA}/\mathrm{A^-}$, or any buffer pair $\mathrm{B}/\mathrm{BH^+}$, can move passively across the plasma membrane of a prototype cell. Then, they applied these two general equations to their specific experimental conditions, namely exposure of a cell (a squid giant axon) to equilibrated extracellular $\mathrm{CO_2}/\mathrm{HCO_3^-}$ or to equilibrated extracellular $\mathrm{NH_3}/\mathrm{NH_4^+}$.

Here, following BDW’s approach, we begin by deriving the equations for $\mathrm{HA}/\mathrm{A^-}$. In the next section, we apply the same general formalism to $\mathrm{B}/\mathrm{BH^+}$.

Derivation for weak acids.¶

Imagine that a cell is exposed to a solution containing equilibrated $\mathrm{HA}/\mathrm{A^-}$ and that both $\mathrm{HA}$ and $\mathrm{A^-}$ initially move into the cell — because of the chemical gradient in the case of $\mathrm{HA}$, and because of the electrochemical gradient in the case of $\mathrm{A^-}$.

An integrated form of Fick’s first law of diffusion describes the net passive influx^[7] of $\mathrm{HA}$ ($J_\mathrm{HA}$)

where $J_\mathrm{HA}$ is flux ($\mathrm{mol\cdot m^{-2}\cdot s^{-1}}$) and $P_\mathrm{HA}$ ($\mathrm{m\cdot s^{-1}}$) is the membrane permeability to the uncharged weak acid $\mathrm{HA}$. Note that this is a passive diffusion equation because $\mathrm{HA}$ is uncharged.

The constant field equation — also known as the Goldman, Hodgkin, Katz (GHK) Goldman, 1943Hodgkin & Huxley, 1952 equation — describes the net passive influx of $\mathrm{A^-}$ ($J_\mathrm{A^-}$):

where $P_\mathrm{A^-}$ ($\mathrm{m\cdot s^{-1}}$) is the membrane permeability to the charged conjugate base $\mathrm{A^-}$, $V_\mathrm{m}$ is the membrane potential (intracellular relative to extracellular potential), and $\epsilon$ is a shorthand for $e^{{-V_\mathrm{m}F}/{RT}}$. Note that $J_\mathrm{HA}$ and $J_\mathrm{A^-}$ have units of $\mathrm{mol\cdot m^{-2}\cdot s^{-1}}$.

Although $\mathrm{HA}$ and $\mathrm{A^-}$ can interconvert in the cytosol, BDW assumed that the intracellular concentration of total weak acid $\mathrm{[TA]_i}$ only can change due to the transmembrane fluxes of $\mathrm{HA}$ and $\mathrm{A^-}$ (see Figure 3). Thus, the time rate of change of $\mathrm{[TA]_i}$ is

where $\rho$ ($\mathrm{m^{-1}}$) is the area-to-volume ratio for the cell, and converts the transmembrane flux per unit area (in units of $\mathrm{mol\cdot m^{-2}\cdot s^{-1}}$) to a time rate of change per unit cell volume ($\mathrm{mol\cdot m^{-3}\cdot s^{-1}}$ or $\mathrm{mM\cdot s^{-1}}$). (11) is the first of two ODEs of the BDW model for the buffer pair $\mathrm{HA}/\mathrm{A^-}$.

Later, Bevensee and Boron defined the time rate of change per unit volume (e.g., $\mathrm{d[TA]_i/dt}$) as a ‘pseudoflux’ $\phi$, with the area-to-volume ratio folded into the value of $\phi$ Bevensee & Boron, 2013. Physiologists sometimes prefer to present experimental data in terms of pseudoflux because most mammalian cells often have complex geometries that make it difficult to estimate surface area.

In deriving the second ODE of their model, BDW started by noting that the time rate of change of free protons, $\mathrm{d[H^+]_i}/\mathrm{d}t$, depends on the rate at which acids are added into the cytosol per unit volume and per unit time — denoted $\mathrm{d}Q/\mathrm{d}t$ ($\mathrm{mol\cdot m^{-3}\cdot s^{-1}}$) where $Q$ is the total intracellular acid content. Like $\mathrm{d[TA]_i}/\mathrm{d}t$, both $\mathrm{d[H^+]_i}/\mathrm{d}t$ and $\mathrm{d}Q/\mathrm{d}t$ are pseudofluxes.

In their simple system of a squid giant axon exposed to $\mathrm{CO_2}$ (i.e., $\mathrm{HA}$) and $\mathrm{HCO_3^-}$ (i.e., $\mathrm{A^-}$), BDW assumed that only three general processes affect $\mathrm{d}Q/\mathrm{d}t$: (i) the release of $\mathrm{H^+}$ by some fraction ($x$) of the entering $\mathrm{HA}$ (i.e., $xJ_\mathrm{HA}$), (ii) the consumption of $\mathrm{H^+}$ by some fraction $(y)$ of the entering $\mathrm{A^-}$ (i.e., $yJ_\mathrm{A^-}$), and (iii) the “additional” rate of intracellular consumption or active extrusion of $\mathrm{H^+}$ ($J_\mathrm{H^+}$; see Figure 3) above the fixed background rate of $\mathrm{H^+}$ extrusion necessary to balance the fixed background rate of acid loading (i.e., addition of $\mathrm{H^+}$ or equivalent acid, or removal of $\mathrm{OH^-}$ or equivalent base) in the absence of $\mathrm{HA}/\mathrm{A^-}$. Thus,

A critical insight by BDW is that during each infinitesimal increment in time during which a bolus of $\mathrm{HA}$ enters the cell, the entering $\mathrm{HA}$ redistributes itself between $\mathrm{HA}$ (and $\mathrm{H^+}$) vs $\mathrm{A^-}$, according to the pre-existing ratios $\mathrm{[HA]_i}/\mathrm{[TA]_i}$ and $\mathrm{[A^-]_i}/\mathrm{[TA]_i}$. Thus, the fraction $y$ of entering $\mathrm{HA}$ that remains $\mathrm{HA}$ is

This is also the fraction of entering $\mathrm{A^-}$ that combines with $\mathrm{H^+}$ and becomes $\mathrm{HA}$. Combining the above expression with (2),

which BDW defined as $\alpha$. Conversely, the fraction $x$ of $\mathrm{A^-}$ that remains $\mathrm{A^-}$ is

This is also the fraction of entering $\mathrm{HA}$ that dissociates to form $\mathrm{A^-}$ and $\mathrm{H^+}$. Combining the above expression with (2),

In summary, (12) becomes:

BDW modelled $J_\mathrm{H^+}$ ($\mathrm{mol\cdot m^{-2}\cdot s^{-1}}$) as the additional proton-extrusion rate above the fixed background rate

where $k$ ($\mathrm{s^{-1}}$) is the proton-pumping rate constant, $(k/\rho)\mathrm{[H^+]_i}$ is the additional flux of $\mathrm{H^+}$ above the background $\mathrm{H^+}$ flux of $(k/\rho)\mathrm{[H^+]'_i}$, which occurs at the resting $\mathrm{[H^+]_i}$ of $\mathrm{[H^+]'_i}$ (i.e., resting $\mathrm{pH_i}$ of $\mathrm{{pH}'_i}$). Note that $k/\rho$ has units of ($\mathrm{m\cdot s^{-1}}$), consistent with the membrane permeability terms $P_\mathrm{HA}$ and $P_\mathrm{A^-}$ in (9) and (10).^[8]

The BDW authors used the definition of buffering power, in its infinitesimal form, to derive the relation between $\mathrm{d}{\mathrm{[H^+]_i}}/\mathrm{d}t$ and $\mathrm{d}Q/\mathrm{d}t$, as shown in the following steps.

Our first goal is to obtain an expression for $\mathrm{d}{\mathrm{pH}}/\mathrm{d}t$ in terms of $\mathrm{d}Q/\mathrm{d}t$. According to the chain rule:

By definition (see (8)), $\beta=-{\mathrm{d}Q}/{\mathrm{dpH}}$ ($\mathrm{mol\cdot m^{-3}}$), or equivalently

Combining (19) and (20)

Our next goal is to obtain an expression for $\dfrac{\mathrm{dpH}}{\mathrm{d}t}$ in terms of $\dfrac{\mathrm{d[H^+]_i}}{\mathrm{d}t}$. According to the chain rule:

By definition, $\mathrm{pH}={-\ln{\mathrm{[H^+]_i}}}/{2.303}$, so that:

Combining (22) and (23), we have

or equivalently,

Substituting (21) into (25), we obtain

Finally, substituting (17) into (26),

which is the second equation of the BDW model.

Substituting for $\alpha$ (from (14)), $[\mathrm{HA]_i} = \alpha [\mathrm{TA]_i}, [\mathrm{A^-]_i} = (1-\alpha) [\mathrm{TA]_i}$, in (11) and (27), we obtain the two ODEs of the BDW model in terms of $[\mathrm{TA]_i}$ and $\mathrm{[H^+]_i}$:

where $J_\mathrm{HA}$ (from (9)), and $J_\mathrm{A^-}$ (from (10)) are given by:

and $J_\mathrm{H^+}$ is given by (18).

The numerical solution of the above two equations yields the time courses of $[\mathrm{TA]_i}$ and $\mathrm{[H^+]_i}$, which in turn yield the time-courses of $[\mathrm{HA]_i}$ and $[\mathrm{A^-]_i}$ via:

where

Simulation for $\mathrm{CO_2}/\mathrm{HCO_3^-}$ experiments.¶

BDW employed (28) and (29) to simulate the experiments in which they exposed a squid giant axon to a solution containing equilibrated $\mathrm{CO_2}/\mathrm{HCO_3^-}$. Their simulation protocol was a step change in (a) extracellular $p_\mathrm{CO_2}$ from $0$ to $5\%$ $\mathrm{CO_2}$ ($37$ $\mathrm{mmHg}$ or, with $s=0.0321$ $\mathrm{mM\cdot mmHg^{-1}}$, $\mathrm{[CO_2]_o}=s.p_\mathrm{CO_2}=1.1877~\mathrm{mM}$) and (b) extracellular $\mathrm{[HCO_3^-]}$ from $0$ to $59.5260$ $\mathrm{mM}$ (the value that $\mathrm{[HCO_3^-]_o}$ has in a solution containing $5\%$ $\mathrm{CO_2}$ at $\mathrm{pH_o}$ of $7.70$)^[9]. The step change is applied for $2700$ $\mathrm{s}$ ($45~\mathrm{min}$) at constant $\mathrm{pH_o}=7.70$.

Table 1 and Table 2 report the parameter values used by BDW. Table 1 provides parameter values that are common to both the $\mathrm{CO_2}/\mathrm{HCO_3^-}$ and the $\mathrm{NH_3}/\mathrm{NH_4^+}$ experiments. Table 2 provides parameter values exclusive to the $\mathrm{CO_2}/\mathrm{HCO_3^-}$ experiments only.

In the present work, the differential (28) and (29) — when coded in CellML and solved with OpenCOR — produce the plots in Figure 4. The simulation file Boron-CO2.sedml contains the computational setting for running the model. Open the .sedml file in OpenCOR and click Run Simulation. The initial conditions are $[\mathrm{TA]_i}=0$ $\mathrm{mM}$ and $\mathrm{pH_i}=7.40$. Note that Figure 4 illustrates the time courses not only of $\mathrm{pH_i}$ — as presented by BDW — but also of quantities (e.g., various solute concentrations and fluxes) not displayed in the original paper; these values are useful for understanding the processes that contribute to the $\mathrm{pH_i}$ transient. Moreover, our curated and annotated version of the BDW model also allows one to alter the parameter values from those originally chosen by BDW, thereby extending the ability of the user to investigate the predictive power of the computational model.

4The Boron & De Weer Model for the Permeation by an Uncharged Weak Base and its Conjugate, Cationic Weak Acid¶

Following an approach analogous to the one outlined above for weak acids, BDW derived two time-dependent ODEs. The first describes the time-course of the concentration of total intracellular buffer ($[\mathrm{TB]_i} = [\mathrm{B]_i}+[\mathrm{BH^+]_i}$), and the other the time-course of the intracellular free $\mathrm{[H^+]_i}$, for any buffer pair $\mathrm{B}/\mathrm{BH^+}$.

Derivation for weak bases.¶

Imagine that a cell is exposed to a solution containing equilibrated $\mathrm{B}/\mathrm{BH^+}$, and that both $\mathrm{B}$ and $\mathrm{BH^+}$ initially move into the cell — because of the chemical gradient in the case of $\mathrm{B}$, and because of the electrochemical gradient in the case of $\mathrm{BH^+}$.

Assuming, as in Figure 5, that $\mathrm{[TB]_i}$ only can change due to the transmembrane fluxes of $\mathrm{B}$ ($J_\mathrm{B}$) and $\mathrm{BH^+}$ ($J_{\mathrm{BH^+}}$), the time rate of change of $\mathrm{[TB]_i]}$ — analogous to (11) above — is

where $\rho$ ($\mathrm{m^{-1}}$) is again the area-to-volume ratio for the cell. The equation

is an integrated form of Fick’s first law of diffusion that describes the net passive flux of B, and

describes the net passive influx of $\mathrm{BH^+}$ according to the GHK equation. In the previous two equations, $P_\mathrm{B}$ ($\mathrm{m\cdot s^{-1}}$) is the membrane permeability to the uncharged weak base $\mathrm{B}$, $P_\mathrm{{BH^+}}$ ($\mathrm{m\cdot s^{-1}}$) is the membrane permeability to the charged conjugate weak acid $\mathrm{BH^+}$, and $\epsilon '$ is a shorthand for $e^{{V_\mathrm{m}F}/{RT}}$. (33) is the first of two ODEs of the BDW model for the buffer pair $\mathrm{B}/\mathrm{BH^+}$.

The second equation of the BDW model for a weak base — analogous to (27) above — is

where $J_\mathrm{H+}$ is the same as in (18) and

and

Substituting for $\alpha$, $[\mathrm{BH^+]_i} = \alpha [\mathrm{TB]_i}, [\mathrm{B]_i} = (1-\alpha) [\mathrm{TB]_i}$, $J_\mathrm{BH^+}$ , $J_\mathrm{B}$ in (33) and (36), we obtain the two ODEs of the BDW model in terms of $[\mathrm{TB]_i}$ and $\mathrm{[H^+]_i}$

where

and $J_\mathrm{H^+}$ is given by (18).

Numerically integrating the above two equations yields the time courses of $[\mathrm{TB]_i}$ and $\mathrm{[H^+]_i}$, from which we can compute the time-courses of $[\mathrm{BH^+]_i}$ and $[\mathrm{B]_i}$ from