Abstract
Purpose:
To identify and rank the lamina cribrosa (LC) morphologic factors that influence LC microcapillary hemodynamics and oxygen concentrations using computational fluid dynamics (CFD).
Methods:
We generated 12,000 ‘artificial' LC microcapillary networks and predicted blood flow velocities and oxygen concentrations within the microcapillaries using CFD. Across models, we varied the average pore size of the LC (5500 ± 2400 μm2), the microcapillary arrangement (radial, isotropic, or circumferential), the LC diameter (1.9 ± 0.3 mm), the inferior–superior curvature (340 ± 116 m−1), and the nasal–temporal curvature (−78 ± 130 m−1). We assumed that blood flow originated from the Circle of Zinn-Haller, fed the LC uniformly at its periphery, and was drained into the central retinal vein. Arterial (50 ± 6 mm Hg) and venous (17.7 ± 6 mm Hg) pressures were applied as boundary conditions and were also varied within our simulations. Finally, we performed linear regression analysis to rank the influence of each factor on LC hemodynamics and oxygen concentrations.
Results:
The factors influencing LC hemodynamics and oxygen concentrations the most were: LC diameter, arterial pressure, and venous pressure, and to a lesser extent: the microcapillary arrangement (anisotropy) and nasal–temporal curvature. Lamina cribrosa pore size and superior–inferior curvature had almost no impact. Specifically, we found that LCs with a smaller diameter, a radial arrangement of the microcapillaries, an elevated arterial pressure and a decreased venous pressure had higher oxygen concentrations across their networks.
Conclusion:
This study described LC hemodynamics using a computational modeling approach. Our study may provide clinically relevant information for the understanding of ischemia-induced neuronal cell death in optic neuropathies.
Open angle glaucoma (OAG) is an optic neuropathy disease, characterized by a progressive loss of retinal ganglion cells (RGC), and permanent vision loss. A well-known risk factor for OAG is an elevated IOP and it is currently considered the only treatable factor. However, the link between the pathogenesis of OAG and elevated IOP is not firmly established.
1 Up to one-third of glaucoma patients do not have elevated IOP
2 and a significant number of people with elevated IOP never develop glaucoma.
3
Intraocular pressure is known to deform the posterior structures of the eyes, particularly the structurally weak lamina cribrosa (LC)
4,5—a major site of RGC damage in glaucoma. The LC is a sieve-like network of collagenous connective tissue beams located in the optic nerve head (ONH) containing a comprehensive network of capillaries that provides nutritional and oxygen support to the RGCs.
6 The LC undergoes morphologic changes (e.g., ‘cupping') and significant remodeling during OAG development
7–10 and these changes have become the focus of many studies that have attempted to link changes in IOP to RGC death. A prominent hypothesis suggests that IOP-induced ONH strains damage the RGCs either directly
11 or indirectly through a mechanical disruption of axonal transport in RGCs, thus depriving RGCs from essential trophic factors.
12–14 An alternative hypothesis suggests that ONH strains can alter hemodynamics within the LC, which will in turn reduce the diffusion of nutrients to astrocytes and/or induce ischemia, thus resulting in RGC death.
To date, there is evidence to support the latter vascular hypothesis. Several studies have reported a correlation between poor ocular hemodynamics and glaucoma
15–17; for instance, Kanakamedal et al.
18 have suggested that people of African descent had a much higher prevalence of OAG compared with those of European descent due to their relatively weaker vascular components. The study also revealed a significant correlation between changes in blood flow and ONH morphology. Unfortunately, accurate measurement of blood circulation in vivo, especially at the level of the LC, has remained extremely challenging. This represents a barrier to establish impaired ocular hemodynamics as a contributing cause of OAG.
In view of the technical limitations, computational fluid dynamics (CFD) is a powerful alternative to understand how deep ocular blood flow could be involved in the development and progression of OAG. A CFD study by Carichino et al.
19 revealed that IOP-induced LC deformation could reduce the lumen size of the central retinal artery (CRA), and thus the overall ocular blood flow volume. A computational study by Causin et al.
20 also suggested that an increase in IOP could disrupt LC hemodynamics through biomechanical interactions. However, no CFD studies have yet modeled the LC as a detailed capillary network and studied how changes in ONH morphology could affect hemodynamics and oxygen concentrations within such a capillary network. We believe it is of critical importance to understand how ONH morphologic changes impact the LC blood flow at the microscale, because this could be a major reason for RGC death in OAG pathophysiology.
The goal of this study was to provide a modeling study of LC hemodynamics and oxygen concentration at the microscale. Specifically, we aimed to identify the main LC morphologic parameters that could affect LC blood flow and oxygen distribution. Because a majority of these morphologic parameters can now be measured in vivo with optical coherence tomography (OCT), our study may be relevant for the clinical management of glaucoma.
Methods
Our General Approach to Study LC Hemodynamics
In this study, our aim was to identify and rank the morphologic factors influencing the hemodynamics and oxygen concentrations in the microcapillaries of the LC. To this end, we generated 12,000 artificially-reconstructed LC microcapillary networks with varied morphologic characteristics. A CFD approach was then used to predict blood flow velocities and oxygen concentrations in all LC networks.
Assumptions on LC Morphology and Microcapillary Blood Flow
For all our models, several general assumptions were made: (1) The LC is a curved surface with principal curvatures aligned along the inferior–superior and nasal–temporal directions,
21 (2) each LC collagenous beam contains a microcapillary in its center that supplies nutrients and oxygen to adjacent axons,
22 (3) the LC beams, and thus the microcapillaries, are parallel to the LC surface, (4) the LC microcapillaries have a uniform diameter of 8 μm,
23 (5) oxygenated blood, originating from the Circle of Zinn-Haller or the short posterior arteries, feeds the LC uniformly at its periphery,
24 and (6) venous drainage of the LC microcapillaries occurs through the central retinal vein that passes through the center of the LC.
25
Generation of Artificial LC Microcapillary Networks
Under our general assumptions, we used a custom-written algorithm (Matlab 2014b; The Mathworks, Inc., Natick, MA, USA) to generate artificial LC microcapillary networks that were representative of a healthy (i.e., nonglaucomatous) human LC. Across networks, we varied the pore size of the LC, the microcapillary arrangement, the LC diameter, and the main LC curvatures (
Table 1).
Table 1 Input Parameters Used for the CFD Simulations
Table 1 Input Parameters Used for the CFD Simulations
Table 2 Table of Constants Used in the Governing Equations
Table 2 Table of Constants Used in the Governing Equations
LC Network Generator Algorithm.
Given a set of control parameters (defined in subsequent sections), we numerically generated multiple ellipses in an enclosed two-dimensional (2D) space (
Fig. 1a). Each new ellipse was discarded if it overlapped with its neighbors, but kept if it did not. We iterated the process until the space was densely packed and additional ellipses could not be fitted into the space.
The negative ‘ellipse' image (white region in
Fig. 1a) was then ‘skeletonized', a process whereby a thin line (1-pixel thick) was generated within the white region, equidistant to the region's boundaries. A ‘pruning' algorithm was then used to remove isolated branches to ensure that all LC capillaries were connected except at the periphery where blood flowed in. The central retinal artery and vein (off-centered by 0.05 mm toward the nasal direction) were added to complete each LC capillary network.
26 The size of the artery (0.02 mm
2) and that of the vein (0.03 mm
2) were fixed for all networks
26 independently of the total LC area.
For each LC network, our algorithm received four input parameters (LC pore size, anisotropy, diameter, superior–inferior curvature, and nasal–temporal curvature), each of which was a normally distributed random number.
The final LC capillary image was then converted to a vector file (DXF format) using ImageJ
27 and Matlab in order to generate nodal and connectivity information (i.e., mesh) of each LC capillary network, as required for CFD processing.
LC Pore Size.
Generally, the LC exhibits larger pores and less connective tissue in the superior and inferior regions than its temporal and nasal regions
28–30 (
Table 1). Furthermore, the pores are smaller and more uniform in the central area.
7 To mimic such heterogeneity, all our networks were subdivided into nine regions (
Fig. 1c): the central region (inner 20% of the LC diameter), four regions (superior [S], inferior [I], nasal [N], and temporal [T],) in the middle area (defined between 20% and 60% of the LC diameter), and four additional regions (S, I, N and T) in the peripheral area (defined between 60% and 100% of the LC diameter).
For each LC network, a pore size ‘seed' value δ was selected randomly from a normal distribution with a mean of 0 and a SD of 5.4 μm. From such a value, the minor axis of each ellipse in the central region was randomly chosen between δ + 30 μm and δ + 32 μm; in the N-T peripheral regions between δ + 45 μm and δ + 50 μm; and in the S-I peripheral regions between δ + 58 μm and δ + 81 μm. This method allowed us to generate heterogeneously distributed LC pores with a single control parameter. For instance,
Figure 1b was generated with a δ value of 8 μm and
Figure 1d with a value of 0 mm.
LC Anisotropy.
In its periphery, the LC was either considered: (1) anisotropic with a radial arrangement of the LC beams (and thus capillaries) as observed in humans and primates,
31 (2) isotropic (also observed in humans),
31 or (3) anisotropic with a circumferential arrangement (hypothetical scenario). For our LC network generator algorithm, the anisotropy input parameter (σ) was defined as the ellipse major-to-minor axis ratio minus 1. For each network, such a parameter was chosen randomly from a normal distribution (mean: 0; SD: 0.05). Here, σ less than 0 indicates anisotropy with a circumferential arrangement, σ equal to zero indicates isotropy, and σ greater than 0 indicates anisotropy with a radial arrangement. Note that a higher absolute value corresponds to more elongated ellipses, and thus a higher anisotropy degree. The central and middle regions of the LC were considered isotropic for all LC networks.
LC Diameter.
For each network, the LC diameter was used as an input to define the height and the width of the original 2D space (before ellipses packing). This method ensured that LC pore size was independent of the LC diameter. Once ellipses were generated, the LC diameter was further used to crop the LC to a circular shape. We used the following normal distribution of healthy LCs' diameters: 1.9 ± 0.3 mm (
Table 1).
LC Curvature.
For each LC network, we varied its N-T and S-I curvatures K1 and K2 by ‘morphing' the initially flat geometry onto a hyperbolic three-dimensional (3D) plane (of principal curvatures K1 and K2). K1 and K2 followed normal distributions, K1: 340 ± 116 m
−1 and K2: −78 ± 130 m
−1 (
Table 1), as recently measured by our group in a healthy population.
21 From K1 and K2, and for each LC network, we computed the corresponding global shape index (GSI) – a quantity describing the global shape of the LC and that varied between −1 (cup shape) and +1 (cap shape) (
Fig. 3c).
21
Hemodynamics Boundary Conditions
We applied arterial inlet pressure (50 ± 6 mm Hg) and venous (17.7 ± 6 mm Hg)
32 outlet pressure as boundary conditions to each LC network and each value was randomly selected from normal distributions (
Table 1). Because there was no reported data on the arterial inlet pressure of the LC microcapillary network, we assumed the mean arterial inlet pressure to be 50 mm Hg, which would give an average blood flow velocity of approximately 1.0 mm/s
33—a reasonable assumption for microcapillaries. We also assumed arterial pressure SD to be the same as that of the venous pressure.
Governing Equations for LC Blood Flow Simulation
Blood Flow Equations.
As a first approximation, we assumed that blood was an incompressible Newtonian fluid,
34 whose flow behavior can be described with the Hagen-Poiseuille equation:
Here, η (Pa·s) the blood viscosity, Q the volumetric flow rate (m
3/s), Δp is the pressure drop along the vessel (Pa) and L is the vessel length (m) and R
c is the capillary radius.
Oxygen mass transport inside the microcapillary was solved with the following equations:
where f is the rate of convective oxygen transport within the microcapillary segments (mol/s), C is oxygen concentration in plasma (mol/m
3), Q is the volume flow rate (m
3/s), C
HB is the concentration of hemoglobin-bound oxygen in a fully saturated red blood cell and H is the hematocrit. We assumed H and C
HB to be constant for all capillaries segment. SO
2(C) represents the oxygen-hemoglobin saturation function approximated by the Adair equation
35:
where α is the oxygen solubility coefficient in blood (mol/m
3/mm Hg), α
1 = 0.01524, α
2 = 2.7 × 10
−6, α
3 = 0, and α
4 = 2.7 × 10
−6.
35
Oxygen Diffusion and Consumption Equations.
To predict oxygen diffusion and consumption inside the nervous tissues, we followed the standard Krogh cylinder-type model
36 (
Fig. 2), in which a tissue cylinder is assumed to receive oxygen supply from a single capillary that lies in its center. The following governing equation (Fick's law of diffusion
36) was imposed to the nervous tissues:
where D is the diffusion coefficient of oxygen of the nervous tissues (m
2/s), C is the oxygen concentration inside the nervous tissues (mol/m
3), M(C) is the oxygen consumption rate of the nervous tissues (mol/m
3/s). We assumed that the nervous tissues were homogeneous with the same oxygen diffusion coefficient and consumption regardless of location. M(C) can be approximated by Michaelis-Menten kinetics—the most common model used to model tissue oxygen consumption rate
36:
where M
0 is tissue-dependent maximum oxygen consumption rate and C
MM is the Michaelis-Menten constant corresponding to the oxygen concentration where the consumption rate is 50% of M
0. We assumed M
0 to be uniform throughout the nervous tissues.
We further assumed that the transfer of oxygen from the LC microcapillaries to the nervous tissues occurred at the capillary wall (R
c) as governed by the following equation
37:
where r is the radial distance (m) from the vessel centerline, R
c is the capillary radius (m), and q is the oxygen flux across the vessel wall (mol/m/s). The term
represents the oxygen flux at the capillary wall, which is a function of the average oxygen concentration inside the capillary.
In our study, we simplified the 3D Krogh model to a one-dimensional (1D) line model to reduce the computational time to process a large number of complex vascular networks.
Equation 1 was first used to solve for blood flow rate inside all capillary segments and
Equation 2 was used to approximate plasma oxygen concentration in the segments. The oxygen consumption–diffusion equation (
Eq. 4) was applied to the nervous tissue, and by conservation of mass (
Eq. 6), we determined the flux of oxygen from the capillaries to the nervous tissues and quantified the oxygen concentration drop along the axial length of each capillary (refer to
Table 2 for list of contants used in the equations). Only the average oxygen concentration and average flow velocity were computed in the capillary. These values were representative of the oxygen concentration in the surrounding nervous tissue (
Supplementary Material).
Data Processing and Statistical Analysis
A total of 12,000 simulations were performed in this study. The simulation outputs consisted of three parameters—the oxygen concentration, the blood pressure, and the blood flow velocity inside the microcapillaries. We also identified capillary segments of each LC that had oxygen concentrations of less than 10 mm Hg O
2 as hypoxic regions.
38 Because no study has been conducted on the metabolism of LC's nervous tissue, this value was assumed from the study of oxygen-intensive human tumors and it should serve as an upper threshold for normal cells. We performed multiple linear regression on a normalized mean value of oxygen concentration and flow velocity for each LC network, and used the regression coefficient (β) to rank the importance of each LC characteristic (
Table 1) in the LC hemodynamics.
Verification of Generated LC Networks
Each generated LC network was re-evaluated using a digital image processing algorithm to obtain each individual pore size in μm
2 and its regional distribution. A Fourier transform algorithm
39 was used to compute the principal alignment and the degree of anisotropy of the microcapillaries in each local region.
Results
Verification of Generated LC Networks
We analyzed the various study parameters of our LC networks and found the following. First, the pore size distribution (
Fig. 3a) of our models closely matched that of normal human eyes, with S-I regions having significantly larger pores than N-T regions (pore size by region: Nasal: 2200 ± 1061 μm
2, Temporal: 2351 ± 950 μm
2, Superior: 5093 ± 1750 μm
2, Inferior: 5100 ± 1765 μm
2). Next, for anisotropy, our LC networks had normally distributed LC anisotropy values that described a range of anisotropic scenarios from radial to tangential capillary arrangements (
Fig. 3b). Finally, the distribution of the Global Shape Index (GSI) values (a measure of LC curvature) of our models had an average value of −0.38, corresponding to a saddle-rut shape and similar to that reported in human eyes (
Fig. 3c).
Regional Variations in Oxygen Concentration and Flow Velocity
The average flow velocity (m/s) and average oxygen concentration (mm Hg O
2) across the LC network are shown in
Figure 4d. Notably, the average oxygen concentration decreased significantly (
P < 0.001) from the peripheral region (64.8 mm Hg O
2) toward the mid (35.8 mm Hg O
2) and central regions (25.5 mm Hg O
2). Variations in oxygen concentration between N-T-S-I regions were not significant. The overall oxygen concentration across all LC networks was 50.4 ± 9.5 mm Hg O
2.
Flow velocity increased significantly (P < 0.001) from the peripheral region (0.7 × 10−3 m×/s) toward the mid (1.2 × 10−3 m/s) and central region (2.9 × 10−3 m/s). No significant variations were observed between N-T-S-I regions and the LC networks had an average flow velocity of 1.2 × 10−3 ± 0.3 × 10−3 m/s.
Influence of the LC Parameters on Oxygen Concentration, Flow Velocity, and Hypoxic Area
The association between LC parameters and the resulting LC blood oxygen concentration and flow velocity are plotted in
Figures 5 and
6, respectively. In addition, we also quantified a third characteristic—hypoxic regions (<10 mm Hg in oxygen) as a proportion of total LC microcapillary area—and looked at its association with all LC parameters (
Supplementary Figure S1).
We ranked the influence of each LC parameter on the three characteristics based on the magnitude of the coefficient (β) obtained from multiple linear regression analysis (a positive coefficient suggests a positive correlation between the LC parameter and the hemodynamics characteristics and vice versa). First, the factors influencing average oxygen concentration the most were, in order of importance, LC diameter (β = +0.72), arterial pressure (β = +0.58), venous pressure (β = −0.37), anisotropy (β = +0.14), K1 (β = −0.07), pore size (β = −0.01), and K2 (β = +0.002). Second, for average flow velocity, the factors were ranked as follows, arterial pressure (β = +0.69), diameter (β = −0.57), venous pressure (β = −0.44), pore size (β = +0.82), K1 (β = −0.06), anisotropy (β = +0.047), and K2 (β = +0.006). Finally, the ranking of factors associated with the proportion of hypoxic areas was similar to that of average oxygen distribution. The rankings based on regression results are summarized in
Figure 7.
Discussion
In this study, we developed a novel method to generate LC microcapillary networks that were incorporated into a CFD model for blood flow analysis. We were able to finely control the important LC morphologic parameters (pore size, anisotropy, diameter, and curvature) and study their effects on LC hemodynamics and oxygen concentration. Our models predicted that LC hemodynamics was highly influenced by the LC diameter, the arterial pressure, and the venous pressure. Notably, LCs with large diameters (>2.0 mm) had significantly poorer oxygen concentrations. Other parameters, such as pore size, anisotropy, and curvature had much less influence on LC hemodynamics.
We found that the regional variations in both the flow velocity and oxygen concentration were significant along the radial direction, but not in the circumferential N-T-S-I direction (
Fig. 4d). On average, oxygen concentrations decreased significantly from 64 mm Hg O
2 in the peripheral region to 25 mm Hg O
2 in the central region. While there is no literature to date that measured human ONH oxygen concentrations in vivo, our range of values is consistent with that observed in mice and porcine models.
40,41 The absolute velocity values obtained from our simulations also expand our knowledge of ocular blood flow in the human LC, because conventional techniques (e.g., laser speckle flowgraphy) yield only relative velocities and cannot provide deep flow measurements within the LC.
42
Lamina cribosa diameter was found to be the most influential parameter affecting both oxygen concentration and hemodynamics. Lamina cribosas with large diameters were found to exhibit a lower average oxygen concentration, a lower average flow velocity and a higher percentage of hypoxic area (up to 80% for LCs larger than 2.0 mm). Previous studies have reported an increase in optic disc diameter in glaucoma patients.
43,44 Combined with our findings, this suggests a plausible mechanism for RGC axonal death in glaucoma: elevated IOP could result in an increase in LC diameter (acutely and chronically),
45–47 which may in turn result in lower blood flow velocities and oxygen concentrations within the LC, starving RGC axons of oxygen, and possibly other nutrients. Further research needs to be done on the role of optic disc enlargement as a triggering cause for glaucoma. Our finding may also explain why certain races have a greater predisposition to OAG.
18 For instance, people of African descent have been shown to exhibit significantly larger optic discs
48 and poorer vascular components compared with those of European descent.
18
Not surprisingly, arterial and venous blood pressures were the other two most significant parameters that affected oxygen concentration and flow velocity. A decrease in arterial pressure and an increase in venous pressure both resulted in lower average oxygen concentration and flow velocity. This result is reasonable, as the aforementioned changes in pressure will reduce the overall perfusion pressure of the LC. It has also been shown that both an increase in central venous pressure and a decrease in arterial pressure are risk factors for the development of OAG,
49–52 and our result suggests that changes in blood pressure can significantly impact LC hemodynamics and potentially contribute to OAG progression.
We found LC pore size to have a significantly smaller influence on LC hemodynamics when compared with other factors. Changes in pore size had almost no effect on the oxygen concentration, while larger pores slightly increased the average blood flow velocity. Several monkey model studies have shown that there are no significant differences in pore size between normal and glaucomatous eyes.
53,54 A histologic study by Roberts et al.
8 found an increase in connective tissue beams in experimental glaucoma eyes, but the connective tissue volume fraction (analogous to pore area) remained relatively constant. Our findings further suggest that changes in collagenous beams and the LC microcapillary pore size may not significantly impact LC hemodynamics.
Anisotropy of the microcapillary network was shown to have greater impact on average oxygen concentration than pore size. Lamina cribosas with radially oriented capillaries at the peripheral region had higher average oxygen concentration than those with isotropic or circumferential configurations. Note that human LCs naturally exhibit a radial arrangement of the microcapillaries at the periphery
31 and any deviations from such an arrangement (possible through focal defects
55) may contribute to the development of OAG.
We found that an increase in curvature along the N-T direction (K1) resulted in both lower average oxygen concentration and blood flow velocity, whereas changes in S-I curvature (K2) had almost no impact on the two output parameters. The increase in N-T curvature corresponds to an increase in LC ‘cupping' and LC depth,
21 which are well-known to be associated with glaucoma development.
56,57 It was notable in our study that K1 was significantly more influential than K2. This difference may arise from the heterogeneity in LC capillary density. The N-T region has a greater density of microcapillaries than the S-I region (due to smaller pore area). This may mean that changes to curvature along the N-T direction have a greater influence on microcapillary morphology and hemodynamics.
It should be mentioned that a computational study on LC hemodynamics was previously conducted by Causin et al.
20 The authors used a poroelastic model, whereby blood capillaries were modeled as homogeneous pores in an elastic matrix. The work accounted for the biomechanical action of IOP, retrolaminar tissue pressure, and scleral tension on hemodynamics, but it did not take into account complex capillary networks and oxygen transport. On the other hand, our work took those latter into account, but ignored the effects of mechanical stress. We believe both works should be regarded as complementary and contribute to the current understanding of LC hemodynamics.
The followings are the main limitations of our study. First, we modelled the LC as a surface. This is not biologically accurate as the LC has a thickness of approximately 270 μm.
58 Considering that the LC may exhibit highly heterogeneous stresses (induced by IOP, the retrolaminar pressure, or during eye movements
59), it would be logical to observe variations in LC hemodynamics and oxygen concentrations through the thickness of the LC, as predicted by Causin et al.
20 LC thickness has also been shown to change with glaucoma progression, which may in turn compromise blood flow further in glaucoma subjects.
8 While using ‘surface-like' LCs (as performed herein) should provide a first degree of understanding of LC hemodynamics, future improvements to our models should prioritize on incorporating LC thickness and interactions with LC stresses.
Another limitation is that we simplified the study by reducing our 3D models into 1D ‘line' models in which only the average blood flow velocity and oxygen concentration were solved for each cross section of the capillary. According to our 3D preliminary work (
Supplementary Material), the average oxygen concentration inside the capillary segment was found representative of the surrounding neural tissue's oxygenation state, because the oxygen concentration dropped at most by 10% from the capillary to the adjacent LC nervous tissue. However, this assumption undermined the effect of very large pores in our findings as there could be a significant drop in oxygen concentration in the middle of the pore that was not taken into account in our models. Further models, which include nervous tissue between the LC pores, should be built, to better understand the effect of pore size.
We also imposed several constraints on the CFD simulation. The constraints associated with blood flow equations and oxygen diffusion in tissue are mostly justified based on similar computational studies conducted on capillary blood flow.
35,36 However, the blood viscosity, especially in microcirculation, may not be constant as we assumed. Blood becomes more particulate in small capillaries and its bulk viscosity can be a function of hematocrit, flow velocity, and capillary diameter.
60 One could try to account for this particulate nature by increasing the bulk viscosity of blood by 30% as proposed by several studies.
61,62 Interestingly, we found that increasing the viscosity did not considerably change the linear regression coefficients of all parameters (coefficients values changed by 5% at most when the viscosity was increased by 30% to 5.2 × 10
−3 Pa/s) and the ranking of factors remained the same. Changes in average flow velocity and oxygen concentration were observed but remained relatively small (20% and 12%, respectively). However, increasing viscosity may not necessarily make our models more biofidelic as blood viscosity in the microcapillary is determined by many factors, including the hematocrit. Furthermore, variations in viscosity could also happen locally in each capillary segment.
In addition, our assumption of constant capillary diameter may not be biologically accurate. A computational study by Causin et al.
20 also suggested that an increases in IOP could induce changes in the capillaries' diameter, which could compromise the lamina's perfusion. Autoregulation process can dynamically change the LC capillary diameter and autoregulation impairment at the ONH has been considered as a possible cause for glaucoma.
63,64 Unfortunately, no studies have yet provided a comprehensive map of capillary size within the LC. It is plausible that capillary radius could change regionally, and this may be further influenced by autoregulation, and IOP-induced LC stresses (these latter being highly heterogeneous). Note that all our other parameters were derived from experimental and clinical observations. Due to a lack of information, we have considered not to vary capillary radius in our models. Further studies are necessary to quantify the size of capillaries in the human LC. Such information could considerably enhance our models.
Our findings suggest that the central region of the LC may be more susceptible to ischemia than other regions. This is in contrast with other studies, which indicate that glaucomatous axonal loss starts from the peripheral region.
65 Also, an hour-glass pattern of axonal loss, in which S-I nerve fibers are selectively more damaged than the N-T nerve fibers,
65 was not evident from this study. These discrepancies may arise from the inhomogeneity of both the microcapillary characteristics (diameter, autoregulation capabilities, breakage, etc.) and nerve fiber density in each region of the LC. It is also possible that mechanical insult plays a greater role than hemodynamics in creating these damage patterns. Further studies that include the regional variations of both nerve fibers and microcapillaries may reconcile the observed axonal loss pattern and our hemodynamics result.
Lastly, accurate in vivo measurement of oxygen concentration and blood flow velocity are not yet possible with current techniques. This makes the validation of our model difficult. Nevertheless, our output lies within normal biological values in capillary networks of comparable size.
In summary, our study revealed, for the first time, the relationship between clinically relevant LC parameters and the hemodynamics and oxygen distribution of the LC. The clear ranking of each parameter's importance presented in the study has the potential to be clinically useful in glaucoma management, diagnosis, and risk profiling. Finally, our approach could serve as a numerical testbed for hypotheses related to LC hemodynamics.
Acknowledgements
Supported by grants from the Ministry of Education, Academic Research Funds, Tier 1 (MJAG; R-397-000-181-112; Singapore) and from an NUS Young Investigator Award (MJAG; NUSYIA_FY13_P03, R-397-000-174-133; Singapore).
Disclosure: T. Chuangsuwanich, None; K.E. Birgersson, None; A. Thiery, None; S.G. Thakku, None; H.L. Leo, None; M.J.A. Girard, None
References
Rieck
J.
The pathogenesis of glaucoma in the interplay with the immune system: a summary of the current knowledge.
Invest Ophthalmol Vis Sci.
2013;
54:
2393–2409.
Hendrickx
KH,
van den Enden
A,
Rasker
MT,
Hoyng
PF.
Cumulative incidence of patients with disc hemorrhages in glaucoma and the effect of therapy.
Ophthalmology.
1994;
101:
1165–1172.
Racette
L,
Wilson
MR,
Zangwill
LM,
Weinreb
RN,
Sample
PA.
Primary open-angle glaucoma in blacks: a review.
Surv Ophthalmol.
2003;
48:
295–313.
Weinreb
RN,
Aung
T,
Medeiros
FA.
The pathophysiology and treatment of glaucoma: a review.
JAMA.
2014;
311:
1901–1911.
Girard
MJ,
Beotra
MR,
Chin
KS,
et al.
In vivo 3-dimensional strain mapping of the optic nerve head following intraocular pressure lowering by trabeculectomy.
Ophthalmology.
2016;
123:
1190–1200.
Burgoyne
CF.
A biomechanical paradigm for axonal insult within the optic nerve head in aging and glaucoma.
Exp Eye Res.
2011;
93:
120–132.
Brown
DJ,
Morishige
N,
Neekhra
A,
Minckler
DS,
Jester
JV.
Application of second harmonic imaging microscopy to assess structural changes in optic nerve head structure ex vivo.
J Biomed Optics.
2007;
12:
024029.
Roberts
MD,
Grau
V,
Grimm
J,
et al.
Remodeling of the connective tissue microarchitecture of the lamina cribrosa in early experimental glaucoma.
Invest Ophthalmol Vis Sci.
2009;
50:
681–690.
Shoji
T,
Kuroda
H,
Suzuki
M,
et al.
Correlation between lamina cribrosa tilt angles, myopia and glaucoma using OCT with a wide bandwidth femtosecond mode-locked laser.
PLoS One.
2014;
9:
e116305.
Kimura
Y,
Akagi
T,
Hangai
M,
et al.
Lamina cribrosa defects and optic disc morphology in primary open angle glaucoma with high myopia.
PLoS One.
2014;
9:
e115313.
Strouthidis
NG,
Girard
MJ.
Altering the way the optic nerve head responds to intraocular pressure—a potential approach to glaucoma therapy.
Curr Opin Pharmacol.
2013;
13:
83–89.
Fechtner
RD,
Weinreb
RN.
Mechanisms of optic nerve damage in primary open angle glaucoma.
Surv Ophthalmol.
1994;
39:
23–42.
Burgoyne
CF,
Downs
JC,
Bellezza
AJ,
Suh
J-KF,
Hart
RT.
The optic nerve head as a biomechanical structure: a new paradigm for understanding the role of IOP-related stress and strain in the pathophysiology of glaucomatous optic nerve head damage.
Prog Retin Eye Res.
2005;
24:
39–73.
Quigley
HA,
McKinnon
SJ,
Zack
DJ,
et al.
Retrograde axonal transport of BDNF in retinal ganglion cells is blocked by acute IOP elevation in rats.
Invest Ophthalmol Vis Sci.
2000;
41:
3460–3466.
Harris
A,
Kagemann
L,
Ehrlich
R,
Rospigliosi
C,
Moore
D,
Siesky
B.
Measuring and interpreting ocular blood flow and metabolism in glaucoma.
Can J Ophthalmol.
2008;
43:
328–336.
Moore
D,
Harris
A,
Wudunn
D,
Kheradiya
N,
Siesky
B.
Dysfunctional regulation of ocular blood flow: a risk factor for glaucoma?
Clin Ophthalmol.
2008;
2:
849–861.
Galassi
F,
Giambene
B,
Varriale
R.
Systemic vascular dysregulation and retrobulbar hemodynamics in normal-tension glaucoma.
Invest Ophthalmol Vis Sci.
2011;
52:
4467–4471.
Kanakamedala
P,
Harris
A,
Siesky
B,
et al.
Optic nerve head morphology in glaucoma patients of African descent is strongly correlated to retinal blood flow.
Br J Ophthalmol.
2014;
98:
1551–1554.
Carichino
L,
Guidoboni
G,
Arieli
Y,
Siesky
BA,
Harris
A.
Effect of lamina cribrosa deformation on the hemodynamics of the central retinal artery: a mathematical model.
Invest Ophthalmol Vis Sci.
2012;
53:
2836–2836.
Causin
P,
Guidoboni
G,
Harris
A,
Prada
D,
Sacco
R,
Terragni
S.
A poroelastic model for the perfusion of the lamina cribrosa in the optic nerve head.
Math Biosci.
2014;
257:
33–41.
Thakku
SG,
Tham
Y-C,
Baskaran
M,
et al.
A global shape index to characterize anterior lamina cribrosa morphology and its determinants in healthy indian eyes LC-GSI to characterize LC morphology.
Invest Ophthalmol Vis Sci.
2015;
56:
3604–3614.
Schacknow
PN,
Samples
JR.
The Glaucoma Book: A Practical, Evidence-Based Approach to Patient Care.
New York:
Springer Science & Business Media;
2010.
Freitas
RA.
Nanomedicine Volume I: Basic Capabilities.
Georgetown:
Landes Bioscience;
1999.
Flammer
J,
Orgül
S.
Optic nerve blood-flow abnormalities in glaucoma.
Prog Retin Eye Res.
1998;
17:
267–289.
Flammer
J,
Orgül
S,
Costa
VP,
et al.
The impact of ocular blood flow in glaucoma.
Prog Retin Eye Res.
2002;
21:
359–393.
Jonas
J,
Mardin
CY,
Schlötzer-Schrehardt
U,
Naumann
G.
Morphometry of the human lamina cribrosa surface.
Invest Ophthalmol Vis Sci.
1991;
32:
401–405.
Schneider
CA,
Rasband
WS,
Eliceiri
KW.
NIH image to ImageJ: 25 years of image analysis.
Nat Methods.
2012;
9:
671–675.
Quigley
HA,
Addicks
EM,
Green
W,
Maumenee
AE.
Optic nerve damage in human glaucoma: II. The site of injury and susceptibility to damage.
Arch Ophthalmol.
1981;
99:
635–649.
Quigley
HA,
Addicks
EM.
Regional differences in the structure of the lamina cribrosa and their relation to glaucomatous optic nerve damage.
Arch Ophthalmol.
1981;
99:
137–143.
Dandona
L,
Quigley
HA,
Brown
AE,
Enger
C.
Quantitative regional structure of the normal human lamina cribrosa. A racial comparison.
Arch Ophthalmol.
1990;
108:
393–398.
Zhang
L,
Albon
J,
Jones
H,
et al.
Collagen microstructural factors influencing optic nerve head biomechanics collagen microstructural factors.
Invest Ophthalmol Vis Sci.
2015;
56:
2031–2042.
Mozaffarieh
M,
Bärtschi
M,
Henrich
P,
Schoetzau
A,
Flammer
J.
Retinal venous pressure in the non-affected eye of patients with retinal vein occlusions.
Graefes Arch Clin Exp Ophthalmol.
2014;
252:
1569–1571.
Gasser
P,
Flammer
J.
Blood-cell velocity in the nailfold capillaries of patients with normal-tension and high-tension glaucoma.
Am J Ophthalmol.
1991;
111:
585–588.
Guidoboni
G,
Harris
A,
Arciero
JC,
et al.
Mathematical modeling approaches in the study of glaucoma disparities among people of African and European descents.
J Coupled Syst Multiscale Dyn.
2013;
1:
1.
Beard
DA,
Bassingthwaighte
JB.
Modeling advection and diffusion of oxygen in complex vascular networks.
Ann Biomed Eng.
2001;
29:
298–310.
Goldman
D.
Theoretical models of microvascular oxygen transport to tissue.
Microcirculation.
2008;
15:
795–811.
McGuire
BJ,
Secomb
T.
A theoretical model for oxygen transport in skeletal muscle under conditions of high oxygen demand.
J Applied Physiol.
2001;
91:
2255–2265.
Carreau
A,
Hafny-Rahbi
BE,
Matejuk
A,
Grillon
C,
Kieda
C.
Why is the partial oxygen pressure of human tissues a crucial parameter? Small molecules and hypoxia.
J Cell Mol Med.
2011;
15:
1239–1253.
Sander
EA,
Barocas
VH.
Comparison of 2D fiber network orientation measurement methods.
J Biomed Mat Res A.
2009;
88:
322–331.
Mezu-Ndubuisi
OJ,
Teng
P-y,
Wanek
J.
In vivo retinal vascular oxygen tension imaging and fluorescein angiography in the mouse model of oxygen-induced retinopathy.
Invest Ophthalmol Vis Sci.
2013;
54:
6968–6972.
La Cour
M,
Kiilgaard
JF,
Eysteinsson
T,
et al.
Optic nerve oxygen tension: effects of intraocular pressure and dorzolamide.
Br J Ophthalmol.
2000;
84:
1045–1049.
Iwase
T,
Yamamoto
K,
Ra
E,
Murotani
K,
Matsui
S,
Terasaki
H.
Diurnal variations in blood flow at optic nerve head and choroid in healthy eyes: diurnal variations in blood flow.
Medicine.
2015;
94.
Downs
JC,
Roberts
MD,
Sigal
IA.
Glaucomatous cupping of the lamina cribrosa: a review of the evidence for active progressive remodeling as a mechanism.
Exp Eye Res.
2011;
93:
133–140.
Poostchi
A,
Wong
T,
Chan
KC,
et al.
Optic disc diameter increases during acute elevations of intraocular pressure.
Invest Ophthalmol Vis Sci.
2010;
51:
2313–2316.
Roberts
MD,
Liang
Y,
Sigal
IA,
et al.
Correlation between local stress and strain and lamina cribrosa connective tissue volume fraction in normal monkey eyes.
Invest Ophthalmol Vis Sci.
2010;
51:
295–307.
Sigal
IA,
Flanagan
JG,
Tertinegg
I,
Ethier
CR.
Modeling individual-specific human optic nerve head biomechanics. Part I: IOP-induced deformations and influence of geometry.
Biomech Model Mechanobiol.
2009;
8:
85–98.
Roberts
MD,
Sigal
IA,
Liang
Y,
Burgoyne
CF,
Downs
JC.
Changes in the biomechanical response of the optic nerve head in early experimental glaucoma.
Invest Ophthalmol Vis Sci.
2010;
51:
5675–5684.
O'Rese
JK,
Girkin
CA,
Budenz
DL,
Durbin
MK,
Feuer
WJ;
Cirrus OCT Normative Database Study Group.
Effect of race, age and axial length on optic nerve head parameters and retinal nerve fiber layer thickness measured by Cirrus HD-OCT.
Arch Ophthalmol.
2012;
130:
312–318.
Morgan
WH,
House
PH,
Hazelton
ML,
et al.
Intraocular pressure reduction is associated with reduced venous pulsation pressure.
PLoS One.
2016;
11:
e0147915.
Flammer
J,
Konieczka
K.
Retinal venous pressure: the role of endothelin.
EPMA J.
2015;
6:
1–12.
Chen
H-F,
Chen
M-C,
Lai
C-C,
et al.
Neovascular glaucoma after central retinal vein occlusion in pre-existing glaucoma.
BMC Ophthalmol.
2014;
14:
1.
Charlson
ME,
de Moraes
CG,
Link
A,
et al.
Nocturnal systemic hypotension increases the risk of glaucoma progression.
Ophthalmology.
2014;
121:
2004–2012.
Sredar
N,
Ivers
KM,
Queener
HM,
Zouridakis
G,
Porter
J.
3D modeling to characterize lamina cribrosa surface and pore geometries using in vivo images from normal and glaucomatous eyes.
Biomed Opt Express.
2013;
4:
1153–1165.
Ivers
KM,
Sredar
N,
Patel
NB,
et al.
In vivo changes in lamina cribrosa microarchitecture and optic nerve head structure in early experimental glaucoma.
PLoS One.
2015;
10:
e0134223.
Faridi
OS,
Park
SC,
Kabadi
R,
et al.
Effect of focal lamina cribrosa defect on glaucomatous visual field progression.
Ophthalmology.
2014;
121:
1524–1530.
Jung
KI,
Park
C.
Lamina cribrosa depth is associated with cup-to-disc ratio in eyes with large optic disc cupping and cup-to-disc ratio asymmetry.
Invest Ophthalmol Vis Sci.
2015;
56:
1004–1004.
Sawada
Y,
Hangai
M,
Murata
K,
Ishikawa
M,
Yoshitomi
T.
Lamina cribrosa depth variation measured by spectral-domain optical coherence tomography within and between four glaucomatous optic disc phenotypeslamina cribrosa depth in glaucomatous phenotypes.
Invest Ophthalmol Vis Sci.
2015;
56:
5777–5784.
Lee
KM,
Kim
T-W,
Weinreb
RN,
Lee
EJ,
Girard
MJ,
Mari
JM.
Anterior lamina cribrosa insertion in primary open-angle glaucoma patients and healthy subjects.
PLoS One.
2014;
9:
e114935.
Wang
X,
Rumpel
H,
Lim
WEH,
et al.
Finite element analysis predicts large optic nerve head strains during horizontal eye movementseye movements induce optic nerve head strains.
Invest Ophthalmol Vis Sci.
2016;
57:
2452–2462.
Pries
A,
Secomb
TW.
Rheology of the microcirculation.
Clin Hemorheol Microcirc.
2003;
29:
143148.
Pries
AR,
Secomb
TW.
Blood flow in microvascular networks.
Comp Physiol.
2011;
1:
3–36.
Cristini
V,
Kassab
GS.
Computer modeling of red blood cell rheology in the microcirculation: a brief overview.
Ann Biomed Eng.
2005;
33:
1724–1727.
Wang
L,
Cull
G,
Burgoyne
CF,
Thompson
S,
Fortune
B.
Longitudinal alterations in the dynamic autoregulation of optic nerve head blood flow revealed in experimental glaucoma dynamic blood flow autoregulation in glaucoma.
Invest Ophthalmol Vis Sci.
2014;
55:
3509–3516.
Wang
L,
Burgoyne
CF,
Cull
G,
Thompson
S,
Fortune
B.
Static blood flow autoregulation in the optic nerve head in normal and experimental glaucomaoptic nerve autoregulation in experimental glaucoma.
Invest Ophthalmol Vis Sci.
2014;
55:
873–880.
Quigley
HA,
Addicks
EM,
Green
WR.
Optic nerve damage in human glaucoma: III. Quantitative correlation of nerve fiber loss and visual field defect in glaucoma, ischemic neuropathy, papilledema and toxic neuropathy.
Arch Ophthalmol.
1982;
100:
135–146.
Rosenberg
KI,
Park
SC,
Su
D,
Banik
R,
Liebmann
JM,
Ritch
R.
Dimensions of the neural canal at the optic nerve head in non-arteritic ischemic optic neuropathy compared to normal subjects.
Invest Ophthalmol Vis Sci.
2012;
53:
3928–3928.
Hinghofer-Szalkay
H,
Greenleaf
J.
Continuous monitoring of blood volume changes in humans.
J Applied Physiol.
1987;
63:
1003–1007.
Buchwald
P.
FEM-based oxygen consumption and cell viability models for avascular pancreatic islets.
Theor Biol Med Model.
2009;
6:
5.
Secomb
TW,
Hsu
R,
Park
EY,
Dewhirst
MW.
Green's function methods for analysis of oxygen delivery to tissue by microvascular networks.
Ann Biomed Eng.
2004;
32:
1519–1529.