Computational modeling of peritoneal dialysis: An overview

1.   Peritoneal dialysis

Over 850 million people worldwide, i.e., ~1 in every 10, suffer to some degree from chronic kidney disease ^[1], among which over 3.8 million are on dialysis, either peritoneal dialysis (PD) or hemodialysis (HD) ^[2]. During PD, dialysis fluid is instilled into the abdominal cavity via a permanent catheter. The lining of the abdominal cavity (the peritoneum) acts as a semi-permeable membrane for solute and water transport. PD removes waste products from blood plasma by diffusion and convection across the peritoneal membrane into the dialysis fluid in the abdominal cavity (see Figure 1). Excess water is removed via osmosis using glucose or -less commonly- amino acids as crystalloid osmotic agents in the dialysis fluid and icodextrin, a mixture of glucose polymers, as colloid osmotic agent ^[3,4]. While glucose and amino acids are rapidly absorbed across the peritoneal membranes, icodextrin is absorbed only to a limited extent due to its large molecular weight (on average 12–20 kDa ^[5]) exerting a sustained osmotic effect suitable for sustained fluid removal during longer dwell times ^[3,6]. Three types of pores are distinguished to properly describe solute and fluid transport across the peritoneal membrane: Small pores, which facilitate small solute transport; large pores, which are responsible for movement of larger solutes such as proteins; and ultra-small pores, which enable the passage of water molecules ^[7]. Selective water transport through the latter is responsible for an initial decrease in dialysate sodium concentration during the early dwell phase and it may explain glucose's osmotic efficacy despite its high diffusibility. Continuous ambulatory PD (CAPD) and automated PD (APD) are the two major modalities of PD used in routine practice. With CAPD, the dialysis fluid (1–2.5 L) is exchanged (drained and instilled) 3–4 times daily manually through the catheter ^[8]. During daytime, glucose-based solutions are typically used with dwells lasting 4–6 hours. During the night, a long dwell with an icodextrin solution is commonly used. With APD, a cycler with patient specific treatment parameters is connected to the abdominal catheter throughout the session (preferably at night) to automatically perform exchanges with (glucose-based) dialysis fluid over a period of time (usually 4–5 per night) ^[9]. This is often combined with a long dwell with an icodextrin solution (14–15 h) during the day. Both CAPD and APD are usually conducted every day of the week.

PD has several advantages compared to HD: it enables continuous gradual removal of waste (instead of intermittent HD, which is characterized by a 'saw-tooth pattern'), does not require blood access, provides more patient autonomy as the treatment is performed at home, and is less expensive. Residual kidney function is also better preserved ^[10] compared to HD. However, PD has important shortcomings. Technique survival is limited (median 3.7 years ^[11]) due to recurrent peritonitis (inflammation of the peritoneal membrane), catheter dysfunction or membrane failure (due to exposure to high (harmful) dialysate glucose concentrations required for osmotic water removal), and small solute clearance is relatively low. Due to technique failure or low small solute clearance (with disappearance of residual diuresis) patients often have to switch to HD after several years ^[11,12].

Despite 60 years of progress in PD, it is used only by ~11% of dialysis patients with considerable variations across countries, mostly due to non-medical reasons ^[13,14,15]. Strikingly, Hong Kong has ~80% end stage kidney disease (ESKD) patients on PD to reduce expenditure ^[16] after an implementation of "PD first" policies. Mexico, Thailand, and Singapore also have high rates of PD due to availability of medical personnel to assist with PD. To increase PD usage, improve outcomes, and better personalize PD regimens, there are several novel techniques under development such as in home production of peritoneal dialysate ^[13], novel PD solutions with improved biocompatibility and/or glucose-sparing, continuous flow PD ^[17,18], and sorbent assisted PD ^[19,20,21] (Figure 1). Glucose-sparing strategies include the use of alternative osmotic ingredients such as (commercially available) icodextrin and amino acids or alternatives that are under development such as taurine, polyglycerol, carnitine, xylitol, and alanyl-glutamine ^[3,22]. Moreover, combinations of crystalloid and colloid osmosis (so-called bi-modal osmosis) have been explored by mixing glucose with icodextrin, This approach leverages rapid initial ultrafiltration provided by glucose and sustained ultrafiltration by icodextrin and was shown to enhance fluid removal while reducing overall glucose exposure ^[4,23]. A study by Stachowska-Pietka et al. showed that comparing glucose 2.27% and icodextrin 7.5% in a patient group with average transfer rate, there was a 204% increase in ultrafiltration, 3.6% increase in glucose and carbohydrates absorbed, 163% increase in sodium removal, 13% increase in urea removal, and 17% increase in creatinine removal ^[6]. However, potential consequences are accumulation of icodextrin metabolites. In addition, using icodextrin and bimodal solutions increases the sodium removal leading to hyponatremia ^[24]. Since icodextrin modeling is relatively new, it has been briefly included in Sections 4 and 5 but has not been discussed extensively in this review.

In continuous flow PD (CFPD), there is a continuous flow of dialysate through an inflow catheter and outflow catheter ^[17,18,25]. The continuous flow of fresh dialysate through the abdominal cavity maintains a large plasma–dialysate concentration gradient, increasing solute transfer across the peritoneal membrane. Dialysate glucose concentration can be kept nearly constant, thereby maintaining a constant osmotic gradient and ultrafiltration rate at lower (peak) glucose concentrations than in conventional PD ^[26], possibly slowing functional deterioration of the peritoneal membrane and reducing the rate of peritonitis (both associated with high glucose concentrations). Sorbent-assisted PD regenerates the dialysis fluid using a sorbent cartridge ^[19,27,28].

The development of new technologies often takes a considerable amount of time owing to the learning curve, ensuring transparency, ethical requirements, trial regulations et cetera ^[29]. The new (post-2020) European Union Medical Device Regulation (EU-MDR) rules related to manufacturing of medical devices are more stringent, leading to delays in certification. As such, there is a need to supplement the pre-clinical research and clinical trials with other methodologies to accelerate the design and manufacture and market novel PD technologies. Computational models are a revolutionary tool in the field of health, medicine, and life sciences due to the ease of optimization, non-intrusiveness, and most importantly, the interpretation of complex interdependent physical processes. They are being increasingly used to study the influence of biological processes in medical devices, including, for example, thrombogenic reaction to biomaterials ^[30], controlled drug release ^[31], and effect of implant surface roughness on protein adsorption ^[32]. Computational modeling of PD may help clinicians make informed decisions, such as those involving catheter placement and treatment optimization ^[33,34]. Moreover, computational modeling is often an inexpensive and effective way to simulate complex natural phenomena.

In this review, we focus on computational models for PD, especially ordinary differential equation (ODE) models. ODE models represent homogeneous well-mixed compartmental systems that describe the evolution of the (solute) concentrations over time in the various compartments. In contrast, PDE models capture systems that are not well-mixed (accounting for catheter obstructions, recirculation, etc.), enabling the solute concentrations to be defined as a function of both time and space. As PDE models are beyond the scope of this discussion, they will not be addressed here. We follow the history of PD modeling to highlight the advances and remaining lacunas as well as learn which level of simplicity and physiological detail is necessary to reach a research aim. This review, containing general concepts of PD models, is specifically intended for experimental scientists or clinicians with little or no computational modeling experience, but who are interested in introducing modeling into their research practice. By providing an overview of computational models, the individual aspects of the models, and the scope of the model, we aim to simplify, enhance, and accelerate the integration of modeling into clinical practice to promote better understanding of device-patient interaction. In Section 6, we show an example of a scenario to illustrate how one can set up a mathematical model of PD using the models highlighted in this review.

2.   Compartmental models for peritoneal dialysis

There are many types of computational models, including linear and non-linear, deterministic and stochastic, discrete and continuous and spatial and non-spatial ones. For general reviews on computational modeling we refer to Yates et al. ^[35], King et al. ^[36], and Brown et al. ^[37]. The most common design for models of PD is the compartmental model. They are simple in nature for physiological, kinetic and dynamic modeling. In compartmental models, the body is divided into theoretical compartments such as the peritoneal membrane, peritoneal cavity, and total body water (examples in Figure 2). The general assumption is that a tissue or organ can be represented as a homogeneous compartment, governed by conservation of mass and other properties of interest (e.g., charge). A compartmental model consists of volumes connected by fluxes of some entity. Each compartment can be a volume representation (in case of PD) and it characterizes the essential physics and chemistry of the biological environment. The flow rates and interactions between the compartments are described by the parameters of the compartmental model. Some PD models use a single compartment model with just the peritoneal cavity and body considered to be a constant source of solutes ^[38] while some models include multiple compartments (e.g., the distributed model, ^[39]), with the peritoneal tissue, peritoneal cavity, and interstitium as different compartments. Compartment models are lumped models but despite their simplicity, they usually capture the underlying physical and biological phenomena quite well, which is why they are commonly used for PD modeling.

The models usually employ fluid and solute mass balance equations, which are described in Sections 3 and 4. This depends on the type of fluid flow (static or continuous), complexity (single compartment or distributed), and membrane model (continuous semi-permeable membrane or 3 pore) to decide which and how the fluid and solute flows should be properly modeled. In case of PD, physiologically, the barrier is a semi-permeable peritoneal membrane. The compartments can exchange the fluid and solutes without any barrier or with a barrier.

The simplest interpretation of the peritoneal membrane is a simple membrane that enables passive diffusion of solutes ^[40,41,42]. The diffusive flow of molecules (J_s in mol/s) will be dependent on the concentration gradient $\frac{\mathrm{d}c}{\mathrm{d}x}$ across the barrier:

where $D$ is the diffusion coefficient (${\mathrm{m}}^{2}{\mathrm{s}}^{-1}$), A is the effective surface area available for diffusion (m²), c is the concentration $\left(\mathrm{m}\mathrm{o}\mathrm{l}/\mathrm{l}\right),$ and $x$ is the distance across the membrane (m). The negative sign is because the solute flux is in the opposite direction of the concentration gradient. Integrating the above equation (Eqs (1) and (2)) from x = 0 to x = Δx (across the membrane),

where $PS = DA/\mathrm{\Delta }x$, and ${c}_{1}, {c}_{2}$ are the concentrations on either side of the membrane, at x = 0 and x = Δx. ${c}_{1}$ and ${c}_{2}$ are usually chosen so that the sign of ${J}_{s}$ is positive when transport is directed into the peritoneal cavity. Eq (3) can be used to calculate the diffusive plasma clearance of a particular solute from c₁ to c₂ (J_s/c₁) or vice versa (J_s/c₂). The factor PS is known as the diffusion capacity of the solute species and is defined as the maximal absolute diffusive plasma clearance (when the concentration is zero on one side of the membrane). Many abbreviations and terms are in use for this parameter such as mass transfer area coefficient (MTAC), k_BD etc.

2.1.   Three-pore model of peritoneal dialysis

A perplexing thing about the peritoneal membrane is that it is a semi-permeable membrane that enables passage of albumin and other large proteins to a limited extent but also restricts the bulk movement of electrolytes ^{[19,43,44,45,46]}. In 1981, Nolph et al. stated, "It is a system that displays characteristics of some very large-pore radii when assessed by diffusion studies, and some very small-pore radii when assessed by ultrafiltration and solute reflection coefficients" ^[47]. The reflection coefficient, $\sigma,$ describes the convective hindrance of a molecule (a value of 100% means convective clearance does not occur). Nolph et al. were the first to hypothesize that the peritoneal membrane is a heteroporous membrane with both small and large pores. Small pores with low solute permeability would be located in the proximal capillary and facilitate ultrafiltration due to the relatively high hydraulic pressure and high osmotic gradient in the proximal capillary. Toward the venous capillaries, the hydraulic pressure drops, and the oncotic pressure increases as a result of increased protein concentration caused by capillary ultrafiltration. A predominance of large "pores" in the venular capillaries with high permeability would facilitate diffusive solute exchange. Glucose would be more readily absorbed, resulting in lower osmotic pressure, and the ultrafiltration rate would be reduced at the venular capillaries also because of the relatively low hydrostatic and high oncotic pressure. In summary, most ultrafiltration would occur through the "small" pores with low hydraulic permeability and most solute exchange through highly permeable "large" pores (Figure 3A)). Aside from these two pathways, there is also a separate channel for movement of water that is inaccessible to solutes. To account for all three pathways, the three-pore model was developed.

In 1991 Rippe et al. ^[7] proposed the three pore model (Figure 3B)). This is the most common representation of the peritoneal membrane. It divides the body into two compartments – the distribution volume (specific to each solute) and the peritoneal cavity. The solute and volume flows through the different pores that are defined by the Starling and Patlak equations, respectively (explained in Eqs (4) and (29)). The most abundant (99% of the pore fraction) is the protein restrictive water-soluble pathway (15 to 36 $Å$), which is responsible for 90% of the hydraulic conductance. The "large pores" of 250 $Å$ constitute 0.01% of the pore population and represent 8% of the hydraulic conductance, and the final pore fraction is the ultrasmall pores (0.99%) and is responsible for only 2% of the hydraulic conductance. The "small pores" correspond to the gaps between the endothelial cells and the "large pores" correspond to the venular interendothelial pathways. These ultrasmall pores (2.3 to 15 $Å$) are permeable only to water (reflection coefficient is unity for all solutes) and were later shown to represent the aquaporin (AQP-1) water channels. AQP-1 was first identified in animal biopsies one year after Rippe proposed its existence ^[48]. Yang et al. showed that knocking out the aquaporin channels in mice reduced the water transport by 50%, thus proving the importance of the transcellular pores in the three-pore model ^[49]. The Rippe model can also predict the transport of intermediate solutes as well as large solutes with reasonable accuracy.

Venturoli et al. proposed a series of two porous membranes (Figure 3C)) to model the bi-directional clearances of macromolecules ^[50]. The first layer is the three-pore membrane by Rippe, identified as the capillary endothelium. This is followed by a second three pore membrane consisting mostly of large pores (95%), transcellular pores (2%), and small pores (3%), which is a lumped representation of extracellular interstitium. This model is able to mathematically explain the build-up of tracer albumin ^[51] seen in rats, which the one layer three pore model cannot.

Another approach to modeling the peritoneal membrane is the distributed model by Flessner et al. ^[39] (Figure 3D)). They model the peritoneal tissue space as a tube-and-shell exchanger with a constant void space (blood capillary or the plasma space) and the interstitium. Water movement occurs throughout the tissue while solute convection occurs only across the plasma capillaries. All plasma capillaries are assumed to have the same shape and size, which makes this a single pore model. This model is also complementary to the Venturoli model ^[50] in that it enables calculation of accumulation and release of certain substances dissolved in dialysis fluid. Because of the additional modeling of the surrounding tissue space, the Flessner model can give insight into the mechanisms ongoing in the peritoneal tissue which could, for example, be useful for modeling transport of drugs being administered intra-peritoneally.

3.   Volume flow calculations

The osmotic gradient created by a hypertonic dialysis fluid drives the flow of water from the plasma by osmosis, commonly referred to as ultrafiltration ^[52]. In APD, underestimation or overestimation of the ultrafiltration volume may result in setting a wrong drain volume of the peritoneal cavity. The use of wrong drain settings may cause discomfort to the patient and reduce the efficiency of the dialysis session. As such, it is important to properly calculate all the flows occurring during PD, i.e., the amount of net ultrafiltrate as determined by

- the glucose concentration used in the dialysis fluid and subsequent water flow due to the osmotic gradient

- the flow of dialysate into the peritoneal cavity at a certain flow rate (APD, CFPD) or by gravity (static dwell, CAPD)

- the lymphatic flow of fluid from the intraperitoneal space towards the lymphatic space

- other physical (and patient-specific) attributes such as the condition of the peritoneal membrane, intraperitoneal volume and pressure and the lymphatic system of the patient etc.

There have been limited efforts to model volume change in the intraperitoneal cavity. We describe the different ways of estimating peritoneal volume change and the ways to obtain the parameters for the flux equations in further subsections.

To describe the water flow (volume flow), the Starling equation is often used, which is essentially a modified version of Ohms law, $\mathrm{F}\mathrm{l}\mathrm{o}\mathrm{w} = \mathrm{H}\mathrm{y}\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{u}\mathrm{l}\mathrm{i}\mathrm{c}\mathrm{ }\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{d}\mathrm{u}\mathrm{c}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}\cdot \mathrm{\Delta }\mathrm{P}\mathrm{r}\mathrm{e}\mathrm{s}\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{e}$. The original Starling equation was written to describe the fluid flow from a capillary to the interstitial space ^[53]. This equation was later corroborated by Kedem and Katchalsky from a thermodynamics point of view ^[54]. As fluid builds up in the peritoneal cavity, it exerts an outward interstitial force on the peritoneal membrane (${P}_{\mathrm{I}},$ Figure 4). Due to the fluid moving into the interstitial space, the corresponding fluid pressure rises, and this opposing pressure is the capillary hydrostatic pressure (${P}_{\mathrm{V}}$, Figure 4). The plasma proteins exert a colloid osmotic pressure or oncotic pressure (${\pi }_{\mathrm{V}},$ Figure 4) to draw water back from the peritoneal cavity to the plasma (opposite to glucose). The osmotic agent in the dialysate exerts an outward osmotic pressure (${\pi }_{\mathrm{I}},$ Figure 4). To model the above, Eq (4) consists of two inward forces, i.e., interstitial hydrostatic (${P}_{\mathrm{I}})$ and plasma colloid osmotic pressure (${\pi }_{\mathrm{V}}$), and two outward forces, i.e., hydrostatic pressure (${P}_{\mathrm{V}})$ and interstitial osmotic pressure (${\pi }_{\mathrm{I}}$) (see Figure 4). Here, we will also take into account that the peritoneal membrane is semi-permeable and will have different reflection coefficients for solutes of different sizes. Note that by inward or outward, we mean into or out of the peritoneal cavity. We refer the reader to the appendix for an overview of all notations.

The volume flow ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$) across the peritoneal membrane is thus modeled as follows:

The volume in the peritoneal cavity is $V,$ the hydraulic conductance of the peritoneal membrane is given by ${L}_{p}S$ ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}\mathrm{m}\mathrm{m}\mathrm{H}{\mathrm{g}}^{-1}$), and the reflection coefficient of the solute is given by $\sigma$ (dimensionless). I denotes interstitial parameters and V denotes the plasma parameters. We have to make certain modifications to Eq (4) as the fluid is not in full contact with the peritoneal membrane, which means that the whole membrane surface cannot be taken into account for the calculation of the volume flow but rather the peritoneal surface area in contact with the dialysis fluid ($af$, fraction of the peritoneal surface area in use, dimensionless). This gives us the modified Starling equation (Eq (5)).

In addition to the fractional peritoneal surface area in use, we can divide the volume flow among the various types of pores present in the peritoneal membrane. If we use the classical 3-pore model here, we can divide the net volume flow into flows across the ultrasmall, small, and large pores by multiplying their contribution to the ultrafiltration (α, dimensionless).

where ${J}_{\mathrm{v}\mathrm{C}}, {J}_{\mathrm{v}\mathrm{S}}$ and ${J}_{\mathrm{v}\mathrm{L}}$ are the volume flows across the ultrasmall, small, and large pores, respectively.

The change in volume is then calculated as a total of all the volume flows minus the net lymphatic flow ^[38],

where $L$ is the lymphatic flow rate ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$) (a sum of all lymphatic flows that drain the peritoneal cavity, i.e., interstitial, diaphragmatic, pelvic and omental ^[55]), considered to be around 0.15–0.5 ml min^-1 (0.216–0.72 L/day).

Öberg et al. extended Eq (9) for APD ^[56] and CFPD ^[38],

where ${J}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}$ and ${J}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{n}}$ are the fill and drain flow rates during CFPD ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$).

The initial fill volume and drain volume are in general known. The residual volume, ${V}_{\mathrm{r}},$ may be estimated from the dilution of a solute concentration (e.g., albumin, creatinine, or total protein) measured in the drained effluent, by measuring concentration just after instillation of a known fill volume V_fill ^[57] (Eq (11)):

where ${c}_{0}$ is the measured concentration of the solute just after filling the peritoneal cavity (rapid mixing is assumed), and ${c}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{n}}$ is the concentration of the solute in the drain bag collected in the previous session. Knowing these volumes, the net ultrafiltration volume, UF, can be calculated as,

where $t1$ is the time of instillment of the dialysis fluid and $t2$ is the time of draining.

Lymphatic flow may or may not be considered depending on the patient and the model. Volume flows are also estimated assuming that the rate of ultrafiltration is known or can be calculated ^[58].

3.1.   How to estimate volume flow parameters?

Few parameters are derived from experimental observations and some are obtained from optimization of models.

3.1.1.   Residual Volume, ${V}_{r}$

Residual volume, ${V}_{\mathrm{r}}$, is a critical factor for accurate assessments of peritoneal membrane parameters such as osmotic conductance to glucose (Section 3.1.4). Clinical conditions such as constipation, catheter position or dysfunction, intraperitoneal adhesions or changes in peritoneal cavity compliance can affect residual volume. The residual volume may be estimated from the dilution of a solute (e.g., albumin, creatinine, urea, total protein, dextran 70 or inulin ^[59]) after a fresh fill calculated with the concentrations measured in the drained effluent and in the fresh fill just after the instillation of a known fill volume V_fill ^[57] (Eq (13)):

where ${c}_{0}$ is the measured concentration of the solute just after filling the peritoneal cavity (rapid mixing is assumed) and ${c}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{n}}$ is the concentration of the solute in the drained effluent. Urea and creatinine-based calculations lead to an overestimation of ${V}_{r}$ and can be handled with a correction factor ^[60]. Albumin-based volumes conformed strongly with three-pore model estimates ^[60].

3.1.2.   Difference in hydrostatic pressure, $\mathrm{\Delta }P$

The hydrostatic pressure difference between the interstitium and the capillaries in most organs is usually around 10–17 mmHg (13.6–23.1 cm H₂O) ^[61]. Durand et al. measured the mean intraperitoneal hydrostatic pressure (IPP) to be 9.56 mmHg (13 cm H₂O) ^[62]. However, in the case of a static PD dwell, the pressure builds up inside the peritoneal cavity due to ultrafiltration so that $IPP$ is a function of time for a patient in a sitting position. Twardowski et al. derived the empirical formula for dependence of IPP as a function of intra-peritoneal volume as follows ^[63],

where $\Delta {P}_{0}$ is the baseline hydrostatic pressure (without dialysate), ${V}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}, {V}_{\mathrm{r}},$ and ${V}_{t}$ are the initial fill volume, residual volume, and the intraperitoneal volume at time $t$, respectively. ${V}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}$ is known (1.0–2.5 l depending on the patient), and the residual volume ${V}_{\mathrm{r}}$ is calculated from eq 13. The peritoneal volume at different times, ${V}_{\mathrm{t}}$, can, for example, be calculated from the dilution of radioactive ¹²⁵I serum albumin (RISA) ^[64,65] during the dwell or from a direct volume recovery technique ^[66]. 490 comes from the slope of the fitting equation of intra-abdominal pressure versus intra-peritoneal volume, as measured by Twardowski et al. ^[18].

For continuous flow PD, the hydrostatic pressure changes depend on drain and fill flow rates, and will generally decrease (drain) or increase (fill) in line with intra-peritoneal volume.

3.1.3.   Reflection coefficients, $\sigma$

Sieving or reflection coefficients for solutes can be determined from experiments ^[67,68] or from analytical equations such as that of Drake et al. ^[69]

where $\lambda$ (dimensionless) is the ratio of the solute radius to the membrane pore radius. The reflection coefficients are calculated separately for different pore sizes and solutes.

Another analytic solution widely used for the three pore model is ^[70],

Theoretically, one can also obtain the (steric) reflection coefficient from the equation by Anderson et al. if one knows the equilibrium concentration of the solute on both sides of the peritoneum, i.e., the blood plasma and the peritoneal cavity ^[71]

where $\varPhi$ is the equilibrium partition coefficient (dimensionless), given by $\varPhi = {\left(\frac{{C}_{1}}{{C}_{2}}\right)}_{\mathrm{e}\mathrm{q}\mathrm{l}\mathrm{b}}$ and ${c}_{1}$ and ${c}_{2}$ are equilibrium concentrations on both sides of the peritoneal membrane.

3.1.4.   Osmotic conductance to glucose

The osmotic conductance to glucose (OCG) of the peritoneal membrane is usually calculated from the double mini-peritoneal equilibration test (dm-PET) described by La Milia et al. ^[72]. OCG is used to measure the ability of the peritoneal membrane to transport water in response to a crystalloid osmotic gradient. It may hold particular value for the long-term follow-up of PD patients to evaluate the integrity of the peritoneal membrane over time and helps to identify patients at risk for encapsulating peritoneal sclerosis ^[73]. In the clinic, a 60 min dwell 1.5% (~83 mmol/l) glucose is followed by 60 min 4.25% glucose (~236 mmol/l), and the difference between the drained volume is used to calculate the osmotic conductance to glucose (OCG or ${\sigma }_{\mathrm{g}}{L}_{p}S, ml\;{min}^{-1}\;mmH{g}^{-1}$). Assuming that the initial volume flow is due to osmosis only, from Eq (5), we get

where $R$ is the gas constant (82 $\mathrm{m}\mathrm{l}\cdot \mathrm{a}\mathrm{t}\mathrm{m}{\mathrm{K}}^{-1}\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$) and $T$ is the absolute temperature ($\mathrm{K}$). Typical values for OCG are between 0.003–0.004 $ml\;{min}^{-1}\;mmH{g}^{-1}$. Knowing the reflection coefficient of glucose from Eqs (15) or (16) or from experiments, ${L}_{p}S$ can be calculated to be around 0.08 ml/min per mmHg ^[7,74,75,76]. OCG may also be derived from a single dwell, as described by Martus et al. ^[61]. However, it must be noted that a high residual volume limits the reliability of OCG estimation ^[73].

3.1.5.   Fractional peritoneal surface area in contact with the dialysate, $af$

Keshaviah et al. calculated the empirical relationship between fill volume and peritoneal surface area in contact ^[77,78],

where V is the intra-peritoneal volume and 0.00077 is an empirical constant in ml^-1.

The fractional peritoneal surface area ($af,$ dimensionless) in contact with dialysate can also be obtained from stereological experiments by superimposing a grid over the CT scan of the peritoneum ^[79,80] or by magnetic resonance imaging ^[81]. The fractional surface area in contact with dialysate is an important determinant of the mass transfer coefficient of solutes. It can be improved by increasing the fill volume ^[82], agitating the dialysate, or adding surfactants ^[83] but the latter has not been tried in humans. Eq (20) has been shown to be identical to the cube-square law for intra-peritoneal volumes up to ~2300 ml ^[38].

3.1.6.   Contribution to hydraulic conductance, $\alpha$

Fractional hydraulic conductances (α) may be estimated by fitting the experimental results to the model ^[50,78]. The usual values used in the three-pore model are given in Figure 3.

3.1.7.   Ultrafiltration rate, $UFR$

There have been multiple ways of determining ultrafiltration rate $\left(\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}\right)$ throughout the literature. The average ultrafiltration rate may be estimated from the drain, ${V}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{n}}$, and fill volume, ${V}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}$ during a static dwell ^[58],

The UF rate can also be theoretically calculated ^[84,85],

Lymphatic flow is often considered to be constant throughout the dwell.

The UF rate is sometimes calculated as a function of time from the interpolated intraperitoneal volumes ^{[86,87,88,89,90,91]}, as the UFR usually decreases during a static dwell ^[84,92,93].

A simple empirical exponential model was used by Randerson et al. to capture the time dependent UFR ^[94],

where $A$ is a fitting constant and $\beta$ is the time constant (a value of 0.0192 min^-1 was used).

Total ultrafiltrate volume can also be calculated from the dilution of initial dialysate albumin concentration, as done by Krediet et al. ^[95],

where ${c}_{D, 0}$ and ${c}_{D, t}$ are the dialysate albumin concentrations at time 0 and time $t,$ which is the end of the dwell, ${c}_{D, g}$ is the geometric mean of the dialysate albumin concentration, $\sqrt{{c}_{\mathrm{D}, 0}*{c}_{\mathrm{D}, t}}$, and ${V}_{0}$ and ${V}_{t}$ are the intraperitoneal volumes at time 0 and $t$.

For CFPD, Öberg et al. theoretically derived the following relation for the UF rate ^[38],

where ${J}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}$ is the fill volume ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$), $MTA{C}_{\mathrm{g}\mathrm{l}\mathrm{u}}$ is the diffusion capacity ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$), and ${U}_{\mathrm{m}\mathrm{a}\mathrm{x}}$ is a function of glucose concentration ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$). The equation slightly overestimates the ultrafiltration.

where $RT$ is the product of gas constant and the temperature in degree Kelvin $\left(mmHg\;mmo{l}^{-1}L\right)$, OCG is the osmotic conductance to glucose $\left(ml\;{min}^{-1}mmH{g}^{-1}\right),$ and 3.1 is a constant to account for the lymphatics, opposing forces, and hydrostatic pressure gradient for other solutes in plasma.

Gotch also provided an empirical formula for CAPD using a dextrose solution ^[96],

where $\%D$ is the w/w dextrose solution (dimensionless) and $t$ is the time elapsed. For other dialysate solutions, similar curves for total UF versus time at different concentrations can be drawn to derive an analytical equation. The equation can then be used to determine the dialysate solution for the desired ultrafiltration.

Depending on the importance of ultrafiltration in the objective of the modeling efforts, one can opt for the simplified UFR values such as Eqs (21) and (23) or for patient specific effects Eqs (24) and (28). Eqs (25) or (26) can provide UFR estimates with sufficient precision but require measuring the concentration of albumin at different time-points, or the estimation of OCG and diffusion capacity of glucose, respectively.

4.   Solute flow calculations

Solute flow across a semipermeable peritoneal membrane occurs because of two simultaneous processes. The first is diffusion due to the (electro)chemical gradient between the peritoneal cavity and the blood plasma. Glucose and bicarbonate (and/ or lactate), which are usually present in high concentrations in the dialysate, move to the plasma. Other solutes such as potassium, phosphate, and toxins move from the plasma to the peritoneal cavity. The second is "convection" due to the water flow ("ultrafiltration") that drags solutes across the membrane (from the plasma to the peritoneal cavity). There have been many ways of describing solute flow in the literature since small solute clearance was the main priority of many early PD models. In the following sections, we discuss how the solute flows can be calculated depending on the complexity required. In subsections, we also discuss how to obtain parameters for solute flow equations.

The Patlak equation captures such a two-part transport of solutes across the thick inhomogeneous peritoneal membrane (see Eq (29), Figure 5) ^[97]. The Patlak equation is essentially an extension of the second equation proposed by Kedem and Katchalsky's thermodynamics-driven volume and solute flux across a thin semi-permeable membrane ^[54]. First, there is a term that models the diffusive flux across the membrane due to the concentration gradient, given by Fick's first law. The second part is the convective transfer arising from ultrafiltration.

where $D$ is the diffusion coefficient (${\mathrm{m}}^{2}{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$), and $A$ is the effective surface area available for diffusion (${\mathrm{m}}^{2}$), $C$ is the intramembrane concentration ($\mathrm{m}\mathrm{o}\mathrm{l}\mathrm{ }\;\mathrm{m}{\mathrm{l}}^{-1}$), and ${J}_{\mathrm{v}}$ is the volume flow ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$).

Rearranging and integrating the ordinary differential Eq (29) gives

where ${c}_{\mathrm{p}}$ and ${c}_{\mathrm{D}}$ are the plasma and dialysate solute concentration, respectively, $Pe$ is the Péclet number, which is a ratio of diffusional and convectional mass transfer (dimensionless), $MTAC$ is the diffusion capacity ("mass transfer area coefficient") ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$), $\sigma$ is the reflection coefficient (dimensionless), and ${J}_{\mathrm{v}}$ is the volume flow ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$). For details of the derivation, refer to ^[74].

There are other representations of diffusion and convection in various models. Pure diffusion models are sometimes used to define solute transport such as those by Leypoldt et al. and Villarroel et al. ^[42,98,99].

where $V$ is the peritoneal dialysate volume.

Babb et al. used the following equation to represent the solute flux ^[58,84] (Eq (33)),

where $SiCo$ is the sieving or transmittance coefficient or membrane selectivity of the particular solute and $UFR$ is the ultrafiltration rate ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$, for calculation of parameter, see Section 3.1.7). Different interpretations of $\boldsymbol{c}$ (mol ml^-1) are reported in the literature:

             $\boldsymbol{c} = {c}_{p}$                   when UF is considered from only one direction ^[58,100] or UF is large ^[84]

$\boldsymbol{c} = \left(1-f\right){c}_{p}+f{c}_{D}$      when $\boldsymbol{c}$ is considered as the intramembrane solute concentration ^[40,85,94]

         $\boldsymbol{c} = \frac{{c}_{p}+{c}_{D}}{2}$                 when UF is low ^[84]

here, $f$ is a function of the Péclet number, which is a ratio of convection to diffusion, given by

Graff and Fugleberg et al. compared six models of the peritoneal solute transport of urea, creatinine, glucose, potassium, and phosphate ^{[86,87,88,89,90,91]}, using Eq (35)

where $fct$ is the equilibrium ratio for solute concentration in dialysate and plasma concentration ($\frac{{C}_{\mathrm{D}}}{{C}_{\mathrm{p}}}$), $SiCo$ is the sieving coefficient, $\boldsymbol{c}$ is the intramembrane solute concentration to account for the non-lymphatic convective transfer across the peritoneal membrane, and ${C}^{*}$ is the concentration in the lymph vessels, which depends on the direction of the lymphatic flow (${C}^{*} = {c}_{\mathrm{D}}$ if flow is from peritoneal cavity to lymph (vessels) or else ${C}^{*} = {c}_{\mathrm{p}}$). They found that glucose transport is purely diffusive while lymphatic flow was important in urea and creatinine transport. Non-lymphatic convective transport is important for urea, creatinine, potassium, phosphate, and sodium with a sieving coefficient close to 1 (all molecules passing through). Their models showcase that one general model cannot be applied across all body solutes but differs depending on the size, concentration, and diffusivity of solutes.

The lymphatic flow is denoted by $L$. Other studies that have included lymphatics flow are Waniewski et al. ^[101], who show that MTAC values are underestimated for total protein if lymphatic flow is neglected (Eq (36)).

The fitting constant $f$ is usually to be fitted for all solutes separately, but an assumption of 0.5 works for most small solutes excluding sodium ^[92].

So far, we have mainly discussed compartmental ODE models. Some efforts have included PDE descriptions of the PD process. For example, multi-compartment models such as those by Flessner et al. ^[39,102], consider the mass transfer from the capillaries into the surrounding peritoneal tissue and lymphatics to be contributors for solute flux in the peritoneal cavity.

where $af$ is the fractional peritoneal surface area in contact with the fluid (Section 3.1.5), ${c}_{t}$ is the concentration of the solute in the surrounding tissue, $S$ is the surface area of the peritoneal membrane (${\mathrm{m}}^{2}$) and, ${D}_{t}$ is the diffusion coefficient in the tissue (${\mathrm{m}}^{2}{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$), $\tau$ is the tortuosity (dimensionless), $r$ is the retardation factor (dimensionless), $L$ is the lymphatic flow rate ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$), ${c}_{p}$ Is the blood plasma concentration ($\mathrm{m}\mathrm{o}\mathrm{l}\;\mathrm{ }\mathrm{m}{\mathrm{l}}^{-1}$), and ${J}_{v}'$ is the local volume flux at distance $x$ into the tissue ($\mathrm{m}\mathrm{l}\;{\mathrm{m}}^{-2}{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$).

where $a$ is the capillary surface area per unit tissue (${\mathrm{m}}^{-1}$) and ${x}_{t}$ is the thickness of the tissue ($\mathrm{m}$). The distributed model is used to study fluid absorption and swelling of peritoneal tissue during peritoneal dwell ^[102,103]. Spatial models can also be utilized to monitor the penetration of medications, such as chemotherapeutics, administered into the peritoneal cavity. Other spatial models can be found in ^[104,105].

Gotch added the drain flow rate ($= {J}_{fill}+UFR$) to derive the solute flux into the peritoneal cavity for single pass CFPD as ^[96] (Eq (39)),

Öberg et al. used the three-pore model for both volume and solute flux for CFPD ^[38]. Volume fluxes are represented in Eqs (6)–(8). Using the volume flux and Patlak equation (Eq (30)), they calculated the solute fluxes as Eq (40),

where ${J}_{\mathrm{s}\mathrm{S}}$ and ${J}_{\mathrm{s}\mathrm{L}}$ are the solute fluxes over the small and large pores (solute flux over the ultrasmall pores is non-existent).

Icodextrin is a mixture of different dextrin molecule fractions (molecular weight ranging from 12–20 kDa). Briefly, for icodextrin modeling, an additional term is added that refers to the hydrolysis of the higher molecular weight icodextrin fraction, caused by the enzyme $\alpha$ -amylase.

where ${r}_{i}^{j}$ is the appearance rate of icodextrin fraction $i$ due to hydrolysis of fraction $j$ ($ml\;{U}^{-1}{min}^{-1}$), ${c}_{D, amylase}$ is the amount of amylase (U), and ${c}_{D, Ic{o}_{j}}$ is the fraction $j$ of icodextrin ^[6,106]. 1 unit (U) is the amount of enzyme that catalyzes the reaction of 1 µmol of substrate per minute.

From Eq (30),

where $Pe$ is calculated from Eq (31).

${J}_{\mathrm{v}\mathrm{C}}$, ${J}_{\mathrm{v}\mathrm{S}}$, and ${J}_{\mathrm{v}\mathrm{L}}$ are the volume fluxes over the three types of pores, ${J}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}$ is the fill flow rate and ${c}_{\mathrm{D}0}$ is the solute concentration in the fresh dialysate.

4.1.   How to calculate the solute flow parameters

After calculating the volume flow, one can use it to calculate the solute flux. Standard peritoneal permeability analysis (SPA) is a standardized tool used to assess the membrane transport properties of a specific patient. These assessments are then used to determine the PD prescription for the specific patient ^[107]. For a SPA test, a static dwell is used with regular dialysate sampling and blood sampling ^[108].

4.1.1.   Diffusion Coefficient, $D$

The diffusion coefficient, $D$ (${\mathrm{m}}^{2}{\mathrm{s}}^{-1}$) can be calculated for charged and uncharged particles through liquid flow at low flow rates,

where $\mu$ is the electric mobility of the solute $\left(\frac{{\mathrm{m}}^{2}}{\mathrm{V}\mathrm{s}}\right)$, ${k}_{B}$ is the Boltzmann constant (Joule, J per Kelvin), $T$ is the temperature ($\mathrm{K}$), $q$ is the charge of the solute (Coulomb), $\eta$ is the dynamic viscosity $\left(\mathrm{P}\mathrm{a}\cdot \mathrm{s}\right),$ and ${r}_{\mathrm{s}}$ is the solute radius ($\mathrm{m}$). Equation (44) is for charged solutes and Eq (45) is for uncharged solutes.

4.1.2.   Mass transfer area coefficients, $MTAC$

The capacity for diffusion MTAC ($\mathrm{m}\mathrm{l}\;{\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$) is the maximal diffusive clearance and can be calculated from

where A₀/Δx is the unrestricted surface area to diffusion length ratio (in cm). Typically, a patient with an average peritoneal solute transfer rate has an A₀/Δx of 25,000 cm, whereas a value < 16,000 cm or > 40,000 cm may indicate slow- and fast peritoneal transport, respectively. The factor A/A₀ represents the diffusive hindrance factor and is usually estimated using the equation by Mason, Wendt, and Bresler ^[70], as follows

where λ is the solute to membrane pore radius ratio. MTAC may also be estimated from experimental data ^[92,100,109]

where ${V}_{t}$ and ${V}_{0}$ is the intraperitoneal volume at time $t$ and 0, $f$ is from Eq (34) and ${c}_{D, 0}$ and ${c}_{D, t}$ is the dialysate solute concentration at time 0 and t, respectively, and ${c}_{p}$ is the plasma solute concentration.

Keshaviah et al. compared different functions for MTAC and found that the parabolic and negative exponential functions for urea, creatinine, and glucose best fit the dialysate volume profile ^[77],

where ${a}_{1}, {a}_{2}, {a}_{3}$ are fitting constants that are different for each of the three solutes.

4.1.3.   Sieving coefficients, $SiCo$

The sieving coefficient, $SiCo$ (dimensionless), is usually fitted in a model ^{[58,86,87,88,89,90,91]}, microscopy study ^[110] or calculated from experimental observations such as the formula derived by Chen et al. ^[67],

where ${V}_{t}$ and ${V}_{0}$ is the intraperitoneal volume at time $t$ and 0, respectively, ${c}_{\mathrm{D}, 0}$ and ${c}_{\mathrm{D}, t}$ are the dialysate solute concentration at time 0 and t, respectively, and ${c}_{\mathrm{p}}$ is the plasma solute concentration, ${V}_{\mathrm{U}\mathrm{F}}$ is the net ultrafiltration volume, and $C{l}_{\mathrm{a}\mathrm{l}\mathrm{b}\mathrm{u}\mathrm{m}\mathrm{i}\mathrm{n}}$ is the clearance of albumin.

Rippe et al. also calculated the sieving coefficient incorporating hematocrit ^[111],

where ${H}_{0}, {H}_{t}$ are the initial and final hematocrit and ${c}_{0}, {c}_{t}$, and ${c}^{*}$ are the initial, final and average plasma proteins concentrations.

It can also be simply calculated from the reflection coefficients calculated in Section 3.1.3,

5.   Summary of models for peritoneal dialysis

In this section, we provide a comprehensive look at the different models of PD that have been developed throughout the years (1966–2019) (Figure 6). We looked for publications that mentioned "kinetic modeling", "mathematical model", "computational model" with "peritoneal dialysis". In total, we found 26 distinct models of PD. We have mentioned what kind of compartments were considered along with what were the main aim/hypotheses of the modeling effort. Modeling was used as a tool to establish the efficiency of a PD model (5 out of 26) ^{[38,41,78,94,112]} while four models were designed to understand the fundamentals of solute and volume transport ^{[42,58,98,113]}. We have also characterized the types of transport processes were included in the volume and solute transport equations of the model. Three of the 26 models were purely diffusive ^[41,42,112] while other models (23 out of 26) included convection and diffusion into the surrounding tissue ^[83,114,115]. Some models were generalized for any type of solute ^[38,84,115] while some were prepared for specific solutes only ^{[86,87,88,89,90,91,99]}. Different PD modalities can be found in ^{[7,38,51,78,93,116]}, and models for glucose-sparring solutions are in ^[6,106,117]. Table 1 shows a consecutive overview of 26 PD models that have been published, and figure 6 shows the evolution of PD models. Depending on the solute(s) or drug and which part of the PD process one is interested in, one can choose any of these models.

6.   Test case for developing a mathematical model of PD

Many PD-related questions can be answered with mathematical modeling, for example:

- What is the optimal glucose concentration for PD, balancing ultrafiltration (volume and efficiency) versus adverse effects (e.g., increased peritonitis risk)?

- What is the influence of two catheters vs one single lumen catheter on solute clearance in continuous flow vs tidal PD?

- How does the intraperitoneal dialysate volume affect solute clearance in CAPD or CFPD?

- Which flow rates are ideal in CFPD?

- What is the clearance of a particular solute/drug for a particular PD modality (CAPD, APD, or CFPD)?

- How does tidal PD with partial drainage of the solution, leaving a residual volume in the peritoneal cavity, affect solute clearance as compared to complete drainage?

In Table 2, we demonstrate how we can create and implement our mathematical model if we know which research question we want to answer, what experimental data or literature data we possess, and what kind of math is necessary for the model. With the realized model, we can then play with the parameters to determine an optimal treatment scenario for the patient, as we can analyze both the short-term and long-term effects of a session.

7.   Discussion

In this review, we give an overview of the physical principles that govern peritoneal dialysis and summarize the essential (differential) equations for volume and solute flux that are required to model peritoneal dialysis. These models are sometimes based on simple principles and parameters are lumped together to study a compartmentalized version of the body ^{[41,58,86,87,88,89,90,91]}, while other models are very complex and require many parameters to capture the physical processes in depth ^[38,102,115]. In this review, we chose to focus more on the time dimensional models, and thus a spatial model discussion is out of scope. However, spatial models can be found in ^{[39,83,103,104,105,114,115]}. We have also not included other glucose-sparring solutions that have been modelled, and valuable information on this can be found in ^[6,106,117]. We also list different ways the in which various simulation parameters are obtained. Most often, parameters are derived directly from patient data (effective peritoneal surface area, dialysate clearance, etc.), some of the parameters are derived from analytical formulae (reflection coefficient, hydrostatic, and osmotic pressure difference), others are derived from the model itself (ultrafiltration, pore contribution to ultrafiltration, and fractional hydraulic conductance), and some parameters may be fitted to obtain the best (patient-specific) interpretation of the dialysate and plasma solute concentration. Each model can be modified to fit the required patient data and get outputs such as ultrafiltration rate, mass transfer area coefficients, and residual volumes. From Eqs (15) and (16), we do not see a significant difference in the estimation of sieving coefficients (Supplementary Figure 1). This underlines the fact that multiple studies have been done to converge on the most optimal values of reflection coefficients. However, we must note that reflection coefficients and thus, sieving coefficients may differ with prolonged exposure of the membrane to glucose and glucose-degradation products. This highlights a critical gap in the literature, necessitating further investigation to refine these parameters for enhanced model accuracy.

The human body is a complex system, which makes it difficult to model. Modeling its various parts should eventually lead us to understand the whole. Over the years, PD models have generally included more and more transport processes. The first PD models were purely diffusive which is representative of small molecular transport such as that of urea and creatinine ^[41,42,112]. As the importance of middle-sized molecules became appreciated, convection was added to the models, as this is an important route to eliminate these compounds ^[58,99]. Lymphatic absorption ^[38,96] and peritoneal tissue surface area ^[115] have since been added to the models. As per the new ISPD guidelines, ultrafiltration is a crucial parameter to determine the efficacy of a PD treatment ^[132]. Thus, we see that recent models are also trying to model time-dependent ultrafiltration in patients ^[38]. This requires rigorous evaluations and comparisons of ultrafiltration rates. In Section 3.1.7, we discuss multiple ways to obtain ultrafiltration rate and which efforts can be made to make the best estimate.

With an increase in computational power and physical understanding, the model complexity can be increased to create a better picture of the underlying mechanisms in the patient's peritoneal cavity. Herein lies a caveat. The increased complexity makes the model space hyper-parameterized. Thus, finding efficient parameter combinations that satisfy physiological relations is difficult to achieve without multiple assumptions. This can be achieved by fitting parameters to clinical data (for example, in ^[133], Stachowska-Pietka et al. fitted 16 parameters while assuming 9 parameters to be fixed) and calibrating the model in a different scenario (they showed that the predictions of interstitial concentration of mannitol in rat abdominal wall were in agreement with experiments ^[133]). Future models can be improved to include metabolism of glucose in the peritoneum, and vascular inflammation ^[134] and how these affect the peritoneal membrane, positioning of catheters in patients and cellular contributions to the solute clearances (for example for glucose ^[135]). Future PD models can also be made spatial (using for example partial differential equations) to understand the spatial influence of flow and the solute gradient. How does the continuous flow change the boundary layer of the peritoneal membrane? What happens to the residual volume? What fraction of recirculation occurs just at the tip of the catheter? Can patient-reported outcome measures be linked to PD (efficacy) parameters? Partial differential equation modeling (used for space localisation) of the peritoneal cavity may help us answer these questions. Agent based modeling (ABM) can also be helpful in this aspect. ABM models take into account different agents (e.g., patient characteristics, residual kidney function, smoking, dietary habits, transport status, sex, weight, and osmotic agent, etc.) to assess the problem and make decisions based on a complex behavioral pattern. ABM models have already been used to determine the optimum treatment pathway for HD patients based on patients', nephrologists', and surgeons' attributes ^[136]. From recent developments in PD, we definitely see that patient membrane characteristics and preferences play a huge role in managing the patient. Efforts have already been made in the past to suggest the best mode of PD treatment for different patient characteristics ^{[38,85,99,112,137]} but more work needs to be done on personalizing a single PD mode based on the sex, age, weight, peritoneal characteristics, and residual kidney function of the patient. Nevertheless, there are many aspects that a virtual PD system cannot cover, such as catheter dysfunction and nonadherence to medication.

We can picture the growth of mathematical modeling of PD in two complementary directions: Fundamental understanding and personalization. More technical components can be introduced to the existing models, and parameters can be analyzed to facilitate a faster understanding of the new devices and optimize them before market entry. One aspect that has been little explored is to employ computational modeling to evaluate different bio-compatible osmotic agents replacing glucose to avoid the adverse effects of glucose and glucose by-products in PD patients. Recent advances in modeling includes using PD models to create virtual clinical trials and personalize treatments specific to the patient with the help of their clinicians.

Virtual Patients

The interest in virtual patients is growing, as we have seen with PBPK (physiologically based pharmacokinetic) and PKPD (pharmacokinetic and pharmacodynamic) modeling studies becoming mandatory in drug studies ^{[138,139,140]}. Similarly, in silico modeling could become essential to medical device development. To realize this, we need digital patient twins. Enabling the use of a virtual platform to test and develop medical devices, would reduce risk and regulatory burden. The FDA also has directed funds towards the establishment of computational modeling as a regulatory tool ^[141]. A growing virtual database (including patients and animals) could mean that with a computational model, we could identify the risk groups and drug side effects, reduce animal testing, and establish safe protocols.

Complexity of PD devices

Currently, there are several CFPD devices under development using continuous sorbent-based dialysate regeneration (preclinical or clinical developmental phase) ^[142,143]. To realize these devices in mathematical terms, future work needs to focus on adding flow through a sorbent chamber. There needs to be an evaluation of dead volume- and recirculation-related loss of efficiency in the system due to application of rapid flow cycling via a single lumen catheter. Considering all the potential parameter settings, which likely influence each other, investigating all parameter combinations is very time consuming, requiring large clinical trials, with many resources and tremendous planning ^[144]. Computational modeling can help in rapid optimization and avoid unnecessary scenario testing through careful calibration and validation.

Personalization

Other problems with clinical trials include the lack of randomization of patients (inclusivity in terms of sex, age, etc.), difficulty in studying long term effects (due to patient drop out e.g., due to kidney transplantation), lack of blinding, the learning curve (an already existing system may be easier to handle), and general reluctance (by both doctors and patients) to try new technologies ^[29]. Chronic kidney disease patients also use various medications, which may interfere with the efficiency of the PD session. With mathematical modeling, we can optimize PD treatments ^[56,145] and conduct virtual PD trials with an inclusive patient spread (age, sex, stage of renal disease, transport parameters, intra-abdominal volume, etc.) and study the short- and long-term effects of different modes of PD. Moreover, PD prescription models such as PatientOnLine ^[34,146] and PD ADEQUEST ^[147,148], developed to optimize and personalize PD prescription, have also been validated in multicentre studies. However, the focus of these prescription models is on optimizing the dialysis dose and ultrafiltration and not on patient quality of life. Incorporation of more patient specific factors ^{[149,150,151]} into the model may further personalize PD prescription and contribute to patient well-being.

8.   Conclusions

In silico modeling is a powerful tool that can be used to understand gaps in knowledge of both new and old PD devices and facilitates their quick and efficient transition from in vitro to in vivo to patient care. In this review, we looked at the different modeling approaches explored through the years to model PD and the essential components in a PD compartmental model. We are optimistic that a joint effort of computational modelers and clinicians can help in the technical improvement of PD technology and in personalizing PD treatment to ensure a higher quality of life for patients.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgements

This work is supported by the partners of Regenerative Medicine Crossing Borders (RegMedXB), a public-private partnership that uses regenerative medicine strategies to cure common chronic diseases, by the Dutch Kidney Foundation and Dutch Ministry of Economic Affairs by means of the PPP Allowance made available by the Top Sector Life Sciences & Health to stimulate public-private partnerships (DKF project code PPS08), by a grant from the Dutch Kidney Foundation (22OK1018) and by the European Union (CORDIAL, Horizon 2020 research and innovation program, grant agreement no. 945207).

Conflict of interest

The authors declare that there is no conflict of interest.

Acronyms

ACE: Angiotensin-Converting Enzyme

APD: Automated Peritoneal Dialysis

AQP1: Aquaporin-1

CAPD: Continuous Ambulatory Peritoneal Dialysis

CFPD: Continuous Flow Peritoneal Dialysis

CCPD: Continuous Cycling Peritoneal Dialysis

CT: Computed Tomography

ESKD: End-Stage Kidney Disease

HD: Hemodialysis

IPV: Intraperitoneal Volume

MTAC: Mass Transfer Area Coefficients

OCG: Osmotic Conductance to Glucose

ODE: Ordinary Differential Equation

PDE: Partial Differential Equation

PD: Peritoneal Dialysis

PK: Pharmacokinetic Model

SAPD: Standard Automated Peritoneal Dialysis

SPA: Standard Peritoneal Assessment

TPM: Three-Pore Model

UF: Ultrafiltration Volume

UFR: Ultrafiltration Rate

Glossary

A = fitting constant for exponential decrease of UF over time (dimensionless)

$af$ = Fraction of peritoneal membrane in contact with fluid (dimensionless); $af = 16.18*\frac{1-{e}^{-0.00077*V}}{13.3187}$ ^[78]

${c}_{\mathrm{D}}$ = dialysate solute concentration (mmol/l or mmol/m³), proteins are displayed in g/L or g/dL.

${c}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{n}}$ = drain solute concentration (mmol/l)

${c}_{\mathrm{p}}$ = peripheral vein plasma water solute concentration (mmol/l), proteins are displayed in g/L or g/dL.

$\boldsymbol{c}$ = intramembrane solute concentration (mmol/l)

$D$ = diffusion coefficient (dimensionless)

$f$ = function of Pe (dimensionless)

$fct$ = equilibrium ${c}_{p}/{c}_{D}$ (dimensionless)

${J}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{i}\mathrm{n}}$ = drain flow rate (l/min)

${J}_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}}$ = fill flow rate (l/min)

${J}_{\mathrm{v}}$ = volume flux (l/min)

${J}_{\mathrm{s}}$ = solute flux (mmol/min)

$K$ = permeability coefficient (dimensionless)

$L$ = lymphatic flow (l/min)

${L}_{p}$ = hydraulic conductivity (l/(min.cm².mmHg))

$MTAC$ = mass transfer area coefficients (ml/min or m³/min)

$P$ = hydrostatic pressure (mmHg)

$Pe$ = Peclet number (dimensionless) to determine the importance of convection over diffusion

$r$ = retardation factor (dimensionless)

$S$ = peritoneal surface area (m²)

$SiCo$ = sieving coefficient (dimensionless)

$t$ = time of session (hr)

$UFR$ = ultrafiltration rate (l/min)

$V$ = intraperitoneal volume (l)

${V}_{r}$ = residual volume (l)

${V}_{fill}$ = initial fill volume (l)

${V}_{drain}$ = drain volume (l)

$W$ = body weight (kg)

Greek symbols

$\alpha$ = contribution to ultrafiltration coefficient (dimensionless); for ultrasmall pores ${\alpha }_{C}$ = 0.02, small pores ${\alpha }_{S}$ = 0.9 and for large pores ${\alpha }_{L}$ = 0.08

$\beta$ = time constant for exponential decrease of UF over time (dimensionless)

$\lambda$ = solute radius/membrane pore radius (dimensionless)

$\sigma$ = reflection coefficient (dimensionless); calculated for different solutes.

$\pi$ = oncotic pressure (mmHg)

$\mathrm{\varPhi }$ = equilibrium partition coefficient (dimensionless)

$\tau$ = tortuosity factor (dimensionless)

Subscripts

$glu$ = glucose

$t$ = at time t

$fill$ = filling the peritoneal cavity at time 0 (static dwell) or during the session (CFPD, APD)

$drain$ = drain of the peritoneal cavity after (static dwell) or during the dwell (CFPD, APD)

$s$ = solute flux

$v$ = volume flux

$C$ = ultrasmall pores

$S$ = small pores

$L$ = large pores

$r$ = residual