New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Wound angiogenesis as a function of tissue oxygen tension: A mathematical model

Contributed by Avner Friedman, December 10, 2007 (received for review October 18, 2007)
Abstract
Wound healing represents a well orchestrated reparative response that is induced by injuries. Angiogenesis plays a central role in wound healing. In this work, we sought to develop the first mathematical model directed at addressing the role of tissue oxygen tension on cutaneous wound healing. Key components of the developed model include capillary tips, capillary sprouts, fibroblasts, inflammatory cells, chemoattractants, oxygen, and the extracellular matrix. The model consists of a system of nonlinear partial differential equations describing the interactions in space and time of these variables. The simulated results agree with the reported literature on the biology of wound healing. The proposed model represents a useful tool to analyze strategies for improved healing and generate a hypothesis for experimental testing.
Injury incites the most rapid and definable neovascularization (1). Healing wounds make excellent experimental models that can provide otherwise unavailable insights into the nature of angiogenesis (2). However, a series of events must occur for a vascular network to form in the wound site. The process begins with bleeding, when angiogenic blood products are spilled into the injured site. Platelets and proteins come into the wound to form a clot. Inflammatory cells, like neutrophils and macrophages, migrate into the wound to degrade necrotic tissue and bacteria. The coagulation proteins, inflammatory cells, and metabolic byproducts provide a flood of growth and angiogenic factors into the wound environment, which is conditioned to a moderate hypoxia, an accumulation of lactate, a slight lowering of pH, and an elevation of pCO_{2} (3, 4). These angiogenic factors recruit fibroblasts and endothelial cells, which accumulate in the wound region. The accumulation of cells increases the requirement for a blood supply, usually well beyond the capacity of the normal vessels that once fed the site. The fibroblasts respond to the metabolic need of new vessels by laying the extracellular matrix in which the capillary tubes form, whereas the endothelial cells line the capillary tubes that allow for blood flow into the repaired wound site. Once the metabolic needs of healing are met and the wound closes, the new vessels regress, leaving a vascularized fibrous scar containing a collagenous matrix that persists somewhat in proportion to the degree to which inflammation has occurred and the time that has been required for healing. In short, wound angiogenesis represents a response to coagulation, inflammation, and temporary metabolic need (3, 4).
Disruption of the vasculature at the wound site is responsible for woundsite hypoxia (5, 6). Although hypoxia is generally recognized as a physiological cue to induce angiogenesis (7–10), severe hypoxia cannot sustain the growth of functional blood vessels (2, 11–14). Interestingly, like hypoxia, hyperoxia also induces the expression of angiogenic factors and supports wound angiogenesis and healing (2, 15–17). Oxygen therapy is commonly used in the wound clinics to treat wound hypoxia (20–24). Extreme hyperoxia, on the other hand, can cause oxygen toxicity and derail tissue repair (16, 17, 19).
One of the first mathematical models on dermal wound angiogenesis was a threevariable system of partial differential equations that describes the evolution of the capillarytip endothelial cells, macrophagederived chemoattractants, and the new blood vessels during the tissue repair process (25). This model was later expanded to incorporate the roles of fibroblasts, oxygen, and extracellular matrix in the woundhealing process (26). More recently, a threevariable model, similar to that in ref. (25), was shown to be in agreement with clinical data on wound surface area for both normal and chronic wounds (27). Later on, mathematical models of angiogenic networks, such as through the induction of vascular networks by vascular endothelial growth factors (VEGFs) (28, 29), were developed by McDougall and coworkers (30, 31), based in part on the work of Anderson and Chaplain (32), in connection with chemotherapeutic strategies. In this work, our objective was to develop a mathematical model directed at addressing excisional cutaneous wound angiogenesis as a function of tissue oxygen tension. The present work builds on the mathematical model developed in ref. (26) by also incorporating the effects of macrophages and the use of oxygen treatment for the wound.
Results and Discussion
A system of nonlinear partial differential equations describing the interactions in space and time of select biological factors that are known to regulate wound angiogenesis was used to develop the proposed model. The process of wound healing involves complex interactions among a number of factors. In the early inflammation phase, the bloodborne inflammatory cells migrate into the wound site releasing a variety of angiogenic growth factors, such as VEGF. These factors interact chemotactically with resident cells such as fibroblasts. The fibroblasts produce collagen and other components of the ECM. Successful laying of the ECM facilitates the migration of fibroblasts, smooth muscle cells, marrowderived cells, and endothelial cells into the wound. Endothelial cells proliferate and align on the threedimensional matrix to form a network of capillaries. As soon as the anastomoses start to occur, blood begins to flow through the capillaries. In this study, for mathematical simplicity, the angiogenic response was estimated by the density of endothelial cells lining the capillary walls and was referred to as the capillary sprouts. Finally, oxygen is centrally essential for wound disinfection and to support the metabolism and signaling underlying wound repair (5, 33).
The proposed mathematical model is based on the following variables: n, density of capillary tips; b, density of capillary sprouts; w, concentration of oxygen; m, density of inflammatory cells (e.g., macrophages and neutrophils); a, concentration of chemoattractants (e.g., VEGF and other growth factors); f, density of fibroblasts; and ρ, density of ECM.
The primary aim of developing the proposed model was to produce simulation results that qualitatively agree with experimental findings reported in the current literature. Where applicable, parameters were obtained from experimental data reported in the literature. Other parameters were estimated based on related observations in the literature. The proposed model also includes a control term corresponding to therapeutic oxygen supplementation regimens. Hence, the model developed herein can be used to generate hypotheses on optimal treatment of wounds by varying the oxygen therapy regimen for wound care. Assumptions related to the variables considered for developing the proposed mathematical model are described below.
Capillary Tips.
It was assumed that the flux of the capillary tips is given by J _{n} = −(D _{n}∇n − χ _{n}(n,ρ)∇a) where D _{n} is a diffusion coefficient and χ _{n}(n,ρ) = χ _{n} (ρ/ρ_{0}) nH(n _{0} − n) is a chemotactic term. H(n) is the Heaviside function, where H(n) = 1 for n ≥ 0 and H(n) = 0, otherwise. The Heaviside function is introduced here to ensure that n does not exceed n _{0}. Note that there is no chemotactic response in the absence of extracellular matrix (34, 35). The new capillary tips are assumed to be produced from the tips themselves, as well as from sprouts, and that capillary tip production depends on the presence of the chemoattractant. Loss terms caused by sprout–tip and tip–tip anastomoses are also included (36, 37). From these assumptions, the equation for the capillary tips is:
Capillary Sprouts.
It was assumed that the capillary sprouts undergo diffusion, and, at the same time, are dragged along the flux of the tips (38–40). The velocity of the tips is v _{n} = (1/n) J _{n} , so that the sprouts move in the same direction as the tips with velocity v _{b} = g _{b}(n) v _{n} . Note that g _{b}(n) depends on n because the drag velocity will be very small if n is small; it is assumed that g _{b}(n) = An/n _{0}. A logistic term is also included, where it is assumed that the maximum capillary sprout density depends on the tissue oxygen tension G _{b}(w/w _{0}) and G _{b}(1) = 1. Because the material derivative of b is given by ∂b/∂t + ∇· (bv _{b}), the equation for the capillary sprouts is:
Oxygen.
It was assumed that the oxygen undergoes diffusion through the capillary vasculature (41, 42), and oxygen is lost through consumption by fibroblasts and macrophages (43, 44). An intervention term G(w/w _{0}) that represents the input of oxygen by using therapeutic approaches was also included. Thus, the equation for the concentration of oxygen is: An appropriate choice of G(w/w _{0}) to examine how the healing process can be improved is discussed later in this article.
Inflammatory Cells.
Inflammatory cells (i.e., macrophages and neutrophils) assist in wound healing by migrating to the injury site where they release angiogenic growth factors that promote the ingrowth of new blood vessels to the wound space. It was assumed that the flux of macrophages is given by J _{m} = − (D _{m}∇m + χ _{m}(m)∇w), where χ_{m}(m) = χ _{m}mH(m _{0} − m) is a chemotactic term. As in the case of the capillary tips, 0 ≤ m ≤ m _{0} is expected. It was also assumed that macrophages decay or lose functionality at a linear rate. Hence, the equation for inflammatory cells is:
Chemoattractant.
Many different families of inflammatory growth factors are known to regulate fibroblasts and endothelial cells. In this work, a single, generic growth factor whose properties simulate the combined effect of in vivo promoters of healing, such as VEGF, was considered. It was assumed that the concentration of chemoattractant undergoes random motion, and that it decreases through decay and through uptake by tips and sprouts. Finally, the concentration rate was assumed to increase in proportion to both the concentration of macrophages and the density of oxygen G _{a}(w/w _{0}), where G _{a}(1) = 1. Hence, the equation for the chemoattractant is:
Fibroblasts.
It was assumed that the flux of fibroblasts is given by J_{f} = − (D _{f}∇f − χ _{f}(f)∇a), where χ _{f}(f) = χ _{f}H(f _{0} − f). A logistic term, which depends on the tissue oxygen concentration G _{f}(w/w _{0}), was also included, where G _{f}(1) = 1. With the inclusion of a death term, the equation for the fibroblasts is:
ECM.
It was assumed that the buildup of ECM depends on fibroblast function (45, 46) in the same manner that capillary sprout motion depends on the motion of the capillary tips; that is, the velocity v _{ρ} of the ECM density is given by v _{ρ} = (g _{ρ}(f)f) J_{f} , where g _{ρ}(f) = B (f/f _{0}) is prescribed. A term that drives the ECM to its repaired value ρ_{0} and a diffusion term with a small coefficient D _{ρ} are also included. Because the material derivative of ρ is given by ∂ρ/∂t + ∇ · (ρv _{ρ}), the equation for the extracellular matrix is:
Boundary and Initial Conditions.
For simplicity, a radial crosssection of the wound was taken as an interval 0 ≤ x ≤ L, where L is the radius of the wound, and the partial differential equation system Eqs. 1 – 7 is solved in this interval. The endpoint x = L corresponds to the edge of the healthy tissue, whereas x = 0 represents the center of the radial wound. Along the wound edge, it was assumed that the densities of capillary sprouts, fibroblasts, and ECM remain constant at the level of the unwounded state, and the concentration of oxygen depends on G(w/w _{0}), which represents the concentration of oxygen in the blood in the case of systemic oxygen therapy (e.g., hyperbaric oxygen). The densities of the capillary tips and inflammatory cells are assumed to decrease at the boundary of the healthy tissue as the healing process progresses, and the flux of the chemoattractant across this boundary is negligible. Hence, the boundary conditions at x = L are: At x = 0, the following noflux boundary conditions are used: To complete the system, a set of initial conditions is needed at t = 0, the time when the inflammatory cells begin to infiltrate into the wound space. The initial conditions are taken to be: These initial conditions agree with boundary conditions 8 and 9 at x = 0 and x = L, and ε is selected such that all the densities are initially near zero away from the wound edge.
Table 1 lists the numerical values of the parameters for the PDE system (Eqs. 1 – 7 ) and the boundary and initial conditions 8 – 10 . Where applicable, parameters are taken from the cited literature and are based on firm experimental data. Whenever such data were not available, the order of magnitude of the parameters was estimated and choices were made that gave the best experimentally expected results regarding angiogenesis outcome. The parameters were then tuned to make sure that they are robust, namely, that the simulation results are essentially unaffected if the changed parameters are within the same order of magnitude.
Table 1 lists the parameters in both dimensional [centimetergramsecond (cgs) units] and dimensionless form. It is the latter values that are used in the simulations. The conversion from the dimensional to dimensionless units was done as follows: With the above substitutions, the nondimensionalized Eqs. 1 – 7 and boundary and initial conditions 8 – 10 become (with removal of the * for notational simplicity): where H(·) is the Heaviside function, G(·) is the input of hyperbaric oxygen, and
Scaling Coefficients ( n _{0}, b _{0}, w _{0}, m _{0}, a _{0}, f _{0}, ρ_{0}).
Estimates for the number of endothelial cells ranges from 1 × 10^{2} cells per mm^{3} (47) to 2 × 10^{3} cells per mm^{2} in a monolayer (48); which corresponds to 1 × 10^{5} cells per mm^{3}. If a cube of size 1 mm^{3} was saturated with endothelial cells with volume 10 μm × 10 μm × 1 μm, there would be 10^{7} cells. Because the cells weigh the same as water, the endothelial cells in 1 cm^{3} weigh between 10^{−5} and 10^{−2} g. Hence, b _{0} = 1 × 10^{−3} g·cm^{−3} and correspondingly n _{0} = 1 × 10^{−4} g·cm^{−3} were taken. Similarly, the average values of m _{0} = f _{0} = 1 × 10^{−3} g·cm^{−3} were chosen, because macrophages are 10 μm × 10 μm × 10 μm and about 65 cells can occupy a volume of 0.07 mm^{3} for a 1cmthick wound (49), fibroblasts are 100 μm × 100 μm × 10 μm and about 10^{4} cells are found in 1 cm^{3} of the dermis (50). Finally, w _{0} = 168 μM × 32 g·mol^{−1}, where 168 μM represents the mouse reference oxygen concentration, which is taken to be 80% of the human concentration (51), and 32 g·mol^{−1} is the molecular mass of oxygen, O_{2}.
Reaction Coefficients (λ_{i}s).
The parameter λ_{2} was expected to be of the same order of magnitude as reported (52), although the precise value was not known. λ_{2} was chosen to be smaller than reported (52). λ_{1} was expected to be larger than λ_{2}, and in the absence of experimental results, λ_{1} = 10 λ_{2} was chosen. The order of magnitude of the parameter λ_{8} was estimated by solving the equation dw̄dt = λ_{8} m _{0} w _{0} for λ_{8}, where dw̄dt = k × 32 g/mol, k is the oxygen consumption rate of macrophages in a cubic centimeter, and 32 g/mol is the molecular mass of O_{2}. The oxygen consumption rate (OCR) for macrophages is 1.87 nmol·min^{−1}·(10^{6} cells)^{−1} (51) and estimating that 10^{6} cells occupy a volume of 1 cm^{3}, k = 1.87 nmol·min^{−1}·cm^{−3}. Because the OCR of endothelial cells is about onefifteenth of the rate for macrophages (51) and the OCR of fibroblasts is similar to that of endothelial cells, λ_{8}* = λ_{7}*/15. The wound being inherently hypoxic, λ_{6}* < λ_{7}* under the assumptions that b ∼ 1, f ∼ 1, and m ∼ 0 for t large. Hence, λ_{6}* = λ_{7}*/2 was chosen. The parameter λ_{9}, the halflife of macrophages in a mouse, was estimated to be one third of that for macrophages in HIV studies (53). Uptake of the chemoattractant by endothelial cell tips was assumed to be similar to the uptake by the capillary sprouts, so λ_{11} = λ_{12}. By using estimates for the VEGFC halflife to be between 3.5 and 8 h (54, 55), λ_{13} = 4 × 10^{−5} s^{−1}. In the absence of parameter estimates from the literature, λ_{10}*, λ_{14}*, and λ_{16}* were chosen such that the healing of the mouse wound should be completed in ∼10 days (1).
The Functions (G_{i}s).
It was assumed that G _{b} (w/w _{0}) = 2w/w _{0} − 1 for w/w _{0} ≥ 1/2, 0 for 0 ≤ w/w _{0} < 1/2, and that G _{f} (w/w _{0}) = G _{b} (w/w _{0}). G _{a} (w/w _{0}) will depend on the oxygen input during healing. Because a slight deviation from normoxia to hypoxia and increased hyperoxia both induce growth factors such as VEGF (2, 16, 56, 57); where the internal minimum corresponds to normoxia.
Boundary Coefficients (γ _{i}S).
It was assumed that the capillary tips along the fixed boundary are at ≈10% of the original concentration after about 9 days, and that macrophages along the fixed boundary are at ≈10% of the original concentration after 6 days. Hence, γ_{1} = (ln 10)/9 day^{−1} and γ_{2} = (ln 10)/6 day^{−1}.
For the simulation in Fig. 1, time profiles of the average density of endothelial cells in the entire wound, normalized to their density in the nonwounded state, are shown for different levels of tissue oxygen tension as therapeutically managed for 10 days. The oxygen concentrations are normalized such that concentration = 1 represents a normoxic environment, <1 is considered a hypoxic setting, and >1 characterizes a hyperoxic scenario. The endothelial cell density increased as the tissue oxygen concentration changed from 0.25 in the hypoxic state to 2 in the hyperoxic state. The endothelial cell density decreased as the tissue oxygen tension increased from 2 (moderately hyperoxic) to 4 (extremely hyperoxic). Furthermore, cases of extreme hypoxia (when the external oxygen concentration ranges from 0.25 to 0.5) could barely recover half of the endothelial cell density that is found in the unwounded state.
Fig. 2 illustrates the difference between endothelial cell densities contrasting the normoxic case of one to administering oxygen therapy for 90 min per day for 10 days for different hyperoxic concentrations by using an oxygensupplementing device such as a hyperbaric oxygen chamber. A positive endothelial cell density difference would indicate that oxygen treatment generated an improved outcome compared with conditions of no oxygen therapy to the inherently ischemic wound. A negative difference would imply the adverse effect of oxygen treatment. For the most part, by increasing the normalized oxygen treatment levels to 2, 3, and 4, an increasingly higher density of endothelial cells was achieved. Simulations (Fig. 3) in the case when the base level represented the hyperoxic concentration of 2 showed no improvement in endothelial cell density when the oxygen level is increased to 3 or 4 for 90 min per day. This agrees with the result that high levels of oxygen are known to cause growth arrest and senescence of cells (19, 65). Fig. 4 suggests that, for the most part, using intermittent oxygen treatment for 90 min twice per day is slightly better than 90 min one time per day.
These results predicted by the developed model agree with the experimental literature on wound healing. First, the model suggests that an extremely hypoxic wound environment cannot sustain vascular growth (10–12, 14). Second, hyperoxia promotes wound angiogenesis and healing (2, 16, 18, 57). Furthermore, the model simulations indicate that there is an optimal level of hyperoxia beyond which the beneficial effects of hyperoxia may be reversed. Third, the use of intermittent oxygen treatment may stimulate angiogenic response (66).
It should be noted that the mathematical model is clearly a simplification of the complex process of wound healing. For example, it lumps all the growth factors into one species, introduces functions G _{i}(w) whose precise form is unknown, and describes the formation of new blood vessels in terms of endothelial cell densities. Nonetheless, the proposed mathematical model represents a starting point in the modeling of one aspect of wound angiogenesis, that is, the effect of tissue oxygen tension. Our effort was directed at focusing on multiple well established biological factors that are known to regulate wound angiogenesis. The list of factors considered herein clearly has room for expansion as the biological complexities continue to be experimentally unveiled. This work provides a mathematical model that enables simulations of how wound tissue oxygen tension may influence angiogenic outcomes. The simulated results agree with a variety of experimental findings. As such, the model represents a useful tool to analyze strategies for improved healing and generate hypotheses for biological testing. In addition, the current model lends itself to refinement to the next level of sophistication addressing additional factors implicated in the systems biology of wound healing.
Methods
In the process of wound healing a very large number of proteins and smaller molecules are involved. For example, there may be as many as 10^{5} cells in 1 mm^{3}. To simulate the complex interactions among these molecules, a standard approximation was used, whereby individual cells were replaced by their density. The density of endothelial cells at location x is defined as the number of cells in a small volume containing x divided by this volume, as this volume becomes “arbitrarily” small. The densities of fibroblasts, VEGF, etc., are similarly defined. These densities depend on space and time, and they undergo continuous changes as a result of their mutual complex interactions. These changes are expressed by a system of nonlinear partial differential equations.
The wound occupies a threedimensional space, but to simplify the computations, only its twodimensional projection was considered. It was further assumed that this projection is circular, so that the boundary of the wound is a circle and all the densities depend only on time and on the distance from the center of the wound. The problem was solved on a radial crosssection of the circular wound. The system of partial differential equations was solved by using the Matlab command pdepe, which uses the method in ref. 68 for the spatial discretization and the Matlab routine of ode15s for the time integration. At each time, t, the profile of the endothelial cell density was computed and then averaged over the wound area. These averages were displayed in Fig. 1 for different choices of the oxygen intervention term G(w/w _{0}).
Acknowledgments
This work was supported by National Science Foundation Award 0112050 and National Institutes of Health Grants RO1 HL 073087, GM 077185, and GM 069589.
Footnotes
 ^{§}To whom correspondence should be addressed. Email: afriedman{at}mbi.osu.edu

Author contributions: R.C.S., A.F., and C.K.S. designed research; R.C.S., A.F., and R.Z. performed research; R.C.S., A.F., and R.Z. analyzed data; and R.C.S., A.F., and C.K.S. wrote the paper.

The authors declare no conflict of interest.
 © 2008 by The National Academy of Sciences of the USA
References
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Gibson JJ ,
 Angeles AP ,
 Hunt TK

↵
 Oberringer M ,
 Jennewein M ,
 Motsch SE ,
 Pohlemann T ,
 Seekamp A
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 AlWaili NS ,
 Butler GJ

↵
 Gallagher KA ,
 Goldstein LJ ,
 Thom SR ,
 Velazquez OC

↵
 Mills C ,
 Bryson P
 ↵
 ↵
 ↵

↵
 Pettet GJ ,
 Chaplain MA ,
 McElwain DL ,
 Byrne HM

↵
 Byrne HM ,
 Chaplain MA ,
 Evans DL ,
 Hopkinson J
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Kearns DB ,
 Bonner PJ ,
 Smith DR ,
 Shimkets LJ
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Olaso E ,
 et al.
 ↵

↵
 Montminy J
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Ristimäki A ,
 Narko K ,
 Enholm B ,
 Joukov V ,
 Alitalo K
 ↵

↵
 Sheikh AY ,
 et al.

 Sherratt JA ,
 Murray JD

 Stokes CL ,
 Lauffenburger DA ,
 Williams SK

 Baldwin ME ,
 et al.

↵
 Morgan CJ ,
 Pledger WJ
 Cohen IK ,
 Diegelmann RF ,
 Lindblad WJ

↵
 Balasubramaniam V ,
 Mervis CF ,
 Maxey AM ,
 Markham NE ,
 Abman SH

 Rakusan K ,
 et al.
 ↵