1.IntroductionHemodynamic changes in the brain are closely coupled with neuronal activity by a relationship referred to as neurovascular coupling (NVC).1,2 Studying NVC is important due to its impact on neurological disorders as well as to understand its mechanisms.2,3 Hemodynamic recordings from the brain can be performed through many modalities such as widefield systems, head-mountable microscopes, and fiber-based instruments.4 The latter are minimally invasive and enable deep brain investigations during preclinical and clinical studies.5–9 During these fiber-based recordings, one or more fibers are used to deliver and collect narrowband or broadband light from the sample of interest, e.g., tissue. The collected light signals are then processed to extract the optical properties of the tissue sample.10 Motion artifacts during fiber-based recordings have been reported to give rise to problems such as changes in the transmitted/collected light leading to misinterpretation of data.6,11–13 A limited number of strategies to mitigate motion artifacts from signals acquired from fiber-based instruments have been previously implemented. These include trial rejection, development of instrument design for rejecting motion artifacts, and trial averaging.14 Motion artifact correction (MAC) algorithms have been developed for near-infrared spectroscopy (NIRS) and fiber-based fluorescence measurements.13,15 Motion artifacts in the widely used gCaMP-based fluorescence fiber photometry for monitoring intracellular calcium () are estimated and corrected using to 410 nm light, which is an isosbestic wavelength of the fluorophor that is responsive only to the total concentration and does not vary with concentration.13 This isosbestic channel removes the artifacts present in the 465 nm fluorescence excitation channel. During NIRS measurements, the motion of the scalp relative to the fiber results in motion artifacts.14 Signal processing methods, such as principal component analysis, and spline interpolation, have been established to correct for motion artifacts in NIRS signals.15 However, these techniques are mostly modality/paradigm specific. In this work, we investigate motion artifacts during recordings from a class of optical fiber-based instruments known as single fiber systems (SFS). Development of SFS for diffuse optical measurements has been previously reported in literature.6,8,16–18 Our research group has established the use of an SFS to enable hyperspectral reflectance recordings from deep brain regions of freely moving rodents.6,17 It has been shown to reliably measure changes in oxygen saturation (% ) and blood perfusion changes in the past.6,17,19 Data recorded using SFS showed an increase in variations in the % as a mouse transitioned from the anesthetized state to a freely moving state.6 It was unclear whether the observed variations were part of hemodynamics or due to noise, such as motion artifacts. In this work, first, we identify a wavelength region in the reflectance data acquired by the SFS that is less sensitive to blood absorption-induced changes (i.e., hemodynamics). Using this wavelength region, we developed a dissimilarity metric (DM) to quantify the effect of motion artifacts during animal experiments. Animal experiments were performed specifically to produce episodes of large motion during SFS recordings from a deep brain region known as the periventricular nucleus (PVN) of the hypothalamus. In experiment-1, we transitioned a mouse between two arenas, which is a scenario representative of spontaneous brain activity and large motion artifacts during transfer. Footshocks (FS) are known to increase brain activity in multiple brain regions including the PVN, while leading to jumping behavior due to short periods of panic.20,21 In experiment-2, FS were applied to the animal, which is representative of a scenario of large motion artifacts during FS due to jumping behavior. We used a combination of experimental data and Monte Carlo (MC) simulations to develop a mathematical modelling framework for the motion artifacts in the SFS signals. This framework was then extended to a MAC algorithm that was validated on simulation and experimental data. 2.MethodsThe animal experiments carried out were in accordance with the guidelines of the Canadian Council on Animal Care and were approved by the Institutional Animal Care and Use Committee of the University of Calgary. The optical setup of SFS has lenses - L1, L2, and L3, an aperture (A), and a beamsplitter (BS), as shown in Fig. 1(a).6 The collected light is detected using a spectrometer (S, MayaPro 2000, Ocean Optics). A microcontroller was used to synchronize SFS with behavioral recording cameras and FS generator through infrared LED pulses (peak wavelength of 940 nm) recorded by the spectrometer. Six Ai-148 mice were implanted with custom-made 16-deg angle polished ferrules (Doric Lenses, Canada) at the PVN region of the brain through a stereotactic procedure.6,20,22 Based on our previous study, the 16-deg angle was determined to reduce the back reflection during the coupling of light into the fiber while minimizing the tip sharpness that can cause brain damage.23 During the surgery, mice are first anesthetized and maintained under anesthesia using an isoflurane-oxygen mixture. The region between the ears is shaved and disinfected using 70% alcohol and Betadine. A hole is made in the skull to target the PVN region (anterior–posterior (AP) = −0.7 mm, dorsal–ventral (DV) = −4.7 mm from the dura and medial-lateral (ML) = −0.2 mm from the bregma), following which the implant ferrule is lowered.20,22 Dental cement and meta bond are used to secure the ferrules on the skull of the animal. After post-operative care, animals were then used for recordings. During our experiments, animals were 6 to 8 weeks old at the time of surgery and 12 to 16 weeks old at the time of recordings. The recordings were performed using optical fiber patch cords of core, 0.53 NA, and 2 m in length with 16-deg angle polished ends (Doric Lenses, Canada). A cross section of the brain where light is interacting with the tissue is shown in Fig. 1(b). Prior to recording in all experiments, the animals were left to habituate in the recording room for at least 20 min. First, spectra were measured using the spectrometer with the fiber in air () and in glycerine (). Then the optical fiber patch cord was connected to the implant ferrule using a ceramic sleeve, and spectra were measured from the brain (). refers to the backscatter from the tissue. Before connecting the patch cord to the implant ferrule, a drop of glycerine is used as an index-matching medium to maximize light coupling into the implant ferrule. During recording spectra, we set the typical exposure times (e.g., 20 to 25 ms) and number of spectra averaged of the spectrometer (e.g., 4 to 5) to maintain a sampling rate of about 10 Hz. After recording from an animal, the fiber tip is cleaned using water. and are then remeasured before recording from next animal. The measured reflectance spectrum () is then calculated as Similar to previous work from our group, the spectra were fit to an empirical model of photon transport in tissue to obtain % .6,17,19,24–26 Briefly, the method involves fitting the using a non-linear least square fitting approach in MATLAB to extract % , which is one of the parameters of the model. During the calculation, we select and which results in a lower fitting error. All scripts used for data processing have been provided in the online repository mentioned in Sec. 7. 2.1.Experiment-1Certain experiments in neuroscience involve transitioning animals from one arena to another, where large motion artifacts can be possible due to the transfer process between arenas.20 We recorded hemodynamics data while we transitioned animals from arena-1 after 5 min to arena-2 for 5 min and back to arena-1 for a duration of 5 min. Three animals were used for this experiment. 2.2.Experiment-2Footshocks (FS) are applied to mice in experiments where stress-based behavior is studied or in fear conditioning studies.21 Earlier studies have shown robust responses recorded during fluorescence photometry from the periventricular nucleus (PVN region) of the hypothalamus.20,27 During an FS, the mouse jumps aggressively, which causes movement of the optical fiber patch cord (F), which creates a situation where large effect of motion artifact can be speculated. FS were generated using ANY-Maze system and were synchronized to the SFS recording using the IR LED [shown in Fig. 1(a)]. Ten FS of 2-s long duration were delivered to the animal in the foot-shock cage with a 30-s inter-stimulus period. An FS trial was defined as 8 s prior to the onset of FS to 22 s post onset of FS. This resulted in a total trial duration of 30 s. Three animals were used for this experiment. 2.3.Monte Carlo SimulationTo mimic light-tissue interaction in the SFS, we developed a simulation model in the non-sequential simulation mode of Zemax OpticStudio. In the simulation model, we omitted certain elements that remain static and do not interfere with the signals we are concerned with in these simulations. A rectangular volume was used to simulate the brain tissue. An optical geometry consisting of a source and detector was set up to couple and detect light respectively from the brain tissue [as shown in Fig. 1(c)]. The far end of the optical fiber had a 16 deg end face and was embedded in the tissue to mimic the experimental setting closely (See Sec. 2).23 A light source was configured to launch photons into the tissue through the optical fiber. A justification for the number of photons is included in the supplementary data (Fig. S1 in the Supplementary Material). The refractive index of the brain tissue and optical fiber was modelled using Cauchy’s formula using data from the literature.28,29 Tissue absorption of the brain tissue was modeled using the equation , where , are the absorption coefficients of oxyhemoglobin and deoxyhemoglobin, and BV is the blood volume fraction.10 Scattering was modelled using an intralipid scattering model with the Henyey–Greenstein (HG) scattering phase function.30,31 The HG function is given by where represents the probability of scattering in the direction of the scattering angle () and is the anisotropy factor.30 We modelled the anisotropy factor () of intralipid using the equation leading to values between for the wavelength range .31 Using these values in the HG phase function reveals forward scattering characteristics, which agree with the scattering properties of biological tissues.32 The simulation model, along with a sample photon trajectory, is shown in Fig. 1(c). During simulations, the % and BV were varied to calculate different absorption coefficients. The absorption coefficient and the scattering coefficients were used to calculate the mean free path (MFP) and transmission parameter, which along with are inputs to the Zemax OpticStudio software. Changes in optical transmission during motion artifact simulations were modeled by forcing a transmission coating at the fiber/brain interface. A macro written in Zemax programming language was used to automate the simulation process. 2.3.1.Static Monte Carlo simulationsTen MC simulations were run for constant tissue properties (35% , and ) for various wavelengths (specified in Secs. 3.1 and 3.3).6,33 Simulations were initially run to model the light collected with the fiber in air () and glycerine () followed by in tissue (). Data generated in Zemax Opticstudio were processed in MATLAB to yield simulated using the following equation: The calculated was used to calculate % using the empirical model as described in Sec. 2. Prior to fitting, we additionally interpolated the simulated data in MATLAB to improve the fitting process.25 2.3.2.Dynamic Monte Carlo simulationsDuring dynamic hemodynamic simulation, we modelled blood absorption changes across time at two different wavelengths (i.e., 680 and 545 nm). The 545 nm wavelength was chosen to evaluate changes due to hemodynamics due to its relatively large sensitivity to hemoglobin absorption. The 680 nm wavelength was used to observe effects due to non-hemodynamic changes such as motion artifacts due to its lower sensitivity to hemoglobin absorption. We performed simulation with and without including motion artifacts. Hemodynamic signals have been previously modeled using hemodynamic response functions (HRFs).34 We used definitions provided by Kamran et al. to generate the dynamics of the temporal variation of the signal.35 The hemodynamic response function [HRF(t), Fig. S2(a) in the Supplementary Material] was produced by convolving a canonical hemodynamic impulse response function [] with the stimulus function [, Fig. S2(a) in the Supplementary Material], i.e., . The canonical hemodynamic response function is described as where represents the time points; , are parameters; and refers to the gamma function. Next, the final physiological hemodynamic response () is written as a weighted superposition of the , sinusoidal terms and a DC bias () as described in Eq. (4). The sinusoidal terms in Eq. (4) are representative of cardiac pulsations, respiration, and Mayer waves with frequencies , , and . We used , , and for our simulations [Fig. S2(b) in the Supplementary Material].35 The terms , and , and are weighting factors. In our work, we have used to model the changes in light intensity recorded by the SFS (i.e., )The normalized , represented as [shown in Fig. S2(c) in the Supplementary Material], is appropriately scaled to simulate changes in the concentration of oxy () and deoxyhemoglobin (), as shown in Fig. S2(d) in the Supplementary Material. Following this, the absorption coefficient of tissue is calculated as where is the concentration of total hemoglobin. We used , and and for the simulations. These concentrations are in agreement with experimentally observed values.36,37 at 545 nm and 680 nm was calculated using MC simulations at a time resolution of 10 Hz with a stimulus () from 2 to 4 s. To simulate motion artifacts, a second series of simulations were done with a 96% transmission coating that was applied on the fiber/brain interface in the model between 2 and 3.5 s.2.4.Quantifying Changes in Signals2.4.1.Quantifying shifts in reflectance spectra () between experimental data and simulationWe quantified the offset-like changes observed in the using root mean square deviation (RMSD) given by the following equation: where is the RMSD calculated for a test curve () with respect to a reference curve (), and is the number of points in .2.4.2.Dissimilarity metricTo quantify the magnitude of changes seen in due to motion artifacts, we developed a DM25 defined by the following equation: DM represents the normalized instantaneous difference in the intensity of collected light () in a wavelength range from to from the average intensity over a window of length centered around a sample . We set and which represents the hemodynamic insensitive region (Sec. 3.1). We then calculate the standard deviation of DM () to evaluate the effect of motion artifacts. A higher corresponds to a higher variation in the wavelength insensitive region of , indicating a higher contribution from motion artifacts. 2.4.3.Change in backscatter intensity ()545 nm wavelength is an isosbestic wavelength of hemoglobin with a relatively large sensitivity to hemoglobin due to the large absorption coefficient of hemoglobin. We refer to the measured intensity of light backscattered from the brain as the backscatter intensity (). We calculated the % change in at 545 nm relative to the mean of the baseline period as an indicator of blood perfusion for analyzing both experimental and simulation data.19 For experiment 2 (Sec. 2.1), the baseline period was defined as the 10 s prior to the onset of FS, while for the dynamic simulations (Sec. 2.3.2), it was defined as the 2 s period prior to stimulus onset. 3.ResultsWe first identify a wavelength region that is relatively independent of hemodynamics (Sec. 3.1) and then quantify the relative magnitude of motion artifacts (Sec. 3.2). We next simulate motion artifacts using MC simulations to understand experimental findings (Sec. 3.3). Next, a model for light propagation through SFS and the tissue is proposed and extended to a MAC algorithm (Sec. 3.4). The MAC is first validated using simulation data (Sec. 3.5), then its application is illustrated in experimental data (Sec. 3.6). 3.1.Identifying a Hemodynamics-Independent Wavelength RegionThe spectra of light collected through SFS are dependent on the incident spectra, tissue properties (such as scattering and absorption), and motion artifacts. Hemodynamic changes are strongly affected by blood absorption changes, which are wavelength-dependent. Motion artifacts can cause wavelength-independent changes in contrast to wavelength-specific changes during hemodynamics.6,25,36 To quantify motion artifacts, the objective in this section was to identify a wavelength region where we collected a reasonable amount of light that was minimally sensitive to blood absorption changes. We performed 10 static MC simulation runs for two cases—(i) 1% BV, 35% with intralipid scattering model, and (ii) only intralipid scattering model (i.e., ) using the simulation model described in Sec. 2.3.1 for 55 wavelengths between 450 and 790 nm. We used photons for the MC simulations to minimize the inherent MC simulation noise [Fig. S1(a) in the Supplementary Material]. We averaged the 10 simulation runs to further reduce the noise. Figure 2(a) shows each wavelength’s sensitivity to hemodynamics, defined as the absolute % change in the light collected between the scattering with absorption and the scattering-only cases. The figure also shows the normalized spectrum of the light incident on tissue with our experimental setup. Our simulations showed that while wavelengths beyond 700 nm showed the least sensitivity to changes in blood absorption, we collect very little light beyond 700 nm [Fig. 2(a)]. As a tradeoff between the quality of the signal and its independence from hemodynamics, we chose the 670 to 680 nm region as largely representing the motion artifacts. 3.2.Motion Artifacts in the Measured Backscatter Intensity () SignalsIn this section, we analyze the signals acquired during experiment-1. Figure 2(b) shows the signals at 545 nm for a single animal during experiment-1. We observe that the shows large, abrupt changes at multiple time points. To investigate these changes further, Figs. 2(c) and 2(d) show and across the entire spectrum for the three-time points marked in Fig. 2(b) respectively. We observed wavelength-wide, offset-like changes in during the abrupt changes. Since hemodynamics is a slow process and causes wavelength-specific changes, the rapid broadband changes observed in Figs. 2(b)–2(d) are likely artifactual. We quantified the offset-like changes observed in Fig. 2(d) using RMSD as described in Sec. 2.4. To calculate RMSD, wavelengths ranging from , with 1 nm spacing were used, giving a total points. The RMSD of shown in blue and pink with respect to the one in orange is 0.19 and 0.41, respectively. Figure 2(d) shows the corresponding % extracted for each of the using the empirical model as described in Sec. 2. Given that the % depends on the spectrum’s shape, the % extracted is expected to be relatively immune to motion artifacts.38 The immunity of % to transmission changes has been further addressed in Sec. 3.3. Therefore, in this work, we focus on understanding and mitigating motion artifacts in the perfusion signal. 3.2.1.Using dissimilarity metric to quantify motion artifactsTo illustrate the utility of DM (described in Sec. 2.4.2) in quantifying motion artifacts, we partitioned the data using a sliding window into 1-s chunks, which were then used to calculate DM values. The standard deviation of DM () over those chunks is plotted in Fig. 2(e). It is observed that peaks when there is an abrupt change in the signal as shown in Fig. 2(b). This indicates that can be used as a reliable measure of changes in the signals. Next, we used DM to evaluate the extent of motion artifact during experiment-2. FS in freely moving mice induce short episodes of panic and elicit jumping behavior. The thirty FS trials (3 mice, 10 trials per mouse) were divided into two data groups. The pre-stim condition consisted of data for the first 8 s of the trial, where we expect minimal motion artifact. The stim condition consisted of data 8 s post onset of the FS where large motion artifacts are possible. DM was calculated for the two conditions followed by the calculation of . Figure 2(f) shows the for both conditions displayed in a pair-wise plot. We observed that there was a statistically significant difference (, using Wilcoxon signed rank test) in the amount of motion artifact between the two conditions with larger artifacts after the onset of the FS (i.e., for the stim condition). FS in mice have been shown to increase blood flow in the PVN due to an increase in the activity of PVN neurons of the hypothalamus (unpublished data20,22). We calculated the % change in backscatter intensity based on Sec. 2.4.3. This was used to calculate the average % backscatter change in 300 ms windows prior to the start of and just before the end of FS to represent the pre-stim and stim periods, respectively. Since an increase in perfusion leads to a decrease in backscatter intensity, Fig. 2(g) illustrates the negative of % change in backscatter intensity. We observed an increase in the perfusion response (, using Wilcoxon signed rank test) in the stim condition compared to pre-stim. Since motion artifacts impact [Figs. 2(b)–2(d)], the increase observed in the perfusion in the stim condition can be considered a superposition of effects due to true hemodynamics and motion artifacts. 3.3.Simulating Motion Artifacts Using Monte Carlo SimulationsFrom Sec. 3.2, we observed that motion artifacts caused changes in . As described in Sec. 3.2, such motion-induced changes in optical signals can be modelled as changes in the transmission of the light path.25 Changes in transmission can occur at optical interfaces in the lightpath that are susceptible to motion, e.g., the fiber implant ferrule interface, fiber-brain interface, or movement of the optical fiber. In this section, we present MC simulations with transient changes in the lightpath transmission to model variations in , , and % that can occur during motion artifacts. The simulated tissue back-scatter is represented as , in contrast to the measured tissue backscatter that is referred to as . With photons in the simulation, Fig. S1(b) in the Supplementary Material shows that the hemodynamic features (e.g., hemoglobin absorption valleys at and ) in the simulated were distinctly visible. This shows that the simulated hemodynamic features were well over the noise due to inherent randomness in MC simulations which is important for % estimation. The tissue oxygenation and blood volume were set to 35 % and 1 % BV, respectively, while other simulation properties were the same as described in Sec. 2.3,33 While the transmission can change anywhere along the lightpath, we simulated changes by adding a transmission coating at the fiber-tissue interface. To stay consistent with how is calculated, we first simulated the nominal back-reflection in glycerine () and air (). Next, MC simulations to extract the simulated tissue backscatter were carried out at coating transmission values 96%, 92%, 90%, and 80%. The nominal case representing no artifacts was simulated without any coating. Resulting based on calculations described in Sec. 2.3.1 are shown in Fig. 3(a). To quantify the relative shift in the spectrum, we used the RMSD calculation described in Sec. 2.4.1 using for nominal transmission as the reference. Figure 3(b) shows the boxplot for the % calculated for various transmissions based on methods described in Sec. 2.3.1. The height of the boxplot is representative of the interquartile range (IQR), which is a measure of the spread of the data. As the fitting process is sensitive to noise, each simulation was repeated 10 times resulting in 10 values of % for each transmission. The small variations in the IQR can be attributed to the fitting noise and inherent variability during MC simulations [Fig. S1(b) in the Supplementary Material]. To test if the changes in the transmission affect the % among all the different transmission levels, we performed non-parametric one-way analysis of variance (ANOVA) using the Kruskal–Wallis test. We observed a -value of 0.44 implying that the % for various transmissions was not statistically significantly different from each other. The offset-like shifts observed in the simulated reflectance spectra shown in Fig. 3(a) during transmission changes were also observed in the measured reflectance spectra shown in Fig. 2(d). While there are certain differences in the characteristic hemodynamic features (e.g., shape of the curves and slope), the focus here is to characterize large non-hemodynamic offset like changes. Therefore, to assess the offset-like shifts with respect to the nominal transmission quantitatively, we calculated RMSD [shown in Fig. 3(a)]. The RMSD observed in the experiment (i.e., 0.19 and 0.41, described in Sec. 3.2) is comparable to what we observed in Fig. 3(a) (i.e., 0.21 to 1.1). This provides evidence to support the idea that transmission changes in the lightpath could be used to model offset-like changes during motion artifacts. 3.4.Motion Artifact Correction Algorithm for Reflectance SignalsIn this section, we propose an MAC algorithm to correct changes in the perfusion signal due to motion artifacts caused by transmission changes. Equation (8) describes the intensity of light measured by SFS from the brain () during recording that accounts for Fresnel reflection at the optical interfaces and tissue backscatter. A detailed explanation of the model is provided elsewhere25 In Eq. (8), the incident spectra just before coupling to the fiber is denoted by . and refer to the transmission and reflection coefficients, respectively, as light goes from medium to at location . The media is indicated by superscripts: , , refer to air, fiber, and brain. The locations are indicated by subscripts: and refer to the near end (where light enters the fiber) and the far end (where light enters the brain) of the fiber, respectively. is the fiber transmission. is the actual reflectance spectrum of the tissue. Further, is the light intensity entering the brain at the far end. We showed in Sec. 3.1 that at 670-680 nm wavelength has a little contribution from hemodynamics and therefore largely depends on changes due to Fresnel reflection. This implies that, the term in Eq. (8) contributes negligibly as . Further, in Sec. 3.2 we have uncovered that changes observed in the 670 to 680 nm can be used as a proxy for motion artifacts. To simplify calculations, we use 680 nm to represent the range. Therefore, the light collected at 680 nm becomes can be used to approximate all transmission-based changes (i.e., motion artifacts) in the recorded artifactual backscatter signal (). Therefore, using Eqs. (8) and (9), we can write , where is a wavelength dependent factor and is the back-scatter signal of interest. The parameter is important to scale wavelength-dependent changes observed in the wavelength to the at other wavelengths. We calculate the corrected (i.e., ) as where and represents the mean of the time series. The parameter is extracted using a linear regression-based approach.25 3.5.Motion Artifact Correction Validation through Simulation DataTo further validate the MAC algorithm we simulated the changes in the tissue backscatter intensity () measured by SFS at two wavelengths, i.e., 545 and 680 nm as a function of tissue absorption as described in Sec 2.3.2. at 545 nm represents the perfusion signal, and that at 680 nm represents the correction signal. Figure S1(c) in the Supplementary Material shows that the photons used for the simulation were sufficient to produce a change in the at 545 nm which was well over the noise due to inherent randomness during MC simulations. There was no significant change seen in at 680 nm due to its insensitivity to hemodynamics as per the proposed MC simulation model (Sec. 2.3). For Fig. 4, the was processed to calculate the % change in backscatter intensity with respect to baseline period prior to stimulus-evoked change (Sec. 2.4.3). Figure 4(a) shows the results from the simulation without motion artifacts (MA). As expected, the 680 nm, signal has no features and correcting the perfusion signal does not change it. Figure 4(b) shows the case with artifacts. Large variations can be seen in both channels and using the correction described in Sec. 3.4 recovers the perfusion measurement. Figure 4(c) compares the ideal trace with the corrected artifactual trace. We performed mean squared error (MSE) calculations to evaluate the effectiveness of the MAC algorithm. We calculated the MSE between the traces shown in the bottom panel of Figs. 4(a)–4(c). This was done by segregating the plots to three regions – (i) region-1 ( to 1.8 s, prior to motion artifact), (ii) region-2 ( to 3.6 s, during motion artifact) and (iii) region-3 ( to 7 s, after motion artifacts) as shown in Fig. 4(c). We observed a small MSE for the region-1 and region-3, as shown in Figs. 4(a)–4(c), indicating that the correction does not significantly change signal which is free of artifacts. During motion artifact [region-2, Fig. 4(b)], the MSE increased by a factor of about 4 to 6 orders of magnitude in comparison to the MSE values in regions-1 and 2, highlighting the large effect of motion. The MSE was reduced by about 2 orders of magnitude after MAC for region 2 in Fig. 4(c). 3.6.Applications in Experimental DataWe first show the application of MAC for data from experiment 1, followed by experiment 2. Figure 5(a) shows the raw () and corrected () backscatter intensity at 545 nm. As described earlier in Sec. 3.2, we observed that motion artifacts cause abrupt changes in the as shown in Fig. 5(a) that are minimized following correction in . Figure 5(b) shows the spectrum of the collected light post-correction (i.e., ) where the shifts observed in the spectrum (i.e., ) in Fig. 2(c) prior to correction have been minimized. Further, the and time series were broken down to 15 60-s-long chunks for each animal. The 45 chunks were used to compute a series of DM values. Figure 5(c) shows the calculations prior to correction, i.e., pre-MAC, and after correction, i.e., post-MAC. We observed that the post-MAC was lower than the pre-MAC . A Wilcoxon signed rank was used to compare the pre-MAC vs post-MAC data, which revealed a statistically significant difference exists between the two (). This shows that the MAC is successfully able to minimize the variations present in the backscatter signal. Next, we applied MAC for experiment-2. Figure 5(d) shows % change in backscatter intensity calculated relative to the 10 s baseline period for a single FS trial. We observe that the peak % change in perfusion is less for the corrected trace (in red) than the raw trace (in blue). This decrease observed during correction can be attributed to the combination of (i) a decrease in motion artifacts due to the correction, and (ii) attenuation in the perfusion signal due to correction using 680 nm wavelength, which is not completely insensitive to absorption due to hemodynamics, i.e., due to the non-zero absorption coefficient of hemoglobin at 680 nm. Identical to the approach in Fig. 2(f), in Fig. 5(e) we used the to calculate for pre-stim and stim data. We observed that both conditions led to comparable . A comparison using a paired Wilcoxon signed rank test showed that the difference is not statistically significant (). Figure 5(f) shows the % change in perfusion between pre-stim and stim conditions calculated after correction using similar to the comparison of the raw data in Fig. 2(g). Post correction plot in Fig. 5(f) shows that the two are still significantly different (, using Wilcoxon signed-rank test), indicating that the correction does not obscure the effect of stimulation. 4.DiscussionAn understanding of motion artifacts is critical to discriminate between signal and noise. Motion artifact mitigation is an active area of research across the biomedical sciences, including light-based techniques. However, knowledge of motion artifacts and how they impact signals during fiber-based measurements, especially SFS is limited and is the focus of our work. The large rapid changes observed during our experiments [e.g., in Figs. 2(b)–2(d)] have been observed by researchers during optical probe-based diffuse optical imaging, and fiber-based NIRS, and have been recognized as motion artifacts.11,12,14,15 There can be several reasons for the shifts observed in the collected light () and the calculated reflectance () such as those shown in Figs. 2(c) and 2(d). First, movement at the fiber patch cord ferrule/implant ferrule or the fiber/brain interface can lead to changes in the optical coupling at the interface. This will change both the light exiting and entering the fiber, leading to artifacts. Changes in pressure on the tissue due to motion can cause alterations to local blood flow leading to changes in and .16,39 Further, movement of the fiber has been known to induce artifacts which can be an additional factor causing changes in light transmission and collection.13 The first objective of this work was to discriminate between true hemodynamics signals and motion artifacts. This was achieved by finding a wavelength region in the collected light, which is relatively independent of tissue backscatter and more sensitive to motion artifacts. The absorption coefficient of whole blood is relatively higher at than beyond 650 nm.36 The scattering coefficient of tissue also decreases beyond 650 nm.10 A decrease in blood absorption and scattering results in higher light penetration into the tissue for longer wavelengths of light.40 Therefore, at longer wavelengths, the relative contribution from changes in Fresnel reflections (e.g., during motion artifacts) at the optical interfaces will be higher than from tissue backscatter in the SFS signal. This is supported by our simulation data where we saw light from 670-680 nm to be about 28 times less sensitive to blood absorption than that at 545 nm (Sec. 3.1). Further, we used this 670 to 680 nm wavelength to quantify the relative degree of motion artifacts during two different experiments. This was achieved by calculating as in Sec. 3.2.1. We performed MC simulations, which revealed that transmission-based changes indeed cause baseline shifts in the collected backscatter () and calculated reflectance (). The RMSD calculations (Sec. 3.3) show similarity in the magnitude of shifts between simulation and experiment. We further show in Sec. 3.3 that the extracted % is relatively immune to such artifacts in due to its dependence on the shape of the reflectance spectrum.38 However, perfusion measurements are susceptible to these artifacts. A mathematical model for light transmission in the SFS is presented in Sec. 3.4 and extended to an MAC algorithm. During MAC, is used to correct for motion artifacts using a linear regression-based method. The MAC algorithm was first validated using dynamic simulation data in Sec. 3.5. The advantage of the simulation is that it enables the computation of data without artifacts which is difficult to do with confidence during experiments. Simulations showed optimal performance by the 680 nm wavelength in reducing the artifactual changes. Further, we used mean square error (MSE) to evaluate MAC performance. While we observed a reduction of MSE after the correction of motion artifacts by two orders of magnitude, the minor variations observed in the MSE between Figs. 4(a)–4(c) in the artifact-free region-1 and region-3 can be attributed to the variations produced during MC simulations as supported in the supplementary data [Fig. S1(c) in the Supplementary Material]. To show the importance of selecting a hemodynamic insensitive wavelength, we have provided data where we test the MAC algorithm with an arbitrarily picked 590 nm wavelength. Results show suboptimal performance for MAC using the 590 nm wavelength due to its increased sensitivity to blood absorption compared to that of 680 nm (Fig. S3 in the Supplementary Material). Testing the MAC algorithm using in-vivo data showed a dramatic reduction of motion artifact-based changes as seen comparing Figs. 2(b) and 5(a). The DM [e.g., as seen by higher in Fig. 2(f)] showed a higher variability prior to correction, which was minimized post-correction as seen in Figs. 5(c) and 5(e). While the MAC algorithm was successful in reducing the variation in the signals, it also appears that it reduces the perfusion signal [Figs. 5(d) and 5(f)]. Prior to the correction, the perfusion during pre-stim condition was significantly lower than that during stim [Fig. 2(g)], but the relative contributions of hemodynamics and motion artifacts to this change were unclear. After correction, there is still a significant difference [Figs. 5(d) and 5(f)] which can be attributed substantially to hemodynamics as motion artifact-based changes are minimized by the MAC algorithm. The corrected trace in Fig. 5(a) shows slow changes, which are characteristic of hemodynamic signals. 4.1.LimitationsThere are some limitations to this study. Firstly, the differences in the hemodynamic features between simulations [Fig. 3(a)] and experiment [Fig. 2(d)] can be attributed to the approximations used in the MC simulation to model the in-vivo situation. For example, tissue was assumed to be homogeneous and the data to estimate the various MC model parameters (e.g., BV, scattering, and absorption properties) from deep brain regions is limited, which can result in the observed differences in the hemodynamic features. Future studies could be aimed at identifying more reliable estimates of various MC model parameters for SFS-based measurements. Second, the small sample size of three animals per experiment, is a limitation of our study. However, to emphasize results, we performed statistical tests on the data by partitioning it as described in Secs. 3.2.1 and 3.6. A similar approach has been used previously to compare peak versus baseline activity in the PVN.22 The Wilcoxon sign test used to test significance in Figs. 2(f)–2(g), Figs. 5(c), and 5(e)–5(f) models the different samples as independent. This assumption is reasonable in our data due to the random and uncorrelated nature of the variations due to motion artifacts in the data. Thirdly, due to the non-zero sensitivity of 680 nm light to hemodynamics, the MAC process leads to some attenuation of the corrected trace. For example, in Fig. 5(d), we observe that the raw trace shows a larger perfusion level than the corrected trace. The faster return to the baseline of the corrected trace than the raw trace can also be attributed to the correction. More studies will be required to study the difference in the temporal features between the corrected and raw traces. While our setup was limited to significant light ( of maximum) from 425 to 710 nm, other wavelengths may work better. Future work aims to identify further longer wavelengths in the near-infrared regions, which have lower sensitivity to hemodynamics than 680 nm to further improve the MAC process. Nevertheless, the methods presented in this paper can be used to visualize and analyze corrected signals, provided we apply identical MAC to all signals being analyzed. Certain spike-like artifacts that were present in Fig. 2(b) are minimized but not completely removed after applying the MAC algorithm. Testing revealed that these abrupt spikes in the post-MAC are due to outliers limiting the regression-based MAC to minimize such spikes. Future work can investigate constrained optimization using a smoothing cost function and interpolation-based smoothing as an additional step for fine-tuning the processed signal after applying the MAC algorithm described in Sec. 3.4. 5.ConclusionsIn conclusion, we found that the wavelength region from 670 to 680 nm of the reflectance data is a good indicator for the amount of motion artifacts present in the hemodynamic signals acquired by the SFS. We observed that the metric can be used to effectively represent changes in the acquired signals during episodes of large speculated motion. This work, as supported with the literature reveals the fact that the motion artifacts can be modeled as transmission changes in the illumination/collection path of light. Such transmission changes can cause shifts in the reflectance data in the experiment. This observation is in agreement with the MC simulations where we varied the transmission in the optical path of light. Secondly, the mathematical modeling of light intensity has been done at various optical interfaces during light propagation through the optical fiber. This framework was extended to a MAC algorithm. The MAC algorithm has been developed and tested using data generated from MC simulations. It is observed that the MAC can minimize motion artifact-based changes while conserving the stimulus-evoked hemodynamic change. Finally, we have tested the MAC algorithm on experimental data where we observed a decrease in the metric post-correction augmenting the importance of the proposed MAC algorithm. We envisage that the methods to address and account for motion artifacts that have been presented in this work will act as stepping stone for future studies dealing with motion artifacts for SFS and other optical probe-based modalities. Code and Data AvailabilityDetails about materials and methods are given in the paper. All code and data in support of the findings of this paper are available at https://figshare.com/s/9f0fe4b39272bb315a87. AcknowledgmentsAnupam Bisht acknowledges Dr. Zelma Kiss and Dr. Patrick Whelan for their support and comments during the initial part of the work. We also thank Govind Peringod, Dr. Ning Cheng, and Dr. Grant Gordon for help during preliminary data collection. We thank the Natural Sciences and Engineering Research Council (NSERC) of Canada (Grant No. NSERC RGPIN/04323) and CMC Microsystems for the funding. Anupam was partially supported by scholarships from the Department of Biomedical Engineering, University of Calgary, and by the Alberta Graduate Excellence Scholarship. We also thank the Hotchkiss brain institute (HBI), and the Cumming School of Medicine (CSM) optogenetics facility and their staff for their help and support during the experiment. ReferencesC. Iadecola,
“The neurovascular unit coming of age: a journey through neurovascular coupling in health and disease,”
Neuron, 96
(1), 17
–42 https://doi.org/10.1016/j.neuron.2017.07.030 NERNET 0896-6273
(2017).
Google Scholar
Y.-R. Gao et al.,
“Time to wake up: Studying neurovascular coupling and brain-wide circuit function in the un-anesthetized animal,”
NeuroImage, 153 382
–398 https://doi.org/10.1016/j.neuroimage.2016.11.069 NEIMEF 1053-8119
(2017).
Google Scholar
H. L. McConnell et al.,
“Astrocyte dysfunction and neurovascular impairment in neurological disorders: correlation or causation?,”
Neurochem. Int., 128 70
–84 https://doi.org/10.1016/j.neuint.2019.04.005 NEUIDS 0197-0186
(2019).
Google Scholar
H. Yu et al.,
“Miniaturized optical neuroimaging in unrestrained animals,”
NeuroImage, 113 397
–406 https://doi.org/10.1016/j.neuroimage.2015.02.070 NEIMEF 1053-8119
(2015).
Google Scholar
P. Rejmstad, P. Zsigmond and K. Wårdell,
“Oxygen saturation estimation in brain tissue using diffuse reflectance spectroscopy along stereotactic trajectories,”
Opt. Express, 25
(7), 8192 https://doi.org/10.1364/OE.25.008192 OPEXFF 1094-4087
(2017).
Google Scholar
L. Yu et al.,
“Fiber photometry for monitoring cerebral oxygen saturation in freely-moving rodents,”
Biomed. Opt. Express, 11
(7), 3491 https://doi.org/10.1364/BOE.393295 BOEICL 2156-7085
(2020).
Google Scholar
S. Skyrman et al.,
“Diffuse reflectance spectroscopy sensor to differentiate between glial tumor and healthy brain tissue: a proof-of-concept study,”
Biomed. Opt. Express, 13
(12), 6470 https://doi.org/10.1364/BOE.474344 BOEICL 2156-7085
(2022).
Google Scholar
F. van Leeuwen–van Zaane et al.,
“In vivo quantification of the scattering properties of tissue using multi-diameter single fiber reflectance spectroscopy,”
Biomed. Opt. Express, 4
(5), 696 https://doi.org/10.1364/BOE.4.000696 BOEICL 2156-7085
(2013).
Google Scholar
X. U. Zhang et al.,
“Refractive index measurement using single fiber reflectance spectroscopy,”
J. Biophotonics, 12
(7), 1
–11 https://doi.org/10.1002/jbio.201900019
(2019).
Google Scholar
S. L. Jacques,
“Optical properties of biological tissues: a review,”
Phys. Med. Biol., 58
(11), R37
–R61 https://doi.org/10.1088/0031-9155/58/11/R37 PHMBA7 0031-9155
(2013).
Google Scholar
W.-C. Lin et al.,
“Diffuse reflectance spectroscopy for in vivo pediatric brain tumor detection,”
J. Biomed. Opt., 15
(6), 061709 https://doi.org/10.1117/1.3505012 JBOPFO 1083-3668
(2010).
Google Scholar
N. Yadav et al.,
“Evaluating and improving the quality of time-dependent, diffuse reflectance spectroscopic signals measured from in vivo brain during craniotomy,”
Med. Eng. Phys., 35
(11), 1551
–1557 https://doi.org/10.1016/j.medengphy.2013.04.006 MEPHEO 1350-4533
(2013).
Google Scholar
E. Martianova, S. Aronson and C. D. Proulx,
“Multi-fiber photometry to record neural activity in freely-moving animals,”
J. Visual. Exp., 2019
(152), 1
–9 https://doi.org/10.3791/60278
(2019).
Google Scholar
R. Huang et al.,
“Motion artifacts removal and evaluation techniques for functional near-infrared spectroscopy signals: a review,”
Front. Neurosci., 16 878750 https://doi.org/10.3389/fnins.2022.878750 1662-453X
(2022).
Google Scholar
R. J. Cooper et al.,
“A systematic comparison of motion artifact correction techniques for functional near-infrared spectroscopy,”
Front. Neurosci., 6 1
–10 https://doi.org/10.3389/fnins.2012.00147 1662-453X
(2012).
Google Scholar
X. U. Zhang et al.,
“Effect of probe pressure on skin tissue optical properties measurement using multi-diameter single fiber reflectance spectroscopy,”
J. Phys. Photonics, 2
(3), 034008 https://doi.org/10.1088/2515-7647/ab9071
(2020).
Google Scholar
L. Yu et al.,
“In-vivo monitoring of tissue oxygen saturation in deep brain structures using a single fiber optical system,”
Biomed. Opt. Express, 7
(11), 4685 https://doi.org/10.1364/BOE.7.004685 BOEICL 2156-7085
(2016).
Google Scholar
W. T. Zhang et al.,
“Spectral fiber photometry derives hemoglobin concentration changes for accurate measurement of fluorescent sensor activity,”
Cell Rep. Methods, 2
(7), 100243 https://doi.org/10.1016/j.crmeth.2022.100243
(2022).
Google Scholar
L. Yu et al.,
“Monitoring stimulus-evoked hemodynamic response during deep brain stimulation with single fiber spectroscopy,”
J. Biophotonics, 15
(11), 1
–11 https://doi.org/10.1002/jbio.202200076
(2022).
Google Scholar
T. Füzesi et al.,
“Hypothalamic CRH neurons orchestrate complex behaviours after stress,”
Nat. Commun., 7 1
–14 https://doi.org/10.1038/ncomms11937 NCAOBW 2041-1723
(2016).
Google Scholar
S. Maren and M. S. Fanselow,
“The amygdala and fear conditioning: has the nut been cracked?,”
Neuron, 16
(2), 237
–240 https://doi.org/10.1016/S0896-6273(00)80041-0 NERNET 0896-6273
(1996).
Google Scholar
N. Daviu et al.,
“Paraventricular nucleus CRH neurons encode stress controllability and regulate defensive behavior selection,”
Nat. Neurosci., 23
(3), 398
–410 https://doi.org/10.1038/s41593-020-0591-0 NANEFN 1097-6256
(2020).
Google Scholar
L. Yu et al.,
“Hemodynamic monitoring in different cortical layers with a single fiber optical system,”
Proc. SPIE, 10481 104811H https://doi.org/10.1117/12.2289138 PSISDG 0277-786X
(2018).
Google Scholar
S. C. Kanick et al.,
“Monte Carlo analysis of single fiber reflectance spectroscopy: photon path length and sampling depth,”
Phys. Med. Biol., 54
(22), 6991
–7008 https://doi.org/10.1088/0031-9155/54/22/016 PHMBA7 0031-9155
(2009).
Google Scholar
A. Bisht et al.,
“Investigating presence of motion artifacts in the oxygen saturation signal during in-vivo fiber photometry,”
Proc. SPIE, 11946 1194603 https://doi.org/10.1117/12.2609877 PSISDG 0277-786X
(2022).
Google Scholar
A. Bisht et al.,
“System for simultaneous measurement of fluorescence and hemodynamics (SSMFH) from deep brain regions of freely moving rodents,”
Proc. SPIE, 12828 128280A https://doi.org/10.1117/12.3002776 PSISDG 0277-786X
(2024).
Google Scholar
K. Simone et al.,
“Open-source, cost-effective system for low-light in vivo fiber photometry,”
Neurophotonics, 5
(2), 025006 https://doi.org/10.1117/1.NPh.5.2.025006
(2018).
Google Scholar
T. M. Gonçalves et al.,
“Spectral optical properties of rabbit brain cortex between 200 and 1000 nm,”
Photochem, 1
(2), 190
–208 https://doi.org/10.3390/photochem1020011
(2021).
Google Scholar
J. Rheims, J. Köser and T. Wriedt,
“Refractive-index measurements in the near-IR using an Abbe refractometer,”
Meas. Sci. Technol., 8
(6), 601
–605 https://doi.org/10.1088/0957-0233/8/6/003 MSTCEP 0957-0233
(1997).
Google Scholar
F. Vaudelle, J. P. L’Huillier and M. L. Askoura,
“Light source distribution and scattering phase function influence light transport in diffuse multi-layered media,”
Opt. Commun., 392 268
–281 https://doi.org/10.1016/j.optcom.2017.02.001 OPCOB8 0030-4018
(2017).
Google Scholar
H. J. van Staveren et al.,
“Light scattering in lntralipid-10% in the wavelength range of 400–1100 nm,”
Appl. Opt., 30
(31), 4507 https://doi.org/10.1364/AO.30.004507 APOPAI 0003-6935
(1991).
Google Scholar
A. Kienle, F. K. Forster and R. Hibst,
“Anisotropy of light propagation in biological tissue,”
Opt. Lett., 29
(22), 2617 https://doi.org/10.1364/OL.29.002617 OPLEDP 0146-9592
(2004).
Google Scholar
W. A. Copen, M. H. Lev and O. Rapalino, Brain Perfusion: Computed Tomography and Magnetic Resonance Techniques, 135 1st ed.Elsevier B.V.(
(2016). Google Scholar
A. F. Abdelnour and T. Huppert,
“Real-time imaging of human brain function by near-infrared spectroscopy using an adaptive general linear model,”
NeuroImage, 46
(1), 133
–143 https://doi.org/10.1016/j.neuroimage.2009.01.033 NEIMEF 1053-8119
(2009).
Google Scholar
M. A. Kamran, M. Y. Jeong and M. M. Mannan,
“Optimal hemodynamic response model for functional near-infrared spectroscopy,”
Front. Behav. Neurosci., 9 1
–11 https://doi.org/10.3389/fnbeh.2015.00151
(2015).
Google Scholar
Y. Ma et al.,
“Wide-field optical mapping of neural activity and brain haemodynamics: considerations and novel approaches,”
Philos. Trans. R. Soc. B: Biol. Sci., 371
(1705), 20150360 https://doi.org/10.1098/rstb.2015.0360 PTRBAE 0962-8436
(2016).
Google Scholar
D. J. F. Cohen et al.,
“Fast estimation of adult cerebral blood content and oxygenation with hyperspectral time-resolved near-infrared spectroscopy,”
Front. Neurosci., 17 1
–12 https://doi.org/10.3389/fnins.2023.1020151 1662-453X
(2023).
Google Scholar
W. Bachir and O. Hamadah,
“Second derivative diffuse reflectance spectroscopy for estimating tissue hypoxia,”
OSA Contin., 4
(2), 650 https://doi.org/10.1364/OSAC.410807
(2021).
Google Scholar
Y. Ti and W. C. Lin,
“Probe contact pressure effects on in vivo diffuse reflectance and fluorescence spectroscopy,”
Biomed. Opt., 16
(6), 4250
–4262 https://doi.org/10.1364/oe.16.004250 BOEICL 2156-7085
(2008).
Google Scholar
C. Ash et al.,
“Effect of wavelength and beam width on penetration in light-tissue interaction using computational methods,”
Lasers Med. Sci., 32
(8), 1909
–1918 https://doi.org/10.1007/s10103-017-2317-4
(2017).
Google Scholar
|
Simulations
Hemodynamics
Optical fibers
Brain
Motion measurement
Backscatter
Animals