2.1.

Animal Preparations and Experimental Protocols

All experimental protocols in animals were approved by the Institutional Animal Care and Use Committee at the University of Kentucky. Two male neonatal piglets (age: 32-days old) were purchased through the Division of Laboratory Animal Resources from the College of Agriculture Swine Center at the University of Kentucky. Their genetic background included 62.5% Landrace, 25% Large White, and 12.5% Yorkshire. The two piglets were housed together in one dedicated room in a pigpen, with the support of a heat lamp and towel-covered floor. They were fed ad libitum with milk replacer (Birthright, Ralco) for the first postnatal week and gradually transitioned to pellet food (Ignite Pre-Starter, Ralco) at 2 to 3 weeks.

To study LFO variations, transient global cerebral ischemia was induced in the two neonatal piglets (Fig. 1). The piglets were fasted for at least 4 h before surgery to avoid reflux of stomach contents. Piglets were intubated with a pediatric endotracheal tube (internal diameter: 3 mm) and maintained under anesthesia with 1.5% to 2% isoflurane. Following previous studies in neonatal piglets,⁴⁸ the anesthetized animal was placed prone and its head secured on a customized stereotaxic frame. The animal’s body was wrapped with an electric heated blanket under the surveillance of a rectal thermometer to avoid hypothermia. A multichannel respirator-oximeter sensor (8400, Smiths Medical) was used to record heartbeat per minute (BPM), respiration per minute (RPM), end-tidal carbon dioxide (${\mathrm{EtCO}}_2$), and peripheral artery blood oxygen saturation (${\mathrm{SpO}}_2$).

Fig. 1

Experimental setup for continuous CBF monitoring of neonatal piglets during transient global ischemia. (a) The anesthetized piglets were secured on a stereotaxic frame. A respirator-oximeter probe was clipped on the ear to acquire BPM, RPM, ${\mathrm{EtCO}}_2$, and ${\mathrm{SpO}}_2$. The noncontact scDCT was set up overhead and used for 3D imaging of rCBF distributions. (b) For scDCT, an ROI of $40\times 40{\mathrm{mm}}^2$ was projected on the intact skull of piglet head without the need for an invasive transcranial window. (c) The ROI was covered by evenly distributed $5\times 5$ sources (scanned) and $21\times 21$ detectors (grouped pixels).

Both piglets underwent transient unilateral and bilateral CCA ligations to create sequential decreases in CBF within the left and right cerebral hemispheres. After hairs on the head and ventral cervical skin were shaved with a clipper and cleaned with hair cream, the skin was disinfected with betadine and wiped with 70% ethanol. The animal was laid supine, and a 5-cm midline incision was made in the cervical area. The skin was retracted with four customized surgical hooks for better exposure of the surgical area. The underlying tissue was separated with blunt dissection until both CCAs were exposed. A 6-0 braided nylon suture was placed around each CCA, and a loose knot was made on each suture without disrupting the blood flow to the brain. The piglet was then laid prone. A $50\times 50{\mathrm{mm}}^2$ C-shape incision was made on the scalp to expose the skull. After a baseline measurement of scDCT for 5 min, the right knot was tightened to occlude the right CCA for $\sim 5\min$, followed by tightening the left knot for $\sim 2\min$ to induce a transient global cerebral ischemia (i.e., both right and left CCA occlusions). The left and right knots were then released sequentially, allowing for the restoration of the blood flow to the brain (i.e., recovery phase). After the restoration of CBF from transient global cerebral ischemia, the piglet was euthanized with 100% ${\mathrm{CO}}_2$ for $\sim 4\min$.

2.2.

scDCT System for Data Acquisition

Details about scDCT techniques for blood flow measurements are available in our previous publications.^{48,52,54,55,57,60–62} Briefly, a high-speed scanning galvanometer mirror positioning system (maximum scan angle: $\pm 12.5\mathrm{deg}$, switching time $<1\mathrm{ms}$, GVS002, Thorlabs) remotely and sequentially delivered coherent focused-point NIR light (780 nm, CrystaLaser) to multiple source positions on the tissue surface for boundary data acquisition [Fig. 1(a)]. An adjustable iris diaphragm (SM05D5, Thorlab) was placed into the source path to optimize the intensity and spot size of the incident light. A scientific CMOS (sCMOS) camera (pixels: $2048\times 2048$; frame rate: 30/s; quantum-efficiency: 50% at 800 nm, ORCA-Flash4.0, Hamamatsu Photonics) was used as a high-density 2D detector array to collect re-emitted light from the tissue for measures of spatial diffuse speckle contrasts in a selected ROI. The movement of red blood cells in the measured tissue volume produced spatial fluctuations of laser speckles captured by the sCMOS camera with a typical exposure time of 2 ms. A zoom lens (Zoom 7000, Navitar) was connected to the camera for adjusting the size of the ROI [Fig. 1(b)]. A long-pass filter ($>750\mathrm{nm}$, #84-761, EdmundOptics) was installed in front of the zoom lens to reduce the impact of ambient light. A pair of polarizers were added across the source and detection paths to reduce specular reflection directly from the scanning light sources on the tissue surface.

The noncontact scDCT system was adjusted above the piglet to cover an ROI of $40\times 40{\mathrm{mm}}^2$ on the head after skull exposure and before CCA occlusion commencement. The continuous data acquisition procedure consisted of sequentially scanning a set of $5\times 5$ source locations on the ROI in a cycling fashion [Fig. 1(c)]. The total sampling time for scanning over 25 source positions was 5 s (0.2 s per source position). For each source position, the sCMOS camera, as a 2D detector array, output a raw intensity image.

The power of incident point light that reached tissue surface from the scDCT was $<0.5\mathrm{mW}$, according to the measurement by a power meter. This low level of power is considered safe according to the American National Standards Institute standard under the Accessible Emission Limit Class 3R classification.⁵⁵ Moreover, the eyes of the anesthetized piglet were always closed during scDCT measurements. Thus the point light scanned over the selected ROI on the head and never shined on the eyes. Nonetheless, further eye protection will be made in the future by placing a protective eyewear on their eyes.

2.3.

Boundary Data Processing for CBF

For the 3D reconstruction of CBF, source and detector positions on the mesh surface must be supplied. To precisely find the actual source location, the intensity image [Fig. 2(a)] was first converted to a binary image [Fig. 2(b)]. This was done by replacing all intensity values above a globally determined threshold with “1” and setting all other values to “0.” The threshold value was determined by Otsu’s method to minimize the intraclass variance of the thresholded black and white pixels.⁶³ After thresholding, those thresholded pixels with a value of “1” that connected to 8 pixels in local neighborhoods (so called 8-connected) were found and grouped together. In this automatic process, the “bwlabel” built-in function in MATLAB was employed to label the connected components in the binary image. By setting the “conn” parameter in the bwlabel function to 8, eight neighboring pixels were determined. Then the “regionprops” built-in function with the “area” property was applied to compute the connected area, labeled as pixel numbers. The area with the maximum number of pixels was identified using the max function in MATLAB and other areas were then removed from the 2D binary image. Morphological operators, such as opening, closing, and filling, were then applied to the binary image for removing thin protrusions and isolated areas to reduce noises [Fig. 2(c)]. Finally, the centroid and radius of the point source were determined by estimating the center of the remaining largest area using the regionprops function with the “centroid” property [Fig. 2(d)]. This automatic process was repeated for all source locations [Fig. 2(e)].

Fig. 2

scDCT image processing for determining source locations and reconstruction steps for generating 3D CBF images. (a) A raw image captured by the camera. (b) The raw image from (a) converted to a binary image, derived by replacing all values above a globally determined threshold with “1” and setting all other values to “0.” (c) The object from (b) with the maximum number of connected pixels where morphological operators were applied to remove the noise from the image. (d) The centroid and radius of the point source was determined by estimating the center of mass for the remaining object in (c). This process was performed on all raw images to determine each source location. (e) Twenty-five raw images with their centroid location were added for presentation purposes. (f) Conventional reconstruction paradigm: each 3D CBF reconstruction was performed after fully updated boundary data were available (every 5 s). (g) Moving window reconstruction method: the first 3D image reconstruction took 5 s. However, a new 3D CBF image was then produced as soon as the raw image from one new source position was available, which occurred approximately every 0.2 s.

After finding these source locations on the $X-Y$ plane, 441 detectors ($21\times 21$) were then evenly distributed between the 25 sources in the selected ROI of $40\times 40{\mathrm{mm}}^2$ to cover both hemispheres [Fig. 1(c)]. With this S–D configuration, the distance between two adjacent detectors was 2 mm. This configuration allowed for arranging 441 detectors on the ROI without any overlap and with balance between the spatial resolution and computation time. Moreover, S–D pairs with signal-to-noise ratios (SNRs) below 3 (linear scale) or above 90% of saturated intensity were excluded from the data analysis. The SNR was calculated by the ratio of $I/I_{\text{noise}}$, where $I$ and $I_{\text{noise}}$ denote the detected light intensity and dark noise of the camera, respectively.

Another processing step was to extract boundary blood flow indices (BFIs) from measured speckle contrasts in the raw intensity images. The speckle contrast relies on the statistical property of spatial signal intensity distribution and was quantified by calculating the intensity ratio of standard deviation ($\sigma$) to mean ($⟨I⟩$), i.e., $K=\sigma /⟨I⟩$ in a pixel window of $7\times 7\text{pixels}$ on the sCMOS camera. We then averaged $K$ values over $3\times 3$ adjacent pixel windows with an area of $0.2\times 0.2{\mathrm{mm}}^2$ (as a single detector out of total 441 detectors) to improve SNRs. These speckle contrast calculations were carried out at the effective S–D distances of 6 to 16 mm, which corresponded to penetration depths up to $\sim 8\mathrm{mm}$ (i.e., one half of the maximal S–D separation). Boundary BFI data were then obtained through a nonlinear relation with $K$. The derivation of this nonlinear relation can be found in our previous publication.⁵² The boundary BFI at one S–D pair was extracted by iteratively minimizing the difference between the measured and theoretical speckle contrast, with a total computation time of $\sim 100\mathrm{ms}$.

2.4.

3D Image Reconstruction of CBF

The S–D pair coordinates, their boundary BFIs, and a solid volume tetrahedral mesh with a 30 mm height, $60\times 60{\mathrm{mm}}^2$ bottom area, and 20k nodes were input into our FEM-based reconstruction program, i.e., modified NIR fluorescence and spectral tomography, to generate a complete 3D image of CBF distributions in the measured tissue volume.^51–53 The relative change in CBF (rCBF) was calculated by normalizing the reconstructed BFI to its corresponding baseline before physiological manipulations.

In the conventional reconstruction method, blood flow distribution in the measured tissue volume is reconstructed once a full source scanning cycle is completed. In other words, the intensity images for all source locations (i.e., $5\times 5$ sources) should be available to generate a complete tomographic flow image [Fig. 2(f)]. Then the process is repeated for subsequent scanning cycles (i.e., another $5\times 5$ sources). As a result, the total scanning time in previous studies using $5\times 5$ sources was 5 s (0.2 Hz). Decreasing the number of source positions would improve the temporal resolution. However, this would undesirably reduce the sampling density and SNRs. Alternatively, we implemented a new reconstruction method based on a moving window, in which a tomographic reconstruction was obtained each time one new source position became available [Fig. 2(g)].⁵⁸ The moving window method updated data at each new source position every 0.2 s. As a result, the sampling rate of the moving window method (5 Hz) was 25 times higher than the conventional reconstruction method (0.2 Hz). Using this higher sampling rate (5 Hz), carrier frequencies from cardiac and respiratory variations are unlikely to be aliased into LFO bands and can be removed with simple low-pass filters. Typical heartbeat and respiratory rates for piglets depend on their age and the anesthetic used. In our study, the heartbeat and respiratory rates for piglets were around 2.5 and 0.5 Hz, respectively.

The reconstruction procedure itself is often an expensive task with a high computational overhead. The higher sampling rate achieved in this study carried a greater burden of reconstructing for all time points. With original reconstruction procedures, processing such a magnified dataset size could prohibit the practical use and longitudinal applications. To improve the speed of reconstructing multiple images, we parallelized the for-loop iterations using MATLAB’s parfor control flow statement available in the Parallel Computing Toolbox™. This was possible because the 3D reconstructions of multiple rCBF images were independent, thus allowing for the execution of iterative calculations simultaneously. Our approach utilized powerful computers at the University of Kentucky High-Performance Computing (HPC) Center to achieve fast parallel computations of Jacobian matrices across multiple S–D pairs and iterative solvers.^51–53 The computer configuration was Intel^® Xeon^® Gold 6252 high-performance server microprocessor, 16 CPU cores at 3.70 GHz (max turbo frequency). For a single 3D image reconstruction of CBF distribution in the slab mesh with total mesh nodes of 20k as used in this study, the sublimated speed and memory efficiency led to 18 times reduction of computation time (from 15 min to 50 s).

2.5.

Power Spectral Density Analysis for LFO Intensity

The PSD analysis for LFO intensity quantification has been established in diffuse optical technologies including NIRS^16,64 and DCS.^7,35,36 Following these methods, LFO intensities of CBF under different physiological conditions (i.e., resting baseline, right CCA occlusion, transient global ischemia, releasing left CCA occlusion, and sequentially releasing right CCA occlusion for recovery) were extracted from PSDs calculated by Welch’s method in a nonparametric approach.³⁵ Specifically, the PSD was calculated as follows: $\mathrm{PSD}=1/(F_s\times N_{\mathrm{seg}})\sum |{\mathrm{FFT}}_{X_{\mathrm{seg}}}(f){|}^2$, where $F_s$ is the sampling frequency, $N_{\mathrm{seg}}$ is the number of segments, and ${\mathrm{FFT}}_{X_{\mathrm{seg}}}$ is the fast Fourier transform of the segmented signals.⁶⁵ Briefly, the time-course rCBF data under each physiological condition at the sampling rate of $\sim 5\mathrm{Hz}$ (frame per second) were first detrended using a built-in “detrend” function in MATLAB to obtain the best straight-line fit linear trend. The detrended rCBF data were then transferred into a Butterworth filter with a passband at the LFO range from 0.0 to 0.1 Hz. Finally, the built-in “pwelch” function in MATLAB was used to generate the PSD over the LFO bandwidth.

To compare PSD levels across different phases of transient global cerebral ischemia, PSD at each phase was normalized to the AUC of the baseline PSD. The trapezoidal integration method in MATLAB built-in “trapz” function was used to calculate the AUC of PSD. The half-width-at-half-max (HWHM) value was then determined from the normalized PSD by quantifying the half width of the frequency band, in which the PSD was above half of its maximum value.

3.1.

rCBF and Vital Variations During Transient Global Ischemia

Figure 3 shows physiological changes measured by the scDCT and multichannel respirator-oximeter during transient global ischemia in two neonatal piglet brains following the CCA occlusion protocol. Figure 3(a) shows rCBF data averaged over the ROI at the depth of 5-mm beneath the skull of piglet #1, quantified by the moving window (blue curve) and conventional (red curve) reconstruction methods, respectively. Based on neonatal piglet head anatomy, the depth of 5 mm reaches the brain cortex.⁴⁸ During right CCA occlusion, rCBF declined, reaching $75\%\pm 11\%$ and $79\%\pm 4\%$ (mean ± standard deviation in timeseries of rCBF values) of the baseline from the moving window and conventional reconstruction methods, respectively. Bilateral ligation further decreased rCBF to $63\%\pm 11\%$ (moving window) and $63\%\pm 6\%$ (conventional reconstruction), respectively. After sequentially releasing left and right ligations, rCBF values were $68\%\pm 11\%$ (moving window)/$70\%\pm 5\%$ (conventional reconstruction) and $70\%\pm 12\%$ (moving window)/$74\%\pm 5\%$ (conventional reconstruction), respectively. At the end of the experiment, euthanasia was induced by $100\%{\mathrm{CO}}_2$ inhalation until the heartbeat stopped. rCBF had an immediate but brief reactive hyperemic response and then dropped to a minimum of 7% (moving window)/9% (conventional reconstruction) by the conclusion. As expected, both the moving window and conventional reconstruction methods portrayed similar trends over all protocol phases.

Fig. 3

Measurement results from piglet #1 and piglet #2 during transient global ischemia. (a), (d) Time-course changes in rCBF measured by scDCT at 5-mm depth during sequential unilateral and bilateral CCA ligations. The experimental protocol included the baseline, right CCA ligation, bilateral ligation (transient global ischemia), releasing of left ligation, releasing of both left and right ligations (recovery), and euthanasia (100% ${\mathrm{CO}}_2$). Red curves indicate the time points at which a 3D CBF reconstruction was available with the conventional reconstruction method, and blue curves indicate the time points with the moving window method. Two zoom-in magnified areas are shown to highlight the increased temporal resolution of the rCBF by the moving window method (blue). Time-course changes in (b) ${\mathrm{EtCO}}_2$, (e) ${\mathrm{SpO}}_2$, (c) RPM, and (f) BPM.

Figure 3(d) shows time-course changes in rCBF detected by scDCT during transient global ischemia in piglet #2. Compared with piglet #1, more pronounced decreases occurred during bilateral ligation compared with unilateral ligation. Specifically, rCBF declined to $90\%\pm 12\%$ (moving window)/$88\%\pm 2\%$ (conventional reconstruction) of the baseline during right CCA ligation and to $45\%\pm 7\%$ (moving window)/$53\%\pm 8\%$ (conventional reconstruction) during bilateral ligation. Then rCBF increased to $68\%\pm 11\%$ (moving window)/$72\%\pm 9\%$ (conventional reconstruction) during release of left ligation and $91\%\pm 11\%$ (moving window)/$90\%\pm 3\%$ (conventional reconstruction) during release of right ligation, respectively. Transient reactive hyperemic responses were detected during the induction of euthanasia by $100\%{\mathrm{CO}}_2$. rCBF dropped to 30% (moving window)/42% (conventional reconstruction) during the euthanasia, which were less dramatic than those observed in piglet #1.

To summarize rCBF results from the two piglets, the model of sequential CCA ligations enabled the prompt manipulation of rCBF. Right CCA ligation induced an instant reduction in rCBF and sequential bilateral occlusions of CCA resulted in a transient global cerebral ischemia. Sequential release of left and right CCAs allowed for the restoration of rCBF. These changes in rCBF are expected and consistent with previous findings in rodents and neonatal piglets during transient global cerebral ischemia.^48,54,57 Although both methods showed similar trends, the moving window method with a higher sampling rate demonstrated larger variations in rCBF, which are likely associated with other physiological changes, such as heartbeat and respiration rates. Importantly, the high sampling rate of the moving window method enabled the detection of LFOs without any interference from heartbeat and respiration.

Fluctuations of ${\mathrm{EtCO}}_2$, ${\mathrm{SpO}}_2$, RPM, and BPM during transient global ischemia are shown in Figs. 3(b) and 3(c) for piglet #1 and Figs. 3(e) and 3(f) for piglet #2. The increased ${\mathrm{CO}}_2$ concentration (${\mathrm{EtCO}}_2$, blue dots in Figs. 3(b) and 3(e)] along with the reduced ${\mathrm{O}}_2$ concentration during euthanasia inhibited respiratory and cardiac functions, thus leading to large drops in rCBF [Figs. 3(a) and 3(d)], ${\mathrm{SpO}}_2$ [red dots in Figs. 3(b) and 3(e)], and BPM [red dots in Figs. 3(c) and 3(f)]. Piglet #2 was found to have a marginal change in ${\mathrm{EtCO}}_2$, which is consistent with the higher observed rCBF minimum upon euthanasia compared with piglet #1. Unfortunately, data were not collected from the respirator-oximeter device during euthanasia of piglet # 2 due to a technical problem with the device, but a BPM decline was confirmed by measuring the heartbeat with a stethoscope. The RPM [blue dots in Figs. 3(c) and 3(f)] response in piglet #1 was negligible, whereas that in piglet #2 displayed a slow decrease throughout. These physiological variations resulting from the hypoxic challenge agreed with pathophysiological processes.

3.2.

LFOs During Transient Global Ischemia

LFOs were isolated from moving window rCBF curves [Figs. 3(a) and 3(d) blue curves] for each session of physiological conditions. The LFO range ($<0.1\mathrm{Hz}$) was identified through the power spectra of rCBF as shown in Fig. 4. rCBF signals were taken from the baseline measurement interval before physiological manipulations. With the new moving window method, scDCT achieved a sampling rate of 5 Hz, which was sufficient for capturing the heartbeat (2.5 Hz) and respiration rate (0.5 Hz). Apparently, high-frequency signals arising from respiration and cardiac cycles were aliased into the desired LFOs. Using a low-pass filter, rCBF signals were reanalyzed for LFO isolation while filtering the artifactual signals unrelated to neural activity.

Fig. 4

Power spectra of rCBF signals during baseline measurements. Power spectra were determined from reconstructed rCBF at 5 mm beneath the skull as acquired using scDCT with the moving window reconstruction method. Baseline measurements preceded physiological manipulations (i.e., CCA ligations). LFOs ($<0.1\mathrm{Hz}$), heartbeat (2.5 Hz), and respiratory rate (0.5 Hz) were captured.

Figure 5 shows the time course for rCBF during the baseline, right CCA occlusion, transient global ischemia, releasing the left CCA occlusion, releasing both left and right CCA occlusions, and euthanasia periods after detrending and LFO bandpass filtering (0.01 to 0.1 Hz). Hemodynamic fluctuations during transient global ischemia and euthanasia periods exhibited marked differences. This positive finding verified that LFOs as measured by the optimized scDCT have the potential for indicating abnormal cerebrovascular conditions.

Fig. 5

Detrended and bandpass filtered time course changes in rCBF over each interval for both piglets. (a) Piglet #1 and (b) piglet #2: time course changes in reanalyzed rCBF, from 5 mm beneath the skull, during the periods of: baseline, right CCA ligation, bilateral ligation (transient global ischemia), releasing of left ligation, releasing of both left and right ligations (recovery), and euthanasia (100% ${\mathrm{CO}}_2$). The detrended and bandpass filtered (0.01 to 0.1 Hz) rCBF represents spontaneous activity that occurred within the piglet’s brain.

3.3.

PSD Changes in LFOs for Identification of Neurological Variations

To better characterize the potential of LFOs as a measure of abnormal neurological status, we calculated changes in PSD within the LFO range. Figures 6(a) and 6(b) show normalized PSD distributions for both piglets. During the period of bilateral ligation, PSD decreased to $\sim 35\%$ of the baseline for both piglet #1 and piglet #2 [Fig. 6(c)]. Correspondingly, HWHM values increased from 0.026 Hz at the baseline to 0.042 Hz (1.6 times increase) and 0.031 Hz (1.2 times increase) in piglet #1 and piglet #2, respectively [Fig. 6(d)]. Alterations in PSD and HWHM are likely associated with impairments of autonomic nervous activities during transient global cerebral ischemia. After release of bilateral ligation, both PSD and HWHM recovered toward their baseline levels.

Fig. 6

PSD analysis of LFOs in both piglets. (a), (b) Normalized PSD distributions from piglet #1 and piglet #2, respectively. PSDs were normalized to the AUC of the baseline PSD. (c), (d) PSD and HWHM changes during different phases of transient global cerebral ischemia from piglet #1 and piglet #2, respectively.