An integrated lumped-parameter model of the cardiovascular system for the simulation of acute ischemic stroke: description of instantaneous changes in hemodynamics

Acute Ischemic Stroke (AIS) is defined as the acute condition of occlusion of a cerebral artery and is often caused by a Hypertensive Condition (HC). Due to its sudden occurrence, AIS is not observable the right moment it occurs, thus information about instantaneous changes in hemodynamics is limited. This study aimed to propose an integrated Lumped Parameter (LP) model of the cardiovascular system to simulate an AIS and describe instantaneous changes in hemodynamics. In the integrated LP model of the cardiovascular system, heart chambers have been modelled with elastance systems with controlled pressure inputs; heart valves have been modelled with static open/closed pressure-controlled valves; eventually, the vasculature has been modelled with resistor-inductorcapacitor (RLC) direct circuits and have been linked to the rest of the system through a series connection. After simulating physiological conditions, HC has been simulated by changing pressure inputs and constant RLC parameters. Then, AIS occurring in arteries of different sizes have been simulated by considering time-dependent RLC parameters due to the elimination from the model of the occluding artery; instantaneous changes in hemodynamics have been evaluated by Systemic Arteriolar Flow (Qa) and Systemic Arteriolar Pressure (Pa) drop with respect to those measured in HC. Occlusion of arteries of different sizes leaded to an average Qa drop of 0.38 ml/s per cardiac cycle (with minimum and maximum values of 0.04 ml/s and 1.93 ml/s) and average Pa drop of 0.39 mmHg, (with minimum and maximum values of 0.04 mmHg and 1.98 mmHg). In conclusion, hemodynamic variations due to AIS are very small with respect to HC. A direct relation between the inverse of the length of the artery in which the occlusion occurs and the hemodynamic variations has been highlighted; this may allow to link the severity of AIS to the length of the interested artery.


Introduction
Acute Ischemic Stroke (AIS) is defined as the acute condition of occlusion of a cerebral artery. Since blood can flow no more through that artery, a region of the brain suffers of blood perfusion loss. Many risk factors for AIS have been identified, the main of which is the Hypertensive Condition (HC), present in most of the patients experiencing AIS [1]. HC is defined as a chronic situation in which systolic Blood Pressure (BP) is higher than about 130 mmHg, depending on age, sex, smoke, weight and general personal health conditions [2]. However, 56% of AIS patients presented systolic BP higher than 140 mmHg, while small fractions presented high diastolic BP and/or high mean BP [1]. HC on long terms leads to a modification in arterial radius, vessel thickness, Young modulus and thus changes vessel impedance [3]. Knowing under which initial conditions mechanisms of pressure control work is important to evaluate their action and then to improve therapy against AIS risk factors and AIS itself.
In pre-hospital care, AIS causes unique effects that allow to identify it with simple protocols [4], while quantitative hemodynamics remains unobserved due to the lack of specialistic instrumentation the right moment an AIS usually occurs. Even in hospital care, those observations can be difficult due to a time delay of seconds up to some minutes between the triggering of AIS and its assessment through scanning and imaging techniques [5]. Hemodynamic instantaneous changes in occluded human arteries when AIS takes place are then quantitatively assessable only via animal models or, in the last decade, via computational models.
Computational models of the cardiovascular system are particularly useful in simulating non directly observable events and models with different levels of detail have been proposed in relation to the intended purpose; 3D models, also known as Computational Fluid Dynamics (CFD) model, allow a highly-detailed description of a precise vessel trait [6]; less detailed description may be obtained through 2D, 1D or 0D models, the latter also defined as Lumped-Parameter (LP) models [7][8][9][10][11][12]. Of note, LP models may be suitable to describe the cardiovascular system during AIS in which the time window of interest is limited to the very first seconds after the occlusion, and thus wave reflection and transmission effects are not considered. Moreover, characteristic dimension of the occluding vessel allows to consider the flow in that vessel steady due to the fact that the pulsatility of blood induced by the periodic pumping of the heart is largely damped by the compliance of the large arteries [13].
A previously proposed LP model has been used to study cardiac flow in certain healthy and diseased heart valves [10] but it has never been used to simulate AIS. Thus, the aim of this study was to simulate AIS and describe instantaneous changes in hemodynamics through an integrated LP model of the cardiovascular system, derived from the model proposed in [10]. To this aim, three different regimens have been considered: starting from the physiological regimen, HC regimen considered modifications introduced by this status and then evolution towards AIS regimen has been described.

The integrated lumped-parameter model of the cardiovascular system
The integrated LP model of the cardiovascular system is composed by four heart chambers (left ventricle; right ventricle; left atrium; right atrium), four heart valves (mitral, aortic, tricuspid, pulmonary), Systemic Circulation (SC), Pulmonary Circulation (PC). Heart chambers have been described with an elastance model integrating approaches proposed in [10,12]; valves have been described as static systems with binary (open/closed) states according to [10]; SC and PC have been described with a LP model that accounts for vessel viscous resistance, fluid inertance and vessel wall elasticity.

Heart chambers
The heart chambers have been described with an elastance model integrating approaches proposed in [10,12]. The mathematical relation between volume of the chamber and flow through it is described in Eq (1) expressing mass conservation for incompressible fluids: where V(t) is the internal volume of the chamber, Qin is the inlet volumetric blood flow and Qout is the flow out. The volume of the chamber can be linked to internal pressure by Eq (2): where P is the pressure inside the chamber [14]. ℰ(t) is an elasticity time-dependent function that can be described by Eq (3) [14]: where ℰmax and ℰmin are characteristic parameters of each chamber. The function e(t) is the activation function that varies between 0 (total diastole) and 1 (systolic peak) and has two different descriptions depending on the chamber and the beginning of cardiac cycle, here assumed at time equal to 0. For the ventricles, the elasticity function ev corresponds to Eq (4), while for the atria ea is given by Eq (5) [10,15].
Equation (2) can be integrated between V0 (i.e., generic initial volume) and V(t) obtaining the first term of Eq (6). A term accounting for the resistive behavior of the chambers, indicated with Rc(t), and expressing the dependence of P(t) on Q(t) is also considered in Eq (6), as already done in [12]: The sign of the resistive term is negative because it dissipates energy from the system. Rc(t) is a viscous resistance shown to be proportional to the first term of Eq (6) [16] as described in Eq (7): The proportionality constant Kc is characteristic for each chamber and proportional to the pressure achieved with isometric contraction [16].

Cardiac valves
Cardiac valves have been described as static pressure-controlled systems with binary (open/closed) states according to [10]. Because of this, the retrograde flow characteristic of the closing process can't be studied. They are controlled by pressure difference across the leaflets, here indicated with ΔP. The blood flow Q through a valve is here described with Eq (8): where Q0 is a parameter dependent on valve dimensions and A(ΔP) is described by Eq (9):

Systemic circulation and pulmonary circulation
As described in similar works about LP models [17][18][19], SC and the PC have been described with different traits of LP models that account for vessel viscous resistance, fluid inertance and vessel wall elasticity. The SC has been divided into: systemic aortic sinus, systemic arteries, systemic arterioles, systemic capillaries, systemic veins. The PC has been divided into: pulmonary aortic sinus, pulmonary arteries, pulmonary arterioles, pulmonary capillaries, pulmonary veins.
Blood flow in larger vessels (i.e., systemic aortic sinus, systemic arteries, pulmonary aortic sinus and pulmonary arteries) is pulsatile, so vessel viscous resistance, fluid inertance and vessel wall elasticity were taken into account [17]. In arterioles and capillaries, the flow can be considered steady and the vessel viscous resistance is dominant on fluid inertance and vessel wall elasticity [13]; thus, the latter two can be neglected. The magnitude of flow oscillations is much smaller in systemic veins and pulmonary veins and thus, it can be considered steady [17]; as a consequence, fluid inertance is neglectable and only viscous resistance and vessel wall elasticity have been taken into account.

Model implementation
The integrated LP model of the cardiovascular system has been implemented in a closed-loop electrical circuit in MATLAB SIMULINK R2019b (MathWorks, Natick, MA, USA) and the Simscape Electrical™ package has been used to implement the electrical analogue as in [17,18]. Model implementation has been obtained by linking in series heart chambers, valves, SC and PC. Units of variables and parameters used in the model have been described in Table 1. Moreover, since the model uses medical units while its implementation uses electric units, conversion between the two and the reference units in the International System has been described in Table 2. Heart chambers and cardiac valves have been implemented according to their governing equations as described in Sections 2.1.1 and 2.1.2, respectively. Electrical analogue for SC and PC has been implemented as a direct resistorinductor-capacitor (RLC) circuit that takes into account vessel viscous resistance, fluid inertance and vessel wall elasticity. On the basis of physiological considerations done in Section 2.1.3, larger vessels (i.e., systemic aortic sinus, systemic arteries, pulmonary aortic sinus and pulmonary arteries) have been implemented as complete direct RLC circuits; arterioles and capillaries have been implemented with a simple resistor each; systemic veins and pulmonary veins have been implemented with a resistor and a capacitor each. The general model implementation is shown in Figure 1 (panel a).

Physiological, hypertensive condition and acute ischemic stroke regimens
The model has been considered working in three different regimens: 1) physiological; 2) HC; 3) AIS. The three regimens differ for the implementation of systemic arteries and for pressure inputs. In the physiological regimen and HC regimen, systemic arteries have been implemented as complete direct RLC circuits ( Figure 1, panel b), but differing in the set of constant parameters used. In the AIS regimen, systemic arteries have been implemented with complete direct RLC circuits with timedependent parameters since, when AIS occurs, the branches of the occluding artery must be eliminated from the model (i.e., open circuit), thus the RLC implementation of systemic arteries changes in time, as shown in Figure 1 (panel b); because of the complex net of cerebral vessels, the specific vasculature of the brain has not been followed. A parallel connection between occluding artery and the rest of the systemic arteries has been considered to model the variation of RLC implementation when the occluding artery is eliminated.
The three regimens have the same activation functions for ventricles and atria, as shown in Figure 2. Moreover, they have also the same pressure input for left atrium, right atrium and right ventricle, while the pressure for left ventricle changes from physiological regimen with respect to HC and AIS regimens, for which it is the same (Figure 3).

Model parameters for physiological, hypertensive and acute ischemic stroke regimens
Different values for model parameters have been used for the three regimens. After comparing different studies [7][8][9][10]12,20,21], the final parameters for the physiological regimen were derived from [10,12]. Starting from parameters values of the physiological regimen and considering anatomical similarity between human and rat, parameters for HC and AIS have been derived considering progressive anatomical pathological reshaping and related R, L and C values ha ve been computed using Eqs (10), (11) and (12) [22]: in which l is the length of the lumped artery/vein, r is its radius, E is its elastic modulus, h is its thickness, A is the cross section, μ is the cinematic viscosity of blood, considered as 4·10 -3 (Pa·s), and ρ is blood density, considered constant and equal to 1060 (kg/m 3 ).
In particular, HC in rats has been estimated to increase vessel viscous resistance of 50% with respect to the physiological condition [23]. This increase corresponds to a reduction of arterial radius of 10%, as calculated with Eq (10). Because of this, RLC parameters of systemic arteries of SC were recalculated in HC according to Eq (10) with a 10% reduced radius than the physiological regimen; since in HC regimen arteries change their properties adapting their internal radius, wall thickness and elasticity [24], it has been supposed an overall 5% reduction of compliance with respect to physiological regimen. For the PC no modification has been applied since the pressure is very low also in HC in this part of the vasculature. It has also been assumed that the parameters related to the cardiac valves and the heart chambers are not modified in HC. The dimension of the occluding artery can be different based on the dimension of the occluding agent. In order to take into account this variability, ten independent occluding arteries (A1-A10) of different lengths have been considered based on the anatomical values derived from [20]. The set of considered arteries were generated under the assumption that the length of the occluding artery l is the parameter that most influences anatomical variability with respect to internal radius, vessel wall thickness and parameters present in Eqs (10), (11) and (12). As such, length l has been assumed as the main parameter able to determine the gravity of an AIS. The remaining parameters of internal radius, vessel wall thickness and elastic modulus have been assumed constant All the numerical values of the parameters are listed in Table 1 (for the physiological and HC regimens) and in Table 3 (for AIS regimen).  For each lumped occluding artery: l is the length; r is the radius; E is elastic modulus; h is the thickness; R is the viscous resistance; L is the fluid inertance; C is the vessel wall elasticity.

Calculation
A variable step solver has been used to solve model equations. Simulation of the physiological regimen has been performed for 5000 s, corresponding to 5000 cardiac cycles. The transitive response was limited to the first 3 cycles, after which the response became periodic. Then, simulation of HC and AIS regimens has been done for 20 s; in particular, HC has been simulated during the first 10 s whereas AIS during the last 10 s. Change from HC to AIS regimen has been simulated as an instantaneous event and change in model parameters values has been simulated with a step function at time equal to 10 s. Systolic, diastolic and mean BP have been measured downstream to every element of the SC only for physiological regimen by means of a voltage sensor (Figure 1, panel a). For both physiological and HC regimens, Aortic Pressure (AP) and Cardiac Output (CO) have been computed at the outlet of the aortic valve by means of a voltage and a current sensor. Systemic Arteriolar Flow (Qa) and Systemic Arteriolar Pressure (Pa) are assumed as output of the model, and the difference between Qa at HC at the last peak before AIS occurrence and Qa at the first peak after AIS occurrence has been indicated as ΔQa. Analogously, the difference between Pa at HC at the last peak before AIS occurrence and Pa at the first peak after AIS occurrence has been indicated as ΔPa. The values of Qa and Pa have been measured only in AIS regimen with current and voltage sensors downstream with respect to systemic arteries ( Figure 1, panel a). Those measurements were then used to compute ΔQa and ΔPa.

Model validation and relation of systemic arteriolar flow and pressure with length of artery interested by acute ischemic stroke
The measured values of BP in the SC were then used to validate the model by comparing the BP trend in the connected elements with physiological data [25]. The values of AP and CO in physiological regimen have also been compared to physiological data [26] as supporting evaluation of model validation. To evaluate association of ΔQa and ΔPa with anatomical variability introduced in Section 2.2.1, a linear regression analysis has been performed between the reciprocal of ΔQa and the length of the occluding artery, l. The same analysis has been done between the reciprocal of ΔPa and the length of the occluding artery, l.

Results
Systolic, diastolic and mean BP are shown in Figure 4. The values of systolic and diastolic BP get closer to the mean value as smaller vessels are considered, following a physiological trend.
AP and CO of the physiological and HC regimens are displayed in Figure 5 with the data used for validation. AP and CO in HC are higher than those in the physiological regimen. AP in both physiological and HC regimens is lower than pressure inputs because of the pressure drop in the left ventricle. CO at closed valve is maintained at 0 ml/s per cycle in both physiological and HC regimens; the peak it reaches in physiological regimen is 980.50 ml/s, while in HC it is higher and reaches 1049.70 ml/s. In physiological condition systolic and diastolic AP reach 118.90 mmHg and 76.30 mmHg respectively; in HC systolic and diastolic AP reach 137.90 mmHg and 87.70 mmHg, respectively.
ΔQa and ΔPa are reported for every occluding artery in Table 4. The average ΔQa after the ischemia is 0.38 ml/s per cycle, with maximum and minimum variations respectively of 0.04 ml/s and 1.93 ml/s. These values correspond to a 0.3% loss in Qa peak from before to after occlusion. For the Pa the results are similar, with a mean value of variation of 0.39 mmHg per cycle and minimum and maximum variations of 0.04 mmHg and 1.98 mmHg per cycle. When considering only arteries shorter than 10 cm, the average loss for Qa and Pa is respectively of 0.59 ml/s and 0.61 mmHg per cycle, this one corresponding to 0.5% of the simulated hypertensive peak. Results of the linear regression analysis of ΔQa and ΔPa with length of the occluding artery l are shown in Figure 6.

Discussion
This study proposed an integrated LP model of the cardiovascular system for the simulation of AIS and of hemodynamics instantaneous changes. For the sake of simplicity of the model, and for the complexity of cerebral vasculature anatomy, the model has been kept as simple as possible, so that anatomical variability has been simulated in a feasible manner.
In this study, control mechanisms have been neglected even though, physiologically, they act in a short and long term to modify pressure in front of body requirements. The main control mechanism playing a role in controlling pressure, in the order of seconds after modification, is the baroreflex control, for which the physiological latency varies between 0.5 s and 3-6 s [27,28]. Therefore, in physiological conditions the first peak after AIS for both Pa and Qa should be modified, especially in our model, where peak interval was of 1 s, as shown by [29]. However, baroreflex latency has been studied in some works applying different protocols and different species, that in principle limit the application of their results to the same testing conditions of the study; consequently, in recent years a clear general definition of the smallest baroreflex latency has been questioned [29]. It has been proven that a reliable value of baroreflex latency is in the order of 1 s or more, and that it varies with testing protocol, heart rate and other control mechanisms [30,31]. Moreover, it has also been proven that baroreflex activity changes with pathologies, for example with HC [32], and with sex [33]; in particular it has been shown that those factors mainly influence baroreflex sensitivity and its control gain, not reporting extensive details about baroreflex latency [33]. Other pathological conditions have also been shown to have a negative effect on baroreflex latency, increasing it up to more than 2 s [34]. Some studies have been conducted about the impact of pathologies on pressure control mechanisms, but still adequate control in HC is poorly defined [35]. Considering pressure control mechanisms is important for the validity of a model simulating physiological conditions to study their effect on pressure peaks [28]. Even if our model neglects pressure control mechanisms to study the first peak after AIS, it is supposed not to be a major limit to the validity of the model, since the simulation in pathological HC regimen considers slow control mechanisms as discussed above, that are not meant to affect to a great extent the first peak after AIS.
The small dimensions of the occluding artery during AIS lead to a small loss of Qa and Pa with respect to the HC. However, reliable values can be variations of 0.60 both ml/s or mmHg per cycle, or 0.5% of the normal physiological values. These results explain why it is so hard to detect an AIS simply by measuring drops of pressure or blood flow in order to use them to evaluate single AIS cases. Other techniques based on diagnostic imaging give better results and allow to evaluate the complexity of the occurring AIS in relation to the other brain vessels. From carotid arteries to deep brain vessels there is a big variety of dimensions. In this work a paradigmatic anatomy has been used considering different occluding arteries and a correlation has been highlighted between the length of the occluding artery and the inverse of arteriolar pressure drop; it can be seen in Figure 6 that shorter arteries lead to higher drops in arteriolar pressure, but since their magnitude is very low, the reliability of that linear regression has to be verified. In fact, a longer occluding artery leaves more cells without perfusion than a smaller one, and then makes more damages. That could be used to classify AIS in intensity levels, since one of the most spread ways to do it is by evaluating macroscopic symptoms. In the end, experimental observations are needed to validate results of Figure 6 and Table 4.
The choice of a LP rather than a CFD model is based on the fact that a CFD model needs a very specific geometry to solve Navier Stokes equations, so the results obtained are usually patient specific. Furthermore, CFD models hardly allow to do a macroscopic analysis of the whole elastic cardiovascular system, that in the proposed LP model is represented as the direct RLC. On one hand LP models need less parameters and are based on a simpler geometry, while on the other their simplifications are source of error. Moreover, for both typologies of model the definition of the optimal value of each parameter is challenging, because of the high number of parameters required. For this reason, an alternative approach based on the simulation of occluding arteries of variable length has been proposed in this study. The results obtained with the proposed model for the physiological regimen are coherent with physiological reality: the values of systolic, diastolic and mean BP follow well reality for the whole model from the left ventricle to the right atrium and the rest of the circulation as compared with data from literature [25]. Specifically, the trend of AP is similar to the trend proposed in literature [26], even if the characteristic dicrotic notch and the minor oscillations have not been replicated ( Figure 5). The trend of CO is slightly different in shape ( Figure 5) since it does not replicate the retrograde incision in flow being the motion of the valves modelled as instantaneous; moreover, the peak value of CO obtained with simulation is higher than the experimental value, probably due to flow velocity reduction from aortic valve to the ascending aorta, where the pressure has been measured [26]. However, the trend of CO is similar to the one obtained from other models [10]. In the HC regimen, a higher systolic BP has been imposed, and as a result a higher CO can be observed. The main assumptions for the simulation of HC is that heart and valves anatomy do not change with an higher pressure, since HC is most of the cases related to a reduction of the elastance of the arteries, increasing their pressure conductivity [3]. Parameter values assumed in the HC may be affected by error because the reshaping effect of HC has been simulated here based on considerations partially based on animal models, especially for the computation of the modified arterial internal radius.

Future perspectives
One of the main trend in cardiovascular modelling is the development of models as bed-side tools [36]. The reduced complexity of LP models with respect to CFD models make the former candidates for a bed-side tool because of the dramatically reduced computational costs, and then higher speed. Many LP models implement a high detailed part of the cardiovascular system, while the rest is modelled with low complexity [7,10,12,18,37,38]. A more complex model would lead to more accurate results, but would even increase computational cost. In literature, cardiovascular models that integrate a more or less detailed description of SC and PC with physiological and pathological control mechanisms are extensively studied [27]. For many applications a model without control mechanisms would be very limited, but for the evaluation intended to be made in this work control mechanisms can be considered neglectable. Moreover, such a bed-side model should also consider patient specificity, such as height, age and others. Model personalization is useful for both patient-specification and clinical training, allowing clinicians to understand the variation of some parameters on the entire system and reproducing pathological conditions [39]. This kind of model may be helpful in supporting the surgery field, where the clinician could design the surgery in advance and assess cardiovascular variations in real time, even during the surgery itself. This kind of model could be applied even in decision support in clinical practice in the field of AIS, as well as in training programs.

Conclusions
An integrated LP model of the cardiovascular system which simulated AIS and hemodynamics instantaneous changes has been proposed. The results about the magnitude of ΔQa and ΔPa highlight the difficulty of experimentally evaluate hemodynamics when an AIS occurs. The numerical values obtained are small with respect to systolic and diastolic BP and usually comparable with physiological pressure variations. However, a strong relation between the length of the occluding artery and the reciprocal values of ΔQa and ΔPa has been highlighted, even though its reliability requires further validation.