A Mathematical Model for Wound Contraction and Angiogenesis

Cutaneous wounds, ulcers or burns, as a result of external damage, such as intensive solar exposure or accidents, occur in high numbers and may cause ever lasting traumas. In some cases, the wounds are very painful and impair an individual’s daily life in terms of health, and social interaction with his or her surroundings, but also in terms of mobility and ability to work or to perform leisure activities. Furthermore, the aesthetic feeling of patients is highly improved by the treatment of burns and other wounds since otherwise the wounds canmake a detrimental impact on the appearance of the patient’s skin. Besides aesthetic appearance, also the functioning of the damaged skin is different since infections are possibly more likely to occur near the (scar tissue of the) burn. Furthermore, the conception of temperature and pain by a patient is also altered by the burn, ulcer or cutaneous wound. Therefore, the treatment of these defects, which occurs in very high numbers, is important and hence it is crucial to improve treatments and to make treatments more efficient.


Introduction
Cutaneous wounds, ulcers or burns, as a result of external damage, such as intensive solar exposure or accidents, occur in high numbers and may cause ever lasting traumas.In some cases, the wounds are very painful and impair an individual's daily life in terms of health, and social interaction with his or her surroundings, but also in terms of mobility and ability to work or to perform leisure activities.Furthermore, the aesthetic feeling of patients is highly improved by the treatment of burns and other wounds since otherwise the wounds can make a detrimental impact on the appearance of the patient's skin.Besides aesthetic appearance, also the functioning of the damaged skin is different since infections are possibly more likely to occur near the (scar tissue of the) burn.Furthermore, the conception of temperature and pain by a patient is also altered by the burn, ulcer or cutaneous wound.Therefore, the treatment of these defects, which occurs in very high numbers, is important and hence it is crucial to improve treatments and to make treatments more efficient.
Nowadays many different treatments, such as, in more serious cases, skin transplantation, injection of cell cultures and the implant of artificial skin, are commonly carried out by plastical surgeons in order to relieve the patient and to minimize the aesthetic effects of scar tissue as much as possible.In order to be able to improve medical treatment of scars, it is important to understand the biological mechanism behind scar formation in terms of plastic deformation and a regenerated excess of collageneous matrix.More information about the treatment of scars and the biology behind the scar formation can, among others, be found in Bloemen (2011); van der Veer (2009).The deformations are caused by the contractile mechanism due to the pulling behavior of the (myo-)fibroblasts.Dallon (2008) experimentally found that the nature of contraction is predominantly determined by the following two mechanisms: -The degree and orientation of motility of fibroblasts.This motility mechanism of the fibroblasts is caused by the attachment and movement of the lamellipodia of the fibroblasts onto the collageneous fibres; -The attachment of the philapods of the less motile myofibroblasts onto the extra cellular matrix.
These mechanisms were simulated by Dallon (2010) by the use of a spring model.Since Dallon's two-dimensional model tracks the forces exerted by each individual philapod or lamellipodium of the individual (myo)fibroblast, the formalism becomes very expensive from a computational point of view when applied to burns of realistic dimensions with millions of fibroblasts and collageneous fibres.However, his observations and calculations are very useful for a thorough understanding of the implications of the fundamental processes and for the calibration of continuum models that are based on partial differential equations.
The severity of a wound is determined by the thickness of the damaged skin layer of the patient.A minor injury concerns the damage on the corneum, which is replaced with the patient's tissues relatively quickly.This type of injury is referred to as a first degree skin burn or wound.A more severe damage concerns the impairment of the epidermis.In this case reepithialialization has to repair the epidermis by migration and proliferation of keratinocytes.This damage is classified as a second degree burn or wound, but this type of wound will heal without many problems for healthy patients.Contraction does not take place and hence the burn or wound will recover entirely without any plastic deformations or significant scars.Third degree burns or wounds pose a more serious damage in the sense that also (part of) the dermis, and even the subcutis could be disrupted.Here the fibroblast rich dermis, consisting of the collageneous fibrous tissue has to be repaired.This takes place by a sequence of signaling processes to initiate proliferation and movement of fibroblasts and the restoration of the vascular network.As the fibroblasts produce an excess of collageneous fibres, on which they, and the myofibroblasts to a larger extent, exert contractile forces, wound contraction takes place.The process of wound contraction is a useful mechanism for a rapid minimization of exposure of the underlying tissues to hazardous external environments, when the wound is caused by mechanical damage.In particular, in skins of rabbits, percentages of up to 80 percent of the initial wound area have been reported in various experimental studies.
In humans, experimental studies evidence that wound contraction is much less significant, being in the order of 5-10 percent.However, in the case of healing of burns, this contraction phenomenon is undesirable, as it gives rise to significant deformations, which can be plastic and furthermore, the patient is left with an excess of collageneous fibres and unpleasantly looking scars.Remodeling is a very slow process which cannot remove all wrinkles as a result of the plastic deformations due to wound contraction.
In order to be able to improve surgical treatments, in terms of effectiveness and minimization of invasion into the patient, it is crucial to know which behavior of the (myo-)fibroblasts cause the plastic deformation of the skin.In the case of a third degree burn or wound, fibroblasts enter the wound area during the proliferative stages of the healing process of the dermis.Subsequently, the fibroblasts start to proliferate and to produce extracellular matrix, referred to as fibrous tissue.Furthermore, as a result of exposure to high strains, fibroblasts differentiate into myofibroblasts, which are considered as weak muscle cells.Myofibroblasts increase the measure of contraction of the wound area of the dermis.This differentiation process also takes place under the influence of a growth factor TF-β,whichisproducedbythe fibroblasts and myofibroblasts depending on the amount of extracellular matrix present.The growth factor also stimulates the production of fibroblasts and its regeneration of collagen.
Next to the enhancement of the production of several biological entities, the growth factor increases the chemotactic transport of the fibroblasts towards the wound area.To get more quantitative insight into the influences of the growth factors and other agents on the amount of wound contraction caused by the pulling forces of the (myo)fibroblasts, detailed mathematical models with parameter sensitivity analysis are indispensable.This sensitivity analysis can give insight into the quantification of the influence of all biological parameters involved 490 Tissue Regeneration -From Basic Biology to Clinical Application www.intechopen.com in this process on contraction and healing, as well as values such that healing of burns is optimized with respect to a minimal level of wound contraction and maximum healing speed.These results are helpful to surgeons to design innovative minimal invasive treatments that are applicable to patients.Furthermore, the mathematical analysis is helpful to make the treatment patient specific.
In this manuscript, we address two partial processes during the proliferative phase: wound contraction and angiogenesis.When the tissue is provided with enough oxygen and nutrients the process of wound closure starts.Cells in the epidermis, which mainly consist of keratinocytes, start regenerating the upperlayer of the wound.Usually the skin can not be replaced fully and some scars are left where the wound was located, due to excess of regenerated collagen.The second and third stage of wound healing do not take place at the same location in the wound.The former is located in the dermis, the latter is limited to the epidermis.The epidermis and the dermis consist of different type of cells and are separated by a so-called basal membrane, see also Figure 1.
In this manuscript, we first present a selection of the currently available mathematical models that seek to describe the biological processes of wound healing as well as possible.The healing process is very complex and many factors contribute to it, therefore simplifications have to be made.The main topics in this thesis are combining two models for angiogenesis and coupling models for the different stages of wound healing.Currently models exist for angiogenesis and dermal regeneration (fibroblasts and collagen) separately, however, only scarcely have there been attempts to couple the models.This is vital as the various stages of wound healing overlap and hence influence each other.The models that we consider in this study are all formulated in terms of partial differential equations (PDEs) that are based on conservation principles.Such coupled models could give more insights into how the process of wound

491
A Mathematical Model for Wound Contraction and Angiogenesis www.intechopen.comhealing evolves.These insights might lead to treatments that reduce healing time, e.g. the use of certain hormones to speed up the healing process.Also scars and other deformations due to incomplete healing might be prevented or reduced.This issue is especially crucial for the treatment of burns, where hypotrophic scars may result after injury.
A lot of simulations have been done in literature to give more insight on how these models behave.Also the dependence of the models on certain parameters is investigated in the present paper.Finally recommendations are made for future research in the topic of mathematically describing wound healing.Many studies are based on PDEs, for instance the work by Murray (2004); Olsen (1995); Maggelakis (2003); Gaffney (2002); Sherratt (1991); Adam (1999); Javierre (2009); Vermolen (2009;2010;2011); Schugart (2008); Xue (2009).Many other modeling studies are based on cell-based or Monte-Carlo methods, such as the cellular Potts models used by Glazier (1992); Merks (2009); Plank (2004) to mention a few.This class of model is lattice-based and minimizes a virtual energy functional.Recently, a semi-stochastic continuous cell-based model was formulated by Vermolen-Gefen (2011).
The present manuscript falls within the class of papers, in which one attempts to combine existing mathematical models for various subprocesses occuring during wound healing.Some recent work by this team are Vermolen (2010), andVermolen (2011), in which the most simplified models for wound contraction, closure and angiogenesis are coupled.The present model is based on a combination of some more advanced models for wound contraction and angiogenesis and gives a more sophisticated formalism for dermal regeneration.The paper is organized as follows.In Section 2, we present several models for processes of angiogenesis and wound contraction.Here, we also present some results of simulations.This is an innovation with respect to other studies.In Section 3, we consider the coupling of the models, including simulations.In Section 4, we describe the numerical methods.Finally, some conclusions are drawn in Section 5.
The character of the current manuscript is rather informal and descriptive about the mathematical concepts used in the present study.The level of abstract mathematics has been reduced tremendously.Further, we note that the simulations shown in this manuscript are still in a preliminary state.

Current models
The mathematical model for the proliferative stage of wound healing is usually separated in three distinct parts representing three stages of wound healing.These three stages are wound contraction, angiogenesis and wound closure.Note that the inflammation stage mentioned in Section 1 is not taken into account.This is due to the fact that inflammation only contains the damage and only after the inflammation stage is finished the real healing process starts.
In this chapter we present some of the currently available models on the above mentioned wound healing stages.In Section 2.1, we first present the model on wound contraction.Next in Section 2.2, a model for angiogenesis is presented.This model is a combination of two accepted models that describe a specific feature in the angiogenesis process.The two models take a very different approach on how to model the growth of new blood vessels in the wound.In Section 2.3, a study that attempts to combine models of dermal regeneration and angiogenesis is briefly discussed.In all the models that we deal with, we consider a bounded simply connected domain Ω ⊂ R 2 .The boundary is denoted by ∂Ω.

The wound contraction models
During the wound contraction stage fibroblasts (connective tissue cells) invade the wound site and contract the extracellular matrix (ECM).The contraction decreases the area of contact between the wound and its surroundings, thereby reducing the chance of contamination and infection.Furthermore this process is vital in assuring that new blood vessels can be formed in the wound during angiogenesis, since the fibroblasts invading the wound form the tissue in which the new capillaries can grow.The wound contraction stage is limited to the dermis, however, the contraction of the ECM also effects the tissue in the epidermis.
We use the model for contraction due to Javierre (2009), which deals with the presence of myofibroblasts.Myofibroblasts are a kind of weak muscle cells.They are nonmotile cells that differentiate from fibroblasts and transmit and amplify the traction forces generated by the fibroblasts, Vermolen (2009).Secondly, the model incorporates the effects of a growth factor that triggers wound contraction.
The equation concerning the fibroblast concentration u fib becomes where c ecm and u myo respectively denote the growth factor and myofibroblast concentration.The first term in the left-hand side corresponds to the total accumulation, the second term follows from passive convection due to deformation of the tissue.The third and fourth terms of the left-hand side account for random walk and chemotaxis (movement towards the gradient of the growth factor).In the right-hand side of the above equation, we respectively have the logistic proliferation term up to an equilibrium, in which proliferation is enhanced by the presence of the growth factor, the differentiation term to myofibroblasts under presence of the growth factor, back differiation from myofibroblasts and finally a term dealing with cell death.The cell death rate is denoted by d fib , the myofibroblast to fibroblast differentiation rate by k 2 and the fibroblast to myofibroblast differentiation rate by k 1 .Furt h ermor eλ 0 fib , C 1/2 and C k are known constants that monitor the growth factor's influence on the contraction process and K is a parameter that regulates the equilibrium concentration.
The production term, first on the right hand side, now also incorporates growth factor stimulated proliferation.The other three terms on the right hand side respectively account for differentiation to and from myofibroblasts and cell death.
The PDE for the myofibroblast concentration u myo is similar to equation (1).However, since myofibroblast are nonmotile cells, they will only move due to passive convection.The myofibroblast concentration thus obeys where ε ecm is a proportionality constant, further d myo and u 0 myo denote the myofibroblasts death rate and the myofibroblast equilibrium concentration respectively.The terms on the left-hand side account for accumulation and passive convection due to deformation of the tissue.The right-hand side contains proliferation of myofibroblasts under presence of the growth factor, differentiation of fibroblasts to myofibroblasts, back differentiation to fibroblasts, and (programmed) cell death (apoptosis).
Both fibroblasts and myofibroblasts contribute to the production of the ECM, furthermore the production is chemically enhanced by the growth factor, Vermolen (2009).The PDE for the ECM density ρ then is given by where λ ρ and d ρ are the ECM production and death rate respectively.The terms on the left-hand side describe accumulation and passive convection, and the right-hand side describes growth factor enhanced production of collagen, as well as decay of collagen.Furthermore λ 0 ρ and C ρ are known constants that monitor the growth factor's influence on the contraction process.Further, R ρ is a parameter that quantifies how the ECM production rate depends on the ECM density itself and η b and η d are proportionality constants.
The dynamics of the growth factor concentration c ecm are mainly determined by the fibroblasts and myofibroblasts as they produce the growth factor.Also the growth factor is motile, so it is subject to active convection.This leads to the following PDE for the growth factor concentration Here the second and third term on the left hand side account for passive convection and ordinary diffusion.The terms on the right hand side account for growth factor production and growth factor decay respectively.Furthermore D c is the growth factor diffusion coefficient, k c denotes the growth factor production rate and d c the natural decay rate.Also Γ is a parameter that quantifies how the growth factor production rate depends on the growth factor concentration itself and ζ is a proportionality constant.
The initial and boundary conditions for the fibroblast concentration and the ECM density remain the same as in the model due to Tranquillo.For the myofibroblast and growth factor concentration the following initial conditions are imposed for x ∈ Ω u .H e r ec 0 ecm denotes the growth factor equilibrium concentration.Furthermore the growth factor concentration satisfy a no flux boundary condition on all boundaries.For the myofibroblast concentration we can apply the same reasoning as for the ECM density ρ and hence no boundary conditions have to be imposed.The mechanical part of the wound contraction model is based on the linear viscoelastic equations, i.e. −∇•σ = f ext .(5) Here σ = σ ecm + σ cell , the stresstensor, accounts for the ECM related stress, σ ecm ,andthecell stress, σ cell .F u r t h e r m o r ef ext represents the external forces acting on the tissue.
In all the below discussed models the ECM related stress tensor σ ecm is given as Here the first two terms on the right hand side respresent the viscous effects and the last term the elastic effects.If we let u = u(x, t) denote the displacement of the ECM, then the strain tensor ǫ and the dilation θ in equation ( 6) are respectively given by Furthermore, in relation ( 6), I denotes the identity tensor and μ 1 , μ 2 , E and ν respectively represent the dynamic and kinematic viscosity, Young's modulus and Poisson's ratio.
The external forces acting on the tissue, f ext , are modelled similarly in all three models, i.e.
Here ρ = ρ(x, t) denotes the ECM density and s is the tethering elasticity coefficient.
At time t = 0 it is assumed that there is no displacement of the ECM, i.e. u(x,0)=0.A l s o we assume that u vanishes at the boundary far away from the wound, i.e. u(x, t)=0 for x ∈ ∂Ω.This can be justified by taking the computational domain Ω sufficiently large, so that the boundary effects can be ignored.
Since the myofibroblasts transmit and amplify the traction forces generated by the fibroblasts this is also visible in the cell traction term σ cell .For the model due to Olsen et al. this term is given by where ξ is a constant of proportionality and R τ quantifies how the cell traction depends on the ECM density.
In Javierre ( 2009) an extension of the model due to Olsen et al. is presented.Javierre et al.
propose the mechanical stress to act as a factor that effects the differentiation from fibroblasts to myofibroblasts.They introduce an estimation of the mechanical stimulus that depends on the dilation θ = ∇•u as

495
A Mathematical Model for Wound Contraction and Angiogenesis
and the first two terms on the right hand side account for the contractile stress generated internally by the myosin machinery and transmitted through the actin bundles, Javierre ( 2009).The third term on the right hand side establishes the contractile stress supported by the passive resistance of the cell.See also Figure 2 for a plot of p cell (θ).
Furthermore, in equation ( 11), the compression and traction strain limits are respectively denoted by θ 1 and θ 2 , p max respresents the maximal contractile force exerted by the actomyosin machinery and K max and K pas the volumetric stiffness moduli of the active and passive components of the cell.Also the parameter θ * can be computed from K act and p max as θ * = p max K act , Javierre (2009).To incorporate the effects of the mechanical stimulus on the fibroblast to myofibroblast differentiation an extra factor is found before the differentiation term (the second term on the right hand side) in equation ( 1).The fibroblast to myofibroblast differentiation term changes to p cell (θ) where τ d is a parameter that quantifies how the differentiation rate depends on the mechanical stimulus.Note that this term is also present in the PDE for the myofibroblast concentration and that thus the second term on the right hand side of equation ( 2) also changes to the above.
The mechanical stimulus also effects the the cell stresses and thus the cell traction term oe cell .
In Javierre (2009) it is assumed that oe cell depends linearly on p cell (θ) and thus that In order to give more insight in the process of wound contraction we did some simulations with the model due to Javierre (2009).We will describe the numerical techniques in Section 4. Further, the input data can be found in the appendix.As a computational domain we use the unit square, i.e.Ω = {x =(x, y)| 0 ≤ x, y ≤ 1}, and the initial wound is given by Ω w = x||x|≤ 1 2 .The results are given in Figures 3 to 5, where we show the solution two days after injury.The computations have been done using parameter values taken from Javierre (2009), see also the appendix.
In Figure 3 we show the fibroblast and myofibroblast concentration two days after injury.The myofibroblast concentration is highly concentrated around the wound edge, whereas the fibroblast have invaded the wound.Figure 4 shows the ECM density and the growth factor concentration two days after injury.We see that the ECM density is slightly elevated at the wound edge, which is to be expected since the wound is healing there.Furthermore the growth factor concentration has spread throughout the computational domain, but is still concentrated inside the wound.In Figure 5, the displacement pattern due to contractile forces exerted by fibroblasts and myofibroblasts is plotted.If we compare Figure 5 with Figure 3 we see that the ECM displacement is largest at places of high myofibroblast concentration.Here  the effect of the myofibroblasts, as they transmit and amplify the traction forces generated by the fibroblasts, can clearly be seen.
The normalized solutions in time at x = 0, furthest in the wound, are shown in Figure 6.This gives an indication of the development of the solutions in time.We see that the fibroblast concentration increases towards it equilibrium, whereas the myofibroblast concentration first rises and then falls again.This is due to the fact that, as the wound heals, the myofibroblasts either differentiate back to fibroblasts again or undergo apoptosis and hence eventually will disappear completely.Also the growth factor concentration eventually goes to zero as the wound is healed.The ECM density slowly rises in time and will eventually reach its equilibrium, allthough it can take some time before the ECM is completely restored.We finally  remark that all the simulations in this subsection were obtained without any coupling with the model for angiogenesis.

Angiogenesis
In this section, we present the model for angiogenesis.These two models each take a very different approach on how to model this process.Where the model due to Maggelakis (2003) focusses on the relation between a lack of oxygen and capillary growth, the model due to Gaffney (2002), attempts to model the migration of endothelial cells into the wound.Both models adress an important aspect of angiogenesis, but on the other hand both also miss an important aspect.
In this section we present a novel angiogenesis model based on the models of Maggelakis and Gaffney et al.In this model we attempt to unite the two approaches and thus create a model that is suitable for dealing with both aspects of angiogenesis.This gives a new model for angiogenesis.By combining these two aspects in one model we hope to create a more accurate model of the process of angiogenesis.As a basis we take the model of Gaffney et al. and extend it with a negative feedback mechanism for oxygen.This assures that both the endothelial cell migration and the oxygen shortage aspects are covered in this novel angiogenesis model.where the left-hand sides account for accumulation and transport of tips and endothelial cells by means of a biased random walk process.Further, the functions f and g are given by

The model due to Gaffney et al. is given by
These functions were taken from Gaffney (2002).The first term of the function f models tip-branching, as each tip has a likelihood to bifurcate.The second term of the function f accounts for tip-tip anastomosis, which is the process that two tips have a probability to join and hence form a loop, by which a tip is no longer a tip, and hence the number of tips decreases then.The third term of the function f stands for tip-sprout anastomosis, which is the process that a tip may converge into a branch, which closes the network as well and thereby also decreasing the number of tips.For the function g, we distinguish the following terms: The first one contains ordinary logistic proliferation under standard conditions, the second term accounts for an increased logistic growth due to the presence of tips.Finally, the third term accounts for the extent of anastomosis.To incorporate the effects of the macrophage derived growth factor (MDGF) concentration c md on the growth of capillary tips and the proliferation of endothelial cells we assume that λ 2 and λ 6 are functions of the macrophage derived growth factor c md , which is released as a result of a shortage of oxygen.These functions must be zero if there is no MDGF present and rise as the MDGF concentration rises.Furthermore the effects of additional MDGF must be lower if there is already a lot of MDGF present.Therefore we let where τ tip and τ end qualify how respectively λ 2 and λ 6 , accounting for tip regeneration and logistic proliferation of endothelial cells respectively, depend on the MDGF concentration.The functions λ 2 (c md ) and λ 6 (c md ) are given in Figure 7 for 0 ≤ c md ≤ 1.
This settles the trigger mechanism that the MDGF concentration fullfills in the growth of new capillaries.To get a good overview of what the model looks like we give the PDEs that drive the model, i.e.
∂u oxy ∂t = D oxy ∆u oxy − λ oxy u oxy + λ 13 u tip + λ oxy u θ u end where we explain the terms occuring in the above PDEs.-the first equation: the first, second and third term in the right-hand side, respectively, stand for diffusion of oxygen (u oxy ), consumption of oxygen, and supply of oxygen due to the blood flow in the tips and capillaries (via diffusion through the capillary walls); -the second equation: the first, second and third term in the right-hand side, respectively, model diffusion of macrophage derived growth factor (VEGF), natural decay of VEGF, and secretion by the macrophages as a result of a shortage of oxygen; The last two equations were described earlier.
Further, the functions f and g are given by To illustrate how this new model behaves we show some results in Figures 8 to 10.The input-data can be found in the appendix.As a computational domain we use the unit square, i.e.Ω = {x =(x, y)| 0 ≤ x, y ≤ 1}, and the initial wound is given by Ω w = x||x|≤ 1 2 .
In Figure 8 we see the oxygen and MDGF concentration seven days after injury.The oxygen tension decreases as one moves into the center of the wound.The relation between a lack of oxygen and the production of the macrophage derived growth factor can clearly be seen.In areas where the oxygen concentration is low, i.e. inside the wound, the MDGF concentration is at its peak and vice versa.Hence, the MDGF concentration is maximal at the center of the wound.The time for ingress of macrophages has been assumed to be negligible.Also in Figure 10 we see that as the oxygen concentration rises the MDGF concentration drops at approximately the same rate.
Furthermore in Figure 9 we see a front of capillary tips, which slowly moves towards the center of the wound.Also the endothelial cell density is slighty elevated at the wound edge, since this is where the new capillaries are formed.This can also be seen in Figure 10, where the endothelial cell density first peaks and then drops again towards an equilibrium.The capillary tip concentration also peaks, but then drops to zero again.This is to be expected since eventually all capillary tips join together in the newly formed capillary network.The local maximum of the oxygen concentration at t ≈ 18 days, follows from the peak in capillary tips.As the shortage of oxygen decreases, the concentration of macrophage derived growth factors decreases down to zero.

The coupled model
For all the stages during the proliferative phase of wound healing mathematical models have been presented that attempt to describe the processes involved in these stages.These models only focus on one wound healing stage particularly and do not take into account the interactions between these stages.One example of interaction is the need of oxygen for the growth of cells and thus the dependency of the profileration terms on the oxygen concentration.But there are several more interactions that we can think of, including some that will not be taken into account here either.
In Vermolen (2010) and Vermolen (2011) an attempt is made to combine models of the three wound healing stages into one mathematical model.In this chapter we take the novel model on angiogenesis, presented the previous subsection, and make an attempt to couple it with the wound contraction model due to Javierre (2009), given in Section 2.1.We will investigate both the influence of the oxygen concentration on the wound contraction mechanism and the influence of the fibroblast concentration on capillary growth.Furthermore the displacement of the extracellular matrix (ECM) will also effect the angiogenesis process, this will also be taken into account.For now we will not consider wound closure.
In the wound contraction stage several processes require energy to work.The growth of new tissue consumes energy as cells divide and form new cells, i.e. mitosis.Also the active movement of the fibroblasts requires energy since the cells crawl over each other.This energy needed for wound contraction comes primarly from the food that an individual consumes, however, oxygen is needed to put the energy to use, i.e. the cells require oxygen to get the most out of the energy.This means that the oxygen concentration is good indicator for the amount of energy that is available.

www.intechopen.com
If we look at the partial differential equation (PDE) for the fibroblast concentration in the Javierre model, the third term on the left hand side denotes the active mobility due to random walk and the first term on the right hand side represents fibroblast proliferation.To incorporate the effects of the oxygen level on these two processes we replace these terms with the following, for the mobility and for the proliferation of fibroblasts.The extra term in front, u oxy τ oxy +u oxy , assures that if there is no oxygen present, i.e. no energy can be consumed, the two processes are stopped.When the tissue is saturated with oxygen the processes continue at their normal rate (the term then approaches one).Note that due to the PDE for the oxygen concentration the oxygen level will never rise to dangerous levels (the concentration moves towards an equilibrium), i.e. oxygen poisoning will never take place.
The proliferation of myofibroblasts also consumes energy and thus the same dependency on the oxygen concentration can be found in the PDE for the myofibroblasts.This concludes the influence of the oxygen level on wound contraction covered in this model.
The other way around wound contraction also effects the growth of new capillaries.First the displacement of the ECM causes passive convection of the variables of the angiogenesis model.This means that in each of the four equations an extra term is found that describes this passive convection, i.e.
where u i denote the four variables of the novel angiogenesis model.
Furthermore the new capillaries need tissue to grow in.At first all tissue in the wound has been destroyed by the injury and so the capillaries can not grow.Gradually fibroblasts invade the wound and new tissue is formed.Only then can the recovery of the capillary network start.This is why it is reasonable to let the growth of capillaries depend on the fibroblast concentration.
The growth terms in ( 17) and ( 18 respectively.Similarly to how the oxygen concentration is coupled with the fibroblast PDE we couple the fibroblasts to the capillary tip and endothelial cell density PDEs.Therefor we introduce the factor u fib /u 0 fib τ fib + u fib /u 0 fib in front of both production terms.This results in no capillary growth if there are no fibroblasts present, i.e. there is no appropriate dermal tissue for them to grow in.As the fibroblast concentration rises the rate of capillary growth also rises towards its normal rate, since the term approaches one.
This concludes the coupling between the wound contraction model due to Javierre and the novel angiogenesis model covered in this thesis.Of course there are several other interaction between them.One can for instance think about the differentation from fibroblasts to myofibroblasts and back, which surely also consumes energy.Alternatively, we could think of also linking the chemotaxis term to the oxygen content.This could be a subject for further study.
Next we present some results of the computation done on the coupled model.The values of the parameters used can be found in appendix.We vary the diffusion speeds of the fibroblasts and the oxygen, since we consider these to be of great importance (the problem seems to be diffusion dominated) and we would like to study the effect of these parameters on the solution.
In Figures 11 to 13 the results of the simulations can be found.The figures show the normalized solutions in time furthest in the wound, i.e. at x = 0.The difference in diffusion speeds can clearly be seen in the graphs.
We see that with slow fibroblast diffusion the growth of capillaries start off at a later point in time.This is to be expected since the fibroblasts provide the tissue in which the capillaries can The growth of the capillaries does follow a similar course.First the capillary tips find their way to the center of the wound and after that the capillary network is restored, which decreases the number of tips.As the number of tips decreases down to zero, the oxygen content decreases for a short while due to the lack of this tip-source.Later, the oxygen tension increases to its undamaged equilibrium.
Further, having slower oxygen diffusion, the oxygen concentration depends more on the growth of new capillaries.In Figure 13 we see that the oxygen concentration rises far slower than 11 and 12.The oxygen concentration grows due to capillaries transporting oxygen to the wound side and not primarly due to diffusion.
We see that the varying the diffusion speeds has a major effect on the solutions.This confirms our idea that the problem is diffusion dominated.Although the solutions follow similar patterns in all cases, the differences can also be seen clearly.
Furthermore the coupling between the two models can also be seen.Especially with low diffusion speeds we see that the fibroblast concentration clearly depends on the oxygen concentration.Also the growth of new capillaries starts off later than in the uncoupled model, see Section 2.2.This can also be explained by the coupling, since the growth of capillaries now  depends on the fibroblast concentration.So the growth can not begin untill there is enough tissue available, i.e. the fibroblast concentration is non-zero.

Numerical methods
We use a Galerkin finite-element method with linear triangular elements to solve the resulting system of PDEs.The reaction-transport equations are solved using a second-order Trapezoidal Time Numerical Integration Method.The nonlinear problem at each time step is solved using a Newton's method, which gives quadratic convergence.For the numerical integration of the spatial integrals in the weak formulation, we use the second order accurate Newton-Cotes integration rule.The numerical solution will thereby have a quadratic accuracy in both element size and time step.For an efficient implementation, the reactive-transport equations are solved as a fully coupled problem in which the degrees of freedom are numbered in a nodewise fashion.The mechanical parameters were taken at the current time step, since at each time-step the mechanical problem is solved before the biological reactive-transport equations are solved numerically.
The visco-elastic equations are solved using a Galerkin finite-element method with linear triangular elements.The system is solved as a fully coupled problem with a nodewise arrangement of the degrees of freedom, being the vertical and horizontal displacement.The concentrations and cellular densities from the previous time step are used in the mechanical balance equations from visco-elasticity.For the time-integration, a Trapezoidal Numerical Integration Method is used as well to have second order accuracy in time as well.The linear elements give a second-order accuracy for the element size.

Discussion
In this paper, we described and coupled some mathematical models for wound contraction and angiogenesis as subprocesses for dermal regeneration post-wounding.It can be seen in Figure 10 that in the course of healing, the endothelial cell density, which is the most important measure for the vascular network density, first increases up to value above the equilibrium value, which corresponds to the fully healed state.This increase coincides with the temporary high capillary tip density at these times.After this intermediate stage, the capillary tip density decreases down to zero, and the endothelial cell density decreases to the equilibrium as time proceeds.This qualitatively agrees with the rabbit ear chamber images from the experiments that were carried out by Komori (2005), among others.
Furthermore, the numerical solutions converge to the undamaged state at the longest scales.This is consistent with clinical situations if the initial wound is very small.However, scar formation due to an excessive regeneration of collagen is not taken into account in the present model.Hence, the present model is incapable of predicting the occurrance of ulcers or hypertrophic scars that can result for large wounds or burns.This issue will be dealt with in future, as well as the extension of the present model by the incorporation of the epidermal layer which starts healing as a result of triggering of keratinocytes by the signaling agents that are secreted by the fibroblasts in the dermis.We are also interested in the development of the basal membrane, which separates the dermis from epidermis, in the course of time.
We also note that the model is entirely based on the continuum hypothesis, and that all densities and concentrations should be considered as averaged quantities over a representative volume element.At the smaller scale, one has to use a more cellular approach, such as the cellular automata (cellular Potts) models or (semi-) stochastic approaches.A translation between these classes of models is highly desirable and hardly found in literature.
Finally, we note that modeling of wound healing is very complicated, due to the enormous complexity of the healing process.The present model can already help in obtaining insight into the relation between several parameters, such as the influence of various parameters on the speed and quality of wound healing.In future studies, we intend to include the effects of several treatment of wounds, such as pressure or drug treatments, on the kinetics of wound healing.The models can then be used as a tool for medical doctors to evaluate the influence of certain treatments on wound healing, as well as to determine the circumstances for optimal healing with respect to minimization of hypertrophic scars in the case of burns.Here, we also remark that all data are patient-dependent and hence all parameters in the model are exposed to a stochastic nature.In order to deal with the issues of uncertainties, we intend to employ stochastic finite element methods in near future.

Conclusions
We presented an advanced model for the coupling of angiogenesis and wound contraction.The model gives a more complete picture than most of the presently existing formalisms in literature.However, still some room for improvement remains.The current model could also be extended with the healing of the epidermis located on top of the dermis.This will be done in a future study.The first important innovations in the present study are the formulation of a new angiogenesis model, which consists of a driving force due to lack of oxygen (Maggelakis) and a balance of both the capillary tips and capillaries by tracking these quantities as individual parameters.
The second important innovation concerns the combination of the wound contraction model due to Javierre (2009) with the newly developed angiogenesis model.As an important factor of communication between both processes, we use the oxygen tension (on fibroblast mobility and proliferation), the presence of macrophage derived growth factor influencing the logistic proliferation rate of endothelial cells and the fibroblast density on the production of endothelial cells.Of course, these coupling terms have been assumed from educated guesses.
An experimental validation could be decisive on the validity of the assumptions used in the present model.The calculations in the present manuscript are preliminary and a parameter sensitivity analysis is to be carried out.
We also note that there is a large uncertainty for the values of all the parameters used.We intend to use stochastic methods to better deal with the occurrance of patient dependencies and other uncertainties in the data.This will be tackled in future studies.The present paper was more descriptive, rather than mathematical.

Appendix: Parameter values
In this Appendix, we show the data that we used for the calculations.The first two tables show the parameter values for the skin regeneration models ((myo)fibroblasts and collagen).

Fig. 1 .
Fig. 1.A schematic of the events during wound healing.The dermis and epidermis are illustrated.The picture was taken with permission from http://www.bioscience.org/2006/v11/af/1843/figures.htm

493A
Mathematical Model for Wound Contraction and Angiogenesis www.intechopen.com

494
Tissue Regeneration -From Basic Biology to Clinical Application www.intechopen.comA Mathematical Model for Wound Contraction and Angiogenesis 7

) 499 A
Mathematical Model for Wound Contraction and Angiogenesis www.intechopen.com

507A
Mathematical Model for Wound Contraction and Angiogenesis www.intechopen.com

508
Tissue Regeneration -From Basic Biology to Clinical Application www.intechopen.comA Mathematical Model for Wound Contraction and Angiogenesis 21

Table 1 .
The third table gives the data of the visco-elastic model for wound contraction.The fourth table gives the input for the angiogenesis model.Further, the values of the coupling parameters are scattered throughout the tables as well.Parameters related to the fibroblast and myofibroblast equations.509 A Mathematical Model for Wound Contraction and Angiogenesis www.intechopen.com

Table 2 .
Parameters related to the collagen and growth factor equations.

Table 3 .
Parameters related to the mechanical behaviour of cells and the ECM.

Table 4 .
Parameters related to angiogenesis.