# Hemodynamic Features in Stenosed Coronary Arteries: CFD Analysis Based on Histological Images

Authors: Mahsa Dabagh; Wakako Takabe; Payman Jalali; Stephen White; Hanjoong Jo
Journal: Journal of Applied Mathematics
 (2013)

---

## Abstract

Histological images from the longitudinal section of four diseased coronary arteries were used, for the first time, to study the pulsatile blood flow distribution within the lumen of the arteries by means of computational fluid dynamics (CFD). Results indicate a strong dependence of the hemodynamics on the morphology of atherosclerotic lesion. Distinctive flow patterns appear in different stenosed regions corresponding to the specific geometry of any artery. Results show that the stenosis affects the wall shear stress (WSS) locally along the diseased arterial wall as well as other adjacent walls. The maximum magnitude of WSS is observed in the throat of stenosis. Moreover, high oscillatory shear index (OSI) is observed along the stenosed wall and the high curvature regions. The present study is capable of providing information on the shear environment in the longitudinal section of the diseased coronary arteries, based on the models created from histological images. The computational method may be used as an effective way to predict plaque forming regions in healthy arterial walls.


---

## Body

## 1. Introduction

Atherosclerosis in the coronary arteries tends to be localized in regions with curvature and branching associated with complex, unsteady, and turbulent pattern of the flow [1]. In these regions, the fluid shear stress deviates from its normal spatial and temporal distribution patterns in straight vessels. The role of WSS in the localization of the atherosclerosis has been widely accepted [1–6]. Shahcheraghi et al. [3] explained clinically the fact that a low shear stress region is generally found in the downstream of developing stenotic plaques, leading to an increased density of leaky junctions associated with an increased level of low-density-lipoprotein (LDL) flux into the arterial wall. Moreover, the OSI is known as the predictor of formation of atherosclerosis and vulnerability for plaque in coronary arteries [7–9]. High OSI has been widely used as indicators of atherosclerosis because the regions with high OSI are likely to experience stagnation or backflow [7–11]. This is strongly related to the progression of atherosclerosis because it affects the arrangement of endothelial cells in adjacent tissues. It has been shown that the areas with high values of OSI are usually located in the regions where wall shear stress is low [7]. Taking all together, the creation and development of atherosclerosis depend not only on the biological features but also on the biomechanical factors, including WSS and the cyclic force caused by the pulsatile blood pressure. Therefore, a full description of the flow pattern through the diseased coronaries, with the emphasis on the patient-specific data, is essential for quantifying the wall shear stress that may occur in a stenosed artery during the cardiac cycle. The knowledge will also lead to find the most vulnerable sites for the plaque rupture.Previous studies have shown the flow pattern through the idealized geometry [6–8] or the reconstructed geometry of arteries from a series of slices obtained using computerized tomography (CT) scan images [9, 10], magnetic resonance images (MRI) [11], the spiral CT scan [12], the phase-contrast MRI method [13], or intravascular ultrasound (IVUS) [4]. However, simplifications were applied in these studies; for example, the cross-section of the artery was approximated as circles which fit the centerline to a cubic spline curve [7, 13], or a constant diameter was applied throughout the model [13, 14].In the present study, the histological images from longitudinal section of four diseased coronary arteries were used, for the first time, to investigate the local hemodynamic more precisely through the stenosed artery models. The geometry is reconstructed directly by graphical reading of points on the luminal surface of the arteries. Any simplifications have been avoided in the construction of the model. The pulsatile blood flow pattern was applied through the reconstructed stenosed coronary arteries, and the blood flow characteristics along the arterial wall was analyzed by means of computational fluid dynamics (CFD).


## 2. Methods

### 2.1. Reconstruction of 2D Geometry of Coronary Artery from Histological Image

Two-dimensional (2D) geometries of the human coronaries are reconstructed from the histological images which were acquired from coronary arteries of four donor hearts (Faculty of Medicine and Density, University of Bristol, Bristol, UK). The samples contained regions with mild to severe atherosclerosis. The diseased coronary arteries were cut in longitudinal section, stained with H&E, and imaged on a nanozoomer, namely, NanoZoomer 2.0HT (Hamamatsu photonics, Hamamatsu, Japan) with 40x magnification at the pathology core lab in the Winship Cancer Institute, Emory University.To construct the coronary artery geometries, all the images were processed individually. The images were converted first to the boundary plots within MATLAB v. 7.11.0 (R2010b). The lumen area was identified via image segmentation. Points representing the surface of the endothelium were obtained as the final output of this step. The vertices acquired from the processing of the images were connected to generate a 2D surface representing the lumen of coronary. The 2D model was trimmed with the CFD preprocessing and meshing code Gambit V2.4.6 to generate the flat inlet and outlet boundaries. The computational mesh grid consisted of about 50,000 grid elements. Thin boundary layers were used for the precise calculations near the walls. Several steady-state simulations were performed with different grid sizes to ensure the mesh independence of the results. The boundary layer mesh used in these simulations consists of 6 layers with the first layer attached to the luminal surface with the thickness of 0.01 mm. The thickness of consecutive layers grows with the factor of 1.2. Figures1(a)–1(d) demonstrate the original histological image of each coronary artery. Figure 1(e) shows the grids generated for the models.(a)–(d) LS images of four diseased coronary arteries, (e) the computational meshes at the stenosis location, and (f) velocity pulse at the inlet.
(a)
(b)
(c)
(d)
(e)
(f)


### 2.2. Governing Equations and Boundary Conditions

The blood flow is assumed to beincompressible with uniform properties including the viscosity of3.5×10-3 Pa·s (Newtonian fluid) and the density of 1067 kgm−3. The assumption of Newtonian fluid for the blood is valid wherever the shear rate is greater than 100 S−1 [9, 10, 12, 13]. This assumption is feasible for coronary artery. The time-dependent Navier-Stokes equations govern the blood flow in the coronary artery, as
(1)∇·U→=0,ρ∂U→∂t+ρ(U→·∇)U→=-∇p+μ∇2U→,
where U→ is the velocity vector, while p, ρ, and μ denote pressure, density, and viscosity, respectively. The maximum Reynolds number, Re=ρU-D/μ (where D is the inlet diameter of coronary and U- is the maximum pulsating inflow velocity), is calculated as 232. The inlet diameter of the distal left anterior descending coronary artery is 1.9±0.4 mm, and the maximum pulsating inflow velocity is taken as 0.4 m/s [14]. Figure 1(e) demonstrates the pulse for inlet velocity. The critical Reynolds number for an unsteady flow ranges from 5000 to 20,000 [15–17]. Therefore, the laminar flow assumption is reasonable in the present simulation.The pulse of the inlet velocity (Figure1(f)) applied in the proximal site of the right coronary artery (RCA) is decomposed into a Fourier series of trigonometric basis functions with 3 harmonics. The input velocity shown in Figure 1(f) is based on the velocity waveforms measured with a commercial ultrasound Doppler probe FloWire (Volcano TM Corporation) with 1 kHz sampling rate at the proximal site of RCA, provided by Torii et al. [18]. The inlet velocity was kept constant across the inlet surface at each time instant of the cardiac cycle. The points indicated by T1, T2, and T3 in Figure 1(f) correspond to the diastole, the peak systole, and accelerating phase within the cardiac cycle, respectively. The outflow boundary conditions were imposed at the outlets. The coronary arterial walls are assumed to be rigid with no-slip boundary condition on the luminal surface. The governing equations along with the given boundary conditions are solved with the finite-volume (FV) solver package, Fluent V12.1.4. The SIMPLE algorithm is used for the coupling of the pressure-velocity terms. The simulations are carried out for four inlet cycles to ensure that the spatial and temporal flow dynamics are periodic. The total time of 1 s with the time step size of 1 ms is taken for each cycle. The simulations are performed for four cardiac cycles, of which only the results of the forth cycle are demonstrated in the following section at different times to ensure that all the phenomena associated with the initial conditions are damped out during the first three cycles. The oscillatory shear index for the pulsatile flow simulation is calculated as [4]
(2)OSI=0.5[1-|∫0Tτ→wdt|∫0T|τ→w|dt].


## 2.1. Reconstruction of 2D Geometry of Coronary Artery from Histological Image

Two-dimensional (2D) geometries of the human coronaries are reconstructed from the histological images which were acquired from coronary arteries of four donor hearts (Faculty of Medicine and Density, University of Bristol, Bristol, UK). The samples contained regions with mild to severe atherosclerosis. The diseased coronary arteries were cut in longitudinal section, stained with H&E, and imaged on a nanozoomer, namely, NanoZoomer 2.0HT (Hamamatsu photonics, Hamamatsu, Japan) with 40x magnification at the pathology core lab in the Winship Cancer Institute, Emory University.To construct the coronary artery geometries, all the images were processed individually. The images were converted first to the boundary plots within MATLAB v. 7.11.0 (R2010b). The lumen area was identified via image segmentation. Points representing the surface of the endothelium were obtained as the final output of this step. The vertices acquired from the processing of the images were connected to generate a 2D surface representing the lumen of coronary. The 2D model was trimmed with the CFD preprocessing and meshing code Gambit V2.4.6 to generate the flat inlet and outlet boundaries. The computational mesh grid consisted of about 50,000 grid elements. Thin boundary layers were used for the precise calculations near the walls. Several steady-state simulations were performed with different grid sizes to ensure the mesh independence of the results. The boundary layer mesh used in these simulations consists of 6 layers with the first layer attached to the luminal surface with the thickness of 0.01 mm. The thickness of consecutive layers grows with the factor of 1.2. Figures1(a)–1(d) demonstrate the original histological image of each coronary artery. Figure 1(e) shows the grids generated for the models.(a)–(d) LS images of four diseased coronary arteries, (e) the computational meshes at the stenosis location, and (f) velocity pulse at the inlet.
(a)
(b)
(c)
(d)
(e)
(f)


## 2.2. Governing Equations and Boundary Conditions

The blood flow is assumed to beincompressible with uniform properties including the viscosity of3.5×10-3 Pa·s (Newtonian fluid) and the density of 1067 kgm−3. The assumption of Newtonian fluid for the blood is valid wherever the shear rate is greater than 100 S−1 [9, 10, 12, 13]. This assumption is feasible for coronary artery. The time-dependent Navier-Stokes equations govern the blood flow in the coronary artery, as
(1)∇·U→=0,ρ∂U→∂t+ρ(U→·∇)U→=-∇p+μ∇2U→,
where U→ is the velocity vector, while p, ρ, and μ denote pressure, density, and viscosity, respectively. The maximum Reynolds number, Re=ρU-D/μ (where D is the inlet diameter of coronary and U- is the maximum pulsating inflow velocity), is calculated as 232. The inlet diameter of the distal left anterior descending coronary artery is 1.9±0.4 mm, and the maximum pulsating inflow velocity is taken as 0.4 m/s [14]. Figure 1(e) demonstrates the pulse for inlet velocity. The critical Reynolds number for an unsteady flow ranges from 5000 to 20,000 [15–17]. Therefore, the laminar flow assumption is reasonable in the present simulation.The pulse of the inlet velocity (Figure1(f)) applied in the proximal site of the right coronary artery (RCA) is decomposed into a Fourier series of trigonometric basis functions with 3 harmonics. The input velocity shown in Figure 1(f) is based on the velocity waveforms measured with a commercial ultrasound Doppler probe FloWire (Volcano TM Corporation) with 1 kHz sampling rate at the proximal site of RCA, provided by Torii et al. [18]. The inlet velocity was kept constant across the inlet surface at each time instant of the cardiac cycle. The points indicated by T1, T2, and T3 in Figure 1(f) correspond to the diastole, the peak systole, and accelerating phase within the cardiac cycle, respectively. The outflow boundary conditions were imposed at the outlets. The coronary arterial walls are assumed to be rigid with no-slip boundary condition on the luminal surface. The governing equations along with the given boundary conditions are solved with the finite-volume (FV) solver package, Fluent V12.1.4. The SIMPLE algorithm is used for the coupling of the pressure-velocity terms. The simulations are carried out for four inlet cycles to ensure that the spatial and temporal flow dynamics are periodic. The total time of 1 s with the time step size of 1 ms is taken for each cycle. The simulations are performed for four cardiac cycles, of which only the results of the forth cycle are demonstrated in the following section at different times to ensure that all the phenomena associated with the initial conditions are damped out during the first three cycles. The oscillatory shear index for the pulsatile flow simulation is calculated as [4]
(2)OSI=0.5[1-|∫0Tτ→wdt|∫0T|τ→w|dt].


## 3. Results and Discussion

### 3.1. Wall Shear Stress

The wall shear stress is among the first mechanisms proposed to relate the blood flow to the localization of atherosclerosis [18] and the plaque rupture [19, 20]. Figures 2–5 demonstrate the WSS distributions along the coronary walls versus the curve length (S), at three different time instants, for four different diseased coronary arteries. In the present study, the WSS is defined as the magnitude of all three components at any location on the wall.(a) Boundary plots of the coronary arterial walls 1, 2, and 3 shown in red, green, and blue, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)(a) Boundary plots of the coronary arterial walls 1 and 2 shown in green and blue, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)(a) Boundary plots of the coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)(a) Boundary plots of coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)Figure2(a) demonstrates the boundary plots for the images of the first stenosed coronary artery (Figure 1(a)). The WSS distributions along the coronary walls versus curve length at diastole, the peak systole, and accelerating phase are shown in Figures 2(b)–2(d), respectively. The areas with maximum WSS can be observed along the distal side and on the throat of the stenosis at diastole time point. Note that a minimum value of WSS appears at the top of the atherosclerotic plaque at the peak systole and accelerating time points. It is worth mentioning that the surfaces of the coronary walls, in the present study, are generated based on the original images of patients. Surfaces are not smoothed after being read in MATLAB. The reason is associated with the aim of the current study at the investigation of the connection of atherosclerotic plaque to the local characteristics of the luminal surface. As Figure 2(a) shows, the wall on the top of the plaque is not smooth which is according to the original histological image (Figure 1(a)). This is the region with a minimum magnitude of WSS.Figures3(a)–3(d) show the WSS versus curve length along the walls of the second stenosed coronary artery (Figure 1(b)). The maximum magnitude of WSS occurs at the top of the plaque. Another stenosed region shows the second peak of WSS values at systole and accelerating time points.Figures4(a)–4(d) demonstrate the WSS along the walls of the third diseased coronary artery (Figure 1(c)), while Figures 5(a)–5(d) demonstrate the WSS along the walls of the forth diseased coronary artery (Figure 1(d)). Maximum of the WSS magnitude occurs at the throat of the stenoses in both models, at diastole time point. As seen in Figures 4(c) and 4(d), another maximum value of WSS is observed at systole and accelerating time points in the region with high curvature. In Figures 5(c) and 5(d), the proximal side of the plaque shows relatively high magnitude of WSS at systole and accelerating time points. In this artery the stenosis has been created at the bifurcation.


### 3.2. Oscillatory Shear Index (OSI)

OSI is a known predictor of formation of atherosclerosis and vulnerability for plaque in coronary arteries. Figures6(b)–6(d) demonstrates OSI for the first stenosed coronary artery. Figure 6(b) shows that OSI has high values at the high curvature region. Therefore, this location may be predicted to experience lesion in the future. As shown in Figure 6(c), OSI is high along the stenosed wall of coronary artery (shown in green in Figure 6(a)). OSI distribution reveals that the whole length of the green wall is vulnerable to atherosclerosis. Figure 6(d) demonstrates a relatively high OSI values at the curvature region along the wall.(a) Boundary plots of the coronary arterial walls 1, 2 and 3 shown in red, green, and blue, respectively. (b)–(d) OSI distribution along coronary walls 1, 2, and 3, respectively.
(a)
(b)
(c)
(d)Figures7, 8, and 9 show OSI for the next stenosed coronary arteries. Maximum values of OSI occur in stenosed regions, in all arteries. The next area with maximum of OSI is the region with high curvature. Note that this behavior is observed even in healthy curves.(a) Boundary plots of the coronary arterial walls 1 and 2 shown in green and blue, respectively. (b)-(c) OSI distribution along coronary walls 1 and 2.
(a)
(b)
(c)(a) Boundary plots of coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(e) OSI distribution along coronary walls 1, 2, 3, and 4.
(a)
(b)
(c)
(d)
(e)(a) Boundary plots of coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(e) OSI distribution along coronary walls 1, 2, 3, and 4.
(a)
(b)
(c)
(d)
(e)


## 3.1. Wall Shear Stress

The wall shear stress is among the first mechanisms proposed to relate the blood flow to the localization of atherosclerosis [18] and the plaque rupture [19, 20]. Figures 2–5 demonstrate the WSS distributions along the coronary walls versus the curve length (S), at three different time instants, for four different diseased coronary arteries. In the present study, the WSS is defined as the magnitude of all three components at any location on the wall.(a) Boundary plots of the coronary arterial walls 1, 2, and 3 shown in red, green, and blue, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)(a) Boundary plots of the coronary arterial walls 1 and 2 shown in green and blue, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)(a) Boundary plots of the coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)(a) Boundary plots of coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(d) WSS distribution along coronary walls at timeT1, T2, and T3, respectively.
(a)
(b)
(c)
(d)Figure2(a) demonstrates the boundary plots for the images of the first stenosed coronary artery (Figure 1(a)). The WSS distributions along the coronary walls versus curve length at diastole, the peak systole, and accelerating phase are shown in Figures 2(b)–2(d), respectively. The areas with maximum WSS can be observed along the distal side and on the throat of the stenosis at diastole time point. Note that a minimum value of WSS appears at the top of the atherosclerotic plaque at the peak systole and accelerating time points. It is worth mentioning that the surfaces of the coronary walls, in the present study, are generated based on the original images of patients. Surfaces are not smoothed after being read in MATLAB. The reason is associated with the aim of the current study at the investigation of the connection of atherosclerotic plaque to the local characteristics of the luminal surface. As Figure 2(a) shows, the wall on the top of the plaque is not smooth which is according to the original histological image (Figure 1(a)). This is the region with a minimum magnitude of WSS.Figures3(a)–3(d) show the WSS versus curve length along the walls of the second stenosed coronary artery (Figure 1(b)). The maximum magnitude of WSS occurs at the top of the plaque. Another stenosed region shows the second peak of WSS values at systole and accelerating time points.Figures4(a)–4(d) demonstrate the WSS along the walls of the third diseased coronary artery (Figure 1(c)), while Figures 5(a)–5(d) demonstrate the WSS along the walls of the forth diseased coronary artery (Figure 1(d)). Maximum of the WSS magnitude occurs at the throat of the stenoses in both models, at diastole time point. As seen in Figures 4(c) and 4(d), another maximum value of WSS is observed at systole and accelerating time points in the region with high curvature. In Figures 5(c) and 5(d), the proximal side of the plaque shows relatively high magnitude of WSS at systole and accelerating time points. In this artery the stenosis has been created at the bifurcation.


## 3.2. Oscillatory Shear Index (OSI)

OSI is a known predictor of formation of atherosclerosis and vulnerability for plaque in coronary arteries. Figures6(b)–6(d) demonstrates OSI for the first stenosed coronary artery. Figure 6(b) shows that OSI has high values at the high curvature region. Therefore, this location may be predicted to experience lesion in the future. As shown in Figure 6(c), OSI is high along the stenosed wall of coronary artery (shown in green in Figure 6(a)). OSI distribution reveals that the whole length of the green wall is vulnerable to atherosclerosis. Figure 6(d) demonstrates a relatively high OSI values at the curvature region along the wall.(a) Boundary plots of the coronary arterial walls 1, 2 and 3 shown in red, green, and blue, respectively. (b)–(d) OSI distribution along coronary walls 1, 2, and 3, respectively.
(a)
(b)
(c)
(d)Figures7, 8, and 9 show OSI for the next stenosed coronary arteries. Maximum values of OSI occur in stenosed regions, in all arteries. The next area with maximum of OSI is the region with high curvature. Note that this behavior is observed even in healthy curves.(a) Boundary plots of the coronary arterial walls 1 and 2 shown in green and blue, respectively. (b)-(c) OSI distribution along coronary walls 1 and 2.
(a)
(b)
(c)(a) Boundary plots of coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(e) OSI distribution along coronary walls 1, 2, 3, and 4.
(a)
(b)
(c)
(d)
(e)(a) Boundary plots of coronary arterial walls 1, 2, 3, and 4 shown in red, green, blue, and brown, respectively. (b)–(e) OSI distribution along coronary walls 1, 2, 3, and 4.
(a)
(b)
(c)
(d)
(e)


## 4. Conclusions

In this study, pulsatile blood flow is studied in four reconstructed models of the diseased coronary arteries from histological images. Surface smoothing and approximation of the vessels with regular shapes such as straight or curved cylinders are avoided unlike most of earlier studies. The real geometric model and the conditions of pulsatile flow allow to spot vulnerable areas to the development of atherosclerotic lesions. The results indicate a strong dependence of the hemodynamics on the morphology of atherosclerotic lesion. Distinctive flow patterns appear in different stenosed regions corresponding to the specific geometry of the arteries. The OSI is a well-known predictor of formation of atherosclerosis and vulnerability for plaque in coronary arteries. High OSI has been widely used as indicators of atherosclerosis because the region with high OSI is likely to experience stagnation or backflow. This is strongly related to the progression of atherosclerosis because it affects the arrangement of endothelial cells in adjacent tissues. It has been shown that the areas with high values of OSI are usually located in the regions where wall shear stress is low. In the present study, OSI is high in diseased region as well as healthy regions with low values of WSS which can be predicted to experience lesion in the future. The findings support the low shear stress and oscillatory shear stress theory which relates abnormal biological effects inside arteries to high OSI. Moreover, the present study demonstrates the shear environment in longitudinal section of the diseased coronary. The stenosis affects the WSS locally along the diseased arterial wall as well as adjacent walls. The maximum value of WSS is observed in the throat of stenosis which highlights the risk of plaque rupture. The arterial plaque can break apart as a result of high stresses acting on or within the arterial wall. The severity and the peripheral distribution of stenosis are the crucial geometric factors influencing the risk of breakage.It is hoped that computational studies of this type will enable scientists to directly relate the hemodynamics to the formation and progression of atherosclerosis in the coronary arteries. More than 60% of myocardial infarction is caused by the rupture of vulnerable plaques [17]. A significant contribution of the present study is the demonstration of robust and highly efficient computational methods for studying the pulsatile blood flow field in different stenosed coronary arteries. The current study is the first detailed numerical investigation of the flow field within the diseased coronary arteries which are reconstructed from histological images. The spatial and temporal variation of wall shear stress from the present study can be fed into other models to investigate the effect of wall shear stress on stresses and strains that individual endothelial cells may tolerate locally.



---

*Source: 715407-2013-06-17.xml*