Development of Novel Fractal Method for Characterizing the Distribution of Blood Flow in Multi-Scale Vascular Tree

Blood perfusion is an important index for the function of the cardiovascular system and it can be indicated by the blood flow distribution in the vascular tree. As the blood flow in a vascular tree varies in a large range of scales and fractal analysis owns the ability to describe multi-scale properties, it is reasonable to apply fractal analysis to depict the blood flow distribution. The objective of this study is to establish fractal methods for analyzing the blood flow distribution which can be applied to real vascular trees. For this purpose, the modified methods in fractal geometry were applied and a special strategy was raised to make sure that these methods are applicable to an arbitrary vascular tree. The validation of the proposed methods on real arterial trees verified the ability of the produced parameters (fractal dimension and multifractal spectrum) in distinguishing the blood flow distribution under different physiological states. Furthermore, the physiological significance of the fractal parameters was investigated in two situations. For the first situation, the vascular tree was set as a perfect binary tree and the blood flow distribution was adjusted by the split ratio. As the split ratio of the vascular tree decreases, the fractal dimension decreases and the multifractal spectrum expands. The results indicate that both fractal parameters can quantify the degree of blood flow heterogeneity. While for the second situation, artificial vascular trees with different structures were constructed and the hemodynamics in these vascular trees was simulated. The results suggest that both the vascular structure and the blood flow distribution affect the fractal parameters for blood flow. The fractal dimension declares the integrated information about the heterogeneity of vascular structure and blood flow distribution. In contrast, the multifractal spectrum identifies the heterogeneity features in blood flow distribution or vascular structure by its width and height. The results verified that the proposed methods are capable of depicting the multi-scale features of the blood flow distribution in the vascular tree and further are potential for investigating vascular physiology.


INTRODUCTION
The microcirculation is the end destination of the cardiovascular system and the patency of microvascular perfusion is essential for the maintenance of tissue metabolism (Ince, 2005;Guven et al., 2020). Various cardiovascular diseases influence the blood perfusion and thus impair the physiological function of organs (Efimova et al., 2008;Kitagawa et al., 2009;Alosco et al., 2013Alosco et al., , 2014. These findings imply that the blood perfusion may act as an important index for the physiological states of living bodies. The blood perfusion can be indicated by the blood flow distribution in the vascular tree. The blood flow in a vascular tree is distributed at different generations, varying in a large range of scales, and the blood flow distribution at a certain generation is directly affected by the superior generation. In the meantime, there is a huge difference between the magnitude of the blood flow at different generations. However, the conventional statistical parameters for characterizing the blood flow distribution, like the coefficient of variation (CV) (Bassingthwaighte et al., 2001;Pries and Secomb, 2009), ignored the connection of the blood flow among multiple scales. To develop a unified description of the blood flow distribution covering all scales remains a big challenge.
To depict the scale-independent characteristic of objects, the fractal theory provides an efficient approach for multi-scale analysis (Mandelbrot, 1982). Presently, a few studies have made an effort on investigating the fractal characteristics of blood flow distribution in the vascular tree directly or indirectly. Van Beek et al. (1989) uncovered the fractality of the relative dispersion of blood flow distribution. Zamir (2001) defined the fractal dimension based on the relationship between the vessel diameter and blood flow according to Murray's law. Grasman et al. (2003) described that the distribution of blood flow at the same generation is multifractal. In all these studies, the unified description for a vascular tree by the fractal parameter all demands that the vascular tree should be a perfect binary tree, in which all interior branch nodes have two daughter branches and all terminals have the same depth or generation. However, the structures of real vascular trees are diverse which limits the physiological application of the methods above.
The fractal analysis has been widely used to investigate the geometrical characteristics of the vasculatures (Cheng and Huang, 2003;Stosic and Stosic, 2006;Lorthois and Cassot, 2010;Gould et al., 2011;Nadal et al., 2020). It inspires us to introduce the established fractal methods for geometrical architecture analysis into the hemodynamic study, and further develop a universal fractal depiction for blood flow distribution. The conservation law is common during the emergence of fractal and multifractal (Hassan, 2019). In fractal geometry, this law presents as the conservation of the number of signal pixels in an image. On the other hand, the total volume of blood flow in the vascular tree also obeys the conservation law. This consistency makes it possible to apply the principle of the fractal method for geometry to the analysis of blood flow by appropriate modification.
In this study, the primary aim is to establish fractal methods for analyzing the blood flow distribution which is potential to be applied to real vascular trees. To achieve this goal, we firstly modified the fractal methods in geometry to accommodate the situation of blood flow and then applied the established methods on experimental data to test the validity. Further, to explore the physiological significance of the yielded fractal parameters, the blood flow distribution in vascular trees with fixed structure or with varying structures were examined in which the hemodynamics was simulated based on a hemodynamic model (Yang and Wang, 2013) and a rheological model (Pries and Secomb, 2005).

Establishment of Fractal Methods for Blood Flow
The fractal dimension is the most important parameter to quantify the fractality of objects. And measuring the information dimension is an efficient way to estimate the fractal dimension in geometry (Pitsianis et al., 1989;Liu et al., 2018). For the calculation of information dimension, non-overlapping boxes are adopted to cover the image of the object and the mass probability of each box, which is defined as the ratio of the number of signal pixels in the box to that of the whole image, is obtained. And the information dimension D I (Pitsianis et al., 1989) is estimated as: where N(L) i=1 −P i log P i is the total entropy of mass according to the information theory, P i is the mass probability of the ith box and N(L) is the number of boxes needed to cover the image with size L.
The total mass, which is the number of signal pixels in the whole image, obeys the law of conservation regardless of the box size. And so does the blood flow. As shown in Figure 1, the total volume of blood flow at the same generation also follows the law of conservation no matter how many times the vascular tree bifurcates. Thus, with appropriate modification, the fractal methods in fractal geometry can be introduced to investigate the fractality of blood flow. By replacing the mass probability P i in Eq. 1 with the flow probability p i and 1/L with the number of vessel segments N(g) at generation g, the expression of the fractal dimension for the blood flow D Q is derived as: Practically, D Q is estimated as the slope of the linear fitting curve of the discrete data pair (log N g , −p i log p i ) (Wang et al., 2019). The precondition of the fractal analysis is scale invariance. And the strong linearity over three orders of magnitude of the fitting curve can be taken as the criteria of the existence of scale invariance (Halley et al., 2004).
For the fractal analysis in geometry, the box-counting dimension (So et al., 2017;Nayak et al., 2019), also known as capacity dimension, is the most popular. This method also requires the image to be covered by non-overlapping boxes. However, the box-counting dimension method only considers the existence of signal pixels in the box but ignores the number of pixels. When applied to blood flow analysis, it produces fractal dimension about the vascular structure other than the distribution of blood flow. In contrast, the information dimension method takes the quantity of blood flow into account and thus can reflect the blood flow distribution.
Very few objects possess perfect mono-fractality exhibiting a single fractal dimension (Gould et al., 2011). In reality, objects with the subsets having different scaling properties are much more common and the estimation of multifractality is more desirable. For the multifractal measure of blood flow, the multifractal spectrum f (α)∼α of the blood flow is adopted (Chhabra and Jensen, 1989) and modified as: in which where q is the moment order. And the range of the spectrum α = α max − α min can be used to measure the degree of multifractality (Halsey et al., 1986).

Generalization of the Established Methods
The methods given in Section "Establishment of Fractal Methods for Blood Flow" are based on the premise that the total blood flow at the same generation in a vascular tree obeys the law of conservation. This premise is valid for a perfect binary tree, as shown in Figure 1, but not for real vascular trees as shown in Figure 2A. If the branch which stops bifurcating before reaching the maximal generation is regarded as a branch covering multiple generations, the vascular tree in Figure 2A can be thought of as the perfect binary tree as shown in Figure 2B. In this case, the proposed methods can be applied and a strategy for the calculation is raised: If a vessel segment stops bifurcating at generation n (n is smaller than the maximal generation of vascular tree), it will be involved in the calculation at all the generations greater than n.
That is to say, for the vascular tree in Figure 2A, V 2,1 will be included in the calculation of total entropy at generation 3 and V 1,2 will be included in the calculation at both generation 2 and generation 3. It should be noted that the number of vessel segments N in Eq. 2 should be 2 n at generation n but not the number of vessel segments at this generation. The reason is explained below.  The blood flow distribution in the vascular tree in Figure 2B can be analogous to a one-dimension graph as shown in Figure 2C. The gray segments indicate the range of the graph and the black segments reflect the blood flow in vessel segments. Each bisection corresponds to a bifurcation of the vascular tree and divides this geometrical structure into smaller subsections. When adopting Eq. 1 for the estimation of fractal dimension for this one-dimension graph, L is corresponding to the length of the smallest subsection. And each bisection of the graph halves L no matter whether there is always a black segment in each subsection. The number of subsections should be 2 n after n bisections. Similarly, the number of vessel segments N(g) in Eq. 2 should be 2 n at generation n.

Vascular Tree Construction and Hemodynamic Simulation
By now, the methods established in "Establishment of Fractal Methods for Blood Flow" and "Generalization of the Established Methods" can be used to characterize the blood flow distribution in arbitrary vascular trees. For validation, the established fractal methods were tested in a real arterial tree (Reglin et al., 2009(Reglin et al., , 2017Wang et al., 2019) under normal and ischemic state. Besides, the capability of the derived fractal parameters was examined in two situations.
For the first situation, the vascular tree was fixed to be a perfect binary tree and the blood flow distribution was adjusted by the split ratio. For a bifurcation with a parent vessel segment and two daughter branches, the split ratio r is defined as the ratio of the smaller blood flow to the larger one in the two daughter branches, ranging from 0 to 1. Assuming that r is constant throughout the perfect binary tree, the fractal dimension for blood flow can be obtained as shown in Eq. 6. The detailed derivation of the equation is given in the Appendix.
For the second situation, a series of vascular trees were constructed. The structures of these vascular trees were diverse while the blood flow distribution was estimated under the same boundary condition. The successive dichotomous division is the most common branching pattern of the vascular tree, in which a parent vessel segment is divided into two daughter branches (Zamir, 2001).
Based on this pattern, the construction of a vascular tree calls for the determination of the vessel diameter and length. For a bifurcation with the diameter of the parent vessel being d 0 and those of the two daughter branches being d 1 and d 2 , a powerlaw relationship between the diameters is given as shown in Eq. 7 according to Murray's law.
where d 0 > d 1 ≥ d 2 , k is bifurcation exponent and λ is asymmetry ratio. It is reported that the k value varies from 2.33 to 3.0 (Gabrys et al., 2005). And a value above 0.6 is most commonly observed for λ (Schmidt et al., 2004;Cheung et al., 2011;Takahashi, 2014). Based on the power-law relationship, the diameters of all vessel segments in a vascular tree can be estimated with the given root diameter and cut-off diameter. In this study, the root diameters of all constructed vascular trees are set as 300 µm and all the terminal diameters are 10 µm, the size of capillaries. An empirical formula (Takahashi et al., 2009) is adopted to obtain the length l (µm) from the diameter d (µm): With the estimated vessel diameters and lengths, a vascular tree can be constructed, and serves for the hemodynamic simulation. According to Hagen-Poiseuille's law as shown in Eq. 9, the blood flow Q (µm 3 /s) in a vessel segment is proportional to the pressure drop P (Pa) between the inlet and outlet. And the flow resistance R is determined by vessel diameter d (µm), vessel length l (µm), and blood viscosity µ (Pa · s).
Fahraeus and Lindqvist (1931) reported a decline in apparent blood viscosity with decreasing tube diameter, the so-called Fahraeus-Lindqvist effect. Among the models to describe the relationship between the blood viscosity µ (mPa · s) and the diameter of the vessel segment d (µm), the model proposed by Pries and Secomb (2005) matches well with the in vivo experimental data. In this model, the in vitro viscosity is firstly estimated as: in which H d (%) is the hematocrit and: C = 0.8 + e −0.075d · 1 1 + 10 −11 · d 12 − 1 + 1 1 + 10 −11 · d 12 Frontiers in Physiology | www.frontiersin.org The estimation of the in vivo viscosity should take account of the effect of the endothelial surface layer. This will involve two parameters, the effective diameter d eff = d − 2W eff and the physical diameter d ph = d − 2W ph . The effective thickness of the layer W eff and physical thickness of the layer W ph are estimated as: W ph = W as + W peak · E peak (14) Based on the experimental data, E HD = 1.18, E peak = 0.6, E amp = 1.1, E width = 0.03, D off = 2.4 µm, D crit = 10.5 µm, D 50 = 100 µm, and W max = 2.6µm. By replacing the d in Eqs 11 and 12 with d ph , we can get the in vitro viscosity µ vitro . And the in vivo viscosity µ vivo can be obtained as below.
For the bifurcation as shown in Figure 3, the relationship between the hemodynamic parameters can be given by Eq. 18 based on Hagen-Poiseuille's law and the conservation law of flow.
Without losing generality, in the hemodynamic simulation we prescribed the inlet pressure as 1 mmHg and the outlet pressure at all terminal branches as 0 mmHg (Yang and Wang, 2013).
With each bifurcation of the vascular tree assigned an equation like Eq. 18, the blood pressure at each branch node can be obtained by solving these equations. Further, the blood flow in each vessel segment is estimated by Eq. 9 and finally the blood flow distribution in a tree can be acquired.

Numerical Solution
In this study, all the calculations and simulations were programmed by MATLAB R2019a (MathWorks Co., MA, United States). Firstly, the node information for each constructed vascular tree was obtained. Then, the blood flow in each vessel segment of the constructed vascular tree was captured by solving the equations set. Ultimately, the fractal dimension and multifractal spectrum were calculated.
All the results about the fractal parameters were presented as Mean ± SD.

The Validation of the Proposed Methods
To test the validity of the proposed methods, the fractal, and multifractal analysis were conducted on a real arterial tree under normal and ischemic state. The blood flow distribution in these two states is as shown in Figures 4A,B. The fractal dimension for the normal and ischemic state are 0.53 and 0.40, respectively. As for the multifractal spectrum for the blood flow, the results are shown in Figures 4C,D. We can observe that the multifractal spectrums for both states appear as curves indicating the existence of multifractality. While the maximal values of the multifractal spectrum for the two states are the same, the range of the multifractal spectrum α for the ischemic state is wider than that for the normal state which indicates a higher degree of multifractality for the ischemic state. The results verified the ability of the proposed methods in distinguishing different blood flow distribution in a real vascular tree.

The Fractality of Blood Flow in the Perfect Binary Vascular Tree
The blood flow distribution in the perfect binary vascular tree was evaluated by both the fractal parameters and the CV, which is defined as the standard deviation divided by the mean value. For a perfect binary vascular tree with the identical split ratio for each bifurcation, the fractal dimension for blood flow was obtained based on Eq. 6. Figure 5A shows the trends of fractal dimension and the coefficient of variation with the change of split ratio. It is noticed that the fractal dimension increases monotonically from 0 to 1 with the increment of the split ratio. In the case r = 1, the distribution of blood flow has the highest value of fractal dimension. And the CV decreases from 10.1 to 4.4 with the increment of the split ratio. As for the multifractal characteristic, the multifractal spectrums of the blood flow are shown in Figure 5B. With the decrement of the split ratio, the range of the multifractal spectrum α expands while the maximal value remains unchanged. It is worth noting that the multifractal spectrum is presented as a point when r = 1, implying the absence of multifractality.
We also examined the fractal dimension and the CV of the blood flow in perfect binary vascular trees with different maximal generations. As shown in Figure 6, the fractal dimension holds steady with the change of maximal generation while the CV varies greatly.

The Fractality of Blood Flow in Constructed Vascular Trees
Vascular trees with diverse structures were constructed. To reflect the heterogeneity in the real vascular tree, the bifurcation exponent k, and asymmetry ratio λ of the bifurcations in each constructed vascular tree were set following the normal distribution. And ten vascular trees were constructed for each pair of k and λ. The statistical characteristics of these vascular trees are shown in Table 1. In these constructed vascular trees, the hemodynamics was simulated and the results are shown  in Figure 7. We can see that there exists a strong linear relationship between the logarithmic values of the diameter and the blood flow rate.
The fractal dimensions for blood flow in the constructed vascular trees are shown in Figure 8. For the vascular trees with the mean k value of 2.7, the fractal dimension increases monotonically from 0.74 ± 0.01 to 1.00 ± 0.00 with λ rising from 0.60 to 1. While for the vascular trees with the mean λ value of 0.8, the fractal dimension fluctuates between 0.91 and 0.95 with the increment of k.  The multifractal spectrums of the blood flow for these vascular trees are presented in Figures 9A,C. It is observed that the multifractality of blood flow exists in all vascular trees but the multifractal spectrums fluctuate. With the increment of λ, the range of the multifractal spectrum α narrows from 0.89 ± 0.04 to 0.04 ± 0.01 and the maximal value grows from 0.93 ± 0.01 to 1.00 ± 0.00. And with the increment of k, α expands from 0.32 ± 0.04 to 0.65 ± 0.03 and the maximal value grows from 0.95 ± 0.01 to 0.99 ± 0.00. The features of the multifractal spectrums are shown in Figures 9B,D. Compared with the fractal dimensions, the difference of multifractal spectrums among these vascular trees is more striking.

Validity of the Hemodynamic Simulation
The hemodynamic simulation was conducted to investigate the variation of fractal parameters with varying blood flow distribution. To make sure that the obtained blood flow distribution is reasonable, a quantitative comparison of the hemodynamic simulation with the existing physiological studies is necessary. For avoiding losing the generality, the boundary condition in the present work was prescribed with an inlet pressure of 1 mmHg and an outlet pressure of 0 mmHg. However, the pressure drop between the inlet and outlet may vary in different studies. Thus, for quantitative comparison, it is more appropriate to examine the relative indices.
As shown in Figure 7, a strong linear relationship exists between the logarithmic values of the diameter and the blood flow rate. This is consistent with the assertion in Huo and Kassab (2016) that there is a scaling law between the blood flow rate and diameter. The slope of the fitting line indicates the relation between the blood flow rate and vascular diameter thus can be an indicator for quantitative comparison. Within a similar diameter range, the slope with a value of 2.40 produced in our work is comparable with the reported work of 1.97 (Wang et al., 2009), 2.0 ± 0.2 (Pijewska et al., 2020), 2.33 (Huo and Kassab, 2012), and 2.49 ± 0.09 (Haindl et al., 2016). Thus, it can be concluded that the simulation results of our work are reasonable.

Physiological Significance of the Fractal Parameters
Two fractal parameters, i.e., fractal dimension and multifractal spectrum, were obtained in this study to investigate the fractality and multifractality of blood flow. By definition, the fractal dimension is determined by the total entropy of blood flow. And the total entropy is calculated by considering the existence as well as the quantity of blood flow in the vessel segment. The existence and the quantity of blood flow are corresponding to the vascular structure and the blood flow distribution, respectively. Thus, the fractal dimension characterizes the combination of the features of vascular structure and blood flow distribution. When the vascular structure is fixed, the lower entropy is obtained from the more heterogeneous distribution according to the information theory. That is to say, the fractal dimension reflects the degree of the blood flow heterogeneity for a specific vascular tree and the lower fractal dimension comes from the blood flow distribution with a higher degree of heterogeneity. The results in Figure 5A that the lower fractal dimension is corresponding to the lower split ratio also support this conclusion.
For a fractal object, the multifractal spectrum describes the scaling properties in different subsets. And when multifractality presents, the subsets of this object will be scaled by different multiples at the same q order moment during the calculation of the multifractal spectrum. Therefore, the degree of multifractality, which is measured by the width of the multifractal spectrum, essentially describes the degree of heterogeneity within the fractal object and it rises with the increment of the blood flow heterogeneity. This judgment is consistent with the results as shown in Figure 5B. The multifractal spectrum f (α(q)) reaches its maximal value when q = 0. In this case, the quantity of the blood flow volume no longer has an effect on the value of f (α(0)). In other words, the height of the multifractal spectrum reflects the heterogeneity or asymmetry of the vascular structure. The higher the multifractal spectrum the closer the vascular tree is to the perfect binary tree. And this makes clear why the height of the multifractal spectrums in Figures 4C,D or Figure 5B is the same. By means of the width and height, the multifractal spectrum separates the information about the vascular structure and blood flow distribution. In this sense, the multifractal spectrum makes the evaluation of the blood flow distribution in different vascular trees possible. As shown in Figure 9D, the heterogeneity of the vascular structure decreases with the increment of the bifurcation exponent while the heterogeneity of the blood flow distribution increases. The interaction of these two opposite trends may explain why the fractal dimension changes slightly with the bifurcation exponent as shown in Figure 8B.
Both the fractal dimension and the multifractal spectrum reflect the blood flow heterogeneity. Physiologically speaking, the change of blood flow heterogeneity is usually associated with pathological conditions. For microcirculation, the increment of blood flow heterogeneity can be an early indicator of diseases, such as sepsis and shock (Ince, 2005;Dubin et al., 2018;Ince et al., 2018) as well as peripheral vascular disease (Butcher et al., 2013). And the increase of blood flow heterogeneity can be depicted by the decrease of fractal dimension and the broadening of the multifractal spectrum. The multifractal spectrum can separate the information about the vascular structure and blood flow distribution. Thus, the multifractal spectrum is also able to distinguish the causes responsible for the change in blood flow heterogeneity, either by hemodynamic problems or by structural alteration due to the diseases such as large vessel stenosis (Kharche et al., 2018).
There are also some other quantitative or semi-quantitative methods for characterizing the blood flow heterogeneity Pries and Secomb, 2009;Ince et al., 2018). In this study, the blood flow heterogeneity is also evaluated by CV. The results in Figure 5A confirm the availability of this statistical parameter. However, the CV treats different vessel segments in a vascular tree as independent components ignoring the connection of blood flow along the whole tree. This would make this parameter less accurate in some cases as discussed below. When the split ratio for each bifurcation in the perfect binary vascular tree is 1, the blood flow is evenly distributed for each generation in the vascular tree. The degree of blood flow heterogeneity should remain unchanged no matter what the maximal generation of the vascular tree is. And this property holds for the other values of the split ratio. As shown in Figure 6A, the fractal dimension stays the same with the change of the maximal generation. However, the values of CV for the vascular trees with different maximal generations are quite different. Considering both the connection and difference of blood flow in different vessel segments, the fractal parameters can provide a more accurate description of the blood flow heterogeneity for the tree-like vasculature.

Limitations
It should be pointed out the hemodynamic simulation in this study was simplified. Nowadays, the RCL model has been developed for hemodynamic simulation in which the resistance (R), capacitance (C), and inductance (L) elements were used to mimic the effects of vessel resistance, vessel compliance, and blood inertia, respectively (Muller and Toro, 2014;Zhang et al., 2014). And models from 0D to 3D were established (Arciero et al., 2017;Liu et al., 2020). In the present hemodynamic model of  microcirculation, only the resistance element for a vessel segment was considered. Although in the microcirculation the resistance element plays a dominant role in hemodynamics (Katanov et al., 2015;Nichols et al., 2015;Secomb, 2017) and the results show that the model is sufficient for produce varying blood flow distribution in a tree, a comprehensive model is worth being introduced in the future study.

CONCLUSION
In this study, the fractal methods were introduced, with appropriate modification, to characterize the multi-scale properties of blood flow. The application of the methods to the real physiological data verified its ability in distinguishing the variety of blood flow distribution. The yielded parameters, as the fractal dimension and the multifractal spectrum for blood flow, can quantify the degree of blood flow heterogeneity. With the increase of blood flow heterogeneity, the fractal dimension decreases and the multifractal spectrum expands. And the investigation on various constructed vascular trees suggests that both the vascular structure and the blood flow distribution influence the fractal parameters. With the aid of the fractal dimension, it is possible to look into the change of blood flow heterogeneity in a specific vascular tree. While the multifractal spectrum can be utilized to assess the blood flow heterogeneity for different vascular trees by considering the blood flow distribution and the structure of vascular trees separately. It can be concluded that the proposed methods provide efficient tools to describe the multi-scale properties of the blood flow distribution and has the potential to assist the study of multi-scale vascular physiology.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding authors.

AUTHOR CONTRIBUTIONS
PL and QP: conceptualization and methodology. PL: algorithm and writing (original draft). SJ and GN: writing (review and editing). MY and JY: discussion of the results and their relevance. GN and JY: supervision and project administration. All authors approved the manuscript.

FUNDING
This study was supported by the National Natural Science Foundation of China (Grant Numbers 81871454 and 31870938).