# Morphometric Reconstruction of Coronary Vasculature Incorporating Uniformity of Flow Dispersion

^{1}California Medical Innovations Institute Inc., San Diego, CA, United States^{2}Faculty of Biomedical Engineering, Technion, Haifa, Israel

Experimental limitations in measurements of coronary flow in the beating heart have led to the development of *in silico* models of reconstructed coronary trees. Previous coronary reconstructions relied primarily on anatomical data, including statistical morphometry (e.g., diameters, length, connectivity, longitudinal position). Such reconstructions are non-unique, however, often leading to unrealistic predicted flow features. Thus, it is necessary to impose physiological flow constraints to ensure realistic tree reconstruction. Since a vessel flow depends on its diameter to fourth power, diameters are the logical candidates to guide vascular reconstructions to achieve realistic flows. Here, a diameter assignment method was developed where each vessel diameter was determined depending on its downstream tree size, aimed to reduce flow dispersion to within measured range. Since the coronary micro-vessels are responsible for a major portion of the flow resistance, the auto regulated coronary flow was analyzed in a morphometry-based reconstructed 400 vessel arterial microvascular sub-tree spanning vessel orders 1–6. Diameters in this subtree were re-assigned based on the flow criteria. The results revealed that diameter re-assignment, while adhering to measured morphometry, significantly reduced the flow dispersion to realistic levels while adhering to measured morphometry. The resulting network flow has longitudinal pressure distribution, flow fractal nature, and near-neighboring flow autocorrelation, which agree with measured coronary flow characteristics. Collectively, these results suggest that a realistic coronary tree reconstruction should impose not only morphometric data but also flow considerations. The work is of broad significance in providing a novel computational framework in the field of coronary microcirculation. It is essential for the study of coronary circulation by model simulation, based on a realistic network structure.

## Introduction

The coronary network supplies oxygen-carrying blood to each myocyte in the myocardium to meet its metabolic demand. This task is assigned to the coronary network consisting of millions of vessels of different diameters which are embedded within the myocardium. Flow in the network is determined by three major factors: (1) The network tree-like structure of asymmetrically bifurcating arteries and respective collecting venous one; (2) The vessels' compliance which responds to external loading by the periodically contracting myocardium to affect the coronary lumen area and the blood flow; (3) The three major flow control mechanisms which act in concert to actively regulate the coronary arterial diameters, aimed to meet the myocardium metabolic demand (Feigl, 1983). The flow control mechanisms are myogenic one which responds to the local trans-vascular pressure (Johnson, 1980; Kuo et al., 1988, 1990a; Davis, 2012); shear (flow) mechanism which responds to the local wall shear of the flowing blood (Kuo et al., 1990b, 1995; Jones et al., 1995b); and the metabolic mechanism which responds to imbalance in oxygen supply/demand in the myocardium (Feigl, 1983; Kanatsuka et al., 1989).

Despite considerable progress in recent decades in understanding coronary flow, further advances are impeded primarily due to considerable difficulties in experimental measurement of flow features in the heart which are due to the significant movement during the cardiac cycle, the inaccessibility of deep myocardial layer vessels, the dynamic myocardium-vessel interaction (MVI), and the difficulties in controlling the levels of each flow regulation mechanism under *in vivo* conditions. Mathematical modeling can be a powerful tool to supplement experiments. Previous modeling studies considered either the network structure (Karch et al., 1999; Beard and Bassingthwaighte, 2000) or the vessel interaction with the contracting myocardium (Downey and Kirk, 1975; Spaan et al., 1981; Krams et al., 1989; Rabbany et al., 1989; Kresh et al., 1990; Manor et al., 1994; Zinemanas et al., 1994; Vis et al., 1997; Smith, 2004; Westerhof et al., 2006; Algranati et al., 2010), or only the flow regulation (Liao and Kuo, 1997; Cornelissen et al., 2002; Carlson et al., 2008).

On the network structure, the very large number of coronary vessels precludes the possibility of deterministic reconstruction of the network. Statistical data-base is thus the method of choice. Such morphometric data was collected for the pig coronary tree (Kassab et al., 1993) and served as a basis for a number of flow studies in reconstructed coronary trees (Beard and Bassingthwaighte, 2000; Smith, 2004; Kaimovitz et al., 2005, 2010). Reconstruction based on statistical morphometric data is, however, non-unique. It may, in addition, lead to unrealistic flow features. Our preliminary flow analysis in morphometry-based reconstructed trees revealed that the predicted flow did not match measured perfusion dispersion in the terminal arterioles, being significantly higher compared to data measured by radioactive microspheres (King et al., 1985; Austin et al., 1994; Mori et al., 1995) and molecular microspheres (Stapleton et al., 1995; Matsumoto and Kajiya, 2001). In the present study, it was hypothesized that additional flow related constraints should be applied in order to produce realistic simulated coronary trees. The need to add flow constraints arises since the statistical morphometric data-base (Kassab et al., 1993) contains vessel order connectivity and statistics of order dependent diameter and length, but no data on diameter and length connectivity (correlation). Since flow is predominantly affected by the vessels' diameters, being dependent on their fourth power (Poiseuille equation), diameter re-assignment was used for flow optimization.

In the present study, a flow-based diameter re-assignment method was developed. In addition to the measured network morphometric data, the model considers the vessel interaction with surrounding myocardium and active diameter regulation. The procedure developed here was previously applied in the analysis of coronary flow regulation (Namani et al., 2018). It was used without details of rationale, the method procedure, no comparison with other studies or proof of validity. Here, we present analysis of the methodology and rationale, proof of validity compared with measured flow characteristics, and supplemented by detailed sensitivity analysis to the model parameters.

In view of the excessive computational cost of simulating dynamic regulated coronary flow in the entire arterial tree (few million vessels), and given the dominance of the microvasculature in flow resistance and control (Feigl, 1983; Jones et al., 1995a; Zamir et al., 2015), the present analysis focused on a microvascular sub-tree spanning vessel orders 1–6. This network was pruned from our previously reconstructed whole arterial LCX tree (Kaimovitz et al., 2005). To test the study hypothesis, the new diameter re-assignment was implemented, and the resulting flow was analyzed and compared with measured features of coronary flow.

## Methods

The tree reconstruction procedure consisted of three stages. First, an arterial microvascular subtree was pruned from the reconstructed LCX tree (Kaimovitz et al., 2005). Next, the subtree auto regulated flow was numerically analyzed. Finally, based on the predicted terminal order 1 arteriole flows, vessel diameters were iteratively re-assigned. The aim was to achieve measured physiological level of flow dispersion which is defined as the variance around the mean level of the flows at pre-capillary vessels. The flow dispersion is quantified by the flow coefficient of variation (CV = Mean/SD), where Mean is the mean flow across all pre-capillary arterioles, and SD is the respective standard deviation.

### Reconstruction of the Coronary Tree

The measured characteristics of the pig coronary tree are connectivity of vessel orders, vessel length and diameter, and myocardial depth. The LCX arterial tree was stochastically reconstructed (Kaimovitz et al., 2005) based on such morphometric data (Kassab et al., 1993). A 400-vessel arterial microvascular subtree was pruned from the reconstructed LCX tree. The stem vessel is of order 6 and the tree contains all vessels down to the pre-capillary (order 1) arterioles. Based on the morphometric data (Kassab et al., 1993), the subtree has both bifurcations and trifurcations. The conductivity matrix (Table 1), an essential input to the flow analysis, was derived from the network structure as elaborated in the Appendix.

### Flow in a Single Coronary Vessel

In the coronary microvasculature, flow is determined predominantly by viscosity forces under given driving pressures, while the effect of inertia is negligible. Under this condition, conservation of momentum in a single cylindrical vessel (Figure 1A) yields the Poiseuille relation:

where *q*(*t*) is the flow rate, *r*(*t*) is the vessel dynamic radius, Δ*P*_{L}(*t*) is the longitudinal pressure drop, μ(*r*) is the dynamic viscosity which is a function of the vessel radius, and *L* is the vessel length. An equivalent validated non-linear 3-element Windkessel model (Jacobs et al., 2008) was used here for each vessel in lieu of Equation 1 (Figure 1B). The network model was assembled from such vessel models.

**Figure 1. (A)** Schematic of a single coronary vessel and the applied dynamic pressure boundary conditions. *P*^{in} is the inlet pressure, *P*^{out} is the outlet pressure, *P*^{T} is the tissue pressure exerted by the myocardium, *P*^{mid} is the unknown intravascular pressure. **(B)** The vessel 3-element Windkessel model comprising two resistors R_{1}, R_{2} and a capacitive element C that accounts for the vessel volume changes.

### Boundary Conditions

Apart from the inlet and outlet pressures [*P*^{in}*(t)* and *P*^{out}*(t)*], the boundary conditions affecting the coronary flow are the myocardial vessel interaction (MVI) which is the myocardial loading on each vessel during the heartbeat, *P*^{T}(*t*), consisting of the sum of three contributions. The first two are the intra-myocardial pressure (IMP) and the shortening-induced intra-myocyte pressure, *P*^{SIP}*(t)* (Algranati et al., 2010). The third myocardial loading reflects the vessel radial traction by its tethering to the surrounding myocardium (Borg and Caulfield, 1979; Caulfield and Borg, 1979). IMP was taken to vary linearly with the myocardial depth (MRD) from zero at the epicardium to the dynamic LV cavity pressure *P*^{LV}*(t)* in the endocardium.

The inlet pressure boundary condition *P*^{in}*(t)*, adopted from measured pressures in order 6 arterioles (Chilian, 1991), was 87/55 mmHg in systole/diastole (average of 66 mmHg). The outlet pressure, *P*^{out}*(t)* at a terminal order 1 arteriole was taken to depend on each vessel's myocardial relative depth (MRD) as predicted in a previous study under non-regulated conditions (Algranati et al., 2010). It was thus 82/−32 mmHg (average = 22 mmHg) in the endocardium and 34/0 mmHg (average = 13 mmHg) in the epicardium for systole/diastole. LV pressure signal, *P*^{LV}(t), was taken to be 99/7 mmHg (average = 36 mmHg), as predicted from a distributive LV mechanical model (Nevo and Lanir, 1989) under a resting heart rate of 75 BPM. More details on the boundary conditions are presented in the Appendix.

### Diameter Re-assignment

The basis for the diameter reassignment is the notion that given the vessel length, and assuming uniform flow in all terminal order 1 vessels, diameter of an upstream vessel should be compatible with a flow rate which is proportional to the number of terminal arterioles it feeds. Hence, from Poiseuille's formula (Equation 1), the vessel diameter *d* should follow as:

where *n*_{pc} is the number of terminal order 1 arterioles in the crown of that vessel, *L* is the vessel length, and the pressure drop along the vessel length (Δ*P*_{L}) is determined from the network flow analysis. The constant parameter *B* represents the target terminal arteriole flow *q*_{target} needed to balance the myocardium metabolic demand. Equations 1 and 2 yield for *B*

The target flow rate for the terminal arterioles *q*_{target} was set to be 0.4 × 10^{−3} mm^{3}/s, within the measured range of resting diastolic flow which is 0.4–0.6 × 10^{−3} mm^{3}/s for vessels of order 1 in dog hearts (Tillmanns et al., 1974; Ashikawa et al., 1986; Stepp et al., 1999). The dog heart level was adopted due to the absence of such data in the swine, and based on the remarkably similar relationship between left ventricle myocardial blood flow (LVMBF) and heart rate in dog and swine hearts (Figure 4 in Duncker and Bache, 2008). Diameter re-assignment was carried out under flow control by myogenic and shear regulations. Metabolic regulation was not included here based on the rationale that the low resting target flow rate represents a state in which the metabolic regulation is nearly inactive. The metabolic regulation plays a crucial role under higher metabolic demands. Diameter reassignment was an iterative procedure since the pressure drop along each vessel, Δ*P*_{L}, changes each time diameters are modified. The re-assignment was subject to the constraints that the diameter level remains within the measured range of that vessel order (Kassab et al., 1993). If the diameter fell outside that range, the diameter was set to be at the limit of its order range. This was done to ensure the reconstructed network statistics to be compatible to that of the database (Kassab et al., 1993). In addition, re-assignment was restricted to comply with the constraint that a vessel diameter be smaller than its mother vessel (hydraulic continuity). If this condition was violated, diameter was assigned to be similar (difference was 0.1 μm) to that of the mother vessel. For order 0 vessels (capillaries), no correction was applied. Convergence of this iterative method was determined when the coefficient of variation in terminal arteriolar flows reduced to < 25%, i.e., within the range of measured data (Austin et al., 1994; Matsumoto and Kajiya, 2001).

#### Network Flow Solution and the Diameters Re-assignment

For a given state of the network vessels unloaded diameters, flow was iteratively solved under the dynamic pressure boundary conditions by imposing mass conservation in all the tree nodes at each time along the cardiac cycle. This yields a system of ordinary differential equations (ODEs). The method of solution is presented in detail in the Appendix. Briefly, for the given vessel unloaded radius, the loaded radius was taken to vary with the pressure difference across the vessel wall, which depends on the vessel's passive and active pressure-diameter relationships (PDR). The passive properties are the non-linear pressure-diameter relationships (Young et al., 2012) constrained by the vessel tethering to its surrounding myocardium. The active properties reflect the contractile and dilatory effects of the vessel's smooth muscle by myogenic (pressure-driven) and shear regulation mechanisms. Their contributions were taken from *in vitro* data on single isolated coronary micro-vessels (Liao and Kuo, 1997). Details are provided in the Appendix.

Once the network flow was solved for a given diameter state, each vessel unloaded diameter was iteratively re-assigned (by application of Equations 2, 3), and the network flow was solved following each diameter re-assignment until convergence was attained.

### Model Predictions

In presenting the model results, emphasis will be on those predictions which can be compared with observed characteristics of the coronary network flow.

#### Longitudinal Pressure Distribution

Numerical solution of the network flow yields the intravascular pressures at each node and vessel in the tree. The associated pressure distribution can be contrasted with the few *in vivo* data that exist on the intravascular pressure in the coronary micro-vessels (Tillmanns et al., 1981; Chilian, 1991). In view of the paucity of coronary data, the predicted pressure distribution was also compared with parallel data in the cat mesentery (Fronek and Zweifach, 1975; Zweifach and Lipowsky, 1977), and in the hamster cheek pouch (Davis et al., 1986).

#### Spatial Flow Heterogeneity

As another test of validity of the diameter re-assignment methodology presented here, the model predicted flow dispersion in terminal arterioles was compared with data. This validation has three related but different aspects (Bassingthwaighte and Beyer, 1991; Yipintsoi et al., 2016). The flow dispersion level (CV), the flow fractal nature (self-similarity), and the flow spatial auto-correlation.

##### Flow dispersion as a function of sample mass

Experimental evidence showed that flow dispersion decreases with increasing sample size (Bassingthwaighte et al., 1989; Matsumoto and Kajiya, 2001). In addition, the dispersion was found to have a fractal (self-similar) nature (Bassingthwaighte et al., 1989). When plotted on a log-log scale, self-similarity is seen as a linear dispersion vs. sample size relationship where the slope is 1–*D, D* being the fractal dimension (Bassingthwaighte et al., 1989). A value of *D* = 1.0 indicates highly correlated near neighbor flows and *D* = 1.5 indicates that neighboring flows are uncorrelated having spatial randomness. These observed characteristics of the coronary flow were examined in the predicted flow as another test of validity of the proposed diameter re-assignment. The tissue volume fed by a terminal arteriole was previously estimated to be 0.2–1.0 mm^{3} (Bassingthwaighte et al., 1974; Matsumoto and Kajiya, 2001). Since there are interspecies variations, the sample mass fed by a terminal arteriole was estimated in an independent manner. Based on the morphometric data of Kassab et al. (1993), approximately 10^{6} order 1 arterioles feed the 150 *gm* heart, which implies that each arteriole supplies approximately 0.15 *mg*. From this and the above listed estimates, the mass perfused by a terminal arteriole was set to be 0.2 *mg* in the present study. To estimate the dispersion vs. sample mass relationship, sample size was increased by grouping together tissue units supplied each by a single terminal arteriole. Sample size was increased up to 4.4 *mg*. The flow dispersion was represented by the coefficient of variation (CV = SD/Mean), SD being the standard deviation.

##### Spatial autocorrelation of pre-capillary flows

Flow spatial autocorrelation is a measure of the flow continuity. The idea is to estimate the correlation of near-neighboring flows between regions as a function of the distance between them. The spatial location of vessels in 3-D space given by its x, y, z coordinates were used to calculate the Euclidean distance, *e* $=\sqrt{{({x}_{1}-{x}_{2})}^{2}+{({y}_{1}-{y}_{2})}^{2}+{({z}_{1}-{z}_{2})}^{2}}$, between any pair of terminal arterioles with spatial co-ordinates, (*x*_{1}*, y*_{1}*, z*_{1}) and (*x*_{2}*, y*_{2}*, z*_{2}). Per definition, the spatial autocorrelation function, *f(e)*, is the sum of co-variances in flow over all pairs of terminal arterioles separated by a distance *e*, divided by the flow variance over all terminal arterioles of the subtree. It is thus expressed as

where *n*_{e} are the number of regions that are separated by a distance *e*. The regional flow was normalized by the total mean flow to derive relative flow *Q*. The regional relative flows are *Q*_{i}(*e*) and *Q*_{j}(*e*) in the *i**th* and *j**th* regions which are separated by a distance *e*, and *Q*_{term} is the flow in the pre-capillaries. Two aspects of the flow autocorrelation have been measured and modeled in previous studies. Accordingly, both the model predicted flow autocorrelation, and the predicted nearest-neighbor autocorrelation coefficient were estimated before and after diameter re-assignment and compared with data (Matsumoto and Kajiya, 2001) and with a previous model (Beard and Bassingthwaighte, 2000).

##### Fractal structure of pre-capillary flows

The predicted autocorrelation function (Equation 4) was compared with a theoretical prediction of this function for terminal arteriole flows in a self-similar fractal tree structure (Beard and Bassingthwaighte, 2000). The spatial autocorrelation function of a fractal network with fractal dimension *D* as a function of the spatial distance *e* is given by

The spatial autocorrelation function *f(e)* that was determined from Equation 4 was plugged into Equation 5. Linear regression was carried out between *f(e)* and *e* and the exponent 4–2*D* was estimated, and hence *D*. This estimate of the fractal dimension *D* was compared with the *D* estimated from the dispersion-sample size analysis described above.

### Sensitivity Analysis

The model parameters used for diameter re-assignment are those of the vessels' myogenic and flow regulations ρ_{m}, ϕ_{m}, *C*_{m}, *F*_{τmax} (Appendix), the target terminal flow *q*_{target} and the network outlet pressure *P*^{out}(Figure 1). Their values were adopted from published experimental studies (Liao and Kuo, 1997) and from simulation of non-active full sized network (Algranati et al., 2010). Since biological systems exhibit inherent variability, a sensitivity analysis was carried out by perturbing each parameter uniformly in all the network vessels. Each parameter was varied from the reference level. The perturbed levels of *q*_{target} were 0.2 and 0.6 × 10^{−3} mm^{3}/s (the reference being 0.4 × 10^{−3} mm^{3}/s). The outlet pressures *P*^{out} were perturbed by ±20 mmHg from their depth dependent reference levels. The activation parameters ρ_{m}, ϕ_{m}, *C*_{m}, *F*_{τmax} were perturbed by ±25% from their reference levels. Upon each perturbation, diameter re-assignment and flow analysis were redone with the new set of parameters.

The sensitivity of the model predictions were assessed in relation to vessel diameters, terminal flow level and its dispersion (CV), flow fractal dimension estimated from both the dispersion and spatial autocorrelation of terminal flows, and network longitudinal pressure distribution.

## Results

### Flow Effects of Diameter Re-assignment

Diameter re-assignment significantly changed the frequency distribution of the terminal arteriole flow, *q*_{term} to be narrower (more uniform) with a higher peak close to the target flow of 0.4 × 10^{−3} mm^{3}/s (Figure 2). The range of flows in the unmodified reconstructed tree was 0.04–2.57 × 10^{−3} mm^{3}/s which changed after diameter re-assignment to be considerably narrower in the range of 0.14–0.73 × 10^{−3} mm^{3}/s. In addition, in the native reconstructed network (prior to diameter re-assignment) the mean terminal flow was 0.57 × 10^{−3} mm^{3}/s, while in the re-assigned network it reduced to 0.39 × 10^{−3} mm^{3}/s which is within the physiological range. Terminal flow dispersion (CV) in the native reconstructed network was 0.97, and reduced in the diameter re-assigned network to be 0.22, which is within the physiological measured range (Table 2).

**Figure 2**. Comparison of the flow distribution in the terminal pre-capillary arterioles in the 400-vessel subtree under shear and myogenic controls prior vs. after diameter re-assignment. Each bin on the X-axis represents a 0.05 × 10^{−3} mm^{3}/s range of terminal pre-capillary flows per cardiac cycle. The Y-axis is the pre-capillary frequency for each bin.

**Table 2**. Flow in the pre-capillary vessels under a target terminal flow of 0.4 × 10^{−3} mm^{3}/s, before and after diameter re-assignment, in the passive state and under myogenic and shear flow regulations.

### Longitudinal Pressure Distribution

The predicted pressure distribution (Figure 3) decreased from 62 mmHg in vessel order 6 to 27 mmHg in vessel order 1. The curve shows elevated pressure drop over the low vessel orders, with the highest drop in vessel orders 1–4. The pressures in these small vessels have the highest standard deviations. The current model predictions are compared with measured data. The comparison shows very good agreement with Chilian (1991) pig coronary data for order 6 vessels, and likewise with Tillmanns et al. (1981) rat and cat coronary data for a range of vessel sizes in the intermediate (orders 3,4,5) and order 6 vessels, as well as with Fronek and Zweifach (1975) cat mesentery data for intermediate vessel range. The cheek pouch data (Davis et al., 1986) differ substantially from the model predictions.

**Figure 3**. Comparison of the model-predicted time-averaged intravascular pressure distribution in the subtree after diameter re-assignment with the data sets of Chilian (1991) for pigs coronaries, of Tillmanns et al. (1981) for the epi-myocardium of rats and cats, of Fronek and Zweifach (1975) and Zweifach and Lipowsky (1977) for the cat mesentery, and of Davis et al. (1986), for the hamster cheek pouch. Model predictions and mesentery data show significant pressure drop in small micro-vessels.

### Spatial Flow Heterogeneity

#### Flow Dispersion as a Function of Sample Mass

The flow dispersion (CV) in myocardial regions was found to decrease exponentially (linearly in a log-log plot) with increase in sample tissue mass (Figure 4). In the native reconstructed tree prior to diameter re-assignment, the highest CV of 0.97 was found at the lowest sample mass of 0.2 *mg*, corresponding to the perfused mass by one terminal arteriole in the present study. The lowest CV of 0.16 was found at a sampled tissue mass of 4.4 *mg*. These values are one half to one order of magnitude higher than in the measured data of the rabbit ventricular free walls (Matsumoto et al., 1999). After diameter re-assignment, the CV reduced to 0.24 for a sample mass of 0.2 *mg*, which decreases linearly in the log-log plot to 0.06 for a sample mass of 4.4 *mg*. The re-assigned network flow dispersions vs. sample size relationship is in close agreement (Figure 4) with the data of Matsumoto et al. from rabbit's free ventricular wall (Matsumoto et al., 1999; Matsumoto and Kajiya, 2001). The regression of log(CV) vs. log (sample mass) is linear with a slope of −0.62 in the original reconstructed tree which changed to −0.44 after diameter re-assignment, in close agreement with the experimental data. These slope values yield fractal dimension *D* of 1.62 in the native reconstructed network and changed to 1.44 after diameter re-assignment.

**Figure 4**. Log-log plot of model prediction and data of flow dispersion (represented by CV = SD/Mean) vs. perfused sample mass. The slope of the linear fit equals 1–D, where D being the fractal dimension (Bassingthwaighte et al., 1989). The results yield *D* = 1.44 for the diameter-reassigned tree, *D* = 1.62 for the original reconstructed tree, *D* = 1.24 for the data of Matsumoto et al. (1999) in the sub-endocardial region. There is a close agreement between the data of Matsumoto et al. and the present model prediction for the diameter re-assigned tree, while predictions of the original reconstructed tree yields a CV of up to an order of magnitude higher than the data.

#### Spatial Autocorrelation of Pre-capillary Flows

Diameter re-assignment was found to significantly increase the spatial flow auto-correlation *f*(*e*) (Equation 4) at distances above ca. 0.4 mm (Figure 5) from an asymptotic average of 0.02 ± 0.01 for an unmodified tree to a level of 0.08 ± 0.01 after re-assignment. The function *f*(*e*) has an exponential decay with increasing spatial distance between terminal arterioles (Figure 5). A single term exponential of the form *a* · *exp*(−*e*/*e*_{0}) + *c* was fit to *f*(*e*) for the subtrees prior and after re-assignment (Figure 5). The fit yielded parameter estimates (Table 3) showing that the most significant change induced by diameter re-assignment is a 300% increase of the long-range autocorrelation (from 0.02 to 0.08).

**Figure 5**. Flow spatial autocorrelation function between pairs of pre-capillary arterioles as a function of their spatial distance (Equation 4). Symbols are model estimates. The curves are single term exponentials of the form *a* · *exp* (−*e*/*e*_{0}) + *c* fit to *f*(*e*) of the two subtrees. The parameter estimates are listed in Table 3.

**Table 3**. Parameter of the spatial autocorrelation model $\widehat{f}\left(e\right)=a\xb7exp\left(-e/{e}_{0}\right)+c\text{}$ estimated to fit the calculated spatial autocorrelation *f*(*e*) for the original reconstructed and the diameter re-assigned subtrees.

#### Fractal Structure of Pre-capillary Flows

Estimates of the fractal dimension, *D* obtained from the flow dispersion analysis and from its auto-correlation are shown in Table 4. The difference in the estimates of *D* between the two methods is 4% for the unmodified reconstructed tree and 9% in the diameter reassigned tree. Both dispersion and spatial autocorrelation analyses gave a *D* > 1.5 in the original reconstructed sub-tree which means that flows are negatively correlated. In the diameter reassigned tree, *D* = 1.44 from the dispersion analysis and *D* = 1.31 from the spatial autocorrelation. Thus, the flow fractal nature of the network is maintained after diameter reassignment.

**Table 4**. Fractal dimension, D predicted from flow dispersion and from spatial autocorrelation analyses in a 400 vessel reconstructed coronary tree before and after diameter reassignment.

### Sensitivity Analysis

Two features are noticeable from the results of the diameter sensitivity analysis (Table 6). First, the diameter changes resulting from perturbations of the model parameters are small, <5%. The exceptions are the substantial diameter changes in order 6 vessels under changes in the following: (1) Outlet pressure *P*^{out}, (2) increasing terminal target flow *q*_{target}, (3) decreased myogenic parameter ρ_{m}, and (4) under changes in the myogenic parameter *C*_{m}. Likewise, the diameter changes are high in orders 5 and 6 vessels under reduction of the myogenic parameter ϕ_{m}. The changes under decreased terminal target flow *q*_{target} are also > 5% even when consideration is given to the higher perturbation set for that flow (± 50% of the reference level) compared to the ±25% perturbations of the other parameters. Second, the trend of diameter change under a specific perturbation is different for the different vessel orders; i.e., diameter in some orders may increase while in others they decrease. The exceptions are the diameter responses to changes in the terminal target flow rate *q*_{target} which are positive in all orders under increased *q*_{target} and negative when it decreases, and corresponding changes under reduced *P*^{out} which are negative for all vessel orders.

Parameter perturbations affect the terminal flow rate *q*_{term} and its dispersion CV (Table 7). The terminal target flow rate *q*_{target} has the highest effect on the terminal flow, followed by that of the myogenic parameter ρ_{m}. Any perturbation of *F*_{τmax} (i.e., both increase and decrease) changes the CV of the terminal flow in the same direction while for the other parameter opposite changes have opposite effects. The flow CV increases from its reference under all parameter perturbations except when *C*_{m} is reduced, where the highest effect on CV is that of the shear parameter *F*_{τmax}, followed by those of the myogenic parameter ρ_{m} and the target flow rate *q*_{target}.

Sensitivity analysis of the terminal flow fractal dimension *D* (Equation 5) to the model parameters (Table 8) shows that the two approaches of estimating the fractal dimension not only produce different values, but in addition predict at times the same trend of change in *D* while in other cases they predict a mutually opposite change. The dispersion analysis approach (left panel) predicts, in most cases, an increase in *D* values under perturbations, which in few cases crosses the boundary *D* = 1.50 between positive and negative correlation of near neighbor flows (Bassingthwaighte and Beyer, 1991). In contrast, the spatial autocorrelation estimated *D* (right panel in Table 8) is mostly lower than its dispersion estimated value except in two cases–reduced ρ_{m} and increased *F*_{τmax}. In the latter case, the estimated (*D* = 1.65) is well within the negative correlated near neighbor flows.

Perturbations in the model parameters seem to have a small influence on the longitudinal (downstream) pressure distribution (Figure 6). For the effects of the myogenic and shear parameters ρ_{m}, ϕ_{m}, *C*_{m} and of *F*_{τmax} (Figures 6A–D), these insensitivities may stem from the uniform percent perturbations of the respective parameters across all vessel orders. The outlet pressure *P*^{out} has a small effect on the pressures in the smallest vessels (order 1 and 2) which are close to the network outlet, but has insignificant influence on the higher order vessels (Figure 6E). The target terminal flow rate *q*_{target} has a very small effect on the pressure distribution (Figure 6F).

**Figure 6**. The network longitudinal (downstream) pressure distribution under ±25% perturbations of the myogenic and flow parameters ρ_{m} **(A)**, ϕ_{m} **(B)**, *C*_{m} **(C)**, and *F*_{τmax}**(D)**, under ± 20 mmHg variations of the network outlet pressure *P*^{out} **(E)**, and under ±50% variations of the terminal target flow rate *q*_{target} **(F)**. Depicted in each sub-figure are the pressure distribution in the reference case (solid line), and under positive (broken line) and negative (dash-dot line) perturbation.

## Discussion

This study addresses a longstanding difficulty in modeling coronary flow: how to reconstruct realistic coronary trees which mimic native trees in both structure and flow. Such reconstructions are indispensable for analyzing the coronary circulation. To this end, physiological constraints were proposed in conjunction with measured morphometric data. The present study used previous reconstructed coronary trees (Kaimovitz et al., 2005) which were built based on morphometric data (Kassab et al., 1993), and incorporated a new diameter re-assignment procedure designed to reduce flow dispersion to within measured ranges, without affecting the vessel orders or their lengths. The results show that application of this flow constraint provided reconstructed trees that conform with the morphometric data and at the same time, yield realistic flows which correspond with measured data in terms of the longitudinal pressure distribution, flow dispersion, spatial autocorrelation, and fractal structure (Austin et al., 1994; Matsumoto and Kajiya, 2001).

The current network reconstruction is the first which considers a comprehensive set of physiological aspects of the coronary circulation: measured morphometric-based tree-like structure of asymmetrically bifurcating vessels (Kassab et al., 1993); vessel compliance which allows for the external loading by the periodically contracting myocardium to affect the vessel lumen area and the blood flow; dynamic coronary flow driven by the likewise dynamic inlet and outlet pressures, myocardium-vessel interaction (Algranati et al., 2010); and the three major flow control mechanisms which act in concert to actively regulate the vessel diameters to meet the myocardium metabolic demand (Feigl, 1983).

A number of previous network flow studies considered morphometric data of vascular dimensions in symmetric vascular trees. Liao and Kuo (1997) applied a symmetric four element arterial network representing small arteries and large, intermediate, and small arterioles. Cornelissen et al. (2002) considered a network which consisted of 10 compartments in series where the first nine represented the arterial tree with resistance depending on the local pressure, and the distal one lumped the capillaries and the venules. Based on data of Van Bavel and Spaan (1992) and (Kassab et al., 1993, 1994; Kassab and Fung, 1994), Vis et al. (1997) applied a model network of symmetric bifurcating arterial and venous networks of 19 generations each, connected by a capillary network of parallel vessels. Arciero et al. (2008) and Carlson et al. (2008) network includes seven compartments connected in series, each consisting of identical parallel segments subjected to the same hemodynamic and metabolic conditions. They represent the upstream artery, the large arterioles, small arterioles, capillaries, small venules, large venules, and veins. The models of Schreiner's group (Schreiner and Buxbaum, 1993; Karch et al., 1999, 2003) generated 2D and 3D coronary networks by the method of constrained constructive optimization by successively adding new terminal segments while maintaining a set of physiological constraints. The actual network structure of asymmetrical bifurcating tree was considered by very few studies. Smith (2004) analyzed the flow in a discrete anatomically accurate model of the largest six generations of the coronary arterial tree. Beard and Bassingthwaighte (2000) network was reconstructed based on Kassab et al. (1993) morphometric data supplemented by self-similarity and avoidance algorithm. Kaimovitz et al. (2005, 2010) developed a large scale 3D reconstruction of the entire porcine coronary vasculature (arterial, capillaries, and veins) based on the morphometric data of Kassab et al. (1993, 1994) and subject to both global constraints relating to the location of the larger veins, and to local constraints of measured morphological features. No previous reconstruction considered flow constraints associated with the coronary flow dynamics, with the associated myocardium-vessel interaction (MVI), and included the active flow regulation.

As expected, flow dispersion in the myocardial regions was modified by diameter re-assignment (Table 2, Figure 4). A novel finding is that diameter re-assignment based on perfusion homogeneity consideration has dramatic effect on the dispersion and fractal nature of the flow, endowing it with dispersion/sample mass relationship and fractal nature that are similar to those of native *in vivo* coronary networks. These results suggest that the flow fractal nature is a consequence of not only the network asymmetry (Bassingthwaighte et al., 1989), but is also determined by the vessel diameters.

The results in Table 2 show that flow dispersion in the terminal pre-capillary arterioles is highly sensitive to diameter re-assignment. This sensitivity is attributed to the fact that by Poiseuille relationship (Equation 1), the rate of flow depends on the fourth power of diameter. Hence, diameter re-assignment is a powerful scheme for achieving target flow features. In addition, the developed scheme can be readily extended to incorporate other physiological features such as different profiles of non-uniform pre-capillary pressure distribution (Keelan et al., 2016) and a transmural non-uniform myocyte metabolic demand (Namani et al., 2018).

The pressure gradient along the subtree drives the flow. *In vivo* pressure data in the intermediate and small arterioles below order 6 that could serve for comparison with the model results are not available for the swine coronary micro-vessels. Data in these vessels are available for the cat mesenteric arterioles (Fronek and Zweifach, 1975; Zweifach and Lipowsky, 1977). Model comparison with these data shows close match with the intravascular pressures predicted by the model (Figure 3). The small differences between the slopes of the longitudinal pressure distribution curves between the model and mesenteric data may result from the difference in the vessel boundary conditions, in particular, the dynamic myocardial-vessel interaction (MVI) which is not present in the mesenteric microcirculation.

Flow dispersion in the terminal arterioles of the reconstructed and diameter re-assigned trees were found to increase with a decrease in sample size, similar to the trend observed experimentally by microsphere deposition (Bassingthwaighte et al., 1989). The log-log relation between flow dispersion and sample mass is linear (Figure 4) in both the model and measured data (Bassingthwaighte et al., 1989; Mori et al., 1995; Matsumoto and Kajiya, 2001). At sample sizes similar to the model ones, flow dispersion data are available (Matsumoto and Kajiya, 2001). There is close agreement between these data and the model predicted dispersions in the terminal arterioles. This stems from the realistic tree structure and from inclusion of the effects of myocardial vessel interaction and flow regulation in the present model, while not all of these factors were considered in previous models. This favorable comparison of the present model with data provides a rationale for using the diameter reassignment for reconstruction of coronary trees in future studies and offers flow dispersion as a new criterion for validation of these coronary trees.

The model predicted spatial auto-correlation of near neighboring flows, was found to decrease exponentially with the distance (Figure 5), a trend similar to that observed by microsphere experiments (Austin et al., 1994; Mori et al., 1995; Bassingthwaighte et al., 2012; Yipintsoi et al., 2016). Diameter re-assignment induced a significant increase of the long-range autocorrelation (Table 3) which suggests that distant flows are more correlated. This can be accounted for by recalling that re-assignment produces more uniform flows with a narrower dispersion range (Figure 2).

Our model estimates for nearest-neighbor autocorrelation coefficient of terminal arterioles flows are 0.47 and 0.52 (Table 5) for the original reconstructed and for the re-assigned sub-trees, respectively. These figures are comparable to the experimental value of 0.4 using a similar sample mass (0.2 mg) (Matsumoto and Kajiya, 2001), and close to some experimental estimates using higher sample masses (Mori et al., 1995; Yipintsoi et al., 2016).

**Table 5**. The nearest neighbor autocorrelation coefficients in the present model, and in several experimental data and other models.

Although, sample sizes used in previous studies were mostly different from the present model, the closeness of most reported coefficients supports the fractal nature of the network flow. In fractal networks the nearest-neighbor autocorrelation coefficient is independent of the sample size (Beard and Bassingthwaighte, 2000).

In models of the coronary tree, Beard and Bassingthwaighte (2000) found the nearest neighbor correlation coefficient to be in the range of 0.38–0.43 which is close to the range predicted by the present model.

The fractal dimension *D* is an established parameter for quantifying the global heterogeneity of coronary blood flow. It can be estimated either from the flow dispersion or from the spatial auto-correlation. Dispersion analysis estimated *D* to be 1.44 after diameter re-assignment. Spatial autocorrelation analysis gave an estimate for *D* of 1.31 which is closer to most of the measured data (Table 4). In contrast, prior to re-assignment, *D* was estimated to be 1.62 (> 1.5) which represents a tree with negatively correlated flows (*f*(*e*) < 0). This was rarely observed, and mainly in abnormal physiological conditions such as under reduced perfusion pressure (Mori et al., 1995) and stenotic conditions (Kleen et al., 1997). These results suggest that diameter re-assignment produces coronary flows that are fractal in nature and are comparable to *in vivo* observations; i.e., a fractal dimension of *D* < 1.5 with positive correlated near neighboring flows.

In comparison with previous models, the present model predicted *D* = 1.31 from spatial autocorrelation analysis of the diameter re-assigned tree is close to the model estimated values by Karch et al. (2003) (*D* = 1.38 – 1.45), by Van Bavel and Spaan (1992) (*D* = 1.23) and by Beard and Bassingthwaighte (2000) (*D* = 1.24) (Table 4).

Diameter sensitivity analysis of the model parameters (Table 6) reveals that the effects are mostly small. This counter-intuitive outcome is likely due to the large effect of the diameter on the vessel resistance to flow (Equation 1). The other finding (Table 6), is that the direction of diameter changes (increasing vs. decreasing) is mostly inconsistent, which seems to be also counter-intuitive. This is likely due to the diameter change by the re-assignment results from a complex non-linear interaction between the network asymmetric structure, the order-dependent passive and active mechanical properties, and the blood flow.

**Table 6**. Sensitivity of the network re-assigned vessel diameters in each order to perturbations in the model active parameters ρ_{m}, ϕ_{m}, *C*_{m}, *F*_{τ max} (±25% of their reference), network output pressure *P*_{out} (±20 mmHg of its references), and the target flow rate q_{target} (±50% of its reference 0.4 × 10^{−3} mm^{3}/s).

The results presented in Table 6 allow comparison of the effects of the myogenic and flow regulations on the vessel diameters. The results show that the effects of their perturbations (ρ_{m}, *F*_{max}, respectively) on the network vessel diameters are small and of approximately comparable levels. In addition, it is seen that the effects of these perturbations on order 1 vessel diameter is generally opposite to their effects on other vessel orders. This result may be accounted for if one recalls that at *in vivo* conditions, the highest effect of the metabolic regulation is on order 1 vessels since the metabolic conducted response is highest at the interface of capillaries with these vessels and drops exponentially toward higher order vessels (Namani et al., 2018).

An interesting conclusion from the results in Table 7 is that under virtually all perturbations, the flow dispersion (CV) is higher than the reference level (= 0.22). The basis and implications of this prediction requires further inquiry.

**Table 7**. Sensitivity of the terminal flow rate and its dispersion over the terminal order 1 vessels to perturbations in the model active parameters ρ_{m}, ϕ_{m}, *C*_{m}, *F*_{τmax} (±25% of their reference), network output pressure *P*^{out} (±20 mmHg of its reference), and the target flow rate q_{target} (±50% of its reference 0.4 × 10^{−3} mm^{3}/s).

It is not surprising that with the increase in the terminal flow dispersion (CV), the flow fractal dimension *D* changes as well (Table 8). Most affected is the *D* estimates from the dispersion vs. sample mass relationship, which in several cases changes the flow from being spatially positively correlated (*D* < 1.50) to being negatively correlated (*D* > 1.50).

**Table 8**. Sensitivity of the terminal flow fractal dimension D to perturbations in the model active parameters ρ_{m}, ϕ_{m}, *C*_{m}, *F*_{τmax} (±25% of their reference), network output pressure *P*_{out} (±20 mmHg of its reference), and the target flow rate q_{target} (±50% of its reference 0.4 × 10^{−3} mm^{3}/s).

## Limitations of the Study

A limitation of the study is that the vessel re-assigned diameters and resulting flow in them depend on the level of the target terminal flow (Equation 3). The applied target flow level was selected from measured data of resting metabolic demand. Another limitation is the size of the network analyzed, which although consisting of 400 vessels is relatively small compared with the full-sized network. Applying the current diameter re-assignment methodology to the entire coronary tree is needed since it will allow for more accurate pressure boundary conditions based on measured data. Otherwise, the network size should have no effect on the diameter re-assignment, nor on the flow in each vessel, provided the appropriate pressure boundary conditions of inlet, outlet and tissue pressure (MVI) are applied. The rationale is that the flow in any subtree chosen is determined by the vessels dimensions and by the inlet, outlet and tissue pressure, and is otherwise independent of flow events in other branches. Finally, consideration of additional physiological criteria to the flow dispersion, such as myocardial metabolic heterogeneity, may provide a more accurate diameter re-assignment. Such considerations can be readily incorporated to the methodology developed here.

## Conclusion

A diameter re-assignment method was developed and applied on morphometry-based reconstructed coronary networks, with the goal to obtain realistic coronary flow characteristics. This method focuses on the size of each vessel crown and the number of terminal pre-capillary vessels it supplies. The results show that whereas flow features in the original reconstructed subtree (without diameter re-assignment) significantly differ from measured data, application of diameter re-assignment was found to improve the physiological realism of the flow field in terms of correspondence with measured data of the longitudinal pressure distribution and of the terminal flow dispersion, spatial autocorrelation, and fractal structure. These results suggest that future coronary flow studies using reconstructed vascular networks should apply both morphometric data and flow considerations to reconstruct the network. This conclusion likely applies (with appropriate modifications) to circulation studies in skeletal muscle, and perhaps in other tissues as well.

## Author Contributions

RN, YL, and GK made significant contributions to creation and design of the study, interpretation of data, and writing the manuscript. RN and YL developed the algorithms. RN wrote the code, ran the simulations and collected the data. RN, YL, and GK approve the final copy of the manuscript, agree to be responsible for all aspects of the research in the manuscript and confirm the role of authorship.

## Funding

This research was funded by the US-Israel Binational Science Foundation (BSF grant 2009029) and National Institutes of Health (U01HL118738).

## Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphys.2018.01069/full#supplementary-material

## References

Algranati, D., Kassab, G. S., and Lanir, Y. (2010). Mechanisms of myocardium-coronary vessel interaction. *Am. J. Physiol. Heart Circ. Physiol.* 298, H861–873. doi: 10.1152/ajpheart.00925.2009

Arciero, J. C., Carlson, B. E., and Secomb, T. W. (2008). Theoretical model of metabolic blood flow regulation: roles of ATP release by red blood cells and conducted responses. *Am. J. Physiol. Heart. Circ. Physiol.* 295, H1562–H1571. doi: 10.1152/ajpheart.00261.2008

Ashikawa, K., Kanatsuka, H., Suzuki, T., and Takishima, T. (1986). Phasic blood flow velocity pattern in epimyocardial microvessels in the beating canine left ventricle. *Circ. Res.* 59, 704–711. doi: 10.1161/01.RES.59.6.704

Austin, R. E., Smedira, N. G., Squiers, T. M., and Hoffman, J. I. (1994). Influence of cardiac contraction and coronary vasomotor tone on regional myocardial blood flow. *Am. J. Physiol.* 266, H2542–H2553. doi: 10.1152/ajpheart.1994.266.6.H2542

Bassingthwaighte, J. B., Beard, D. A., Carlson, B. E., Dash, R. K., and Vinnakota, K. (2012). Modeling to link regional myocardial work, metabolism and blood flows. *Ann. Biomed. Eng.* 40, 2379–2398. doi: 10.1007/s10439-012-0613-5

Bassingthwaighte, J. B., and Beyer, R. P. (1991). Fractal Correlation in Heterogeneous Systems. *Physica D* 53, 71–84. doi: 10.1016/0167-2789(91)90165-6

Bassingthwaighte, J. B., King, R. B., and Roger, S. A. (1989). Fractal nature of regional myocardial blood flow heterogeneity. *Circ. Res.* 65, 578–590. doi: 10.1161/01.RES.65.3.578

Bassingthwaighte, J. B., Yipintsoi, T., and Harvey, R. B. (1974). Microvasculature of the dog left ventricular myocardium. *Microvasc. Res.* 7, 229–249. doi: 10.1016/0026-2862(74)90008-9

Beard, D. A., and Bassingthwaighte, J. B. (2000). The fractal nature of myocardial blood flow emerges from a whole-organ model of arterial network. *J. Vasc. Res.* 37, 282–296. doi: 10.1159/000025742

Carlson, B. E., Arciero, J. C., and Secomb, T. W. (2008). Theoretical model of blood flow autoregulation: roles of myogenic, shear-dependent, and metabolic responses. *Am. J. Physiol. Heart Circ. Physiol.* 295, H1572–H1579. doi: 10.1152/ajpheart.00262.2008

Caulfield, J. B., and Borg, T. K. (1979). The collagen network of the heart. *Lab. Invest.* 40, 364–372.

Chilian, W. M. (1991). Microvascular pressures and resistances in the left ventricular subepicardium and subendocardium. *Circ. Res.* 69, 561–570. doi: 10.1161/01.RES.69.3.561

Cornelissen, A. J., Dankelman, J., Vanbavel, E., and Spaan, J. A. (2002). Balance between myogenic, flow-dependent, and metabolic flow control in coronary arterial tree: a model study. *Am. J. Physiol. Heart Circ. Physiol.* 282, H2224–H2237. doi: 10.1152/ajpheart.00491.2001

Davis, M. J. (2012). Perspective: physiological role(s) of the vascular myogenic response. *Microcirculation* 19, 99–114. doi: 10.1111/j.1549-8719.2011.00131.x

Davis, M. J., Ferrer, P. N., and Gore, R. W. (1986). Vascular anatomy and hydrostatic pressure profile in the hamster cheek pouch. *Am. J. Physiol.* 250, H291–H303. doi: 10.1152/ajpheart.1986.250.2.H291

Downey, J. M., and Kirk, E. S. (1975). Inhibition of coronary blood flow by a vascular waterfall mechanism. *Circ. Res.* 36, 753–760. doi: 10.1161/01.RES.36.6.753

Duncker, D. J., and Bache, R. J. (2008). Regulation of coronary blood flow during exercise. *Physiol. Rev.* 88, 1009–1086. doi: 10.1152/physrev.00045.2006

Fronek, K., and Zweifach, B. W. (1975). Microvascular pressure distribution in skeletal muscle and the effect of vasodilation. *Am. J. Physiol.* 228, 791–796. doi: 10.1152/ajplegacy.1975.228.3.791

Iversen, P. O., and Nicolaysen, G. (1995). Fractals describe blood flow heterogeneity within skeletal muscle and within myocardium. *Am. J. Physiol.* 268, H112–H116. doi: 10.1152/ajpheart.1995.268.1.H112

Jacobs, J., Algranati, D., and Lanir, Y. (2008). Lumped modeling of dynamic flow in coronary vessels *J. Biomechanical Eng*. 130:054504. doi: 10.1115/1.2979877

Johnson, P. C. (1980). “The myogenic response” in *Handbook of Physiology, Section 2: The Cardiovascular System. Vol. II: Vascular Smooth Muscle*, eds D. F. Bohr, A. P. Somlyo, and H. V. Sparks Jr. (Bethesda, MD: American Physiological Society), 409–442.

Jones, C. J., Kuo, L., Davis, M. J., and Chilian, W. M. (1995a). Regulation of coronary blood flow: coordination of heterogeneous control mechanisms in vascular microdomains. *Cardiovasc. Res.* 29, 585–596.

Jones, C. J., Kuo, L., Davis, M. J., Defily, D. V., and Chilian, W. M. (1995b). Role of nitric oxide in the coronary microvascular responses to adenosine and increased metabolic demand. *Circulation* 91, 1807–1813.

Kaimovitz, B., Lanir, Y., and Kassab, G. S. (2005). Large-scale 3-D geometric reconstruction of the porcine coronary arterial vasculature based on detailed anatomical data. *Ann. Biomed. Eng.* 33, 1517–1535. doi: 10.1007/s10439-005-7544-3

Kaimovitz, B., Lanir, Y., and Kassab, G. S. (2010). A full 3-D reconstruction of the entire porcine coronary vasculature. *Am. J. Physiol. Heart Circ. Physiol.* 299, H1064–H1076. doi: 10.1152/ajpheart.00151.2010

Kanatsuka, H., Lamping, K. G., Eastham, C. L., Dellsperger, K. C., and Marcus, M. L. (1989). Comparison of the effects of increased myocardial oxygen consumption and adenosine on the coronary microvascular resistance. *Circ. Res.* 65, 1296–1305. doi: 10.1161/01.RES.65.5.1296

Karch, R., Neumann, F., Neumann, M., and Schreiner, W. (1999). A three-dimensional model for arterial tree representation, generated by constrained constructive optimization. *Comput. Biol. Med.* 29, 19–38. doi: 10.1016/S0010-4825(98)00045-6

Karch, R., Neumann, F., Podesser, B. K., Neumann, M., Szawlowski, P., and Schreiner, W. (2003). Fractal properties of perfusion heterogeneity in optimized arterial trees: a model study. *J. Gen. Physiol.* 122, 307–321. doi: 10.1085/jgp.200208747

Kassab, G. S., and Fung, Y. C. (1994). Topology and dimensions of pig coronary capillary network. *Am. J. Physiol.* 267, H319–H325.

Kassab, G. S., Lin, D. H., and Fung, Y. C. (1994). Morphometry of pig coronary venous system. *Am. J. Physiol.* 267, H2100–H2113.

Kassab, G. S., Rider, C. A., Tang, N. J., and Fung, Y. C. (1993). Morphometry of pig coronary arterial trees. *Am. J. Physiol.* 265, H350–H365. doi: 10.1152/ajpheart.1993.265.1.H350

Keelan, J., Chung, E. M., and Hague, J. P. (2016). Simulated annealing approach to vascular structure with application to the coronary arteries. *R. Soc. Open. Sci.* 3:150431. doi: 10.1098/rsos.150431

King, R. B., Bassingthwaighte, J. B., Hales, J. R., and Rowell, L. B. (1985). Stability of heterogeneity of myocardial blood flow in normal awake baboons. *Circ. Res.* 57, 285–295. doi: 10.1161/01.RES.57.2.285

Kleen, M., Welte, M., Lackermeier, P., Habler, O., Kemming, G., and Messmer, K. (1997). Myocardial blood flow heterogeneity in shock and small-volume resuscitation in pigs with coronary stenosis. *J. Appl. Physiol.* 83, 1832–1841. doi: 10.1152/jappl.1997.83.6.1832

Krams, R., Sipkema, P., Zegers, J., and Westerhof, N. (1989). Contractility is the main determinant of coronary systolic flow impediment. *Am. J. Physiol.* 257, H1936–H1944. doi: 10.1152/ajpheart.1989.257.6.H1936

Kresh, J. Y., Fox, M., Brockman, S. K., and Noordergraaf, A. (1990). Model-based analysis of transmural vessel impedance and myocardial circulation dynamics. *Am. J. Physiol.* 258, H262–H276. doi: 10.1152/ajpheart.1990.258.1.H262

Kuo, L., Chilian, W. M., and Davis, M. J. (1990a). Coronary arteriolar myogenic response is independent of endothelium. *Circ. Res.* 66, 860–866.

Kuo, L., Davis, M. J., and Chilian, W. M. (1988). Myogenic activity in isolated subepicardial and subendocardial coronary arterioles. *Am. J. Physiol.* 255, H1558–H1562. doi: 10.1152/ajpheart.1988.255.6.H1558

Kuo, L., Davis, M. J., and Chilian, W. M. (1990b). Endothelium-dependent, flow-induced dilation of isolated coronary arterioles. *Am. J. Physiol.* 259, H1063–H1070.

Kuo, L., Davis, M. J., and Chilian, W. M. (1995). Longitudinal gradients for endothelium-dependent and -independent vascular responses in the coronary microcirculation. *Circulation* 92, 518–525. doi: 10.1161/01.CIR.92.3.518

Liao, J. C., and Kuo, L. (1997). Interaction between adenosine and flow-induced dilation in coronary microvascular network. *Am. J. Physiol.* 272, H1571–H1581. doi: 10.1152/ajpheart.1997.272.4.H1571

Manor, D., Sideman, S., Dinnar, U., and Beyar, R. (1994). Analysis of flow in coronary epicardial arterial tree and intramyocardial circulation. *Med. Biol. Eng. Comput.* 32, S133–S143. doi: 10.1007/BF02523339

Matsumoto, T., Ebata, J., Tachibana, H., Goto, M., and Kajiya, F. (1999). Transmural microcirculatory blood flow distribution in right and left ventricular free walls of rabbits. *Am. J. Physiol.* 277, H183–H191. doi: 10.1152/ajpheart.1999.277.1.H183

Matsumoto, T., and Kajiya, F. (2001). Microheterogeneity of myocardial blood flow. *Basic Res. Cardiol.* 96, 547–552. doi: 10.1007/s003950170005

Mori, H., Chujo, M., Haruyama, S., Sakamoto, H., Shinozaki, Y., Uddin-Mohammed, M., et al. (1995). Local continuity of myocardial blood flow studied by monochromatic synchrotron radiation-excited x-ray fluorescence spectrometry. *Circ. Res.* 76, 1088–1100. doi: 10.1161/01.RES.76.6.1088

Namani, R., Kassab, G. S., and Lanir, Y. (2018). Integrative model of coronary flow in anatomically based vasculature under myogenic, shear, and metabolic regulation. *J. Gen. Physiol.* 150, 145–168. doi: 10.1085/jgp.201711795

Nevo, E., and Lanir, Y. (1989). Structural finite deformation model of the left ventricle during diastole and systole. *J. Biomech. Eng.* 111, 342–349. doi: 10.1115/1.3168389

Rabbany, S. Y., Kresh, J. Y., and Noordergraaf, A. (1989). Intramyocardial pressure: interaction of myocardial fluid pressure and fiber stress. *Am. J. Physiol.* 257, H357–H364. doi: 10.1152/ajpheart.1989.257.2.H357

Schreiner, W., and Buxbaum, P. F. (1993). Computer-optimization of vascular trees. *IEEE Trans. Biomed. Eng.* 40, 482–491. doi: 10.1109/10.243413

Smith, N. P. (2004). A computational study of the interaction between coronary blood flow and myocardial mechanics. *Physiol. Meas.* 25, 863–877. doi: 10.1088/0967-3334/25/4/007

Spaan, J. A., Breuls, N. P., and Laird, J. D. (1981). Diastolic-systolic coronary flow differences are caused by intramyocardial pump action in the anesthetized dog. *Circ. Res.* 49, 584–593. doi: 10.1161/01.RES.49.3.584

Stapleton, D. D., Moffett, T. C., Baskin, D. G., and Bassingthwaighte, J. B. (1995). Autoradiographic assessment of blood flow heterogeneity in the hamster heart. *Microcirculation* 2, 277–282. doi: 10.3109/10739689509146773

Stepp, D. W., Nishikawa, Y., and Chilian, W. M. (1999). Regulation of shear stress in the canine coronary microcirculation. *Circulation* 100, 1555–1561. doi: 10.1161/01.CIR.100.14.1555

Tillmanns, H., Ikeda, S., Hansen, H., Sarma, J. S., Fauvel, J. M., and Bing, R. J. (1974). Microcirculation in the ventricle of the dog and turtle. *Circ. Res.* 34, 561–569.

Tillmanns, H., Steinhausen, M., Leinberger, H., Thederan, H., and Kübler, W. (1981). Pressure measurements in the terminal vascular bed of the epimyocardium of rats and cats. *Circ. Res.* 49, 1202–1211. doi: 10.1161/01.RES.49.5.1202

Van Bavel, E., and Spaan, J. A. (1992). Branching patterns in the porcine coronary arterial tree. Estimation of flow heterogeneity. *Circ. Res.* 71, 1200–1212.

Vis, M. A., Sipkema, P., and Westerhof, N. (1997). Modeling pressure-flow relations in cardiac muscle in diastole and systole. *Am. J. Physiol.* 272, H1516–H1526. doi: 10.1152/ajpheart.1997.272.3.H1516

Westerhof, N., Boer, C., Lamberts, R. R., and Sipkema, P. (2006). Cross-talk between cardiac muscle and coronary vasculature. *Physiol. Rev.* 86, 1263–1308. doi: 10.1152/physrev.00029.2005

Yipintsoi, T., Kroll, K., and Bassingthwaighte, J. B. (2016). Fractal regional myocardial blood flows pattern according to metabolism, not vascular anatomy. *Am. J. Physiol. Heart Circ. Physiol.* 310, H351–H364. doi: 10.1152/ajpheart.00632.2015

Young, J. M., Choy, J. S., Kassab, G. S., and Lanir, Y. (2012). Slackness between vessel and myocardium is necessary for coronary flow reserve. *Am. J. Physiol. Heart Circ. Physiol.* 302, H2230–H2242. doi: 10.1152/ajpheart.01184.2011

Zamir, M., Vercnocke, A. J., Edwards, P. K., Anderson, J. L., Jorgensen, S. M., and Ritman, E. L. (2015). Myocardial perfusion: characteristics of distal intramyocardial arteriolar trees. *Ann. Biomed. Eng.* 43, 2771–2779. doi: 10.1007/s10439-015-1325-4

Zinemanas, D., Beyar, R., and Sideman, S. (1994). Effects of myocardial contraction on coronary blood flow: an integrated model. *Ann. Biomed. Eng.* 22, 638–652. doi: 10.1007/BF02368289

Keywords: coronary circulation, computer simulation, network reconstruction, morphometry, diameter assignment, perfusion dispersion

Citation: Namani R, Kassab GS and Lanir Y (2018) Morphometric Reconstruction of Coronary Vasculature Incorporating Uniformity of Flow Dispersion. *Front. Physiol.* 9:1069. doi: 10.3389/fphys.2018.01069

Received: 03 March 2018; Accepted: 17 July 2018;

Published: 29 August 2018.

Edited by:

Timothy W. Secomb, University of Arizona, United StatesReviewed by:

Brendan Fry, Metropolitan State University of Denver, United StatesWalter Lee Murfee, Tulane University, United States

Copyright © 2018 Namani, Kassab and Lanir. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Yoram Lanir, yoramlanir@yahoo.com