Up-down biphasic volume response of human red blood cells to PIEZO1 activation during capillary transits

In this paper we apply a novel JAVA version of a model on the homeostasis of human red blood cells (RBCs) to investigate the changes RBCs experience during single capillary transits. In the companion paper we apply a model extension to investigate the changes in RBC homeostasis over the approximately 200000 capillary transits during the ~120 days lifespan of the cells. These are topics inaccessible to direct experimentation but rendered mature for a computational modelling approach by the large body of recent and early experimental results which robustly constrain the range of parameter values and model outcomes, offering a unique opportunity for an in depth study of the mechanisms involved. Capillary transit times vary between 0.5 and 1.5s during which the red blood cells squeeze and deform in the capillary stream transiently opening stress-gated PIEZO1 channels allowing ion gradient dissipation and creating minuscule quantal changes in RBC ion contents and volume. Widely accepted views, based on the effects of experimental shear stress on human RBCs, suggested that quantal changes generated during capillary transits add up over time to develop the documented changes in RBC density and composition during their long circulatory lifespan, the quantal hypothesis. Applying the new red cell model (RCM) we investigated here the changes in homeostatic variables that may be expected during single capillary transits resulting from transient PIEZO1 channel activation. The predicted quantal volume changes were infinitesimal in magnitude, biphasic in nature, and essentially irreversible within inter-transit periods. A sub-second transient PIEZO1 activation triggered a sharp swelling peak followed by a much slower recovery period towards lower-than-baseline volumes. The peak response was caused by net CaCl2 and fluid gain via PIEZO1 channels driven by the steep electrochemical inward Ca2+ gradient. The ensuing dehydration followed a complex time-course with sequential, but partially overlapping contributions by KCl loss via Ca2+-activated Gardos channels, restorative Ca2+ extrusion by the plasma membrane calcium pump, and chloride efflux by the Jacobs-Steward mechanism. The change in relative cell volume predicted for single capillary transits was around 10−5, an infinitesimal volume change incompatible with a functional role in capillary flow. The biphasic response predicted by the RCM appears to conform to the quantal hypothesis, but whether its cumulative effects could account for the documented changes in density during RBC senescence required an investigation of the effects of myriad transits over the full four months circulatory lifespan of the cells, the subject of the next paper.

Osmotic equilibrium in the reference steady state: M Os = COs (8)

Cytoplasmic buffering of protons, calcium and magnesium.
Heamoglobin is the major cytoplasmic buffer for protons (eq 4) and for calcium (α-buffer in eq 9c). The main magnesium buffers are ATP and 2,3-DPG, compounds integrated within the X phenomenology. Because the bound forms of Ca and Mg are contained within CX − , they are not included as separate osmolarity contributors in eq 6, leaving only the free forms of Ca 2+ and M g 2+ as osmotic contributors. Cytoplamic Ca 2+ and M g 2+ buffering have been measured with precision in intact RBCs [5][6][7][8] enabling accurate representations in the model. The total Ca and Mg content of the cells, QCa and QMg, is reported in units of mmol/(340g Hb) (or mmol/Loc) whereas concentrations of the free forms, CCa 2+ and CM g 2+ , are expressed in units of mmol/Lcw, a conversion requiring translation for operational reasons in the model. Equation 9a translates QCa in units of mmol/Loc to CCa in units of mmol/Lcw using: The total calcium concentration is the sum of free and bound forms: There are two buffer systems for binding calcium in the RBC cytoplasm, α (mostly haemoglobin), and the BCa/KBCa buffer [6]. The concentration of bound calcium, CCaB, at each total calcium concentration, CCa, is represented by: CCa 2+ is solved from the implicit equation: by the Newton-Raphson routine in the RS and at the end of the computations in each iteration cycle. The measured values of the calcium binding parameters are α = 0.30, CB Ca = 0.026 mmol/Loc, and K BCa = 0.014 mM [6].
The corresponding equations for cytoplasmic magnesium buffering and CM g 2+ are: CM g 2+ is solved from the implicit equation: The measured values of the Mg bufferes [8] are: CB M g 1 = 1.2 mmol/Loc, KB M g 1 = 0.08 mM; CB M g 2 = 7.5 mmol/Loc (15 mEq/Loc), KB M g 2 = 3.6 mM; CB M g 3 = 0.05 mmol/Loc. B M g 1 represents ATP, BM g2 represents 2,3-DPG and miscellaneous phosphate groups, and BM g3 is an unidentified high affinity magnesium buffer.
1.1.6 Effects of deoxygenation on cytoplasmic M g 2+ buffering and pHi.
Deoxygenation increases haemoglobin binding of ATP and 2,3-DPG thus reducing their availability for buffering intracellular magnesium. CB M g 1 is reduced by half and CB M g 2 by 1.7 [8]. This is particularly relevant for simulating accurately the effects of changing the oxygenation condition of RBCs, a process in which changes in CM g 2+ become enmeshed with effects arising from changes in the isoelectric point of haemoglobin, pI (eq 4). In vivo, RBCs are continuously changing between oxy and deoxy states as they flow between the arterial and venous vasculature causing well documented alternating changes in pHi, CM g 2+ , CA and cell volume which the model accurately reproduces. Although these changes are fully reversible in physiological conditions, deoxygenation of sickle RBCs can lead to hyperdense collapse, as shown in the next paper [2].
Hb is assumed to be in a oxy-state by default, the most frequent experimental condition. Deoxygenation of Hb (Deoxy) changes its pI(0oC) from 7.2 to 7.5. The model automatically adjusts the actual pI change for the temperature of the experiment. The pI shifts during oxy-deoxy transitions cause sudden changes in the protonization condition of Hb with secondary changes in pHi and CM g 2+ , changes which the model predicts with verified accuracy [9][10][11]. Electroneutrality preservation during oxy-deoxy transitions requires constancy of nHb values (eq 4) when pI changes, from which the compensatory changes in pHi can be derived according to [10]: On deoxygenation: On reoxygenation: pHioxy = pHideoxy + pIoxy − pIdeoxy (4.b)

The dynamic state (DS)
A first requirement at the start of simulations is to define the relative volume occupied by cells in the cell suspension system, the cell volume fraction, CVF. Perturbations alter the flux of transported solutes and water across the plasma membrane of the cell thus initiating a cascade of downstream changes in the compositions of cell and suspending medium. It is therefore important to start by listing the membrane transport component of the cell and of the equations describing their basic kinetic properties.

Flux equations of the model, F i and F j
All flux equations are defined as products between permeabilities (P) or rate constants (k) and driving forces. The substrates of the RBC membrane transporters are Na, K, A, H, Ca, Mg and water, the "i" in F i . The sign-convention applied in the equations is for positive fluxes into the cell (influx) and for negative fluxes into the medium (efflux). The name convention adopted here for the transport of substrate X by the different membrane transporters is as follows: FPX = pump-mediated flux of X, with P = NaP for the Na/K pump or CaP for the calcium pump (PMCA); FGX = X-flux through electrodiffusional channel defined with constant field kinetics; FXA = electroneutral carrier-mediated cotransport of cation X and anion A defined with low-saturation kinetics; FzX = electrodiffusional flux of X through PIEZO1 channel; FCoX = electroneutral cotransport of X mediated by the Na:K:2Cl symport, of minimal expression and activity in human RBCs; FA23X = electroneutral M 2+ : 2H + exchange flux through the divalent cation ionophore A23187, the only exogenous membrane transporter included in the model; Fw = water flux mediated mainly by aquaporins and partly by partition diffusion through the plasma membrane.
1.2.2 Flux pathways for each transported substrate, F i : There are no data on PIEZO1-mediated M g 2+ fluxes in RBCs. Although PzMg most certainly has a small finite value, F z M g is likely to be very small under the usually low electrochemical M g 2+ gradients across the RBC membrane. With this level of uncertainty, F z M g was not included in the current model version.

Kinetic descriptions of individual transporters
Certain transporter kinetics are reported in the equations with the default numerical values used for dissociation and rate constants in the model, based on well established values in the literature and on the good semiquantitative fits to experimental data provided in the past [12][13][14]. The default values of all permeabilities and rate-constants used in the model (P or k) correspond to experimentally measured values at 37oC, P = P (37) or k = k(37). Temperature-changed values of P or k, P = P (T ) or k(T ), are computed relative to P (37) or k(37) using the Q10-derived formalism for temperature coefficients. In eq (11) we use P to represent both Por k-defined values: 1.2.4 Na/K pump mediated fluxes of Na and K (f = forward; r = reverse) [15,16] 1.2.6 Electrodiffusional fluxes of i (Na, K, Ca, H and A) through endogenous channels, FGi , Gardos channels, FGGardos, and PIEZO1 channels, Fzi, are represented with constant field kinetics [19]: with PGi representing the Goldmanian i-permeability in h −1 units
Note that k HA , the rate constant of the H:A cotransport phenomenology representing the operation of the Jacob-Stewart mechanism (JS) is between five and six orders of magnitude faster than that of any of the other ion transporters in the membrane (see User Guide for details and references).

Electroneutral Na:K:2A cotransport
The CX and MX values in eq 17b are those set for the RS d is a wildcard factor introduced to set F Co = 0 only in the RS. Its value is set by the initial Na, K and A concentrations in the RS. d remains as a fixed-value parameter during dynamic state computations.

1.2.11
Electroneutral M 2+ : 2H + exchange fluxes of Ca 2+ and M g 2+ mediated by the divalent cation ionophore A23187 The divalent cation ionophore A23187 mediates an electroneutral M 2+ : 2H + exchange when incorporated into cell membranes [23]. Divalent cation ionophores became essential and extensively used tools in research on calcium and magnesium function and dysfunction in RBCs [5,7,24,25] and in many other cell types. To emulate experimental protocols with the use of divalent cation ionophores it became necessary to represent their transport properties in the model as an optional exogenous transporter of the RBC membrane.
In albumin-free RBC suspensions, the RBC/medium partition ratio of the lipophilic ionophore A23187 is 60/1, 20 to 50% of it confined to the cell membrane [26]. The transport kinetics of the ionophore was modeled with symmetric binding (Km) and inhibitory (KI) dissociation constants for Ca 2+ and M g 2+ on each membrane side, as follows: Following extensive preliminary tests [27], default values of 10 mM for the four Km and KI parameter set were found to deliver excellent agreement between predicted and measured ionophore-mediated fluxes, and to ensure adequate compliance with the measured equilibrium distribution of the transported ions when ionophoremediated net fluxes approach zero [23]: Combining the Ca 2+ , M g 2+ and H + driving gradients we obtain: The ionophore-mediated fluxes of Ca 2+ , M g 2+ and H + , F A23 Ca, F A23 M g and F A23 H, respectively, can now be computed from: Where P A23 is the ionophore-mediated permeability. P A23 is a power function of the RBC ionophore concentration, P A23 = 0.22 × [I] 1.45 , when P A23 is expressed in units of 10 −6 cm/s, and [I] in µmol/Loc [26,28,29]. Within the units-set in the model, numerical values of P A23 in the range 10 17 to 2 × 10 18 offered a perfectly adequate minimalist emulation of the effects of different ionophore concentrations on the fluxes and distributions of Ca 2+ , M g 2+ and H + ions in RBCs in a large variety of experimental conditions [5,27,[30][31][32][33].

Equation sequence for the computations of dynamic states.
Following perturbations, sustained charge conservation and electroneutrality is implemented by: where I j represents the current carried by each of the j-membrane transporters. I j is therefore the fist equation that has to be solved at the start of each iteration in the computational sequence of dynamic states. Capacitative currents (Ic = C(dV /dt)) are ignored because their magnitude and time-course decay are orders of magnitude smaller than those of the homeostatic relevant currents. The relation between currents and fluxes, F j , is given by With the electrogenic flux components in the model (z j = 0), I j = 0 renders: With Em, the new z j F j values for each of the electrodiffusional terms in eq 18c can be computed. We can now add up the absolute values of the new computed fluxes to the values of the electroneutral fluxes in the previous iteration |F j | to assign a new ∆t duration to each iteration interval, as follows: The value of a, under user control, optimises ∆t scales for different simulations ("frequencyfactor" in the RCM); b is a small zero-avoidance parameter in the denominator. The advantage of this strategy over using regular iteration intervals is that by setting a constant value for the cycles per outcome ("cyclesperprint(epochs)" in the RCM) the density of data output points automatically adjusts to the overall rate of change in the system, emulating the way good experimental practice seeks to sample for data at the bench, thus optimizing comparisons between predicted and experimental results.
With the new F t i and ∆t the new Qi t may be computed using the values of F N a t , F K t , F A t , F H t , F Ca t and F A23 M g t from equations (10a-f) as follows: ∆H is a special case because ∆H adds to the only titratable proton buffer nHb × QHb, so that: From which we can now compute the new cell pHi from eq 4 by solving for pH t : The new intracellular H + concentration in molar units is: With the new Qi t , we need the new cell water volume, V w t in order to compute the new cell concentrations, Ci t = Qi t /V w t . The water flux across the RBC membrane, F w, is driven by the osmotic gradient across the RBC membrane (eqs 5 and 6): COs t can be computed from the altered osmotic load resulting from the ∆Qi changes during ∆t operating on the cell volume at the start of the each iteration interval: The new cell water volume, V w t , and volume-associated variables, RCV t , M CHC t , Density t and Hct t , can now be computed from: With V w t we proceed to compute next the new intracellular concentrations of Na, K, A, H, Ca, Hb, and X: The new osmotic coefficient of Hb, f t Hb , can now be calculated from eq 7 and the new CHb t : 1.2.13 Computation of the medium concentrations at time = t.
Medium concentration changes arise from independent solute and water transfers between cells and medium under mass conservation. At constant suspension volume, water transfers between cells and medium generate self-compensating changes in cell and medium volume fractions, CVF and (1-CVF), respectively, according to: By mass conservation, the Qi changes during ∆t, ∆Qi, are transferred to the medium, ∆Qim, so that: ∆Qim + ∆Qi = 0 (26b) ∆Qim can be expressed in terms of M i changes during ∆t as follows: Replacing ∆Qim by −∆Qi (eq 26b) in equation 26c and solving for M i t , we obtain: With eq 26d we can now compute the new medium concentrations at time = t for transported solutes, eqs 27a-f, and for impermeant solutes (∆Qi = 0) whose concentration changes only because of water shifts, eqs 27g-j: