Abstract

We deal with a mathematical model of atherosclerosis plaque formation, which describes the early formation of atherosclerotic lesions. The model assumes that the inflammatory process starts with the penetration of low-density lipoproteins cholesterol in the intima, and that penetration will occur in the area of lower shear stress. Using a system of reaction-diffusion equations, we first provide a one-dimensional model of lesion growth. Then we perform numerical simulations on an idealized two-dimensional geometry of the carotid artery bifurcation before and after the formation of the atherosclerotic plaque. For that purpose, we consider the blood as an incompressible non-Newtonian fluid with shear-thinning viscosity. We also present a study of the wall shear stress and blood velocity behavior in a geometry with one plaque and also with two plaques in different positions.

1. Introduction

Cardiovascular diseases (CVD), like coronary artery disease and stroke, are among the largest causes of death and disability in the industrialized world. Each year CVD causes over 4 million deaths in Europe (47%) (data from European CVD statistics of 2012) [1]. For example, in Cape Verde, according to data of 2010, from statistical report of Health Ministry [2], CVD are at the top of the list of mortality causes (25.3%). Therefore, a better understanding of the mechanisms underlying atherosclerosis is essential for the development of new therapeutic approaches, and the progress in mathematical modeling and numerical simulations of the associated phenomena have a key role in this framework.

Coronary heart disease and stroke are clinical symptoms of atherosclerosis which is an inflammatory disease, in which high plasma concentrations of cholesterol, in particular those of low-density lipoproteins (LDL) cholesterol, are among the principal risk factors [3]. More princely, atherosclerosis is a pathological process promoted by an inflammation of the inner arterial wall (intima), initiated by an excess of LDL in the blood. The LDL particles as well as the high-density lipoproteins (HDL) become oxidized by ongoing chemical reactions within the body and form ox-LDL, contributing to the atherosclerotic process when ox-LDL concentration exceeds a threshold. At this stage, endothelial cells activate the immune system (monocytes and T cells) to deal with the problem. Once in the intima, the incoming immune cells instantaneously differentiate into active macrophages, which absorb ox-LDL by phagocytosis. This reaction transforms macrophages into foam cells that should be removed by the immune system, yielding the secretion of a proinflammatory signal contributing to the recruitment of new monocytes, and consequently this starts a chronic inflammatory reaction (autoamplification phenomenon). The inflammation process involves the proliferation and the migration of smooth muscle cells to create a fibrous cap over the lipid deposit, isolating it from the blood flow. The fibrous cap changes the geometry of the vessel and modifies the blood flow. The occurrence of a heart attack or stroke is thus attributed to the degradation and rupture of the cap, the formation of a blood clot in the lumen, and the subsequent obstruction of the artery.

Mathematical modeling of these processes in time and space leads to complex systems of flow, transport, chemical reactions, interactions of fluid and elastic structures, movement of cells, coagulation and growth processes, and the additional complex dynamics of the vessel walls.

The starting point of this work is to reproduce the results obtained in [4], using a system of reaction-diffusion equations to describe the formation of the atherosclerotic plaque and numerical simulations to understand the behavior of plaque growth. Moreover, the main goals are (i)to perform numerical simulations on an idealized two-dimensional geometry of the carotid artery bifurcation before and after the formation of the atherosclerotic plaque, considering blood as an incompressible non-Newtonian fluid with shear-thinning viscosity; (ii)to present a study of potential regions of risks to the formation of new plaques, based on blood velocity behavior and WSS values.

This paper is organized into 7 sections. In the second one we explain all biological processes of atherosclerosis, separated in four main steps. In section three we present a mathematical model of atherosclerosis plaque formation, consisting of a system of 2D reaction-diffusion equations. In the fourth section we model the blood flow as an incompressible non-Newtonian fluid with shear-thinning viscosity in a two-dimensional domain. In section five we present numerical simulations to understand better the atherosclerotic plaque growth and the behavior of fluid flow in an idealized carotid bifurcation before and after the atherosclerotic plaque. In section six we consider a deformed idealized geometry and study the possible zones of risk of the formation of new plaques, based on blood velocity behavior and WSS values. Finally, section seven is devoted to conclusions and plans for future work.

2. The Biological Process

Before explaining the biological process of atherosclerosis we start with a brief presentation of the protagonists of this process, namely, the blood components and the blood vessel wall structure.

Blood is a suspension of particles (blood formed elements, representing approximately 45% by volume of the normal human blood) in an aqueous ionic protein solution called plasma; see Figure 1.

The most important blood particles are (i)red blood cells (RBCs) or erythrocytes, responsible for the exchange of oxygen and carbon dioxide with the cells; (ii)white blood cells (WBCs) or leukocytes, which play a vital role in the immune system fighting infection; (iii)platelets or thrombocytes, the main particles responsible for blood coagulation.

The vessel wall is divided into three main layers: tunica intima, tunica media, and tunica adventitia, represented in Figure 2.

The intima layer is made essentially of the endothelial cells (EC) and a few smooth muscle cells (in human arteries). The internal elastic membrane separates the intima from the media that comprises numerous layers of smooth muscle cells (SMC). An external elastic lamina separates the media from the adventitia layer.

The resident cells of the arterial wall (EC and SMC) in concert with cells emigrated from the blood (in particular T lymphocytes and monocytes) and their secretory products (chemokines, cytokines, enzymes), through ample cross talk and signaling, contribute to the initiation, evolution, and fate of the atherosclerotic plaque.

2.1. Atherosclerosis Stages

According to Ross [3], atherosclerosis is an inflammatory disease, in which high concentrations of low-density lipoprotein (LDL) cholesterol in the blood is one of the principal risk factors.

Atherosclerotic lesions occur mainly in large and medium size arteries such as the abdominal aorta, the coronary arteries, or the carotid bifurcation.

The atherosclerosis process can be divided into four main stages.

(i) Inflammation and Initiation. The initiation of atherosclerosis begins with the penetration of LDL into the intima of the blood vessel, where they are oxidized (ox-LDL), [3]. Oxidized LDL is considered as a dangerous agent, hence an anti-inflammatory reaction is launched: monocytes adhere to the endothelium; then they penetrate into the intima, where they are transformed into active macrophages (Figure 3).

(ii) Lipid Accumulation. Active macrophages absorb ox-LDL in the intima by the phagocytosis process. This reaction transforms macrophages into foam cells (lipid-laden cells) which in turn have to be removed by the immune system. Hence, they set up a chronic inflammatory reaction (autoamplification phenomenon) by secreting proinflammatory cytokines that promote the recruitment of new monocytes and the production of new proinflammatory cytokines (Figure 4).

(iii) Growth and Cap Formation. The inflammation process involves the proliferation (growth or production of cells by multiplication of their parts) and the migration of smooth muscle cells (SMC) to create a fibrous cap over the lipid deposit.

The fibrous cap changes the geometry of the vessel and modifies the blood flow (Figure 5).

(iv) Plaque Rupture. Due to several reasons, such as local blood flow behavior, SMC apoptosis, and endothelial cells apoptosis, the plaque can be degraded and can rupture (Figure 6), leading to the formation of a blood clot in the lumen and the subsequent obstruction of the artery, and finally to a heart attack or stroke.

3. Mathematical Modeling of Atherosclerosis Plaque Formation

Stages (i) and (ii) of atherosclerotic plaque formation can be mathematically described using reaction-diffusion partial differential equations as in [4].

First of all we define the geometry of the intima layer.

The intima layer will be considered as a 2D band of height , with and (see Figure 7). The coordinate denotes the height of the intima layer with values between and . The value corresponds to the interface with blood and to the interface with tunica media.

To simplify, we consider that the intima's lower boundary is fixed.

3.1. Oxidized LDL ()

The concentration of oxidized LDL (ox-LDL) can be modeled by the following equations:

The second left-hand side term represents the lesion growth, stating that the molecules of ox-LDL are transported with the tissue displacement having velocity . On the right-hand side, the first term is the diffusion term, where is the diffusion coefficient.

Here we are assuming the law of mass action, which states that the rate of an elementary reaction (a reaction that proceeds through only one transition state) is proportional to the product of the concentrations of the participating molecules. Thus, the second term on the right-hand side represents the reaction between ox-LDL () and macrophages (), where is the constant of proportionality.

The function is the permeability of the blood vessel, and is the LDL-cholesterol concentration. Here we are assuming that depends on the wall shear stress (WSS), which is the mechanical force imposed on the endothelium by the flowing blood. As we know, the WSS favours the penetration of both LDL and macrophages.

3.2. Macrophage ()

For the macrophages density we use the following equations:

Here we are assuming that the incoming monocytes immediately differentiate into macrophages and the recruitment of new monocytes depends on a general proinflammatory signal which gathers both chemokines and cytokines. This signal acts through the function , which is defined, in order to impose a limit in the macrophages recruitment, as

3.3. Foam Cells ()

Assuming that foam cells () cannot diffuse, the description of the foam cells concentration is given by

Under a local incompressibility assumption, when foam cells are created the intima volume is locally increasing.

3.4. Signal ()

For the signal () emission we consider the equations

The signal can diffuse, which is taken into account by the first term on the right-hand side, where is the diffusion coefficient. The second term denotes the natural death of the cells, and is a degradation rate.

The starting point of the signal emission is assumed to be a too high oxidized LDL concentration. Thus, corresponds to a given ox-LDL quantity, and is an activation rate.

3.5. New Intima's Height—

Denoting by the biomass which constitutes the rest of the intimal medium (extracellular matrix, smooth cells, and othesr) and assuming that does not contribute to the inflammatory process, we consider the equation

We also assume a local matter incompressibility stating that there exists a constant such that

We introduce the vector field which is the velocity of lesion growth and assume that the matter growth induced by foam cells accumulation is in the direction. In other words, by assumption .

Adding (4) and (6) and using the velocity, , and the local incompressibility assumptions (7), we get where, for any concentration ,

4. Mathematical Model for the Blood Flow

The blood flow is governed by the conservation of linear momentum and the conservation of mass. Assuming that the fluid is incompressible, we have where is the blood density, represents the velocity field, and is the pressure.

Here, we consider the blood modeled as a generalized Newtonian fluid: where represents the Cauchy stress tensor, is the dynamic viscosity, is the symmetric part of the velocity gradient, and is the shear rate.

The viscosity law, , is defined by Carreau-Yasuda Model: where are the asymptotic viscosities.

Here, we consider physiological data according to [5]:

We assume no slip boundary conditions at the wall, , zero normal stress at the outlet, and a Poiseuille velocity profile at the inlet: where is the spatial mean velocity (taken from [6]) and is the radius of the corresponding artery.

The initial conditions are set to .

5. Numerical Simulations

The numerical simulations were obtained using Comsol Multiphysics, [7] with the followings steps.(1)Definition of the geometry. (2)Computation of the velocity of the fluid in this geometry and analysis of its behavior. (3)Computation of the wall shear stress and choice of the region where the plaque will form, using the function .(4)Computation of , , , , and .(5)In the new geometry (deformed), computation of the fluid velocity, the wall shear stress, and analysis of the behavior of the fluid after deformation.

5.1. Geometry

Atherosclerotic lesions occur mainly in large and medium size arteries such as the abdominal aorta, the coronary arteries, or the carotid arteries. In this paper, we choose to study the formation of the atherosclerotic plaque in carotid arteries.

Carotid arteries are major blood vessels in the neck that supply blood to the brain, neck, and face. There are two carotid arteries, one on the right and one on the left. Each carotid artery can refer to common carotid artery (CCA), internal carotid artery (ICA), and external carotid artery (ECA) (see Figure 8).

We perform numerical simulations on an idealized two-dimensional geometry mimicking the carotid artery bifurcation (Figure 9).

We use the following values as approximations of physiological data, according to [6]:

Here , , and are, respectively, the lengths of CCA, ICA, and ECA, and , , and are, respectively, the diameters of CCA, ICA, and ECA.

The initial height of the intima is , based on [4].

5.2. Inflammatory Process, Plaque Growth, and Fluid Flow Behavior

In order to simplify the numerical simulations, we consider a one-dimensional model for the intima. We neglect the effect of the transport due to the lesion growth on the ox-LDL and on the macrophage.

Thus, the reaction-diffusion equations take the form with , satisfying

We consider the following physical and biological parameters, taken from [4]: cm/s (diffusion coefficient), g/cm3 (initial concentration of LDL in intima), g/cm (a given oxidized LDL quantity), g/cm3, s−1 (degradation rate), cm (intima initial height), cm/(g·s), s−1, g/cm3,and the initial values considered are and .

After computing the velocity field we compute the wall shear stress (WSS). Figure 10 shows the region of lower WSS values.

The permeability function, , will be defined as a characteristic function, such that in the area of low WSS we have penetration of LDL, and otherwise we have no penetration.

Figure 11 shows the plaque already formed near the bifurcation.

The deformed geometry, represented in Figure 12, shows the blood velocity profile and the streamlines.

We can observe that the velocity is higher in the stenotic region, and we can also verify a small region of recirculating flow ion in the right and left corners of the plaque.

Figure 13 shows that in the stenotic region the WSS is higher than before the plaque formation.

6. Study of Possible Regions of Formation of New Plaques

Considering the fact that multiple atherosclerotic plaques can exist, next we study blood velocity profile, WSS values and, consequently, the risky regions of new plaque formation in the following cases: (1)an idealized geometry with one plaque in three different positions (upper wall of ICA, upper wall of CCA, or lower wall of CCA), (2)an idealized geometry with two plaques in two different positions (one in the lower wall of CCA and the other in the lower wall of ICA or one in the lower wall of CCA and the other in the upper wall of ICA ).

In case (1), we first consider a plaque in the upper wall of ICA, next in the upper wall of CCA, and finally in the lower wall of CCA (Figure 14) and observe the blood velocity changes with the position of the atherosclerotic plaque.

Figures 15, 16, and 17 show the WSS computed at the lower walls (lower wall of CCA and ICA), upper walls (upper wall of CCA and ECA), and middle walls (lower wall of ECA and upper wall of ICA), respectively. In each figure, starting from left to right, we considered the following cases: the first one corresponds to the idealized geometry without deformation and the three others to the geometry deformed by an atherosclerotic plaque, at the upper wall of ICA, at the upper wall of CCA, and at the lower wall of CCA, respectively.

Looking to Figures 15, 16, and 17, we can also compare the WSS for the initial geometry, without deformation, with the WSS in the different situations exposed above and observe that(a)the WSS is higher in the stenotic region and lower in the corners of the plaque, where we can also observe the existence of a recirculation flow; (b)the atherosclerotic plaque in CCA does not change the WSS at the lower wall of ECA and upper wall of ICA; (c)when the plaque is at the lower wall, the behavior of WSS at the upper wall does not change much; (d)the atherosclerotic plaque at the upper wall of ICA has a strong influence on the WSS behavior; (e)based on the WSS values we can eventually conclude the possible formation of new plaques on the corners of the older plaques.

Figures 18 and 19 correspond to the second case, showing the behavior of the blood velocity field in the presence of two atherosclerotic plaques in the carotid artery.

Remaining in the second case, Figures 20, 21, and 22 present the study of WSS computed at the lower walls (lower wall of CCA and lower wall of ICA), upper walls (upper wall of CCA and ECA), and middle walls (lower wall of ECA and upper wall of ICA), respectively. In all these Figures 20, 21, and 22, the case of the idealized geometry without deformation is on the left, the geometry deformed by an atherosclerotic plaque in the lower wall of CCA and another in the lower wall of ICA is in the center, and, finally, the geometry deformed by an atherosclerotic plaque in the lower wall of CCA and another in the upper wall of ICA is on the right.

Comparing the WSS depicted in Figures 20, 21, and 22, we observe the same conclusions as those obtained for the different situations studied in the first case.

7. Conclusions and Future Perspectives

Atherosclerosis stages are very complex, and the rupture of the atherosclerotic plaque is an alarming stage. Through the numerical simulations we could better understand this process and how it changes completely the behavior of the blood flow. By computing the WSS in an idealized geometry with and without deformation, we could also detect potential regions of plaque formation.

This work is only the starting point of a larger project. The next steps will be to work with biological parameters and a more realistic geometry (3D geometry reconstructed from medical images). Also, we would like to consider the lower boundary of the interface with the media as an elastic material, including the action of smooth muscle cells, and to use a particle method to simulate the earliest atherosclerosis stages. Finally, we pretend to simulate all stages of atherosclerosis and to use realistic data, obtained from experiments.

Acknowledgments

This work has been partially funded by the University of Cape Verde and by FCT (Fundação para a Ciência e a Tecnologia, Portugal) through the Grants SFRH/BPD/66638/2009 and the Project EXCL/MAT-NAN/0114/2012.