Normalized field autocorrelation function-based optical coherence tomography three-dimensional angiography

Abstract. Optical coherence tomography angiography (OCTA) has been widely used for en face visualization of the microvasculature, but is challenged for real three-dimensional (3-D) topologic imaging due to the “tail” artifacts that appear below large vessels. Further, OCTA is generally incapable of differentiating descending arterioles from ascending venules. We introduce a normalized field autocorrelation function-based OCTA (g1-OCTA), which minimizes the tail artifacts and is capable of distinguishing penetrating arterioles from venules in the 3-D image. g1 (τ) is calculated from repeated optical coherence tomography (OCT) acquisitions for each spatial location. The decay amplitude of g1 (τ) is retrieved to represent the dynamics for each voxel. To account for the small g1 (τ) decay in capillaries where red blood cells are flowing slowly and discontinuously, Intralipid is injected to enhance the OCT signal. We demonstrate that the proposed technique realizes 3-D OCTA with negligible tail projections and the penetrating arteries are readily identified. In addition, compared to regular OCTA, the proposed g1-OCTA largely increased the depth-of-field. This technique provides a more accurate rendering of the vascular 3-D anatomy and has the potential for more quantitative characterization of vascular networks.


Introduction
Optical coherence tomography angiography (OCTA) is an implementation of optical coherence tomography (OCT) that enables noninvasive high-resolution volumetric imaging of perfused vessels. [1][2][3] The motion of scattering cells within blood vessels [the majority are red blood cells (RBCs)] causes signal fluctuation, which produces OCT signal decorrelation. OCTA arises from the motion contrast of these dynamic scatters by comparing signals at the same spatial location at different times. The perfused vessels are then visualized by amplitude/intensitybased, [4][5][6] phase-based, 7 or complex signal-based 8,9 OCTA data processing algorithms.
Compared to Doppler OCT, which is limited to measuring axial flow, [10][11][12] OCTA provides comprehensive images of both penetrating and transverse vessels since OCTA measures the intensity and phase fluctuation arising from RBC dynamics. However, unlike Doppler OCT, OCTA does not distinguish descending vessels (flowing into brain) from ascending vessels (flowing toward brain surface). In comparison to two-photon microscopy angiography, 13 which offers topological threedimensional (3-D) image of blood vessel networks by labeling the plasma, OCTA is often limited to en face two-dimensional projections of vascular networks due to the "tail" artifacts along the axial direction that likely arise from the multiple scattering of light within the vessel. 14 To minimize the tail artifacts, Leahy et al. 15 used a high numerical aperture objective to reject the signal from the multiple scattering tail, but this required dense depth scanning to obtain a volumetric image; Vakoc et al. 16 introduced a method to attenuate the tail signal by applying a step-down exponential filtering beneath the vessels; and Zhang et al. 14 proposed a layer-based subtraction approach, which considers the sample structural information to minimize the retinal vessel tail in the choroid layer.
In this paper, we introduce a OCTA processing algorithm (g 1 -OCTA) based on the normalized first-order field autocorrelation function (g 1 ) calculation to address the tail artifacts and identify penetrating arteries in the 3-D image. Generally, the g 1 analysis of OCT data is used to quantify particle/blood flow velocity [17][18][19][20][21] and diffusion [22][23][24] within individual voxels, providing quantitative 3-D dynamic maps by taking the advantage of 3-D imaging ability of OCT. This requires a relatively longer data acquisition time for each location (>3 to 4 ms for OCTbased methods 20 ) to produce a sufficient ensemble averaged estimate of the experimental g 1 . In this work, we use g 1 to qualitatively estimate the dynamics (denoted as dynamic index) of each OCT resolution voxel by calculating the maximum decorrelation of |g 1 | in a shorter data acquisition time (<0.6 ms for each location). The use of the shorter data acquisition time is the key to minimize the decorrelation in regions beneath the vessels (i.e., tail artifacts), where the OCT signal is a combination of RBC multiple scattering (dynamic) and tissue back scattering (static). However, due to the slow flow velocity and single-file nature of RBCs flowing in capillaries, it is highly possible that either no or very little decorrelation will be detected in some segments of capillaries. To address this issue, Intralipid solution is administrated to enhance signal detection in capillary networks, as utilized previously for OCTA imaging. [25][26][27] In addition, using the high-frequency phase information of g 1 , we are able to identify the flow direction of penetrating vessels. We compared the proposed g 1 -OCTA with regular complex signal-based OCTA, 8 in terms of mitigating vessel tail artifacts and improving the imaging depth-of-field, by imaging an anesthetized mouse brain.

OCTA Tail Artifacts and Contrast Enhancement with Intralipid
OCT imaging uses interferometric detection of scattered light, more precisely the backward scattered light, to resolve the depth-dependent scattering contrast. RBCs are the dominant scattering objects inside blood vessels as the refractive index of the RBC is larger than the surrounding plasma by ∼3% to 6%. 28,29 As illustrated in Fig. 1(a), multiple forward scattering events followed by a backward scattering event are likely to happen within or through a macrovessel (top), whereas RBC backward scattering or RBC forward-tissue backward scattering may happen in a capillary. The multiple forward scattering events result in dynamic contrast appearing a path length deeper than the vessel in the otherwise static tissue, which is commonly observed in OCTA images as shown in Fig. 1(b). Signals measuring photons that have experienced forward scattering through the vessel, and then backscattering from the static tissue, will decorrelate slower than signals measuring photons that experience backscattering from dynamic RBCs within the vessel. Further, these dynamically scattered photons will mix with static photons that have traveled the same pathlength but have only scattered from the static tissue, resulting in a smaller decay amplitude. Thus we expect the tails beneath large vessels to have a smaller and slower decorrelation. As OCTA detects the decorrelation of the OCT signal at the same spatial location occurring within a certain time interval, it is expected that a relatively longer time interval will result in a larger decorrelation and thus higher signal contrast. However, with the prolonged time interval used in OCTA of generally 6 ms or longer, the decorrelation beneath vessels increases and produces signal contrast comparable to that observed inside vessels. Figure 1(c) illustrates the decorrelation of jg 1 ðτÞj [Eq. (1)] with time lag spanning 4 ms at voxels above, inside, and beneath the vessels. Note that the fastest decay happens inside the blood vessel (red curves), tails have moderate decorrelation (magenta curves), whereas little decorrelation is observed in static tissue above the large vessel (black curves). For regular OCTA using Bscan repeat acquisition, the time interval for the same location is usually >6 ms, which results in a comparable decorrelation in regions beneath the vessel to that within the vessel. Thus tails are widely observed in an OCTA 3-D image.
To minimize the decay of vessel tails, one feasible solution is to reduce the time interval between serial measurements. Based on the observation that the OCT signal decorrelation time inside a blood vessel is usually on the scale of 1 to 4 ms, a submillisecond calculation of g 1 ðτÞ would be appropriate to get sufficient decorrelation inside the vessel while minimizing the decorrelation from the tails. However, as shown in the top row of Figs. 2(a) and 2(c), we noticed that it is difficult to obtain sufficient decorrelation in some capillaries in such a short time [e.g., v2 and vessels marked with white arrows in the bottom panel of Fig. 2(a)] where RBCs are flowing slowly and in single file. To address this issue, we employed Intralipid, which is an FDA approved nutritional supplement, to enhance OCT contrast in small vessels following the protocols previously described in Refs. [25][26][27]. The Intralipid solution is a highly scattering lipid emulsion and with a mean particle diameter of 226 nm. 30 When injected into blood vessels, Intralipid behaves like a blood plasma tracer and exhibits translational and diffusive motion, both of which contribute to the OCT signal fluctuation. As shown in the bottom row of Figs. 2(a) and 2(c), Intralipid greatly enhances the OCT signal contrast in blood vessels. In addition, we noticed that the OCT signal with Intralipid decorrelates faster than that with intrinsic contrast (RBCs), as illustrated in Fig. 2(b). Therefore, the administration of Intralipid solution not only enhances dynamic contrast in blood vessels in terms of flow continuity (compared to discrete RBCs) but also reduces the acquisition time required to achieve an equivalent decorrelation.

g 1 -OCTA
As noted in Sec. 2.1, a shorter time interval is required to minimize the tail artifacts and Intralipid administration is needed to enhance the signal in the capillaries. Therefore, we applied a repeat A-scan scanning protocol (i.e., the imaging beam was fixed at one location to collect a number of consecutive A-scans then swept laterally for the whole volume) for data acquisition instead of the typical repeat B-scan scanning protocol used by regular OCTA, and we utilized Intralipid to enhance the capillary signals. The normalized field autocorrelation function g 1 ðτÞ for each voxel was first obtained via E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 1 ; 6 3 ; 4 1 5 g 1 ðτÞ ¼ where RðtÞ is the complex OCT signal at time t, R Ã ðtÞ is the complex conjugate, h i indicates ensemble averaging, and τ is the time lag. We calculate a dynamic contrast index I d , as the maximum decay of the correlation function since the first time lag over a submillisecond time period E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 2 ; 3 2 6 ; 5 0 3 I d ðx; y; zÞ ¼ jg 1 ð1Þ ðx;y;zÞ j − min jg 1 ðτÞ ðx;y;zÞ j (2) Figure 3 illustrates the calculation of the dynamic contrast index I d and the determination of flow direction. For example, the dynamic index is  stands for Intralipid enhanced, normalized field autocorrelation function-based OCTA unless otherwise stated.

Image Processing
The signal-to-noise ratio (SNR) of where I d;vessel and I d;bkgd are the dynamic index of blood vessels and background, respectively. I d;vessel and I d;bkgd were distinguished using a threshold that was selected to be the average signal in the volumetric image plus the standard deviation.
As the majority of the image volume is background signal, this provides a threshold below which most voxels are background and above which most voxels are vessels.
In Figs. 5(d) and 5(e), blood vessel cross section is calculated in each XZ plane using the MATLAB label connected component function (bwlabel function in MATLAB) based on the binarized image of I d;vessel . The number of pixels for each vessel cross section was counted and vessels with <15 pixels were omitted since a minimum XZ cross-sectional area of ∼60 μm 2 was assumed for a capillary, noting that the imaging point spread function is ∼3.5 μm isoptropic for our OCT system.

Animal Preparation and Imaging Protocol
To validate the advantages of g 1 -OCTA, we imaged an anesthetized mouse brain through a cranial window using a Thorlabs OCT system. During surgery, the animal was anesthetized with isoflurane (2% to 3% induction, 1% to 2% maintenance, in 100% oxygen) and body temperature was maintained with a homeothermic unit (Harvard Apparatus). A craniotomy (3 × 3 mm) was performed over the left somatosensory cortex and the dura was removed. The cortex was covered with agarose (1% in saline), then protected with a 5-mm-diameter glass. The window was sealed with dental cement. The animal was then placed under the OCT imaging system. A spectral-domain OCT system (1310-nm center wavelength with a bandwidth of 170 nm, Thorlabs Inc.) 21,31 was used for this experiment. The axial resolution of the system was 3.5 μm and imaging speed was 47,000 A-scan/s. A 10× objective was used for this study allowing a transverse resolution of 3.5 mm. Both regular OCTA acquisition (i.e., twice-repeated B-scans with 10 times whole volume repeats) and the proposed g 1 -OCTA acquisition (i.e., repeated A-scans) were performed before and after Intralipid administrated (20% Intralipid solution, 6 ml/kg of body weight, and femoral artery injection). The imaged region of interest (ROI) was 600 × 600 μm (400 × 400 pixels in X and Y) with light focused ∼150 μm beneath the brain surface. In order to compare the effect of decorrelation time upon g 1 -OCTA image quality, A-scans were acquired 100 times at each spatial location (lasting ∼2 ms), then g 1 ðτÞ was calculated with different autocorrelation time periods (i.e., n τ ) over different observation times (i.e., n t ) using Eq. (1).
The largest decorrelation of g 1 ðτÞ is retrieved to represent the dynamic contrast index (I d ) for each (x; y; z) voxel through Eq. (2). All animal experimental procedures were reviewed and approved by the Massachusetts General Hospital Subcommittee on Research Animal Care. Figure 4 compares the results obtained with regular OCTA and g 1 -OCTA before versus after Intralipid administration. Regular OCTA were obtained with our routine practice that 10 volume repeats were acquired for averaging (total acquisition time for a 400 × 400 XY pixel ROI with 10 times averaging: ∼70 s) and g 1 -OCTA were calculated from g 1 ðτÞ with n τ ¼ 25 (time lag) over n t ¼ 30 (sampling points) without averaging (equivalent total acquisition time for the same ROI: ∼105 s). We denote g 1 -OCTA calculated with n τ time lags (i.e., autocorrelation time period) over n t sampling points (i.e., observation time) as n τ ∕n t , e.g., 25/30. Note that the physically relevant parameters are the time delays over which the autocorrelation function is calculated (i.e., autocorrelation time period) and the total observation time. These physically relevant parameters are linearly related, respectively, to our reported variables of n τ and n t by the sampling time between consecutive A-lines which depends on the particular OCT system used in the measurement. For a fair comparison, regular OCTA voxels were 3-D normalized to the range of [0 1] and all regular OCTA and g 1 -OCTA images were presented with a threshold of mean minus standard deviation of the 3-D volume. The administration of Intralipid greatly enhanced the dynamic contrast for both regular OCTA and g 1 -OCTA, especially helpful for g 1 -OCTA, as shown by the en face maximum intensity projection (MIP) images in Figs. 4(a)-4(d). Compared to regular OCTA, g 1 -OCTA largely mitigates the tail artifacts as shown by the XZ cross plane images [Figs. 4(c) and 4(d)] and the single-depth en face images in Fig. 4(e). The white arrows in Fig. 4(e) clearly suggests the suppression of tail artifacts by g 1 -OCTA. Since both images in Fig. 4(e) were XZ crosssectional images obtained at the same depth (thickness: ∼3 μm), the more network-like feature of the left figure further indicates the tail artifacts issue of regular OCTA and the tail artifacts suppression by g 1 -OCTA. The ability to accurately identify vessel cross sections of g 1 -OCTA enables 3-D topologic imaging as shown by supplemental movies for 3-D comparison. In addition, we noticed that g 1 -OCTA largely increased the depth-of-field compared to regular OCTA as shown by the XZ single/MIP images in Figs. 4(c) and 4(d) and the single-depth en face image of Fig. 4(f), in which the color bar was scaled to [0.02, 0.1] for the left figure of regular OCTA. This greater depth of field arises in part because of the noise suppression provided by the g 1 calculation and the use of g 1 ð1Þ (the first time lag) minus the minimum value of g 1 ðτÞ to get the dynamic contrast index (I d ). In the out-of-focus region, the SNR is reduced due to the decrease of number of photons backscattered from these depths. The change in phase between photons scattered over time at the same depth, however, should not change.

Results
Therefore, phase-based imaging techniques like g 1 -OCTA are less sensitive to the decrease in number of backscattered photons thus providing a larger depth of field. In addition, the calculation of I d using g 1 ð1Þ instead of g 1 ð0Þ to minus the minimum of g 1 ðτÞ helps to cancel out voxels with low SNR. For example, the minimum value of g 1 ðτÞ will be close to 0 and g 1 ð1Þ will also approach 0 when SNR is decreasing, which results in I d reducing to 0. Therefore, we used g 1 ð1Þ instead of g 1 ð0Þ [g 1 ð0Þ is independent of noise and equals to 1] minus the minimum value of g 1 ðτÞ to get the dynamic contrast index I d .
We investigated the effect of autocorrelation time period (the number of time lags n τ ) and observation time (the number of sampling points n t ) on g 1 -OCTA. Qualitatively, increasing n t (i.e., longer observation time) results in a higher SNR as indicated by comparing the g 1 -OCTA results of 25/30, 25/50, and 25/100 in Fig. 5(a); whereas increasing n τ provides a longer observation time period of signal decorrelation resulting in a larger dynamic contrast index (I d ) as indicated by the increased image intensity of the g 1 -OCTA results of 25/100, 50/100, and 75/100 in Fig. 5(a). However, increasing n τ also leads to larger tail artifacts as indicated by the green rectangles in the middle row of Fig. 5(a). The reason for a larger I d in penetrating vessels [e.g., the one indicated by dark blue arrow in the en face MIP of 25/50 in Fig. 5(a)] compared to that of the parent transverse vessels [e.g., the one indicated by light blue arrow in the en face MIP of 25/50 in Fig. 5(a)] is due to the fact that g 1 decays faster from axial speed, as shown in Fig. 5(b). Figure 5(c) compares the SNR and the background noise for different n τ /n t pairs. The SNR increased from 3.5 to 8 and the background noise reduced from 0.08 to 0.03 when n t increased from 30 to 100 while keeping n τ unchanged (n τ ¼ 25) suggesting that the SNR of g 1 -OCTA is improved with longer observation time. We further compared the vessel detectability of g 1 -OCTA and regular OCTA by counting the number of vessel cross sections in each XZ plane and the total number of pixels occupied by these vessel cross sections (600X × 400Z μm 2 ; threshold: mean + std.) as described in Sec. 2.3. As shown in Fig. 5(d), g 1 -OCTA can detect almost twice the number of vessels and more total number of "vessel pixels" in XZ plane compared to regular OCTA (either with intrinsic RBC contrast or with Intralipid enhanced contrast). This arises because g 1 -OCTA has a better imaging ability for dynamic contrast in the outof-focus regions so that more vessel cross sections could be detected. Figure 5(e) shows the average number of pixels per vessel cross section. The larger value arises because vessels detected by regular OCTA have a larger size because of larger vessel tails reconstructed by regular OCTA. In addition, we noticed that the number of vessels, number of vessel pixels, and number of pixels per vessel obtained with g 1 -OCTA increased slowly with increasing n τ (e.g., 25/100, 50/100, and 75/100) as the pixels with a slower decaying autocorrelation function, including the tail pixels, become detectable. Figure 6(a) compares the cross-sectional g 1 -OCTA images and the profiles along the center of the two vessels obtained with g 1 -OCTA with different pairs of n τ ∕n t and with regular OCTA after Intralipid administered. Note that while the contrast-to-background ratio of vessels is enhanced, the artifact signal from the tail decorrelation becomes stronger when increasing the time lag. Figure 6(b) shows the normalized axial profiles averaged across 22 vessels randomly selected from the en face MIP (with no prior viewing of the XZ projection) with diameters ranging from 21 to 25 μm (full width at half maximum). Clearly, g 1 -OCTA provides a better tail artifact suppression and higher SNR compared to regular OCTA. We also noticed that the tail signal increases when increasing the time lag for g 1 -OCTA. For instance, the g1-40/80 and g1-50/100 profiles show larger tails compared to the profiles obtained with shorter time lags, e.g., g1-15/30, g1-20/30, g1-25/30, and g1-30/60. In our opinion, g 1 -OCTA obtained with an autocorrelation time period of ∼0.53 ms (i.e., n τ ¼ 25 for our OCT system) over an observation time of ∼0.636 ms (i.e., n t ¼ 30 for our OCT system) would be the appropriate choice considering the image quality and data acquisition time. Note that this is the recommendation for the physically relevant parameters of autocorrelation time period and observation time and the specific n τ and n t values are subject to change depending on the sampling time for consecutive A-lines for the given OCT system. For example, n τ should be around 50 and n t should be around 60 for an OCT system with a 100-kHz A-line rate. Therefore, the total data acquisition time for a 400 × 400 pixel image would be ∼105 s. Finally, Fig. 6(c) shows the en face MIP and en face images at depths = 24, 48, and 71 μm of the ROI marked by the white rectangle.

Conclusion and Discussion
In this paper, we introduced g 1 -OCTA that uses the normalized field autocorrelation function to obtain a microvascular angiogram that effectively suppress vessel tail artifacts, maintains a high SNR in the out-of-focus region, and distinguishes axial flow direction. The normalized field autocorrelation function ½g 1 ðτÞ was first calculated from a time series signal acquired by repeating A-scans at the same spatial location, then the maximum decay of g 1 ðτÞ was retrieved to represent the dynamic contrast in this voxel. To enhance dynamic contrast in capillaries where RBCs may flow slowly and discontinuously, Intralipid was administered through a femoral artery injection. Compared to regular OCTA repeating B-scans (∼5 to 10 ms time interval), 4,5-9 g 1 -OCTA detects dynamics over a much shorter time scale (∼0.02-ms time interval with total acquisition time <1 ms). The decorrelation of the vessel tail signal is negligible during such a short decorrelation time compared to that within the vessels, thus minimizing the tail artifacts. Compared to regular OCTA, g 1 -OCTA effectively suppresses the tail artifacts and is sensitive enough to detect dynamics in the out-offocus region, providing a greater depth-of-field. In our opinion, g 1 -OCTA obtained with an autocorrelation time period of ∼0.53 ms (i.e., n τ ¼ 25 for our OCT system) and an observation time of ∼0.636 ms (i.e., n t ¼ 30 for our OCT system) would be the appropriate choice to provide high image quality with an acceptable total data acquisition time (∼105 s for an ROI consisting of 400 × 400 pixels). Increasing the observation time (i.e., n t ) while keeping the autocorrelation time period (i.e., n τ ) constant will enhance the SNR of g 1 -OCTA. As expected, increasing the autocorrelation time period (i.e., n τ ) results in larger tail artifacts. We note that the selection of the autocorrelation time period is preferably done adaptively as g 1 decays faster for large vessels due to higher particle flow speed while it takes longer for g 1 to decay to the same level for small vessels or capillaries with slower flow. A shorter autocorrelation time period is preferred for large vessels to more effectively suppress the tail artifacts, whereas a longer autocorrelation time period is preferred for slow flowing vessels to enhance the dynamic index, I d . One drawback of g 1 -OCTA compared to regular OCTA is the longer acquisition time required for a volume image. However, it is worth noting that multiple averaging (usually 10 times averaging) is performed in practice for regular OCTA imaging, which takes ∼70 s for the same ROI. One approach to reduce the imaging time of g 1 -OCTA is to apply it with a fast full field OCT, 32,33 which will significantly reduce the imaging time to a few seconds. In summary, we have demonstrated that g 1 -OCTA has the potential to provide accurate 3-D maps of the vasculature. In the future, this method may be further applied to mapping the RBC transient time in the microvasculature with OCT velocimetry, such as phase resolved Doppler OCT 12 and dynamic light scattering-optical coherence tomography, 21 which are capable of measuring blood flow velocity in both large vessels and small capillaries.

Disclosures
The authors declare no competing financial interest.