Understanding the Parameter Influence on Lesion Growth for a Mechanobiology Model of Atherosclerosis

Abstract

In this work, we analyse the influence of the parameters of a mathematical model, previously proposed by the authors, for reproducing atheroma plaque in arteries. The model uses Navier–Stokes equations to calculate the blood flow along the lumen in a transient mode. It also uses Darcy’s law, Kedem–Katchalsky equations, and the three-pore model to simulate plasma and substance flows across the endothelium. The behaviours of all substances in the arterial wall are modelled with convection–diffusion–reaction equations, and finally, plaque growth is calculated. We consider a 2D geometry of a carotid artery, but the model can be extrapolated to other geometries or arteries, such as the coronaries or the aorta. A mono-variant sensitivity analysis of the model parameters was performed, with values of $\pm 25\%$ and $\pm 10\%$, with respect to the values of the previous model. The results were analysed with respect to the volume in the plaque of foam cells (FC), synthetic smooth muscle cells (SSMC), and collagen fibre. It was observed that the volume in the plaque of the different substances (FC, SSMC, and collagen) has a strong influence on the results, so it could be used to analyse the vulnerability of plaque. The stenosis ratio of the plaque was also analysed, showing a strong influence on the results as well. Parameters that influence all the results considered when ranged $\pm 10\%$ are the rate of LDL degradation and the diffusion coefficients of LDL and monocytes in the arterial wall. Furthermore, it was observed that the change in the volume of foam cells in the plaque has a greater influence on the stenosis ratio than the change of synthetic smooth muscle cells or collagen fibre.

1. Introduction

Cardiovascular diseases, including atherosclerosis, are some of the main causes of mortality in developed countries [1]. They consist of the development of atheroma plaque in the arterial wall, which is caused by lipid deposition. This leads to an increase in the thickness of the arterial wall and, therefore, to a decrease in the area where blood flows, called the lumen. In addition, plaque can break and cause a blood clot that can travel through the circulatory system. For these reasons, it can lead to several consequences that are dependent on the artery affected, e.g., myocardial attacks, ischaemia, or ictus, among others.

The process of atheroma plaque formation begins with the deposition of low-density lipoproteins (LDLs) in the arterial wall and, once there, oxidise. Then, monocytes are deposited from the lumen into the arterial wall and differentiate into macrophages, which ingest oxidised LDL. Once they cannot ingest more oxidised LDL, they become foam cells. These foam cells compose the lipid core of the atheroma plaque. Concurrently, at the beginning of the inflammatory process, all muscular cells in the arterial wall are contractile, i.e., they cannot move or react with any other substance. However, due to the presence of cytokines in the arterial wall (segregated by macrophages), these muscle cells change their phenotype and become synthetic smooth muscle cells, which can move and interact with other substances in the arterial wall. Then synthetic smooth muscle cells surround foam cells, macrophages, and oxidised LDL, and segregate collagen fibre to isolate the lipid core, forming the fibrous layer of the plaque [2,3,4].

Therefore, atheroma plaques are developed as a consequence of an increase in endothelial permeability. This variation in permeability was accepted to be due to a change in the shape of endothelial cells, depending on some mechanical stimuli caused by blood flow [5]. These mechanical stimuli could be, among others, the wall shear stress of blood with the endothelium ($WSS$), time-averaged wall shear stress ($TAWSS$), and oscillatory shear index ($OSI$) in a cardiac cycle [6,7], or a combination of some of them [8].

There are some studies in the literature related to the influence of blood hemodynamics on the location of atheroma plaque for patient-specific geometries. Malvè et al. [9] modelled a carotid artery bifurcation and analysed the influence of $TAWSS$ and $OSI$ on the shape index ($SI$) of endothelial cells. Sáez et al. [6] analysed $WSS$ on eight different geometries of coronary bifurcations. Alimohammadi et al. [7] analysed $TAWSS$, $OSI$, and a new index that they previously proposed [10], Highly Oscillatory Low Magnitude Shear (HOLMES), in the bifurcation of the abdominal aorta.

Some continuum mathematical models of the development of atherosclerosis have been developed. There is a huge quantity of these models developed as two-dimensional axisymmetric models. Olgac et al. [11] consider the LDL flow from the lumen to the arterial wall with the three-pore model, depending on the $WSS$, while Tomaso et al. [12] considered the flow of monocytes as well as LDL flow. Calvez et al. [13] used a two-dimensional model which, in addition to LDL and monocytes, also considered macrophages, cytokines, and foam cells. Finally, they modelled the fluid–structure interaction between the mesh displacement due to plaque formation and blood flow. Cilla et al. [4] modelled a two-dimensional axisymmetric geometry with LDL and monocyte flows into the arterial wall depending on $WSS$ and the three-pore model, and the interactions between all the substances commented on before, as well as smooth muscle cells and collagen fibre. In addition, Shahzad et al. [14] studied the influence of hemodynamics in a two-dimensional geometry of an artery bifurcation with an obstacle plaque that disturbs blood flow. They used the fluid–structure interaction (FSI) with an elastic wall and computational fluid dynamics (CFD) considering a rigid wall.

There are also mechanobiological models with patient-specific geometries, such as by Siogkas et al. [15], who consider LDL that becomes oxidised, monocytes that differentiate into macrophages, and cytokines in the arterial wall. Filipovic et al. [16] modelled oxidised LDL, macrophages, and cytokines, depending on the $WSS$ on coronary arteries, and adjusted the parameters of their model in order to obtain results based on experimental data. Hernández-López et al. [8] included also foam cells, smooth muscle cells, and collagen fibre, and analysed the effect of mechanical stimuli on the growth of atheroma plaque in carotid arteries. Kenjereš and De Loor [17] developed a computational model to simulate the flow of LDL through the arterial wall in a realistic geometry of a carotid bifurcation with a multilayered wall.

On the other hand, there are some agent-based models that study plaque development in the arteries [18,19,20]. The main advantage of continuum models against agent-based models is that they allow the simulation of plaque development in real geometries, whereas the advantage of agent-based models is that they can take into account the random behaviour of cells, which is not possible with continuum models. Olivares et al. [18] developed a model for early-stage atherosclerosis, in which they considered substances such as LDL, oxidised LDL, macrophage, foam cells, smooth muscle cells, endothelial cells, and autoantibodies. Bhui and Hayenga [19] developed a three-dimensional model of transendothelial migration during atherogenesis, in which they modelled endothelial cells, monocyte, macrophage, lymphocyte, neutrophil, and foam cells. In addition, Corti et al. [20] modelled atheroma plaque growth depending on $WSS$ in three-dimensional geometry. They considered a multilayered wall, composed of the intima, media, and adventitia layers, modelling also the internal and external elastic laminae.

There is a large number of substances and parameters involved in the process of formation of atheroma plaque, so it is important to understand how each substance and parameter influences plaque growth.

There are also some studies focused on parameter influence in previous related models, such as Tomaso et al. [12], who analysed the influence of different LDL concentrations in blood to determine how it affects plaque growth. Olivares et al. [18] analysed, in 8 cases, the influence in the plaque of the oxidation rate of LDL, migration speed of agents, and auto-antibodies maximum concentration. Cilla et al. [21] analysed the effect of transmural transport properties on atheroma plaque formation, attending to LDL and oxidised LDL diffusion coefficient anisotropy. Finally, Corti et al. [20] made both, mono-parametric and multi-parametric sensitivity analyses to evaluate the changes in their results caused by changes in the parameters of an agent-based model of plaque growth.

The aim of this study is to analyse the influence of the parameters of a previous computational model on the formation of atheroma plaque in the arteries [8], to determine the effect of these parameters on the growth and composition of atheroma plaque.

The mathematical model used has a total of 52 parameters, which come from different studies and can be related to experimental or computational analysis. Among the experimental ones, there are some differences between the analysis conditions. For example, some of the parameters have been determined from in vitro tests, while others are from in vivo tests. Moreover, they come from studies of different species and arteries (coronary, carotid, or aorta). Finally, some of them have been estimated based on computational results. Therefore, the values of the parameters can have a large variation, so it is relevant to perform a sensitivity analysis of the parameters of the model.

Although the study of parameter variation was performed in a geometry of a carotid artery, the influence of the parameters in the model would be the same in the case of other arteries or geometries.

2. Geometry

Due to the high computational requirements of the model and because the objective of this study was to analyse the influences of the different parameters on plaque growth, a 2D axisymmetric approximation was considered enough for this study.

The geometry was developed based on that of Olgac et al. [11] as it reproduces the mechanical stimuli to which real patients are subjected, leading to the appearance of plaque ($TAWSS$ and $OSI$). The geometry was adapted to carotid arteries by modifying the diameter and thickness of the vessel to the corresponding values for carotid arteries: 3.63 mm and 0.7 mm, respectively [22]. Figure 1 shows the geometry, which is composed of the lumen and the arterial wall, considered a monolayer. It also has an obstacle plaque that disturbs blood flow, causing low $TAWSS$ and high $OSI$ downstream and favouring the appearance of a new plaque in that area. It is a phenomenon that was observed in real patients [23,24]. This second plaque is the one that is going to be analysed.

The model was discretised using the finite element method. A sensitivity analysis of the length of the geometry and the mesh was developed, obtaining that the minimal length of the artery to ensure a fully developed blood flow is 90 mm. The domains of the model were meshed using quadrilateral elements, with three boundary layers in the lumen near the endothelium. A sensitivity analysis of the mesh was performed in both lumen and arterial wall meshes, in order to determine the most suitable mesh to solve the model with minimal computational costs, but without affecting the results. The total number of elements is equal to 36,731 for the lumen and 36,036 for the arterial wall, resulting in a total number of elements in the geometry of 72,767. A discretisation with P1-P1 elements was made for the calculation of blood flow. Therefore, linear interpolation was considered for both the velocity and pressure of blood flow in the lumen of the artery. To perform the sensitivity analysis of the mesh, a coarse mesh was selected in the first place and progressively refined, until the results obtained did not vary more than $5\%$. The same procedure was done with the number of boundary layers of the lumen mesh. In addition, linear and quadratic serendipity elements were used for the inflammatory process in the arterial wall and for plasma flow through the endothelium and the growth of the plaque, respectively. Two details of the selected mesh are shown in Figure 1. The solution to the transient problems was obtained using the Backwards differentiation formula (BDF) method, which is implicit, following Newton’s method for nonlinear problems.

The total number of degrees of freedom of the lumen can be determined by multiplying the components of the velocity of blood and the number of elements in the lumen mesh. For the case of the arterial wall, the number of degrees of freedom would be equal to the multiplication of the components of the plasma flow velocity, the number of substance concentrations that we calculate and the number of elements of the arterial wall mesh.

3. Mathematical Model

In this section, the mathematical model is described, showing the equations in different subsections, referring to blood flow in the lumen, plasma flow across the endothelium, the inflammatory process in the arterial wall, and plaque growth.

3.1. Blood Flow in the Lumen

Blood flow along the arterial lumen is considered laminar due to the Reynolds number in the arteries under physiological conditions ($Re\approx 950$ for the mean diameter of the considered artery) [25]. Furthermore, blood is modelled as a Newtonian incompressible fluid, since the diameter of the lumen in the arteries is greater than 0.5 mm [26,27]. Moreover, it is modelled as a homogeneous fluid due to the small size of the particles with respect to the lumen diameter of the arteries [25].

The Navier–Stokes and continuity equations govern blood flow along the lumen:

$${\rho }_b(u_l·\nabla )u_l=\nabla ·[-P_{l}I+{\mu }_b(\nabla u_l+{(\nabla u_l)}^T)]+F_l$$

(1)

$${\rho }_b\nabla ·u_l=0,$$

(2)

where subscripts b and l refer to blood and lumen, respectively. The parameters ${\rho }_b$ and ${\mu }_b$ are, respectively, the density and dynamic viscosity of blood, while $u_l$ and $P_l$ are the velocity and pressure of blood flow. $F_l$ is the internal force of blood, which is negligible compared to the friction of blood flow with the arterial wall. All parameters referring to blood flow along the lumen are in

Table 1. Parameters to calculate blood flow.

Blood Flow Parameters
Parameter	Description	Value	Reference
${\rho }_b$	Blood density	$1050\frac{kg}{m^3}$	Milnor [28]
${\mu }_b$	Blood viscosity	$0.0035Pa·s$	Milnor [28]
T	Cardiac cycle period	$0.85$s	Malvè et al. [25]

.

Blood flow is modelled in a transient step. An analysis of the number of cardiac cycles needed to model blood flow was performed, determining that three cardiac cycles are enough to develop blood flow in the lumen.

Transient mass flow and pressure are imposed, respectively, at the inlet and outlet of the geometry as boundary conditions [25], as can be seen in Figure 2. These boundary conditions influence the location and size of the generated plaque as they directly affect the recirculation that the obstacle plaque produces downstream. Therefore, they directly affect the mechanical stimuli that initiate plaque growth. For example, when increasing the inlet mass flow, the recirculation is bigger and, therefore, the areas with the mechanical stimuli that initiate plaque growth move to the end of the geometry. Finally, a no-slip condition is imposed at the endothelium.

3.2. Plasma Flow across the Endothelium

As a result of the porous nature of the arterial wall, plasma can enter the bloodstream, causing a plasma flow through the endothelium. Plasma flow across the arterial wall is modelled with Darcy’s law and continuity equations [11,29,30,31,32,33]. Darcy’s law relates the velocity of the plasma flow with the pressure drop in the arterial wall and considers some parameters related to the porous nature of the arterial wall (permeability and porosity).

$$u_w=-\frac{k_w}{{\mu }_p}\nabla p_w,$$

(3)

$$\frac{\partial \left({\rho }_p{ϵ}_w\right)}{\partial t}+{\rho }_p\nabla ·u_w=J_v$$

(4)

Parameters related to the arterial wall are $u_w$, $k_w$, ${ϵ}_w$, and $\nabla p_w$, and represent the velocity of the plasma flow, the permeability and porosity of the arterial wall, and the pressure gradient, respectively. ${\rho }_p$ and ${\mu }_p$ are the dynamic viscosity and density of the plasma. Finally, $J_v$ is the plasma flow across the endothelium, which is calculated with the Kedem–Katchalsky equations [34] and the three-pore model [11]:

$$J_v=J{v}_{nj}+J{v}_{lj}+\cancel{J{v}_{vp}}$$

(5)

As can be seen in Equation (5), the three-pore model considers three different types of plasma flow through the endothelium: plasma flow through normal junctions ($J{v}_{nj}$), leaky junctions ($J{v}_{lj}$), and vesicular pathways ($J{v}_{vp}$). Nevertheless, vesicular transport is intended for molecular transport, so plasma flow is negligible, relative to that across normal and leaky junctions.

Plasma flow through normal and vesicular pathways can be defined as:

$$J{v}_{nj}=L{p}_{nj}·\mathsf{\Delta }{P}_{End}·(1-{\mathsf{\Phi }}_{lj})$$

(6)

$$J{v}_{lj}=L{p}_{lj}·\mathsf{\Delta }{P}_{End}$$

(7)

$L{p}_{nj}$ and $L{p}_{lj}$ depend on the considered artery and are the hydraulic conductivity of normal and leaky junctions [35], while $\mathsf{\Delta }P$ is the pressure drop across the endothelium, which depends on the intraluminal pressure [35]. Finally, ${\mathsf{\Phi }}_{lj}$ is the endothelial fraction occupied by leaky cells [36,37,38].

The hydraulic conductivity of leaky junctions can be calculated as [38,39]:

$$L{p}_{lj}=\frac{Ap}{S}·L{p}_{slj},$$

(8)

with $\frac{Ap}{S}$ being the fraction of the total area occupied by leaky junctions, and $L{p}_{slj}$ being the hydraulic conductivity of a single leaky junction. To determine these parameters, we considered that leaky cells are enclosed by leaky junctions, which have ring shapes and are randomly distributed with a distance equal to ${ϵ}_{lj}$, being $2{ϵ}_{lj}$ of the permeability of a leaky junction. Moreover, the spaces between endothelial cells are cylindrical [38,39]. Therefore, $\frac{Ap}{S}$ can be calculated as:

$$\frac{A_p}{S}=\frac{A_{slj}}{\pi ·{{ϵ}_{lj}}^2},$$

(9)

where $\pi ·{{ϵ}_{lj}}^2$ is the total area occupied by leaky junctions, and $A_{slj}$ is the area of a single leaky junction that can be calculated, assuming reduced thickness, as:

$$A_{slj}=2\pi ·R_{cell}·2w_l,$$

(10)

where $R_{cell}$ is the radius of endothelial cells and $w_l$ is the half-width of a leaky junction [36,38,39].

On the other hand, ${\mathsf{\Phi }}_{lj}$ can be defined as [36,38]:

$${\mathsf{\Phi }}_{lj}=\frac{{R_{cell}}^2}{{{ϵ}_{lj}}^2}$$

(11)

Combining Equations (9)–(11), the fraction of area occupied by leaky junctions can be rewritten as (12):

$$\frac{A_p}{S}=\frac{4w_l}{R_{cell}}·{\mathsf{\Phi }}_{lj}$$

(12)

The ratio ${\mathsf{\Phi }}_{lj}$ is known to be a function of some mechanical stimuli, such as the wall shear stress ($WSS$), the time-averaged wall shear stress ($TAWSS$), and the oscillation shear index ($OSI$) [40]. As blood flow is modelled in transient mode, $WSS$ was not considered due to its changes according to the cardiac cycle. Instead of $WSS$, $TAWSS$ was defined. This variable takes into account the tangential stresses that blood produces in the endothelium during a cardiac cycle and that causes changes in the shape index of the endothelial cells. Therefore, in areas with low $TAWSS$, endothelial cells will be rounder, allowing plasma and substances to flow between their pores [16,40,41,42]:

$$TAWSS=\frac{1}{T}{\int }_o^T\left|\tau \left(t\right)\right|·dt,$$

(13)

where T is the period of a cardiac cycle and $\left|\tau \right(t\left)\right|$ the magnitude of $WSS$ dependent on time. On the other hand, high values of $OSI$ imply that blood is highly oscillatory there. Thus, the tangential stresses of blood flow in the endothelium are also oscillatory and, due to their variation, they cause a change in the shape index of the endothelial cells, making them rounder and allowing the flow of plasma and substances through the endothelium, being areas of high $OSI$ considered as atheroprones [8,40]. $OSI$ is defined as:

$$OSI=0.5·\left(1-\frac{|\frac{1}{T}{\int }_o^T\tau \left(t\right)·dt|}{TAWSS}\right)$$

(14)

Moreover, the ratio ${\mathsf{\Phi }}_{lj}$, $TAWSS$, and $OSI$ are related by experimental correlations [6,7,8,11,40]. The first correlation is a parameter named shape index ($SI$), which can be calculated as:

$$SI=\frac{4\pi ·Area}{Perimete{r}^2},$$

(15)

where $Area$ and $Perimeter$ correspond to a single endothelial cell.

To determine the behaviour of endothelial cells with $TAWSS$ and $OSI$, we proposed in a previous work [8] a numerical correlation based on the experimental results of Levesque et al. [40].

$$SI=0.0264·e^{5.647·OSI}+0.5513·e^{-0.1815·{\left(TAWSS\right)}^2}$$

(16)

The total number of mitotic cells ($MC$) in the endothelium depends on the $SI$ [43], so the next experimental correlation was developed, with a unit area of 0.64 mm${}^2$:

$$MC=0.003797·e^{(14.75·SI)}$$

(17)

Finally, the number of leaky cells ($LC$) and mitotic cells ($MC$) is correlated by the following empirical equation [11,44]:

$$LC=0.307+\frac{0.805·MC}{0.453}$$

(18)

On the other hand, ${\mathsf{\Phi }}_{lj}$ can be defined as the ratio between the area of leaky cells and the total area of cells. It can be calculated as:

$${\mathsf{\Phi }}_{lj}=\frac{LC·\pi ·R_{cell}^2}{A_{unit}},$$

(19)

considering $A_{unit}$ equal to the unit area in all the experimental correlations before (0.64 mm${}^2$). By calculating blood flow in the model, it is now possible to determine ${\mathsf{\Phi }}_{lj}$ and therefore $\frac{A_p}{S}$.

Finally, the hydraulic conductivity of a single leaky junction, $L{p}_{slj}$, is defined according to Olgac et al. [11]:

$$L{p}_{slj}=\frac{w_l^2}{3·{\mu }_p·l_{lj}},$$

(20)

where $w_l$ and $l_{lj}$ are, respectively, the width and the length of a leaky junction.

So, the plasma flow across the endothelium is completely defined, and the parameters needed to calculate it are shown in

Table 2. Parameters to model plasma flow across the endothelium.

Plasma Flow Parameters
Param.	Description	Value	Reference
$R_{cell}$	Endothelial cell radius	$15\mathsf{\mu }$m	Weinbaum et al. [38]
$w_l$	Half-width of a leaky junction	$20$nm	Weinbaum et al. [38]
$l_{lj}$	Leaky junction length	$2\mathsf{\mu }$m	Weinbaum et al. [38]
${\rho }_p$	Plasma density	$1050\frac{kg}{m^3}$	Milnor [28]
${\mu }_p$	Plasma viscosity	$0.001Pa·s$	Milnor [28]
$k_w$	Darcian artery permeability	1.2 $·{10}^{-18}$m${}^2$	Huang et al. [36], Prosi et al. [45]
${ϵ}_p$	Intima porosity	0.96	Ai and Vafai [46]
$L_{p,nj}$	Normal junction conductivity	$1.984·{10}^{-12}\frac{m}{s·Pa}$	Tedgui and Lever [35]
$\mathsf{\Delta }{P}_{End}$	Endothelial pressure difference	20.727 mmHg	Tedgui and Lever [35]
$A_{unit}$	Unit area for the exp. correlations	$0.64$mm${}^2$	Chien [43]
$P_{adv}$	Pressure of the adventitia	$17.5$mmHg	Olgac et al. [11]

.

The boundary conditions of the plasma flow across the endothelium are the normal velocity of the plasma flow ($J_v$), which has already been defined, and the pressure in the adventitia, which [11] must be defined. Both an increase in the normal velocity of the plasma flow and a decrease in pressure at the adventitia would produce an increase in plasma flow, which will also cause an increase in substances flow, with the consequent increase in plaque size [8].

3.3. Inflammatory Process in the Arterial Wall

There are many substances involved in the inflammatory process of the arterial wall. In the model, we consider molecules, such as low-density lipoproteins (LDL) and oxidised LDL (oxLDL). We also consider cells, such as monocytes (m), macrophages (M), foam cells (FC), and contractile and synthetic smooth muscle cells (CSMC and SSMC). We also model cytokines (C) and collagen fibre (Cg).

All of these substances can have convection and diffusion and interact with the others, so the inflammatory process of the arterial wall can be described with convection–diffusion–reaction equations.

In addition, the flow of substances across the arterial wall can be defined as:

$$N=-D_{X_i}\nabla X_i+u_w{X}_i,$$

(21)

where $D_{X_i}$, $X_i$, and $u_w$ are the diffusion coefficient of the considered substance, its concentration and its convection velocity in the arterial wall. Due to the structure of the arterial wall, the diffusion is anisotropic, as the longitudinal diffusion is 3 times higher than radial diffusion [39].

At the beginning of the process, in the healthy artery, all muscle cells in the arterial wall are contractile. Thus, the only substance that has an initial concentration in the arterial wall is CSMC.

A high level of LDL concentration in blood is considered, according to a pathological model ($2.7\frac{mg}{ml}$ [47], which corresponds to $6.98\frac{mol}{m^3}$). Moreover, we impose a physiological monocyte inlet concentration in the lumen of $550·{10}^{-9}\frac{cells}{m^3}$ [48]. An increase in both inlet concentrations would cause an increase in the growth of the resultant plaque.

3.3.1. LDL

Once the LDL molecules cross the endothelial barrier, they can suffer both convection and diffusion, because of their small sizes. When an LDL molecule enters the arterial wall, it is oxidised. Therefore, the reactive term of LDL molecules represents their oxidation.

$$\frac{\partial C_{LDL,w}}{\partial t}+\nabla ·(-D_{C_{LDL,w}}\nabla C_{LDL,w})+K_{lag}·u_w·\nabla C_{LDL,w}=-d_{LDL}{C}_{LDL,w}$$

(22)

$d_{LDL}$ is the degradation ratio of LDL in the arterial wall and $C_{LDL,w}$ is its concentration at each time.

It is necessary to define as a boundary condition the LDL concentration in the adventitia ($C_{LDL,adv}$) [49]. The Kedem–Katchalsky equation is used to determine the LDL flow across the endothelium [11]:

$$J_{S,LDL}=C_{LDL,l}·LD{L}_{dep}·P_{app},$$

(23)

where $C_{LDL,l}$ is the LDL concentration in blood and $LD{L}_{dep}$ is the amount of LDL molecules that are deposited from the lumen to the arterial wall. Finally, $P_{app}$ is the apparent permeability of the arterial wall. LDL can flow from the lumen to the arterial wall through normal and leaky junctions and vesicular pathways, so the apparent permeability of the arterial wall is composed of three types of permeabilities [50,51]:

$$P_{app}=P_{app,lj}+\cancel{P_{app,nj}}+P_{app,vp},$$

(24)

with $P_{app,nj}$, $P_{app,lj}$, and $P_{app,vp}$ being the permeabilities of normal and leaky junctions and vesicular pathways, respectively.

The transport of molecules through the endothelium is dependent on the size of the particles. For the case of LDL, which has a radius of 11 nm [50]), transport across normal junctions is not possible due to their small size. Therefore, LDL transport through the endothelium can only occur through leaky junctions and vesicular pathways [52].

In addition, LDL transport through vesicular pathways is $10\%$ of the flux through leaky junctions [11]:

$$P_{app,vp}=0.1·P_{app,lj}$$

(25)

The apparent permeability of leaky junctions can be defined as:

$$P_{app,lj}=P_{lj}{Z}_{lj}+J_{v,lj}·(1-{\sigma }_{f,lj}),$$

(26)

where $P_{lj}$, $Z_{lj}$, and ${\sigma }_{f,lj}$ are the diffusive permeabilities of the leaky junctions, factors of reduction of the LDL concentration gradient at the inlet of the flow and the solvent-drag coefficient of leaky junctions. So, the LDL flux across the endothelium can be written as:

$$J_{S,LDL}=1.1·C_{LDL,l}·LD{L}_{dep}·(P_{lj}{Z}_{lj}+J_{v,lj}(1-{\sigma }_{f,lj}))$$

(27)

The diffusive permeability of leaky junctions is defined as:

$$P_{lj}=\frac{A_p}{S}\chi P_{slj},$$

(28)

where $\chi$ is the difference between the total area of endothelial cells and the area of cells separated by leaky junctions:

$$\chi =1-{\alpha }_{lj},$$

(29)

with ${\alpha }_{lj}$ being the ratio between the radius of an LDL molecule ($a_m$) and the half-width of a leaky junction ($w_l$):

$${\alpha }_{lj}=\frac{a_m}{w_l}$$

(30)

$P_{slj}$ is the permeability of a single leaky junction, and can be determined by knowing the LDL diffusion coefficient in a leaky junction ($D_{lj}$) and the length of a leaky junction ($l_{lj}$):

$$P_{slj}=\frac{D_{lj}}{l_{lj}}$$

(31)

LDL diffusion coefficient in a leaky junction is related to the LDL diffusion coefficient by [11]:

$$\frac{D_{lj}}{D_l}=F\left({\alpha }_{lj}\right)=1-1.004{\alpha }_{lj}+0.418{\alpha }_{lj}^3-0.16{\alpha }_{lj}^5$$

(32)

On the other hand, $Z_{lj}$ depends on a modified Péclet number:

$$Z_{lj}=\frac{P{e}_{lj}}{e^{\left(P{e}_{lj}\right)}-1}$$

(33)

This modified Péclet number can be defined as:

$$P{e}_{lj}=\frac{Jv,lj·(1-{\sigma }_{f,lj})}{P_{lj}}$$

(34)

Finally, the solvent-drag coefficient of leaky junctions is given by [11]:

$${\sigma }_{f,lj}=1-\frac{2}{3}{\alpha }_{lj}^2(1-{\alpha }_{lj})·F\left({\alpha }_{lj}\right)-(1-{\alpha }_{lj})\left(\frac{2}{3}+\frac{2{\alpha }_{lj}}{3}-\frac{7{\alpha }_{lj}^2}{12}\right)$$

(35)

3.3.2. Oxidised LDL

Once LDL is oxidised, it is considered to not experience convection [4]. However, due to its similar size to LDL, oxidised LDL shows diffusion in the arterial wall. Once LDL enters the arterial wall, it becomes oxidised, so one of its reactive terms refers to this oxidation. On the other hand, oxidised LDL is phagocytosed by macrophages, which corresponds to its second reactive term.

$$\frac{\partial C_{oxLDL,w}}{\partial t}+\nabla ·(-D_{C_{oxLDL,w}}\nabla C_{oxLDL,w})=d_{LDL}{C}_{LDL,w}-LD{L}_{ox,r}{C}_{oxLD{L}_w}{C}_{M,w},$$

(36)

where $C_{oxLD{L}_w}$ and $C_{M,w}$ are the oxidised LDL and macrophage concentrations in the arterial wall. In addition, $LD{L}_{ox,r}$ is a ratio of the quantity of oxidised LDL that a macrophage can ingest.

3.3.3. Monocytes

Due to their sizes, we hypothesise that monocytes do not have convection, but they experience diffusion in the arterial wall. Once monocytes are in the arterial wall, they differentiate into macrophages, and also, they suffer apoptosis. These two phenomena are represented by their reactive terms.

$$\frac{\partial C_{m,w}}{\partial t}+\nabla ·(-D_{C_{m,w}}\nabla C_{m,w})=-d_m{C}_{m,w}-m_d{C}_{m,w},$$

(37)

where $C_{m,w}$ is the monocyte concentration in the arterial wall. $d_m$ is the rate of monocytes that differentiate into macrophages and $m_d$ is the apoptosis rate of monocytes.

In addition, a physiological monocyte inlet concentration in the lumen of $550·{10}^{-9}\frac{Monocyte}{m^3}$ is imposed [48].

Monocyte flow through the endothelium is dependent on hemodynamics. Areas of low $TAWSS$ are known to be atheroprone; particularly for carotid arteries, areas of $TAWSS$ lower than 2 Pa can develop atheroma plaque [16,41,42]. For the case of $OSI$, according to our previous work [8], it can be determined that areas of $OSI$ greater than 0.1910 are atheroprones. Monocyte flow across the endothelium can be defined as [8]:

$$J_{s,m}(TAWSS,OSI)=m_r·(0.8588·e^{-0.6301·TAWSS}+0.1295·e^{3.963·OSI})·C_{LDL,ox,w}{C}_{m,l},$$

(38)

with $m_r$ being the monocyte recruitments from the lumen.

3.3.4. Macrophages

Similar to monocytes, we hypothesise that macrophages experience diffusion but not convection. Once monocytes are in the arterial wall, they differentiate into macrophages. Furthermore, macrophage phagocytes oxidised LDL, and once they cannot ingest a greater amount of oxidised LDL, they suffer apoptosis and become foam cells. These two phenomena are represented in the macrophage reactive terms:

$$\frac{\partial C_{M,w}}{\partial t}+\nabla ·(-D_{C_{M,w}}\nabla C_{M,w})=d_m{C}_{m,w}-\frac{LD{L}_{ox,r}}{n_{FC}}·C_{oxLD{L}_w}{C}_{M,w},$$

(39)

In Equation (39), $LD{L}_{ox,r}$ is the constant rate of oxidised LDL uptaken by macrophages, and $n_{FC}$ is the maximum amount of oxidised LDL that a single macrophage must ingest to transform into a foam cell.

3.3.5. Cytokines

Due to their sizes, cytokines do not have convection. In addition, they are enclosed by macrophages, so they do not have diffusion either. Cytokines in the arterial wall are segregated by macrophages and also experience degradation, phenomena that can be seen in their reactive terms:

$$\frac{\partial C_{c,w}}{\partial t}=C_r{C}_{oxLD{L}_w}{C}_{M,w}-d_c{C}_{c,w},$$

(40)

where $C_{c,w}$ is cytokine concentration in the arterial wall. $C_r$ and $d_c$ are, respectively, the ratios of production and degradation of cytokines.

3.3.6. Contractile Smooth Muscle Cells

They have neither convection nor diffusion because of their large size. At the beginning of the process, all smooth muscle cells of the arterial wall have a contractile phenotype. Then, due to the presence of cytokines, they experience a change of phenotype into synthetic smooth muscle cells, which is represented in their reactive term:

$$\frac{\partial C_{CSMC,w}}{\partial t}=-C_{CSMC,w}·S_r·\left(\frac{C_{c,w}}{k_c·C_{c,w}^{th}+C_{c,w}}\right)$$

(41)

$C_{CSMC,w}$ is CSMC concentration in the arterial wall. $S_r$ is the CSMC differentiation rate into SSMC and $C_{c,w}^{th}$ is the maximum cytokine concentration allowed in the arterial wall. Finally, $k_c$ is a constant for the saturation equation.

3.3.7. Synthetic Smooth Muscle Cells

Equal to CSMC and because of their sizes, SSMCs neither have convection nor diffusion. CSMCs change their phenotypes into SSMCs due to the presence of cytokines in the arterial wall. They also experience proliferation and apoptosis. These three phenomena are represented in their reactive terms:

$$\begin{array}{c}\hfill \frac{\partial C_{SSMC,w}}{\partial t}=C_{CSMC,w}·S_r·\left(\frac{C_{c,w}}{k_c·C_{c,w}^{th}+C_{c,w}}\right)+\\ \hfill \left(\frac{p_{ss}{C}_{c,w}}{C_{c,w/2}^{th}+C_{c,w}}\right)C_{SSMC,w}-r_{Apop}·C_{SSMC,w}\end{array}$$

(42)

$C_{SSMC,w}$ represents the SSMC concentration in the arterial wall. In addition, $p_{ss}$ is the SSMC proliferation rate and $r_{Apop}$ is the SSMC apoptosis rate.

3.3.8. Foam Cells

Due to their large sizes, foam cells neither have convection nor diffusion. Once a macrophage cannot ingest a greater quantity of oxidised LDL, it becomes a foam cell, which is represented in the reactive term of foam cells:

$$\frac{\partial C_{FC,w}}{\partial t}=\frac{LD{L}_{ox,r}}{n_{FC}}·C_{LDLo{x}_w}{C}_{M,w}$$

(43)

3.3.9. Collagen Fibre

Collagen fibre do not experience convection or diffusion because they are composed of many molecules. Collagen fibre experience segregation by SSMC, which can be seen in its reactive term.

$$\frac{\partial C_{Cg,w}}{\partial t}=G_r·C_{SSMC,w},$$

(44)

where $C_{Cg,w}$ is collagen concentration in the arterial wall, and $G_r$ is its secretion rate due to plaque formation. Natural segregation and degradation of collagen in the arterial wall were not considered, as they occur in healthy areas of the arterial wall, not related to plaque generation.

Table 3. Parameters to calculate the inflammatory process in the arterial wall.

Inflammatory Process Parameters
Parameter	Description	Value	Reference
$D{r}_{LDL,w}$	LDL radial diff. coeff.	$8·{10}^{-13}\frac{m^2}{s}$	Prosi et al. [45]
$D{r}_{m,w}$	Monocyte radial diff. coeff.	$8·{10}^{-15}\frac{m^2}{s}$	Cilla et al. [4]
$D{r}_{oxLDL,w}$	Ox. LDL radial diff. coeff.	$8·{10}^{-13}\frac{m^2}{s}$	Prosi et al. [45]
$D{r}_{M,w}$	Macroph. radial diff. coeff.	$8·{10}^{-15}\frac{m^2}{s}$	Cilla et al. [4]
$d_{LDL}$	LDL oxidation	$2.85·{10}^{-4}$s${}^{-1}$	Ai and Vafai [46]
$d_m$	Monocyte differentiation	$1.15·{10}^{-6}$s${}^{-1}$	Bulelzai and Dubbeldam [53]
$m_d$	Monocyte natural death	$\frac{1}{60}$d${}^{-1}$	Krstic [54]
$LD{L}_{ox,r}$	OxLDL uptake	$2.45·{10}^{-23}\frac{m^3}{Macrophage·s}$	Zhao et al. [55]
$n_{FC}$	Max. oxLDL uptake	$2.72·{10}^{-11}\frac{mo{l}_{oxLDL}}{Macrophage}$	Hernández-López et al. [8]
$C_r$	Cytokine production	$3·{10}^{-10}\frac{mo{l}_C·m^3}{mo{l}_{oxLDL}·Macrophage·s}$	Cilla et al. [4]
$d_c$	Cytokine degradation	$2.3148·{10}^{-5}$s${}^{-1}$	Zhao et al. [56]
$S_r$	SMC’s differentiation	$0.0036$d${}^{-1}$	Cilla et al. [4]
$p_{ss}$	SSMC’s proliferation	$0.0202$d${}^{-1}$	Budu-Grajdeanu et al. [57]
$G_r$	Collagen production	$2.49·{10}^{-21}\frac{k{g}_{Cg}}{SSMC·s}$	Zahedmanesh et al. [58]
$w_l$	Half-width of a leaky junct.	$20$nm	Weinbaum et al. [38]
$l_{lj}$	Leaky junction length	$2$$\mathsf{\mu }$m	Weinbaum et al. [38]
$R_{LDL}$	LDL radius	$11$nm	Tarbell [50]
$C_{c,w}^{th}$	Cytokine threshold	$1.235·{10}^{13}\frac{mo{l}_C}{m^3}$	Hernández-López et al. [8]
$C_{0,LDL}$	LDL initial conc.	$6.98\frac{mo{l}_{LDL}}{m^3}$	Cilla [47]
$C_{0,m}$	Monocyte initial conc.	$550·{10}^9\frac{Monocytes}{m^3}$	[48]
$C_{0,CSMC}$	CSMC initial conc.	$3.16·{10}^{13}\frac{CSMC}{m^3}$	Escuer et al. [59]
$C_{LDL,adv}$	LDL conc. at adventitia	$11.6$‰ $·C_{LDL,l}$	Meyer et al. [49]
$k_c$	Cytokine threshold factor	$0.65093$	Hernández-López et al. [8]
$LD{L}_{dep}$	LDL at the endothelium	${10}^{-2}·C_{LDL,l}$	Meyer et al. [49]
$m_r$	Monocyte recruitment	$6.636·{10}^{-4}\frac{m^4}{mo{l}_{oxLDL}·day}$	Steinberg et al. [60]
$r_{apop}$	SSMC apoptosis rate	$8.011·{10}^{-8}$s${}^{-1}$	Kockx et al. [61]
${\rho }_{LDL}$	LDL density	$1063\frac{kg}{m^3}$	Ivanova et al. [62]
$M{w}_{LDL}$	LDL molecular weight	$386.65\frac{g}{mo{l}_{LDL}}$	Cilla [47]
$k_{lag}$	Solute lag coefficient of LDL	$0.893$	Dabagh et al. [63]

contains all the parameters to calculate the inflammatory process in the arterial wall.

3.4. Plaque Growth

Finally, with Equation (45), it is possible to calculate the growth of plaque in the arterial wall. We consider plaque to be composed of a lipid nucleus of foam cells and a fibrous layer of synthetic smooth muscle cells and collagen fibre. Therefore, considering the isotropic growth of plaque, it is possible to determine the change in volume in the arterial wall due to plaque appearance:

$$\nabla ·v=\frac{\partial C_{FC,w}}{\partial t}·Vo{l}_{FC}+\frac{\partial \mathsf{\Delta }{C}_{SSMC,w}}{\partial t}·Vo{l}_{SSMC}+\frac{\partial C_{Cg,w}}{\partial t}·\frac{1}{{\rho }_{Cg}},$$

(45)

where $\frac{\partial C_{i,w}}{\partial t}$ is the variation of concentration with respect to the initial concentration of the considered substance. $Vo{l}_{FC}$ and $Vo{l}_{SSMC}$ are the volumes of a foam cell and a synthetic smooth muscle cell, respectively. Finally, ${\rho }_{Cg}$ is the collagen density.

To calculate the volume of foam cells, they have been approximated as spherical geometries, while synthetic smooth muscle cells are modelled as ellipsoids, so their volumes can be calculated with Equations (46) and (47).

$$Vo{l}_{FC}=\frac{4}{3}\pi {R_{FC}}^3$$

(46)

$$Vo{l}_{SSMC}=\frac{4}{3}\pi {R_{SSMC}}^2·l_{SSMC},$$

(47)

with $R_{FC}$ and $R_{SSMC}$ being foam cells and the synthetic smooth muscle cell radius, and $l_{SSMC}$ being the lengths of synthetic smooth muscle cells. Parameters for calculating plaque growth in the arterial wall are given in

Table 4. Parameters to calculate plaque growth in the arterial wall.

Plaque Growth Parameters
Parameter	Description	Value	Reference
$R_{FC}$	Foam cell radius	$15.264\mathsf{\mu }$m	Krombach et al. [64], Cannon and Swanson [65]
$R_{SSMC}$	SMC radius	$3.75\mathsf{\mu }$m	Cilla et al. [4]
$l_{SSMC}$	SMC length	$115\mathsf{\mu }$m	Cilla et al. [4]
${\rho }_{Cg}$	Collagen density	$1000\frac{kg}{m^3}$	Sáez et al. [66]

.

4. Numerical Methods

The software COMSOL Multiphysics (COMSOL AB, Burlington, MA, USA) was used to computationally solve the model. It was modelled using four consecutive steps, which can be seen in Figure 3. In the first step, three cardiac cycles and hemodynamic stimuli were calculated in a transient step, determining the values of $TAWSS$ and $OSI$ for the third cardiac cycle. Then, in the second step, which is stationary, the plasma flow through the endothelium is calculated. In the third step, the inflammatory process is calculated in a transient mode, determining the concentrations of all substances in the arterial wall for a total of 30 years. Finally, a last stationary step was developed to calculate the growth of the plaque, knowing the concentrations of all substances for a period of 30 years.

A direct solver (multi-frontal massively parallel sparse direct solver, MUMPS) was used to calculate transient blood flow along the lumen, as well as plasma flow across the endothelium, the inflammatory process in the arterial wall, and plaque growth. The inflammatory process was calculated iteratively, using groups of segregated steps for the different substances.

5. Sensitivity Analysis

As can be seen in Section 3, the computational model has a large number of parameters that can affect the composition and growth of the plaque. There is a wide range of values of these parameters in the literature that also come from studies under different conditions, such as in vivo or in vitro experiments, or different species or arteries. Therefore, a sensitivity analysis can help to understand the role of every parameter in the generated plaque.

The objective of this study is to analyse plaque growth and its composition. To do so, a previous selection of the parameters to analyse was made. Therefore, geometric parameters, initial concentrations, material properties, flow properties, and other factors that are well known have not been considered in the analysis because their influence on the model seems to be clear (for example, an increase in the radius of foam cells would result in an increase in the stenosis ratio due to an increase in the volume of the plaque). Thus, in this study, only the parameters related to the reactive terms of the convection–diffusion–reaction equations of the inflammatory process are analysed, as well as the diffusion coefficients of substances in the arterial wall.

Table 5. Analysed parameters and values.

Parameter	Description	Value −10%	Mean Value	Value +10%
$D_{LDL,w}=D_{LDLox,w}\left(\frac{m^2}{s}\right)$	LDL and oxLDL diff. coeff.	$7.2·{10}^{-13}$	$8·{10}^{-13}$	$8.8·{10}^{-13}$
$D_{m,w}=D_{M,w}(\frac{m^2}{s}$)	Monocyte and macrophage diff. coeff.	$7.2·{10}^{-15}$	$8·{10}^{-15}$	$8.8·{10}^{-15}$
$d_{LDL}\left(s^{-1}\right)$	LDL oxidation	$2.565·{10}^{-4}$	$2.85·{10}^{-4}$	$3.135·{10}^{-4}$
$d_m\left(s^{-1}\right)$	Monocyte differentiation	$1.035·{10}^{-6}$	$1.15·{10}^{-6}$	$1.265·{10}^{-6}$
$m_d\left(d^{-1}\right)$	Monocyte death	$0.015$	$0.016$	$0.018$
$LD{L}_{ox,r}\left(\frac{m^3}{Macrophage·s}\right)$	OxLDL uptake	$2.205·{10}^{-23}$	$2.45·{10}^{-23}$	$2.695·{10}^{-23}$
$n_{FC}\left(\frac{mo{l}_{oxLDL}}{Macrophage}\right)$	Max. oxLDL uptake	$2.448·{10}^{-11}$	$2.72·{10}^{-11}$	$2.992·{10}^{-11}$
$C_r\left(\frac{mo{l}_C·m^3}{mo{l}_{oxLDL}·Macrophage·s}\right)$	Cytokine production	$2.7·{10}^{-10}$	$3·{10}^{-10}$	$3.3·{10}^{-10}$
$d_c\left(s^{-1}\right)$	Cytokine degradation	$2.082·{10}^{-5}$	$2.314·{10}^{-5}$	$2.545·{10}^{-5}$
$S_r\left(d^{-1}\right)$	SMC differentiation	$3.24·{10}^{-3}$	$3.6·{10}^{-3}$	$3.96·{10}^{-3}$
$p_{ss}\left(d^{-1}\right)$	SSMC proliferation	$0.01818$	$0.0202$	$0.0199$
$G_r\left(\frac{k{g}_{Cg}}{SSMC·s}\right)$	Collagen production	$2.241·{10}^{-21}$	$2.49·{10}^{-21}$	$2.739·{10}^{-21}$
$C_{c,w}^{th}\left(\frac{mo{l}_C}{m^3}\right)$	Cytokine threshold	$1.111·{10}^{13}$	$1.235·{10}^{13}$	$1.358·{10}^{13}$
$m_r\left(\frac{m^4}{mo{l}_{oxLDL}·day}\right)$	Monocyte recruitment	$5.972·{10}^{-4}$	$6.636·{10}^{-4}$	$7.299·{10}^{-4}$
$r_{apop}\left(s^{-1}\right)$	SSMC apoptosis rate	$7.209·{10}^{-8}$	$8.011·{10}^{-8}$	$8.812·{10}^{-8}$

contains all the parameters whose variation was analysed and their descriptions. The analysed parameters have been increased and reduced by 25 % and 10 % in 61 different simulations in a mono-variant sensitivity analysis. The values of the parameters for the case of ± 10% variation are included in Table 5, as an example:

The percentage change in the volume of the plaque due to each of the substances involved in its growth (foam cells that compose the lipidic core of the plaque, and synthetic smooth muscle cells and collagen fibre, which correspond to the fibrous layer of the plaque) was analysed, as well as the variation in the stenosis ratio, which is defined as:

$$SR(\%)=\left(1-\frac{Areaoflumenwithplaque}{Areaofhealthylumen}\right)·100$$

(48)

6. Results

Table 6. Results of foam cells, synthetic smooth muscle cells, and collagen volume variations (second, third, and fourth columns, respectively), and stenosis ratio variation (fifth column), by reducing and increasing the values of the parameters of the first column by 25%.

Parameter	FC Volume Variation (%)		SSMC Volume Variation (%)		Cg Volume Variation (%)		SR Variation (%)
	−25%	+25%	−25%	+25%	−25%	+25%	−25%	+25%
$D_{LDL,w}=$ $D_{LDLox,w}$	7.35	−5.65	>100	−92.81	>100	−74.38	-	−36.61
$D_{m,w}$ $=D_{M,w}$	4.82	−4.07	>100	−43.43	>100	−39.57	-	−19.71
$d_{LDL}$	−20.29	18.51	−99.89	>100	−93.98	>100	−82.68	-
$d_m$	−0.98	0.63	−42.43	21.00	−31.32	15.33	−7.28	4.69
$m_d$	0.23	−1.09	15.32	−34.09	12.95	−29.32	7.90	−8.16
$LD{L}_{ox,r}$	−3.73	2.90	>100	−99.89	>100	−93.27	-	−12.99
$n_{FC}$	48.43	−36.00	−0.21	0.13	−0.10	0.24	>100	−77.22
$C_r$	0	0	−99.89	>100	−95.52	>100	−22.05	-
$d_c$	0	0	>100	−99.89	>100	−94.14	-	−21.95
$S_r$	0	0	−5.50	6.46	−7.48	7.71	−0.64	0.032
$p_{ss}$	0	0	−99.89	>100	−96.88	>100	−22.17	-
$G_r$	0	0	0	0	−24.54	48.32	−1.08	2.09
$C_{c,w}^{th}$	0	0	>100	−99.89	>100	−94.14	-	−21.96
$m_r$	−3.77	3.75	−41.23	>100	−35.02	>100	−20.83	-
$r_{apop}$	0	0	>100	−99.91	>100	−96.72	-	−22.19

and

Table 7. Results of the foam cells, synthetic smooth muscle cells, and collagen volume variations (second, third, and fourth columns, respectively), and variation in the stenosis ratio (fifth column), reducing and increasing the values of the parameters of the first column by 10%.

Parameter	FC Volume Variation (%)		SSMC Volume Variation (%)		Cg Volume Variation (%)		SR Variation (%)
	−10%	+10%	−10%	+10%	−10%	+10%	−10%	+10%
$D_{LDL,w}=$ $D_{LDLox,w}$	2.69	−2.42	>100	−67.85	>100	−48.60	-	−21.31
$D_{m,w}$ $=D_{M,w}$	2.32	−1.32	>100	−13.45	>100	−9.32	-	−11.23
$d_{LDL}$	−7.89	7.60	−97.67	>100	−81.81	>100	−45.91	-
$d_m$	−0.54	0.32	−25.54	12.55	−19.89	8.35	−2.67	2.10
$m_d$	0.05	−0.43	7.65	−25.43	4.56	−20.88	3.91	−3.06
$LD{L}_{ox,r}$	−1.36	1.23	>100	−97.85	>100	−82.55	−	−9.61
$n_{FC}$	15.50	−11.52	−0.07	0.04	−0.03	0.08	57.57	−26.82
$C_r$	0	0	−99.13	>100	−86.56	>100	−21.31	-
$d_c$	0	0	>100	−98.69	>100	−84.93	-	−21.16
$S_r$	0	0	−1.17	1.38	−1.59	1.64	−0.13	0.007
$p_{ss}$	0	0	−99.89	>100	−91.40	>100	−21.75	-
$G_r$	0	0	0	0	−4.50	8.77	−0.19	0.36
$C_{c,w}^{th}$	0	0	>100	−98.69	>100	−84.93	-	−21.15
$m_r$	−0.26	0.26	−2.89	>100	−2.45	67.59	−1.49	-
$r_{apop}$	0	0	>100	−99.90	>100	−91.14	-	−21.77

show the results obtained by increasing and decreasing the selected parameters by 25% and 10%, respectively.

The first column of Table 6 and Table 7 represents the parameter whose influence is being analysed. The second, third, and fourth double columns are, respectively, changes in the volume of foam cells, synthetic smooth muscle cells, and collagen fibre in the plaque, caused by the considered parameter variation. Finally, the last double column is the change in the stenosis ratio of the artery due to the change in the considered parameter.

As can be seen in Table 6 and Table 7, the trend of the results is the same in the cases of variation of 25% and 10% parameters variation. Therefore, the results will be discussed only with reference to the 10% variation table (Table 7), and can be extrapolated for the 25% variation table (Table 6).

As can be seen in Table 6 and Table 7, the variation of the substance was limited to a maximum of 100%. Therefore, in cases where the variation in the percentage of a substance was greater than 100%, the stenosis ratio was not calculated.

As can be seen in Table 7, an increase in the diffusion parameters, $m_d$, $LD{L}_{ox,r}$, $n_{FC}$, $d_c$, $c_{c,w}^{th}$, and $r_{apop}$ causes a decrease in the stenosis ratio produced by the plaque. In contrast, an increase in $d_{LDL}$, $d_m$, $C_r$, $S_r$, $p_{ss}$, $G_r$, and $m_r$ induces an increase in it. In addition, there are some parameters of the model that have more influence on the results than others, and their variation causes a change greater than 100% in the volume of any of the substances in the plaque.

When considering the parameters that influence the change in the volume of foam cells in the plaque (Table 7), these are, in order of influence, $n_{FC}$, $d_{LDL}$, $D_{LDL,w}=D_{LDLox,w}$ and $D_{m,w}=D_{M,w}$, which are related to foam cells, LDL, and the diffusion properties of substances in the arterial wall, respectively. As can be seen, none of the variations produces a change in the volume of foam cells greater than 100%.

The parameters that cause a higher change in the volume of synthetic smooth muscle cells are, in order of influence: $d_{LDL}$, $C_r$, $p_{ss}$, $m_r$, $r_{apop}$, $d_c$, $C_{c,w}^{th}$, and $LD{L}_{ox,r}$, when increased (the first four produce a change greater than 100%). When their values decrease, the most influential are, in order: $D_{LDL,w}=D_{LDLox,w}$, $D_{m,w}=D_{M,w}$, $LD{L}_{ox,r}$, $d_c$, $C_{c,w}^{th}$, $r_{apop}$, $p_{ss}$, $C_r$, and $d_{LDL}$, the first six of which cause changes greater than 100%. Therefore, for the change in the volume of synthetic smooth muscle cells volume in the plaque, the parameters $d_{LDL}$, $C_r$, $p_{ss}$, $r_{apop}$, $d_c$, $C_{c,w}^{th}$, and $LD{L}_{ox,r}$ have huge influences, regardless if their values are increased or decreased.

For the case of the influence on volume change due to SSMC, $d_{LDL}$, $C_r$, $p_{ss}$, $r_{apop}$, $d_c$, $C_{c,w}^{th}$ and $LD{L}_{ox,r}$ have a great influence on the results obtained.

For the variation of collagen volume in the plaque, the most influential parameters are, when increased: $d_{LDL}$, $C_r$, $p_{ss}$, $r_{apop}$, $d_c$, $C_{c,w}^{th}$ and $LD{L}_{ox,r}$, while when decreased: $D_{LDL,w}=D_{LDLox,w}$, $D_{m,w}=D_{M,w}$, $LD{L}_{ox,r}$, $d_c$, $C_{c,w}^{th}$, and $r_{apop}$.

In Figure 4, the variation in the volume of foam cells, synthetic smooth muscle cells, and collagen fibre is represented in a graph of parallel bars for variations of $\pm 10\%$. In cases of parameters that cause a volume variation higher than 100% in any of the considered substances, the bar of this substance is represented in red, i.e., the cases of the diffusion coefficients $d_{LDL}$, $LD{L}_{ox,r}$, $C_r$, $d_c$, $p_{ss}$, $c_{c,w}^{th}$, $m_r$, and $r_{apop}$.

Figure 5 represents the change in the stenosis ratio due to the variation of the considered parameters when varied ±10%. The red bars refer to cases in which at least one of the substances that adds volume to the plaque has a volume variation greater than 100%.

7. Discussion

In this work, an analysis of the influence of some parameters of a previous mathematical model of atherosclerosis development in arteries was performed [8]. The mathematical model has a large number of parameters that can affect the growth of the plaque. However, some of them are considered well-known due to, for example, corresponding to geometrical properties of arteries or substances. Therefore, the parameters whose influence on plaque growth was analysed are related to the reactive terms of the equations referred to substances in the arterial wall. These parameters have been modified by increasing and decreasing their value in different cases by $10\%$ and, to determine how they affect plaque growth, variations in the volume of substances that add volume to the plaque have been analysed (foam cells, synthetic smooth muscle cells, and collagen fibre). In addition, the variation in the plaque stenosis ratio was analysed.

As can be seen in the results, on the one hand, the parameters whose variations are directly proportional to the stenosis ratio are the degradation rate of LDL ($d_{LDL}$), the monocyte differentiation rate ($d_m$), the parameters referring to the production and degradation of cytokines ($C_r$ and $S_r$), the proliferation rate of synthetic smooth muscle cells ($p_{ss}$), the segregation rate of collagen ($G_r$) and the parameter related to the amount of monocyte recruited by the endothelium ($m_r$). On the other hand, some of the analysed parameters are inversely proportional to the growth of the plaque, and an increase in their values will cause a decrease in the volume of the plaque and, therefore, in the stenosis ratio. It is the case of the diffusion parameters of substances on the arterial wall ($D_{LDL,w}=D_{LDLox,w}$ and $D_{m,w}=D_{M,w}$), the rate of death of monocytes ($m_d$), the rate of oxidised LDL that is uptaken by macrophages, and the maximum amount of oxidised LDL that a macrophage can ingest ($LD{L}_{ox,r}$ and $n_{FC}$), the cytokine degradation rate, its threshold in the arterial wall ($d_c$ and $C_{c,w}^{th}$), and the rate of apoptosis of synthetic smooth muscle cells ($r_{apop}$).

As can be noticed, the parameters that influence the change in the volume of synthetic smooth muscle cells also influence the change in the volume of collagen. It is due to the segregation of collagen fibre by SSMC, so collagen depends directly on them.

The less influential parameters in the volume change of substances in the plaque are: $d_m$, $m_d$, $S_r$, and $G_r$. It should be noted that $d_m$ and $m_d$ are parameters referring to monocytes, which act at the beginning of the process. Therefore, a great influence on them could be expected. However, monocytes highly affect the results with the parameter $m_r$, which is the monocytes recruitment from the lumen, and its variation has a huge influence on the volume of FC in the plaque and, therefore, in the stenosis ratio.

As can be seen in Figure 4, $r_{apop}$ has a large influence on the variation of volume of synthetic smooth muscle cells and collagen for both cases, when increasing and decreasing its value by 10%. The smaller the $r_{apop}$ value, the more plaque is generated, as it is an apoptosis factor of synthetic smooth muscle cells. However, its influence on the variation of the stenosis ratio is greater in the case of increasing its value than in the case of decreasing it. As it is a parameter related to the apoptosis ratio of synthetic smooth muscle cells, its change does not cause variation in the results of foam cells (Figure 5).

$m_r$ has more influence in both concentrations and stenosis ratio variation when increased (Figure 4 and Figure 5). This is because, when its value is decreased, the amount of monocytes deposited in the arterial wall is reduced.

$C_{c,w}^{th}$ has a large influence in the variations of the results, having more influence when increased (the change of volume and stenosis ratio are greater than 100% in this case, as can be seen in Figure 4 and Figure 5). It is a parameter involved in the differentiation of contractile smooth muscle cells into synthetic ones due to the presence of cytokines in the arterial wall. Thus, it does not influence the volume of foam cells.

As $G_r$ is a parameter of collagen fibre segregation, it only has influence on the change of the volume of collagen in the plaque (Figure 4). Therefore, its influence on the stenosis ratio is limited (Figure 5).

$p_{ss}$ is related to the proliferation of synthetic smooth muscle cells, so it does not influence the results of foam cells. On the contrary, as can be seen in Figure 4, increasing it by 10% produces a change greater than 100% in the variation of synthetic smooth muscle cells and collagen fibre.

$S_r$ changes do not cause a large variation in the volume of any substance or the plaque stenosis ratio (Figure 5).

$d_c$ is a parameter that also has a large influence on the volume variations (Figure 4). This parameter represents the cytokine degradation; thus, the higher it is, the more cytokines are degraded and, thus, the volume and the ratio of stenosis are lower (Figure 5). The same occurs with $C_r$, which represents the cytokine production.

$n_{FC}$ represents the maximum amount of oxidised LDL that a macrophage can ingest before becoming a foam cell. Therefore, an increase in its value produces a reduction in the volume of foam cells and the stenosis ratio (Figure 4 and Figure 5). Its influence on the variation of the volume of the substances is not very large, but it produces an important variation in the stenosis ratio. When these results are contrasted with those of a substance that produces a large variation in synthetic smooth muscle cells and collagen volumes (for example, $d_c$), it can be observed that a smaller change in the volume of foam cells produces a larger change in the stenosis ratio. It can be explained by attending to Equation (45): The volume of a foam cell is equal to $1.489·{10}^{-14}{m}^3$, while the volume of a synthetic smooth muscle cell is equal to $6.774·{10}^{-15}{m}^3$. Therefore, due to their size difference, less change in foam cell volume is needed to produce a similar stenosis ratio variation.

$LD{L}_{ox,r}$ is related to the oxidised LDL uptaken by macrophages, so it affects the volume of each of the considered substances. An increase in its value produces a reduction in the volume of substances (Figure 4) and therefore of the stenosis ratio (Figure 5).

$d_m$ and $m_d$ are both parameters referring to monocytes. The first one is related to their differentiation, while the second one refers to their apoptosis. Therefore, their influences are opposite. Their influence is more notable for synthetic smooth muscle cells and collagen volumes (Figure 4).

$d_{LDL}$ is the degradation rate of LDL, so it has an influence on foam cells, synthetic smooth muscle cells, and collagen fibre and, therefore, in the stenosis ratio of the plaque (Figure 4 and Figure 5). So, it is one of the most influential parameters of the model and the most influential in the stenosis ratio when it is reduced.

$D_{LDL,w}=D_{LDLox,w}$ and $D_{m,w}=D_{M,w}$ are related to the diffusion of substances in the arterial wall, so they affect all the processes. Therefore, they influence the results of all the substances, and are some of the most important parameters in the model (Figure 4 and Figure 5).

With all of this information, knowing the influence of all the parameters, they could be adjusted to achieve more or less vulnerable atheroma plaque, according to the percentage volume of foam cells, synthetic smooth muscle cells, and collagen fibre [67,68]. The vulnerability of a plaque is dependent on multiple factors, such as its size or stresses caused by blood flow in it, but it is also dependent on its composition. There is a high risk of rupture of plaque with a large lipidic nucleus and a thin fibrous cap [67,69]. Therefore, a high quantity of foam cells and a small amount of synthetic smooth muscle cells and collagen fibre will be indicators of plaque with a high risk of rupture (and, therefore, less stable) than one with a large quantity of synthetic smooth muscle cells and collagen [67,68,69,70]. Therefore, reducing the maximum amount of oxidised LDL that a macrophage can ingest and the ratio of oxidation of LDL ($n_{FC}$) will cause plaque with bigger lipid nuclei, which can develop into more unstable plaque. However, as can be seen, it has no influence on SSMC and collagen volumes in the plaque. Conversely, increasing the apoptosis ratio of SSMCs, the cytokine threshold in the arterial wall and its degradation rate, and the rate of oxidised LDL uptaken by macrophages ($r_{apop}$, $c_{c,w}^{th}$, $d_c$ and $LD{L}_{ox,r}$, respectively), and reducing SSMC proliferation, cytokine production, and the oxidation LDL ratio ($p_{ss}$, $C_r$, and $d_{L}DL$, respectively) will produce plaque with less fibrotic layer and, thus, a high risk of rupture.

Some studies in the literature focus on determining the most influential parameters in the development of atheroma plaque in different mathematical models. It is the case of Cilla et al. [21], in which the effect of the anisotropy of the arterial wall on the diffusion coefficients of LDL was analysed. Their results have been considered to implement the anisotropy of the diffusion coefficients in the present model. There are also some studies on parameter influence in agent-based models. It is the case of Olivares et al. [18], who focus on determining how the migration of agents, the velocity of oxidation of LDL, and the maximum amount of autoantibodies can affect the plaque. In addition Corti et al. [20] have mono-parametric and multi-parametric studies to determine the influence of the parameters on their model. However, it is not possible to contrast their results with the ones obtained in this article, as each of the models considers different substances in the process of atheroma plaque formation and, therefore, their parameters are referred to other substance values.

The findings of this study should be interpreted in the context of its limitations. For example, the study of the influence of parameters was done only in the geometry of the carotid artery. However, the behaviour of the mathematical model would be the same for other geometrical configurations and arteries, for example, coronary or aorta arteries, adapting the values of the corresponding parameters if necessary. Another limitation of the study is that it was done with a 2D-axisymmetric model instead of a real geometry. This produces a higher plaque growth than in real cases as diffusion in the circumferential direction is not allowed and therefore causes a higher accumulation of substances in the area of the plaque. However, it should not influence the development of the plaque according to the variations in the parameters that were analysed here. In addition, blood flow and the inflammatory process are not coupled, which could influence the shape and stenosis ratio of the developed plaque. In this study, we also do not consider the influence of the mechanics of the arterial wall in the development of the plaque (such as tortuosity or changes in the permeability of the arterial wall due to the thickness variation of it with the cardiac cycle).

8. Conclusions

The mathematical model has a large number of parameters and their values influence the plaque obtained. They can affect the volume of the substances that provide volume to the plaque and also its stenosis ratio. Therefore, the dependence of the model on the variation of its parameters was analysed. For that, a previous selection of the parameters to analyse was done, and those referred to geometrical parameters, initial concentrations, and material properties were discarded.

As can be seen in Section 6, the variation of the selected parameters carries important variations on the results, and in some cases, this variation can be greater than 100%.

$d_{LDL}$, $d_m$, $C_r$, $S_r$, $p_{ss}$, $G_r$, and $m_r$ are directly proportional to the change of substances volume (FC, SSMC, and collagen) and to the stenosis ratio, while $D_{LDL,w}=D_{LDLox,w}$, $D_{m,w}=D_{M,w}$, $m_d$, $LD{L}_{ox,r}$, $n_{FC}$, $d_c$, $C_{c,w}^{th}$, and $r_{apop}$ are inversely proportional.

In addition, it was noticed that a variation in foam cell volume results in more of a change in the plaque stenosis ratio than in the volume of synthetic smooth muscle cells or collagen, due to the larger volume.

For all of this, it could be interesting to study the vulnerability of plaque by changing the analysed parameters, knowing how each one of them affects the volume of foam cells, synthetic smooth muscle cells, and collagen fibre in the plaque.