A modular and reusable model of epithelial transport in the proximal convoluted tubule

We review a collection of published renal epithelial transport models, from which we build a consistent and reusable mathematical model able to reproduce many observations and predictions from the literature. The flexible modular model we present here can be adapted to specific configurations of epithelial transport, and in this work we focus on transport in the proximal convoluted tubule of the renal nephron. Our mathematical model of the epithelial proximal convoluted tubule describes the cellular and subcellular mechanisms of the transporters, intracellular buffering, solute fluxes, and other processes. We provide free and open access to the Python implementation to ensure our multiscale proximal tubule model is accessible; enabling the reader to explore the model through setting their own simulations, reproducibility tests, and sensitivity analyses.


I. Introduction
Kidneys are vital organs and play an essential role in the overall homeostasis of the body in mammals. A human kidney is typically composed of one million nephrons, which are the primary functional unit of the kidney. The nephron consists of the renal corpuscle and the renal tubule. A renal tubule is a tubular structure composed of a single layer of epithelial cells divided into various functional segments. Each segment of the nephron has specific functions in the regulation of blood and urine composition. The proximal convoluted tubule (PCT) is considered one of the most significant functional segments in the nephron and a key contributor to pathologies such as hypertension and diabetes.
To gain a deeper understanding of the mechanisms and investigate any hypothesis regarding the underlying physiopathological conditions, such as hypertension, diabetes, or other kidney diseases, a virtual nephron model is an invaluable tool. A model like this can serve as a virtual laboratory, ideally be inexpensive to run, and should target the minimisation of animal experiments.
Many models of epithelial transport have been published that would be relevant and useful to integrate into a virtual nephron model, but often there is insufficient information in the literature to enable readers to reproduce the published observations and predictions, thus making their reuse in novel studies time consuming and resource intensive. Furthermore, as models evolve over time, any given version of a model is usually designed to investigate specific hypotheses. The various instances of a given model, or a family of models, are therefore inconsistent and require the reader to search the literature to discover or infer the modifications required to integrate the models and successfully reproduce published results. Such modifications, for example, may be as trivial as changes in parameter values or alterations of physical units, or as complex as specific assumptions made or mathematical formulations chosen. These problems hamper efforts to develop a virtual nephron model. To help address this problem, we introduce here an integrative mathematical model of the PCT that reproduces the capabilities of existing renal models. We used a modular approach to build our model, which we believe will improve reusability while ensuring that the segment function predicted by the source models can be reproduced. At the same time, this approach is sufficiently flexible and configurable to be able to adapt to different segments or epithelia. This model has been implemented in Python and is freely available under an open-source and permissive license at https://github.com/iNephron/W-PCT-E. Our long-term goal is to make this resource available in a more reproducible and reusable form via the CellML standard [1,2] and the bond graph approach [3][4][5][6][7][8] that ensures energy conservation as well as mass and charge conservation across different physical systems (biochemical, electrical, mechanical and metabolic). However, to bring this work into a consistent and unified mathematical framework and to understand the various parameter combinations and manipulations that have been used in a series of publications spanning 30+ years, we first implement in a Python scripting environment with clearly defined parameters and consistent physical units throughout.
This work would not have been possible without the comprehensive epithelial and tubular modelling published by Alan Weinstein, and colleagues, examining many aspects of renal function [9-16, for example]. Our implementation is derived from information available in the literature, supplemented with additional equations and parameter values discovered in Weinstein's original implementation available at https://github.com/amweins/kidney-modelsamw and also the work of other researchers who have adopted and confirmed the Weinstein model, such as [17][18][19] (with accompanying implementation available at https://github.com/ Layton-Lab/nephron). As mentioned above, our goal here is to make this resource available to the renal modelling community in a manner that enhances reproducibility and reusability.

II. Design and implementation
In this work, we follow the comprehensive PCT epithelial model collectively presented in Weinstein et al. [9][10][11][12][13] which we refer to here as the W-PCT-E. The W-PCT-E consists of four logical blocks, see Fig 1. The first one provides the geometrical definitions of the system and the equations of the selected components, such as cell volume, lateral intercellular volume, basement membrane area, and epithelial thickness. These parameters are variable, whereas the rest of the geometrical components, such as membrane area for the other regions, are constant (see Section vii). The second block defines five different intra-epithelial fluxes: water fluxes, convective fluxes, passive fluxes, coupled solute fluxes, and active (metabolically dependent) fluxes. The output from the first block is coupled to the input of the second block, plus the specification of the membrane types and solutes that appear in the chemical reactions. The outcome of the second block are the total membrane solute fluxes and membrane water fluxes.
In the third block, the W-PCT-E model enforces the mass conservation of the system through differential equations within each compartment (or each membrane) with the total solute fluxes and water fluxes as input for each membrane.
In the final block, the W-PCT-E model applys more constraints by defining buffer pairs, pH equilibrium, and electroneutrality of the system.

III. Model illustrations
Here, we define the various components in the modelling framework in the W-PCT-E system, including all equations, definitions, and assorted tables of constants and variables that appear in the epithelial model's compartments. In addition, we introduce the ingredients present in our Python code and, whenever possible, the provenance of the parameter values. The comprehensive W-PCT-E model consists of cellular and lateral intercellular compartments between luminal and peritubular solutions.  The first block indicates the geometrical inputs (see Section vii). The second block shows the intra-epithelial fluxes; these fluxes will be involved in the mass conservation equations. The third block indicates the total mass conservation differential equations present in the epithelial model. The differential equations are updated by the system's buffer pairs, pH equilibrium and electroneutrality to add more constraints to the system. The sequential placement of the blocks is not related to a procedural description of the model simulation but helps to understand the model construction. membrane, MI), tight junction (ME), cell-lateral membrane (IE), interspace basement membrane (ES), or cell-basal membrane (IS). The order of the two letters indicates the positive direction of the mass flow. J αβ and J ναβ represent the solute flux and water flux, respectively, through the corresponding membrane; A is the corresponding membrane surface area; V is the volume; E is the trans-epithelial potential difference. Symbols are defined in the following sections as introduced by Weinstein et al. [12] in the epithelial PCT model. Intra-epithelial fluxes are designated J αβ (i), where αβ refers to the different membranes. Models can include many different solutes in various compartments. In this paper, according to [12], the model consists of 15  , H 2 CO 2 , and glucose, as well as two impermeant species within the system; a nonreactive anion and a cytosolic buffer. The solutes considered in a specific simulation experiment can vary, leading to the dynamics of some of the solutes being ignored under certain conditions. There are 14 transporters (symporters, antiporters, complex transporters, and ATPases) that produce electrochemical fluxes in the W-PCT-E model.
Overall, there are six symporters on different membranes; two on the lumen-cell (MI): Sodium-Glucose (SGLT) and Sodium-Phosphate (NPT1), two on the cell-basal (IS): Potassium-Chloride (KCC4) and Sodium-Bicarbonate (NBC1), and two on the cell-lateral (IE): Potassium-Chloride (KCC4) and Sodium-Bicarbonate (NBC1). It is important to note that the cell-basal and cell-lateral membranes share a similar layout.
There are two complex exchangers on the cell-lateral membrane: Sodium-Bicarbonate/ Chloride (NCBE) and one on the cell-basal membrane: Sodium-Bicarbonate/Chloride (NCBE). The Na/K pumps (NaK-ATPase) are located in both cell-basal and cell-lateral membranes. H-pumps (H-ATPase) are located in the lumen-cell membrane, regulating the pH in the W-PCT-E model. We noticed two different approaches to translate the behaviour of sodium/hydrogen antiporter in Weinstein's work. In the first approach, Weinstein et al. used the equivalent definition of two antiporters: Sodium/Hydrogen and Sodium/Ammonium [11]. In the second approach, reported in [12], the authors used a detailed model of the Na + /H + antiporter from [20].

i. State equations
First, we introduce the sets of compartments (C) and membranes (M) indexes Indexes in M might appear interchanged in some flux definitions, which means that the flux sign must be inverted. Now, we introduce the following set of solutes which is divided into non-reacting (NR) and reacting (R) solutes In the W-PCT-E model, the system of state equations represents two different compartments: the cell and the lateral intercellular spaces. To formulate the mass conservation equations within each compartment, the net generation of each species S α (i) is defined as an intermediate variable within the compartment [11]. The generation of multiple reacting solutes is the sum of the net exchange of the flux plus the accumulation of that solute in each compartment, that is where S I (i) and S E (i) indicate the generation for solute i 2 S in the cell space and lateral interspace compartments, respectively. Within the epithelium, the flux of solute i across the membrane αβ is denoted as J αβ (i) [mmol. s −1 . cm −2 ], ab 2 M and V α is the compartment volume per unit surface area [cm 3 /cm 2 . epithelium], a 2 C. The principle of mass conservation relating water fluxes and volume change for different compartments reads where J vαβ [ml. s −1 . cm −2 ] is denoted as the transmembrane volume flux. It is important to mention that for nonreacting solutes we have The mass conservation then defines the change of the concentration of the i 2 S species in the intracellular solution as the transport of solute i into and out of the cell through the apical and basolateral membrane. This is a direct transport of solutes through the membrane, and a contribution from convective transport due to the flow of water through the membranes. For such a nonreacting solute, the combination of (6) and (8) yields This equation holds for each solute i 2 S being considered in a particular instantiation of the model. If required for a particular model, similar equations can be introduced for the solute concentrations in the mucosal and/or serosal solutions. Consequently, from (8), the conservation of cellular water yields the equation below, with the rate of change of cell volume, V I , defined as where each of the total membrane water fluxes, J vαβ , ab 2 M, is scaled (multiplied) by its respective membrane area to take into account the averaged behaviour of the representative membrane. For more detailed information of all various fluxes across the membrane, see Section IV.
To calculate the concentration of the cellular buffer pairs, the following expression is employed [11] ðC Buf À þ C HBuf where C Buf À and C HBuf are cell buffer and protonated cell buffer concentration, respectively. V I and V 0I indicate the cell volume and cell volume reference, respectively. To add more constraints on the proton cellular concentration, and to generate the paired equation, we have the following equilibrium equations where pK HPO 4 2À and pK Buf À are the equilibrium constants for the phosphate and cell buffer pairs. Considering that the concentration of the H + remains constant for all buffer pairs in the model, expressions (15) and (16) are combined to give To include cellular proton concentration in the conservation of proton mass, S I (H + ), the following modification is applied to the Eq (6) which can equivalently be written in the following form

ii. Buffer pairs and pH equilibrium
The W-PCT-E model defines different types of buffer pairs, the mass conservation principle for the phosphate and formate buffer pairs takes the form Similar equations apply to the ammonia pair within the lateral interspace and cell [11], that is For an impermeant cytosolic buffer, Buf − , it is All buffer species are assumed to be at chemical equilibrium. Within each compartment, there are four additional pH equilibrium relations, corresponding to the four buffer pairs. The algebraic relations of the model include the pH equilibria of four buffer pairs, pH ¼ pK þ log 10 Base À HBase The collection of chosen buffer pairs and even the definition of mass conservation equation can be different over various studies depending on the specific focus of the study. As an example of the variability in these equations, one can compare the presentation of the mass conservation equation for NH 3 :NH 4 + buffer pairs within a cell. In Weinstein's reports from 1992 and 2009, see [11,21], the authors introduced Eq (22) to define the mass conservation equation for these buffer pairs within a cell, while in [12] they used the following equation where Q I (NH 4 + ) is defined as an ammoniagenesis factor, see [12]. According to [11,21], conservation of charge among the buffer reactions requires that We should mention that in [12], the authors employed a modified version of expression (26) as the requirement of the charge conservation among the cellular buffer reaction in the following form In turn, the charge conservation for the lateral interspace (E) buffer reaction stays unchanged, as given by (26). Although peritubular PCO 2 is specified, the CO 2 concentrations in cell and interspace are model variables. The mass conservation for HCO 3 − , H 2 CO 3 , CO 2 is expressed as below where k h , k d are the hydration and dehydration rates for CO 2 , respectively. The buffering described above and implemented in our unified W-PCT-E model provides a flexible and modular method to explore various features of the system. This covers the broadest range of buffering included in the Weinstein models while enabling replication of specific instances where different approaches are used. For example [12], in which the impermeant buffer conservation equation was omitted, resulting in the expression (22) being replaced by (25) and for the expression (26) is replaced by (27).

iii. Electroneutrality constraints
When considering the movement of charged solutes, with valence Z i , i 2 S, the proposed system of equations is not sufficient to guarantee that the cell and interspace remain electrically neutral. An electroneutrality relation for the cell compartment is expressed through the following balance equation X i2S Z i C I ðiÞ þ Z I;Imp C I;Imp À C Buf À ¼ 0; where C Imp and C Buf À denote the concentration of cell impermeant solute and cell unprotonated buffer, Z I,Imp is the cell impermeant valence, and Z i is the valence of species i, i 2 S. The law of electroneutrality for the interspace is defined as follows X and for all of the buffer reactions, there is conservation of protons, which implies that the fol- The electroneutrality condition is effectively prescribed by considering that the net charge fluxes into and out of the epithelium are the same. The membrane charge fluxes can be represented as electrical currents using the following relationships where F is Faraday's constant. Balancing the flow of charge into and out of the epithelium therefore results in which must hold true at all times. In the current work, according to [12], expression (3) is only considered for the balance of charge transfer across the membranes such that I In = 0. Eqs (26)-(34) are a collection of different electroneutrality constraints that were introduced in different studies (e.g., [11,12,21]). However, it is important to note that not all these equations were utilised in all different studies. Instead, they were selectively chosen based on the context of each study.

IV. Model specialisation
The basic principles of mass conservation, pH equilibrium of buffer species, and maintenance of electroneutrality described above apply to epithelial transport in general. To instantiate the general model into a mathematical model for a specific epithelium, all that remains is to define the actual membrane solute and water fluxes of interest to create the specialised model.

i. Water fluxes
With respect to water fluxes, the volume conservation equations for lateral interspace and cell are considered to compute the lateral interspace hydrostatic pressure, and cell volume. Across each cell membrane, the transmembrane water fluxes are proportional to the hydrostatic, oncotic, and osmotic driving forces where P α and π α are the hydrostatic and oncotic pressures within compartment a 2 C, L pαβ is the membrane water permeability and σ αβ (i) is the reflection coefficient of membrane ab 2 M to solute i 2 S, and R and T are the gas constant and absolute temperature, respectively. It is important to mention that all L pαβ constants which are represented in the implementation of this model (see the Python code) are multiplied by R and T. The reflection coefficients, σ αβ (i), stay identical in most of Weinstein's body of work, see Table 1 in S1 File. In turn, there are some variations in the model parameters (such as coupled transporter coefficients, cell membrane water permeability, and cell membrane solute permeability) across the different works.

ii. Convective solute fluxes
In the model proposed in [11], it is assumed that there are convective fluxes for all intraepithelial solutes. Taking into account the definition of the water fluxes, see (36), the convective fluxes are defined as follows [12] where � C ab ðiÞ is the logarithmic mean membrane solute concentration described by the expression Studying the reflection coefficient values σ αβ (i) (defined as membrane-solute properties), one can see that the reflection coefficient is mostly one in lumen-cell (MI), cell-lateral (IE) and cell basal membrane (IS). The most effective membrane to produce the convective fluxes is the interspace basement membrane (ES) with the reflection coefficient mostly zero, and then in the second place it is the tight junction (ME).

iii. Passive solute fluxes
In the W-PCT-E model, passive solute fluxes across all membranes are assumed to occur by electrodiffusion and to conform to the Goldman-Hodgkin-Katz constant-field flux equation [22]. Passive solute fluxes of the species i 2 S across the membrane ab 2 M are given by where, for solute i 2 S and for membrane ab 2 M, h αβ (i) [cm.s −1 ] is the permeability and z αβ (i) is a normalised electrical potential difference (dimensionless), defined by where ψ α and ψ β are electrical potentials within compartments α and β, respectively. F is Faraday's constant, R and T are the gas constant and absolute temperature, sequentially. Weinstein et al. did not hold on to one solid definition for the permeability, in some cases permeability was multiplied by the area of the corresponding membrane h αβ (i)A αβ [10 −5 cm 3 .s −1 .cm −2 . epithelium] as an example see Weinstein et al. [12]. For the uncoupled permeation of neutral solutes across membranes, the Fick law is utilised

iv. Coupled solute fluxes
Coupled solute fluxes in the W-PCT-E model include three different categories of transporters: simple cotransporters, simple exchangers, and complex exchangers. All coupled solute transporters in this model have been represented according to linear nonequilibrium thermodynamics, so that the solute permeation rates are proportional to the electrochemical driving force of the aggregate species, with a single permeation coefficient. Simple cotransporters consist of peritubular K + -Cl − , luminal Na + -Gluc and Na + -H 2 PO 4 − , which are in the form of the following equations It is important to mention that all transporters within cell-basal (IS) membrane are also considered for the cell-lateral (IE) membrane. In the equations above, the fluxes of two different species across the cotransporter are equal (1: 1 stoichiometry).
Simple exchangers such as Na + /H + , Na + /NH 4 + , Cl − /HCO 2 − , and Cl − /HCO 3 − are located at the lumen-cell (MI) membrane, and represented by the equations below In [11], the authors introduced the NHE3 exchanger in the luminal membrane through Eqs (45) and (46). However, the NHE3 exchanger has been developed through the kinetic model proposed by Weinstein et al. in 1995 [20]. In the current work, the NHE3 exchanger is modelled by employing the mathematical system introduced in [20]. There are also two more complex transporters at the peritubular membrane: Na + -HCO 3 − and Na

PLOS ONE
are defined by the following equations (see [11]) In the expressions above, L (i,j) is the transporter coupling coefficient, a single proportionality constant which specifies the relative activity of the transporters. In turn, � m ab ðiÞ is the electrochemical potential difference of species i across all membranes Here, we also define the following quantity It is important to mention that all permeability coefficients, L (i,j) , which are represented in Table 1, from [11], are scaled (multiplied) by their respective membrane area to take into account the effective behaviour of the representative membrane.

v. Active solute fluxes
In the W-PCT-E model, there are two ATPases, the apical membrane H + -ATPase and a peritubular Na + /K + -ATPase. To model the H + -ATPase, an expression of the following form is utilised Note that the rate of proton pumping varies as a function of the transmembrane electrochemical potential difference. Here, [J MI (H + )] max is a maximal rate of transport, and � MI is a steepness coefficient. The Na + /K + − ATPase exchanges three cytosolic Na + ions for two peritubular cations, K + or NH 4 + , in a way that competes for binding. The following expressions represent

PLOS ONE
all three different fluxes due to the Na + /K + − ATPase activities in which K Na þ , the half-maximal of Na + concentration, which scales linearly with the cellular concentration of K + ; K K þ , which is the half-maximal of K + concentration, rises linearly with the external concentration of K + . The pump flux of K + plus NH 4 + reflects the 3:2 stoichiometry. Similar expressions are considered for active transport at the cell-lateral membrane (IE), denoted by J A IE . Fig 3 illustrates the proximal PCT cell featuring coupled transport pathways and ion channels within luminal and basolateral membranes, the coupled transport pathways within the cell-lateral are not included.

vi. Total membrane solute fluxes
In summary, in the W-PCT-E model, the intraepithelial solute transport through membrane ab 2 M results from the contribution of convective flux J C ab , passive flux J P ab , electrodiffusive coupled flux J EC ab , and/or metabolically driven flux J A ab , that is where J αβ (i)[mmol.s.cm −2 ] is the total flux of i solute flow from the compartment α to the compartment β (through membrane αβ).

vii. Model compliance parameters
In the W-PCT-E model, the cell is compliant in a manner that there is no hydrostatic pressure difference between the cell and lumen, therefore, we have P I = P M . There is a substantial oncotic force within the cell, π I , that increases with decrements in the cell volume. Here, it is assumed that the total cell protein content C I,Imp V 0I is fixed and that π I is proportional to C I, Imp , for this reason C I,Imp replaces π I as one of the model unknowns (for more information see [11,12]) The lateral interspace volume (V E ) and its basement membrane area (A ES ) are functions of

PLOS ONE
interspace hydrostatic pressure, P E , that is where A ES0 and V E0 are reference values for the outlet area and volume, ν A and ν V are the compliance parameters used to ensure that the suitable pressures are maintained in the system (see Table 1 in S1 File). Table 1 in S1 File includes the model parameters and their definitions either in the current document or in the Python code. You can see the geometric parameters for all compartments in the tables reported below. As can be seen, the luminal and peritubular cell membranes have equal areas, i.e., 36 [cm 2 /cm 2 . epithelium]. The lateral interspace is compliant and distends

PLOS ONE
with transport-associated increments in interspace pressure. In this model, the tight junction properties are fixed and do not vary with transjunctional pressure differences.

V. Model calculations
In this work, the W-PCT-E model is bathed on both luminal and peritubular sides by solutions of equal concentration. Baseline bath and lumen conditions are those reported in Table 1 in S1 File. The choices for the model parameters appear in Table 2 in S1 File. The geometric parameters are completely different from [11], A IE = 35 [cm 2 /cm 2 . epithelium], A IM = 36 [cm 2 /cm 2 . epithelium]. As we mentioned earlier, the parameter values used in the current work are derived from information available in the literature and, in some cases, discovered in Weinstein's original implementation [https://github.com/amweins/kidney-models-amw/blob/ master/epithelia/pct/param.tem].
In the LIS interspace, the LIS basement membrane area and volume are compliant and distend with transport-associated increments in interspace pressure. However, the membrane areas in the cell are fixed and do not vary with transjunctional pressure differences. For a fixed quantity of cellular impermeant solutes, increasing the cell volume will decrease the impermeant concentration. The suitability of these parameters was not tested here. For the W-PCT-E simulations, the 35 nonlinear ordinary differential equations are solved using a finite difference numerical method for time discretisation along the Python solver "scipy.optimize.root". Evaluation of the model involves integrating the mass conservation equations from an initial time to a final time using small time increments. The simulation time is chosen to ensure that a steady-state regime is reached. Here, we should highlight that the solver also provides the eigenvalues of the W-PCT-E model. Studying the eigenvalues demonstrates that the W-PCT-E model is strongly stable with negative real parts and zero imaginary parts (for more information, see [23]). To analyse the W-PCT-E model stability, we extracted the Jacobian matrices from the results provided by Python solver "scipy.optimize.root". Then, we applied the built-in function "scipy.linalg.eig" to calculate the eigenvalues. Our stability results represent a stable system.
Here, there is a list of all the variables in the model: First, all variables appear in the lateral interspace compartment,

VI. Results
We investigated the validity of the WPCT-E model by designing several experiments; the analyses are performed over the steady-state solutions found from numerical simulations. To test the robustness of the W-PCT-E model, we investigated the sensitivity of steady-state solutions to different sets of initial conditions or time steps. Furthermore, we explored reproducibility by replicating some simulation experiments reported in [11,12] using the W-PCT-E. We then investigated the sensitivity to salt intake or luminal salt concentration in the W-PCT-E model based on earlier work [24]. Structural analyses were performed by inhibiting key transporters in different membranes, such as the Na + /K + -ATPase in the peritubular membrane or SGLT, NHE3, and Na + -H 2 PO 4 transporters in the apical membrane, and relating the predicted responses to observed biological phenomena.

i. Model sensitivity analysis
A time step size no greater than Δt = 0.1 s is required to ensure a converged numerical solution for epithelial transport models. In our simulations, we considered a time-step of Δt = 0.1 s.
In exploring the sensitivity of our W-PCT-E model to the initial conditions, we found that as long as the initial conditions were within a reasonable physiological range that the steadystate solution was insensitive to the initial conditions (to check the initial conditions, see the https://github.com/iNephron/W-PCT-E). However, outside that range, the results of the model are highly variable. To ensure consistent behaviour across the range of simulation experiments presented here, we use the same initial conditions. To determine that initial state, we performed a simulation experiment whereby we disabled all active transporters and allowed our W-PCT-E model to reach a steady state using only the passive processes. This steady state was then used as the initial state for all subsequent simulation experiments. When a given simulation experiment requires the addition of active processes to the model, we follow a similar protocol starting from the passive steady state, introducing the active process(es), and allowing the system to reach a new steady state before introducing further perturbation. This ensures that the simulations then begin at a suitable initial state in which the predicted steady state solution is insensitive. Applying the above approach to define the initial conditions does not affect our simulation results or the behaviour of the model. Even when the cell volume substantially increases because of the process of defining the initial conditions, the procedure yields stable numerical simulations.

ii. Model reproducibility
The present W-PCT-E model is built from a collection of mathematical representations reported in the literature. Our efforts have been to compile model parameters and equations from different scientific studies in which the components of PCT have been reported. Moreover, we provide the community with a freely available implementation to speed up research. In doing so, we aim to ensure that the behaviour of the W-PCT-E model reproduces that of the source model. Here, we provide exemplars exhibiting the flexibility and reproducibility of the W-PCT-E model.
Flow-dependent transport in the PCT. The mathematical model of the rat proximal tubule [12] was designed to include the calculation of microvillous torque and to incorporate torque-dependent solute transport in a compliant tubule. Here, we aim to reproduce some of the results reported in [12] by tuning the parameters according to [12] (see Tables 1-2 in that article). For those constant parameters or boundary conditions which were not defined in [12], we infer them from earlier works, specifically [10,11]. Reproducing these simulation experiments, we found no significant discrepancies, as shown in Table 1, for steady-state solutions such as solute concentrations, luminal voltage, and luminal pressure.
Chloride transport in the PCT. Here, we aim to reproduce some of the results from chloride transport in the proximal tubule [11] to explore the possible interactions between the individual transporter pathways and their contribution to overall chloride reabsorption in a proximal tubule. At the apical membrane, Cl − /HCO 2 − and Cl − /HCO 3 − exchangers are the main pathways for Cl − entry, and across the peritubular membrane Cl − /2HCO 3 − -Na + and Cl − -K + for the exit of Cl − . At the tight junction, chloride fluxes are both diffusive and convective. We performed the same experiments as reported in [11] and observed that the W-PCT-E model (with parameters modified accordingly) predicted similar behaviour as [11]. tions for the case of the small H 2 CO 2 permeability, in [11] the authors assumed that all luminal Cl − entry is through Cl − /HCO 2 − exchange. The W-PCT-E demonstrates that small luminal and peritubular cell membrane permeability for H 2 CO 2 could not sustain any substantial luminal Cl − flux. These results confirm previously published findings according to [11,25]. We highlight that the ordinate values in Figure 4 are 2-fold greater than Weinstein's Figure 8 [11]; Weinstein represented a tubule model while here we are describing a cellular epithelial model. To maintain the overall membrane reabsorption, the transporter parameters are multiplied by a factor of two (the transport per unit length needed to be constant).
Salt sensitivity. Next, we test the flexibility, reusability, and reproducibility of the W-PCT-E model by reproducing a simple model of Na + transport in the mammalian urinary bladder to study the salt sensitivity [24]. Here, we aim to design the similar experiment as in Figure 8, [24] where all transporters in the W-PCT-E model are disabled except for the Na + / K + -ATPase pump. The permeabilities are modified to values much higher than the default values (increased by a factor of six), and solute types are also matched accordingly. The Na + and Cl − concentrations are increased in a step-wise manner, just as with the values in the primary paper. To preserve electroneutrality, there is an increase of 6.65 mM in both peritubular and luminal sodium and chloride concentrations. In the final step of the current experiment, Latta et al. eliminated the effect of the Na + /K + -ATPase, and here we decreased this effect by a factor of ten. As a consequence, we obtained a step-wise increase in the cellular activities for the primary solutes Na + and Cl − , which is displayed in Fig 5(a). Fig 5(b) shows the original result Table 1 [12]. Electrical potentials are in mV, pressures in mmHg and concentrations in mM.

Variables
Cell Interspace   . Fig 12, [11]). Here, the luminal Cl − entry proceeds exclusively through Cl − /HCO 2 − exchanger, which means the coupled transport coefficient for Cl − /HCO 3 − is considered to be zero. Also, apical and peritubular membrane permeabilities to H 2 CO 2 are set at 10% of original reference in Table 1

PLOS ONE
from the primary paper ( Figure 8 in [24]). It can be seen that there is a good agreement between the W-PCT-E model and the Na + transport model in the mammalian urinary bladder reported in [24,26]. We compared these two models qualitatively as there is not much similarity between model parameters and mathematical formulations, although we defined the same configuration for model structures. But we can see a similar behaviour even though the range of the changes in solute concentration is different. We should mention that in the studies cited above, there are only three solutes (Na + , K + and Cl − ), an impermeant monovalent anion to balance electrical charge in the cell and baths, a nonelectrolyte to produce solution osmolarities, and one active transporter, Na + /K + -ATPase. To see the new set of results, see Fig 5.

iii. Structural analysis
The goal of the present study was also to create and make available a mathematical model of epithelial transport that is sufficiently flexible to accommodate the investigation of different physiological phenomena in the epithelial system. Here, we performed a structural analysis of the W-PCT-E model to both demonstrate this flexibility and to explore the application of this model to a range of physiological perturbations.
To investigate the effect of each transporter in the W-PCT-E model on the overall behaviour, we performed experiments in which we individually inhibited each of the transporters and compared the total epithelial fluxes. We illustrate some of these results in the following section; the first set of simulations addresses the inhibition of both basal and cell-lateral transporters and the second set of simulations addresses the inhibition of the apical cell transporters.
We should highlight that the following experiments aim to investigate the model behaviour corresponding to either biological or, in some cases, extreme assumptions concerning the model configuration. Considering these assumptions, we can test the models' limitations, flexibility, or reliability.
For example, we can see a noticeable increase in cell volume in the case of the inhibition of Na + /K + -ATPase. But, this does not mean the model lacks a self-regulatory system. By inhibiting Na + /K + -ATPase, we disturb the haemostasis in the epithelial model, causing an

PLOS ONE
accumulation of solutes in the cell, causing a noticeable increase in the cell volume. For more information, see Figs 3 and 4 in S1 File.
iii.1 Inhibition of peritubular (IS and IE) transporters. We separately eliminated the Na + /K + -ATPase and two symporters (K + -Cl − and Na + -HCO 3 − ) on both the cell-basal and celllateral membranes and observed the resulting membrane fluxes and cellular concentrations. Inhibition of each transporter was accomplished by setting the coupling transport coefficient to zero.
We present our results in Fig 6. Panel (a) displays the membrane fluxes (ES, IE, IS, ME, MI) and cellular concentrations for the four primary solutes (Na + , K + , Cl − , Glucose) in the case of the original full W-PCT-E model. Panel (b) represents the results when considering the inhibition of the Na + /K + -ATPase, from which one can observe a reduction in all the membrane fluxes and notable changes in the cellular solute concentrations; see bottom row in Fig 6(b). There is a considerable reduction in Na + membrane fluxes, which demonstrates the critical role of the Na + /K + -ATPase in the production of Na + fluxes. Inhibition of the pump stops sodium exit and potassium entry into the cell; thereby, sodium concentration increases (from 19.6 mmol/L to 142.0 mmol/L) while there is a decrease in the cellular potassium concentration (from 137.3 mmol/L to 8.4 mmol/L). Fig 6(c) illustrates the effect of the inhibition of K + -Cl − transporter on the W-PCT-E model activity. One can observe a decrease in total fluxes for both sodium and chloride. Although there is a decrease in the sodium concentration, there is an increase in potassium and chloride concentrations. Fig 6(d) illustrates the response of the W-PCT-E model in the case of elimination of Na + -HCO 3 − transporter. There is a significant decrease in sodium total flux, accompanied by notable growth in glucose fluxes. While there is a decrease in the sodium concentration, one can see an increase in other solute concentrations. To better understand the underlying mechanisms of the W-PCT-E model responses to these structural changes, we narrow our focus to sodium. We then study how the elimination of Na + /K + -ATPase can affect the sodium fluxes across the epithelial membrane. We need to clarify that there are no transporters across the tight junction and interspace basement and the corresponding fluxes across these membranes are either convective or passive. While in the cell-basal, cell-lateral, and apical cell membrane, the passive and convective fluxes are negligible, for the reflection and permeability coefficient values are considered to be very small, see Table 1. Electrochemical fluxes represent the primary source of fluxes across the cell-basal, cell-lateral, and apical cell membranes. Fig 7(a) features all different sodium membrane fluxes. In Fig 7(b), we narrow our focus down to sodium fluxes on the epithelial membrane (summation of the tight junction and apical cell membrane activities) and its components. The epithelial activities subdivide into three components: convective, passive, and electrochemical activities. In Fig 7c, we subdivide the electrodiffusive coupled fluxes into their segments, which are NHE3, SGLT, and Na + -H 2 PO 4 − transporters.
The first row in Fig 7 represents the sodium fluxes for the full W-PCT-E model, considered as the control version. In the second row, we illustrate the sodium fluxes due to the elimination of the Na + /K + -ATPase. By a simple comparison between the first and second rows, one can see that the membrane activities drop considerably after the Na + /K + -ATPase inhibition; as an example, the total epithelial sodium flux decreases from 7.39 nmol.s -1 .cm -2 to 2.19 nmol.s -1 . cm -2 , as seen by comparing the first and second rows in Fig 7(a). To gain insight into the exacerbated decline in epithelial activity, we further divide the epithelial activity into its components, as seen in the second row of Fig 7(b). One can observe a marked drop in the convective fluxes (from 1.348 nmol.s -1 .cm -2 to 0.58 nmol.s -1 .cm -2 ), which is due to the notable drop in the tight junction water fluxes (from 16.7 nmol.s -1 .cm -2 to 7.21 nmol.s -1 .cm -2 ), see Eq (37) and Fig  Fig 6. Changes in the membrane fluxes and cellular concentrations due to the inhibition of transporters on the cell-basal and 7(b) (first and second rows). In contrast, there is an intensification in passive activities (from 0.34 nmol.s -1 .cm -2 to 1.02 nmol.s -1 .cm -2 , which is due to the changes of the normalised electrical potential differences, which appear in the form of linear and exponential expressions, as described by Eq (39), and seen in Fig 7(c) (first and second rows).

PLOS ONE
Here, we aim to highlight the flexibility of the W-PCT-E model, as the user can have a better insight into the system behaviour due to the elimination or cooperation of transporters by simply turning them on or off.
Because of the flexibility of the W-PCT-E model, there is future opportunity for similar analyses to describe the system behaviour due to the elimination of other transporters and their impact on the different solutes.
Here, we applied extreme assumptions to study and understand the performance of the model under these conditions (such as changes in model configurations and model parameters). However, one should not expect to see biological behaviour due to these extreme assumptions. As mentioned in Subsection iii, there is an increase in cell volume due to the elimination of Na + /K + -ATPase; for more detailed information, see Figs 3 and 4 in S1 File.
iii.2 Inhibition of apical membrane (MI) transporters. In this section, we separately eliminate the NHE3 antiporter and apical symporters (SGLT and Na + -H 2 PO 4 − ) and then we study the behaviour of the W-PCT-E model by analysing the results for membrane fluxes and cellular concentrations relative to each scenario. In Fig 8, we present the membrane fluxes in the first row and the cellular concentrations in the second row. By inhibiting NHE3, one can observe a notable decrease not only in sodium membrane fluxes, but also in cellular sodium and chloride concentrations. While there is an increase in glucose membrane fluxes and cellular glucose concentration, as seen in Fig 8(b).
By eliminating SGLT, we consider the total absence of glucose membrane fluxes which is consistent with the model structure as there are no other sources to produce either convective or passive fluxes of glucose. There is a decrease in both Na + and Cl − fluxes; while cellular sodium and chloride concentrations depict decreases, potassium and glucose concentrations demonstrate increases, as can be seen in Fig 8(c).
Removal of Na + -H 2 PO 4 does not show any significant changes neither in epithelial fluxes nor cellular concentrations in any of the four primary solutes, as clearly shown in Fig 8(d).
Here, we focus on the sodium fluxes, as NHE3 is the primary source of sodium fluxes in the epithelial model. NHE3 inhibition stops the exit of hydrogen and ammonia from the cell and the entry of sodium into the cell. Sodium concentration drops from 19.6 mmol/L to 10.5 mmol/L, which demonstrates the role of NHE3 in the production of Na + fluxes.
To have a better understanding of Fig 8(b) and a deeper insight into the underlying mechanisms of NHE3 and its effect on the sodium fluxes, we study all sources of sodium fluxes, and the results are featured in the third row in Fig 7. After inhibition of NHE3 in Fig 7 (third row), sodium activities in epithelial membranes (ME and MI) drop considerably (the total activity drops from 7.39 nmol.s -1 .cm -2 to 3.37 nmol. s -1 .cm -2 , Fig 7(a)). To better understand the origin of these changes, we further divide the epithelial sodium fluxes into convective, passive, and electrochemical component fluxes, as seen in Fig 7(b) (third row).
As we mentioned before, convective and passive fluxes in the epithelial membrane are mainly through the tight junction. In Fig 7(b) (third row), one can observe a notable drop in the convective fluxes (from 1.348 nmol.s -1 .cm -2 to 0.88 nmol.s -1 .cm -2 ), which is mostly due to the reduction in the water fluxes (from 16.7 nmol.s -1 .cm -2 to 11.01 nmol.s -1 .cm -2 ), see Eq (37).
There is a slight reduction in passive fluxes accompanied by changes in direction (from 0.34 nmol.s -1 .cm -2 to −0.15 nmol.s -1 .cm -2 ). To explain the decrease in the passive fluxes, we

PLOS ONE
need to bring to attention that the main driving force for the passive fluxes is not only the normalised electrical potential differences (which are in the form of linear and exponential components, see Eq (39)) but also solute concentrations. The total electrochemical flux across the epithelial membrane falls from 5.67 nmol.s -1 .cm -2 to 2.64 nmol.s -1 .cm -2 , see Fig 7(b). We then investigate the coupled transport fluxes in Fig 7(c) by visualising the components individually: NHE3, SGLT, and Na + -H 2 PO 4 − . Sodium fluxes for NHE3 dropped to zero as a result of inhibition. Due to changes in sodium cellular concentration (from 19.6 mmol/L to 10.7 mmol/L), the sodium concentration gradient between the lumen (140 mmol/L) and cell increases notably, increasing the electrochemical potential difference of sodium across the apical cell membrane, which increases the activity of transporters related to the production of sodium fluxes (SGLT from 2.05 nmol.s -1 .cm -2 to 2.50 nmol.s -1 .cm -2 , and Na + -H 2 PO 4 − from 0.11 nmol.s -1 . cm -2 to 0.143 nmol.s -1 .cm -2 ).

VII. Discussion
Mathematical modelling provides a tool for investigating complex physiological phenomena that is hardly accessible to experimentation. A comprehensive example of these tools is the work done by Alan Weinstein in the field of renal modelling. Weinstein has provided a valuable resource for mechanistic modelling of many aspects of renal function through many publications [9-16, for example]. However, existing mathematical models of epithelial transporters are, in general, not readily findable, accessible, interoperable, or reusable (FAIR) [27] for researchers who are new to this field. Here, we want to provide a well-tested reusable implementation of Weinstein epithelial PCT models, collated into a single consistent implementation that we believe encourages reuse.
We have presented here what we believe to be a comprehensive and FAIR epithelial model for the PCT of the renal nephron. This model encapsulates and recapitulates the seminal work performed by Weinstein and colleagues [9][10][11][12][13] over many years. In the case of model reproducibility, we have demonstrated that the W-PCT-E model reported here can reproduce many different aspects from related or earlier works. We chose three exemplars [10,11,24], the constant parameters and boundary conditions tuned according to the model of interest. In all cases, we observed a close agreement between the W-PCT-E simulation outcomes and the results reported in previous works, see Section ii (Table 1, Figs 4 and 5).
To show the flexibility of the implementation of the W-PCT-E through the application of structural analysis, we investigated the impact of each transporter on the W-PCT-E model responses (total fluxes and cellular concentrations) through the inhibition of each transporter, see Section iii (Figs 6-8).
To further study the limitations and performance of our model, we examined it under some extreme conditions. Due to the non-physiological nature of these conditions, one should not expect to observe real-world physiological behaviours. Here, we aim to demonstrate the mathematical feasibility of the W-PCT-E model configuration by removing or appending different elements. As an example, there is an increase in the cell volume due to the elimination of Na + /K + -ATPase (Fig 6).
To demonstrate the comprehensiveness and flexibility of the W-PCT-E model, we now briefly explore various physiological phenomena using our model.

Clinical studies have shown that excess glucose in the cell and bloodstream is associated with Type II Diabetes (T2D) [28-36]
The inhibition of SGLT for the treatment of T2D [37,38] has shown improvement of glycemic control and T2D [39][40][41]. According to the W-PCT-E model, the inhibition of the SGLT transporters decreases the cellular Na + and glucose fluxes which decrease the cellular concentration of glucose, see Fig 8(c). The model also predicts a decrease in interspace concentration of glucose, which would lead to less glucose available for reabsorbtion into the blood. Thus, demonstrating a clear consistency between the W-PCT-E and the clinical findings regarding the inhibition of SGLT and improvement of T2D.
Experimental reports illustrate that NHE3 residing in the apical membrane mediates transcellular reabsorption of Na + and fluid reabsorption [42,43], and show that NHE3 deficiency can cause a reduction in Na + reabsorption [8,44] These observations put renal Na + reabsorption via NHE3 in a central position in the development and control of salt loading-and volume expansion-mediated hypertension [45]. We performed some simulation experiments to check that our W-PCT-E model is consistent with these findings by inhibiting NHE3. In doing so, we observed that inhibition of the apical NHE3 decreases the intracellular concentration of HCO 3 − and HCO 2 − , which in turn is transported via the apical Cl − /HCO 3 − or Cl − /HCO 2 − exchangers. The inhibited apical NHE3 and the Cl − /base exchangers work in parallel and produce net Na + reduction and Cl − reabsorption in the proximal tubule, see Fig 1(b)-1(e) in S1 File (cellular concentrations) and Fig 1(F) in S1 File (cell volume). These findings are also consistent with other studies demonstrating that the inhibition of either the apical NHE3 and Cl − /base exchangers inhibits net Na + and Cl − absorption in the proximal tubule [11,[46][47][48][49].

Experimental results reveal that a global knockout of NHE3 gene on the proximal tubule reduces water fluxes and HCO 3 − absorption [50]
The W-PCT-E model demonstrates that the apical NHE3 exchanger can inhibit HCO 3 − flexes in the proximal epithelial model. We designed a set of experiments, in which we decreased the coupled transport coefficient for NHE3, see Fig 1 in S1 File. We observed that the inhibition of NHE3 lowers the Cl − /HCO 3 − and reduces the cell volume. These findings are consistent with the earlier findings [11,[50][51][52].
Experimental reports show that there is a coordination between basolateral Na + /K + -ATPase and apical NHE3 activities; they simultaneously regulate the Na + transport which can be the cause of inhibition or activation of transepithelial Na + transport [31,45,[53][54][55] In the W-PCT-E model, one can see that there is strong coordination between basolateral Na + /K + -ATPase and apical NHE3 activities; inhibition of NHE3 (Na + /K + -ATPase) can cause the inhibition of Na + /K + -ATPase (NHE3). We can conclude from this observation that in the W-PCT-E model, sodium is primarily reabsorbed via NHE3 which is regulated by Na + /K + -ATPase, see Fig 2 in S1 File.

The kidney plays a critical role in the regulation of body electrolyte and fluid balance, primarily occurring in the proximal tubule segment of the nephron [56-58]
A balance of body extracellular electrolyte composition and fluid volume is essential for all animals and humans to survive. Either excess or shortage of crucial extracellular electrolytes or overall fluid volume may lead to disturbance of blood pressure and abnormalities in cellular functions, including cell volume [58][59][60][61][62]. The W-PCT-E model can capture the impact of modulation of different transporters on the fluid balance and cellular volume. For example, with the inhibition of basal-cell Na + /K + -ATPase or the activation of apical NHE3 or SGLT, we observed increase in the cell volume. Such observations are consistent with previous modelling approaches and experimental reports [31,60,62,63]. In subsection iii, we mentioned that there is an increase in cell volume due to the elimination of Na + /K + -ATPase; we also added Figs 3 and 4 in S1 File. The implementation of the present epithelial transport model to support the composition and parameterisation of the proximal tubule epithelial system (the W-PCT-E model) is available on Github under a license which allows open and unrestricted reuse. The implementation is in Python (version 3) to make it broadly accessible. As demonstrated above, the implementation is sufficiently flexible and configurable to support the generation of different epithelial models. Users need only to extract the required constituent transporter modules from the available repository and then integrate them into the desired configurations. The boundary conditions need to be changed according to the epithelium of interest. In sharing the model implementation between the authors of this manuscript, we tested for reusability as we each work with W-PCT-E model.

VII. Conclusions
We believe that the time is now right to develop a reproducible and FAIR virtual nephron. Here, we term this the iNephron. It is critical to ensure that the capabilities of published models are captured while being sufficiently flexible and configurable to support the generation of novel models to investigate specific scientific or clinical questions. To achieve this, iNephron will provide a tool to investigate the fundamental mechanisms involved in hypertension, diabetes, and many other kidney diseases. Furthermore, iNephron will be developed and shared following established best practices in open and reproducible science to guarantee that the scientific community can benefit and extend from this work to improve our collective understanding of these diseases.
The work we presented here is our first step toward achieving the iNephron. As we look to grow the repository of available transporter modules, we also look to better follow FAIR principles [27,64], and move toward a standards-based repository of reusable modules (e.g., [65]). To ensure the composed models are meaningful and that model composition can occur reliably, we are planning to migrate our repository of reusable modules to ensure thermodynamic consistency [3][4][5][6][7]. This is a requirement for the arbitrary composition of models needed for iNephron.
As we have shown, the W-PCT-E model is actually a generic epithelial model which is flexible and configurable to support the generation of different epithelial models, where the user is able to meet their design requirements, easing the process for "getting started" with a novel modelling study. All one needs is to provide the model with a set of transporters and boundary conditions appropriate for the epithelial model of interest. Additionally, by establishing a comprehensive ability to perform sensitivity analyses, we provided tools by which future users are able to test their own additions or modifications of this model with confidence. Our testing, summarised in this manuscript with more detailed findings in the supplemental material, shows that our model behaves as expected in physiological terms.