Dynamic Bayesian Network Modeling of Vascularization in Engineered Tissues Caner Komurlu Jinjian Shao Mustafa Bilgic Computer Science Department Computer Science Department Computer Science Department Illinois Institute of Technology Illinois Institute of Technology Illinois Institute of Technology Chicago, IL, 60616 Chicago, IL, 60616 Chicago, IL, 60616 ckomurlu@hawk.iit.edu jshao3@hawk.iit.edu mbilgic@iit.edu Abstract depths of the tissue, and form connections to allow blood circulation. In this paper, we present a dynamic Bayesian The formation of new blood vessels are triggered and network (DBN) approach to modeling vascu- a↵ected by growth factors that are released by dis- larization in engineered tissues. Injuries and tressed cells that are far from the existing blood ves- diseases can cause significant tissue loss to sels. When these growth factors reach existing blood the degree where the body is unable to heal vessels, they sprout new branches and these branches itself. Tissue engineering aims to replace the “search” for the distressed cells by following the gradi- lost tissue through use of stem cells and bio- ent of the growth factor. This process, however, is materials. For tissue cells to multiply and stochastic for at least two reasons: i) even though migrate, they need to be close to blood ves- growth factors are the main ingredients for causing sels, and hence proper vascularization of the sprouts, they are not the only elements that a↵ect vas- tissue is an essential component of the en- cularization, and ii) the growth factors are increasingly gineering process. We model vascularization more uniformly distributed as they go further away through a DBN whose structure and parame- from the distressed cells, and hence the gradient is al- ters are elicited from experts. The DBN pro- most uniform, hindering the capability of the blood vides spatial and temporal probabilistic rea- vessel finding its way correctly. soning, enabling tissue engineers to test sen- This inherent stochasticity in the vascularization pro- sitivity of vascularization to various factors cess, the spatial nature of the tissue, and the temporal and gain useful insights into the vasculariza- aspect of the vascularization make temporal graphical tion process. We present initial results in this models a great fit for reasoning with uncertainty in paper and then discuss a number of future re- vascularization. In this paper, we present a dynamic search problems and challenges. Bayesian network (DBN) for modeling vascularization in engineered tissues. We elicit the structure of the DBN from tissue engineering experts and we experi- 1 INTRODUCTION ment with various parameter settings to provide fur- ther insights into the vascularization process. Because People lose tissue due to accidents, medical opera- this is a first and novel application of DBNs to tissue tions, treatments, and illnesses. While some organs, engineering, it avails itself to many interesting future e.g. liver, can replace the lost tissue most cannot espe- research directions and challenges. cially when the damage is too severe. For these kinds Our contributions in this paper include: of tissue damages, the lost tissue can be replaced by engineering a new tissue through stem cells and bio- materials [18]. • We present a novel application of DBNs to vascu- larization in engineered tissues An essential process for engineering a healthy tissue is the proper vascularization (formation of new blood • We present initial results and insights, where we vessels) of the tissue, as the tissue cells need to be experiment with various parameter settings, and close to the blood vessels both to discharge their waste and to receive nutrition and oxygen. The blood ves- • We discuss several future research challenges and sels need to spread out in the tissue, invade into the opportunities in detail. 89 The rest of the paper is organized as follows: in Sec- tion 2, we provide a brief background on tissue engi- neering and vascularization. In Section 3, we describe our DBN model for vascularization. We present our experimental setup and results in Section 4. In Sec- tion 5, we briefly discuss related work. We then dis- cuss future research directions and challenges in detail in Section 6, and then conclude. 2 BACKGROUND In this section, we first provide a brief background on tissue engineering and vascularization and then discuss briefly why dynamic Bayesian networks (DBNs) are a good fit for modeling vascularization. People lose tissue due to accidents, treatments, and Figure 1: Illustration of vascularization, including the illnesses. Some organs, e.g. liver, can replace the lost tip cells (active cells), the fixed cells (stalk cells), and tissue while others cannot. Sometimes, the damage anastomosis. [19] can be so severe that the body cannot heal itself. For example, bones can heal after smooth fractures. Yet, some fractures damage bone body so severely that the bone cannot regenerate. For these kinds of damages, Vascularization is a key process in tissue development. the lost tissue can be replaced by engineering a new When cells that are emitting VEGF cannot be reached tissue through stem cells and biomaterials. in time by the new blood vessels, the cells first fall in hypoxia (i.e., lack of Oxygen) and then start dy- Stem cells are generic types of cells that have the abil- ing. Hence the formation of healthy tissue depends on ity to replicate and transform to any tissue. Stem appropriate vascularization; the blood vessels need to cells, like all other cells, need to be close enough to the spread out in the newly-formed tissue, invade into the blood vessels so that they can forward their biological depth, and need to form connections to allow blood wastes to the vessels and they can be fed with nutri- circulation. tion and oxygen carried by the blood vessels. When a tissue is engineered through replication and trans- Though it is well-known that the VEGF is a major formation of stem and tissue cells, there is no existing contributor to sprouting of new blood vessels and that blood vessel web in the environment; the only blood the tip of the blood vessel typically follows the gradient vessels available are the original vessels located at the of the VEGF, there are still unknown factors that af- edges, ready to sprout and progress to the depths of fect vascularization. Moreover, the VEGF distribution the newly-formed tissue. becomes more uniform as we get further away from the source of the emission and hence the gradient does The stem cells that do not have access to blood vessels not necessarily point to the distressed cell. Therefore, will not be able to discharge waste and receive nutri- given our knowledge of the VEGF distribution the en- tion and oxygen. In such cases, a cell starts signaling vironment, the blood vessels do not necessarily follow about its needs by means of emitting chemicals called a deterministic path; they also do a bit of exploration. vascular endothelial growth factor (VEGF). VEGF dif- This is where the uncertainty reasoning capabilities of fuses and disperses in the environment. When it con- probabilistic graphical models become handy for mod- tacts a blood vessel, it triggers a new sprout of blood eling vascularization. vessel towards the source of emission. The tip of these new sprouts typically follow the gradient of the VEGF In this paper, we model the vascularization process to find the distressed cell. During this process, the through dynamic Bayesian networks (DBNs) to enable newly-formed blood vessel can also branch and sprout tissue engineering researchers to reason with spatial new blood vessels. When the branches meet with other and temporal growth of blood vessels. With the help branches, they merge (this process is called anastomo- of DBNs, the researchers can formulate and query the sis) and a blood circulation through the new vessel DBNs and try a number of parameter settings, without starts. The blood circulation helps nearby stem and the need to experiment with every one of them in the tissue cells, which then stop emitting growth factors. lab. This process allows the researchers to gain further This event is called angiogenesis or vascularization. insights and formulate new in-vivo (on animals) and Please see Figure 1 for an illustration of this process. in-vitro (on glass) experiments. 90 3 APPROACH In this section, we describe our DBN model for vas- 𝐿 𝐿( ) cularization. We made a number of assumptions to simplify the model. In this model, we assume a 2D structure, whereas in real-life scenarios, the tissue ob- viously has a 3D structure. In this 2D structure, which is illustrated in Figure 2, as also assumed in [1], we as- sume that the blood vessel grows bottom-up towards north. Therefore, the status of a location at time t 𝐿 𝐿 𝐿 depends on: i) its status at time t, and ii) the statuses of its south neighbors at t. t t+1 Figure 3: A two-time slice representation of the DBN. A location at a time t + 1 has four parents: itself at ( ) time t and its lower neighbors at time t. 𝐿 𝐿 AC ( ) ( ) ( ) 𝐿( )( ) 𝐿 ( ) 𝐿( )( ) 𝐿( 𝐿 ( 𝐿( )( ) ) )( ) T F t t+1 PAC PSC PE SC Figure 2: The tissue grid. Each cell of the grid rep- T F resents a location, which can be Empty, or can be oc- cupied with an Active Cell or Stalk Cell. Each PAC PSC PE location is represented as a random variable in DBN. To simplify the notation, when we refer to a generic location Ltxy , we will drop the subscripts and hence simply use Lt , and when we refer to its neighbors Figure 4: The CPD for P (L(t+1) |Lt , LtSW , LtS , LtSE ). at its south Lt(x 1)(y 1) , Ltx(y 1) , and Lt(x+1)(y 1) we will simply use LtSW , LtS , and LtSE , corresponding to neighbors at south west, south, and south east, respec- • The tip of a blood vessel (AC) at time tively. We illustrate the relevant 2-time slice dynamic t becomes the body (SC) at time t + 1. Bayesian network in Figure 3. That is P L(t+1) |Lt = AC, LtSW , LtS , LtSE = Each location on the 2D grid is a random variable, P L(t+1) |Lt = AC = h✏, 1 2✏, ✏i, where ✏ is a representing whether that location is Empty, or occu- small noise parameter. pied by a blood vessel cell. Blood vessel cells are two • A Stalk Cell at time t either continues types: the tip of a blood vessel that has the potential to to remain a Stalk Cell at time t + 1 or grow (henceforth called an Active Cell) or the body it might become Active Cell with probabil- of the blood vessel (henceforth called the Stalk Cell). ity to sprout a new blood vessel branch. Therefore, the domain of random variable is [Active That is, P L(t+1) |Lt = SC, LtSW , LtS , LtSE = Cell, Stalk Cell, Empty], abbreviated henceforth as [AC, SC, E]. P L(t+1) |Lt = SC = h , 1 ✏, ✏i. We refer to as the sprout possibility. We model the conditional probability distribution, (CPD), P L(t+1) |Lt , LtSW , LtS , LtSE as a tree CPD as • An Empty location at time t will remain Empty illustrated in Figure 4. To give a simple overview, at time t + 1 if none of its SW, S, or SE neigh- at each step in time, an Active Cell elongates and bors are Active Cell at time t; if there is an moves into a nearby Empty location, forming the body Active Cell at one or more of those neighbor- of the blood vessel (i.e., Stalk Cell) in the process. ing locations at time t, one of them might elon- The transitions are: gate to this Empty location at time t + 1. The 91 probability of that an Empty location being oc- grid over three time slices. Then, we present results cupied by an Active Cell at time t + 1 is mod- on a bit larger scale, 9 ⇥ 9, over nine time slices. Fi- eled as a Noisy-OR of its neighboring locations. nally, we present a framework where we quantify the That is P L(t+1) = AC|Lt = E, LtSW , LtS , LtSE uncertainty over the predictions on the last time slice is a Noisy-OR of LtSW , LtS , LtSE , with parameters and discuss how it is a↵ected by the growth patterns 0 , SW , S , and SE , where 0 is leak param- and sprout possibilities. eter, and SW , S , and SE corresponds to the For inference, in the 3⇥3 case, we used exact inference. possibility that an Active Cell elongates in the For the 9 ⇥ 9 case, we used forward sampling. Note NE, N, or NW direction.1 The magnitude of SW , that we are able to use forward sampling in our settings S , and SE are determined by the VEGF gradi- because we provide the initial condition (all locations ent. We refer to various configurations of the at time t = 0) as evidence and compute probabilities parameters as the growth patterns. for the remaining time slices. 4 EXPERIMENTAL SETUP, 4.1 Detailed Results for 3 ⇥ 3 RESULTS, AND INSIGHTS In this toy setting, we provide the evidence for the In this section, we describe the experiments we per- initial configuration of the experiment, i.e., we provide formed using various settings for the growth pattern evidence for all locations for time t = 0, and compute ( ) and sprout ( ) parameters. In all the experiments probabilities for all locations for all future time slices. to follow, we set the noise ✏ and the leak 0 parameters That is, we compute P (L1 , L2 |L0 ), where Lt denotes to 0.01. For the growth pattern, we present results for all locations at time t. For t = 0, we provide the two settings: evidence as follows: the middle of the bottom row is set as the tip of the blood vessel (i.e, L0x=1,y=0 = AC) and the rest of the locations are set as Empty. Figure 5 • straight-growth: h SW , S , SE i = illustrates this setting. h0.01, 0.98, 0.01i. For this pattern, the blood vessel follows a straight line, growing towards north. E E E • uniform-growth: h SW , S , SE i = h 13 , 13 , 13 i. For this pattern, the blood vessel has equal chance E E E of growing towards north, north west, or north east. E AC E For the sprout possibility, that is a Stalk Cell turn- ing into an Active Cell, we present results for two Figure 5: The initial configuration for the 3 ⇥ 3 grid. settings: The straight-growth results are presented in Figures 6 and 7, and uniform-growth results are presented in • seldom-sprout: = 0.01. For this setting, the Figures 8 and 9. Stalk Cell has very small chance (probability of 0.01) of becoming an Active Cell in the next The simplest setting where the blood vessel grows in time step. a straight path and that does not sprout at all (Fig- ure 6) is fairly straightforward to analyze. The tip of • always-sprout: = 0.98. For this setting, the the blood vessel migrates one location towards north at Stalk Cell has 0.98 probability of becoming ac- each step, forming the body of the vessel along the pro- tive in the next step. This is quite an unrealistic cess. This setting, therefore, serves as a sanity check. setting; we present it only for didactic purposes. In the next setting, which is presented in Fig- ure 7, we keep the growth pattern the same We present results for four possible configurations: the (straight-growth) but increase the sprout possibility cross-product of the growth patterns and sprout pos- to 0.98 (always-sprout). In this setting, the blood sibilities. We first provide detailed results on a 3 ⇥ 3 vessel grows towards north as expected. Unlike the 1 Note that SW denotes the probability that an Active seldom-sprout case, however, a Stalk Cell at time Cell at the SW of an Empty location will move to this t = 1 became active at time t = 2. Empty location; hence SW denotes the possibility that an Active Cell at SW moves in the NE direction to occupy Next, we present results for the uniform-growth an Empty location. cases. In this setting, the blood vessel has uniform 92 .01 .01 .01 .04 .95 .04 .01 .01 .01 .04 .95 .04 AC .02 .98 .02 .02 .01 .02 AC .02 .98 .02 .02 .01 .02 .01 .01 .01 .01 .01 .01 .01 .01 .01 .02 .96 .02 .00 .00 .00 .01 .01 .01 .00 .00 .00 .01 .01 .01 SC .00 .00 .00 .02 .96 .02 SC .00 .00 .00 .02 .96 .02 .00 .98 .00 .02 .97 .02 .00 .98 .00 .02 .02 .02 t=1 t=2 t=1 t=2 Figure 6: straight-growth, seldom-sprout. Figure 7: straight-growth, always-sprout. The This is the most straightforward setting where the blood vessel grows towards north. A location that blood vessel grows one step at a time towards is a Stalk Cell at time t = 1 (Lx=1,y=0 ) becomes north. Active Cell at time t = 2. .01 .01 .01 .22 .31 .22 .01 .01 .01 .22 .31 .22 AC .34 .34 .34 .02 .02 .02 AC .34 .34 .34 .02 .02 .02 .01 .01 .01 .01 .01 .01 .01 .01 .01 .02 .96 .02 .00 .00 .00 .01 .01 .01 .00 .00 .00 .01 .01 .01 SC .00 .00 .00 .34 .33 .34 SC .00 .00 .00 .33 .33 .33 .00 .98 .00 .02 .97 .02 .00 .98 .00 .02 .02 .02 t=1 t=2 t=1 t=2 Figure 8: uniform-growth, seldom-sprout. The Figure 9: uniform-growth, always-sprout. The blood vessel has equal probability to grow in all blood vessel has equal probability to grow in all three directions. A Stalk Cell at time t will most three directions. A Stalk Cell will most likely likely remain as Stalk Cell at t + 1. become Active Cell in the next time slice. probability of growing towards NW, N, and NE. In becomes an Active Cell at t = 2. the seldom-sprout case (Figure 8), the Active Cell These toy experiments provide insights into how the at t = 0 turned into a Stalk Cell at time t = 1 and process typically works. Next, we present results for remained a Stalk Cell at time t = 2. The Active the 9 ⇥ 9 grid. Cell, unlike the straight-growth case, has equal probability of moving in all three directions. In the last time step, the middle of the top row has higher 4.2 Summary Results for 9 ⇥ 9 probability (.31) than the sides (.22) simply because the middle location can be reached from more locations Similar to the 3 ⇥ 3 grid, we provide evidence for t = 0 compared to the side locations. The always-sprout case and compute probabilities for the remaining eight case (Figure 9) is similar except a Stalk Cell at t = 1 time slices. In the initial configuration, the middle 93 .06 .06 .06 .10 .64 .10 .06 .05 .05 .17 .17 .18 .22 .75 .22 .17 .18 .17 .06 .06 .07 .06 .01 .06 .06 .06 .06 .17 .18 .17 .18 .13 .17 .18 .18 .18 .05 .05 .05 .05 .01 .05 .05 .05 .06 .17 .16 .18 .24 .88 .24 .17 .18 .18 .05 .05 .05 .04 .01 .05 .05 .05 .05 .16 .17 .16 .16 .09 .16 .17 .17 .16 Active .04 .04 .04 .04 .01 .04 .04 .04 .04 .16 .14 .16 .22 .91 .22 .15 .16 .15 Cell .04 .04 .03 .04 .01 .03 .04 .03 .04 .13 .14 .14 .13 .07 .14 .14 .14 .14 .03 .03 .03 .03 .01 .03 .03 .03 .03 .12 .11 .12 .17 .92 .16 .12 .12 .12 .02 .03 .02 .02 .01 .02 .02 .02 .02 .09 .09 .10 .10 .07 .09 .09 .10 .09 .02 .01 .01 .01 .01 .01 .01 .01 .01 .06 .06 .06 .06 .86 .05 .05 .06 .05 .22 .23 .22 .23 .23 .23 .24 .23 .23 .14 .14 .14 .15 .14 .14 .14 .14 .14 .23 .23 .23 .27 .86 .27 .23 .23 .22 .14 .14 .15 .18 .77 .18 .14 .15 .14 .23 .23 .23 .28 .87 .26 .24 .23 .22 .14 .15 .14 .14 .11 .14 .15 .14 .15 .23 .22 .23 .26 .88 .26 .23 .23 .22 .14 .13 .14 .20 .89 .20 .14 .15 .14 Stalk .21 .21 .22 .25 .89 .25 .22 .21 .21 .13 .14 .13 .13 .08 .13 .14 .14 .13 Cell .20 .20 .20 .23 .89 .23 .20 .20 .20 .12 .11 .12 .18 .91 .18 .12 .12 .12 .18 .17 .17 .20 .90 .19 .18 .17 .18 .10 .11 .11 .10 .07 .11 .11 .10 .10 .14 .14 .14 .15 .91 .15 .14 .14 .14 .08 .08 .08 .11 .92 .10 .08 .08 .08 .10 .10 .10 .10 .93 .10 .10 .10 .10 .06 .06 .06 .06 .07 .06 .05 .06 .06 straight-growth – seldom-sprout straight-growth – always-sprout Figure 10: AC and SC probabilities for the 9 ⇥ 9 grid at the last time slice (t=8) for straight-growth. Left: seldom-sprout. Right: always-sprout. On the right, it is apparent that the Stalk Cells and Active Cells alternate. of the bottom row is set as an Active Cell and the side of Figure 11), the blood vessel can be anywhere remaining locations are set as Empty. Due to space on the grid, except, as expected, the middle locations limitations, we present results for only the last time have higher probability. In the always-sprout case slice, t = 8. The straight-growth case is shown in (the right side of Figure 11), the Stalk Cells and Figure 10 and the uniform-growth case is shown in Active Cells alternate, as expected. Additionally, Figure 11. the probabilities for locations being a blood vessel (ei- ther Stalk to Active) are higher in the always-sprout In the straight-growth seldom-sprout case (the left case compared to the seldom-sprout case, again as side of Figure 10), we see a straight blood vessel for expected. the middle of the grid, where every cell of the blood vessel except the tip is a Stalk Cell and the tip is The results so far have been nothing surprising, but an Active Cell, as expected. In the always-sprout only confirming our expectations. The value of the case (the right side of Figure 10), the Stalk Cells and DBNs, however, lies at their capability to reason with Active Cells alternate, again as expected. spatial and temporal uncertainty as well as their po- tential for future directions. We discuss one of the In the uniform-growth seldom-sprout case (the left 94 .04 .06 .07 .08 .09 .08 .07 .06 .04 .17 .22 .24 .26 .26 .25 .24 .22 .17 .03 .04 .03 .04 .04 .04 .04 .04 .03 .16 .20 .21 .20 .20 .22 .21 .20 .16 .03 .04 .03 .04 .03 .04 .04 .03 .03 .21 .30 .34 .40 .40 .38 .34 .29 .21 .03 .04 .03 .03 .03 .03 .03 .03 .03 .15 .19 .18 .17 .17 .17 .19 .19 .15 Active .03 .03 .03 .03 .03 .03 .03 .03 .03 .22 .37 .51 .60 .63 .59 .50 .36 .21 Cell .02 .03 .03 .03 .02 .03 .03 .03 .02 .15 .17 .15 .13 .12 .13 .15 .17 .15 .02 .03 .02 .02 .02 .02 .02 .03 .02 .13 .16 .52 .71 .80 .70 .53 .16 .13 .02 .02 .02 .02 .02 .02 .02 .02 .02 .10 .13 .13 .09 .09 .09 .12 .13 .10 .01 .01 .01 .01 .01 .01 .02 .02 .02 .06 .06 .06 .05 .86 .06 .06 .06 .06 .17 .19 .21 .20 .21 .21 .20 .20 .16 .13 .16 .17 .17 .16 .17 .17 .16 .13 .17 .23 .23 .26 .27 .26 .25 .23 .18 .14 .18 .21 .22 .23 .21 .21 .18 .14 .18 .22 .25 .27 .27 .27 .26 .23 .17 .13 .16 .16 .15 .17 .16 .17 .16 .13 .17 .22 .26 .28 .29 .29 .26 .22 .17 .18 .26 .34 .40 .41 .39 .33 .25 .18 Stalk .17 .22 .26 .30 .31 .30 .26 .21 .17 .12 .15 .15 .14 .13 .14 .15 .15 .12 Cell .15 .21 .26 .32 .34 .31 .26 .20 .15 .12 .29 .47 .62 .66 .62 .46 .28 .11 .15 .16 .25 .32 .39 .33 .24 .17 .14 .11 .14 .12 .10 .09 .10 .12 .14 .11 .13 .15 .14 .40 .41 .41 .14 .14 .13 .09 .10 .10 .71 .71 .71 .10 .10 .08 .10 .10 .10 .10 .91 .10 .10 .10 .10 .06 .06 .06 .06 .07 .06 .06 .06 .05 uniform-growth – seldom-sprout uniform-growth – always-sprout Figure 11: AC and SC probabilities for the 9 ⇥ 9 grid at the last time slice (t=8) for uniform-growth. Left: seldom-sprout. Right: always-sprout. In both cases, the blood vessel can grow uniformly in each direction and the middle locations have higher probability of having a blood vessel simply because they can be reached from multiple locations. In the always-sprout case, Stalk Cells and Active Cells alternate. future directions here supplemented with some pre- periment at time t, what is the uncertainty over LT ? liminary results, and discuss more future directions in More practically: when is the earliest time we can stop Section 6. an experiment so that the uncertainty over the last time slice is below a pre-specified target ?. It is im- 4.3 Quantifying Uncertainty portant to note that when an experiment is stopped, the researchers dissect the tissue to analyze its proper- Given an initial condition, L0 , the tissue engineers are ties, such as vascularization, and hence the experiment interested in the final status of the tissue, LT , where cannot continue beyond that point. T denotes the final step of the experiment. Because Given an uncertainty measure, this question can be real-world experiment take a long time, mostly weeks, formulated rather straightforwardly using DBNs. Let they would like to be able to stop an experiment at U N C P LT |l0 , lt denote the uncertainty over the time t < T and still be able to reason about time T . predictions over the last time slice, given the initial Therefore, they are interested in the following ques- condition L0 = l0 and the status of the experiment at tion: given an initial condition L0 , if we stop the ex- 95 time t, Lt = lt . Then, we simply need to find Uncertainties for Two Growth Patterns 8.00 argmin U N C P LT |l0 , lt < (1) t 0 unless we stop the experiment. Therefore, we 4.00 need to take an expectation over all possible outcomes at time t: 3.00 2.00 X 1.00 argmin P Lt = lt |l0 U N C P LT |l0 , lt < t