Numerical Modeling and Simulation of Blood Flow in a Rat Kidney: Coupling of the Myogenic Response and the Vascular Structure

Abstract

A numerical simulation was carried out to investigate the blood flow behavior (i.e., flow rate and pressure) and coupling of a renal vascular network and the myogenic response to various conditions. A vascular segment and an entire kidney vascular network were modeled by assuming one single vessel as a straight pipe whose diameter was determined by Murray’s law. The myogenic response was tested on individual AA (afferent artery)–GC (glomerular capillaries)–EA (efferent artery) systems, thereby regulating blood flow throughout the vascular network. Blood flow in the vascular structure was calculated by network analysis based on Hagen–Poiseuille’s law to various boundary conditions. Simulation results demonstrated that, in the vascular segment, the inlet pressure $P_{\mathrm{inlet}}$ and the vascular structure act together on the myogenic response of each individual AA–GC–EA subsystem, such that the early-branching subsystems in the vascular network reached the well-regulated state first, with an interval of the inlet as $P_{\mathrm{inlet}}$ = 10.5–21.0 kPa, whereas the one that branched last exhibited a later interval with $P_{\mathrm{inlet}}$ = 13.0–24.0 kPa. In the entire vascular network, in contrast to the $P_{\mathrm{inlet}}$ interval (13.0–20.0 kPa) of the unified well-regulated state for all AA–GC–EA subsystems of the symmetric model, the asymmetric model exhibited the differences among subsystems with $P_{\mathrm{inlet}}$ ranging from 12.0–17.0 to 16.0–20.0 kPa, eventually achieving a well-regulated state of 13.0–18.5 kPa for the entire kidney. Furthermore, when $P_{\mathrm{inlet}}$ continued to rise (e.g., 21.0 kPa) beyond the vasoconstriction range of the myogenic response, high glomerular pressure was also related to vascular structure, where $P_{\mathrm{GC}}$ of early-branching subsystems was 9.0 kPa and of late-branching one was 7.5 kPa. These findings demonstrate how the myogenic response regulates renal blood flow in vascular network system that comprises a large number of vessel elements.

1. Introduction

Kidney blood flow plays a vital role due to its supply to the basic units of the organ, the nephrons, to carry out organ functions (i.e., ultrafiltration, reabsorption, and urine secretion) [1]. In a kidney, the distribution of renal blood flow is uneven in temporal and spatial terms [2,3,4,5,6,7,8,9,10,11]. First, renal blood flow is temporally changed by the so-called autoregulation process, in which the diameters of small arteries respond to changes in flow conditions [7,8], among which the myogenic response is an important regulatory mechanism for the elevation of glomerular pressure [9,10]. Its failure or dysregulation will cause hypertension, which plays a major role in chronic kidney disease [10,11]. After this, most of the blood enters the renal cortex [2], whereas the medulla receives the remaining 10%. Even within the cortex, it still exhibits superficial and juxtamedullary flow differences depending on its location [3,4,5,6]. Consequently, in experimental studies, flow rate and pressure on different sites, as well as under various inlet conditions, were measured in order to fully describe them [2,4,5,6].

Mathematical modeling and numerical simulations are useful tools for investigating renal function and flow. Following the excellent simulation work of the circulation of the entire body [12], multiple numerical studies of the kidney were conducted. They include the numerical modeling of renal flow ultrafiltration and solute reabsorption [13], CKD development by functional loss of nephrons and hypertension [14], blood flow regulation by the myogenic response and tubuloglomerular feedback [15,16,17], the neighboring nephrons interactions through a vascular tree [18,19], and fluid dynamics studies of the renal flow [16,18,20].

As they form the framework of blood flow, vascular structures play an important role in determining it, which has already been demonstrated via microcirculation [21,22]. Several numerical simulations of renal blood flow tried to use a realistic vascular model, among which micro-computed tomography data by Nordsletten et al. provided mean values of vessel diameter and length by Strahler orders. [23]. Based on those data, Kleinstreuer et al. [20] investigated autoregulation by modeling the symmetrical vascular network, and Postnov et al. [24] simulated the effects of arterial structure on the interactive behaviors between nephrons. In addition, the asymmetric branching style of small vessels was experimentally observed and numerically simulated [25,26].

However, to date, it is still not possible to fully characterize renal flow. Few numerical studies on renal flow and function were based on complex vascular structures, instead tending to use simplified symmetric models [24,27]. Our recent study introduced a detailed renal vascular network [28]. However, considering the sensitivity of the kidney to high pressure [10,11,13], renal blood flow cannot be fully elaborated without the employment of multiple pressure input and regulation systems.

In this study, we numerically simulated the coupling of the myogenic response and the renal vascular structures. With the regulation of the myogenic response, blood flow simulations were performed via network analysis on the basis of Hagen–Poiseuille’s law, and a single blood vessel was modeled as a straight pipe for the purpose of an electric circuit investigation. Our simulation assessed the flow rate and pressure of the renal vascular network coupled with the myogenic response to various inlet conditions.

2. Materials and Methods

2.1. Model of the Vascular Network

Two vascular networks were utilized in this study; the first was a typical asymmetrical vascular segment based on Marsh et al.’s measurements [24], where fourteen afferent arterioles (AAs) originated on an interlobular artery (ILA) from the nonterminal sites and six AAs branched at the end of the ILA through two intermediate vessels, as shown in Figure 1. The other one was a model of the entire renal vascular network, the EP2 model referred to in our previous study [28].

The individual vessels were modeled as straight pipes. at a branching point, all of the vessels’ diameters conformed to Murray’s law [29] which is described as

$$d_{\mathrm{P}}^3=d_{\mathrm{D}1}^3+d_{\mathrm{D}2}^3,$$

(1)

where $d_{\mathrm{P}}$ [m] is the diameter of the parent vessel, and $d_{\mathrm{D}1}$ [m] and $d_{\mathrm{D}2}$ [m] are the diameters of the two daughter vessels, respectively. In the vascular segment, lengths of the ILA and the two intermediate vessels were obtained via measurement [19] as 556 µm, 206 µm and 96 µm, respectively. The length and diameter were set as 200 µm and 20 µm, respectively, for both the AA and EA [19,23]. For the entire renal vascular network, we employed symmetrical (Sy.) and asymmetrical (Asy.) styles to reproduce the renal vasculature. Detailed information is listed in Appendix A 

Table A1. Morphological parameters of the arterial structure in the asymmetrical (Asy.) and symmetrical (Sy.) models. BS stands for the branching style of the vessels, ranging from renal artery to AAs, where “L” and “S” show the symmetrical and asymmetrical branching styles. $D$, $N_{\mathrm{total}}$ and St are vessel diameter, vessel number and Strahler order, respectively. “Measurement” indicates the experimental data from Nordsletten et al. [23] in Table 2 and Table 5.

Asy. model	BS	L2-L2-L2		L2-L2-L2		L2-L2-L2-S6		S10
	D [µm]	628	314		157		43		20
	$N_{\mathrm{total}}$	1	8		64		3072		30,720
	St	11	8–10		5–7		1–4		0–1
Sy. model	BS	L2-L2-L2		L2-L2-L2		L2-L2-L9		L6
	D [µm]	628	314		157		38		20
	$N_{\mathrm{total}}$	1	8		256		4608		27,648
	St	10	7–9		4–6		1–3		0–1
Measurement	$D$ [µm]	432	278–382		106–172		20–88		20
	$\mathrm{St}$	9–10	8–9		6–7		2–6		0–1

.

The resistance $R_I$ [$\mathrm{Pa}·\mathrm{s}·{\mathrm{m}}^{-3}$] of one single vessel I modeled as a straight pipe with a given diameter $d_I$ [m] and a length of $L_I$ [m] was calculated as

$$R_I=\frac{128{\mu }_I{L}_I}{\pi d_I{}^4}$$

(2)

where the blood viscosity ${\mu }_{\mathrm{I}}$ [Pa·s] of one single vessel was set as a function of the vessel diameter $d_I^{\prime }=1.0\times {10}^6{d}_I$ μm, as follows:

$${\mu }_I={\mu }_0\left[1+\left\{6\exp \left(-0.085d_I^{\prime }\right)+2.2-2.44\exp \left(-0.06d_I^{\prime }{}^{0.645}\right)\right\}{\left(\frac{d_I^{\prime }}{d_I^{\prime }-1.1}\right)}^2\right]{\left(\frac{d_I^{\prime }}{d_I^{\prime }-1.1}\right)}^2$$

(3)

by assuming a Hct (hematocrit) level as 45% [24,29] and ${\mu }_0$ as a constant value of 1 $\mathrm{mPa}·\mathrm{s}$ [28].

2.2. AA-GC (Glomerular Capillaries)-EA (Efferent Artery) Subsystem

The myogenic response model acted at each AA to constrict the vessel [30,31,32], demonstrating the regulation of pressure at individual subsystems [13,27] as

$$R_{\mathrm{M}}=G_d\left(\frac{P_{\mathrm{CRA}}}{P_{\mathrm{aso}}}-1\right)\left(R_{\mathrm{B}}+R_{\mathrm{MA}}+R_{\mathrm{TGF}}+R_{\mathrm{G}}+R_{\mathrm{Eff}}\right)=G_d\left(\frac{P_{\mathrm{GC}}}{P_0}-1\right)R_{\mathrm{AA}}$$

(4)

where $R_{\mathrm{M}}$ is the AA resistance under regulation of the myogenic response; $G_d$ is a scaling coefficient that expresses the regulatory ability of the myogenic response, and was set at 0.3 [13]; $P_{\mathrm{CRA}} \mathrm{and} P_{\mathrm{aso}}$ are the calculated value and criterion of pressure P at the ILAs; $P_{\mathrm{GC}}, P_0$ are the calculated value and criterion of pressure P at the AAs; and $R_{\mathrm{B}}, R_{\mathrm{MA}}, R_{\mathrm{TGF}}, R_{\mathrm{G}} R_{\mathrm{Eff}}$ and $R_{\mathrm{AA}}$ are the resistances of the vessels, the ascending myogenic response, TGF, a glomerulus, the EA and the AA. Moore et al. [14] introduced the first equation of Equation (4) to predict resistance change with a simplified vascular network by employing $P_{\mathrm{aso}}$ as the criterion. However, a wider range of interlobular pressure $P_{\mathrm{CRA}}$ would be obtained in this study with the use of detailed vascular structures, which did not exclude the small vascular structures such as interlobular arteries (Figure 1a). Furthermore, the AAs were designated as the only vessel capable of constriction [30,31,32]. Thus, the glomerular pressure $P_0=6.7$ kPa [33,34,35] was assumed more suitable as a pressure criterion to accurately represent the response of the AA to $P_{\mathrm{GC}}$ changes. And then, $R_{\mathrm{M}}$ was approximated to the resistance caused by the constriction of the AA, as described by the second equation of Equation (4).

2.3. Blood Flow Simulation with Network Analysis

Using the resistance value obtained from Equation (2), we carried out a network analysis based on the Hagen-Poiseuille law [21]. Flow rate $Q_I$ was represented with the pressure values $P_i$ and $P_j$ of sites i and j as

$$Q_I=\frac{P_i-P_j}{R_I}$$

(5)

A continuity equation was described at a branching point i shared by N vessels using a flow rate $Q_I$ [m³/s] of vessel I, as follows:

$${\displaystyle \sum }_I^N{Q}_I=\left\{\begin{array}{c}0 \left(\mathrm{branching} \mathrm{point}\right)\\ Q_{\mathrm{inlet}} \left(\mathrm{inlet}\right)\\ Q_{\mathrm{outlet}} \left(\mathrm{outlet}\right)\end{array}\right.$$

(6)

where $Q_{\mathrm{inlet}}$ and $Q_{\mathrm{outlet}}$ represent the inlet and outlet boundary conditions of the simulation, respectively. By solving Equations (5) and (6), the unknown pressure $P_i$ was calculated and then the flow rate $Q_I$ of a single vessel I was derived by substituting $P_i$ into Equation (5).

We used two sets of boundary conditions for the two vascular networks. There was one inlet and seventeen outlets in the segmental vascular network model. The pressure value $P_{\mathrm{inlet}}$ at the inlet was set as a range of 6.5–28 kPa, the pressure value $P_{\mathrm{outlet}}$ at all the outlets was set as a constant (2.3 kPa) [33]. For the entire vascular network, the range of $P_{\mathrm{inlet}}$ was set as a range from 8.0 to 24.0 kPa, with a single outlet $P_{\mathrm{outlet}}=$ 0.3 kPa [33,34]. The flow rate at the inlet, $Q_{\mathrm{inlet}}$, and the flow rates at the outlets, $Q_{\mathrm{outlet}}$, were obtained by solving Equation (5) with the given $P_{\mathrm{inlet}}$ and $P_{\mathrm{outlet}}$. For a given pressure $P_{\mathrm{inlet}}$, flow rate and pressure values were solved (Equation 6); based on the new-calculated $P_{\mathrm{GC}}$ (Equation (4)), AA diameter was constricted accordingly and the myogenic responses were updated. These processes were repeated until every $P_{\mathrm{GC}}$ converged.

3. Results

3.1. The Myogenic Response in a Vascular Segment with 17 AA-GC-EA Subsystems

To investigate how the asymmetrical vascular structures affect the glomerular flow, we divided all 17 AA–GC–EA subsystems into three groups according to the branching sequences of the AA–GC–EA subsystems on the main vessel ILA in Figure 1: G1, the early branching N1–N4; G2, N5–N11, with asymmetric branching; G3, N12–N17 with symmetrical branching at the end of the ILA.

3.1.1. Flow Responses to a Steady Flow Condition

In this section, $P_{\mathrm{inlet}}=$ 12.7 kPa [27,28] was set as the inlet pressure boundary condition at the beginning of the ILA, while 2.3 kPa [33] was set as an outlet at the end of each EA to simulate the normal blood pressure condition of the kidney without constriction from the AAs. Under these circumstances, the blood flow of the vascular network is shown below.

The blood pressure and blood flow values of the vascular segment decreased in the descending direction on the ILA, as shown in Figure 2. The average $P_{\mathrm{GC}}$ of its 17 glomeruli in the AA-GC-EA subsystem was $6.8\pm 0.4$ kPa, and the subsystem′s blood pressure decreased in the same way as its connected ILA, with the early-branching N1 being the largest and the last symmetrical branching N12–N17 being the smallest. Likewise, among the six symmetric branches, N12, N13 and N14 have showed lower glomerular flow pressure values than N15–N17 due to longer intermediate vessels. The distribution of blood pressure resulting from this asymmetric branching style indicates that this approach can generate variations. The current model is validated by the results of the replicated experimental glomerular pressure measurements [33,34,35].

3.1.2. Flow Responses to a Step Condition

In this section, an increasing inlet pressure value ranging from 6.0 to 28.0 kPa was set as the inlet pressure boundary condition at the beginning of the ILA to investigate the geometry-induced effect on the blood flow regulated by the myogenic response. The outlet pressure was set at 2.3 kPa at the end of each EA.

As shown in Figure 3a, the glomerular pressure was well-regulated at a $P_{\mathrm{inlet}}$ value between 13.0 and 20.0 kPa. During this time, blood pressure in the 17 glomeruli was consistently steady at approximately 6.8 kPa, despite $P_{\mathrm{inlet}}$ increasing, indicating a well-regulated state in this study, which is consistent with previous experimental and modeling studies [33,34,35].

However, as shown in Figure 3a, when the input blood pressure was outside of this region, the mean glomerular pressure value began to rise. When $P_{\mathrm{inlet}}$ was in the region of < 13.0 kPa, with increasing the input pressure, the glomerular pressure gradually increased and approached the criterion value (6.7 kPa). In the region of $P_{\mathrm{inlet}}$ > 21.0 kPa, the increase in pressure exceeded the upper limit of vasoconstriction capacity, resulting in a continued increase. It is worth mentioning that, when compared with the curve from other experimental or numerical results [9,32], the values in the two rising regions where $P_{\mathrm{GC}}$ in Figure 3 are located show more variation (error bar in Figure 3), such that the early-branching G1 in the vascular network reached a well-regulated state first with a $P_{\mathrm{inlet}}=$ 10.5–21.0 kPa, whereas G3 that branched last reached a well-regulated state later with $P_{\mathrm{inlet}}=$13.0–24.0 kPa.

Referring to the results in Section 3.1.1, this variation could be attributed to the use of asymmetric vascular networks. As shown in Figure 3b,c, in the elevated region, the numerical value of $P_{\mathrm{GC}}$ and $Q_{\mathrm{AA}}$ of each glomerulus was also arranged in the order of the branches (G1 > G2 > G3). At the same time, the $P_{\mathrm{inlet}}$ of the three groups of glomeruli that reached the well-regulated state were different.

3.2. The Myogenic Response in the Entire Vascular Network

The blood pressure and flow values in the entire renal vascular network were investigated through simulations based on the EP2 model from our previous study [28], of which two generations of asymmetrical branching vessels (six ILAs and ten AAs) generated a vascular network including a set of 60 AA-GC-EA subsystems. Although the myogenic response took effect as a regulator, it still could not completely eliminate the variations in blood pressure and flow rate caused by an asymmetric vascular network. For example, when $P_{\mathrm{inlet}}=$ 13.0 kPa, as shown in Figure 4, during the process when pressure P dropped from arteries to veins, there still were pressure variations among the AAs.

Overall, our model fits the experimental measurements [36,37] well, in which TGF was blocked by furosemide, regulating renal blood flow during the $P_{\mathrm{inlet}}$ rise from 12.0 to 18.5 kPa (Figure 5a), demonstrating the validity of the modeling. However, as mentioned in Section 3.1, the asymmetrical network resulted in pressure P variation. Furthermore, because the number of AAs in one group in this section is 60, which was greater than that in Section 3.1 (18), and the variation ranges of $P_{\mathrm{GC}}$ became correspondingly larger.

As shown in Figure 5b, N1, the first branch, and N60, the last branch, exhibited different pressure values or AA diameter constrictions as the $P_{\mathrm{inlet}}$ values increased. To provide a more complete view, we set up two scenarios to depict the normotensive (13.0 kPa) [33] and hypertensive (21.0 kPa) [38] states of each AA–GC–EA subsystem in this branch situation, respectively.

3.2.1. Normotensive State

As shown in Figure 5, the effect of the myogenic response was obvious when $P_{\mathrm{inlet}}$ reached 13.0 kPa, as the excessive pressure caused by the vascular structure was regulated (Figure 6(a-1,a-2)) by way of constriction of the AA diameter (Figure 6(c-1,c-2)), resulting in $P_{\mathrm{GC}}$ and $Q_{\mathrm{AA}}$ values of 6.7 ± 0.5 kPa and 246.5 ± 40.6 nL/min, respectively. These results are consistent with previous experimental data [4,5,6,33,34,35], indicating the validity of our modeling and simulation.

The values showed a decreasing trend corresponding with AA branching order, exhibiting larger values of both pressure and flow rate for the glomeruli connected to the shorter ILAs in the vascular network (Figure 6(a-2,c-2)). Although the simulated values of both pressure and flow rate were within the range of experimental measurements [4,5,6], the diameter contractions of the AAs were lower than experimental measurements (Figure 5b), and also showed periodic variation from 18.3% (N1) to 0% (N60).

3.2.2. Hypertensive State

As demonstrated in the present simulation (Figure 3 and Figure 5) and experiments [38], a high blood pressure inlet (i.e., 21.0 kPa) would exceed the limit of regulation and cause a rise in $P_{\mathrm{GC}}$. Figure 7 depicts the renal vascular pressure distribution under hypertensive conditions. The asymmetrical vasculature showed pressure differences distribution which was similar to those in the normotensive state.

It is worth mentioning that the degree of variation in value was reduced with the myogenic response compared to without it, regarding both pressure and flow rate, as shown by the smaller error bars in Figure 8a,c. However, since all AAs reached the maximum limit of constriction in diameter (Figure 8(b-2)), a further increase in pressure could not be avoided, resulting in $P_{\mathrm{GC}}=$ 8.1 ± 0.3 kPa. Furthermore, the pressure distribution was also correlated with the vascular structure, where glomeruli connected to shorter ILAs (e.g., N1) tended to acquire higher pressure P.

3.3. Branching Number Sensitivity Analysis and Comparison

In addition to the asymmetrical branching style, variations in vessel diameter also induced flow rate variations within the same-generation vessels which has been demonstrated in our previous study [28]. By coupling with the myogenic response, more comprehensive flow data were available for a wider range of inlet conditions. The relative standard deviations (RSD) of the flow rate and pressure are shown in Figure 9, indicating that, first, for all vascular networks, the degree of variation was greatest at 8 kPa when the myogenic response was not initiated, which demonstrated a significant effect of the myogenic response in regulation. Second, the RSD values of the symmetrical model were the least for all input conditions, which suggested its limited ability to generate variations. Finally, the RSD of the asymmetrical model increased with increased ILA branching number.

4. Discussion and Conclusions

We conducted a numerical simulation to investigate the blood flow behavior (i.e., flow rate and pressure) coupling of a renal vascular network including the myogenic response to various conditions. Individual glomerular pressure is regulated via the myogenic response, thereby regulating blood flow throughout the vascular network. In the vascular segment, the inlet pressure $P_{\mathrm{inlet}}$ and the vascular structure acted together on the myogenic response of each individual AA–GC–EA subsystem, such that the early-branching G1 in the vascular network, reached a well-regulated state first, with an interval of the inlet $P_{\mathrm{inlet}}=$ 10.5–21.0 kPa, whereas the late-branching G3 reached such a state with $P_{\mathrm{inlet}}=$ 13.0–24.0 kPa. In the entire vascular network, in contrast to the $P_{\mathrm{inlet}}$ interval (13.0–20.0 kPa) of the unified well-regulated state for all AA–GC–EA subsystems in the symmetric model, the asymmetric model exhibited differences among subsystems with $P_{\mathrm{inlet}}$ ranging from 12.0–17.0 to 16.0–20.0 kPa, eventually achieving a well-regulated state of 13.0–18.5 kPa for the entire kidney. Moreover, when $P_{\mathrm{inlet}}$ continued to rise (e.g., up to 21.0 kPa) beyond the vasoconstriction range of the myogenic response, high glomerular pressure was also related to vascular structure, where the $P_{\mathrm{GC}}$ of the early-branching N1 was 9.0 kPa and the late-branching N60 was 7.5 kPa.

When compared to measured data, the myogenic response model accurately replicated the experimentally reported pattern of overall blood flow with increasing inlet pressure (Figure 5a) [36,37]. In addition, degree of contraction of AAs approximated the experimental data [39]. Although the values varied, glomerular blood pressure was also regulated within a certain pressure range (Figure 6a and Figure 8a). These replications of the experimental data indicated the present study′s validity. However, there was a limit to the myogenic response′s ability to prevent high glomerular pressure by constricting blood vessels. Under the condition of an excessive value of blood inlet pressure, such as 21.0 kPa, glomerular pressures rose [32,37] causing glomerular damage [10,11,38] due to the presence of an upper limit to vasoconstriction capacity [39]. This phenomenon is illustrated in this study (Figure 9) and in the experimental comparison [11] of normotensive and hypertensive states.

The asymmetric structure of small vessels has been mentioned in previous renal studies [25,41], but either a regulation mechanism such as coupling with the myogenic response was not examined [28,41], or its effect on flow was not the focus of the investigation [24,25]. Due to technical limitations, experimental measurements may not yield sufficient data to completely describe the relationship between the vascular structures and flow. [25]. Furthermore, experiments deep inside the kidney examining a complex vascular structure and flow may affect the proper working condition of the blood vessels and therefore the accuracy of the data, as preparation procedures to access the interior of kidneys and AAs in isolation may alter the AAs′ behavior [31].

Experimental results suggest that there are differences in glomerular blood flow [3,4,5] and glomerular damage risks [37,39] between different regions. This would not be reproduced by a symmetrical model with diameter variations, as depicted in Figure 8, indicating that such differences could be caused by the asymmetrical branching style [24,28]. At the same time, glomerular pressure remains to be consistent. Therefore, it is a problem in modeling how to ensure stable a glomerular pressure while representing the measured variations in flow rate. The coupling of the myogenic response and renal vascular structure in this study may explain this issue. Simultaneously, as the inlet pressure $P_{\mathrm{inlet}}$ rose beyond the regulatory range of the myogenic response, high glomerular pressure was significantly correlated with its branching conditions, as glomeruli connected to a shorter ILA possessed higher blood pressure, resulting in a higher dysfunction risk. This also represents the different risks associated with high glomerular pressure in various regions [10,37,39].

Beyond the shortcomings of a numerical study of a detailed and highly complex biological system, three limitations of this study are especially apparent. First, the two vascular networks in this study represented renal vascular networks. However, the intricacies of how an arcuate artery branches into AA are still not clear based on observations so far [19,23]. Although some studies have provided important information on renal vascular structure [19,23,28], a more complete vascular geometry from the arcuate artery to the AAs should be assessed for numerical modeling in future works. Second, we omitted vasoconstriction in addition to AAs. Compared to previous autoregulation studies, we used a more complex asymmetric model [20,27], whereby the vascular structure could not be simplified into a dichotomous structure, as shown in Figure 2, and more branching nodes of the ILA resulted in a lack of uniform criteria for controlling the contraction of the superior vessels of the AAs. This may lead to a narrowed range of regulation outcome. Accordingly, the present model exceeded the upper $P_{\mathrm{inlet}}$ regulatory limit of around 18.5 kPa, as shown in Figure 4. However, this prediction curve is consistent with the experimental measurements [36,37], indicating the significant role of AA constriction in the myogenic response. In addition, we used a steady-state model of mean pressure input rather than a pulse wave to collect a large amount of result data while retaining a low computational load. In this way, the mean pressure value of each vascular network node can be derived under various pressure input settings, as shown in Figure 4 and Figure 7. However, future studies with a more realistic dynamic flow model including pulsatile pressure inputs are warranted, as the temporal diversity of pressures may further increase with periodically fluctuating pulse inputs [17,20]. Furthermore, under the same mean pressure input conditions, variations in pulsation amplitude can also induce differences in the myogenic response [10,39].

Despite the above limitations, our model investigated the myogenic response coupled with the entire complex renal vascular structure, which would be difficult to perform in experiments. These results demonstrate how the myogenic response regulates renal blood flow in vascular network system that comprises a large number of vessel elements. The findings might contribute a new idea to numerical investigations and provide the measurement requirements for research on renal vascular modeling and flow research.