A 1D computer model of the arterial circulation in horses: An important resource for studying global interactions between heart and vessels under normal and pathological conditions

Arterial rupture in horses has been observed during exercise, after phenylephrine administration or during parturition (uterine artery). In human pathophysiological research, the use of computer models for studying arterial hemodynamics and understanding normal and abnormal characteristics of arterial pressure and flow waveforms is very common. The objective of this research was to develop a computer model of the equine arterial circulation, in order to study local intra-arterial pressures and flow dynamics in horses. Morphologically, large differences exist between human and equine aortic arch and arterial branching patterns. Development of the present model was based on post-mortem obtained anatomical data of the arterial tree (arterial lengths, diameters and branching angles); in vivo collected ultrasonographic flow profiles from the common carotid artery, external iliac artery, median artery and aorta; and invasively collected pressure curves from carotid artery and aorta. These data were used as input for a previously validated (in humans) 1D arterial network model. Data on terminal resistance and arterial compliance parameters were tuned to equine physiology. Given the large arterial diameters, Womersley theory was used to compute friction coefficients, and the input into the arterial system was provided via a scaled time-varying elastance model of the left heart. Outcomes showed plausible predictions of pressure and flow waveforms throughout the considered arterial tree. Simulated flow waveform morphology was in line with measured flow profiles. Consideration of gravity further improved model based predicted waveforms. Derived flow waveform patterns could be explained using wave power analysis. The model offers possibilities as a research tool to predict changes in flow profiles and local pressures as a result of strenuous exercise or altered arterial wall properties related to age, breed or gender.

Introduction A wide range of one dimensional (1D) computer models of the human arterial circulation is available. Such models allow the computation of pressure and flow waveforms throughout the whole arterial network, and hence allow researchers to study the normal and abnormal physiology of the cardiovascular system, without the need of in vivo measurements [1][2][3][4][5][6][7]. 1D models are well-balanced between complexity and computation costs, making them relevant for many (bio)medical applications. Due to their capability of involving extensive arterial segments, 1D models can provide useful information about characteristics of blood flow at the level of individual branches or even in patient-specific situations [8,9]. These models can also be used as a non-invasive diagnostic tool, helping physicians to understand observed changes in routine clinical blood pressure measurements and their possible physiological origin and to predict surgical operation results [9][10][11].
Due to technical limitations, difficult arterial accessibility, and ethical concerns, in-depth pathophysiological research of the equine arterial tree remains challenging and could be facilitated by the application of a model. Because of large differences between human and equine arteries regarding dimensions and branching patterns, especially of the aortic arch, a horse specific 1D model of the arterial circulation is needed.
The aim of this study was therefore to develop a 1D computer model of the equine arterial circulation. Providing reference data on equine arterial hemodynamics and physiology, this model might contribute to a better understanding of some clinical findings, such as the origin of the more oscillatory flow patterns in horses, the higher prevalence of aortic rupture in Friesians compared to Warmblood horses [12,13], the occurrence of sudden death during exercise due to arterial rupture [14][15][16][17], the higher chance on uterine artery rupture in older mares [18,19], or the higher chance on arterial rupture after phenylephrine administration in older horses [20]. In order to develop a reliable model, several anatomical data of the main equine arterial tree were collected ex vivo and combined both with in vivo invasive blood pressure measurements and non-invasively determined ultrasonographic flow profiles.

Materials and methods
All procedures were approved by the Ethical Committee of the Faculty of Veterinary Medicine, Ghent University (EC 2016/104). Five warmblood horses were investigated, mean age 18 ± 3 years and mean body weight 648 ± 47 kg. All horses were scheduled for euthanasia because of non-cardiovascular reasons. One horse was privately owned (informed owner consent was obtained), the remaining four horses were experimental horses owned by the Faculty of Veterinary Medicine, Ghent University. For euthanasia the following protocol was applied: premedication with detomidine 0.02mg/kg; induction with a combination of ketamine 2.2mg/kg I.V. and midazolam 0.04mg/kg I.V.; and finally euthanasia with 6 ml/kg T61 (Intervet International GmbH, Unterschleissheim, Germany), containing 24mg/kg embutramide, 6 mg/kg mebezoniumjodide and 0,6 mg/kg tetracaïnehydrochloride.

In vivo measurements
Ultrasound. Ultrasound imaging was performed (Vivid IQ, GE Healthcare) on all 5 standing, non-sedated horses. Different regions along the arterial tree were examined: the aorta from a left and right parasternal position, the right common carotid artery 15 cm cranial to the thoracic inlet, the right external iliac artery from the inguinal region and the right median artery just proximal to the carpus on the medial side of the leg. 2D B-mode images were collected, using a 9 MHz linear transducer (9L-RS, GE Healthcare) for the common carotid and the median artery, a 6 MHz phased array probe (6S-RS,GE Healthcare) for the external iliac artery and a 5 MHz phased array probe (M5Sc-RS, GE Healthcare) for the aorta. Mean values were obtained from 3 consecutive cardiac cycles at a heart rate between 35-45 beats per minute. Measurements were performed off-line (Echopac version 201, GE Healthcare). Diastolic diameters were measured from a transverse image for the carotid and external iliac artery, and from a longitudinal image for the aorta and the median artery. Pulsed wave Doppler images were collected at every location, using a 6 MHz phased array probe (6S-RS,GE Healthcare) for the common carotid and the external iliac artery, a 5 MHz phased array probe (M5Sc-RS GE Healthcare) for the aorta, and a 9 MHz linear probe (9L-RS,GE Healthcare) for the median artery. Angle correction was set at 45˚for every image at all locations. Using this fixed angle correction, images were optimised to align with the flow direction.
Invasive blood pressure. In all 5 horses, the blood pressure at the level of the common carotid artery was measured invasively in the standing, awake animal. The right common carotid artery was punctured aseptically under ultrasound guidance (Vivid IQ, GE Healthcare; 9L-RS, GE Healthcare), using an 18 gauche 90mm needle (Terumo spinal needle, Terumo) placed in the middle of the lumen and kept in place for at least 20 consecutive heart cycles. The needle was connected with a fluid filled pressure transducer (MLT0699 Disposable BP Transducer, ADInstruments) interfacing with a digital acquisition station (PowerLab 8, ADInstruments), blood pressure curves were recorded for offline analysis (LabChart, ADInstruments). For each horse the systolic, diastolic and mean arterial pressure was calculated automatically as the mean of 20 consecutive heart beats (heart rate between 35-45 bpm).
In one horse, scheduled for euthanasia, invasive blood pressure measurements over the whole length of the thoracic and abdominal aorta were performed under general anaesthesia (pre-medication: detomidine 0.02mg/kg; induction: combination of ketamine 2.2mg/kg I.V. and midazolam 0.04mg/kg I.V.) with the horse in the dorsal recumbent position. After chirurgical exposure of the right common carotid artery, a 72cm steerable 8.5Fr sheath (Zurpaz, Boston Scientific) was introduced using the Seldinger technique. Under transthoracic ultrasound guidance (Vivid IQ, GE Healthcare; M5Sc-RS, GE Healthcare) the sheath was introduced retrogradely through the brachiocephalic trunk, into the ascending part of the aorta. Once in place, a custom-made pressure tip catheter (Gaeltec) was introduced and advanced to the most caudal end of the aorta under transrectal ultrasound guidance. Blood pressures were recorded (PowerLab 8, ADInstruments) at the most caudal site and subsequently at every 10 cm during step-wise pulling back of the catheter until the ascending part of the aorta was reached. An ECG was recorded simultaneously. At each location, an ensemble-averaged waveform was constructed from at least 5 cardiac cycles. The ensemble-average was aligned in time relative to the peak of the R-wave of the ECG. Pulse wave velocity was calculated from the relation between inter-measurement distance and time delay. The time delay was calculated from delays in the peak of the 2 nd derivative of pressure, assumed to represent the foot of the pressure wave. After the procedure the horse was euthanized while still under general anaesthesia with 6 ml/kg bodyweight T61 (Intervet International GmbH, Unterschleissheim, Germany) containing 24mg/kg bodyweight embutramide, 6 mg/kg bodyweight mebezoniumjodide and 0,6 mg/kg bodyweigth tetracaïnehydrochloride.

Ex vivo measurements
Necropsy of 4 out of 5 horses was performed within 12 hours after euthanasia. A dissection was completed on the aorta and the most important (left sided) side branches and morphometric data (length, diameter and branching angle) were recorded. Arterial length was measured using a tape measure and diameters were measured in the middle of each segment by introducing custom-made iron rods of different diameters into the explored arterial lumen.
Subsequently, post-mortem diameters were scaled to the in vivo diameters, using the in vivo ultrasound measurements of common carotid artery, external iliac artery and median artery as a reference. Angles of the different arterial segments in the 1-dimensional plane were measured using anatomical images of the equine arterial tree [21]. Lengths of terminal segments were only measured in one horse. Right sided circulation was assumed to be the same as the left sided arterial circulation, except for the right sided subclavian circulation (as the branching pattern is different), which was measured separately in one of the horses. All investigated arterial branches and their corresponding average length and diameter are displayed in Table 1. Note that in the mathematical model vessel tapering will be included, meaning that average segment diameters will be adjusted to minimize forward reflections (see chapter 'Distal boundary conditions and bifurcations').

Mathematical model
The mathematical model is based on a previously published and validated 1D arterial network model in humans (22,23). Main differences between the human and the adapted equine model are described below. The system of nonlinear equations is solved using in-house MatLab code. For the solution an implicit finite difference scheme was chosen, with second order of accuracy for the temporal and spatial domains. Forward, central and backward difference approximations of the spatial derivative were used for the proximal, middle and distal nodes in each segment, respectively. The arterial tree was initialized with a pressure of 100 mmHg and a flow of 1ml/s. The solution was found over 8 cardiac cycles yielding pressure and flow waveforms over the entire arterial tree. For more details on the modelling aspects and the mathematical equations, we refer the interested reader to Ref. (22).
Governing equations. The main branches of the equine arterial tree were divided into 117 interconnected straight cylindrical arterial segments and inserted in the mathematical model (Fig 1, Table 1). The integrated forms of the continuity and momentum equations of the Navier-Stokes equations were solved in each of these segments for pressure (P), flow (Q) and cross-sectional area (A), where x and t are the spatial and temporal variables, u is the longitudinal velocity component, and R is the lumen radius. Given that measurements in animals, used as input into the model, were acquired in standing, awake animals, as well as in anesthetized and supine animals, we expect effects of gravity to be important given the height of the animal. We therefore accounted for the effects of gravity including the body forces term in the momentum equation, with g = 9.81 m/s 2 , the gravitational acceleration constant, and θ, the projection angle on the vertical axis. Results will be further reported for both situations, considering the effects of gravity and neglecting gravitational body forces. Blood was assumed to be an incompressible Newtonian fluid with density ρ = 1050 kg/m 3 and dynamic viscosity μ = 0.004 Pa�s. Given the large arterial diameters in horses, the Witzig-Womersley correction (see section 'Velocity profile') was used to approximate the convective acceleration term ( @ @x R A u 2 dA À � ) and the wall friction , both present in the momentum equation. Modelling of the arterial wall. A constitutive equation is needed to account for the elastic properties of the arterial wall, relating the dependency in intra-arterial pressure with the cross- sectional area. The nonlinear elastic behaviour of the arterial wall is assumed with an expression of area compliance (C A ) [22], as the product of a location-dependent function, C d ð � d; P ref Þ, and a pressure-dependent function, C p (P).
The function C d ð � d; P ref Þ gives the compliance for a given local mean arterial lumen diameter, � d, at a given reference pressure value, P ref = 100 mmHg. Reymond et al. [22] fitted an empirical inverse power curve for PWV as a function of � d, from human data reported in the literature.
To the best of our knowledge, measurements of the pressure and diameter dependency of the compliance in horses are lacking in literature but it can be assumed that the intrinsic building blocks (elastin, collagen, smooth muscle cells, proteoglycans) and their organization, and hence the overall mechanical behaviour, is similar as in humans. Initially, following Reymond et al. [22], we used the fitted values on humans for the pressure dependency of the compliance. However, these parameters were subsequently adapted to obtain values of pulse pressure close to the value reported by Boegli et al. [23] for healthy horses, resulting in a 1 = 0.76, b 1 = 5, P maxC = 10 mmHg and P width = 21 mmHg. We also applied the fit obtained from human data to describe the relation between local diameter and compliance (for a 2 = 13.3 and b 2 = 0.3). Nonetheless, we did tune the arterial distensibilities to the equine physiology, by multiplying the distensibilities of all vessels of the 1D model by a common factor with value 0.75. That factor was determined such that in resting conditions, we obtained a value of PWV close to the PWV computed from invasive blood pressure measurements over the abdominal aorta of an anesthetized horse; keeping in mind that differences may arise between an anesthetized and a conscious animal. Velocity profile. The Witzig-Womersley theory describes the effect of flow pulsatility and inertia on the velocity profile. This oscillatory flow theory was needed to calculate the terms of convective acceleration and wall shear stress in the momentum equation derived from the Navier-Stokes equations, given that both terms depend on the instantaneous velocity profile. The Witzig-Womersley theory requires the knowledge of the local flow profile across the arterial lumen over the entire heart cycle, which is a priori unknown in the 1D model. This difficulty was overcome assuming that the solution is periodic; the flow waveform from the previous heart cycle was used to calculate the velocity profile and the wall shear stress using the formulations:  Table 1. https://doi.org/10.1371/journal.pone.0221425.g001 The velocity profile (u) and the wall friction term (τ) are calculated as a Fourier series with harmonics (n) and depend on the harmonic-specific Womersley's number a ¼ R ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi r2pf=m p , where R is the artery radius and f the frequency, r/R is the relative radial position, Q n is the n th harmonic of the flow pulse, and J 0 and J 1 are the Bessel functions of first kind of order 0 and 1, respectively. The oscillatory flow theory is taken into account for vessels with α>3 (75% of the total number of arterial segments), where the wall friction and convective acceleration are related as independent terms in the 1D momentum equation; the solution described by Stergiopulos et al. [3] is considered otherwise. The solution is found solving in time a number of repeating cycles until convergence.
Distal boundary conditions and bifurcations. For the terminal nodes a 3-element Windkessel model was used, to account for the cumulative effect of all distal vessels beyond the terminal sites [3,22]. The equation takes the form: where R 1 is the proximal resistance, R 2 is the distal resistance, and C T is the terminal compliance. Total peripheral resistances R T = R 1 +R 2 , were estimated taking into account both the distribution of flow described in the literature [24] and our own measurements of flow at specific locations from ultrasound data. A total resistance parallel combination of about 0.14 mmHg � s/ml was assumed. The values of R 1 were estimated assuming minimal reflection at high frequencies, with the condition R 1 = Z c , where Z c = ρ�PWV/A is the characteristic impedance of the terminal segments. The values of distal resistance were then calculated as R 2 = R T −R 1 .
Terminal compliance was estimated following Reymond et al. [22], where the terminal compliance of each terminal vessel, C T,i , was assumed to be proportional to the area compliance, C A,i , at the distal end of the terminal vessels: with C T = ∑C T,i the part of the total volume compliance attributed to peripheral vessels not included in the arterial tree model, assumed to be in the order of 20% of the total systemic vascular compliance. The total systemic vascular compliance is the sum of the volume compliance of all vessels and compliance of the terminal beds, so that C V ¼ P i n C V;i þ P i m C T;i , where n = 117 is the total number of arterial segments and m = 62 is the number of terminal segments. To obtain volume compliance of each segment, the area compliance given by Eq (3) is integrated over the segment length. All values of terminal resistance and compliance can be found in S1 Table. Continuity of pressure and flow was imposed throughout the arterial network at bifurcations. Forward wave reflections were minimized by adapting the characteristic impedance of the tributaries so that the absolute value of the reflection coefficient was < 0.2 at all bifurcations. The cross-sectional area of the vessels was therefore slightly adjusted resulting in tapered-structure segments. Cross-sectional areas were determined by minimizing the reflection coefficient, subject to three conditions: (i) measured area is the average between the input and output areas (A in and A out respectively); (ii) A in of tributaries is lower or equal to A out of the parent, and (iii) for each segment A out � A in . Initial searching points were obtained by assuming linear tapering. The forward wave reflection coefficient is calculated as: where Z is the characteristic impedance of the parent and daughter vessels.
Heart model. A model of the left ventricle (LV) was used at the proximal end of the arterial tree (root of the ascending aorta), simulating the blood flow pumped out of the LV. The LV model is based on the time-varying elastance model, which describes the variation of LV pressure (P LV ) and volume (V LV ) during a cardiac cycle.
The interaction between the ventricle and the arterial tree is mainly produced during the ejection phase, when the aortic valve is open. To improve the simulation of the wave reflection phenomena that occur during this period between the ventricle and the aorta, an internal resistance of the LV was introduced. Taking this into consideration, the varying elastance model originally suggested by Sagawa [25] was adapted by Reymond et al. [22], leading to the following expression for the varying elastance of an ejecting heart: where E � is the elastance that would be measured during a nonejecting isovolumic contraction, and κ a constant relating the internal resistance of the LV to the ventricular pressure during the same cardiac phase. Assuming that the elastance curve (when normalized with respect to its peak value) is similar in shape for all mammals, a normalized isovolumic elastance, E � , can be derived from Eq (11) using the global normalized elastance curves, E, reported by Senzaki et al. [26]. The constant κ was derived iteratively by minimizing the difference between the elastance curve resultant from the 1D model and the original elastance curve reported by Senzaki et al. [26]. The value of κ obtained was 55E-06 s/ml. Due to the lack of detailed horse data in current literature, assumptions needed to be made to set most of the input parameters necessary for the heart model. The value of end-diastolic pressure was taken as, P end-dias = 16 mmHg, according to the value reported by Brown and Holmes [27] for a normal horse and well within the standard range (12-24 mmHg) reported in the literature for horses [28], whereas the end-systolic pressure was 113 mmHg [27]. The value of the dead volume of the LV was set to V 0 = 0 ml. Initial reference values for stroke volume (SV = 900 ml/min) and ejection fraction (EF = 60%) were used to estimate the end-diastolic volume (EDV) and the endsystolic volume (ESV), wherewith initial values of minimal and maximal elastance were derived considering Eq (10). These values were further tuned to obtain a close match between the simulated aortic flow velocity and the aortic flow velocity measured from ultrasound imaging, which resulted in final values of E min = 0.01 mmHg/ml and E max = 0.26 mmHg/ml. The standard heart rate was set to 40 bpm, a normal physiological value for the horse at rest, whereas systolic duration was set to 478 ms [29].The heart model simulates the four main phases of the cardiac cycle, starting the loop at the onset of the isovolumic contraction phase, where the volume in the LV equals EDV (derived in the simulation from Eq (10), for a LV pressure equal P end-dias ). With the contraction, the pressure in the ventricle rises over the aortic pressure, which causes the opening of the aortic valve and the start of the ejection phase. During this period, the ventricle-arterial interaction is described by the combination of Eqs (10) and (11). When the flow becomes negative, end systole is reached (aortic valve closes) and the relaxation phase takes place. The filling phase is set when the pressure in the ventricle drops below the initially assumed P end-dias , and the filling flow is modelled from the internal resistance of the LV (assumed as 0.003 mmHg�s�ml -1 ). The solution of the heart model is periodic; at the start of every cardiac cycle, EDV takes the value derived from the previous cardiac cycle.
Coronary model. Coronary arteries were modelled following Reymond et al. [22], assuming that changes in compliance, distensibility and resistance are proportional to the local time varying elastance of each vessel. For the right coronary, it was additionally assumed that the effect of the right ventricular contraction is smaller by a factor proportional to the ratio of maximal pressure in the two ventricles (PLV,max/PRV,max � 3 [28]).

Hemodynamic and wave reflection analysis using wave power
Wave power analysis (WPA) [30] was applied in different locations to study the dynamics of the waves. The method defines the wave power, dπ, as the product of the changes in pressure (dP) and flow (dQ), and equals the energy carried by the wave which is conserved at junctions. Wave power can be separated into its forward and backward components (subscripts (+) and (-), respectively), dπ = dπ + + dπ − , where dπ + = 1/(4�Z c )(dP+Z c dQ) 2 and dπ − = −1/(4�Z c )(dP −Z c dQ) 2 . The concept of wave power is analogous to wave intensity analysis [31], with dπ>0 indicating dominant forward waves and dπ<0 indicating dominant backward waves. The distance travelled by the waves to their reflection points can be estimated as the product of the transit time between a wave and its reflection and the local theoretical PWV. Since the wave travels twice the same distance, the final distance can be computed as L = (Δt�PWV theor )/2.

General physiological parameters
Running the model without including gravity revealed a cardiac output of 33 L/min, an ejection fraction of 65% and a stroke volume 820 ml. Systolic/diastolic pressure in the aortic root was 114/70 mmHg, with a pulse pressure of 44 mmHg and a mean arterial pressure of 93 mmHg. Taking gravity into account, cardiac output was reduced to 30 L/min, with an ejection fraction of 59% and a stroke volume of 740 ml, whereas systolic/diastolic pressure increased to 131/88 mmHg, resulting in an almost unaltered pulse pressure of 43 mmHg and an increased mean arterial pressure of 111 mmHg. The distribution of cardiac output derived from the model for both configurations is summarized in Table 2. Table 3 shows the values of the Womersley number, maximum shear stress, mean values of the convective acceleration approximation and the Reynolds number, derived from the model at different locations in the arterial network, including and neglecting gravity. for both configurations of the model. On the other hand, the total vascular resistance resulted in a value of 0.17 mmHg � s/ml for the model without gravity, and 0.22 mmHg � s/ml for the model with gravity. Note that these values differ from the total combination of resistances in parallel, since the total vascular resistance also account for the resistances in series in the vessels.

Effects of gravity
To assess the importance of considering gravity in the model, the distribution of pressure and flow velocity was plotted all over the arterial tree, both with and without taking gravity into account (Fig 2). Considerable differences in pressure were most evident in the limb arteries, carotid arteries and the arteries of the head. Because of pressure amplification, systolic pressures were relatively high in the front and hind legs (>160 mmHg) when gravity was neglected while mean pressures were almost unaltered over the whole arterial tree. Considering gravity, systolic pressures in the front legs became even higher, whereas systolic pressures in the arteries of the head became low (<100 mmHg) and mean pressures increased from the head to the legs. Mean flow velocity distribution was similar for both configurations, with higher values for vessels in the splanchnic region and toward the head. A more direct comparison on the features of waveforms is displayed in Fig 3, with evident discrepancies in flow waveforms of the common carotid artery and in pressure waveforms for the proximal aorta, common carotid artery and median artery. Table 4 presents peak blood flow velocities and pressure data estimated from the 1D model in both situations, including and neglecting gravity, along with available in vivo measurements. Of note is the higher pulse pressure predicted by both models, in combination with the amplification of pulse pressure towards the periphery. Moreover, both models reveal an amplification in pressure along the aorta. An increase in systolic pressure from the ascending aorta to the iliac artery of~54% was found for the model neglecting gravity, while it was reduced tõ 34% when gravitational forces are included. Fig 4 displays the modelled blood flow velocity waveforms, both with and without including gravitational forces, for the ascending aorta, common carotid artery and the main limb arteries (median artery and external iliac artery), compared with the measured waveforms in standing horses, derived using pulsed wave Doppler ultrasound. A relatively good similarity in waveform shape and their amplitude was found at all arterial locations, with relative errors for peak flow velocity between in vivo data, in the standing awake animal, and simulations without including gravitational forces of 29% for the ascending aorta, 82% for the common carotid artery, 15% for the external iliac artery and 40% for the median artery. When gravitational forces are included, relative errors improved to 2%, 48%, 2% and 23%, respectively.  horse), with and without considering gravitational forces. The effect of wave propagation is well captured by the model. Modelled mean, systolic and diastolic arterial pressures in the proximal aorta showed relative errors of 11%, 3% and 25%, respectively, compared to the invasively recorded pressure measurements (recorded under general anaesthesia with the horse in dorsal recumbent position). When gravitational forces are included relative errors changed to 6% for mean pressure, 12% for systolic pressure and 5% for diastolic pressure. For the distal aorta relative errors were slightly higher, 21% for mean pressure, 20% for systolic pressure and 40% for diastolic pressure. Including gravity had almost no influence on the derived pressures and thus the relative errors. Aortic PWV determined from the simulation was 5.3 m/s when gravity is neglected, changing to 5.2 m/s considering gravity, whereas the value obtained from in vivo data was 5.27 m/s (difference of~0.6% when comparing with the model neglecting gravity, and~-1.3% when comparing with the model including gravity).

Wave power analysis (WPA)
WPA is shown in Fig 6 for the proximal aorta, common carotid, median and external iliac artery. Models with and without including gravitational forces are compared. In both situations, the proximal aorta shows the three typical peaks: i) a forward compression (FC) wave generated by the systolic ejection, ii) a backward compression (BC) wave as the result of peripheral reflection, and iii) a forward expansion (FE) wave due to the slowing of ventricular contraction. Besides these normal peaks, a small FE can be visible at the different locations in mid-systole. The timing of occurrence of this wave coincides with the drop in flow signals for the median and iliac arteries. For all the considered locations, wave power was lower when gravitational forces were included. The largest difference between both situations is observed at the carotid artery, where for the configuration without gravity, a mid-systolic FE wave combines with the typical FE from the ventricular contraction, resulting in a wide FE wave; and this in a more drastic drop in the carotid systolic pressure of the model without gravity than the model with gravity (Fig 3).

Discussion
Current literature provides only limited information on arterial hemodynamics in horses. Indepth fundamental research into equine intra-arterial pressures and flows is therefore expected to increase the understanding of aortic and arterial rupture in this species [14][15][16][17]. Today, due to technical limitations it is difficult to assess pressure and flow in the centrally located arteries and therefore peripheral measurements are usually extrapolated to the rest of the vascular network. The present study aimed to develop a 1D computer model for the equine arterial circulation comprising all major vessels of the arterial tree. The original 1D human model on which this equine model is based, has been validated and is used as a research tool in many studies [5,[32][33][34][35]. Recently, it formed the basis for the development of a 1D model in mice [36]. As aortic branching pattern of horses is completely different compared to humans, arterial segments needed to be redefined. Ex vivo measurements (diameter, length and branching angle) of each  (Table 1). Dimensions (lengths and diameters) of the equine arterial tree are much bigger than those of humans. Moreover, equine heart rate is lower at rest (28-45 bpm in horses and 60-100 bpm in humans) and the heart mass is much larger compared to humans (>1% of the body mass in horses and~0.5% of the body mass in humans [24]). Both elements imply a much higher cardiac output in horses (~35 L/min compared with~7 L/min in humans). Of note is the higher Womersley number in the horse, implying a considerably higher impact of inertia on the flow velocity profile, leading to flatter velocity profiles with steeper velocity gradients near the wall than those predicted by Poiseuille flow. As this leads to an improved estimation of viscous friction, it was important to include the Witzig-Womersley correction factor in the momentum equation.

Importance of gravity
Existing 1D models usually ignore the effects of gravity, because most physiological measurements in humans are performed in supine position, with the heart being at the same level as the rest of the arterial tree. In horses most blood flow velocity and pressure measurements are performed on standing patients, and therefore the effect of gravity was incorporated in the present model. The gradient of pressure in the arterio-venous system, which causes flow, does not fluctuate considerably for the body in the supine position compared with that in the standing position; thus overall flow is similar in both situations. The transmural pressures on the contrary, are strongly different. Given the nonlinear pressure dependency of arterial stiffness, gravity was expected to exert an effect on wave propagation and wave morphology. This was indeed confirmed by our results. In order to assess the impact of gravitation, two configurations of the model (with and without gravity) were compared (Figs 2-6). Taking into account gravity generally improved predicted peak flows. Relative errors improved from 15-82% for the model without gravity to 2-48% for the model including gravity. While peak flow predictions improved when including gravity into the model, relative error for mean pressure, diastolic pressure and systolic pressure did not improve when gravity was included, with relative errors ranging from 3-40% without gravity and 5-56% with gravity included. Fig 2 shows that, when gravity is neglected, mean blood pressure is practically the same in all conduit arteries. As was expected, including gravity into the model revealed more variability in mean pressure throughout the arterial tree, with lower values toward the head and higher values toward the limbs. Fig 5 clearly shows the similarities in pressure waveforms of the descending aorta with and without gravity. These similarities in pressure waveforms are due to the assumption of a horizontal motion of blood (angle = 0) in the descending aorta, causing no influence of gravity.

Wave power analysis
As changes in arterial pressure waves are associated with alterations in the contour of the arterial flow profiles [37], the pronounced oscillations during flow wave measurements at different arteries in horses probably indicate pressure waves returning from the periphery at multiple reflection sites, starting during systole and continuing during diastole. Wave power analysis performed at several locations revealed important wave reflections mainly during systole. Wave power patterns were most complex for the common carotid artery, which also displayed the biggest delay between the first FC wave and its peripheral reflection (BC wave) coming from the head. This was somewhat expected considering that the wave has to travel forth and back along the carotid artery. The reflection of the first forward peak at the peripheral site in the head was higher when gravity was included in the model (~56% with gravity vs.~36% without gravity). The presence of a mid-systolic forward expansion wave immediately after the first peak generated by the heart ejection (Fig 6), was responsible for the drop of pressure and flow velocity that resulted in narrow systolic peaks in pressure and flow waveforms of the investigated vessels. In order to understand the origin of this suction wave, we estimated the distance travelled by the wave to the site of re-reflection. By combining the time difference measured from the foot of the waves with the local theoretical PWV in the involved segments, very short distances were derived. These re-reflections might be occurring locally in the network rather than at the heart level, and may be related to a mismatch of junctions in the peripheral branches. Higher efforts to minimize forward wave reflections at these sites are needed.

Limitations and future work
In general flow waveform patterns (morphology) are well captured by the model, especially by the one including gravity while discrepancies in amplitude are quite obvious (Fig 4), this is probably due to several limitations of this model, which will be explained clearly in this section.
Since the present model is based on averaged data, it enables us to predict generic local pressures and flow profiles in all investigated arterial segments. A fully quantitative validation, however, would require a detailed horse-specific approach, tuning of all input parameters that define the model to each specific animal (such as geometry, elastic properties, peripheral resistance, and cardiac parameters) and comparing the outcomes of the model with in vivo measurements in that specific animal. Such an approach is technically challenging and almost impossible. Therefore, even if we consider the present model to be representative for the average healthy horse, waveforms should only be compared qualitatively and not quantitatively to individual measurements.
The lack of literature data on equine hemodynamics was the major challenge to develop an equine 1-D model. Tuning the model parameters was based on plausible assumptions and scaling factors between human and equine patients. Further fine-tuning of the input parameters to equine physiology will be necessary to obtain a closer match with in vivo flow profiles. Branching pattern and dimensions of the arterial tree, both largely defining flow wave patterns, are well integrated in this model. However, branching patterns and arterial dimensions can vary significantly between individual horses. Other important parameters, influencing flow velocity and pressure waves morphology are arterial elasticity and peripheral resistance. Both, arterial elasticity and peripheral resistance largely define diastolic flow [38] and pulse pressure [33], two parameters still showing large differences between measured values and modelled ones, with relative errors for diastolic pressure ranging from 5-60% and for pulse pressure ranging from 81 to 288%. Equine arterial elasticity and peripheral vascular resistance are therefore interesting criteria to further unravel in the future. In addition, in our simulations the compliance distribution and peripheral resistances remain the same for both, the model with gravity and the model without gravity. Neglecting the autoregulation mechanisms that lead to a cardiovascular response to control blood flow and pressure levels during postural changes, can also partially explain the differences found in our simulations and in vivo data.
The limited research possibilities in equines form another restrictive factor in the development of the equine 1D model. Due to ethical concerns, invasive aortic blood pressure was only collected in one horse, which is a limitation of this study. Moreover these aortic blood pressures were collected in anaesthetized, dorsally recumbent animals, while ultrasound and pressure measurements at the carotid artery were performed in non-sedated, standing horses. This complicates comparison between in vivo measurements and modelled findings. The present model does not yet account for the changes in physiological parameters due to anaesthesia or dorsal recumbent position. Anaesthesia tends to slow down heart rate, reduce both cardiac output and pulse pressure and likely leads to modulations in resistance of vascular beds and mechanical properties of arteries. Changes in smooth muscular tone modulate distensibility and stiffness. The recumbent position of the horse has a huge effect on intrathoracic and arterial transmural pressure, modulating stiffness and leading to important volume shifts of blood affecting cardiac filling and preload and, via the Frank-Starling mechanism, cardiac contractility. In addition, the model does not include physiological control mechanisms such as the baroreflex that neurologically modulates cardiovascular function via sympathetic/parasympathetic mechanisms. When a body changes from the supine to the standing position, these baroreceptor control mechanisms are activated to limit mean arterial pressure decrease to only a few mmHg. The restoring mechanisms include an increase in cardiac output (mainly an increase in heart rate), as well as an increase in systemic vascular resistance.
Next to this, discrepancies in pressure and flow waveforms between measured and obtained modelled values may be affected by angle correction. In order to obtain standardised images, angle correction was set at 45˚when collecting flow velocity profiles in the standing horse. Measurements were optimized for alignment with the flow but this alignment was not always perfect and flow was probably not always captured in the centre of the artery, which implies that captured flow velocities are only an approximation of the true flow profiles.
Discrepancies in pressure and flow waveforms may also be affected by the exclusion of the entire systemic and pulmonary circulation, making the model an open loop system that requires boundary conditions: a cardiac time-varying elastance model at the inlet, and the Windkessel model at the terminal ends. Moreover, the cerebral arterial tree is not included in detail in the equine model, only a simplified representation containing the major vessels that supply the cerebral circulation. A more detailed description of the cerebral arterial tree may provide better predictions of pressure and flow waveforms in the carotid artery and smaller vessels of the head circulation, as has been previously shown for the human and murine models [39]. Furthermore, in the present design, the model used for coronary arteries is simplistic. Last but not least, development of the present model was based on data obtained on a limited number of horses, without accounting for cardiovascular variations due to age or gender.
Nevertheless, despite individual differences in absolute values of flow velocities and arterial pressures, flow wave morphology is well captured by the model as shown in Fig 4. Fig 5 evidences the effect of wave propagation along the aorta both in the measured and the modelled pressure profiles. This indicates the added value of the model for studying trends in arterial flow dynamics in horse populations.
In the future this computer model may be useful to predict changes in flow profiles and local pressures under specific circumstances or conditions (age, exercise). During high-intensive exercise, heart rate may rise up to 8 times above the resting rate and total aerobic capacity can reach a 40-fold increase, which is much greater compared to human athletes. By altering input parameters of the horse-specific model, this model might predict local pressures and flow profiles during these extreme circumstances and contribute to the understanding of the relatively high incidence of sudden death during exercise due to arterial rupture [14][15][16][17]. As increasing age increases the risk of arterial disorders (arterial rupture during exercise, during parturition or after phenylephrine administration) [14,[18][19][20], it might also be interesting to use the present model to study the effect of age on arterial hemodynamics. Moreover the development and adjustment of this kind of computer models, could lead to a better understanding of some intensively studied, but poorly understood clinical situations such as exercise induced pulmonary hemorrhage.

Conclusions
A 1D computer model for the equine arterial circulation has been developed and this provided a unique anatomical dataset for horses. Ex vivo anatomical measurements were combined both with literature data and physiological information from ultrasound analysis in order to predict pressure and flow waveforms in the equine arterial tree by means of 1D modelling. The qualitative validation of the model was carried out by comparing the results with average flow velocities and pressures measured in vivo in horses. Despite its generic character and limitations, outcomes from the model showed plausible predictions of pressure and flow waveforms throughout the considered arterial tree. Simulated flow waveforms reproduce important features observed in ultrasound Doppler images, especially the oscillating pattern (most pronounced at the external iliac artery, median artery and common carotid artery). Adapting the model by taking into account gravity further improved predicted waveforms. Thanks to wave power analysis, the contours of the arterial flow profiles could be explained. Despite the shortcomings of in vivo measured pressures (aortic pressures measured under general anaesthesia with the horse in dorsal recumbent position), modelled pressure data seem in line with invasive measurements. We believe that the present model may be useful, not only to explain flow wave patterns in horses, but also to predict changes in flow profiles and local pressures as a result of strenuous exercise or altered arterial wall properties related to age, breed or gender.
Supporting information S1