# Multiscale bio-chemo-mechanical model of intimal hyperplasia

Jérôme Jansen, Xavier Escriva, Fabien Godeferd, Patrick Feugier

► To cite this version:

Jérôme Jansen, Xavier Escriva, Fabien Godeferd, Patrick Feugier. Multiscale bio-chemo-mechanical model of intimal hyperplasia. Biomechanics and Modeling in Mechanobiology, 2022, 21 (2), pp.709-734. 10.1007/s10237-022-01558-5. hal-03712935

HAL Id: hal-03712935

https://hal.science/hal-03712935

Submitted on 4 Jul 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Multiscale bio-chemo-mechanical model of intimal hyperplasia

Jérôme Jansen · Xavier Escriva · Fabien Godeferd · Patrick Feugier

Received July 9, 2021 / Accepted January 6, 2022 / Published January 29, 2022

**Abstract** We consider a computational multiscale framework of a bio-chemo-mechanical model for intimal hyperplasia. With respect to existing models, we investigate the interactions between hemodynamics, cellular dynamics and biochemistry on the development of the pathology. Within the arterial wall, we propose a mathematical model consisting of kinetic differential equations for key vascular cell types, collagen and growth factors. The luminal hemodynamics is modelled with the Navier-Stokes equations. Coupling hypothesis among time and space scales are proposed to build a tractable modelling of such a complex multifactorial and multiscale pathology. A one-dimensional numerical test-case is presented for validation by comparing the results of the framework with experiments at short and long timescales. Our model permits to capture many cellular phenomena which have a central role in the physiopathology of intimal hyperplasia. Results are quantitatively and qualitatively consistent with experimental findings at both short and long timescales.

**Keywords** Intimal hyperplasia · Growth & re-modeling · Hemodynamics · Multiscale modelling · Biochemistry

**Declarations:**

**Funding:** No funding was received to assist with the preparation of this manuscript.

**Competing interests:** The authors have no conflicts of interest to declare that are relevant to the content of this article.

**Code availability:** All simulation has been achieved with in-house Python packages. https://gitlab.com/jrme.jansen/pytg/

**Availability of data and material:** Python codes used in this study are available at https://gitlab.com/jrme.jansen/pytg/.

**Authors' contributions:** X. Escriva and P. Feugier set up the relevant framework for tackling the modelling of the physiological-mechanical pathology. J. Jansen developed the model, carried out simulations of the test-case and wrote the manuscript. F. Godeferd, X. Escriva and P. Feugier participated in the discussion of the results. J. Jansen and F. Godeferd wrote the first draft of manuscript.

**Ethics approval:** This article does not present research with ethical considerations.

**Consent to participate:** All the authors consented to participate in this study.

**Consent for publication:** All authors proof-read and consent for the publication of the manuscript.

<footnote>J. Jansen · X. Escriva · F.S. Godeferd
Laboratoire de Mécanique des Fluides et d'Acoustique, UMR 5509, Univ Lyon, École Centrale de Lyon, INSA Lyon, Université Claude Bernard Lyon I, CNRS, 36 Avenue Guy de Collongue, 69134 Ecully, France
E-mail: jrme.jansen@gmail.com</footnote>

<footnote>P. Feugier
Service de Chirurgie Vasculaire et Endovasculaire, LYVES
Groupement Hospitalo-Universitaire Lyon Sud, Université
Claude Bernard Lyon 1, 69100, Villeurbanne, France</footnote>
# 1 Introduction

Cardiovascular diseases (CVDs) are the leading cause of global mortality worldwide and are a burden for health systems. Among them, there are vascular diseases (VDs) as atherosclerosis, hypertension, intimal hyperplasia (IH), aneurysms or endofibrosis. VDs may cause opposing arterial morphological changes (narrowing versus enlargement of lumen) and develop according to very different risk factors, e.g. intensive sport activities induce endofibrosis. However, they are all due to a deregulation of the homeostasis of the arterial wall which activates several same types of cellular events (migration, deposit or loss of collagen, hyperplasia). It is widely accepted that the cellular processes described in the VDs developments are modulated by mechanically-induced (Wentzel et al., 2003; Humphrey, 2002) and biochemically-induced (Ducasse et al., 2003) mediators. The biochemical dynamics in VDs have been studied intensively which brought to light the importance of several types of regulatory molecules, or so called growth factors (GFs). The main cellular populations acting in the physiopathology of VDs are endothelial cells (ECs), vascular smooth muscle cells (vSMCs) and various types of immune cells (Clowes, 1993). The VDs features are so similar that links have been proposed between certain pathologies as for atherosclerosis and IH at short-times (Goodman et al., 2016) or endofibrosis and IH (Feugier and Rouvière, 2018).

A clean understanding of the interactions of vascular cells, hemodynamics and GFs is essential to completely understand the mechanisms that control vascular cell growth and pathological onset. Over the last two decades, to address this issue, in addition to experimental studies, numerous mathematical models have been proposed for multiple vascular pathologies. This effort was undertaken with the aim of offering an *in silico* approach as a new way of experimentation to address a number of unresolved questions about physiopathology and even to identify new therapeutic targets. In the following, we summarize the state of the art in modelling cellular and biochemical dynamics coupled with hemodynamics in the context of VDs. Partial differential equations (PDEs) of reaction-diffusion type have been used to model the spatio-temporal lesion evolution for atherosclerosis (Calvez et al., 2010; Goodman et al., 2016) coupled to the local dynamics of blood flow modelled as Navier-Stokes equations. Test-cases was used to investigate the evolution of the lesion in pro-atherogenic flow regions. Ordinary differential equations (ODEs) and delay differential equations (DDEs) of kinetic equations type have been used to develop models of tissue growth time evolution for atherosclerosis (Bulelzai and

Dubbeldam, 2012), hypertension (Wilstein et al., 2018), intimal hyperplasia (Schwartz et al., 1996; Donadoni et al., 2017), tissue-engineered vascular grafts (Khosravi et al., 2020) and wound healing inflammation (Nagaraja et al., 2014).

Given the previous modelling efforts, we propose a multiscale and multiphysics model of IH. For a complete review of IH physiopathology, we refer the reader to the works of Kenagy (2011); Ducasse et al. (2003); Clowes (1993); Model and Dardik (2012). Starting from Donadoni et al. (2017), the models are written as a set of ODEs and DDEs within the intima and media layers of the arterial wall. In this compartmental approach, spatial dependency are implicitly expressed as intima and media layers inside which we model cellular and biochemical dynamics. We aim to model the main cell and biochemical mechanisms considering a coupling with hemodynamics within a tractable model. The challenge of proposing tractable models is to simplify as much as possible while pathologies are sums of cascades of complex biological and physical interactions. In this process of simplification, the kinetics of the species and growth factors modelled in this article will be assumed to be of first order as for collagen turnover (Humphrey, 2002). We believe that proposing and testing a ODEs and DDEs model is essential before any increase in model complexity as explicitly integrated spatial dependencies via PDEs.

The article is organized according to the following plan. In Section 2.1, we present the hemodynamical model. Section 2.2 derives equations of arterial wall species for ECs, vSMCs and collagen fibers within the intima and media layers. Section 2.3 is devoted to the biochemical model of bioavailability. Section 2.4 presents the coupling methods between arterial wall species, biochemistry and hemodynamics. Section 3 presents a one-dimensional test-case of the present model aimed at mimicking the conditions of an endothelial denudation. As IH has been extensively studied through animal models, *in vitro* and *in vivo* experiments, experimental data will be used for model validations. In the associated Sections: 3.1, we initialize the test-case, 3.2, we compare results of the simulation at short- and long-time to experimental data available by aiming for validation of the model, 3.3, we discuss results and present a sort sensitivity analysis of important parameters of the model, 3.4, we discuss the possible paths of model improvement. Conclusion is in the final Section 4.
## 2 A multiscale and multiphysics framework

### 2.1 Hemodynamics modelling

The mechanical contents of our model corresponds to the description of hemodynamics, since hemodynamic forces are known to modulate pathological development (Wentzel et al., 2003). In the case of IH, blood flow is described as “best associated with the formation of intimal hyperplasia” (Davies, 2019). We consider a simplified description of arterial hemodynamics as a stationary laminar flow in a rigid duct, in which blood is assumed to be an incompressible Newtonian fluid. Under these hypotheses, hemodynamics is modelled by the Navier-Stokes equations

$$ \frac{\partial \mathbf{u}}{\partial t} + (\mathbf{u} \cdot \nabla) \mathbf{u} = -\frac{1}{\rho} \nabla p + \nu \Delta \mathbf{u} \quad \text{in } \Omega_1, \quad (1a) $$

$$ \nabla \cdot \mathbf{u} = 0 \quad \text{in } \Omega_1, \quad (1b) $$

where $\mathbf{u}$ is the velocity vector, $p$ the pressure, $\rho$ the blood density, $\nu$ the kinematic viscosity and $\mu = \rho\nu$ the dynamical viscosity. The wall shear stress (WSS) vector $\tau_w$ and its magnitude are

$$ \tau_w = \mu (\nabla \mathbf{u} + \nabla \mathbf{u}^T) \cdot \mathbf{n} \quad \text{in } \Gamma_e, \quad (2a) $$

$$ \tau_w = |\tau_w| \quad \text{in } \Gamma_e, \quad (2b) $$

where $\mathbf{n}$ is the unit vector normal to endothelial surface $\Gamma_e$. The WSS represents the tangential friction stress exercised by blood on the endothelium, a surface lined by ECs. Spatial domains of application of equations are shown in Fig. 1.

Assuming the artery as a rigid duct, we ignore the coupling effects between tissue growth-induced artery deformation and hemodynamics. This hypothesis is based on the timescale separation between tissue growth and hemodynamics discussed in Section 2.4. Thus, we consider the arterial wall deformation as a succession of quasi-static arterial configurations at equilibrium called arterial generation. For each arterial generation, we apply the Navier-Stokes equations (1) and consider only the new WSS stimuli modulating the tissue growth.

In addition, we ignore other mechanical stimuli on the present model, such as circumferential/axial intraparietal stresses, in order to focus on the influence of hemodynamics on IH. This assumption ignores the macroscopic stress/strain relation of the arterial wall material, which is theoretically dependent on the structural properties of the wall modulated by the tissue growth Humphrey (2002).

The timescale separation allows to overcome the additional complexity induced by the fluid-structure interactions, in the aspects of modelling and computational

<figure><img src="image_2.png" /><figcaption>Fig. 1 Idealized artery of length $\Delta z$. $\Omega_1$ is the luminal volume within the luminal radius $R_1$. $\Gamma_e$ is the endothelial surface, $\Omega_i$ is the intima layer of thickness $e_i$ between $R_1$ and $R_{IEL}$, the radius of internal elastic lamina. $\Omega_m$ is the media layer of thickness $e_m$ between $R_{IEL}$ and $R_{EEL}$, the radius of external elastic lamina. $\Omega_a$ is the adventitia layer of thickness $e_a$ between $R_{EEL}$ and $R_{ext}$, the external radius of the artery. The whole arterial wall domain is noted $\Omega_w$</figcaption></figure>

cost. By considering these restrictive hypotheses, we aim to experiment numerically the influence of WSS stimuli on the development of intimal hyperplasia.

The two constant independent parameters of the hemodynamical model are summarized in Table 1.

### 2.2 Arterial species modelling

In this section, we obtain the set of coupled equations describing the evolution of cellular and non-cellular species within the arterial wall, as seven ODEs and two DDEs. The term species is used to refer to materials constituting the vessel wall. These equations are designed to model a pathology which “... results from the net increase in cell number and ECM, which are dependent on rates of cellular migration, proliferation, and death and on rates of ECM synthesis and degradation...” (Kegagy, 2011). We therefore model the dynamical evolution of cell’s main populations and extracellular matrix (ECM) within the intima and media layers in the following three sets of equations:

1. Firstly, ECs are taken into consideration because of their primordial role in vascular homeostasis and pathologies developments (Chiu and Chien, 2011). ECs evolution within the monolayer called endothelium is described from Ultman et al. (2016) as the
logistic equation

$$ \frac{dE}{dt} = r_E E \left( 1 - \frac{E}{E_{\max}} \right) \quad \text{in } \Gamma_e, \qquad (3) $$

where E is the count of endothelial cells, $E_{\max}$ the physiological number of ECs to have a confluent monolayer and $r_E$ the proliferation rate of cells.

2. vSMCs have a main role in arterial wall homeostasis through their phenotypical heterogeneity as contractile (cSMCs) or secretory (sSMCs) (Berk, 2001). From Donadoni et al. (2017), we modify the four coupled equations for vSMCs dynamics within the intima and media. Subscripts “i” and “m” are used to refer to variables related to the intima and media layers and domains of the layers are denoted $\Omega_i$ and $\Omega_m$. The equations for cells populations cSMCs and sSMCs are species functional properties (SFPs) summations with positive sign for production and negative sign for removal as follows.

- The evolution of cSMCs within intima and media is modelled with the following modified logistic equations

$$ \frac{dQ_i}{dt} = r_Q Q_i \left( 1 - \frac{Q_i}{Q_{i \max}} \right) - c_i Q_i + m_i Q_m \quad \text{in } \Omega_i, \qquad (4a) $$

$$ \frac{dQ_m}{dt} = r_Q Q_m \left( 1 - \frac{Q_m}{Q_{m \max}} \right) - c_m Q_m - m_i Q_m \quad \text{in } \Omega_m. \qquad (4b) $$

$Q_{i,m}$ is the count of cSMCs, and $Q_{i,m \max}$ is the physiological amounts of cSMCs in layers needed for vascular homeostasis. The proliferation rate of cells is $r_Q$, and $c_{i,m}$ are dedifferentiation rates from cSMCs to sSMCs. Finally, $m_i$ is a migration rate of cSMCs from media towards intima. The rates $c_{i,m}$ and $m_i$ are bioavailability-dependent growth factors (GFs), detailed in Sections 2.3 and 2.4.

- Equations for sSMCs are summations of cSMCs dedifferentiation, proliferation, apoptosis, and migration as

$$ \frac{dS_i}{dt} = c_i Q_i + (p_i - a_i)(S_i - S_i^{\text{eq}}) + m_i S_m \quad \text{in } \Omega_i, \qquad (5a) $$

$$ \frac{dS_m}{dt} = c_m Q_m + (p_m - a_m)(S_m - S_m^{\text{eq}}) - m_i S_m \quad \text{in } \Omega_m, \qquad (5b) $$

with $c_{i,m}$ and $m_i$ used in (4). $p_{i,m}$ are cell proliferation rates, $a_{i,m}$ cell apoptosis rates and $S_{i,m}^{\text{eq}}$ the equilibrium counts of sSMCs in each layer defined in equation (11). As for the rates $c_{i,m}$ and

$m_i$, $p_{i,m}$ and $a_{i,m}$ are bioavailability-dependent growth factors (see Sections 2.3 and 2.4).

3. The evolution of the ECM contents is globally considered as the evolution of the amount of collagen denoted as $C_{y_{i,m}}$ (Donadoni et al., 2017). It is described as a balance between $S_{i,m}$ production, collagen self-degradation, and collagen ageing terms $s_{i,m}$ defined in equations (7), as

$$ \frac{dC_{y_{i,m}}}{dt} = \lambda_{i,m} S_{i,m} - \chi_{i,m} C_{y_{i,m}} - s_{i,m} \quad \text{in } \Omega_{i,m}. \qquad (6) $$

$\lambda_{i,m}$ are the collagen production rates, $\chi_{i,m}$ the self-degradation rates and $s_{i,m}$ the collagen ageing rates. We do not assume cSMCs having the ability to synthesize collagen fibers, as in Donadoni et al. (2017). As in equations (4) and (5), functional properties related to collagen dynamics are bioavailability-dependent. Feugier and Rouvière (2018); Ojha et al. (2000) describe during histological analysis ageing of lesions in different vascular pathologies. Collagen fibers would have the ability to age without being degraded, or switch from normal/young to old fibers. The old fibers have the ability to remain within tissues without self-degradation at very long time (Feugier and Rouvière, 2018). To model this ageing process, two DDEs are used to describe collagen ageing evolution

$$ \frac{dC_{o_{i,m}}}{dt} = s_{i,m} \quad \text{in } \Omega_{i,m}, \qquad (7) $$

with $s_{i,m}$ the collagen ageing rates and $C_{o_{i,m}}$ old collagen amounts. The rates are derived from the temporal derivation of the amount of old collagen fibers which is assumed as

$$ C_{o_{i,m}} = \frac{1}{\tau_C} \int_{-\infty}^{t} \Phi \left[ C_{y_{i,m}}(u), C_{y \ i,m}^{\text{eq}}(u), C_{y_{i,m}}(u - \tau_C), C_{y \ i,m}^{\text{eq}}(u - \tau_C) \right] du \quad \text{in } \Omega_{i,m}, \qquad (8) $$

with $\tau_C$ collagen ageing delay, $\Phi$ a threshold function (introduced in Section 2.4.2), $C_{y_{i,m}}(t-\tau_C)$ past value of $C_{y_{i,m}}$ at $t-\tau_C$ and $C_{y \ i,m}^{\text{eq}}$ collagen equilibrium amounts. Equation (8) models that, if during $\tau_C$ an excessive amount of collagen is not regulated by the normal balances (6), then it will be expressed as old collagen. We assume that old collagen fibers will remain in layers without self-degradation. Temporal derivation of equation (8) gives expressions for collagen ageing rates as

$$ \frac{dC_{o_{i,m}}}{dt} = s_{i,m} = \frac{1}{\tau_C} \Phi \left[ C_{y_{i,m}}(t), C_{y \ i,m}^{\text{eq}}(t), C_{y_{i,m}}(t - \tau_C), C_{y \ i,m}^{\text{eq}}(t - \tau_C) \right] \quad \text{in } \Omega_{i,m}. \qquad (9) $$
The equilibrium quantities introduced in (5) and (9) are the result of the following consideration. The cellular populations of arterial wall act to maintain its physiological functions as the conservation of the wall mechanical properties (stiffness, elasticity). Without explicit account for mechanical modelling of the arterial wall, we suppose that physiological mechanical properties are related to the concentration of collagen within layers, noted $c_{C_{i,m}}^{ph} = \rho_c \varphi_{i,m}^{C_{i,m}}$, with $\varphi_{i,m}^{C_{i,m}}$ the physiological volume fraction (VF) of collagen within layers and $\rho_c$ the collagen density.

We assume that this conservation goes through maintaining a certain amount of collagen within layers, noted $C_{y \ i,m}^{eq}$. Once this goal is achieved, species of the wall should undergo equilibrium. We define the equilibrium amount of collagen as

$$C_{y \ i,m}^{eq} = c_{C_{i,m}}^{ph} V_{i,m}, \qquad (10)$$

with $V_{i,m}$ the current volumes of the layers defined in equation (12). The amount of sSMCs, noted $S_{i,m}^{eq}$, work to compensate collagen self-degradation within layers $\chi_{i,m}C_{y_{i,m}}$ defined in (6) and restore a preferred collagen concentration. This equilibrium count of sSMCs is then defined as

$$S_{i,m}^{eq} = \frac{\chi_{i,m}}{\lambda_{i,m}} C_{y \ i,m}^{eq}. \qquad (11)$$

Relations (10) and (11) define a couple of equilibrium variables $(C_{y \ i,m}^{eq}, S_{i,m}^{eq})$ which has a dual role. The couple tends to restore in the arterial wall a physiological concentration of collagen at any time and it allows the equations to find a state of equilibrium when $S_{i,m} = S_{i,m}^{eq}$ independently of the bioavailability-dependent sSMCs turnover $p - a$.

Finally, the set of nine differential equations (3), (4), (5), (6) (7) model the dynamics of species in layers. Under physiological conditions, which will be clarified in the test-case in Section 3, equations of species remain under a stable equilibrium unlike the previous model of Donadoni et al. (2017) (under the conditions that $p_{i,m} \le a_{i,m}$ and $t_{c,m} > 1$). This equilibrium state model vascular homeostasis as a static process although it is truly dynamic (Widmaier et al., 2016).

According to Donadoni et al. (2017) and assuming constant species densities, the evolution of volumes in layers is related to species dynamics as

$$V_i = \frac{E}{\rho_E} + \frac{S_i + Q_i}{\rho_s} + \frac{C_{y_i} + C_{o_i}}{\rho_c} + V_i^o \quad \text{in } \Omega_i, \qquad (12a)$$

$$V_m = \frac{S_m + Q_m}{\rho_s} + \frac{C_{y_m} + C_{o_m}}{\rho_c} + V_m^o \quad \text{in } \Omega_m, \qquad (12b)$$

$$V_w = V_i + V_m + V_a \quad \text{in } \Omega_w, \qquad (12c)$$

with $\rho_s$ cell density, $\rho_E$ ECs density, $\rho_c$ collagen density, $V_{i,m}^o$ constant volumes of other not modelled species within intima and media, $V_a$ constant volume of the adventitia and $V_w$ arterial wall volume. The ECs density is calculated by assuming its shape as ellipsoidal and it reads

$$\rho_E = 1 / \left( \frac{4}{3} \pi \frac{w_E}{2} \frac{l_E}{2} \frac{e_E}{2} \right), \qquad (13)$$

with $w_E$ and $l_E$ respectively the width and length of EC and $e_E$ the thickness of an endothelium.

Assumption of constant species densities neglects the hypertrophy capacity of cells. The size of cells do not change over time and volume variations of arterial layers are only due to variations in the amount of collagen and count of cells. When cells hypertrophy, the constant density hypothesis can lead to errors (Binder and Simpson, 2016), but in first approximation, we thus model the development of arterial narrowing pathologies induced by hyperplasia, i.e. an increase in cell number. The model ignores biochemical species (introduced in Section 2.4.2) in volumes calculation and assumes species as homogeneously distributed within layers.

Finally, by knowing the volume of a species X within layers, $V_{i,m}^X$, and the total volume of layers, $V_{i,m}$, we define species VFs as $\varphi_{i,m}^X = V_{i,m}^X/V_{i,m}$ in $\Omega_{i,m}$ and with the following properties $\sum_X \varphi_{i,m}^X = \sum_X \varphi_m^X = 1$. Table 1 summarizes numerical values of the fourteen arterial species model constant parameters.

## 2.3 Biochemical modelling

In Section 2.2, we proposed evolution equations for arterial species as summations of SFPs. These properties are shown to depend strongly on the GFs bioavailability which itself is modulated by hemodynamics and species dynamics (Chiu and Chien, 2011; Berk, 2001). Depending on the discipline of the authors, molecules signaling for SFPs have different names as “growth factors” in cellular biology or “cytokines” in immunology. In the present study, we consider the term GFs in an extensive understanding, including all types of biochemical messengers modulating species behaviours. For more than forty years, studies have focused on the influence of GFs on cells and have demonstrated that “... they provide an essential means for a cell to communicate with its immediate environment and to ensure that there is proper local homeostatic balance between the numerous cells that comprise a tissue ...” (Sporn and Roberts, 1991, Chap. 1).

We therefore develop hereafter a model of GFs dynamics within layers focusing on the most impacting families of GFs. Six GFs families and the nitric oxide
<table><caption>**Table 1** List of definitions, values, units and data source of constant parameters of the model. The parameters are grouped in several topics sections. By abuse of language, we write cells as a unit. We differentiate data if they are prescribed from or estimated from studies, chosen or estimated by us if unavailable in the literature, or either calculated</caption><thead><tr><th>Constant</th><th>Description</th><th>Value</th><th>Units</th><th>Reference(s)</th></tr></thead><tbody><tr><td colspan="5"><strong>Geometry</strong></td></tr><tr><td>w<sub>E</sub></td><td>EC width</td><td>13.2 × 10<sup>-6</sup></td><td>m</td><td>Prescribed from Garipcan et al. (2011)</td></tr><tr><td>l<sub>E</sub></td><td>EC length</td><td>25.8 × 10<sup>-6</sup></td><td>m</td><td>Prescribed from Garipcan et al. (2011)</td></tr><tr><td>e<sub>E</sub></td><td>Endothelium thickness</td><td>2.0 × 10<sup>-6</sup></td><td>m</td><td>Prescribed from Karner and Perktold (2000)</td></tr><tr><td>R<sub>10</sub></td><td>Initial luminal radius</td><td>4 × 10<sup>-3</sup></td><td>m</td><td>Chosen</td></tr><tr><td>e<sub>0</sub></td><td>Initial wall thickness</td><td>5 × 10<sup>-4</sup></td><td>m</td><td>Estimated from Karner and Perktold (2000)</td></tr><tr><td>R<sub>ext</sub></td><td>External wall radius</td><td>4.5 × 10<sup>-3</sup></td><td>m</td><td>Calculated as R<sub>ext</sub> = R<sub>10</sub> + e<sub>0</sub></td></tr><tr><td>Δz</td><td>Length of the arterial segment</td><td>1 × 10<sup>-3</sup></td><td>m</td><td>Chosen</td></tr><tr><td colspan="5"><strong>Hemodynamics</strong></td></tr><tr><td>μ</td><td>Dynamic viscosity of blood</td><td>3.45 × 10<sup>-3</sup></td><td>Pa s</td><td>Prescribed from Robertson et al. (2009)</td></tr><tr><td>ρ</td><td>Density of blood</td><td>1.056 × 10<sup>3</sup></td><td>kg m<sup>-3</sup></td><td>Prescribed from Robertson et al. (2009)</td></tr><tr><td>ν</td><td>Kinematic viscosity of blood</td><td>3.27 × 10<sup>-6</sup></td><td>m<sup>2</sup>/s</td><td>Calculated as ν = μ/ρ</td></tr><tr><td>Re<sub>0</sub></td><td>Initial Reynolds number</td><td>3.0 × 10<sup>2</sup></td><td>—</td><td>Chosen</td></tr><tr><td colspan="5"><strong>Arterial species</strong></td></tr><tr><td>r<sub>E</sub></td><td>ECs proliferation rate</td><td>7.2 × 10<sup>-1</sup></td><td>day<sup>-1</sup></td><td>Prescribed from Ultman et al. (2016)</td></tr><tr><td>r<sub>Q</sub></td><td>cSMCs proliferation rate</td><td>1.0 × 10<sup>-2</sup></td><td>day<sup>-1</sup></td><td>Estimated from Khosravi et al. (2020)</td></tr><tr><td>c<sup>0</sup></td><td>cSMCs to sSMCs dedifferentiation rate</td><td>2.5 × 10<sup>-1</sup></td><td>day<sup>-1</sup></td><td>Chosen</td></tr><tr><td>m<sup>0</sup></td><td>vSMCs migration rate</td><td>1.0 × 10<sup>-1</sup></td><td>day<sup>-1</sup></td><td>Chosen</td></tr><tr><td>p<sup>eq</sup></td><td>sSMCs proliferation rate</td><td>2.4 × 10<sup>-1</sup></td><td>day<sup>-1</sup></td><td>Prescribed from Boyle et al. (2010)</td></tr><tr><td>a<sup>eq</sup></td><td>sSMCs apoptosis rate</td><td>2.4 × 10<sup>-1</sup></td><td>day<sup>-1</sup></td><td>Calculated as p<sub>eq</sub> = a<sub>eq</sub></td></tr><tr><td>λ<sub>eq</sub></td><td>Collagen production rate</td><td>2.16 × 10<sup>-13</sup></td><td>g/cell/day</td><td>Prescribed from Donadoni et al. (2017)</td></tr><tr><td>χ<sup>eq</sup></td><td>Collagen self-degradation rate</td><td>3.3 × 10<sup>-2</sup></td><td>day<sup>-1</sup></td><td>Prescribed from Donadoni et al. (2017)</td></tr><tr><td>τ<sub>C</sub></td><td>Collagen ageing delay</td><td>4.5 × 10<sup>1</sup></td><td>day</td><td>Chosen</td></tr><tr><td>ρ<sub>c</sub></td><td>Collagen density</td><td>2 × 10<sup>3</sup></td><td>g/m<sup>3</sup></td><td>Prescribed from Donadoni et al. (2017)</td></tr><tr><td>ρ<sub>s</sub></td><td>Cell density</td><td>2.18 × 10<sup>14</sup></td><td>cell/m<sup>3</sup></td><td>Prescribed from Donadoni et al. (2017)</td></tr><tr><td>ρ<sub>E</sub></td><td>ECs density</td><td>2.80 × 10<sup>15</sup></td><td>cell/m<sup>3</sup></td><td>Calculated as ρ<sub>E</sub> = ((4/24)w<sub>E</sub>l<sub>E</sub>e<sub>E</sub>)<sup>-1</sup></td></tr><tr><td>t<sub>m</sub></td><td>Migration threshold value</td><td>1.8</td><td>—</td><td>Chosen</td></tr><tr><td>t<sub>c</sub></td><td>Dedifferentiation threshold value</td><td>1.1</td><td>—</td><td>Chosen</td></tr><tr><td colspan="5"><strong>Biochemistry</strong></td></tr><tr><td>ζ</td><td>Degradation rate of GFs</td><td>1 × 10<sup>-1</sup></td><td>day<sup>-1</sup></td><td>Estimated from Marino et al. (2017)</td></tr><tr><td>k<sub>c</sub></td><td>GFs coupling rate</td><td>1.0</td><td>day<sup>-1</sup></td><td>Chosen</td></tr><tr><td>β</td><td>Denudation parameter</td><td>3.0</td><td>—</td><td>Estimated from Conklin and Chen (2007)</td></tr><tr><td>R<sup>NO</sup><sub>max</sub></td><td>Maximal NO production rate</td><td>3.358 × 10<sup>8</sup></td><td>ng/(m<sup>3</sup>day)</td><td>Estimated from Andrews et al. (2010)</td></tr><tr><td>R<sup>NO</sup><sub>min</sub></td><td>Minimal NO production rate</td><td>-6.108 × 10<sup>-1</sup></td><td>ng/(m<sup>3</sup>day)</td><td>Estimated from Andrews et al. (2010)</td></tr><tr><td>α<sup>RNO</sup></td><td>NO tendency relation parameter</td><td>2.5</td><td>Pa<sup>-1</sup></td><td>Estimated from Andrews et al. (2010)</td></tr><tr><td>τ<sup>RNO</sup><sub>w</sub></td><td>NO offset parameter</td><td>9.49 × 10<sup>-1</sup></td><td>Pa</td><td>Estimated from Andrews et al. (2010)</td></tr><tr><td>P<sup>max</sup><sub>min</sub></td><td>Maximal apparent permeability</td><td>4.187 × 10<sup>-4</sup></td><td>m/day</td><td>Estimated from<sup>a</sup></td></tr><tr><td>P<sup>a</sup><sub>min</sub></td><td>Minimal apparent permeability</td><td>4.263 × 10<sup>-6</sup></td><td>m/day</td><td>Estimated from<sup>a</sup></td></tr><tr><td>α<sub>P</sub></td><td>Tendency relation parameter for P<sup>a</sup></td><td>-1.92</td><td>Pa<sup>-1</sup></td><td>Estimated from<sup>a</sup></td></tr><tr><td>τ<sup>P</sup><sub>w</sub></td><td>Offset parameter for P<sup>a</sup></td><td>6.845 × 10<sup>-1</sup></td><td>Pa</td><td>Estimated from<sup>a</sup></td></tr><tr><td>c<sup>NO</sup><sub>D</sub></td><td>Plasma NO concentration</td><td>1.88 × 10<sup>4</sup></td><td>ng/m<sup>3</sup></td><td>Calculated from<sup>b</sup></td></tr><tr><td>c<sup>PDGF</sup><sub>p</sub></td><td>Plasma PDGF concentration</td><td>4.5 × 10<sup>5</sup></td><td>ng/m<sup>3</sup></td><td>Prescribed from Rossi et al. (1998)</td></tr><tr><td>c<sup>FGF</sup><sub>p</sub></td><td>Plasma FGF concentration</td><td>6.4 × 10<sup>3</sup></td><td>ng/m<sup>3</sup></td><td>Prescribed from Larsson et al. (2002)</td></tr><tr><td>c<sup>Ag</sup><sub>p</sub></td><td>Plasma Ag concentration</td><td>7.2 × 10<sup>4</sup></td><td>ng/m<sup>3</sup></td><td>Prescribed from Levy et al. (1996)</td></tr><tr><td>c<sup>TGF</sup><sub>p</sub></td><td>Plasma TGF concentration</td><td>1.3 × 10<sup>6</sup></td><td>ng/m<sup>3</sup></td><td>Prescribed from Ivanovic et al. (2003)</td></tr><tr><td>c<sup>TNF</sup><sub>p</sub></td><td>Plasma TNF concentration</td><td>1.9 × 10<sup>4</sup></td><td>ng/m<sup>3</sup></td><td>Prescribed from Feitosa et al. (2013)</td></tr><tr><td>c<sup>MMP</sup><sub>p</sub></td><td>Plasma MMP concentration</td><td>3.0 × 10<sup>8</sup></td><td>ng/m<sup>3</sup></td><td>Estimated from Jonsson et al. (2016)</td></tr><tr><td>k</td><td>Arterial wall GFs consumption rate</td><td>1 × 10<sup>-2</sup></td><td>s<sup>-1</sup></td><td>Prescribed from Liu et al. (2014)</td></tr><tr><td>D</td><td>Arterial wall GFs diffusion coefficient</td><td>8.48 × 10<sup>-10</sup></td><td>m<sup>2</sup>s<sup>-1</sup></td><td>Prescribed from Liu et al. (2014)</td></tr><tr><td colspan="5"><strong>Temporal integration</strong></td></tr><tr><td>t<sub>0</sub></td><td>Initial time</td><td>0.0</td><td>day</td><td>Chosen</td></tr><tr><td>t<sub>f</sub></td><td>Final time</td><td>4.0 × 10<sup>4</sup></td><td>day</td><td>Chosen<sup>c</sup></td></tr><tr><td>d<sub>0</sub></td><td>Initial damage parameter</td><td>9.9 × 10<sup>-1</sup></td><td>—</td><td>Chosen</td></tr><tr><td>y<sub>0</sub></td><td>Initial conditions</td><td>see Table 5</td><td>—</td><td>Chosen (see Section 3.1)</td></tr><tr><td>h<sub>0</sub></td><td>History conditions</td><td>see Table 5</td><td>—</td><td>Chosen (see Section 3.1)</td></tr><tr><td>A<sub>tol</sub></td><td>Absolute numerical tolerance</td><td>1 × 10<sup>-8</sup></td><td>—</td><td>Chosen</td></tr><tr><td>R<sub>tol</sub></td><td>Relative numerical tolerance</td><td>1 × 10<sup>-6</sup></td><td>—</td><td>Chosen</td></tr><tr><td>ε</td><td>Remodeling sensitivity</td><td>1 × 10<sup>-2</sup></td><td>—</td><td>Chosen</td></tr><tr><td>ε</td><td>Equilibrium sensitivity</td><td>1 × 10<sup>-4</sup></td><td>—</td><td>Chosen</td></tr><tr><td colspan="5"><strong>Volume fractions of species within the arterial wall</strong></td></tr><tr><td>φ<sup>el</sup><sub>w</sub></td><td>Elastin</td><td>2.5 × 10<sup>-1</sup></td><td>—</td><td>Estimated from<sup>d</sup></td></tr><tr><td>φ<sup>vSMCs</sup><sub>w</sub></td><td>vSMCs</td><td>2.4 × 10<sup>-1</sup></td><td>—</td><td>Estimated from<sup>d</sup></td></tr><tr><td>φ<sup>C</sup><sub>w</sub></td><td>Collagen</td><td>4.5 × 10<sup>-1</sup></td><td>—</td><td>Estimated from<sup>d</sup></td></tr><tr><td>φ<sup>E</sup><sub>w</sub></td><td>ECs</td><td>3.7656 × 10<sup>-3</sup></td><td>—</td><td>Calculated</td></tr><tr><td>φ<sup>o</sup><sub>w</sub></td><td>Other species</td><td>5.6234 × 10<sup>-2</sup></td><td>—</td><td>Calculated</td></tr></tbody></table>

<footnote><sup>a</sup> Conklin and Chen (2007); Goodman et al. (2016)</footnote>

<footnote><sup>b</sup> We consider the characteristic concentration relation of $c_p^{NO} = 2R_{10} R_{max}^{NO} e_E / D^{NO}$ with the luminal NO diffusion coefficient $D^{NO} = 3.3 \times 10^{-9} \text{ m}^2 \text{s}^{-1}$ (Plata et al., 2010).</footnote>

<footnote><sup>c</sup> An another final time is reached corresponding to a new equilibrium state presented in Table 5.</footnote>

<footnote><sup>d</sup> O'Connell et al. (2008); Glagov et al. (1981); Bellini et al. (2013); Marino et al. (2017)</footnote>
molecule are considered and listed in Table 2. Regulatory activities of GFs go through several types of mechanisms with organization of producers and targets cells. Table 2 presents all the assumptions used for writing the model.

The amount of a GF named $x$ is denoted $\eta_{i,m}^x$ and its evolution follows the general kinetic differential equation of the form

$$ \frac{d\eta_{i,m}^x}{dt} = \bar{M}_{i,m}^x - [\varphi_{i,m}^S k_c N_c^x + \zeta^x] \eta_{i,m}^x + \sum_y \varphi_{i,m}^S k_c \eta_{i,m}^y \quad \text{in } \Omega_{i,m}, \quad (14) $$

with $\zeta^x$ the natural decay rate of $x$, $\bar{M}_{i,m}^x$ the source terms of $x$, and $\sum_y \varphi_{i,m}^S k_c \eta_{i,m}^y$ and $\varphi_{i,m}^S k_c N_c^x$ the related coupling terms between GFs $y$ and $x$ with $k_c$ the coupling rate. The form of equation (14) is based on previous molecular modelling (Donadoni et al., 2017; Marino et al., 2017). We model the inter-GFs coupling mechanisms, named autocrine signaling in Berk (2001), as follows. Coupling terms are scaled by the volume fractions of sSMCs (Marino et al., 2017) and for each coupling GFs action of $y$ on $x$, the same term appears in equation for GFs $y$ with negative sign to model inactivation after the signaling. This procedure results in the term $-\varphi_{i,m}^S k_c N_c^x \eta_{i,m}^x$ which is the proportion of $x$ inactivated by inter-GFs coupling mechanisms, where $N_c^x$ is the number of coupling mechanism considered for $x$. In this biochemical model, we limit ourselves to the minimal spatial scale of cells, we then not model intra-cellular kinetic phenomena of GFs binds, inactivate, and signal pathways acting on cellular fates as in Starbuck and Lauffenburger (1992). However, we consider a GF degradation term within layers $-\zeta^x \eta_{i,m}^x$.

The source terms $\bar{M}_{i,m}^x$ model the endothelium two-fold influence on a GF $x$: first, as a barrier, endothelium protects tissues from molecules and cells that flow in the lumen; second, endothelium is also a source of $x$, both as a producer and as a membrane that lets $x$ permeates from blood to arterial layers. The terms $\bar{M}_{i,m}^x$ are computed by averaging the term $m^x$ through the arterial wall over the thicknesses of layers $e_{i,m}$ as

$$ \bar{M}_{i}^{x} = \frac{1}{e_{i}} \int_{R_{i}}^{R_{IEL}} m^{x}(r) dr \quad \text{in } \Omega_{i}, \quad (15a) $$

$$ \bar{M}_{m}^{x} = \frac{1}{e_{m}} \int_{R_{IEL}}^{R_{EEL}} m^{x}(r) dr \quad \text{in } \Omega_{m}. \quad (15b) $$

The integrated flux over the surface $m^x$ is solution of a one-dimensional steady diffusion-reaction equation in

cylindrical coordinates and reads

$$ m^x(r) = 2\pi r \Delta z B^x \frac{I_1(\kappa r) K_1(\kappa R_{\text{ext}}) - K_1(\kappa r) I_1(\kappa R_{\text{ext}})}{I_1(\kappa R_l) K_1(\kappa R_{\text{ext}}) - K_1(\kappa R_l) I_1(\kappa R_{\text{ext}})} \quad \text{in } \Omega_w. \quad (16) $$

The boundary condition (BC) $B^x$ at $r = R_l$ is specified in Section 2.4, $\Delta z$ is the arterial length, $I_1$ and $K_1$ are respectively modified Bessel functions of first and second kinds and $\kappa^2 = k/D$ is the ratio of the tissue consumption rate $k$ over the molecular diffusion coefficient $D$ of $x$ in the tissue. Details about the analytical solution (16) are provided in appendix A.

We first note that the modelling of GFs releases by endothelium in equation (16) reminds of the radial exponential decrease proposed by Taber (1998) to account for the $\tau_w$ influence on “... biochemical signal released from the endothelium ...”. A second remark is that no term in equation (14) accounts for sSMCs production and degradation on $x$. To limit complexity on coupling relation in the model, we assume here that sSMCs production of $x$ is exactly compensated by sSMCs degradation of $x$. Thus, at the scale of the layer, we neglect this phenomenon. Furthermore, we justify this assumption by noticing that sSMCs are producer and target of all GFs with regard to Table 2.

In order to limit the number of parameters of the model, we suppose same rate $\zeta$ for all GFs and all coupling mechanisms of GFs have the same coupling rate $k_c$. The twenty biochemical model parameters are summarized in Table 1.

Overall, the biochemical model is the set of the seven pairs of equations (17). The present multiphysics and multiscale differential models of arterial wall dynamics is defined by equations (3), (4), (5), (6), (7) and (17). We write it in the following general vectorial form

$$ \frac{dy}{dt} = f(y, y(t - \tau_c), \tau_w), \quad (18) $$

with $y$ vector of variables, $\tau_w$ the magnitude of WSS and $f$ the right-hand-side terms. Equations (18) are solved on $t_0 \le t \le t_f$ with a given history $y(t) = h(t)$ for $t_0 - \tau_c \le t < t_0$ and an initial conditions $y(t_0) = y_0$.

## 2.4 Coupling hemodynamics with biochemical and species models

The growth of arterial tissue is a biomechanical problem coupled over time, space and the multiple factors that develop it. In this section, we propose coupling methods between biochemical environment, hemodynamics
<table><caption>**Table 2** Enumeration of growth factors considered in the biochemical model. Organizational assumptions of growth factors pro- ducer cells and target cells are listed with the corresponding references. For each growth factor, the inter-GFs coupling mechanisms of sSMCs (autocrine signaling) modelled are summarized with references and the total number of coupling mechanisms N<sub>c</sub> is given in the last column.</caption><thead><tr><th>GF family</th><th>Abbreviation</th><th>Producing cell(s)</th><th>Target(s)</th><th>Coupling relation(s)</th><th>N<sub>c</sub></th></tr></thead><tbody><tr><td>Nitric-Oxide</td><td>NO</td><td>ECs Andrews et al., sSMCs Bernhardt et al.</td><td>sSMCs Qiu et al.</td><td>-</td><td>0</td></tr><tr><td>Platelet-Derived-GF isoforms</td><td>PDGF</td><td>ECs San Martin et al., sSMCs San Martin et al.</td><td>sSMCs Qiu et al.</td><td>Ag Berk, TGF Berk, FGF Berk</td><td>3</td></tr><tr><td>Fibroblast-GF</td><td>FGF</td><td>ECs Model and Dardik, sSMCs Model and Dardik</td><td>sSMCs Berk</td><td>Ag Michel; Berk, PDGF Millette et al., MMP Parks et al.</td><td>3</td></tr><tr><td>Angiotensin-II</td><td>Ag</td><td>ECs Berk; Qiu et al., sSMCs Berk; Michel</td><td>sSMCs Berk; Michel</td><td>PDGF Berk, TGF Berk, FGF Berk</td><td>3</td></tr><tr><td>Transforming-GF</td><td>TGF</td><td>ECs Qiu et al., sSMCs Qiu et al.</td><td>sSMCs Qiu et al.</td><td>Ag Michel; Berk, MMP Parks et al.</td><td>2</td></tr><tr><td>Tumor-Necrosis-Factor</td><td>TNF</td><td>ECs Sigg, sSMCs Sigg</td><td>sSMCs Sigg</td><td>MMP Parks et al.</td><td>1</td></tr><tr><td>Matrix-MetalloProteinase-2-9</td><td>MMP</td><td>ECs Newby, sSMCs Qiu et al.; Beamish et al.</td><td>sSMCs Qiu et al., ECMs Newby</td><td>TNF Newby; Tedgui and Mallat, FGF Newby, PDGF Newby</td><td>3</td></tr></tbody></table>

$$ \frac{d\eta_{i,m}^{\text{NO}}}{dt} = \bar{M}_{i,m}^{\text{NO}} - \zeta\eta_{i,m}^{\text{NO}} \qquad \text{in } \Omega_{i,m} \qquad (17a) $$

$$ \frac{d\eta_{i,m}^{\text{PDGF}}}{dt} = \bar{M}_{i,m}^{\text{PDGF}} - (\zeta + N_c^{\text{PDGF}} \varphi_{i,m}^S k_c) \eta_{i,m}^{\text{PDGF}} + \varphi_{i,m}^S k_c (\eta_{i,m}^{\text{Ag}} + \eta_{i,m}^{\text{FGF}} + \eta_{i,m}^{\text{TGF}}) \quad \text{in } \Omega_{i,m} \qquad (17b) $$

$$ \frac{d\eta_{i,m}^{\text{FGF}}}{dt} = \bar{M}_{i,m}^{\text{FGF}} - (\zeta + N_c^{\text{FGF}} \varphi_{i,m}^S k_c) \eta_{i,m}^{\text{FGF}} + \varphi_{i,m}^S k_c (\eta_{i,m}^{\text{Ag}} + \eta_{i,m}^{\text{PDGF}} + \eta_{i,m}^{\text{MMP}}) \quad \text{in } \Omega_{i,m} \qquad (17c) $$

$$ \frac{d\eta_{i,m}^{\text{Ag}}}{dt} = \bar{M}_{i,m}^{\text{Ag}} - (\zeta + N_c^{\text{Ag}} \varphi_{i,m}^S k_c) \eta_{i,m}^{\text{Ag}} + \varphi_{i,m}^S k_c (\eta_{i,m}^{\text{PDGF}} + \eta_{i,m}^{\text{TGF}} + \eta_{i,m}^{\text{FGF}}) \quad \text{in } \Omega_{i,m} \qquad (17d) $$

$$ \frac{d\eta_{i,m}^{\text{TGF}}}{dt} = \bar{M}_{i,m}^{\text{TGF}} - (\zeta + N_c^{\text{TGF}} \varphi_{i,m}^S k_c) \eta_{i,m}^{\text{TGF}} + \varphi_{i,m}^S k_c (\eta_{i,m}^{\text{Ag}} + \eta_{i,m}^{\text{MMP}}) \quad \text{in } \Omega_{i,m} \qquad (17e) $$

$$ \frac{d\eta_{i,m}^{\text{TNF}}}{dt} = \bar{M}_{i,m}^{\text{TNF}} - (\zeta + N_c^{\text{TNF}} \varphi_{i,m}^S k_c) \eta_{i,m}^{\text{TNF}} + \varphi_{i,m}^S k_c \eta_{i,m}^{\text{MMP}} \quad \text{in } \Omega_{i,m} \qquad (17f) $$

$$ \frac{d\eta_{i,m}^{\text{MMP}}}{dt} = \bar{M}_{i,m}^{\text{MMP}} - (\zeta + N_c^{\text{MMP}} \varphi_{i,m}^S k_c) \eta_{i,m}^{\text{MMP}} + \varphi_{i,m}^S k_c (\eta_{i,m}^{\text{TNF}} + \eta_{i,m}^{\text{PDGF}} + \eta_{i,m}^{\text{FGF}}) \quad \text{in } \Omega_{i,m} \qquad (17g) $$

and volumic growth of arterial layers. Several metrics are available to quantify the local hemodynamical state, but since we concentrate on a steady flow regime in Section 2.1, we choose the magnitude of WSS, $\tau_w$, defined in equation (2b).

### 2.4.1 Hemodynamics stimuli modulate GFs releases

ECs have a variety of biological receptors that sense the flow and turn mechanical signals into biochemical signals transmitted to the underlying tissues (Chiu and Chien, 2011). In this section, we develop a coupling method between WSS and ECs functional property, i.e. endothelium permeability and GFs production as these properties are assumed to be important in regulating the local response to vascular injury (Li et al., 2005).

The main difficulty in the coupling process is feeding adequate numerical parameters in models. Humphrey (2008) insists that “... there remains a pressing need for more quantitative data that will enable the formulation of an integrative mathematical theory that describes and eventually predicts vascular adaptations in response to diverse stimuli.”

ECs produce several GFs at a rate which is cor- related to $\tau_w$. In several studies (Andrews et al., 2010; Humphrey, 2008), this dependency is analytically mod- elled as a sigmoid mathematical function. The relation

between GFs production rates and $\tau_w$ can be propor- tional (Andrews et al., 2010) or inversely proportional (Malek and Izumo, 1992). We assume here that pro- duction rates of all GFs is related to $\tau_w$ by the generic hyperbolic tangent function

$$ R^x(\tau_w) = \frac{R_{\text{max}}^x + R_{\text{min}}^x}{2} + \frac{R_{\text{max}}^x - R_{\text{min}}^x}{2} \tanh \left( \alpha R^x (\tau_w - \tau_w^x) \right), \qquad (19) $$

where $R_{\text{max,min}}^x$ are the extremum of production rates of GF $x$, $\alpha_R$ is a slope adjusting parameter of the func- tion, and $\tau_w^x$ is an offset, both being used to match experimental data.

Experimental data on GFs productions rates are scarce except for NO (Andrews et al., 2010). To ac- count for hemodynamics dependency of all GFs mod- elled in Section 2.3, we assume that parameters of $R^{\text{NO}}$ can be transposable to other GFs via a proportional or inversely proportional relationship. In equation (19), if $\alpha_{R^x} > 0$ the relation ($\tau_w, R^x$) is proportional, and if $\alpha_{R^x} < 0$ the relation becomes inversely proportional. This assumption permits to be only qualitatively con- sistent with experimental finding, as the production rates magnitude of each GF could be different. The de- pendence of $R^{\text{NO}}$ on $\tau_w$ is shown in Fig. 2, for both a negative or positive $\alpha_{R^{\text{NO}}}$ corresponding to inversely
proportional or proportional relationship between the couple ($\tau_w$, $R^{NO}$). Table 3 provide the considered relations between $R^{GFS}$ and $\tau_w$.

The apparent permeability $P^a$ of ECs monolayer is also modulated by $\tau_w$ and Buchanan et al. (2014); Conklin and Chen (2007) report a few discrete numerical values for $P^a$ and $\tau_w$. Previous study of Goodman et al. (2016) propose an analytical piecewise polynomial dependence of $P^a$ on $\tau_w$, while Silva et al. (2020) uses a hyperbolic functional dependence. For consistency with equation (19), we retain the form

$$P^a(\tau_w) = \frac{P_{\max}^a + P_{\min}^a}{2} + \frac{P_{\max}^a - P_{\min}^a}{2} \tanh(\alpha_{P^x}(\tau_w - \tau_w^P)). \quad (20)$$

Fig. 2 shows the variation of $P^a$ with $\tau_w$ as proposed by the hyperbolic function and experimental numerical values available.

From equations (19), and (20), we obtain the mass flux BC at the endothelium, in terms of $\tau_w$, as

$$B^x(\tau_w) = R^x(\tau_w)e_E + P^a(\tau_w)[c_p^x - c_w^x] \quad \text{in } \Gamma_e, \quad (21)$$

where $c_p^x$ is the concentration of $x$ in blood plasma, $c_w^x$ is the concentration of $x$ in arterial wall and $e_E$ is the endothelium thickness. This definition is consistent with previous BC treatments: for the production flux (Plata et al., 2010) and the flux due to permeability (Olgac et al., 2008). The general boundary condition for flux due to permeability considers a concentration difference which couples luminal and arterial wall domains. As a decoupling simplification, the endothelium is supposed highly resistant to $x$ transport and Olgac et al. (2008) propose that $c_p^x \gg c_w^x$. Finally, the mass flux BC reads

$$B^x(\tau_w) = R^x(\tau_w)e_E + P^a(\tau_w)c_p^x \quad \text{in } \Gamma_e. \quad (22)$$

The mass flux BC (22) is assumed to be modulated by the state of endothelium layer. For instance, the absence of ECs leads to no production of GFs, and the permeability parameter will be larger than in physiological condition. To quantify the monolayer state, a damage parameter is defined as $d = 1 - E^\dagger$ where $E^\dagger = E/E_{\max}$. In the limit cases, $d = 0$ if the monolayer is confluent and $d = 1$ if the monolayer is fully damaged. Assuming $B^x$ as linearly dependent on $d$, the modified damage-dependent BC (22) is therefore

$$B^x(\tau_w, d) = R^x(\tau_w)e_E(1-d) + (1+\beta d)P^a(\tau_w)c_p^x \quad \text{in } \Gamma_e, \quad (23)$$

where $\beta$ is a denudation parameter. We estimate $\beta$ equal to three because apparent permeability is multiplied by four-fold between an intact and denuded vessel from ECs at physiological shear (Conklin and Chen, 2007).

Those coupling methods between hemodynamics stimuli, ECs functional properties and monolayer state, model shear-dependent mechanotransduction of ECs and is presented in Fig. 3.

### 2.4.2 Biochemical environment modulates species functional properties

GFs are multifunctional on cells activities and, through the complex network of signaling pathways, they positively or negatively regulate species behaviour (Davies, 2019; Berk, 2001). SFPs modelled in Section 2.2 are influenced by their biochemical environment, whereas most existing mathematical models of vascular pathologies' development propose constant and independent parameters for SFPs as Donadoni et al. (2017); Wilstein et al. (2018). Based on a review of studies presented in Table 3 on the influence of GFs bioavailability on the pathologies' development, we propose in this section a coupling method between biochemical modelling and species dynamics modelling. This type of approach has already been proposed, in the context of dermal wound healing (Rognoni et al., 2018), in the modelling of inflammatory response (Nagaraja et al., 2014), or in the context of arterial growth (Irons and Humphrey, 2020; Khosravi et al., 2020). Remaining at the cellular scale, we will not be interested in modelling the underlying mechanisms of intra-cellular transduction pathways resulting from the binding of GFs to surface receptors of cells as in Irons and Humphrey (2020), we rather propose the following causal coupling.

As for GFs production rates (Section 2.4.1), experimental studies suggest that there is a proportional or inversely proportional relation between a SFP and the bioavailability of a GF within layers. As reviewed in Donadoni et al. (2017); Starbuck and Lauffenburger (1992), experimental findings have shown linear and hyperbolic functional dependence between GFs bioavailability and SFPs. As multiple GFs can influence a SFP, we assume overlap of these effects, so that we model them, for simplicity, as a linear relation between GFs amounts within layers and SFPs by computing weighted arithmetic averages. Table 3 gathers the relations for the considered GFs and SFPs.

We quantify the importance of a GF within layers by the following dimensionless parameter

$$\delta_{i,m}^x = \frac{\eta_{i,m}^x}{\tilde{\eta}_{i,m}^x} \quad \text{in } \Omega_{i,m}, \quad (24)$$

where $\eta_{i,m}^x$ is the amount of $x$ within the layer and $\tilde{\eta}_{i,m}^x$ is the physiological amount within layers, which is here assimilated to the initial equilibrium state. The dimensionless parameters $\delta_{i,m}^x$ permit to assess the relative
<figure><img src="image_3.png" /><figcaption>**Fig. 2** Left y-axis: Production rate of NO, RNO, as a function of WSS, τw, from equation (19). We obtain parameter values for RNO from nonlinear regression of sigmoid data of Andrews et al. (2010) (see Table 1). Right y-axis: Apparent permeability Pa from equation (20) where parameter values was obtained from nonlinear regression of data from literature (symbols)</figcaption></figure>

variations around physiological values, independently of their amplitude.

In the following, proportional relationships between GF x and a functional property will be obtained by summing δx_i,m and an inverse proportionality by summing 1/δx_i,m.

In addition to continuous functional dependencies, biological processes often also involve threshold mechanisms (Starbuck and Lauffenburger, 1992). Thus, according to Wilstein et al. (2018); Marino et al. (2017) we use the Macaulay function defined as

$$ \langle x - t_x \rangle = \frac{x - t_x + |x - t_x|}{2}, \qquad (25) $$

where t_x is the threshold value. By extension, multiple-thresholds mechanisms are dealt with the Φ function, used for collagen ageing rates, as

$$ \Phi(x, t_x, y, t_y) = \begin{cases} x - t_x, & x \ge t_x \text{ and } y \ge t_y \\ 0, & \text{otherwise} \end{cases} \qquad (26) $$

with t_x and t_y two threshold values.

With ratios δx_i,m, GFs bioavailability and SFPs relations, threshold mechanisms (25), and weighted arithmetic mean, we propose the relations (27) to model evolution of SFPs coupled with biochemical equations.

In equations (27), a distinction is made between maintenance functional properties and ones that triggers tissue growth, i.e. the migration m, and the differentiation c. The use of Macaulay function permits to

switch on m and c in certain conditions. This is done via threshold values t_{c,m}.

To simplify and because of a lack of data on relative impact of GFs on SFPs, every weighted parameters γ are set to unity in the present study. Moreover, the linear assumption between δ and SFPs assumes that cellular receptors are always available for GFs without saturation effect so that all amounts of GFs in layers act on SFPs at any time and linearly. Other models proposed different types of relations as the non-linear Michaelis–Menten function (Bulelzai and Dubbeldam, 2012; Khosravi et al., 2020).

As seen in Sections 2.4.1 and 2.3, the amounts of GFs within layers are dependent on τw and on the volume fractions of sSMCs. In this section, we connect bioavailability and species dynamics by taking into account the influences of GFs on SFPs thanks to the weighted arithmetic averages.

### 2.4.3 Loose coupling between hemodynamics and tissue growth

The above systems of equations for hemodynamics, evolution of species and biochemistry are linked a priori in a three-way instantaneous coupling: arterial remodeling, due to populations of cells evolution and collagen deposit, modifies the biochemical and the blood flow, which in turn imposes an evolution of WSS, hence biochemical and species exchanges, and so forth. In the-
$$p_{i,m} = p^{\text{eq}} \frac{\gamma_p^{\text{NO}} \delta_{i,m}^{\text{NO}-1} + \gamma_p^{\text{PDGF}} \delta_{i,m}^{\text{PDGF}} + \gamma_p^{\text{FGF}} \delta_{i,m}^{\text{FGF}} + \gamma_p^{\text{Ag}} \delta_{i,m}^{\text{Ag}} + \gamma_p^{\text{TGF}} \delta_{i,m}^{\text{TGF}} + \gamma_p^{\text{TNF}} \delta_{i,m}^{\text{TNF}}}{\sum_y \gamma_p^y} \quad (27a)$$

$$a_{i,m} = a^{\text{eq}} \frac{\gamma_a^{\text{NO}} \delta_{i,m}^{\text{NO}} + \gamma_a^{\text{FGF}} \delta_{i,m}^{\text{FGF}} - 1 + \gamma_a^{\text{Ag}} \delta_{i,m}^{\text{Ag}} + \gamma_a^{\text{TGF}} \delta_{i,m}^{\text{TGF}} - 1 + \gamma_a^{\text{TNF}} \delta_{i,m}^{\text{TNF}}}{\sum_y \gamma_a^y} \quad (27b)$$

$$m_i = m^0 \frac{\gamma_m^{\text{PDGF}} \langle \delta_i^{\text{PDGF}} - t_m \rangle + \gamma_m^{\text{FGF}} \langle \delta_i^{\text{FGF}} - t_m \rangle + \gamma_m^{\text{Ag}} \langle \delta_i^{\text{Ag}} - t_m \rangle + \gamma_m^{\text{TGF}} \langle \delta_i^{\text{TGF}} - t_m \rangle + \gamma_m^{\text{TNF}} \langle \delta_i^{\text{TNF}} - t_m \rangle}{\sum_y \gamma_m^y} \quad (27c)$$

$$\lambda_{i,m} = \lambda^{\text{eq}} \frac{\gamma_{\lambda}^{\text{NO}} \delta_{i,m}^{\text{NO}-1} + \gamma_{\lambda}^{\text{Ag}} \delta_{i,m}^{\text{Ag}} + \gamma_{\lambda}^{\text{TGF}} \delta_{i,m}^{\text{TGF}} + \gamma_{\lambda}^{\text{TNF}} \delta_{i,m}^{\text{TNF}}}{\sum_y \gamma_{\lambda}^y}, \quad (27d)$$

$$\chi_{i,m} = \chi^{\text{eq}} \delta_{i,m}^{\text{MMP}} \quad (27e)$$

$$c_{i,m} = c^0 \frac{\gamma_c^{\text{PDGF}} \langle \delta_{i,m}^{\text{PDGF}} - t_c \rangle + \gamma_c^{\text{FGF}} \langle \delta_{i,m}^{\text{FGF}} - t_c \rangle + \gamma_c^{\text{Ag}} \langle \delta_{i,m}^{\text{Ag}} - t_c \rangle + \gamma_c^{\text{MMP}} \langle \delta_{i,m}^{\text{MMP}} - t_c \rangle}{\sum_y \gamma_c^y} \quad (27f)$$

<table><caption>**Table 3** Influence of GFs on the following species functional properties: proliferation *p*, apoptosis *a*, migration *m*, dedifferentiation *c*, collagen deposition *λ*, and collagen self-degradation *χ*. The type of relationship between *τ*<sub>w</sub> and the ECs's production rate of GFs *R*<sup>*x*</sup> are presented in last column. Abbreviations “p.”, “i.p.” and “t.” indicate respectively proportional, inversely proportional and threshold relationship. References used are specified at each assumption</caption><thead><tr><th></th><th>p</th><th>a</th><th>m</th><th>c</th><th>λ</th><th>χ</th><th>(τ<sub>w</sub>, R<sup>2</sup>) relation</th></tr></thead><tbody><tr><td>NO</td><td>i.p. Qiu et al.</td><td>p. Qiu et al.</td><td>-</td><td>-</td><td>i.p. Wilstein et al.</td><td>-</td><td>p. Andrews et al.</td></tr><tr><td>PDGF</td><td>p. Qiu et al.; Berk; Sigg; San Martin et al.</td><td>p. Qiu et al.</td><td>p. Qiu et al.; Beamish et al.; San Martin et al.; Berk</td><td>p. Model and Dardik; Beamish et al.; Berk; Zhang et al.; Raines and Ross</td><td>-</td><td>-</td><td>p. Qiu et al.; Berk</td></tr><tr><td>FGF</td><td>p. Qiu et al.; Sigg; San Martin et al.</td><td>p. Berk; Rhoads et al.</td><td>p. Beamish et al.; Zhang et al.; Model and Dardik</td><td>p. San Martin et al.; Model and Dardik</td><td>-</td><td>-</td><td>p. Michel et al.</td></tr><tr><td>Ag</td><td>p. Beamish et al.; Sigg</td><td>p. Sigg</td><td>p. Model and Dardik</td><td>p. Beamish et al.</td><td>p. Sigg; Model and Dardik</td><td>-</td><td>p. Chin and Chen</td></tr><tr><td>TGF</td><td>p. San Martin et al.; Qiu et al.; Raines and Ross</td><td>i.p. Li et al.</td><td>p. Model and Dardik</td><td>p. Michel; Sigg; Raines and Ross</td><td>p. Michel; Sigg; Raines and Ross</td><td>-</td><td>i.p. Qiu et al.; Qiu et al.</td></tr><tr><td>TNF</td><td>p. Reetsenwald et al.; Sigg</td><td>p. Tedgui and Mallat</td><td>p. Sigg; Raines and Ross</td><td>p. Tedgui and Mallat</td><td>p. Newby; Beamish et al.</td><td>-</td><td>i.p. Reetsenwald et al.</td></tr><tr><td>MMP</td><td>-</td><td>-</td><td>-</td><td>p. Zhang et al.; Beamish et al.</td><td>-</td><td>-</td><td>i.p. Magid et al.</td></tr></tbody></table>

ory, this imposes a strong constraint on the way all the equations have to be solved. However, considering the assumptions that have been made at different stages of the derivation of the hemodynamical, biochemical, and cell populations models, it seems reasonable to consider a simplified way of coupling the three realms.

The assumption comes from a separation of timescales. The hemodynamic timescale is the cardiac cycle period, which is of the order of a second. The remodeling timescale ranges from days to months. Hence, we neglect the influence of arterial wall deformation due to growth on hemodynamics, as in Calvez et al. (2010). We therefore propose a loose coupling between arterial tissue dynamics and hemodynamics. Among the different types of vascular remodeling proposed by van Varik et al. (2012), we choose the hypothesis of an inward-hypertrophic remodeling. Thus, we assume that the external radius of the artery $R_{\text{ext}}$ is constant, and we propose two remodeling criteria from which blood flow can be modified by tissue dynamics. Remodeling occurs if either of the following conditions is reached

$$\frac{V_{\text{all}} - V_w}{V_l} > 1 + \varepsilon, \quad (28a)$$

$$\frac{V_{\text{all}} - V_w}{V_l} < 1 - \varepsilon. \quad (28b)$$

$\varepsilon$ is a remodeling parameter indicating hemodynamics sensitivity to tissue growth. $V_{\text{all}}$ is the whole constant arterial volume, i.e. the lumen volume and the

arterial wall volume, calculated as $V_{\text{all}} = \Delta z \pi R_{\text{ext}}^2$. $V_w$ is the time-dependent arterial wall volume defined from species dynamic evolution (12c), and $V_l$ is the morphological-generation-dependent luminal volume calculated as $V_l = \Delta z \pi R_l^2$ where $R_l$ is the luminal radius associated with the current hemodynamical state. The criterion (28a) corresponds to a tissue loss with an enlargement of the lumen, and (28b) corresponds to a tissue growth with luminal narrowing.

When these remodeling criteria are reached, e.g. during the development of a pathology, a new hemodynamical resolution is triggered, otherwise the blood flow remains the same. The latter case reflects the fact that cellular growth is not significant enough for the geometry of the artery to be changed. In the former case, a new arterial generation is produced in the form of a new geometrical configuration defined by its updated luminal volume $V_l$.

Another situation is the return to an equilibrium state during the resorption of the pathology, or when a new equilibrium is reached. To test this, we write a stop criterion and suppose that our system of equations (18) returns to equilibrium if

$$\left\| \frac{d\mathbf{y}^{\dagger}}{dt^{\dagger}} \right\| < \epsilon, \quad (29)$$

where $\epsilon$ is the equilibrium sensitivity parameter, $||\cdot||$ the $L_2$ norm and $d\mathbf{y}^{\dagger}/dt^{\dagger}$ the dimensionless time derivative
of y (seen Table 5 for the definitions of dimensionless variables). The two sensitivity parameters used in the loose coupling are given in Table 1.

The whole model and coupling methodology is presented in the flowchart of Fig. 3. The assumption of loose coupling between hemodynamics, biochemical, arterial species dynamics and arterial remodeling, simplifies greatly the modelling of a complex multifactorial and multiscale pathology both at equations' level and for numerical resolution cost, especially considering the fact that parametric studies may be desirable to test various evolution scenarios.

## 2.5 Numerical resolution

The coupled system of equations (18) is solved by a classical time marching algorithm, implemented as an object-oriented code written in Python. The resolution of the system of ODEs and DDEs is done with a Runge-Kutta solver based on Shampine and Thompson (2009). This solver uses an adaptive time step based on error estimation, continuous extensions, discontinuity tracking and events location routines to detect remodeling criteria (28a), (28b) and (29).

The multiscale, multifactorial system of equations involves large differences in the order of magnitude between variables, e.g. the count of cells and the amount of a GF within arterial wall layers. To avoid numerical issues, before numerical resolution, the species variables are respectively rescaled and non-dimensionalized by initial-physiological values (seen Table 5) and time is non-dimensionalized as $t^\dagger = tm_0$. The two tolerance parameters used in the time marching procedure are shown in Table 1.

## 3 Short and long time evolution of intimal hyperplasia

Our model and its numerical implementation is now tested on a first configuration. We consider the case of an idealized artery, shown in Fig. 1, which suffers a mechanically-induced local damage. The artery is assumed to be damaged infinitely on its length, so that the hemodynamics is a steady fully-developed Poiseuille flow all along the lesion evolution. The analytical solution for WSS reads

$$ \tau_w = \frac{4\mu Q_v}{\pi R_1^3} \quad \text{at } \Gamma_e, \qquad (30) $$

with $R_1$ the luminal radius and $Q_v$ the volumic blood flow rate. The hemodynamics is characterized by the following Reynolds number $Re = 2Q_v/(\pi\nu R_1)$. In this

study, we initialize hemodynamics by chosen $Re_0 = 300$. Hemodynamics is easily updatable from one arterial generation to another, by the use of $Q_v$ conservation between each generation. This simplification frees ourself from the use of numerical solution for the fluid dynamics part of the model. The update relations for the hemodynamical generation $k+1$ are respectively for the magnitude of WSS and the luminal volume

$$ \tau_w^{k+1} = \frac{4\mu Q_v}{\pi (R_1^{k+1})^3}, \qquad (31a) $$

$$ V_1^{k+1} = \pi \Delta z (R_1^{k+1})^2. \qquad (31b) $$

We now analyse the behaviour of our model under an initial damage which models denudation of endothelium as in arterial injury experimental models. We aim to validate the numerical results of our model both qualitatively and quantitatively in comparison to experiments.

### 3.1 Initialization with physiological-damage conditions

In order to initialize the test-case, we first determined from literature the physiological volume fractions of the species in the entire arterial wall. As these data were not found in a single source for a type of artery, they are the result of the data synthesis from O'Connell et al. (2008); Bellini et al. (2013); Marino et al. (2017); Fung (1993). The three main species' volume fractions, i.e. elastin, collagen and vSMCS, are set respectively to $\varphi_w^{el} = 0.25$, $\varphi_w^{vSMCs} = 0.24$ and $\varphi_w^C = 0.45$. Additionally to have $\sum_k \varphi_w^k = 1$, we define $\varphi_w^o = 0.06$, the volume fraction of others species present in the artery and not considered in those studies, e.g., fibroblast, ECs and so forth. As we model ECs dynamics, we estimate volume fraction of ECs from endothelium volume and from the assumption that all endothelium volume is composed of ECs. Thereby, the corresponding ECs density is function of endothelium thickness.

To define species physiological conditions in each layer starting from $\varphi_w$, we impose logical histologic observation, e.g. major part for elastin is in media or that main part of collagen is in adventitia. The volume fractions of species within layers are presented in Table 4. To validate our logical histologic observations, we compare the wall layers thicknesses resulting from Table 4 with those from Karner and Perktold (2000).

With Table 4, the entire arterial wall volume and the densities of species, i.e. collagen, vSMCs and ECs, we define the physiological state of species in each layer. Note that to differentiate sSMCs to cSMCs, we assume that balance between production and degradation of young collagen is under equilibrium, i.e. there is enough
<figure><img src="image_4.png" /><figcaption>**Fig. 3** Flowchart of the present multiscale tissue growth and remodeling framework. Coupling ways between models are represented by dashed curved arrows with the variable of interest specified. The decision tree is drawn in diamond with arrows. Criteria variables and criteria results are respectively before and after each diamond</figcaption></figure>

sSMCs to compensate self-degradation of young collagen fibers within each layer. For physiological quantities of GFs, $\tilde{\eta}_{i,m}^x$, knowing the physiological $\phi_{i,m}^{Sph}$, WSS $\tau_w^{ph}$ and assuming an equilibrium state for equations (17), a resolution of this defined system of equations gives them.

Physiologic state for variable of the system (18) are shown in Table 5 as the vector $\mathbf{y}^{\text{ref}}$. The history condition of the arterial wall is assumed physiological as $\mathbf{h}(t) = \mathbf{h}_0 = \mathbf{y}^{\text{ref}}$ for $t_0 - \tau_C \le t < t_0$. The initial condition models a damage to the only endothelial layer and so $\mathbf{h}_0 = \mathbf{y}^{\text{ref}}$ except for the EC variable which is set to $E(t_0) = 1 - d_0$ with $d_0$ a damage parameter. Note that, under the previously defined initial physiological conditions, with history as $\mathbf{y}(t) = \mathbf{y}^{\text{ref}}$ for $t_0 - \tau_C \le t \le t_0$, our model predicts a stable equilibrium state.

## 3.2 Results

We experiment numerically the tissue growth occurring after an initial denudation of the ECs layer. In our test-case, the damage parameter is fixed at initial time to $d_0 = 0.99$ in order to compare with data from literature. As reviewed by Lemson et al. (2000), most research on IH has been conducted with injured animal models by a balloon catheter (Clowes et al., 1983a,b) or wire filament (Lindner et al., 1993). The experiments damaged principally endothelium, but also the underlying layers. By this numerical experiment, we want to simulate the same type of injury. Furthermore, we assume that no damage occurs on intimal and medial species, except on ECs, in order to investigate the model response to a unique perturbation of the endothelium.

In the following, we will consider the response of the artery to this damage, first at short time of the order of several days, then over a longer evolution period over several months.
<table><caption>**Table 4** Proposed volume fractions, $\phi_w$, of species within layers. The volume fraction of species X within intima is equal to $\phi_w^{X,i} = V_i^X/V_w$ with $V_i^X$ the volume of X in intima and $V_w$ the whole arterial wall volume. The “Whole wall” column is a sum on the layers of the volume fractions for each species and prescribed value from Table 1 are recovered. The volume fraction of each layer is given in row $\phi_w^{layer}$. The thicknesses of layers associated to the defined volume fractions are given and compared with a previous study (Karner and Perktold, 2000)</caption><thead><tr><th></th><th colspan="4">Volume fractions in (-)</th></tr><tr><th></th><th>Intima</th><th>Media</th><th>Adventicia</th><th>Whole wall</th></tr></thead><tbody><tr><td>φ<sub>w</sub><sup>el</sup></td><td>2.25 × 10<sup>-3</sup></td><td>2.2275 × 10<sup>-1</sup></td><td>2.5 × 10<sup>-2</sup></td><td>2.5 × 10<sup>-1</sup></td></tr><tr><td>φ<sub>w</sub><sup>vSMCs</sup></td><td>1.08 × 10<sup>-2</sup></td><td>2.052 × 10<sup>-1</sup></td><td>2.4 × 10<sup>-2</sup></td><td>2.4 × 10<sup>-1</sup></td></tr><tr><td>φ<sub>w</sub><sup>C</sup></td><td>2.25 × 10<sup>-3</sup></td><td>4.275 × 10<sup>-2</sup></td><td>4.05 × 10<sup>-1</sup></td><td>4.5 × 10<sup>-1</sup></td></tr><tr><td>φ<sub>w</sub><sup>E</sup></td><td>3.7656 × 10<sup>-3</sup></td><td>0.0</td><td>0.0</td><td>3.7656 × 10<sup>-3</sup></td></tr><tr><td>φ<sub>w</sub><sup>o</sup></td><td>2.8117 × 10<sup>-3</sup></td><td>2.8117 × 10<sup>-3</sup></td><td>5.0611 × 10<sup>-2</sup></td><td>5.6234 × 10<sup>-2</sup></td></tr><tr><td>φ<sub>w</sub><sup>layer</sup></td><td>2.1877 × 10<sup>-2</sup></td><td>4.7351 × 10<sup>-1</sup></td><td>5.0461 × 10<sup>-1</sup></td><td>1.0</td></tr></tbody></table><table><thead><tr><th></th><th colspan="4">Thickness in (μm)</th></tr><tr><th></th><th>Intima</th><th>Media</th><th>Adventicia</th><th>Whole wall</th></tr></thead><tbody><tr><td>Present study</td><td>1.1606 × 10<sup>1</sup></td><td>2.4344 × 10<sup>2</sup></td><td>2.4496 × 10<sup>2</sup></td><td>5.0 × 10<sup>2</sup></td></tr><tr><td>Karner and Perktold (2000)</td><td>1.2 × 10<sup>1</sup></td><td>3.02 × 10<sup>2</sup></td><td>1.86 × 10<sup>2</sup></td><td>5.0 × 10<sup>2</sup></td></tr></tbody></table>

Despite the diversity of experiments (e.g. animal models, type of artery, experimental conditions and objectives of studies), we tried to collect as many results as possible in any form — as moments of occurrence of events or temporal variation of variables — that may be useful for validating the model. Several ways of comparing these data with model results are possible. We choose hereafter to present the short- and long-time temporal variations of relevant variables, mainly found in the experiments of the literature.

### 3.2.1 Short timescale evolution

First, we examine ECs dynamics. Immediately after the injury and for approximately two weeks, the simulated biological phenomena are set off by ECs loss and regeneration. Fig. 4 shows the time evolution of $E^\dagger$ starting from initial value $E^\dagger(t_0) = 1 - d_0$. $E^\dagger$ is the percentage of ECs affected by the initial damage. The figure shows an initial slow regeneration of $E^\dagger$, with a rapid acceleration after one week, before full re-endothelization, corresponding to $E^\dagger \simeq 1$ after 14 days. Overall, this corresponds to a logistic dynamics of endothelial cells. The timescale of restoration of endothelial monolayer varies greatly according to the experimental conditions and the species considered. In animal models the regeneration has been seen as early as 3 days after procedure and ended after 4 weeks (de Vries et al., 2016), and in 3 weeks after injury (Lindner et al., 1993) or 6–10 weeks in several others animal models (Cornelissen and Vogt, 2019). The breakpoint time of our simulated ECs re-

<figure><img src="image_5.png" /><figcaption>**Fig. 4** Re-endothelization of the ECs layer in the spatial domain $\Gamma_e$ during the first fourteen days. The dimensionless variable is $E^\dagger = E/E_{max}$</figcaption></figure>

generation is consistent with timescale from literature.

Second, our model predicts time evolution of functional properties of vSMCs. Fig. 5 (a) presents the evolution of the dedifferentiation rates $c_i$ and $c_m$ in layers (defined in (27f)), from 0 to 25 days. In intima and media, dedifferentiation starts at day 0.16 (~ 4 hours) while it ends in intima after 24.4 days and in media after 22.3 days. Fig. 5 (b) presents the migration rate $m_i$ (equation (27c)) of vSMCs. Authors reported the beginning of migration phenomenon around 2–4 days
after the injury (Raines, 2004; Clowes, 1993). Migration continues of up to one week (Raines, 2004) or one month (Lemson et al., 2000). We have calibrated our model so that migration starts around the fourth day (at $t = 3.5$ days) after the injury, and as this functional property is modelled by a Macaulay function (equation (27c)) the adjusted threshold parameter is fixed to $t_m = 1.8$. The model prediction for the end of migration is after 17.4 days, as shown by Fig. 5 (b). The time intervals for migration and dedifferentiation are defined by a prone biochemical state within layers. In this state, GFs are overexpressed so that they become greater than threshold values, and thus pilot SFPs. The time evolution of GFs is discussed and shown in appendix B on Fig. 12. Fig. 5 (a) and (b) exhibit a very rapid triggering of large rates $c_i$ and $c_m$ from their activations. This indicates that both phenomena of dedifferentiation from cSMCs to sSMCs and of migration of vSMCs from the media toward intima are predominant in the early phase of the neointimal formation. Fig. 5 (c) focuses on the proliferation $p_i$ and apoptosis $a_i$ of sSMCs in the intima, from their initial equilibrium state where $a_i(t_0) = p_i(t_0)$. It also shows the time evolution of the turnover rate $r_i = p_i - a_i$. During the first 21 days, before the cross-over of the curves for $p_i$ and $a_i$, $r_i > 0$ indicating an exponential growth of sSMCs. This imbalance between cellular proliferation and apoptosis leads to IH as seen in Fig. 6. After $t \simeq 21$ days, $r_i$ becomes negative so that sSMCs decrease. Again, this dynamics is related to the evolution of GFs. The latter explains the non-uniform evolution of apoptosis $a_i$ evolution, whereas proliferation $p_i$ exhibits clear growth and decay phases with a maximum reached approximately at day 7.8.

Finally, the other important phenomenon within short timescale is sSMCs dynamics. sSMCs proliferation in early days is a key factor for the arterial narrowing and for the future deposition of ECM within neointimal lesion. Fig. 6 shows the counts $S_{i,m}$ of sSMCs in intima and media divided by their respective physiological counts $S_{i,m}^{ph}$, in addition to their growth rates $r_{S_{i,m}}$. The figure shows an exponential increase of the count of sSMCs within intima and media until day 21, where the maximal count of sSMCs is predicted. In the literature, the vSMCs' maximal count within the intima was observed at day 14 (Davies, 2019). This sSMCs increase within layers is responsible for the artery's rapid narrowing until day 21 (see Fig. 8). We compare now the times $t_{i,m}^p$ at which the proliferation peaks occur in our model. They are found by computing the growth

rates of $S_{i,m}$ equation as

$$r_{S_{i,m}} = \frac{1}{S_{i,m}} \frac{dS_{i,m}}{dt} \qquad (32)$$

and finding its maximum, as proposed by Schwartz et al. (1996).

As seen in Fig. 6, our model predicts maximal proliferation in intima and media respectively at days 5.7 and 4.4. For comparison, the experimental peak of neointima maximal proliferation for rats is at day 2.04, for pigs day 6.3 and for humans in reconstructed data at day 14.7 (Schwartz et al., 1996). Clowes et al. (1983a) indicate, in an injured rat carotid experiment, that sSMCs proliferation is maximal at day 4 in intima and day 2 in media (data reported on Fig. 6). Clearly, there are important cell kinetic differences across species (Schwartz et al., 1996). Nonetheless, the predictions of our model about peak proliferation instants within intima and media are consistent with this literature.

Peaks of $r_S$ are closely related to maximal IH. As reviewed by Kenagy (2011), maximal IH is reached two months after arterial injury. Davies (2019) mentions that this event occurs after 1 month, mainly because of cellular proliferation and ECM production. Variability of these data between experimental results again can be explained by species and experimental differences. In our simulation, Fig. 6 and 8 indicate that the maximal IH occurs at day 21. From this point of view, day 21 appears to be a kind of turning point in the evolution of the pathology, such that the dynamic evolution of species amounts to a luminal enlargement phase, as a transition between early and mature lesion. The longer time evolution is discussed in the coming section.

### 3.2.2 Long timescale evolution

The ECM changes in IH is a key element at long timescale (Kenagy, 2011). Over time, the injured zone undergoes a maturation phase and becomes paucicelluar and more fibrous. After a period that ranges from a few months to about one year, the lesion reaches a steady state wherein intimal lesion composition is in the range 60%-80% of ECM and 20-40% of vSMCs (Clowes, 1993; Kenagy, 2011), depending on the degree of injury (Schwartz et al., 1996; Ducasse et al., 2003). Our model's predictions for the evolution of intimal collagen $\varphi_i^{C_y+C_o}$ and vSMCs $\varphi_i^{Q+S}$ VFs are shown on Fig. 7. The model predicts a composition of intima layer with 73% of ECM and 20% of vSMCs after 12 months, values that are almost already reached after two months of evolution. The lesion dynamic behaviour at long time is thus qualitatively and quantitatively consistent with data from Clowes (1993); Kenagy (2011). The simulated IH ends
<figure><img src="image_6.png" /><figcaption>**Fig. 5** Time evolution of the vSMCs functional properties between 0 and 25 days. (a) De-differentiation of cSMCs in the intimal $c_i$ (line with symbols) and medial $c_m$ (solid line) layers. (b) Intimal migration $m_i$ of vSMCs. (c) Intimal proliferation $p_i$ (dashed line), apoptosis $a_i$ (solid line with symbol), and turnover rate $r_i = p_i - a_i$ (dashed line with symbol) of sSMCs</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Fig. 6** Left y-axis: Intimal (dashed line) and medial (solid line) sSMCs evolution within the first month as $S^\dagger = S/S^{ph}$. Right y-axis: Intimal (dashed line with symbols) and medial (solid line with symbols) growth rate of $S$, $r_S$, with peak proliferation instants of intima and media layers marked on x-axis respectively at $t_i^p = 5.7$ and $t_m^p = 4.4$ days. Inner intima growth rate $r_{S_i}^C$ (diamond symbols) and first layer of media $r_{S_m}^C$ (triangle symbols) are reported from Clowes et al. (1983a)</figcaption></figure>

by reaching the stop criteria (29) after 37 months. At 24 months, the intimal lesion is stable (data not shown), in the sense of the criterion (29), although in the media this event appears after 37 months because of cSMCs regeneration (data not shown). Finally, by observing the evolution of the luminal radius as a function of time in Fig. 8, we identify three phases in the development of the simulated lesion. From initial time to day 21, arterial wall is in an early exponential growth phase. Then from day 21 to day 150, the luminal radius decreases mainly because of sSMCs apoptosis predominance. Finally, from day 150 until day 1107, the maturation phase leads to a new equilibrium state of the arterial wall. This equilibrium state is characterized by values of the species variables that are provided in Table 5.
<figure><img src="image_8.png" /><figcaption>**Fig. 7** Intimal volume fractions of collagen $\varphi_i^{C_y+C_o}$ (solid lines), and vSMCs $\varphi_i^{Q+S}$ (solid lines with symbols) within the first year of the experiment</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 8** Dimensionless luminal radius predicted in the numerical experiment. The inset figure shows a shorter period until month 1.5. Red crosses indicate instants when criteria (28) are reached and hemodynamics is updated (31). The stop criteria (29) instant is shown with a dot</figcaption></figure>

## 3.3 Discussion

### 3.3.1 Lesion evolution

We first discuss the dynamics of sSMCs predicted by our model, showing non-trivial behaviour at different timescales of evolution of the pathology. We then discuss the transition phase between early cellular proliferation and lesion fibrotic steady state.

Early modelled hyperplasia is marked by an exponential growth until day 21 as seen in Fig. 6. This phase is driven by GFs bioavailability that promotes growth

and proliferation, a feature that comes from the coupling assumptions between GFs and SFPs and from ECs denudation. It follows from this configuration that $p_{i,m} - a_{i,m} > 0$ (Fig. 5). The cellular proliferation peak appears in intima at $t_i^p = 5.7$ days, and is mainly driven by the migration and dedifferentiation phenomena. The $t_i^p$ instant is determined by the constraint imposed in our model, that the migration process starts around 2–4 days, provided by experimental data. The strong increase of the intimal growth rates $r_S$ around day 4 is due to the use of Macaulay functions in the modelling of the migration rate of vSMCs. Discussing the values of growth rate shown in Fig. 6, in intima the model is qualitatively and quantitatively consistent with the reported thymidine index from Clowes et al. (1983a) despite that, between day 21 and day 30, the model predicts negative growth rate whereas experiments report $r_S$ tending towards zero. For the growth rate in the media layer, the model is also consistent with data in $t \in [7, 14]$ days and underestimates growth rate at the onset of the lesion in $t \in [0, 4]$ days. This underestimation may be due to the initial medial difference between simulation and experiment as Clowes et al. (1983a) report damages in media layer after the arterial injury. This unmodelled experimental damage may impact medial sSMCs dynamics at short timescale. As seen in Fig. 8, this cellular growth is responsible for the arterial narrowing at short time before day 21. Interestingly, the hemodynamical updates of the WSS during this narrowing phase slow down the exponential growth by modifying the biochemistry, i.e. production rate of GFs and transport of GFs through the endothelium (data not shown).

In the second transition phase, from day 21 until day 150, our model predicts a regression of tissue growth and an enlargement of the arterial lumen. This phase is due to the exponential decrease of sSMCs population within layers. This phase occurs in intima and media respectively, at day 21 and at day 22. It is due to the change of sign of the functional parameter $p_{i,m} - a_{i,m}$ driven by biochemical dynamics. The predicted bioavailability of GFs in this regression phase no longer permits species growth, according to GFs influence assumptions. Neointimal regression has been discussed in Farb et al. (2004). On rat carotid artery in-stent neointimal, the lesion size decreases between days 28 and 60. On canine coronary, a gradual thinning of neointima 12 months after implant is observed, along with reduced cellularity and increased fibrosis. Similar observations were made on pig and human. The numerical results are thus consistent with experiments about this regression phase during which the lesion becomes increasingly fibrous (Farb et al., 2004). The amount of collagen fiber,
young or aged, produced during this second phase can- not compensate the cellular losses. However, after day 150, a re-narrowing phase is seen because of the inti- mal and medial cSMCs regeneration (data not shown). During the second phase which appears to be a tran- sition phase, the neointima becomes less cellular and increasingly fibrous as seen in Fig. 7.

### 3.3.2 Critical discussion of the proposed model and parameters

Parameter estimation in phenomenological modelling is essential for realistic model prediction, since the solu- tions are highly dependent on the choice of parameters (Donadoni et al., 2017). In biomechanical modelling, from one study to another, parameter values may have significant differences in the reported numerical values. The use of different sources — animal models, species, in vivo or in vitro experiments — can lead to a great variability in the values of parameters and in the ensu- ing phenomenological model results. It is the case for the decay rate of GFs that varies from 0.1 day⁻¹ (Es- cuer et al., 2019) to 0.86 day⁻¹ (Marino et al., 2017), or the SMCs collagen secretion rate λeq = 2.16 × 10⁻¹³ g/cells/day from (Donadoni et al., 2017) and λeq = 2–5 × 10⁻¹¹ g/cells/day from (Khosravi et al., 2020), or for the regeneration rate of endothelium discussed below.

In addition, the choice of a type of equation to model a specie's dynamics must be associated correctly with the values of the parameters chosen. To model cSMCs dynamics, we use a logistic equation as Donadoni et al. (2017) but with a different rQ value. Donadoni et al. (2017) propose rQ = 5 × 10⁻⁴ day⁻¹ with a character- istic timescale τQ = 1/rQ around 5.8 years. For this parameter value applied to equation (4b) and assum- ing initial condition Q†m(t0) = 0.5, the ODE predicts a return to the physiological state after more than 30 years. This order of magnitude for cSMCs regeneration appears to be inconsistent with arterial wound heal- ing. We estimated rQ = 1.0 × 10⁻² day⁻¹ with the corresponding characteristic timescale τQ = 0.27 years, which is consistent with the values retained for cSMCs proliferation by Khosravi et al. (2020). In the same dy- namic test with our rQ value, we obtain a regeneration after 2 years, which is in better agreement with physio- logical repair timescales for arterial wall. This remark is also applicable to the ECs regeneration rate rE. Escuer et al. (2019) propose rE = 8.6 × 10⁻² day⁻¹ while Ult- man et al. (2016) give rE = 0.72 day⁻¹. Using param- eter values proposed in Escuer et al. (2019) or in Ult- man et al. (2016), the dynamic test (with E†(t0) = d0) yields regeneration after respectively 150 days and 15

days. We thus conclude that the choice of the numer- ical values of the parameters is decisive in the present phenomenological modelling. As this choice can be dif- ficult, sensitivity analysis can help to choose between values proposed in the literature.

We investigate the sensitivity of the model to a per- turbation of its initial equilibrium state by different lev- els of endothelium damages. Fig. 9 shows the time evo- lution over months of the luminal radius Rl normalized by Rl0, and its sensitivity on the intensity of the initial damage on endothelium d0, and, as such, has a di- rect influence on the perturbation of static homeostasis equilibrium and of tissue growth modulation. The fig- ure shows that the luminal radius does not decrease, or only of a few percent, after a weak damage is applied, that is at d0 < 0.5, but peaks at less than 80% of the nominal radius at larger d0 = 0.99. For d0 = 0.01, the arterial wall does not trigger a tissue growth and return to an equilibrium state after 53 days. The case where d0 = 0.255 does not develop tissue growth either, and a new equilibrium state is found after 64 months reach- ing Rl/Rl0 = 9.98 × 10⁻¹. The absence of tissue growth after light damages is due to the use of threshold mech- anisms for dedifferentiation and migration of vSMCs. In the two cases, the perturbation of the bioavailabil- ity is not sufficient to trigger lesion onset. It can be observed that the more d0 increases, the more the ex- ponential proliferation phase would be important but over a shorter period. This result is consistent with ex- perimental observations about the influence of the ini- tial injury on the extent of tissue growth (Davies, 2019; Ducasse et al., 2003).

Given the key role in lesion set off, we investigate the effect of the variation of rE on the minimal dimen- sionless luminal radius

$$ \text{MdLR} = \min_{t \in [t_0, t_f]} \frac{R_l}{R_{l0}}, \qquad (33) $$

experienced by the artery during the pathology that starts at time t0 and ends at time tf. As shown in Fig. 10, even a small variation of rE greatly modifies the pro- gression of the pathology and for rE < 0.46 the present model predicts an arterial occlusion. This is not sur- prising since phenomenological modelling is sensitive to parameter values and since the duration of endothelial denudation has a key role in the development of intimal lesions (Fingerle et al., 1990).

### 3.4 Possible paths for model improvement

As seen above, our phenomenological model allows to study a complex pathology in a simplified way, at the price of a series of assumptions. The result is a model
<figure><img src="image_10.png" /><figcaption>**Fig. 9** Time evolution of normalized luminal radius $R_1/R_{10}$ colored with the value of the initial damage $d_0 \in [0.01, 0.99]$. The inset figure shows a shorter period until month 4. The dashed arrow with dots in the insert emphasizes the trend that, the more $d_0$ increase, the more the lumen narrows over a shorter period</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 10** Sensitivity of minimal dimensionless luminal radius MdLR, defined in equation (33), to regeneration rate of ECs layer $r_E \in [0.4, 1.0] \text{ day}^{-1}$</figcaption></figure>

validated on a 1D test-case and which may permits to explore parametric cases. In this section, we propose various paths for improving the model and critical discussions on the hypotheses formulated.

Firstly, the assumption of homogeneously distributed species within the arterial wall allows us to strongly reduce the computational time by reducing the dimensionality of the system. Nevertheless, arterial patholo-

gies are heterogeneous (Davies, 2019; Feugier and Rouvière, 2018; Clowes et al., 1983b). A coarse modelling cannot reproduce spatially localized processes which can be strongly coupled with WSS stimuli and mass transport near the luminal surface.

In terms of species, our model focuses on the three main ones responsible for the pathology, i.e. ECs, vSMCs and collagen. Nonetheless, as seen in Section 3.2, the model can still predict some of the phenomena observed experimentally, although the actual pathology is of course more varied in terms of species and in signaling pathways. Enrichment of the model could be obtained via exploration of others relations between GFs bioavailability, SFPs and ECM; e.g. it is known that ECM has feedback mechanism modulating cells SFPs (Rognoni et al., 2018). Another way for improvement could be the addition of species linked to the process of inflammation (monocytes, macrophages), wound repair (Rognoni et al., 2018) (fibroblasts, myofibroblasts) or platelet activation immediately after injury (Nagaraja et al., 2014). As reviewed by Davies (2019); Ducasse et al. (2003); Forte et al. (2010); Kenagy (2011); Lemson et al. (2000), inflammatory cells from blood flow, fibroblast/myofibroblast from adventitia and platelet activation play a determining role in the onset and/or regulation of the pathology via release of diverse GFs. In addition, we have simplified ECM to collagen fibers while the matrix is composed of multiple fibers that may have an important role in vascular pathology as proteoglycan described in endofibrosis pathology (Feugier and Rouvière, 2018).

The precise dependence of production rate dynamics of each GF on WSS is unknown. There is nothing to suggest that they will have the same responses as nitric oxide at the quantity level produced by ECs but also at the sensitivity to the WSS stimuli. Indeed, data from endothelin-1 production rate as a function of WSS (Humphrey, 2008) shows that WSS sensitivity is different from that of NO. By lack of further information, we had to assume that the production rate values of all GFs are matched to that of $R^{\text{NO}}$ (see Section 2.4.1), which makes the biochemical model only qualitative. The biochemical model was designed to be as basic as possible with a reduced number of constant parameters. We thus supposed that GFs share same properties as coupling and decay rates values while each GFs family have specific ones (Khosravi et al., 2020; Escuer et al., 2019). Another hypothesis was to assume that all secreted GFs were in an active states, while as for TGF-$\beta$ a distinction can be achieved between latent and active state (Nagaraja et al., 2014) or as for angiotensin II which is the end product of the highly complex renin-angiotensin system (Michel, 2004). This kind of im-
provements of the biochemical part could thus greatly enrich the model.

Data used for comparison comes from animal models. As reviewed by Bennett and O'Sullivan (2001), animal model can be used, *but with caution over their applicability*". Depending on the species, the cellular processes can be both qualitatively and quantitatively different, as the associated characteristic timescale can also be largely variable. Again, we use data from these references since they are the only one available to our knowledge. From this viewpoint, our model is validated for animal hyperplasia more than human hyperplasia, although this is a common methodology whereby animal models provide insight into the mechanisms of human vascular diseases.

Finally, our modelling does not take into account the states of intra-parietal stresses during the pathological evolution. It would be interesting to incorporate this aspect of the mechanical stresses in the model of intimal hyperplasia because during lesion development, these stimuli could play an important role in the physiopathology. Indeed, as the structural properties of the wall evolve, *e.g.* volume fraction of collagen or vSMCs, these changes must modulate the mechanical characteristics of the wall, *e.g.* elasticity, rigidity, and its macroscopic mechanical stress/strain relation. The multiple experiments describing circumferential and axial mechanical stresses as modulators of the SFPs Berk (2001); Humphrey (2008) and of macroscopic arterial wall homeostatic state Humphrey (2008) suggest to explore this way of enriching the model.

## 4 Conclusion

In the present study, we propose a novel multiscale and multiphysics model of arterial growth. We rely on a modelling approach similar to previous ones (Donadoni et al., 2017; Wilstein et al., 2018; Khosravi et al., 2020), and we thus use a deterministic model for a vascular pathology through phenomenological equations. The present model is successfully applied to a one-dimensional artery test-case by numerically simulating neointima hyperplasia development after endothelial damage. Starting from experimental hypotheses at both microscopic and macroscopic scales, respectively related to cellular dynamics, biochemistry and hemodynamics, our model permits to obtain results consistent with the literature. Our model predicts a complex phasic tissue growth dynamics, strongly coupled with hemodynamics and biochemistry. This result confirms that the development of intimal hyperplasia is a multifactorial and multiscale (with two main timescales observed) pathology of which many phenomena remain to be understood (Davies, 2019).

Our model is simple enough so that it can be used for testing numerically the multiple findings on the physiopathology that come from several types of studies, *e.g.* animal, *in vivo*, *in vitro* experiments. It is therefore possible to collect and test the experimental results coming from the different scales and physics of the pathology while also being able to explore different coupling hypotheses. In particular, this study permits to investigate the key role of GFs in development of intimal hyperplasia.

Two perspectives are considered below as improvements of our present work.

Firstly, our model of IH was applied to one-dimensional test-case but further implementation with a higher dimensionality is called for. We started in the present work the indispensable validation of equations in a simple geometrical configuration, to focus on tissue growth and avoid the burden of CFD computations. However, the pathology development coupled with hemodynamics involves pathology-protective and pathology-promoting zones linked to complex and three-dimensional flow patterns. As a halfway step, in a coming study, we will investigate the behaviour of our model in a two-dimensional axisymmetric idealized artery suffering a spatially localized desendothelization. Keeping the compartmental approach proposed by Donadoni et al. (2017), we will consider a mesh of the arterial wall composed of a number of compartments, in which we will apply the system of equations (18). Using our hypotheses of coupling between hemodynamics and tissue growth, this two-dimensional extension will evaluate the responses of our multiscale model in the high and low wall shear stress zones.

Secondly, improvements of the arterial wall modelling require taking into account elastic deformation into its growth and remodeling. In this context, kinematic growth theory (Rodriguez et al., 1994) extended to the constrained mixture theory (Humphrey, 2002) then more recently to the homogenized constrained mixture theory (Cyron et al., 2016) allow to model growth and remodeling as the contribution of an elastic deformation and an inelastic growth in the arterial material. By applying these concepts to intimal hyperplasia, we aim to propose a more realistic modelling of this pathology in further works.

**Acknowledgements** We acknowledge fruitful discussions with Dr Lara Chalabreysse and Dr Pauline Baillard about the physiopathology of intimal hyperplasia and endofibrosis.
# A Analytical development about diffusion–reaction equation

## A.1 Dimensional analysis

We use cylindrical coordinates $(r, \theta, z)$, with $r$ the radius, $\theta$ the polar angle, and $z$ the axial coordinate. The unsteady diffusion–reaction equation that describes the evolution of the concentration $C$ of a growth factor (GF) $x$ within an idealized arterial wall, which assumes the only radial dependence of $C$, is

$$ \frac{\partial C}{\partial t} = \frac{D}{r} \frac{\partial}{\partial r} \left( r \frac{\partial C}{\partial r} \right) - kC \quad \text{in } \Omega_w. \qquad (34) $$

$D$ is the $x$ diffusion coefficient within the arterial wall, and $k$ is the consumption rate of the medium. The idealized arterial wall and domain $\Omega_w$ are described in Fig. 1.

The dimensionless unsteady diffusion–reaction equation reads

$$ \frac{L^2}{T_{\text{gr}}D} \frac{\partial C^{\dagger}}{\partial t^{\dagger}} = \frac{1}{r^{\dagger}} \frac{\partial}{\partial r^{\dagger}} \left( r^{\dagger} \frac{\partial C^{\dagger}}{\partial r^{\dagger}} \right) - k^{\dagger} C^{\dagger} \quad \text{in } \Omega_w, \qquad (35) $$

with the dimensionless variables

$$ C^{\dagger} = \frac{C}{\tilde{C}}, \quad r^{\dagger} = \frac{r}{L}, \quad t^{\dagger} = \frac{t}{T_{\text{gr}}}. \qquad (36) $$

$\tilde{C}$, $L$ and $T_{\text{gr}}$ are a characteristic concentration, a characteristic diffusion length and a characteristic growth timescale. The dimensionless consumption is $k^{\dagger} = kL^2/D$, and the dimensionless ratio $L^2/(T_{\text{gr}}D)$ expresses diffusion timescale versus growth timescale.

Upon assuming short growth timescale $T_{\text{gr}} \sim 1$ day, thickness of arterial wall $L = 500$ µm and NO diffusion coefficient $D = 8.48 \times 10^{-10}$ m² s⁻¹, the dimensionless number $L^2/(T_{\text{gr}}D) = 3 \times 10^{-3}$ remains much lower than one. According to this, we can remove the time-dependent term in equation (34) — this was also done in Goodman et al. (2016)— and we get the following steady diffusion reaction

$$ \frac{D}{r} \frac{d}{dr} \left( r \frac{dC}{dr} \right) - kC = 0 \quad \text{in } \Omega_w. \qquad (37) $$

Note that even accounting for a transitional timescale of endothelial cells (ECs) adaptation to hemodynamics, which would be much smaller than $T_{\text{gr}}$ as $T_{\text{tr}} \sim 1$ hour (Hahn and Schwartz, 2009), the time derivative term in (34) remains negligible ($L^2/(T_{\text{tr}}D) = 8 \times 10^{-2}$).

## A.2 Development of source term expression

The complete problem to be solved is equation (37) plus the following Neumann boundary conditions

$$ -D \left. \frac{dC}{dr} \right|_{r=R_1} = B \quad \text{at} \quad r=R_1, \qquad (38a) $$

$$ \left. \frac{dC}{dr} \right|_{r=R_{\text{ext}}} = 0 \quad \text{at} \quad r=R_{\text{ext}}, \qquad (38b) $$

to model ECs influence on the $x$ GF bioavailability within arterial layers. Equations (38a) and (38b) are respectively a flux conservation condition (see equation (23)) and a no flux condition at external arterial radius.

<figure><img src="image_12.png" /><figcaption>Fig. 11 $m^x$ from equation (41) for $r \in [R_1, R_{\text{EEL}}]$. The averaged values of $m^x$ (from equation (15)) within intima and media layers, $\tilde{M}_{i,m}^x$, are represented with dashed lines</figcaption></figure>

Equation (37) is a Bessel equation, whose solution with boundary conditions (38) is

$$ C(r) = \frac{B}{\kappa D} \frac{I_1(\kappa R_{\text{ext}})K_0(\kappa r) + K_1(\kappa R_{\text{ext}})I_0(\kappa r)}{K_1(\kappa R_1)I_1(\kappa R_{\text{ext}}) - I_1(\kappa R_1)K_1(\kappa R_{\text{ext}})} \qquad (39) $$

with $\kappa^2 = k/D$, $I_{0,1}$, $K_{0,1}$ modified Bessel functions of first and second kind. The vector mass flux comes from Fick's law $J = -D\nabla C$, and its amplitude in the radial direction $n$ comes from (39) as

$$ \begin{aligned} J &= (-D\nabla C) \cdot n \\ &= B \frac{I_1(\kappa r)K_1(\kappa R_{\text{ext}}) - K_1(\kappa r)I_1(\kappa R_{\text{ext}})}{I_1(\kappa R_1)K_1(\kappa R_{\text{ext}}) - K_1(\kappa R_1)I_1(\kappa R_{\text{ext}})}. \end{aligned} \qquad (40) $$

To define the source term of the generic equation (14), we integrate $J$ over the surface of length $\Delta z$ at radius $r$, which amounts to multiplying by $2\pi r \Delta z$ as

$$ m(r) = 2\pi r \Delta z B \frac{I_1(\kappa r)K_1(\kappa R_{\text{ext}}) - K_1(\kappa r)I_1(\kappa R_{\text{ext}})}{I_1(\kappa R_1)K_1(\kappa R_{\text{ext}}) - K_1(\kappa R_1)I_1(\kappa R_{\text{ext}})}. \qquad (41) $$

$m$ is the production rate of GF as a function of the radial coordinate. To get average production rates of GF $x$ over the layers considered, namely the intima and the media, equation (41) is integrated over the layers thicknesses, as in (15). This averaging procedure and the radial variation of $m^x$ are shown in Fig. 11.

# B Dynamics of growth factors

For reference, we show in Fig. 12 the evolution of all the GFs in both media and intima layers for the test-case presented in Section 3, between 0 and 30 days. The initial damage of ECs layer causes a subexpression of NO amount (Lemson et al., 2000) and overexpression of other GFs at short timescale (Ducasse et al., 2003). There are two main types of time evolution for GFs, that of NO and that of other GFs (PDGF, FGF, Ag, TGF, TNF, MMP). From an analysis of biochemical equations, the
NO is not involved in inter-GFs coupling mechanisms so its dy-
namics is mainly driven by hemodynamics stimuli. This is seen
between 6 and 30 days, where there is a net increase in NO
amount within intima layer corresponding to the luminal nar-
rowing phase (as WSS increases, so does the production rate).
In the second group, a different type of temporal evolution is
seen because this group is strongly coupled (seen Table 2). They
evolve globally in the same way with an increase in the short
time induced by the denudation of the endothelium and a de-
crease below the equilibrium value $\delta^{\text{GFs}} = 1.0$ around day 21
which corresponds to the turning point ($r_i < 0$) discussed in
Section 3.2.

C Final equilibrium state of the test-case

The final equilibrium state of the denudation test-case devel-
oped in Section 3 is presented in Table 5. We provide the vector
of variables **y** of (18) as **y**† = **y**/*y*<sup>ref</sup> where **y**† is the vector of
rescaled or dimensionless variables by its initial–physiological
values **y**<sup>ref</sup>.

References

Andrews A, Jaron D, Buerk D, Kirby P, Barbee K (2010) Direct, real-time measurement of shear stress-induced nitric oxide produced from endothelial cells in vitro. Nitric Oxide: Biology and Chemistry 23(4):335-42, DOI 10.1016/j.niox.2010.08.003

Beamish J, He P, Kottke-Marchant K, Marchant R (2010) Molecular regulation of contractile smooth muscle cell phenotype: Implications for vascular tissue engineering. Tissue engineering Part B, Reviews 16:467-91, DOI 10.1089/ten.TEB.2009.0630

Bellini C, Ferruzzi J, Roccabianca S, Di Martino E, Humphrey J (2013) A microstructurally motivated model of arterial wall mechanics with mechanobiological implications. Annals of biomedical engineering 42, DOI 10.1007/s10439-013-0928-x

Bennett M, O'Sullivan M (2001) Mechanisms of angioplasty and stent restenosis: implications for design of rational therapy. Pharmacology & Therapeutics 91(2):149-166, DOI 10.1016/S0163-7258(01)00153-X

Berk B (2001) Vascular smooth muscle growth: Autocrine growth mechanisms. Physiological Reviews 81(3):999-1030, DOI 10.1152/physrev.2001.81.3.999, pMID: 11427690, https://doi.org/10.1152/physrev.2001.81.3.999

Bernhardt J, Tschudi M, Dohi Y, Gut I, Urwyler Ba (1991) Release of nitric oxide from human vascular smooth muscle cells. Biochemical and Biophysical Research Communications 180(2):907-912, DOI 10.1016/S0006-291X(05)81151-9

Binder B, Simpson M (2016) Cell density and cell size dynamics during in vitro tissue growth experiments: Implications for mathematical models of collective cell behaviour. Applied Mathematical Modelling 40(4):3438-3446, DOI 10.1016/j.apm.2015.01.016

Boyle C, Lennon A, Early M, Kelly D, Lally C, Prendergast P (2010) Computational simulation methodologies for mechanobiological modelling: A cell-centred approach to neointima development in stents. Philosophical transactions Series A, Mathematical, physical, and engineering sciences 368:2919-2935, DOI 10.1098/rsta.2010.0071

Buchanan C, Verbridge S, Vlachos P, Rylander M (2014) Flow shear stress regulates endothelial barrier function and

expression of angiogenic factors in a 3D microfluidic tu-
mor vascular model. Cell Adhesion & Migration 8, DOI
10.4161/19336918.2014.970001

Bulelzai M, Dubbeldam J (2012) Long time evolution of atherosclerotic plaques. Journal of Theoretical Biology 297:1-10, DOI 10.1016/j.jtbi.2011.11.023

Calvez V, Houot J, Meunier N, Raoult A, Rusnakova G (2010) Mathematical and numerical modeling of early atherosclerotic lesions. ESAIM: Proc 30:1-14, DOI 10.1051/proc/2010002

Chiu J, Chien S (2011) Effects of disturbed flow on vascular endothelium: Pathophysiological basis and clinical perspectives. Physiological Reviews 91:327-387, DOI 10.1152/physrev.00047.2009

Clowes A (1993) Chapter 17 intimal hyperplasia and graft failure. Cardiovascular Pathology 2(3, Supplement):179-186, DOI 10.1016/1054-8807(93)90058-A

Clowes A, Reidy M, Clowes M (1983a) Kinetics of cellular proliferation after arterial injury. I. smooth muscle growth in the absence of endothelium. Laboratory Investigation 49(3):327-333

Clowes A, Reidy M, Clowes M (1983b) Mechanisms of stenosis after arterial injury. Laboratory Investigation 49(3):208-215.

Conklin B, Chen C (2007) Effect of low shear stress on per-
meability and occludin expression in porcine artery en-
dothelial cells. World journal of surgery 31:733-743, DOI
10.1007/s00268-006-0735-8

Cornelissen A, Vogt F (2019) The effects of stenting on coronary endothelium from a molecular biological view: Time for improvement? Journal of Cellular and Molecular Medicine 23(1):39-46, DOI 10.1111/jcmm.13936, https://onlinelibrary.wiley.com/doi/pdf/10.1111/jcmm.13936

Cyron C, Aydin R, Humphrey J (2016) A homogenized constrained mixture (and mechanical analog) model for growth and remodeling of soft tissue. Biomechanics and Modeling in Mechanobiology 15(6):1389-1403, DOI 10.1007/s10237-016-0770-9

Davies M (2019) Rutherford's Vascular Surgery and Endovascular Therapy, Elsevier, chap 5 : Intimal Hyperplasia

Donadoni F, Pichardo-Almarza C, Bartlett M, Dardik A, Homer-Vanniasinkam S, Díaz-Zuccarini V (2017) Patient-specific, multi-scale modeling of neointimal hyperplasia in vein grafts. Frontiers in Physiology 8:226, DOI 10.3389/fphys.2017.00226

Ducasse E, Eschwege J, Chevalier J, de Ravignan D, Puppinck P, Lartigau E (2003) Hyperplasie intimale artérielle par prolifération de cellules musculaires lisses dans la paroi : données actuelles, traitements expérimentaux et perspectives. Journal des maladies vasculaires 28:130-144, DOI JMV-06-2003-28-3-0398-0499-101019-ART6

Escuer J, Martínez M, McGinty S, Peña E (2019) Mathematical modelling of the restenosis process after stent implantation. Journal of The Royal Society Interface 16(157):20190313, DOI 10.1098/rsif.2019.0313, https://royalsocietypublishing.org/doi/pdf/10.1098/rsif.2019.0313

Farb A, Kolodgie F, Hwang J, Burke A, Tefera K, Weber D, Wight T, Virmani R (2004) Extracellular matrix changes in stented human coronary arteries. Circulation 110(8):940-947, DOI 10.1161/01.cir.0000139337.56084.30

Feitosa MCP, de Sousa Lima VB, Neto JMM, do Nascimento Marreiro D (2013) Plasma concentration of IL-6 and TNF-alpha and its relationship with zincemia in obese women. Revista da Associação Médica Brasileira 59:429-434, DOI 10.1016/j.ramb.2013.03.003
<figure><img src="image_13.png" /><figcaption>**Fig. 12** Time evolution of dimensionless GFs amounts $\delta^x = \eta^x/\tilde{\eta}$ in (-) within intima (solid lines) and media (solid lines with symbols) layers between 0 and 30 days</figcaption></figure>

Feugier P, Rouvière O (2018) Pathologie vasculaire du sportif de haut niveau : endofibrose artérielle. EMC-Angéiologie 1750(19):1-10, DOI 10.1016/S2211-0364(18)85302-X

Fingerle J, Tina Au Y, Clowes A, Reidy M (1990) Intimal lesion formation in rat carotid arteries after endothelial denudation in absence of medial injury. Arteriosclerosis: An Official Journal of the American Heart Association, Inc 10(6):1082-1087, DOI 10.1161/01.ATV.10.6.1082, https://www.ahajournals.org/doi/pdf/10.1161/01.ATV.10.6.1082

Forte A, Della Corte A, De Feo M, Cerasuolo F, Cipollaro M (2010) Role of myofibroblasts in vascular remodeling: focus on restenosis and aneurysm. Cardiovascular Research 88(3):395-405, DOI 10.1093/cvr/cvq224, https://academic.oup.com/cardiovascres/article-pdf/88/3/395/1129421/cvq224.pdf

Fung Y (1993) Mechanical Properties and Active Remodeling of Blood Vessels, Springer New York, pp 321-391. DOI 10.1007/978-1-4757-2257-4_8

Garipcan B, Maenz S, Pham T, Settmacher U, Jandt K, Zanow J, Bossert J (2011) Image analysis of endothelial microstructure and endothelial cell dimensions of human arteries — a preliminary study. Advanced Engineering Materials 13(1-2):B54-B57, DOI 10.1002/adem.201080076, https://onlinelibrary.wiley.com/doi/pdf/10.1002/adem.201080076

Glagov S, Grande J, Vesselinovitch D, Zarins C (1981) Quantitation of Cells and Fibers in Histologic Sections of Arterial Walls: Advantages of Contour Tracing on a Digitizing Plate, Springer New York, New York, pp 57-93. DOI 10.1007/978-1-4612-5967-1_2

Goodman M, Luo X, Hill N (2016) Mathematical model on the feedback between wall shear stress and intimal hyperplasia. International Journal of Applied Mechanics 8(7)

Hahn C, Schwartz M (2009) Mechanotransduction in vascular physiology and atherogenesis. Nature Reviews Molecular
<table><caption>**Table 5** New state of equilibrium reached after $t_f \approx 1107$, days ($\approx 37$ months) presented in terms of rescaled and dimensionless variables $y^\dagger = y/y^{\text{ref}}$, with $y^{\text{ref}}$ the initial-physiological value</caption><thead><tr><th>y<sup>†</sup> formulation</th><th>y<sup>†</sup></th><th>y<sup>ref</sup></th></tr><tr><td></td><td colspan="2" style="text-align: center;">Value</td></tr></thead><tbody><tr><td>δ<sub>i</sub><sup>NO</sup> = η<sub>i</sub><sup>NO</sup>/η̃<sub>i</sub><sup>NO</sup></td><td>1.05</td><td>1.25 × 10<sup>-2</sup> ng</td></tr><tr><td>δ<sub>m</sub><sup>NO</sup> = η<sub>m</sub><sup>NO</sup>/η̃<sub>m</sub><sup>NO</sup></td><td>8.57 × 10<sup>-1</sup></td><td>8.02 × 10<sup>-3</sup> ng</td></tr><tr><td>δ<sub>i</sub><sup>PDGF</sup> = η<sub>i</sub><sup>PDGF</sup>/η̃<sub>i</sub><sup>PDGF</sup></td><td>8.38 × 10<sup>-1</sup></td><td>2.66 ng</td></tr><tr><td>δ<sub>m</sub><sup>PDGF</sup> = η<sub>m</sub><sup>PDGF</sup>/η̃<sub>m</sub><sup>PDGF</sup></td><td>6.64 × 10<sup>-1</sup></td><td>1.62 ng</td></tr><tr><td>δ<sub>i</sub><sup>FGF</sup> = η<sub>i</sub><sup>FGF</sup>/η̃<sub>i</sub><sup>FGF</sup></td><td>8.42 × 10<sup>-1</sup></td><td>3.35 ng</td></tr><tr><td>δ<sub>m</sub><sup>FGF</sup> = η<sub>m</sub><sup>FGF</sup>/η̃<sub>m</sub><sup>FGF</sup></td><td>6.81 × 10<sup>-1</sup></td><td>2.11 ng</td></tr><tr><td>δ<sub>i</sub><sup>Ag</sup> = η<sub>i</sub><sup>Ag</sup>/η̃<sub>i</sub><sup>Ag</sup></td><td>8.38 × 10<sup>-1</sup></td><td>2.65 ng</td></tr><tr><td>δ<sub>m</sub><sup>Ag</sup> = η<sub>m</sub><sup>Ag</sup>/η̃<sub>m</sub><sup>Ag</sup></td><td>6.64 × 10<sup>-1</sup></td><td>1.62 ng</td></tr><tr><td>δ<sub>i</sub><sup>TGF</sup> = η<sub>i</sub><sup>TGF</sup>/η̃<sub>i</sub><sup>TGF</sup></td><td>8.43 × 10<sup>-1</sup></td><td>3.67 ng</td></tr><tr><td>δ<sub>m</sub><sup>TGF</sup> = η<sub>m</sub><sup>TGF</sup>/η̃<sub>m</sub><sup>TGF</sup></td><td>6.85 × 10<sup>-1</sup></td><td>2.33 ng</td></tr><tr><td>δ<sub>i</sub><sup>TNF</sup> = η<sub>i</sub><sup>TNF</sup>/η̃<sub>i</sub><sup>TNF</sup></td><td>8.43 × 10<sup>-1</sup></td><td>4.24 ng</td></tr><tr><td>δ<sub>m</sub><sup>TNF</sup> = η<sub>m</sub><sup>TNF</sup>/η̃<sub>m</sub><sup>TNF</sup></td><td>6.86 × 10<sup>-1</sup></td><td>2.70 ng</td></tr><tr><td>δ<sub>i</sub><sup>MMP</sup> = η<sub>i</sub><sup>MMP</sup>/η̃<sub>i</sub><sup>MMP</sup></td><td>8.51 × 10<sup>-1</sup></td><td>7.07 ng</td></tr><tr><td>δ<sub>m</sub><sup>MMP</sup> = η<sub>m</sub><sup>MMP</sup>/η̃<sub>m</sub><sup>MMP</sup></td><td>7.19 × 10<sup>-1</sup></td><td>4.75 ng</td></tr><tr><td>Q<sub>i</sub><sup>†</sup> = Q<sub>i</sub>/Q<sub>i max</sub></td><td>1.00</td><td>2.23 × 10<sup>4</sup> cells</td></tr><tr><td>Q<sub>m</sub><sup>†</sup> = Q<sub>m</sub>/Q<sub>m max</sub></td><td>9.99 × 10<sup>-1</sup></td><td>4.23 × 10<sup>5</sup> cells</td></tr><tr><td>S<sub>i</sub><sup>†</sup> = S<sub>i</sub>/S<sub>i ph</sub></td><td>5.49</td><td>9.18 × 10<sup>3</sup> cells</td></tr><tr><td>S<sub>m</sub><sup>†</sup> = S<sub>m</sub>/S<sub>m ph</sub></td><td>1.17</td><td>1.74 × 10<sup>5</sup> cells</td></tr><tr><td>C<sub>y i</sub><sup>†</sup> = C<sub>y i</sub>/C<sub>y i ph</sub></td><td>5.61</td><td>6.01 × 10<sup>-8</sup> g</td></tr><tr><td>C<sub>y m</sub><sup>†</sup> = C<sub>y m</sub>/C<sub>y m ph</sub></td><td>1.31</td><td>1.14 × 10<sup>-6</sup> g</td></tr><tr><td>C<sub>o i</sub><sup>†</sup> = C<sub>o i</sub>/C<sub>y i ph</sub></td><td>3.39 × 10<sup>1</sup></td><td>6.01 × 10<sup>-8</sup> g</td></tr><tr><td>C<sub>o m</sub><sup>†</sup> = C<sub>o m</sub>/C<sub>y m ph</sub></td><td>2.85</td><td>1.14 × 10<sup>-6</sup> g</td></tr><tr><td>E<sup>†</sup> = E/E<sub>max</sub></td><td>1.00</td><td>1.41 × 10<sup>5</sup> cells</td></tr></tbody></table>

Cell Biology 10:53–62, DOI doi.org/10.1038/nrm2596

Humphrey J (2008) Vascular adaptation and mechanical homeostasis at tissue, cellular, and sub-cellular levels. Cell Biochemistry and Biophysics 50:53–78, DOI 10.1007/s12013-007-9002-3

Humphrey JD (2002) Cardiovascular Solid Mechanics. Springer New York, DOI 10.1007/978-0-387-21576-1

Irons L, Humphrey J (2020) Cell signaling model for arterial mechanobiology. PLOS Computational Biology 16(8):1–22, DOI 10.1371/journal.pcbi.1008161

Ivanovic V, Todorovic-Rakovic N, Demajo M, Z NK (2003) Elevated plasma levels of transforming growth factor-b1 (TGF-b1) in patients with advanced breast cancer: association with disease progression. European Journal of Cancer 39(4):454–461, DOI 10.1016/S0959-8049(02)00502-6

Jonsson A, Hjalmarsson C, Falk P (2016) Levels of matrix metalloproteinases differ in plasma and serum-aspects regarding analysis of biological markers in cancer. British Journal of Cancer 67:703–706, DOI doi.org/10.1038/bjc.2016.127

Karner G, Perktold K (2000) Effect of endothelial injury and increased blood pressure on albumin accumulation in the arterial wall: A numerical study. Journal of Biomechanics 33:709–715, DOI 10.1016/S0021-9290(99)00226-2

Kenagy R (2011) Mechanisms of Vascular Disease, University of Adelaide Press, chap 7: Biology of Restenosis and Targets for Intervention

Khosravi R, Ramachandra A, Szafrań J, Schiavazzi D, Breuer C, Humphrey J (2020) A computational biochemo-mechanical model of in vivo tissue-engineered vascular graft development. Integrative Biology 12(3):47–63, DOI 10.1093/intbio/zyaa004, https://academic.oup.com/ib/article-pdf/12/3/47/33043008/zyaa004.pdf

Larsson A, Sköldenberg E, Ericson H (2002) Serum and plasma levels of FGF-2 and VEGF in healthy blood donors. Angiogenesis 5:107–110

Lemson M, Tordoir J, Daemen M, Kitslaar P (2000) Intimal hyperplasia in vascular grafts. European Journal of Vascular and Endovascular Surgery 19(4):336–350, DOI 10.1053/ejvs.1999.1040

Levy B, Benessiano J, Henrion D, Caputo L, Heymes C, Duriez M, Poitevin P, Samuel J (1996) Chronic blockade of AT2-subtype receptors prevents the effect of angiotensin II on the rat vascular structure. The Journal of Clinical Investigation 98(2):418–425, DOI 10.1172/JCI118807

Li L, Zhang X, Li X, Lv C, Yu H, Xu M, Zhang M, Fu Y, Meng H, Zhou J (2016) TGF-beta inhibits the apoptosis of pulmonary arterial smooth muscle cells and contributes to pulmonary vascular medial thickening via the PI3K/Akt pathway. Molecular Medicine Reports 13(3):2751–2756, DOI 10.3892/mmr.2016.4874

Li Y, Haga J, Chien S (2005) Molecular basis of the effects of shear stress on vascular endothelial cells. Journal of Biomechanics 38(10):1949–1971, DOI 10.1016/j.jbiomech.2004.09.030

Lindner V, Fingerle J, Reidy M (1993) Mouse model of arterial injury. Circulation Research 73(5):792–796, DOI 10.1161/01.res.73.5.792
Liu X, Wang Z, Zhao P, Fan Z, Sun A, Zhan F, Fan Y, Deng X (2014) Nitric Oxide Transport in Normal Human Thoracic Aorta: Effects of Hemodynamics and Nitric Oxide Scavengers. PLoS ONE 9(11):e112395

Magid R, Murphy T, Galis Z (2003) Expression of matrix metalloproteinase-9 in endothelial cells is differentially regulated by shear stress role of c-Myc. The Journal of Biological Chemistry 278:32994-32999, DOI 10.1074/jbc.M304799200

Malek A, Izumo S (1992) Physiological fluid shear stress causes downregulation of endothelin-1 mRNA in bovine aortic endothelium. The American Journal of Physiology 263(2):C389-96, DOI 10.1152/ajpcell.1992.263.2.C389

Malek A, Gibbons G, Dzav V, Izumo S (1993) Fluid shear stress differentially modulates expression of genes encoding basic fibroblast growth factor and platelet-derived growth factor B chain in vascular endothelium. The Journal of Clinical Investigation 92:2013-2021, DOI 10.1172/JCI116796

Marino M, Pontrelli G, Vairo G, Wriggers P (2017) A chemo-mechano-biological formulation for the effects of biochemical alterations on arterial mechanics: the role of molecular transport and multiscale tissue remodelling. Journal of The Royal Society Interface 14(136):20170615, DOI 10.1098/rsif.2017.0615

Michel J (2004) Système rénine-angiotensine et remodelage vasculaire. Med Sci 20(4):409-413, DOI 10.1051/medsci/2004204409

Millette E, Rauch B, Defawe O, Kenagy R, Daum G, Clowes A (2005) Platelet-derived growth factor-BB-induced human smooth muscle cell proliferation depends on basic FGF release and FGFR-1 activation. Circulation Research 96:172-179, DOI 10.1161/01 RES.0000154595.87608.db

Model LS, Dardik A (2012) Haimovici's Vascular Surgery, John Wiley and Sons, Ltd, chap 10: Neointimal Hyperplasia: Basic Considerations, pp 178-196. DOI 10.1002/9781118481370.ch10, https://onlinelibrary.wiley.com/doi/pdf/10.1002/9781118481370.ch10

Nagaraja S, Wallqvist A, Reifman J, Mitrophanov A (2014) Computational approach to characterize causative factors and molecular indicators of chronic wound inflammation. The Journal of Immunology 192(4):1824-1834, DOI 10.4049/jimmunol.1302481, https://www.jimmunol.org/content/192/4/1824.full.pdf

Newby A (2006) Matrix metalloproteinases regulate migration, proliferation, and death of vascular smooth muscle cells by degrading matrix and non-matrix substrates. Cardiovascular Research 69:614-624, DOI 10.1016/j.cardiores.2005.08.002

O'Connell M, Murthy S, Phan S, Xu C, Buchanan J, Spilker R, Dalman R, Zarins C, Denk W, Taylor C (2008) The three-dimensional micro- and nanostructure of the aortic medial lamellar unit measured using 3D confocal & electron microscopy imaging. Matrix Biology : Journal of the International Society for Matrix Biology 27:171-181, DOI 10.1016/j.matbio.2007.10.008

Ojha M, Leask R, Johnston K, David T, Butany J (2000) Histology and morphology of 59 internal thoracic artery grafts and their distal anastomoses. The Annals of Thoracic Surgery 70:1338-1344, DOI 10.1016/S0003-4975(00)01975-5

Olgac U, Kurtcuoglu V, Poulikakos D (2008) Computational modeling of coupled blood-wall mass transport of LDL: effects of local wall shear stress. American Journal of Physiology-Heart and Circulatory Physiology 294(2):H909-H919, DOI 10.1152/ajpheart.01082.2007

Parks WC, Wilson CL, López-Boado YS (2004) Matrix metalloproteinases as modulators of inflammation and innate immunity. Nature Reviews Immunology 4(8):617-629, DOI 10.1038/nri1418

Plata A, Sherwin S, Krams R (2010) Endothelial nitric oxide production and transport in flow chambers: The importance of convection. Annals of Biomedical Engineering 38:2805-2816, DOI 10.1007/s10439-010-0039-x

Qi Y, Jiang J, Jiang X, Wang X (2011) PDGF-BB and TGF-beta1 on cross-talk between endothelial and smooth muscle cells in vascular remodeling induced by low shear stress. Proceedings of the National Academy of Sciences of the United States of America 108(5):1908-1913, DOI 10.1073/pnas.1019219108

Qiu J, Zheng Y, Hu J, Liao D, Gregersen H, Deng X, Fan Y, Wang G (2014) Biomechanical regulation of vascular smooth muscle cell functions: from *in vitro* to *in vivo* understanding. Journal of The Royal Society Interface 11(90):20130852, DOI 10.1098/rsif.2013.0852, https://royalsocietypublishing.org/doi/pdf/10.1098/rsif.2013.0852

Raines E (2004) PDGF and cardiovascular disease. Cytokine & Growth Factor Reviews 15(4):237-254, DOI 10.1016/j.cytogfr.2004.03.004

Raines E, Ross R (1996) Multiple growth factors are associated with lesions of atherosclerosis: Specificity or redundancy? BioEssays 18(4):271-282, DOI 10.1002/bies.950180405, https://onlinelibrary.wiley.com/doi/pdf/10.1002/bies.950180405

Rectenwald E, Moldawer L, Huber S, Seeger M, Ozaki K (2000) Direct evidence for cytokine involvement in neointimal hyperplasia. Circulation 102:1697-1702, DOI 10.1161/01.CIR.102.14.1697

Rhoads D, Eskin S, McIntire L (2000) Fluid flow releases fibroblast growth factor-2 from human aortic smooth muscle cells. Arteriosclerosis, Thrombosis, and Vascular Biology 20(2):416-421, DOI 10.1161/01.ATV.20.2.416, https://www.ahajournals.org/doi/pdf/10.1161/01.ATV.20.2.416

Robertson A, Sequeira A, Owens R (2009) Cardiovascular mathematics modeling and simulation of the circulatory system, Springer, chap 6: Rheological models for blood

Rodriguez E, Hoger A, McCulloch A (1994) Stress-dependent finite growth in soft elastic tissues. Journal of Biomechanics 27(4):455-467, DOI 10.1016/0021-9290(94)90021-3

Rognoni E, Pisco A, Hiratsuka T, Sipilä K, Belmonte J, Mobasseri S, Philippeos C, Dilão R, Watt F (2018) Fibroblast state switching orchestrates dermal maturation and wound healing. Molecular Systems Biology 14(8):e8174, DOI 10.15252/msb.20178174

Rossi E, Casali B, Regolisti G, Davoli S, Perazzoli F, Negro A, Sani C, Tumiati B, Nicoli D (1998) Increased plasma levels of platelet-derived growth factor (PDGF-BB + PDGF-AB) in patients with never-treated mild essential hypertension. American Journal of Hypertension 11(10):1239-1243, DOI 10.1016/s0895-7061(98)00124-1

San Martin A, Hilenski L, Griendling K (2012) Muscle 2-Volume Set, Elsevier, chap 96: Molecular Pathways of Smooth Muscle Disease, pp 1279-1287

Schwartz R, Chu A, Edwards W, Srivatsa S, Simari R, Isner J, Holmes D (1996) A proliferation analysis of arterial neointimal hyperplasia: Lessons for antiproliferative restenosis therapies. International Journal of Cardiology 53:71-80, DOI 10.1016/0167-5273(95)02499-9

Shampine L, Thompson S (2009) Numerical Solution of Delay Differential Equations, Springer US, Boston, MA, pp 1-27. DOI 10.1007/978-0-387-85595-0_9

Sigg DC (2005) Handbook of Cardiac Anatomy, Physiology, and Devices, Humana Press, chap Cardiac and Vascular Receptors and Signal Transduction, pp 149-159. DOI 10.1007/978-1-59259-835-9_11
Silva T, Jäger W, Neuss-Radu M, Sequeira A (2020) Modeling of the early stage of atherosclerosis with emphasis on the regulation of the endothelial permeability. Journal of Theoretical Biology 496:110229, DOI 10.1016/j.jtbi.2020.110229

Sporn M, Roberts A (1991) Peptide Growth Factors and Their Receptors I. Springer-Verlag New York, DOI 10.1007/978-1-4612-3210-0

Starbuck C, Lauffenburger D (1992) Mathematical model for the effects of epidermal growth factor receptor trafficking dynamics on fibroblast proliferation responses. Biotechnology Progress 8(2):132-143, DOI 10.1021/bp00014a007

Taber L (1998) A model for aortic growth based on fluid shear and fiber stresses. Journal of biomechanical engineering 120:348-354, DOI 10.1115/1.2798001

Tedgui A, Mallat Z (2006) Cytokines in atherosclerosis: Pathogenic and regulatory pathways. Physiological reviews 86:515-581, DOI 10.1152/physrev.00024.2005

Ultman J, Baskaran H, Saidel G (2016) Biomedical Mass Transport and Chemical Reaction: Physicochemical Principles and Mathematical Modeling, John Wiley & Sons, chap 17: Cell Population Dynamics

van Varik B, Rennenberg J, Reutelingsperger C, Kroon A, de Leeuw P, Schurgers L (2012) Mechanisms of arterial remodeling: Lessons from genetic diseases. Frontiers in genetics 3:290, DOI 10.3389/fgene.2012.00290

de Vries M, Simons K, Jukema J, Braun J, Quax P (2016) Vein graft failure: from pathophysiology to clinical outcomes. Nature Reviews Cardiology 13:451-470, DOI 10.1038/nrcardio.2016.76

Wentzel J, Gijsen F, Stergiopulos N, Serruys P, Slager C, Krams R (2003) Shear stress, vascular remodeling and neointimal formation. Journal of Biomechanics 36(5):681-688, DOI 10.1016/S0021-9290(02)00446-3

Widmaier E, Raff H, Strang K (2016) Vander's Human Physiology, Mc Graw Hill, chap 1. Homeostasis, pp 1-16

Wilstein Z, Alligood D, McLure V, Miller A (2018) Mathematical model of hypertension-induced arterial remodeling: A chemo-mechanical approach. Mathematical biosciences 303:10-25, DOI 10.1016/j.mbs.2018.05.002

Zhang M, Zhou Y, Chen L, Wang Y, Wang X, Pi Y, Gao C, Li J, Zhang L (2016) An overview of potential molecular mechanisms involved in VSMC phenotypic modulation. Histochemistry and Cell Biology 145:119-130, DOI 10.1007/s00418-015-1386-3