Biomed. Eng.-Biomed. Tech. 2014; 59(6): 471–477
Tim A.S. Kaufmann*, Christoph Leisser, Jeannie Gemsa and Ulrich Steinseifer
Analysis of emboli and blood flow in the ophthalmic artery to understand retinal artery occlusion Abstract: Retinal artery occlusion (RAO) is a common ocular vascular occlusive disorder that may lead to partial or complete retinal ischemia with sudden visual deterioration and visual field defects. Although RAO has been investigated since 1859, the main mechanism is still not fully understood. While hypoperfusion of the ophthalmic artery (OA) due to severe stenosis of the internal carotid artery might lead to RAO, emboli are assumed to be the main reason. Intra-arterial thrombolysis is not a sufficient treatment for RAO, and current research is mainly focused on risk factors. In this study, a computational fluid dynamic model is presented to analyse flow conditions and clot behaviour at the junction of the internal carotid artery and OA based on a realistic geometry from a RAO patient. Clot diameters varied between 5 and 200 μm, and the probability of clots reaching the OA or being washed into the brain was analysed. Results show sufficient blood flow and perfusion pressure at the end of OA. The probability that clots from the main blood flow will to be washed into the brain is 7.32 ± 1.08%. A wall shear stress hotspot is observed at the curvature proximal to the internal carotid artery/OA junction. Clots released from this hotspot have a higher probability of causing RAO. The occurrence of such patient-specific pathophysiologies will have to be considered in the future. Keywords: computational fluid dynamics; embolism; retinal artery occlusion. DOI 10.1515/bmt-2014-0002 Received January 13, 2014; accepted June 20, 2014; online first July 16, 2014
*Corresponding author: Tim A.S. Kaufmann, Institute of Applied Medical Engineering, Helmholtz Institute, Department of Cardiovascular Engineering, RWTH Aachen University, Pauwelsstraße 20, 52074 Aachen, Germany, Phone: +492418080540, Fax: +492418082144, E-mail: [email protected]; and enmodes GmbH, Aachen, Germany Christoph Leisser: Department of Ophthalmology, Hanuschkrankenhaus, Vienna, Austria Jeannie Gemsa and Ulrich Steinseifer: Institute of Applied Medical Engineering, Helmholtz Institute, Department of Cardiovascular Engineering, RWTH Aachen University, Aachen, Germany
Introduction Retinal artery occlusion (RAO) is a common ocular vascular occlusive disorder first reported on in 1859 [31] that leads to partial or complete retinal ischemia with sudden visual deterioration and visual field defects [11]. Several forms of RAO exist [8–10], and the mechanism behind them is still not fully understood [6]. RAO is associated with systemic cardiovascular disease [3, 4, 18, 19, 27, 32–34], and one reason for this is hypoperfusion due to severe stenosis of the carotid artery, but only for stenosis > 70% [1, 7, 24, 29]. It is assumed that systemic emboli, mainly from the carotid artery but also possibly from endocarditis, are the main reason for RAO [11]. Mead et al. [23] found that emboli released from the heart valve are more likely to cause stroke than RAO, probably because they are located in the central flow, while emboli from a carotid artery stenosis are more likely to cause occlusion of retinal arteries since the higher flow circulation due to the stenosis causes a higher concentration of emboli in the near wall flow. Since the EAGLE study showed that intra-arterial thrombolysis is not beneficial for patients with RAO [28], current research has mainly focused on risk factors for RAO [1, 3, 4, 11, 18, 19, 21, 23, 27, 32–34]. Thus, understanding of flow conditions and blood clot dynamics in the internal carotid artery (ICA) and ophthalmic artery (OA) may help to improve treatment concepts for RAO. Systematic analysis of emboli with respect to strokes and RAO is limited, however, and can currently only be done in large clinical trials, most of them analysed retrospectively. In the past, we developed an experimentally validated computational model to analyse flow conditions in the human aortic arch [14, 15, 20]. This model has now been adapted to the flow at the junction of ICA and OA. It allows analysis of flow conditions and the behaviour of blood clots at that junction, giving information about the circumstances under which RAO can occur, since only blood clots washed into the OA contribute to the development of RAO.
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
472 T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion
Materials and methods
Simulation
Model creation
Flow simulations
Based on magnetic resonance imaging (MRI) scans of the human head of a RAO patient, a three-dimensional computer-aided design model of the ICA and the branch to the OA was created by Hounsfield units using the DICOM editing program Mimics (Materialise Inc., Leuven, Belgium). Automatic and manual fixing and smoothing iterations were performed in the computer-aided engineering tool 3-matic (Materialise Inc., Leuven, Belgium). The final model was cut before and after the junction of ICA and OA to limit the numerical domain to this section. The process of model creation is shown in Figure 1. The MRI scans were selected retrospectively out of a group of 28 patients with RAO [21], including three patients who had cerebral MRI scans for neurological reasons. Only one of these sets of MRI scans offered adequate data for reconstruction of the ICA/OA junction.
Numerical flow simulations in the ICA were performed for pulsatile flow conditions using commercial software (ANSYS CFX 14.0, Ansys Germany Inc., Otterfing, Germany). A specified blend factor of 0.75 was set as the advection scheme, meaning that at least 75% of all elements are solved second-order-accurate while up to 25% may be solved first-order-accurate. The time step size was set to 0.01 s. Up to ten internal coefficient loops were solved per time step until the average changes in the transport equations were smaller than the specified convergence target of 1e–4. The flow was assumed to be laminar and blood was modelled as a non-Newtonian fluid according to [2], as shown in equations (1)–(3), with a haematocrit of 44% and a density of 1.056 kg/l.
η = λshear γ
n
shear-1
(1)
d
Mesh generation A mesh of tetrahedral elements was generated within the computer-aided design model using commercial software (ANSYS ICEM CFD, Ansys Germany Inc., Otterfing, Germany). Since a no-slip condition was assumed at the vessel walls, high velocity gradients are expected close to the wall boundaries and five layers of prismatic elements were created around the vessel wall to resolve the boundary layers of the fluid. A mesh independency study was initially undertaken to assure that no information was lost due to mesh size, which varied between 250,000 and 2,000,000 elements. Mesh convergence was analysed for blood flow, particle distribution and wall shear stress. After the initial mesh study, a mesh density of 750,000 elements was chosen for further analysis.
Figure 1 Process of model creation. (l) DICOM data and model reconstruction; (m) reconstructed model; (r) final computer-aided design model.
nshear = n∞ -∆n e
γ - 1+ e γ c
(2)
b
λshear = η∞ + ∆η e
γ - 1+ e γ a
(3)
With n∞ = 1, η∞ = 0.00035 Pa s, a = 70.71 s-1, b = 4.24 s-1, c = 70.71 -1 s , d = 5.66 s-1, Δn = 0.45, Δη = 0.025, Pa s The vessel walls were considered rigid with a no-slip condition. The opening to the inferior ICA was modelled as an inlet. Pulsatile inflow conditions were assumed based on a transient flow profile that was normalized to account for an average systolic velocity of 0.9 m/s [13], which corresponds to a systolic pressure of 126 mm Hg. The cycle time was set to 1 s with a systolic time of 0.3 s. A consecutive number of cardiac cycles were simulated for each condition to avoid influences by the initial guess, which was achieved after three cycles. The flow profile is shown in Figure 2. It was additionally modulated ± 50% in steps of 10% to allow statistical analysis and account for patient variations. The results were averaged over the variations in the inlet flow profile. The boundaries to the superior ICA and the OA were modelled as pressure openings at 80 mm Hg. In addition to this, flow-dependent vascular resistances were set for each vessel, assuming a pressure drop ΔP according to the Darcy-Weißbach equation (4), as presented in Kaufmann et al. [17]:
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion 473
Inlet profile
The number of integration steps per element was set to 100, and the maximum number of integration steps was set to 100,000. These values were chosen based on initial results to ensure that particles that are trapped in areas with low flow velocities are still accounted for. The number of particle starting positions was set to 1500 and the path of each particle was tracked. The probability for a clot being washed into the ophthalmic artery was analysed dependent on the clot size. The particle starting points were defined as being spread equally at the blood flow inlet, representing particles coming from the heart valve or left ventricle. In a second series of simulations, particles were injected from a hotspot of wall shear stress (WSS), which was observed during the first simulations.
0.04 0.035
Mass flow (kg s-1)
0.03 0.025 0.02 0.015 0.01 0.005 0 -0.005 0
0.2
0.4 0.6 Time (s) Inlet flow
0.8
1
Results
Figure 2 Inlet flow profile at standard flow rate.
1 ∆P = L ρU n2 2
(4)
In this equation, ρ is the fluid density and Un the normal velocity at each boundary. The free variable L is defined as the dimensionless loss coefficient, which represents the systemic resistance. L is defined in a way that represents the pressure drop from arterial to venous pressure. The static pressure at each opening is calculated using equation (5):
Pstatic = Prelative + ∆P
(5)
In Kaufmann et al. [16, 17] we found that a loss coefficient of 50–100 is feasible for the aortic arch and outgoing vessels. However, the region of interest of this study is well behind the aortic arch, which should result in a lower value for L. It can be estimated by assuming an average flow velocity of 0.8 m/s at the internal carotid artery and a pressure drop of 75 mm Hg, resulting in a value of 25 for the loss coefficient.
Particle tracking In addition to regular fluid simulations, analysis of blood clots in the flow was performed. Particles with the density of platelets were injected in the flow using a fully coupled Lagrangian particle tracking approach with a SchillerNaumann drag force. Discrete particle diameters of 5, 25, 50, 100 and 200 μm were analysed.
Flow fields Of the mass flow assumed at the inlet (ICA), 92% ± 0.5% were distributed to the superior ICA and 8% ± 0.5% to the OA. The averaged flow velocity at the outlet representing OA was 0.6 m/s ± 0.3 m/s. Calculated from the relative pressure (80 mm Hg) and the pressure drop due to resistance, the perfusion pressure at the end of OA is approximately 49 mm Hg. The flow field at mid-systole is shown in Figure 3. Proximal to the junction, a WSS hotspot is observed that is depicted in Figure 4. It correlates with increased velocities of about 2.5 m/s at the curvature. The maximum WSS is 236 Pa at peak systolic flow. The hotspot occurs only during systole, is highly localized and does not change its position. During diastole, the WSS distribution has a different profile. A comparison with the DICOM images revealed reduced Hounsfield units of 180 at the position of the WSS hotspot, which is in the range of mixed or calcified plaques [12, 26].
Particle tracks Two separate sets of simulations were performed. In the first set, 1500 particles were released from the inlet representing the inferior ICA, thus describing microemboli in the blood flow. The particle tracks for this scenario are shown in Figure 5. All particles that end in the OA are seeded in the near wall flow. In the second set of simulations, the same number of particles was released from the WSS hotspot that
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
474 T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion Velocity 1.2
O6 C6
O1
C7
O2 O3
O4
O2
O1
O5
O3
O4
O6
O5
0.9 0.6
C5
0.3 C4 C4
C2
0 (m s-1) Velocity 2.5
C6
C3 C1
C3
C5
C2
C7
1.9 1.3
C1
0.6 0 (m s-1)
Figure 3 Flow distribution at peak systole in the carotid (lower) and ophthalmic (upper) artery.
Wall shear wall 240 180 120 60 0 (Pa)
Figure 4 Location of wall shear stress hotspot during systole proximal to the ophthalmic artery (l) and possible red-coloured plaque in the curvature based on Hounsfield units (r).
was identified during the flow simulations. Close to the hotspot, a secondary flow pattern was observed, which channels emboli into OA. This is depicted in Figure 6. The overall probability for clots to reach the OA is 7.32 ± 1.08% for microemboli from the inferior ICA and 13.27 ± 6.82% for microemboli released from the WSS hotspot. The high standard deviation of the probability of emboli from the WSS hotspot reaching the OA is due to the fact that at low flow rates (50–70% of average), only 3.43 ± 2.42% of the emboli from the WSS hotspot are washed into OA. This percentage is significantly increased to 18.20 ± 2.83% at higher flow rates. The overall distribution of particles dependent on clot size is presented in Table 1. The average travel time for emboli reaching the OA is 164 ms and 313 ms for emboli in the blood stream or
released from the WSS hotspot, respectively. For emboli washed into the ICA, the travel times are slightly shorter, at 116 ms (blood stream) and 241 ms (WSS hotspot). This is mainly due to the shorter distance between both boundaries.
Discussion Results Computational Fluid Dynamics (CFD) simulations were performed in a realistic RAO patient geometry to analyse the blood flow and clot behaviour at the junction of the OA and ICA. The results show a sufficient perfusion of the
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion 475
Table 1 Particles reaching the ophthalmic artery (OA) dependent on particle diameter. Particle diameter ( μm) 5 25 50 100 200
Figure 5 Emboli seeded in the internal carotid artery reaching the ophthalmic artery.
OA; however, this might be affected by possible stenosis of the ICA prior to the modelled geometry, dependent on its position and grade [1, 3, 4, 7, 24, 27, 29]. The pressure at the end of the OA equals the perfusion pressure of the retinal artery. The simulated pressure of 49 mm Hg is in good agreement with data in the literature, which describe a retinal artery perfusion pressure of 45–55 mm Hg [25]. Velocity 2.0 1.5 1.0 0.5 0 (m s-1)
Figure 6 Upper: Secondary flow pattern at a plane at the WSS hotspot. Lower: Flow path of emboli following main flow (blue) and secondary flow (red).
Particles reaching the OA (%) originating from the OA 8.87 ± 1.50 6.20 ± 0.76 5.90 ± 0.89 8.13 ± 0.97 7.52 ± 1.31
Particles reaching the OA (%) originating from hotspot
14.24 ± 7.28 13.98 ± 7.09 13.79 ± 7.15 12.65 ± 6.64 11.69 ± 5.95
In addition to CFD simulations, the behaviour of emboli in the blood stream was analysed using Lagrangian particle tracking. All particles from the main blood flow that are washed into the OA and might therefore cause RAO were seeded in the near wall flow indicating that particles released from the vessel walls might have a higher potential to cause RAO than particles coming from the heart valve, as already assumed in Mead et al. [23]. The simulations revealed varying WSSs throughout the cycle, with a constant localized hotspot of 240 Pa during systole. Such a WSS hotspot does not correspond to plaque formation, but it might lead to mobilization of existing plaques [5, 30]. Since intravascular plaques may alter flow conditions, this is one reason for RAO that has to be taken into consideration. In addition, secondary flow patterns at the junction of the ICA and OA were observed that channel fluid of the near wall flow of the curvature into the OA, especially at flow rates above 70% of average. Given that mobilization of atherosclerotic plaques is more likely at high flow rates, the effect of plaque mobilization might be even more significant. If a patient shows plaque formation proximal to the curvature, it is therefore possible that such a WSS hotspot might cause RAO by releasing emboli. The occurrence of atherosclerotic plaques in RAO patients proximal to the junction is currently being investigated in a retrospective analysis of angiographic patient data. According to the first results from the numerical simulation of clots behaviour at the ICA/OA junction, one can assume that more than 90% of all clots coming from the ICA are washed into the brain. Astonishingly, the rate of patients having had a stroke previous to RAO is low, at 19% or below [4, 11, 21, 22]. This raises questions about the reasons for RAO. It may happen that small particles do not always cause symptoms if washed into the brain, but may cause RAO if washed into the OA. Also, some of the blood clots causing RAO might not pass the ICA/OA junction, but develop behind it in the OA. Another explanation
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
476 T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion is the aforementioned occurrence of WSS hotspots and pathophysiologies causing secondary flows.
Limitations The model itself is a limitation to this study. Simulations were performed in a vascular system that was averaged over systole and diastole. Possible movement of the vessel wall was neglected. Additionally, no further information is available on the patient condition. It is therefore not possible to separate risk factors in this first feasibility study. A future study will focus on the analysis of several RAO patients, including risk screening and an analysis of parameter correlations. Regarding the particle tracks, only particles with the density of platelets were considered. In future studies, different clot materials, especially plaques, should also be included. The tracks of the particles are very sensitive to small changes in the flow field. Thus, more studies with small variations in the flow parameters and an increased particle number should be performed in order to increase the accuracy of the model. This is especially true for larger particles due to the nature of the Lagrangian particle tracking method.
Conclusions The model presented has the potential to provide insight into the mechanisms behind RAO. It allows us to analyse both general and patient-specific risk factors in a numerical framework, whereby the effects of various parameters may be studied. Further simulations and validating experiments are, however, needed to improve the credibility of this method. Ultimately, results will be compared to an ongoing clinical trial of RAO patients. Conflict of interest statement: The authors declare that there are no conflicts of interests.
References [1] Anderson DC, Kappelle LJ, Eliasziw M, et al. Occurrence of hemispheric and retinal ischemia in atrial fibrillation compared with carotid stenosis. Stroke 2002; 33: 1963–1967. [2] Ballyk PD, Steinman DA, Ethier CR. Simulation of non-Newtonian blood flow in an end-to-side anastomosis. Biorheology 1994; 31: 565–586.
[3] Cheung N, Lim L, Wang JJ, et al. Prevalence and risk factors of retinal arteriolar emboli: the Singapore Malay Eye Study. Am J Ophthalmol 2008; 146: 620–624. [4] Dunlap AB, Kosmorsky GS, Kashyap VS. The fate of patients with retinal artery occlusion and Hollenhorst plaque. J Vasc Surg 2007; 46: 1125–1129. [5] Groen HC, Gijsen FJ, van der Lugt A, et al. Plaque rupture in the carotid artery is localized at the high shear stress region: a case report. Stroke 2007; 38: 2379–2381. [6] Hayreh SS. Prevalent misconceptions about acute retinal vascular occlusive disorders. Prog Retin Eye Res 2005; 24: 493–519. [7] Hayreh SS, Podhajsky P. Ocular neovascularization with retinal vascular occlusion. II. Occurrence in central and branch retinal artery occlusion. Arch Ophthalmol 1982; 100: 1585–1596. [8] Hayreh SS, Zimmerman MB. Central retinal artery occlusion: visual outcome. Am J Ophthalmol 2005; 140: 376–391. [9] Hayreh SS, Zimmerman MB. Fundus changes in central retinal artery occlusion. Retina 2007; 27: 276–289. [10] Hayreh SS, Podhajsky PA, Zimmerman MB. Branch retinal artery occlusion: natural history of visual outcome. Ophthalmology 2009; 116: 1188–1194.e1-4. [11] Hayreh SS, Podhajsky PA, Zimmerman MB. Retinal artery occlusion: associated systemic and ophthalmic abnormalities. Ophthalmology 2009; 116: 1928–1936. [12] Hur J, Kim YJ, Lee HJ, et al. Quantification and characterization of obstructive coronary plaques using 64-slice computed tomography: a comparison with intravascular ultrasound. J Comput Assist Tomogr 2009; 33: 186–192. [13] Kaiser HJ, Schotzau A, Flammer J. Blood-flow velocities in the extraocular vessels in normal volunteers. Am J Ophthalmol 1996; 122: 364–370. [14] Kaufmann TAS, Hormes M, Laumen M, et al. Flow distribution during cardiopulmonary bypass in dependency on the outflow cannula positioning. Artificial Organs 2009; 33: 988–992. [15] Kaufmann TAS, Hormes M, Laumen M, et al. The impact of aortic/subclavian outflow cannulation for cardiopulmonary bypass and cardiac support: A computational fluid dynamics study. Artificial Organs 2009; 33: 727–732. [16] Kaufmann TAS, Schmitz-Rode T, Steinseifer U. Implementation of cerebral autoregulation into computational fluid dynamics studies of cardiopulmonary bypass. Artif Organs 2012; 36: 754–758. [17] Kaufmann TAS, Wong KC, Schmitz-Rode T, Steinseifer U. Mimicking of cerebral autoregulation by flow-dependent cerebrovascular resistance: a feasibility study. Artif Organs 2012; 36: E97–E101. [18] Klein R, Klein BE, Jensen SC, et al. Retinal emboli and stroke: the Beaver Dam Eye Study. Arch Ophthalmol 1999; 117: 1063–1068. [19] Klein R, Klein BE, Moss SE, Meuer SM. Retinal emboli and cardiovascular disease: the Beaver Dam Eye Study. Arch Ophthalmol 2003; 121: 1446–1451. [20] Laumen M, Kaufmann T, Timms D, et al. Flow analysis of ventricular assist device inflow and outflow cannula positioning using a naturally shaped ventricle and aortic branch. Artificial Organs 2010; 34: 798–806. [21] Leisser C. Hyperlipidaemia and increased systemic blood pressure – is there a connection to retinal artery occlusion? Spektrum Augenheilkd 2010; 24: 154–156.
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion 477 [22] Leisser C. Risk factors of retinal artery occlusion. Spektrum Augenheilkd 2012; 26: 202–204. [23] Mead GE, Lewis SC, Wardlaw JM, Dennis MS. Comparison of risk factors in patients with transient and prolonged eye and brain ischemic syndromes. Stroke 2002; 33: 2383–2390. [24] Mizener JB, Podhajsky P, Hayreh SS. Ocular ischemic syndrome. Ophthalmology 1997; 104: 859–864. [25] Pillunat LE, Anderson DR, Knighton RW, Joos KM, Feuer WJ. Autoregulation of human optic nerve head circulation in response to increased intraocular pressure. Exp Eye Res 1997; 64: 737–744. [26] Saba L, Mallarin G. Window settings for the study of calcified carotid plaques with multidetector CT angiography. AJNR Am J Neuroradiol 2009; 30: 1445–1450. [27] Schmidt D, Hetzel A, Geibel-Zehender A, Schulte-Mönting J. Systemic diseases in non-inflammatory branch and central retinal artery occlusion–an overview of 416 patients. Eur J Med Res 2007; 12: 595–603. [28] Schumacher M, Schmidt D, Jurklies B, et al. Central retinal artery occlusion: local intra-arterial fibrinolysis versus conservative treatment, a multicenter randomized trial. Ophthalmology 2010; 117: 1367–1375.e1.
[29] Sharma S, Pater JL, Lam M, Cruess AF. Can different types of retinal emboli be reliably differentiated from one another? An inter- and intraobserver agreement study. Can J Ophthalmol 1998; 33: 144–148. [30] Tang D, Teng Z, Canton G, et al. Local critical stress correlates better than global maximum stress with plaque morphological features linked to atherosclerotic plaque vulnerability: an in vivo multi-patient study. Biomed Eng Online, 2009; 8: 15. [31] von Graefe A. Ueber Embolie der Arteriacentralisretinae als Ursache plotzlicher Erblindung. Albrecht von Graefes ArchOphthalmol 1859; 5: 136–157. [32] Weger M, Renner W, Pinter O, et al. Role of factor V Leiden and prothrombin 20210A in patients with retinal artery occlusion. Eye 2003; 17: 731–734. [33] Wong TY, Klein R. Retinal arteriolar emboli: epidemiology and risk of stroke. Curr Opin Ophthalmol 2002; 13: 142–146. [34] Wong TY, Larsen EK, Klein R, et al. Cardiovascular risk factors for retinal vein occlusion and arteriolar emboli: the Atherosclerosis Risk in Communities and Cardiovascular Health studies. Ophthalmology 2005; 112: 540–547.
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
Tim A.S. Kaufmann*, Christoph Leisser, Jeannie Gemsa and Ulrich Steinseifer
Analysis of emboli and blood flow in the ophthalmic artery to understand retinal artery occlusion Abstract: Retinal artery occlusion (RAO) is a common ocular vascular occlusive disorder that may lead to partial or complete retinal ischemia with sudden visual deterioration and visual field defects. Although RAO has been investigated since 1859, the main mechanism is still not fully understood. While hypoperfusion of the ophthalmic artery (OA) due to severe stenosis of the internal carotid artery might lead to RAO, emboli are assumed to be the main reason. Intra-arterial thrombolysis is not a sufficient treatment for RAO, and current research is mainly focused on risk factors. In this study, a computational fluid dynamic model is presented to analyse flow conditions and clot behaviour at the junction of the internal carotid artery and OA based on a realistic geometry from a RAO patient. Clot diameters varied between 5 and 200 μm, and the probability of clots reaching the OA or being washed into the brain was analysed. Results show sufficient blood flow and perfusion pressure at the end of OA. The probability that clots from the main blood flow will to be washed into the brain is 7.32 ± 1.08%. A wall shear stress hotspot is observed at the curvature proximal to the internal carotid artery/OA junction. Clots released from this hotspot have a higher probability of causing RAO. The occurrence of such patient-specific pathophysiologies will have to be considered in the future. Keywords: computational fluid dynamics; embolism; retinal artery occlusion. DOI 10.1515/bmt-2014-0002 Received January 13, 2014; accepted June 20, 2014; online first July 16, 2014
*Corresponding author: Tim A.S. Kaufmann, Institute of Applied Medical Engineering, Helmholtz Institute, Department of Cardiovascular Engineering, RWTH Aachen University, Pauwelsstraße 20, 52074 Aachen, Germany, Phone: +492418080540, Fax: +492418082144, E-mail: [email protected]; and enmodes GmbH, Aachen, Germany Christoph Leisser: Department of Ophthalmology, Hanuschkrankenhaus, Vienna, Austria Jeannie Gemsa and Ulrich Steinseifer: Institute of Applied Medical Engineering, Helmholtz Institute, Department of Cardiovascular Engineering, RWTH Aachen University, Aachen, Germany
Introduction Retinal artery occlusion (RAO) is a common ocular vascular occlusive disorder first reported on in 1859 [31] that leads to partial or complete retinal ischemia with sudden visual deterioration and visual field defects [11]. Several forms of RAO exist [8–10], and the mechanism behind them is still not fully understood [6]. RAO is associated with systemic cardiovascular disease [3, 4, 18, 19, 27, 32–34], and one reason for this is hypoperfusion due to severe stenosis of the carotid artery, but only for stenosis > 70% [1, 7, 24, 29]. It is assumed that systemic emboli, mainly from the carotid artery but also possibly from endocarditis, are the main reason for RAO [11]. Mead et al. [23] found that emboli released from the heart valve are more likely to cause stroke than RAO, probably because they are located in the central flow, while emboli from a carotid artery stenosis are more likely to cause occlusion of retinal arteries since the higher flow circulation due to the stenosis causes a higher concentration of emboli in the near wall flow. Since the EAGLE study showed that intra-arterial thrombolysis is not beneficial for patients with RAO [28], current research has mainly focused on risk factors for RAO [1, 3, 4, 11, 18, 19, 21, 23, 27, 32–34]. Thus, understanding of flow conditions and blood clot dynamics in the internal carotid artery (ICA) and ophthalmic artery (OA) may help to improve treatment concepts for RAO. Systematic analysis of emboli with respect to strokes and RAO is limited, however, and can currently only be done in large clinical trials, most of them analysed retrospectively. In the past, we developed an experimentally validated computational model to analyse flow conditions in the human aortic arch [14, 15, 20]. This model has now been adapted to the flow at the junction of ICA and OA. It allows analysis of flow conditions and the behaviour of blood clots at that junction, giving information about the circumstances under which RAO can occur, since only blood clots washed into the OA contribute to the development of RAO.
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
472 T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion
Materials and methods
Simulation
Model creation
Flow simulations
Based on magnetic resonance imaging (MRI) scans of the human head of a RAO patient, a three-dimensional computer-aided design model of the ICA and the branch to the OA was created by Hounsfield units using the DICOM editing program Mimics (Materialise Inc., Leuven, Belgium). Automatic and manual fixing and smoothing iterations were performed in the computer-aided engineering tool 3-matic (Materialise Inc., Leuven, Belgium). The final model was cut before and after the junction of ICA and OA to limit the numerical domain to this section. The process of model creation is shown in Figure 1. The MRI scans were selected retrospectively out of a group of 28 patients with RAO [21], including three patients who had cerebral MRI scans for neurological reasons. Only one of these sets of MRI scans offered adequate data for reconstruction of the ICA/OA junction.
Numerical flow simulations in the ICA were performed for pulsatile flow conditions using commercial software (ANSYS CFX 14.0, Ansys Germany Inc., Otterfing, Germany). A specified blend factor of 0.75 was set as the advection scheme, meaning that at least 75% of all elements are solved second-order-accurate while up to 25% may be solved first-order-accurate. The time step size was set to 0.01 s. Up to ten internal coefficient loops were solved per time step until the average changes in the transport equations were smaller than the specified convergence target of 1e–4. The flow was assumed to be laminar and blood was modelled as a non-Newtonian fluid according to [2], as shown in equations (1)–(3), with a haematocrit of 44% and a density of 1.056 kg/l.
η = λshear γ
n
shear-1
(1)
d
Mesh generation A mesh of tetrahedral elements was generated within the computer-aided design model using commercial software (ANSYS ICEM CFD, Ansys Germany Inc., Otterfing, Germany). Since a no-slip condition was assumed at the vessel walls, high velocity gradients are expected close to the wall boundaries and five layers of prismatic elements were created around the vessel wall to resolve the boundary layers of the fluid. A mesh independency study was initially undertaken to assure that no information was lost due to mesh size, which varied between 250,000 and 2,000,000 elements. Mesh convergence was analysed for blood flow, particle distribution and wall shear stress. After the initial mesh study, a mesh density of 750,000 elements was chosen for further analysis.
Figure 1 Process of model creation. (l) DICOM data and model reconstruction; (m) reconstructed model; (r) final computer-aided design model.
nshear = n∞ -∆n e
γ - 1+ e γ c
(2)
b
λshear = η∞ + ∆η e
γ - 1+ e γ a
(3)
With n∞ = 1, η∞ = 0.00035 Pa s, a = 70.71 s-1, b = 4.24 s-1, c = 70.71 -1 s , d = 5.66 s-1, Δn = 0.45, Δη = 0.025, Pa s The vessel walls were considered rigid with a no-slip condition. The opening to the inferior ICA was modelled as an inlet. Pulsatile inflow conditions were assumed based on a transient flow profile that was normalized to account for an average systolic velocity of 0.9 m/s [13], which corresponds to a systolic pressure of 126 mm Hg. The cycle time was set to 1 s with a systolic time of 0.3 s. A consecutive number of cardiac cycles were simulated for each condition to avoid influences by the initial guess, which was achieved after three cycles. The flow profile is shown in Figure 2. It was additionally modulated ± 50% in steps of 10% to allow statistical analysis and account for patient variations. The results were averaged over the variations in the inlet flow profile. The boundaries to the superior ICA and the OA were modelled as pressure openings at 80 mm Hg. In addition to this, flow-dependent vascular resistances were set for each vessel, assuming a pressure drop ΔP according to the Darcy-Weißbach equation (4), as presented in Kaufmann et al. [17]:
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion 473
Inlet profile
The number of integration steps per element was set to 100, and the maximum number of integration steps was set to 100,000. These values were chosen based on initial results to ensure that particles that are trapped in areas with low flow velocities are still accounted for. The number of particle starting positions was set to 1500 and the path of each particle was tracked. The probability for a clot being washed into the ophthalmic artery was analysed dependent on the clot size. The particle starting points were defined as being spread equally at the blood flow inlet, representing particles coming from the heart valve or left ventricle. In a second series of simulations, particles were injected from a hotspot of wall shear stress (WSS), which was observed during the first simulations.
0.04 0.035
Mass flow (kg s-1)
0.03 0.025 0.02 0.015 0.01 0.005 0 -0.005 0
0.2
0.4 0.6 Time (s) Inlet flow
0.8
1
Results
Figure 2 Inlet flow profile at standard flow rate.
1 ∆P = L ρU n2 2
(4)
In this equation, ρ is the fluid density and Un the normal velocity at each boundary. The free variable L is defined as the dimensionless loss coefficient, which represents the systemic resistance. L is defined in a way that represents the pressure drop from arterial to venous pressure. The static pressure at each opening is calculated using equation (5):
Pstatic = Prelative + ∆P
(5)
In Kaufmann et al. [16, 17] we found that a loss coefficient of 50–100 is feasible for the aortic arch and outgoing vessels. However, the region of interest of this study is well behind the aortic arch, which should result in a lower value for L. It can be estimated by assuming an average flow velocity of 0.8 m/s at the internal carotid artery and a pressure drop of 75 mm Hg, resulting in a value of 25 for the loss coefficient.
Particle tracking In addition to regular fluid simulations, analysis of blood clots in the flow was performed. Particles with the density of platelets were injected in the flow using a fully coupled Lagrangian particle tracking approach with a SchillerNaumann drag force. Discrete particle diameters of 5, 25, 50, 100 and 200 μm were analysed.
Flow fields Of the mass flow assumed at the inlet (ICA), 92% ± 0.5% were distributed to the superior ICA and 8% ± 0.5% to the OA. The averaged flow velocity at the outlet representing OA was 0.6 m/s ± 0.3 m/s. Calculated from the relative pressure (80 mm Hg) and the pressure drop due to resistance, the perfusion pressure at the end of OA is approximately 49 mm Hg. The flow field at mid-systole is shown in Figure 3. Proximal to the junction, a WSS hotspot is observed that is depicted in Figure 4. It correlates with increased velocities of about 2.5 m/s at the curvature. The maximum WSS is 236 Pa at peak systolic flow. The hotspot occurs only during systole, is highly localized and does not change its position. During diastole, the WSS distribution has a different profile. A comparison with the DICOM images revealed reduced Hounsfield units of 180 at the position of the WSS hotspot, which is in the range of mixed or calcified plaques [12, 26].
Particle tracks Two separate sets of simulations were performed. In the first set, 1500 particles were released from the inlet representing the inferior ICA, thus describing microemboli in the blood flow. The particle tracks for this scenario are shown in Figure 5. All particles that end in the OA are seeded in the near wall flow. In the second set of simulations, the same number of particles was released from the WSS hotspot that
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
474 T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion Velocity 1.2
O6 C6
O1
C7
O2 O3
O4
O2
O1
O5
O3
O4
O6
O5
0.9 0.6
C5
0.3 C4 C4
C2
0 (m s-1) Velocity 2.5
C6
C3 C1
C3
C5
C2
C7
1.9 1.3
C1
0.6 0 (m s-1)
Figure 3 Flow distribution at peak systole in the carotid (lower) and ophthalmic (upper) artery.
Wall shear wall 240 180 120 60 0 (Pa)
Figure 4 Location of wall shear stress hotspot during systole proximal to the ophthalmic artery (l) and possible red-coloured plaque in the curvature based on Hounsfield units (r).
was identified during the flow simulations. Close to the hotspot, a secondary flow pattern was observed, which channels emboli into OA. This is depicted in Figure 6. The overall probability for clots to reach the OA is 7.32 ± 1.08% for microemboli from the inferior ICA and 13.27 ± 6.82% for microemboli released from the WSS hotspot. The high standard deviation of the probability of emboli from the WSS hotspot reaching the OA is due to the fact that at low flow rates (50–70% of average), only 3.43 ± 2.42% of the emboli from the WSS hotspot are washed into OA. This percentage is significantly increased to 18.20 ± 2.83% at higher flow rates. The overall distribution of particles dependent on clot size is presented in Table 1. The average travel time for emboli reaching the OA is 164 ms and 313 ms for emboli in the blood stream or
released from the WSS hotspot, respectively. For emboli washed into the ICA, the travel times are slightly shorter, at 116 ms (blood stream) and 241 ms (WSS hotspot). This is mainly due to the shorter distance between both boundaries.
Discussion Results Computational Fluid Dynamics (CFD) simulations were performed in a realistic RAO patient geometry to analyse the blood flow and clot behaviour at the junction of the OA and ICA. The results show a sufficient perfusion of the
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion 475
Table 1 Particles reaching the ophthalmic artery (OA) dependent on particle diameter. Particle diameter ( μm) 5 25 50 100 200
Figure 5 Emboli seeded in the internal carotid artery reaching the ophthalmic artery.
OA; however, this might be affected by possible stenosis of the ICA prior to the modelled geometry, dependent on its position and grade [1, 3, 4, 7, 24, 27, 29]. The pressure at the end of the OA equals the perfusion pressure of the retinal artery. The simulated pressure of 49 mm Hg is in good agreement with data in the literature, which describe a retinal artery perfusion pressure of 45–55 mm Hg [25]. Velocity 2.0 1.5 1.0 0.5 0 (m s-1)
Figure 6 Upper: Secondary flow pattern at a plane at the WSS hotspot. Lower: Flow path of emboli following main flow (blue) and secondary flow (red).
Particles reaching the OA (%) originating from the OA 8.87 ± 1.50 6.20 ± 0.76 5.90 ± 0.89 8.13 ± 0.97 7.52 ± 1.31
Particles reaching the OA (%) originating from hotspot
14.24 ± 7.28 13.98 ± 7.09 13.79 ± 7.15 12.65 ± 6.64 11.69 ± 5.95
In addition to CFD simulations, the behaviour of emboli in the blood stream was analysed using Lagrangian particle tracking. All particles from the main blood flow that are washed into the OA and might therefore cause RAO were seeded in the near wall flow indicating that particles released from the vessel walls might have a higher potential to cause RAO than particles coming from the heart valve, as already assumed in Mead et al. [23]. The simulations revealed varying WSSs throughout the cycle, with a constant localized hotspot of 240 Pa during systole. Such a WSS hotspot does not correspond to plaque formation, but it might lead to mobilization of existing plaques [5, 30]. Since intravascular plaques may alter flow conditions, this is one reason for RAO that has to be taken into consideration. In addition, secondary flow patterns at the junction of the ICA and OA were observed that channel fluid of the near wall flow of the curvature into the OA, especially at flow rates above 70% of average. Given that mobilization of atherosclerotic plaques is more likely at high flow rates, the effect of plaque mobilization might be even more significant. If a patient shows plaque formation proximal to the curvature, it is therefore possible that such a WSS hotspot might cause RAO by releasing emboli. The occurrence of atherosclerotic plaques in RAO patients proximal to the junction is currently being investigated in a retrospective analysis of angiographic patient data. According to the first results from the numerical simulation of clots behaviour at the ICA/OA junction, one can assume that more than 90% of all clots coming from the ICA are washed into the brain. Astonishingly, the rate of patients having had a stroke previous to RAO is low, at 19% or below [4, 11, 21, 22]. This raises questions about the reasons for RAO. It may happen that small particles do not always cause symptoms if washed into the brain, but may cause RAO if washed into the OA. Also, some of the blood clots causing RAO might not pass the ICA/OA junction, but develop behind it in the OA. Another explanation
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
476 T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion is the aforementioned occurrence of WSS hotspots and pathophysiologies causing secondary flows.
Limitations The model itself is a limitation to this study. Simulations were performed in a vascular system that was averaged over systole and diastole. Possible movement of the vessel wall was neglected. Additionally, no further information is available on the patient condition. It is therefore not possible to separate risk factors in this first feasibility study. A future study will focus on the analysis of several RAO patients, including risk screening and an analysis of parameter correlations. Regarding the particle tracks, only particles with the density of platelets were considered. In future studies, different clot materials, especially plaques, should also be included. The tracks of the particles are very sensitive to small changes in the flow field. Thus, more studies with small variations in the flow parameters and an increased particle number should be performed in order to increase the accuracy of the model. This is especially true for larger particles due to the nature of the Lagrangian particle tracking method.
Conclusions The model presented has the potential to provide insight into the mechanisms behind RAO. It allows us to analyse both general and patient-specific risk factors in a numerical framework, whereby the effects of various parameters may be studied. Further simulations and validating experiments are, however, needed to improve the credibility of this method. Ultimately, results will be compared to an ongoing clinical trial of RAO patients. Conflict of interest statement: The authors declare that there are no conflicts of interests.
References [1] Anderson DC, Kappelle LJ, Eliasziw M, et al. Occurrence of hemispheric and retinal ischemia in atrial fibrillation compared with carotid stenosis. Stroke 2002; 33: 1963–1967. [2] Ballyk PD, Steinman DA, Ethier CR. Simulation of non-Newtonian blood flow in an end-to-side anastomosis. Biorheology 1994; 31: 565–586.
[3] Cheung N, Lim L, Wang JJ, et al. Prevalence and risk factors of retinal arteriolar emboli: the Singapore Malay Eye Study. Am J Ophthalmol 2008; 146: 620–624. [4] Dunlap AB, Kosmorsky GS, Kashyap VS. The fate of patients with retinal artery occlusion and Hollenhorst plaque. J Vasc Surg 2007; 46: 1125–1129. [5] Groen HC, Gijsen FJ, van der Lugt A, et al. Plaque rupture in the carotid artery is localized at the high shear stress region: a case report. Stroke 2007; 38: 2379–2381. [6] Hayreh SS. Prevalent misconceptions about acute retinal vascular occlusive disorders. Prog Retin Eye Res 2005; 24: 493–519. [7] Hayreh SS, Podhajsky P. Ocular neovascularization with retinal vascular occlusion. II. Occurrence in central and branch retinal artery occlusion. Arch Ophthalmol 1982; 100: 1585–1596. [8] Hayreh SS, Zimmerman MB. Central retinal artery occlusion: visual outcome. Am J Ophthalmol 2005; 140: 376–391. [9] Hayreh SS, Zimmerman MB. Fundus changes in central retinal artery occlusion. Retina 2007; 27: 276–289. [10] Hayreh SS, Podhajsky PA, Zimmerman MB. Branch retinal artery occlusion: natural history of visual outcome. Ophthalmology 2009; 116: 1188–1194.e1-4. [11] Hayreh SS, Podhajsky PA, Zimmerman MB. Retinal artery occlusion: associated systemic and ophthalmic abnormalities. Ophthalmology 2009; 116: 1928–1936. [12] Hur J, Kim YJ, Lee HJ, et al. Quantification and characterization of obstructive coronary plaques using 64-slice computed tomography: a comparison with intravascular ultrasound. J Comput Assist Tomogr 2009; 33: 186–192. [13] Kaiser HJ, Schotzau A, Flammer J. Blood-flow velocities in the extraocular vessels in normal volunteers. Am J Ophthalmol 1996; 122: 364–370. [14] Kaufmann TAS, Hormes M, Laumen M, et al. Flow distribution during cardiopulmonary bypass in dependency on the outflow cannula positioning. Artificial Organs 2009; 33: 988–992. [15] Kaufmann TAS, Hormes M, Laumen M, et al. The impact of aortic/subclavian outflow cannulation for cardiopulmonary bypass and cardiac support: A computational fluid dynamics study. Artificial Organs 2009; 33: 727–732. [16] Kaufmann TAS, Schmitz-Rode T, Steinseifer U. Implementation of cerebral autoregulation into computational fluid dynamics studies of cardiopulmonary bypass. Artif Organs 2012; 36: 754–758. [17] Kaufmann TAS, Wong KC, Schmitz-Rode T, Steinseifer U. Mimicking of cerebral autoregulation by flow-dependent cerebrovascular resistance: a feasibility study. Artif Organs 2012; 36: E97–E101. [18] Klein R, Klein BE, Jensen SC, et al. Retinal emboli and stroke: the Beaver Dam Eye Study. Arch Ophthalmol 1999; 117: 1063–1068. [19] Klein R, Klein BE, Moss SE, Meuer SM. Retinal emboli and cardiovascular disease: the Beaver Dam Eye Study. Arch Ophthalmol 2003; 121: 1446–1451. [20] Laumen M, Kaufmann T, Timms D, et al. Flow analysis of ventricular assist device inflow and outflow cannula positioning using a naturally shaped ventricle and aortic branch. Artificial Organs 2010; 34: 798–806. [21] Leisser C. Hyperlipidaemia and increased systemic blood pressure – is there a connection to retinal artery occlusion? Spektrum Augenheilkd 2010; 24: 154–156.
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
T.A.S. Kaufmann et al.: Analysis of retinal artery occlusion 477 [22] Leisser C. Risk factors of retinal artery occlusion. Spektrum Augenheilkd 2012; 26: 202–204. [23] Mead GE, Lewis SC, Wardlaw JM, Dennis MS. Comparison of risk factors in patients with transient and prolonged eye and brain ischemic syndromes. Stroke 2002; 33: 2383–2390. [24] Mizener JB, Podhajsky P, Hayreh SS. Ocular ischemic syndrome. Ophthalmology 1997; 104: 859–864. [25] Pillunat LE, Anderson DR, Knighton RW, Joos KM, Feuer WJ. Autoregulation of human optic nerve head circulation in response to increased intraocular pressure. Exp Eye Res 1997; 64: 737–744. [26] Saba L, Mallarin G. Window settings for the study of calcified carotid plaques with multidetector CT angiography. AJNR Am J Neuroradiol 2009; 30: 1445–1450. [27] Schmidt D, Hetzel A, Geibel-Zehender A, Schulte-Mönting J. Systemic diseases in non-inflammatory branch and central retinal artery occlusion–an overview of 416 patients. Eur J Med Res 2007; 12: 595–603. [28] Schumacher M, Schmidt D, Jurklies B, et al. Central retinal artery occlusion: local intra-arterial fibrinolysis versus conservative treatment, a multicenter randomized trial. Ophthalmology 2010; 117: 1367–1375.e1.
[29] Sharma S, Pater JL, Lam M, Cruess AF. Can different types of retinal emboli be reliably differentiated from one another? An inter- and intraobserver agreement study. Can J Ophthalmol 1998; 33: 144–148. [30] Tang D, Teng Z, Canton G, et al. Local critical stress correlates better than global maximum stress with plaque morphological features linked to atherosclerotic plaque vulnerability: an in vivo multi-patient study. Biomed Eng Online, 2009; 8: 15. [31] von Graefe A. Ueber Embolie der Arteriacentralisretinae als Ursache plotzlicher Erblindung. Albrecht von Graefes ArchOphthalmol 1859; 5: 136–157. [32] Weger M, Renner W, Pinter O, et al. Role of factor V Leiden and prothrombin 20210A in patients with retinal artery occlusion. Eye 2003; 17: 731–734. [33] Wong TY, Klein R. Retinal arteriolar emboli: epidemiology and risk of stroke. Curr Opin Ophthalmol 2002; 13: 142–146. [34] Wong TY, Larsen EK, Klein R, et al. Cardiovascular risk factors for retinal vein occlusion and arteriolar emboli: the Atherosclerosis Risk in Communities and Cardiovascular Health studies. Ophthalmology 2005; 112: 540–547.
Brought to you by | New York University Bobst Library Technical Services Authenticated Download Date | 5/31/15 8:42 AM
