A new model for evaluating pressure-induced vascular tone in small cerebral arteries

The capacity of small cerebral arteries (SCAs) to adapt to pressure fluctuations has a fundamental physiological role and appears to be relevant in different pathological conditions. Here, we present a new computational model for quantifying the link, and its contributors, between luminal pressure and vascular tone generation in SCAs. This is assembled by combining a chemical sub-model, representing pressure-induced smooth muscle cell (SMC) signalling, with a mechanical sub-model for the tone generation and its transduction at tissue level. The devised model can accurately reproduce the impact of luminal pressure on different cytoplasmic components involved in myogenic signalling, both in the control case and when combined with some specific pharmacological interventions. Furthermore, the model is also able to capture and predict experimentally recorded pressure-outer diameter relationships obtained for vessels under control conditions, both in a Ca2+-free bath and under drug inhibition. The modularity of the proposed framework allows the integration of new components for the study of a broad range of processes involved in the vascular function.


Introduction
Several cerebrovascular pathological conditions involve the malfunctioning of the contractile apparatus in resistance blood vessels (Cipolla et al. 2014b;Good et al. 2018;Fan et al. 2022).In the brain, as well as in many other vital organs such as heart, stomach and kidneys, small arteries and arterioles can actively regulate their diameter to maintain almost constant blood flow in spite of sudden pressure changes.Besides guaranteeing adequate organ perfusion across different physiological conditions, this intrinsic regulation also protects the vascular wall against excessive mechanical stress.Pivotal in this autoregulatory mechanism is the pressure-induced generation of 'myogenic' tone in the vascular SMCs which drives the contraction/dilation of vessel wall (Davis 2012).Early studies (Harder 1984;Knot and Nelson 1998;Osol et al. 2002;Mufti et al. 2010) in small cerebral vessels showed that an augmentation in luminal pressure causes SMC membrane depolarization with consequent intracellular Ca 2+ elevation.Applied mechanical forces can be sensed by the vascular wall through a group of ion channels, whose composition depends on the physiological function of the vessel (Mederos and Schnitzler 2008;Pires et al. 2017;Cui et al. 2022;Harraz et al. 2022).The pressure-induced increase of Ca 2+ in the cytosol enhances the production of myosin light chain kinase (MLCK) that, by phosphorylating the myosin light chains (LC 20 ), favours the actin-myosin interaction and consequently the generation of vascular active tone, leading to eventual constriction.However, experimental evidence (Chrissobolis and Sobey 2001;Cipolla et al. 2002) indicates that Ca 2+ signalling is not the only cellular pathway through which pressure modulates the cell's contractile machinery.It was shown that above certain pressure levels ( ≈ 40-60 mmHg) the cell tone is also regulated by pathways that do not directly affect the level of intracellular Ca 2+ but act rather as parallel mechanisms which can 'change the sensitivity' or 'sensitize' the contractile machinery to Ca 2+ variations and re-organize the cytoskeleton.Pressureactivation of the RhoA/Rho-associated kinase (ROCK) pathway increases the cell contractility by inhibiting the activity of myosin light chain phosphatases (MLCP) which in turn increases the myosin fraction available for the cross-bridges (XBs) formation (Gokina et al. 2005;Johnson et al. 2009;El-Yazbi et al. 2010;Cole and Welsh 2011).High pressure levels also promote actin polymerization via Protein Kinase C (PKC).This cytoskeleton remodelling strengthens the connections between plasma membrane, actin and extracellular matrix which ultimately augments the transmission of contractile force generated by cross-bridges cycling (Gunst and Zhang 2008;Walsh and William 2013;Moreno-Dominguez et al. 2014;El-Yazbi et al. 2015).
To quantify processes underlying disease onset and development, and to virtually assess the efficacy of medical interventions it is mandatory to develop accurate in silico models describing cerebrovascular function.The modelling work by Yang et al. (2003a, b) represents reference studies for describing tone generation and contractility in SCAs.Due to its importance, several authors have devised complex and accurate computational methodologies (Yang et al. 2003a;Kapela et al. 2008;Edwards and Layton 2014) describing the Ca 2+ dynamics regulating SMC tone under different physiological conditions and/or pharmacological inhibition.On the other hand, there is no established SMC contractility model accounting for the ROCK and PKC pathways, which in cerebral arteries are dominant at middle-high pressure levels.Through this work, we lay the foundation towards the realization of a holistic computational framework integrating all known signalling pathways involved in the pressure-induced tone generation.This is achieved by introducing a multi-scale model that is able, for the first time, to quantitatively describe how luminal pressure impacts all major SMC signalling pathways and ultimately the level of tissue contractility.Moreover, the developed framework can be used as an in silico platform for predicting the effect of different types of pharmacological interventions at both intracellular and tissue levels.

Methods
The proposed model allows to evaluate how luminal pressure modulates the vascular tissue tone and consequently the inner diameter of SCAs.This framework is constituted by chemical and mechanical sub-models which are sequentially coupled.The chemical sub-model describes the pressure-induced SMC signalling, with some of its variables representing the input for the mechanical submodel, which translates the active tone information from cellular to tissue level, with consequent diameter adaptation.Model parameter estimation is carried out in two steps by using data under control conditions and pharmacological interventions.

Chemical sub-model
According to previous studies (Lagaud et al. 2002;Lidington et al. 2013), we assume that the pressure-induced pathways are not sequentially activated but they can independently and parallelly contribute to tone generation across wide and overlapping pressure load ranges.Inspired by the recent logic-based model by Irons and Humphrey (2020), we describe the activity level of intracellular ions and proteins involved in the pressuredependent pathways through a signalling network.Activation/inhibition processes (connecting the network nodes) are defined by means of logistic functions.Here, each network node represents either the concentration/ content of an ion/protein or the phosphorylation level of a protein, normalized with respect to the maximum value.With respect to the model by Irons and Humphrey, we introduce some changes that enable us to better capture certain features of the considered SMC signalling.First of all, a generic logistic function ( ) ∈ [0,1] (rep- resenting a intracellular process) depends on the activity level of a node (intracellular variable) ∈ [0,1] via where B 0 represents the basal activity level (which does not need to be necessarily equal to 0, but still < 1 ), whilst n and K are positive coefficients of the logistic function.These three parameters ( B 0 , n and K) need to satisfy the following constraints: (1) and are identified during the model parameter fitting against experimental measurements.In general, if a node ( i ) is affected by multiple independent activation/inhibition processes (i.e.l and m ) of other nodes (i.e.j and k ), the conditional operators 'AND' ( ∧ ), 'OR' ( ∨ ) and 'AND NOT' ( ¬ ) are used: However, if enough information is available, the dynamics of a network node can be described through a specific 'ad hoc' equation that does not necessarily follow the logic-based formalism.The architecture of the intracellular network (i.e. the connections between nodes) is devised based on the current experimental knowledge/evidence on pressure-induced contractility in SCAs (Johnson et al. 2009;Cole and Welsh 2011;Moreno-Dominguez et al. 2014).Hence, we assume (2) (0) = B 0 , (1) = 1, (3) that the stimulus normalized pressure P (with respect to a maximum pressure P max beyond which the generated tone does not increase) is able to modulate the SMC signalling (and consequently the vascular tone) through three distinct pathways, involving cytosolic Ca 2+ , ROCK and PKC signalling (Fig. 1).An increase in luminal pressure (and consequently in transmural stress) leads to a complex Ca 2+ dynamics involving an orchestra of ionic channels which ultimately results in an intracellular Ca 2+ concentration elevation (Jackson 2021).It is worth noticing that cytosolic Ca 2+ increase is required for tone generation since upon extracellular Ca 2+ removal vessels do not exhibit active behaviour (Knot and Nelson 1998;Schubert et al. 2008).Here, for the sake of simplicity, we assume that the effect of pressure on the normalized intracellular Ca 2+ level 0 can be resumed through an activation function  0 ( P) , which is defined within the maximum ranges reported among different experimental studies on cerebral vessels (Knot The symbols ← and ⊢ represent activation and inhibition processes, respectively.The symbol ⋆ indicates that the dynamics of 5 is described through an ad hoc equation that does not follow the logic-based rules of the remaining part of the network, whilst ⇐ implies that the associated variable serves as an input for the mechanical sub-model and Nelson 1998;Osol et al. 2002;Cipolla et al. 2014a;Li and Brayden 2017).The intracellular Ca 2+ level 0 is directly associated with the concentration level of the Calcium-Calmodulin complex, which in turn promotes through MLCK the phosphorylation of myosin LC 20 .This cascade of events/processes is condensed into the Ca 2+ -dependent function 5 ( 0 ) , representing the normalized MLCK capacity to phosphorylate LC 20 .The fraction of phosphorylated LC 20 (pLC 20 ), normalized with respect to the maximum phosphorylatable fraction, is indicated with 5 , whilst 1 − 5 represents the normalized fraction of non- phosphorylated LC 20 that can be phosphorylated.
Also, the pressure-induced ROCK and PKC pathways involve a myriad of proteins and messengers, and here we represent them only through few intracellular variables.The ROCK activity level 1 is up-regulated through pressure-dependent function  1 ( P) and is expected to impact the phosphorylation levels of the proteins MLCP (pMLCP) and Cofilin (pCofilin), indicated with 3 and 4 , through the logistic functions 3 ( 1 ) and 4 ( 1 ) .The phosphorylation of MLCP plays a pivotal role as Ca 2+ -sensitization mechanism as it reduces the MLCP's capacity to de-phosphorylate LC 20 , which is described via 6 ( 3 ) .For describing the combined effect of MLCK and MLCP on LC 20 we do not use the standard logic-based formalism as in Irons and Humphrey (2020).Since MLCK acts on the non-phosphorylated LC 20 fraction, we can write the mass balance of phosphorylated LC 20 fraction (coinciding with the ratio of activated myosin filaments over the total, 5 ) as: where 5 is the characteristic time constant for the dynamic regulating 5 .On the other hand, the phosphorylation of Cofilin impacts cytoskeleton remodelling by down-regulating the level of G-actin content/activity 6 via the logistic function 7 ( 4 ) .Similarly, the pressure-activated PKC path- way is able to modify the cytoskeletal stiffness by promoting (through  2 ( P) ) heat shock protein 27 (HSP27) phosphoryla- tion (pHSP27), indicated with the variable 2 , which in turn down-regulates 6 via 8 ( 2 ) .To recap, the considered intra- cellular variables for this network are the normalized values (with respect to their maximum) of Ca 2+ concentration ( 0 ), ROCK activity level ( 1 ), HSP27 phosphorylation level ( 2 ), MLCP phosphorylation level ( 3 ), Cofilin phosphorylation level ( 4 ), LC 20 phosphorylation level ( 5), and G-actin content ( 6 ).The resulting system of ODEs is reported below: (4) (5) Here, we assume that the myosin pLC 20 coincides with the proportion of myosin filaments able to form cross-bridges with the actin filaments and produce power-stroke.The dimensional fraction of phosphorylated/activated myosin cross-bridges ( n XB ) is retrieved as where n XBmax is the maximum fraction of cross-bridges that can be activated (coinciding with the maximum myosin fraction that can be phosphorylated).The fraction n XB , together with the level of actin polymerization (F-actin content, simply expressed as 7 = 1 − 6 ) serve as input for the mechani- cal sub-model for tissue tone generation.

Mechanical sub-model
The information from the signalling network drives the intracellular contractile machinery, which in turn is translated at tissue level through an active strain-energy function.

Contractile units and actin cortex
For computing tone generation, we follow one of the latest contributions by Murtada et al. (2016), based on previous filament-sliding models (Murtada et al. 2010(Murtada et al. , 2012) ) which have been experimentally validated by using SMCs from blood vessels.With respect to its predecessors, this model incorporates cytoskeleton contribution into the contractile machinery.
Passive surrounding elements (mainly F-actin) constituting the actin cortex are connected with contractile units (CUs) arranged in series to form tissue contractile fibres (CFs).Here, ( 6) where k ACmax is the maximum stiffness under loading condi- tions, whilst n AC and K AC are the coefficients of the associ- ated activation function.With respect to the work by Murtada et al. (2016), we assume that even if cross-bridges in the latched state exist, they do not contribute to the active tone generation.In this way, the total stiffness of all attached cross-bridges in half of a CU ( k tCB ) can be simply evaluated as where L m is the average length of myosin filaments, m is the average distance between myosin monomers heads and k XB represents the elastic stiffness of the cross-bridge.As in Murtada et al. (2016) we define the total stiffness associated with a number of CUs ( N CU ) as where Lfo is a normalized function representing the overlap between thin and thick filaments, to be introduced later on.The length of the CU in the deformed configuration ( l SMC ) defined as where L SMC , u CB , u fs and u FA are all quantities evaluated in the reference configuration.L SMC is the length of the CU, u CB is the average XB elongation, u fs is the relative filament sliding and u FA is the average surrounding passive ele- ments elongation.From ( 16) the circumferential stretch = l SMC ∕L SMC can be also defined as where ūCB = u CB ∕L SMC , ūfs = u fs ∕L SMC and ūFA = u FA ∕L SMC .
According to Murtada et al. (2016), the resulting reaction force due to the resistance from N CU in series with a F-actin element at each extremity is given as total displacement by the resulting stiffness of the system: whilst the average driving force generated from the cycling XBs results (13) where u PS is the average displacement associated with power-stroke.The normalized filament overlap is described through a Gaussian function (Murtada et al. 2016) depending on the (normalized) relative filament sliding ūfs : where ūopt fs is the optimal filament sliding and s f0 is a scaling parameter.A simple evolution law is used for describing the kinetics of the filament sliding caused by the XB cycling: where c is the time constant associated with the process.By considering (17), the active force F a can be re-written as a function of and ūfs : The generated tone from CFs can be quantified through the active First Piola-Kirchhoff stress P a of the tissue material, which can be evaluated as where N CF is the surface density of the contractile fibres.

Vessel wall
The considered vascular tissue, consisting of different SMCs layers, is modelled as an axisymmetric homogeneous hyperelastic thick-walled tube.The luminal pressure P can be defined through the momentum conservation principle along the radial direction (Coccarelli et al. 2021) where P ext is the external pressure acting on the outer surface of the vessel (which in this study is assumed negligible), Ψ is the material strain-energy function, R i and H are the inner radius and the thickness in the reference configuration Ω R , r is the radial coordinate in the deformed configuration Ω D , whilst z and r = R rk z are the axial and radial stretches, respectively, with k accounting for the residual circumfer- ential strain.Moreover, we assume that there is no difference in thickness between the load-free and reference configurations and the material is incompressible.This latter constraint allows the following mapping ( 19) , where r i is the deformed luminal radius.As usual, strain- energy function is decomposed into a passive and active components ( Ψ = Ψ p + Ψ a ), with the former defined as (Gasser et al. 2006) where c 0 , c 1 and c 2 are the media constitutive parameters, z sin 2 with being the orientation angle (with respect to the circumferential direction) of a collagen fibres family.In line with Murtada and Holzapfel (2014), we assume that these collagen fibres are aligned with the contractile fibres CFs, which are oriented along the circumferential direction of the vessel ( = 0 • ).By assuming that collagen fibres cannot sustain compression loads, the anisotropic contribution ( I 4 -1) is set to zero when I 4 < 1.On the other hand, the strain-energy function component representing active tone generation in the SMC tissue can be obtained as in Murtada et al. (2010); Murtada and Holzapfel (2014)

Model parameter identification and solution procedure
The parameters of the chemical and mechanical sub-models were separately identified through a two-tiered approach by solving two optimization problems, which were defined in line with Coccarelli and Pant (2023).For the mechanical sub-model optimization problem, the best-fit parameters of the chemical sub-model were used.For each optimization problem, the solution is the set of model parameters that minimizes the discrepancy (through a cost function) between the simulated data and the corresponding measurements across a vessel population.The cost function is made of two components ( L tot =L data +L sum con ), with L data is the L2 norm of the simulated data with respect to average values from experiments, while L sum con is a penalty-based term accounting for the model constraints that variables/parameters must satisfy.These are introduced to refine and accelerate the search by discarding solutions which do not lie within the expected physical value ranges.Each constraint contribution is included into the cost function as an individual penalty L i con = Δy con , where is the penalty param- eter (1.0e3) and Δy con defines the absolute distance of the variable with respect to the allowed interval boundary.The parameter set is identified by minimizing the cost function through the covariance matrix adaptation evolution strategy (CMA-ES) (Hansen et al. 2009;Hansen 2016).In the CMA-ES settings, the maximum number of iterations 'maxiter' per search is set equal to 1.0e3, with a standard deviation 'sigma0' not exceeding 1. Due to the stochastic nature of the search algorithm (and the high number of model parameters), multiple runs were carried out to find best model fits.
In this study, we focus on the steady-state behaviour of the model.Hence, the associated time constants ( 0 , 1 , 2 , 3 , 4 , 5 , 6 and c ) are all set equal to 1.For any pressure level, the static solution of the chemical submodel nodal variables was obtained by solving the system of ODEs ( 5)-( 11).Regarding the mechanical sub-model, we only consider the ūfs evaluated with a circumferential stretch averaged over the vessel thickness, which is representative of the mean filament sliding occurring across the wall.Furthermore, ( 24) is re-arranged as a function of the luminal radius in the deformed configuration r i = r i (P, ūfs ) , which can be internally evaluated within the solution procedure for computing dū fs dt = dū fs dt (ū fs , r i (P, ūfs )) via ( 21).Hence, the latter equation is solved at steady-state for obtaining r i and ūfs at each pressure level.For the integral in (24) Simpson's rule is used.Both sub-models were implemented in Python3.9 and all equations were solved by means of the function 'root' (method ='lm') from the SciPy library (version 1.6.0),while for the optimizer we used the function 'fmin2' from the cma library (version 3.2.2).

Cell signalling sub-model
Here, we report how the chemical sub-model is able to capture the experimental data by Johnson et al. (2009), El-Yazbi et al. (2010), Moreno-Dominguez et al. (2014) obtained by employing pressure myography on middle cerebral arteries from Sprague Dawley rats.These studies report the levels of pHSP27, pMLCP, pCofilin, pLC 20 and G-actin for different pressures (ranging from 10 to 140 mmHg), under control and pharmacologically manipulated conditions.For this purpose the selective ROCK inhibitor H-1152 dihydrochloride (H1152) and PKC inhibitor GF109203X (GF) were used.For the signalling network variables normalization, we considered a maximum pressure P max equal to 150 mmHg, beyond which all the considered intra-cellular processes reach saturation level.H1152 doses of 0.3 M and 0.5 M are expected to inhibit the ROCK pathway and were simulated by setting the network variable 1 equal to 0.12 and 0.07, respectively (in line with the H1152 dose-inhibition curve reported in Johnson et al. (2009)).On the other hand, 3.0 M GF is expected to block the PKC pathway, and therefore, this condition was simulated by setting the network variable 2 equal to 0.1.The experimental data used for the model parameters identification are summarized in the Appendix.
Regarding the optimization procedure constraints, we ensured that the pressure-cytosolic Ca 2+ curve (  0 ( P) ) lies within the maximum ranges observed among the experimental data reported in Knot and Nelson (1998), Osol et al. (2002), Cipolla et al. (2014a), Li and Brayden (2017) (see Appendix).These experimentally obtained curves exhibit different features but for each of them the EC 50 appears to be significantly higher than 15 mmHg.Here to refine and accelerate the search we add the constraint which penalizes parameters sets providing an EC 50 for  0 ( P) that is lower than 0.1.We also impose that, in case of extracellular Ca 2+ removal (simulated by setting the network variable 0 equal to 0.1), the pLC 20 cannot exceed the value in the control case at 10 mmHg, for pressure increasing up to 120 mmHg.In line with the recordings by Johnson et al. (2009), n XBmax was assumed to be 0.55 (this is necessary for passing the information from the chemical sub-model to the mechanical one).
These settings enabled us to identify the coefficients ( B 0 , K and n) of the nine logistic functions (composing the signalling network) which best capture the above-mentioned experimental data (Johnson et al. 2009;El-Yazbi et al. 2010;Moreno-Dominguez et al. 2014) (see Table 1).
The logistic functions constituting the cell signalling network are reported in Fig. 2. Figure 3 shows how the cell signalling sub-model is able to capture the influence of pressure on the pMLCP and pLC 20 levels in the control case, under 0.3 M H1152 and 3.0 M GF.While for pMLCP values are normalized with respect to the 10 mmHg control case, absolute values are reported for pLC 20 fractions.The experimental pLC 20 value under 3.0 M GF was re-scaled (multiplied by a factor ≈ 1.1 ) by considering the difference between the two control values at 100 mmHg reported in Johnson et al. (2009) (Figs. 8B and 10E therein).
The effect of pressure on the intracellular pathways governing the cytoskeleton remodelling is reported in Fig. 4. The defined chemical sub-model is able to reproduce the pressure-induced pHSP27 elevation in the control case and under 0.5 M H1152.Similarly, the relationships between pressure and pCofilin, and pressure and G-actin levels are well captured in the control case, under 0.5 M H1152 and 3.0 M GF.

Vascular tissue sub-model
Here, we considered other vessels data from Johnson et al. (2009) to fit the mechanical sub-model and validate the All the reported values are normalized with respect to their control case at 10 mmHg.In the 3.0 M GF case, pHSP27 ( 2 ) was set equal to 0.1 across the whole pressure range to simulate the drug effect global methodology.The simulated protocol consists of tube inflation tests at different pressure levels under fixed axial length, with the addition/removal of vasoactive substances.
A first group of five vessels (here renamed 'vessels group 1') was used for comparing the H1152 effect on diameter with respect to the control and extracellular Ca 2+ removal cases.Similarly, GF effect was tested on another group of six arteries (here referred as 'vessels group 2').From the reported experimental data, the only main difference between these two vessels groups seems to be the average size, with the group 1 having a slighter smaller radius near to the load-free condition ( R ′ o =0.011 cm for vessels group 1, R ′ o =0.0126 cm for vessels group 2).For the vessels used in Johnson et al. (2009), the load-free thickness-medium radius ratio ( h w ), the axial stretch ( z ) and the circumferen- tial pre-stretch ( k ) were not reported.We used the pressure- diameter data derived from the assessment on vessels group 1 for identifying (through the methodology reported in Sect.2.3) these geometric factors alongside some passive and active response parameters ( c 0 , c 2 , u PS , Other model parameters are expected to vary less among SMCs, and therefore, they were assumed equal to what is reported in the literature (Murtada and Holzapfel 2014;Murtada et al. 2016).The chemical submodel parameters (as well as the way in which drug interventions are reproduced) are the same as in Sect.3.1.Extracellular Ca 2+ removal dramatically reduces the intracellular Ca 2+ level due to Ca 2+ influx abolition.This condition is modelled by setting the network node 0 equal to 0.1.To accelerate and refine the search, we penalized (as reported in Sect.2.3) any set of parameters for which i) |ū fs | exceeds 0.1 and/or ii) the derivative of the outer diameter (in the deformed configuration) D o with respect to P under control conditions is positive across the range 60-100 mmHg.All the parameters of the mechanical sub-model (including the ones obtained though fitting against the experimental curves) are summarized in Table 2.
It is worth noticing that with this set of parameters the system reaches steady-state conditions for all the considered pressure levels (see Appendix).Geometric factors such as axial stretch ( z ) and circumferential pre-stretch ( k ) may vary quite a lot from vessel to vessel in the brain and depend of course on the experimental procedure (Monson et al. 2008).Furthermore, for resistance cerebral vessels, the role of axial stretch on their pressure-induced vasoreactivity is  not fully clarified (Coats and Hillier 1999).The identified parameters need therefore to be interpreted as representative of the mean behaviour of the vessels group.
Figure 5 shows how the model captures the pressureinduced tone under control, extracellular Ca 2+ removal and H1152 inhibition conditions.In the control case, vessels are able to limit variations in their diameter despite pressure elevation.On the other hand, by removing extracellular Ca 2+ , we obtain a situation very close to the pure passive vessel behaviour.The H1152 case shows the impact of ROCK inhibition on the SMC contractile machinery.In the latter case tone generation is compromised similarly to the case of extracellular Ca 2+ absence.Then, we used the pressure-outer diameter data concerning the vessels group 2 for validating the chemo-mechanical model.This group of vessels appeared to be slightly bigger (different load-free radius) and was subjected to 3.0 M GF (simulated as in Sect.3.1).Overall, for these vessels, predictions fall within the experimentally observed ranges (Fig. 6).It is important to notice that these arteries may have different geometric and passive material properties with respect to the vessel group 1.Furthermore, we assumed that GF impacts the contractile machinery only through the PKC pathway, without considering potential secondary effects on the remaining signalling network.These may be some of the reasons for the observed discrepancies between simulated behaviour and experimental mean values in the pharmacologically manipulated cases.
In Fig. 7, we show how the (normalized) relative filament sliding changes upon pressure elevation across all the considered conditions and vessels groups.Interestingly, in the control cases, the sliding between filaments is limited (slight contraction) the whole considered pressure range.Extracellular Ca 2+ removal hinders the formation of new phosphorylated XBs, leaving the vessel to expand upon pressure loading.On the contrary, alteration of the cytoskeleton function through drugs seems to have the opposite effect (filaments highly contracted).
Pressure-induced signalling pathways ultimately alter the stiffness of contractile units and actin-cortex (Fig. 8).The extracellular Ca 2+ removal case, as expected, has a profound limiting effect on pressure-induced CU stiffness, while the actin cortex behaviour coincides with the control one.Both drug interventions almost nullify the stiffness of the actin cortex across the considered pressure range.On the other hand, while inhibition of the ROCK pathway limits pressure-induced CU stiffness, the GF case exhibits a stiffness increase with pressure that is even slightly higher than the control one.This may be due to the fact that CU stiffness also depends on ūfs , which is significantly differ- ent between the GF and the control case.Although some of the reported numerical findings may appear intuitive, these still need to be corroborated by future experimental studies.

Discussion and concluding remarks
Through this work we introduced a new multi-scale model for computing pressure-induced tone and resulting diameter in SCAs.Luminal pressure is converted into vessel tone through Fig. 5 Mechanical model fitting for vessels group 1 under control, extracellular Ca 2+ removal and 0.3 M H1152 conditions.D o is the outer diameter in the deformed configuration.The experimental data are reported as mean value ± standard deviation a framework made of a signalling network, accounting for the intracellular Ca 2+ , ROCK and PKC pathways, which feeds a mechanical sub-model combining actin-myosin interaction with cytoskeleton remodelling to evaluate tissue contractility.These sub-models were fitted and validated against experimental data by considering control conditions and with selective intracellular pathway manipulations.The model was also used for assessing how stiffness of the contractile machinery components vary with pressure and with/without pharmacological intervention.So far, the vast majority of vascular SMC contractility models (Yang et al. 2003a;Murtada et al. 2010;Carlson and Beard 2011;Murtada et al. 2012;Edwards and Layton 2014;Murtada et al. 2016;Coccarelli et al. 2018) describe the XB kinetics by using the so-called 'Hai and Murphy model' (Hai and Murphy 1988), which hypothesizes a contribution to force generation by XBs in a latched state.Through this study, we showed that this four-state kinetic model is not strictly necessary for reproducing the behaviours observed in experiments, as long as MLCP and cytoskeleton regulation are accounted for.Upon a more detailed experimental characterization of the XB dynamics, complex models such as (Hai and Murphy 1988) could also be integrated in the proposed framework.
The main limitations of the study are intrinsically associated with the uncertainties regarding different cellular processes and cell-tissue signal transduction.Some alleged secondary signalling connections between cellular components which so far have been only partially explored may play a role in pressureinduced contractility and therefore require further elucidation.An earlier study (Gokina et al. 2002) reported indeed that, at intermediate pressure levels, drugs which prevent/inhibit actin polymerization had also an effect on the cellular Ca 2+ influx by modifying the activity of voltage gated Ca 2+ channels.Mufti et al. (2010) proposed that Ca 2+ waves originated from sarcoplasmic reticulum can also modulate MLCP activity.However, the lack of further experimental evidence together with the scarce mention by reference works (Cole and Welsh 2011;Walsh and William 2013) made us consider these connections of secondary importance (and therefore omitted).Ca 2+ dynamics is pivotal for the initiation of tone generation (Jackson 2021) and certainly its modelling representation can be enhanced.A more complete experimental characterization of these processes will certainly help us to model with higher accuracy selective interventions targeting cell contractility (such as GF).Furthermore, we mention that the proposed model does not include yet a detailed link between integrins and Ca 2+ sensitization processes (Colinas et al. 2015), which could be extremely useful in the future for simulating drugs that impact the signalling between SMC and extra-cellular matrix.
Given the high number of model parameters to be identified, we conducted numerous CMA-ES runs, which yielded multiple solutions in the parameter space.Among these, we identified the set of parameters which minimizes the associated cost function and that is most likely close to the global minima of the optimization problem.While a systematic identifiability study can be performed in the future, additional experimental data may be necessary to guarantee well-posedness of the inverse problem and uniqueness of parameters.The chemical sub-model was developed by considering only the data from one laboratory (Johnson et al. 2009;El-Yazbi et al. 2010;Moreno-Dominguez et al. 2014) and therefore a broader validation by using new data from different research groups will be necessary for reinforcing its predictive power and robustness.On the other hand, the mechanical sub-model of our framework is based on the very established work by Murtada et al. (2010Murtada et al. ( , 2012Murtada et al. ( , 2016) ) and, in the future, it could be further enhanced by using displacement-active tension recordings from wire myograph settings.All the presented results were at steady-state conditions, but we are currently investigating how the model is able to reproduce and predict time-dependent behaviours.
Being one of the earliest attempts to quantify the major known contributors to vascular SMC contractile machinery in SCAs, the proposed minimalistic modelling approach aims to provide reliable predictions without sacrificing its interpretability.However, as per the methodology proposed in Irons and Humphrey (2020), Irons et al. (2021), new components can be easily incorporated into the framework once new experimental data are made available.In addition, differences in pressure-induced vasoreactivity between cerebral resistance arteries and arterioles seem mainly due to different Ca 2+ activity rather than Ca 2+ sensitivity (Cipolla et al. 2014a).Hence, upon adequate characterization of the pressure-Ca 2+ relationship, the current model could be also used for modelling the behaviour of small cerebral vessels.Finally, since alterations in some of the considered cellular pathways (such as ROCK) are associated with different vascular pathological Fig. 8 Simulation results-role of luminal pressure on total CU stiffness and actin cortex stiffness.ktCU is the total CU stiffness normalized with respect to the corresponding control value.kAC is the actin cortex stiffness normalized with respect to the corresponding control value conditions (Chrissobolis and Sobey 2001;El-Yazbi and Abd-Elrahman 2017), the proposed methodology may also constitute an in silico basis for the design and testing of novel therapies.

Experimental data for chemical sub-model parameters identification
Experimental data for the identification of the chemical sub-model parameters are reported in

Boundaries of pressure-Ca 2+ relationship
For the chemical sub-model parameter identification, we imposed that the curve representing pressure-Ca 2+ relationship (  0 ( P) ) must be confined within some defined bounda- ries.Here, we considered the experimentally recorded traces from different works (Knot and Nelson 1998;Osol et al. 2002;Cipolla et al. 2014a;Li and Brayden 2017) and we normalized them with respect to their value at 10 mmHg.Among these, the experimentally derived curve by Osol et al. (2002) had the highest values across the considered pressure range.For our normalized pressure-Ca 2+ function, we defined as upper boundary the data (mean value plus standard deviation) from the latter study (Table 6).
On the other hand, since pressure causes an increase in cytosolic Ca 2+ , we imposed that across the considered pressure range the normalized pressure-Ca 2+ function cannot go below 1.0.

Mechanical sub-model steady-states
In Fig. 9 the dependency of dū fs dt on ūfs for different pressure levels is illustrated.

Fig. 1
Fig. 1 Chemical sub-model represented as signalling network.The symbols ← and ⊢ represent activation and inhibition processes, respectively.The symbol ⋆ indicates that the dynamics of 5 is described through an ad hoc equation that does not follow the logic-based rules of the remaining part of the network, whilst ⇐ implies that the associated variable serves as an input for the mechanical sub-model

Fig. 3 Fig. 4
Fig.3Chemical model fitting -pMLCP and pLC 20 .The experimental data are reported as mean value ± standard deviation.The reported pMLCP values are normalized with respect to their control case at 10 mmHg.The experimental measurements of pMLCP and pLC 20 under 3.0 M GF inJohnson et al. (2009) were reported with respect to

Fig. 6 Fig. 7
Fig.6Mechanical model validation for vessels group 2 under control, extracellular Ca 2+ removal and 3.0 M GF conditions.D o is the outer diameter in the deformed configuration.The experimental data are reported as mean value ± standard deviation

Table 2
Mechanical model parameters.It is noticed that in the current model formulation only the ratio (and not the absolute values) between k XB and m is required.Some of the fitted parameters were restrained between specific ranges of values: 0.13 ≤ h w Table 3 (control conditions), Table 4 (H1152 stimulation) and Table 5 (GF stimulation).

Table 3
Experimentally recorded mean values of normalized (with respect to 10 mmHg control case, indicated with superscript '0') pHSP27, pMLCP, pCofilin, pLC 20 and G-actin under control conditions used for model fitting

Table 4
Experimentally recorded mean values of normalized (with respect to 10 mmHg control case, indicated with superscript '0') pHSP27, pMLCP, pCofilin, pLC 20 and G-actin under 0.3 M H1152 inhibition used for model fitting

Table 5
Experimentally recorded mean values of normalized (with respect to 10 mmHg control case, indicated with superscript '0') pHSP27, pMLCP, pCofilin, pLC 20 and G-actin under 3.0 M GF inhibition used for model fitting.*Obtained by re-normalizing the original experimental data with respect to 2 @ 10 mmHg control case.**obtained by re-scaling the original experimental data with respect to the reference 5 @100 mmHg control case (scaling factor ≈ 1.1)

Table 6
Osol et al. (2002)fromOsol et al. (2002).Ca 2+ data, reported as mean value ± standard deviation, was normalized with respect to the value at 10 mmHg