Relative retinal flow velocity detection using optical coherence tomography angiography imaging.

Optical coherence tomography angiography (OCTA) imaging is a valuable tool for the visualization of retinal vasculature at an unprecedented level of details. However, due to relatively long time-interval between repeated scans in the conventional OCTA scanning protocol, the OCTA flow signal suffers from low dynamic range and loss of velocity-intensity correlation. The ability to distinguish fast and slow flow in the retina may provide a powerful tool for the assessment of early-stage retinal diseases such as vein occlusion. Here, we report a method to detect relative flow velocity in human retina using a 67.5 kHz spectral-domain OCTA device. By adapting the selection of A-scan time-intervals within a single OCTA acquisition and combining the resulting OCTA images, we expand the detectable velocity range. After a quantitative validation of this method performing microchannel flow experiments with varying flow velocities, we demonstrate this approach on human eyes using CIRRUS HD-OCT 5000 with AngioPlex (ZEISS, Dublin, CA) through a prototype scanning pattern.


Introduction
Until the invention of optical coherence tomography (OCT) [1,2], ophthalmologists were limited to surface visualization of the retina e.g. by fundus photography. The development of optical coherence tomography angiography (OCTA) techniques made non-invasive imaging of the retinal blood flow possible [3]. Since then, scientists have shown that apart from structural changes (e.g. tractional retinal detachment from fibrovascular proliferation in diabetic retinopathy [4]) dynamical changes within the retinal blood flow come along with impairment of vision. For example, Burgansky-Eliash et al. found that age-related macular degeneration (AMD) is accompanied by decreased blood flow rates [5]. Similarly, choroidal blood flow changes temporarily upon the surgical treatment of rhegmatogenous retinal detachment by scleral buckling [6]. Therefore, tracking the retinal vascular blood flow during the recovery process may be essential for a successful therapy.
Considering the variety of possible complications accompanied with abnormal retinal blood flow, it is therefore desirable to study and diagnose the blood flow behavior. Conventional OCTA does not provide quantitative information of the flow due to relatively long time-intervals between repeated scans. However, researchers have developed advanced methods such as Doppler-based OCTA [7,8] and OCTA capillary velocimetry based on Eigendecomposition [9,10] to map the blood flow. Recently, Ploner et al. implemented a variable time-interval analysis (VISTA) algorithm in an OCTA acquisition using a high A-scan rate SSOCT system to visualize blood flow speeds [11]. Most techniques face at least one of the following limitations: long acquisition times, technical complexity, restriction to visualize only a certain flow flux range, restriction to specific sample orientation. For instance, while traditional Doppler-based OCT (DOCT) acquires quantitative information, its sensitivity is limited to flow perpendicular to the observation direction. However, recent developments in DOCT imaging introduced alternative algorithms, such as Doppler variance method [12,13] and relative standard deviation OCT [14], to overcome the angle dependency of the Doppler signal. The VISTA-algorithm on the other hand provides flow information over a large dynamical range and is not orientation dependent. However, it also requires a high-speed scanning (400 kHz A-scan rate) system which has not been yet realized in ophthalmic clinical devices. Other methods which quantify the retinal blood flow measure the laser-speckle pattern [15] or a fluorescence signal such as in Fluorescein angiography (FA) [16]. FA is capable of imaging a wider field but comes with other drawbacks such as implementation difficulties and invasiveness which can lead to potential adverse events such as vomiting. In this manuscript, we present a fast and simple method using a commercially available OCTA device capable of expanding the resolvable dynamic range of flow velocities. This method combines MB-mode scanning (multiple repeated A-scans at one location contribute to final B-frame image) and adaptive selection of A-scan time-intervals to enable simultaneous representation of multiple flow regimes in a vasculature network such as present in the retina.

Methods
Following OCT scanning terminology, we implemented an MB-scan pattern, which implies multiple M-scans taken over different locations to create a 2D OCT-data set of M-scans, on a commercial CIRRUS HD-OCT 5000 with AngioPlex. The field of view (FOV) of the scanned area was 2.14 × 1 mm. A full scan consisted of 82 B-scans with 175 A-scans per B-scan and two additional B-scans performed perpendicular to the OCT cube for motion correction resulting in a total scanning time of 3.5 s. The scanning laser was operating at a central wavelength of 840 nm with a bandwidth of 60 nm. Each A-scan location was repeated 10 times (n=10) before moving to the next adjacent location ( Fig. 1(A)). The A-scan rate of the system was 67.5 kHz corresponding to a time-interval ∆T of ∼15 µs between two consecutive A-scans. This scanning strategy theoretically enables access to OCTA acquisition intervals of 15 µs -135 µs in increment of 15 µs. Practically, all n A-scans per scanning location could be utilized to perform the OCTA signal calculation. In order to increase the interval-time ∆T we considered every second (i.e. 1 → 2, 2 → 3, . . . ), third (i.e. 1 → 3, 2 → 4, . . . ), fourth (i.e. 1 → 4, 2 → 5, . . . ) and fifth (i.e. 1 → 5, 2 → 6, . . . ) A-scan for the decorrelation calculation. In this way we obtained decorrelation maps at ∆T=15 µs, ∆T=30 µs, ∆T=45 µs and ∆T=60 µs, respectively. According to Woo June Choi et. al. [17], the OCTA signal correlates with the flow velocity v and the time-interval ∆T by where N is the number of scattering particles within an imaging voxel, α the backscattering coefficient of the sample, λ 0 the central wavelength of the OCT light source and ∆φ the difference between the phases, respectively. N can be expressed in terms of the scatter concentration c with where V Voxel is the volume of the imaging voxel and V Particle is the volume of the scattering particle, respectively. Assuming a time independent backscattering coefficient α, a mean phase-difference of ∆φ=0 and a constant number of particles within an imaging voxel (fixed concentration of scattering particles in a voxel volume) we can express a proportionality of the OCTA with the flow velocity by For the microfluidic experiments we used a custom fabricated multi-microchannel phantom ( Fig. 1(B)). The main channel had a rectangular cross-section with a width of 240 µm and a constant height of 40 µm for all channels. The material of the channel body was Polydimethylsiloxane (PDMS) mixed with TiO 2 powder in order to mimic the characteristics of the static retinal tissue. To mimic the retinal blood flow, we used a syringe pump (SyringeOne, New Era Instruments) to inject a mixture of whole milk and distilled water (mixing ratio 1:10) into the phantom at controlled pumping rates. Knowing the channels cross section and the applied pump rate we could calculate the flow velocity v inside the phantom. A focusing scanning lens was mounted on an adaptor replacing the focusing function of the cornea and lens. The numerical aperture of this lens was 0.1 resulting in a transverse resolution of 20 µm. Axial resolution was determined as 5 µm. One imaging voxel of the OCT scan had the dimensions of 20 µm x 20 µm x 2 µm (imaging depth=2 mm, pixels in axial direction=1024). A custom segmentation algorithm was used to detect the rear surface of the phantom ( Fig. 1(C)). Subsequently, 20 voxel layers (channel height/lateral voxel size=40 µm/2 µm) were selected for averaging and image reconstruction.
Next, we scanned three normal human eyes under an IRB-approved study using a commercial angiography 3´3 mm scan and the novel prototype MB-mode 2.14´1 mm scan pattern on CIRRUS 5000 HD-OCT with AngioPlex OCT Angiography. Prior to calculating decorrelation motion correction algorithms were applied as demonstrated previously by An et al. [18]. To correct for axial motion the system acquires two B-scans in perpendicular direction to the OCT cube B-scans. Lateral eye movement was corrected by incorporating the retina tracking in CIRRUS 5000 HD-OCT. If too much motion was detected, a rescan was performed. OCT datasets were processed using Optical Micro Angiography (OMAG C ) among frames with different time-intervals (15 µs, 45 µs, 60 µs). We summed together flow data from one flow acquisition but different time-intervals to create one flow volume per scan. Subsequently, vascular en face images were generated from the inner limiting membrane (ILM) to retinal pigment epithelium (RPE) thereby avoiding signal contribution from the choroidal layer ( Fig. 1(D)).

Results
As shown previously by Choi et al. the OCT decorrelation signal strength depends on the flow velocity v and the time-interval ∆T [19]. High flow flux results in a stronger decorrelation value which is represented by a brighter pixel. Static tissue, on the other hand, is represented by dark pixels due to a low decorrelation value. Based on this we calibrated our system by recording the OCTA decorrelation signal of the microchannel flow for various flow velocities v and time-intervals ∆T (Fig. 2). Note, that we observed a stronger signal contrast at the channel wall location which we assign to an artifact due to dried out fluid. From this we obtained calibration graphs (Fig. 3(A)) which we used to quantify the flow at various flow regimes. To demonstrate the reproducibility of the calibration process we obtained two additional calibration graphs at ∆T=45 µs and ∆T=60 µs by flushing through different channels with a channel width of 240 µm and 60 µm, respectively, (Fig. 3(B)). Each time-interval ∆T limited the OCTA system to a certain flow range with a minimal detectable flux (detection limit) and a maximum distinguishable flux (saturation). For instance, in Fig. 3(C), at ∆T=60 µs, we only started to detect the OCTA decorrelation signal of the influx at flow velocities v>3 mm/s. Similarly, we were not able to distinguish channel flow velocities beyond v>16 mm/s. Note, that due to phantom fabrication, the channel shown in Fig. 3(C) expressed a narrowing. Having no flow inside the channel, the signal remained at the level of the background noise. However, by increasing the pump pressure and thereby the flow rate the smallest detectable flow (∼3 mm/s) was surpassed and the OCTA signal within the microchannel became stronger. Considering the continuity principle for fluids, we expected higher flow velocities at the narrowing compared to the wider channel inlet and outlet. Indeed, the higher OCTA signals indicated larger flow speeds at all flow rate settings at the location of the narrowing (Fig. 3(C)). By varying the interval-time ∆T the slope of the 'signal intensity vs. flow velocity'-graph, which indicates the flow sensitivity of the OCTA signal, changed too (graph in Fig. 3(A) and Fig. 4(A)). This implies that the increase of ∆T-settings results in increased sensitivity for certain flow velocity ranges (Fig. 4(B)). Looking at the behavior of the OCTA signal at different flow speeds and ∆T settings (graph in Fig. 3(C) and Fig. 4(A)), we noticed that for larger ∆T values and higher flow velocities the OCTA signal was saturated while it had greater sensitivity at low flow velocities. In this case, within the same imaging voxel of two A-scans separated by ∆T, all scattering particles were fast enough to be replaced within ∆T by new incoming scattering particles, leading to the maximal decorrelation value. Slow particles on the other hand have not left the beam completely and thereby did not cause a saturated signal.
Therefore, this setting is very efficient when investigating low flow velocities. When choosing a short interval-time ∆T instead we found that this particular choice of experimental parameters is suited for the investigation of high flow velocities. However, due to the presence of the minimal detectable flow, slow flows could be hardly measured and therefore could not be resolved fully. For small interval-times ∆T between two A-scans, scatters which were too slow to leave the beam within ∆T caused a low signal decorrelation. Similar to VISTA [11], we overcame this restriction of resolvable flow ranges by combining the OCTA data of a single measurement in such a way that the most sensitive flow regime in each ∆T-settings ( Fig. 4(B)) was contributing to the overall flow representation (Fig. 4(C)). This allowed us to study cases where several flow velocities were present in one acquisition.
Therefore, to validate our strategy we kept the syringe pump settings fixed (flow velocity v=8 mm/s) and reconstructed the OCTA images at four different interval-times using only OCTA signal data from a single acquisition (Fig. 5). As described before, higher flow speed comes with a stronger OCTA signal. Analyzing the data with an interval-time of ∆T=15 µs, we detected stronger flow with the OCTA signal reaching ∼0.5x the level of saturation at the most narrowing location. Reading the corresponding flow velocity value from Fig. 3(A), we get a value of v 25 mm/s. Considering the continuity principle and the geometrical dimension of the channel at this location (∼100´40 µm 2 ) the measured flow velocity agrees with the prediction. However, we could not resolve the flow at the inlet and outlet area using ∆T=15 µs. Next, we considered the reconstructed OCTA decorrelation signal image at ∆T=30 µs. While the area around the narrowing was fully saturated indicating flow velocities v>25 mm/s we obtained additional flow information at the neighboring region. However, flow information from the wider inlet and outlet of the channel were missing. Therefore, we added the OCTA signal data corresponding to an interval-time of ∆T=60 µs. The decorrelation values at the inlet site were constant (∼0.8x level of saturation), which indicated an even flow flux due to the constant pressure provided by the syringe pump across the whole channel. At the outlet side we detected a decrease of the relative flow velocity with increasing distance from the narrowing. This was due to the geometry of the microchannel phantom. At the outlet the channel was broader. Hence, according to the continuity principle the flow got weaker while the outlet channel became wider. The decorrelation signal at the narrowing of the channel was fully saturated at this setting due to high flow velocities present, which could not be distinguished at ∆T=60 µs. Following our calibration data in Fig. 3(A) the OCTA decorrelation signal indicated flow velocities at the inlet and outlet area of the channel of v 10 mm/s. Choosing a long time-interval such as ∆T=120 µs did not provide us any relative or quantitative flow information (Fig. 5) as the signal across the whole channel was mostly saturated.
The last image in Fig. 5 visualizes the obtained enhanced dynamic range: taking the flow resolving power of first three interval-time-settings (∆T=15 µs, ∆T=30 µs and ∆T=60 µs) we Finally, we demonstrated our concept on human subjects performing the prototype MB-mode scan and reconstructing flow maps choosing different interval-time settings in order to resolve and visualize different flow regimes within the retinal vascular system by color coding (Fig. 6). We observed that small capillaries and 3 rd branch vessels exhibited lower flow signals compared to that of the main arterioles. Moreover, for thick vessels we could visualize the characteristic Poiseuille flow (Fig. 6(A), insert). Superimposing a standard OCT angiography image and the corresponding flow map (Fig. 6(B)) showed that despite a small scanning area of the presented method, it was easy to assign the flow map to the larger image acquired by regular OCT-angiography. In clinical practice an ophthalmologist can therefore acquire a regular angiography/fundus image first and perform an MB-mode scan to obtain a flow map at the area of interest thereafter. Finally, in Fig. 6(C) we show the decorrelation value maps acquired at time-intervals ∆T=15 µs, ∆T=45 µs and ∆T=60 µs, respectively. Here, we were only able to obtain relative flow information since an accurate system calibration as performed with the microfluidic phantom is challenging. High flow flux regions contributed to a detectable decorrelation signal within the central area of the main retinal vessel at ∆T=15 µs whereas signal coming from 2 nd and 3 rd order branches fell into the undetectable flow regime. On the other hand, by increasing the interval-time the flow present in higher order branches became visible. Combining the low flow flux map with the high flow flux map resulted in the simultaneous visualization of all flows (Fig. 6(C), right).

Discussion
Using a relatively slow A-scan repetition-rate device such as spectral domain OCT, we showed a simple strategy to map relative flow dynamics and enhance the detectable flow range. However, a few considerations have to be made when utilizing the presented method.
The theoretical framework of our method was based on the mathematical model described in by Woo June Choi et. al. [17]. However, the proportionality in Eq. (2) assumes a constant time-independent backscattering coefficient α. Neglecting multiple scattering events this simple assumption can be used for scattering of light from spheres such as lipid globules in milk. It is not true for red blood cells (RBCs) due to their disk-like geometrical shape. Nevertheless, for a large number of scatters N within the imaging voxel, it is reasonable to assume a Gaussian distribution for the backscattering coefficient α as well as for the phase-difference ∆φ and treat the mean values as constants.
Following the theoretical discussion, we then validated the analytical model experimentally by performing measurements on a microchannel phantom. As an example, we used a milk-water mixture to simulate retinal blood flow assuming the mean size and optical properties of the scattering particles in milk (lipid globules) [20,21] resemble those of red blood cells [22][23][24]. To mimic the concentration of RBCs in blood [25] and avoid clogging of the microchannels due to drying we diluted the milk with water with a mixing ratio 1:10 which results in a scatter number per sample volume similar to blood [26]. However, while the size distribution of RBCs is narrow (mean diameter ∼7.8 ± 0.6 µm [27]) the size distribution of lipid globules can be quite broad with particle sizes below 1 µm or above 10 µm [28] which affects the scattering of light. Additionally, the average size of whole milk globules is ∼4.5 µm which despite being in the same order is still below the value of RBCs. As we aim to establish our method for biomedical use, we recognize that diluted milk is not a true substitute for blood. However, it is a reasonable sample candidate for validating studies. Optimally, one would perform therefore microchannel experiments with blood samples which may however have limitations due to clotting and coagulation. Alternatively, one could also use commercially available micro-spheres to mimic the RBCs within the blood flow.
During the measurements we controlled the flow velocity by changing the pump rate of the syringe pump. Knowing the geometry of the microfluidic phantom, we were able to calculate the exact flow speeds within the channel and use one experimental dataset which contains all the information about pump rate, time interval ∆T and corresponding OCTA signal to calibrate the whole system (Fig. 2). This allowed us to quantify the enhancement of the resolvable dynamic range and thereby validate our method.
However, to gain access to quantitative data for clinical use one would need to do these measurements within the human eye. Performing a calibration in-vivo is extremely challenging since it is difficult to locate reference points of known flow velocity within the retinal vasculature. Therefore, while our proposed method is capable of visualizing relative changes in the retinal blood flow in-vivo, quantification of the obtained results has yet still to be done. Additionally, as with any optical method, intensity-based correlations will highly depend on the opacity of the scanned material and are therefore naturally limited in their calibration accuracy. One way to calibrate the device for clinical use would include the single use of a second method such as canon laser blood flowmetry [29] in order to obtain the initial calibration graph.
Compared to alternative and recently developed methods of acquiring velocity dependent flow signals, such as Doppler variance (DV) or intensity based Doppler variance (IBDV) signals, our method is less sensitive to relative phase variations than DV. Moreover, although DV and IBDV are less sensitive to the incident angle it may be difficult to provide accurate measurement at very slow flows. The sensitivity of IBVD can be improved increasing the time interval between the A-scans. However, this would come at a cost of the total acquisition time. Therefore, while our presented method allows the representation of multiple flow speed ranges other approaches are potentially limited to a two-level representation of the flow (i.e. high or low) which in clinical practice might not be sufficient to characterize the progress or severity of a disease such as vein occlusions. Moreover, using variable B-scan delays such as the VISTA approach demonstrated by Ploner et al. on a slow A-scan rate system would saturate the flow signal and thereby prohibit the access to any flow information. Therefore, our strategy shows an alternative for visualizing various flow regimes in slow A-scan repetition-rate devices.
On the instrumentation site the range of available interval-times ∆T was restricted by the scanning mechanism. Using an OCT device with A-scan rates larger than 67.5 kHz would not only expand the accessible dynamic range of the flow towards higher flow speeds even further but also increase the number of A-scans n per location which in return would give access to more A-scan pairs of the same time interval ∆T. This is especially beneficial for high time-interval information as more data pairs will contribute to the averaging process and slow flow branches could be resolved more efficiently. Here, we used interval-times below 100 µs as higher time-intervals would provide only 1-2 datasets for these ∆T -values which resulted in noisy flow maps due to lack of averaging options. Another limitation which our system was facing is the small field of view (FOV) which can be expanded at the cost of longer acquisition times. Another aspect of the demonstrated method refers to motion artefacts due to eye movement, which can be divided twofold, axial and lateral motion, respectively. Axial motion correction was accomplished by performing two B-scans oriented perpendicular to the OCT B-scan cube. Lateral eye movements were detected by incorporated retina tracking which operates independently of the OCT scan. In the event of very strong eye movement a rescan was performed. A potential limitation of the MB-scanning is the long acquisition time of a single B-scan. Therefore, in order to efficiently correct for patient eye movement, the B-scan time length should not be longer than the retina tracking refresh period. However, one can overcome this limitation by acquiring multiple B-scans and montaging them together via post-processing.
As we performed measurements on a commercially available device our technique can be readily directly translated into the clinic and easily implemented as part of a more complex clinical diagnostic procedure, therewith enriching the set of tools for a detailed medical evaluation.

Conclusion
We demonstrated the detection of relative flow velocities with a commercial OCT device by applying a modified MB-scanning pattern at an A-scan rate multiple times slower than in previously reported studies. By adaptive selection of A-scan time-intervals, we enhanced the detection range of flow velocities within the retinal blood flow, thus enabling the simultaneous representation of various flows within one acquisition. Experiments on a microchannel phantom validated our method. The implementation of this method in a commercial device makes it an easy-to-use tool among existing diagnostic procedures. Before cost-effective and faster scanning systems capable of representing high dynamic flow-range information are commercially available, our method offers a strategy to achieve similar results. We expect that our results will contribute to further development of quantitative methods to assess the retinal blood flow.