# Nonlinear neural patterns are revealed in high frequency functional near infrared spectroscopy analysis

Ameer Ghouse, Diego Candia-Rivera, Gaetano Valenza

► To cite this version:

Ameer Ghouse, Diego Candia-Rivera, Gaetano Valenza. Nonlinear neural patterns are revealed in high frequency functional near infrared spectroscopy analysis. Brain Research Bulletin, 2023, 203, pp.110759. 10.1016/j.brainresbull.2023.110759. hal-04855384

HAL Id: hal-04855384

https://hal.science/hal-04855384v1

Submitted on 24 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Nonlinear Neural Patterns Are Revealed In High Frequency Functional Near Infrared Spectroscopy Analysis

Article in Brain Research Bulletin · September 2023

DOI: 10.1016/j.brainresbull.2023.110759

CITATIONS

3

3 authors:

<figure><img src="image_2.png" /></figure>

Ameer Ghouse

Jaume I University

17 PUBLICATIONS

55 CITATIONS

SEE PROFILE

<figure><img src="image_3.png" /></figure>

Gaetano Valenza

University of Pisa

358 PUBLICATIONS

6,736 CITATIONS

SEE PROFILE

READS

54

<figure><img src="image_4.png" /></figure>

Diego Candia-Rivera

Paris Brain Institute - ICM

47 PUBLICATIONS

454 CITATIONS

SEE PROFILE
Research report

# Nonlinear neural patterns are revealed in high frequency functional near infrared spectroscopy analysis

Ameer Ghouse<sup>b, *</sup>, Diego Candia-Rivera<sup>c</sup>, Gaetano Valenza<sup>a, b, **</sup>

<figure><img src="image_5.png" /></figure>

Check for
updates

ARTICLE INFO

Keywords:

Neurovascular coupling

Nonlinear analysis

Hemodynamics

Mathematical modeling

Functional near infrared spectroscopy

A B S T R A C T

Functional Near Infrared Spectroscopy (fNIRS) is a useful tool for measuring hemoglobin concentration. Linear theory of the hemodynamic response function supports low frequency analysis (<0.2 Hz). However, we hypothesized that nonlinearities, arising from the complex neurovascular interactions sustaining vasomotor tone, may be revealed in higher frequency components of fNIRS signals. To test this hypothesis, we simulated nonlinear hemodynamic models to explore how blood flow autoregulation changes may alter evoked neurovascular signals in high frequencies. Next, we analyzed experimental fNIRS data to compare neural representations between fast (0.2–0.6 Hz) and slow (<0.2 Hz) waves, demonstrating that only nonlinear representations quantified by sample entropy are distinct between these frequency bands. Finally, we performed group-level distance correlation analysis to show that the cortical distribution of activity is independent only in the nonlinear analysis of fast and slow waves. Our study highlights the importance of analyzing nonlinear higher frequency effects seen in fNIRS for a comprehensive analysis of cortical neurovascular activity. Furthermore, it motivates further exploration of the nonlinear dynamics driving regional blood flow and hemoglobin concentrations

1. Introduction

Functional near-infrared spectroscopy (fNIRS) is a non-invasive method used to measure brain tissue hemodynamics as a proxy for neural activity (Strangman et al., 2002b). Slow hemodynamic response signals in fNIRS oscillations (< 0.2 Hz) are well-suited for linear (i.e., spectral) analysis, particularly when the time intervals between stimuli are long (Glover, 1999). However, such a conservative cut-off frequency is aimed at avoiding overlap with frequency components in vascular dynamics that are independent of brain function, such as cardiac pulsatility or respiratory waves (Pinti et al., 2019). High frequency signals are often filtered out because they are assumed to correspond to systemic physiological or instrumentation noise (Huppert et al., 2009). However, disentangling fNIRS artifacts from systemic brain-heart and, more in general, brain-body interaction effects may be difficult (Vikner et al., 2021; Candia-Rivera et al., 2022a,b).

Recent studies suggest that high frequency (> 0.2 Hz) signals in
fNIRS may provide additional information about neurovascular activity

(Yücel et al., 2021; Santosa et al., 2018; Ghouse et al., 2020). Indeed, the nonlinear nature of the autonomic nervous system (Goldberger et al., 2002; Marmarelis, 2004; Sunagawa et al., 1998; Barbieri et al., 2017) may impact the hemodynamic response in higher frequencies. Recent research by (Ghouse et al., 2020) assessed entropy estimates of fNIRS signals that contained an upper bound of 0.6 Hz and demonstrated complementary areas of activity when compared to neural correlates observed using linear analysis of fNIRS. This finding supports recent suggestions to use cutoff filters as high as 0.5 Hz (Yücel et al., 2021). Some studies even suggest minimal fNIRS preprocessing, avoiding filtering and instead opting for robust statistics (Santosa et al., 2018). These developments underscore the importance of investigating the potential contribution of higher frequency signals in fNIRS analysis to improve on our understanding of neural and neurovascular activity.

More specifically, to investigate the potential contribution of non-
linearities in higher frequency of fNIRS signals, this study assesses dif-
ferences between neural representations of fNIRS in the traditional slow
wave (< 0.2 Hz) and the proposed fast wave (0.2–0.6 Hz) frequencies

<footnote><sup>a</sup> Department of Information Engineering, School of Engineering, University of Pisa, Italy</footnote>

<footnote><sup>b</sup> Bioengineering and Robotics Research Center “E. Piaggio”, School of Engineering, University of Pisa, Italy</footnote>

<footnote><sup>c</sup> Sorbonne Université, Paris Brain Institute (ICM), INRIA, CNRS, INSERM, AP-HP, Hôpital Pitié-Salpêtrière, Paris, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>** Corresponding author at: Department of Information Engineering, School of Engineering, University of Pisa, Italy
E-mail addresses: a.ghouse@studenti.unipi.it (A. Ghouse), gaetano.valenza@unipi.it (G. Valenza).</footnote>
using nonlinear information-theoretic methods. We hypothesize that to fully characterize the cognitive phenomena reflected in fNIRS signals, we need to assess its full spectrum with additional nonlinear analysis. This study uses sample entropy (SampEn) (Richman and Moorman, 2000) to assess fNIRS signals irregularity (and so predictability) based on embeddings of the fNIRS signals (Sauer et al., 1991).

First, we simulated plausible modulations of autoregulatory feedback on blood flow control with stochastic dynamics in the hemodynamic model (Friston et al., 2000) to motivate fast wave nonlinearity analysis. Then, using a mental arithmetic paradigm and a control motor imagery paradigm (Berntson et al., 1996), we investigated the additional information that high frequency components may provide in experimental data when a subject is under cognitive stress. We performed a representational similarity analysis (RSA) (Kriegeskorte et al., 2008; Carlin et al., 2011) on data obtained from a public database (Shin et al., 2017), comparing fast and slow wave fNIRS modalities. Our hypothesis was that fast wave activity would add significant complementary regions, particularly in nonlinear analysis, and that the correlations between fast and slow wave fNIRS activity would be lower during mental arithmetic compared to motor imagery due to regional nonlinear interactions. We assessed multivariate correlation at each detector and hypothesized that mental arithmetic topoplots of activity would be less correlated between fast and slow wave fNIRS.

## 2. Methods and materials

### 2.1. Simulations

For initial validation of plausible effects of autoregulation dynamics on hemoglobin, we first simulate a mechanistic balloon model of hemodynamic activity (Cui et al., 2010; Friston et al., 2000). First, a flow inducing signal (s), i.e. the response of the vascular system to a neural metabolic demand, is linearly described for the sustenance of incoming blood flow ($f_{in}$). Particularly:

$$ \dot{s} = \epsilon u(t) - \frac{s}{\tau_s} - \frac{f_{in} - 1}{\tau_f} \quad (1) $$

$$ \dot{j}_{in} = s $$

u is the control signal which is the neural signal generating the flow. $\epsilon$ is the efficacy by which the neural signal can sustain a flow-inducing signal, i.e. a response that may dilate vessels to increase blood flow into. This is related to the vascular resistance (r), whose dynamics are $\dot{r} = -r^2s$ (Friston et al., 2000). $\tau_s$ is the time constant for this flow-inducing signal and $\tau_f$ is the time constant describing the effects of autoregulatory feedback from blood flow. Considering $\tau_f$ comprises information on autoregulatory effects of blood flow (the time constant for returning back to baseline), this is the parameter we later model to assess its nonlinear effects on hemoglobin concentrations.

From the sustaining blood flow generating signal, a so-called “balloon model” describes the dynamics of the volume of a blood vessel (the balloon), and the permeation of hemoglobin inside and outside the vessel (Buxton et al., 1998). Explicitly, the rate of decay of blood volume is related to blood flowing in and blood flowing out of a vessel:

$$ \begin{aligned} \tau_o \dot{v} &= f_{in} - f_{out}(v) \\ f_{out} &= v^{1/a} \end{aligned} \quad (2) $$

$a$ describes the capacity at which a balloon can expel water, having been distended by the surge of inflow. $\tau_o$ is the time constant which governs the rate of change of the volume, which is similarly intertwined with the rate of change of deoxyhemoglobin in the venous compartment (q):

$$ \begin{aligned} \tau_o \dot{q} &= f_{in} \frac{E(f_{in}, E_o)}{E_o} - f_{out}(v) \frac{q}{v} \\ E(f_{in}, E_o) &= 1 - (1 - E_o)^{1/f_{in}} \end{aligned} \quad (3) $$

The function E describes the oxygen extraction coefficient, or efficacy of the tissue in extracting oxygen from the incoming blood, while $E_o$ is the resting oxygen extraction coefficient. In (Cui et al., 2010), an extension was proposed to relate the blood volume and total hemoglobin concentration as:

$$ \tau_o \dot{p} = (f_{in} - f_{out}(v)) \frac{p}{v} \quad (4) $$

Then, oxyhemoglobin is merely $o = p - q$.

We reiterate the observation in eq. (1) that the parameter $\tau_f$ is particularly related to the autoregulatory feedback. Exactly how its value relates to the true autoregulatory changes is little understood, though literature states that responses to autoregulatory changes occur over a period of 1–2 min (or less than 0.02 Hz) (Lemkuil et al., 2013). To maintain this expected periodicity, while allowing random deviations due to uncertain dynamics, we designed a 2nd order stochastic dynamical model whose power spectral density (PSD) on average has a peak at the frequency of hypothesized autoregulatory changes.

$$ \ddot{t}_f = A(-2(\tau_f - .8) - 10\dot{t}_f)dt + 10dW \quad (5) $$

The steady state value of $\tau_f$ is at 0.8 (normalized units), while the steady state of its change is 0. A is a parameter that modulates how fast it returns to steady state, and the diffusion term is random variations outside the potential normal oscillatory behavior of $\tau_f$. When there is no autoregulatory activity, it quickly returns to steady state, with $A=10$. When there is autoregulatory activity, it is more free to change with $A = .1$. black Eq. (5) was devised to model the realizations of $\tau_f$ through simple, non-trivial stochastic differential equations. Such a generative model embeds interpretable 2nd order linear dynamics showing emergent nonlinear oscillatory properties in the observed variable.black We integrate the stochastic differential equations with the Euler-Maruyama method (Platen and Bruti-Liberati, 2010).

We numerically simulated 100 $\tau_f$ time courses to illustrate expected spectral properties of the stochastic dynamics. This was achieved by obtaining an estimate of the power spectral density (PSD) using the Welch method (Welch, 1967). The units of the PSD for $\tau_f$ is in arbitrary units, considering its the time constant of a normalized blood flow signal seen in eq. (1). We simulated a block design experiment using these generated $\tau_f$ time courses, where each block is 40 s long with a minute long rest, and 5 stimulus blocks would either induce $\tau_f$ modulations (i.e. where $A = 0.1$) or 5 stimulus blocks would not induce the modulations. The final simulated hemoglobin concentration time courses are corrupted with a signal-to-noise (SNR) ratio of 0 dB.

### 2.2. Experiment data

A publicly available dataset was used to obtain fNIRS signals with the desired experimental protocol for this study, as reported in (Shin et al., 2017). In summary, the experiment recruited twenty-nine healthy subjects (aged 28.5 ± 3.7), fifteen of which were female and fourteen male. Of these twenty-nine healthy subjects, all were right-handed expect for one. As according to Shin et al. (2017), all participants were free from neurological, psychiatric, or any brain-related disorders. They were fully informed about the experimental process, and written consent was obtained from each volunteer. After the experiment, they received financial compensation. The study adhered to the guidelines of the declaration of Helsinki and received approval from the Ethics Committee of the Institute of Psychology and Ergonomics, Technical University of Berlin (approval number: SH_01_20150330). Three trials were performed with ten repetitions of mental arithmetic and baseline events, and three trials were performed with ten repetitions of right-hand and
left-hand motor imagery for each subject. We note that the motor imagery tasks and mental arithmetic tasks were done independently, not concurrently. Thirty-six fNIRS series were acquired for each subject with 10 Hz sampling rate.

The experiment design had sixty seconds of resting state to start data acquisition from a subject, after which an instruction was shown on the screen indicating which task was to be performed–either an arithmetic problem, a “–” for a baseline, or a left or right arrow for motor imagery. The subject performed the indicated task for ten seconds, with a subsequent fifteen second resting phase before the next instruction. After twenty repetitions of these instructions and tasks (ten repetitions for each task in an experiment run), a sixty second rest was performed. A total of three trials were performed, for a total of thirty repetitions per event.

## 2.3. fNIRS signals

Thirty six channels of optical densities (OD) were resolved from source detector pairs comprising 760 nm and 850 nm wavelengths at distances of 3 cm covering the frontal, lateral parietal and posterior cortical regions as seen in Fig. 1a. The modified Beer Lambert law was used to convert the ODs to deoxyhemoglobin (Hb) and oxyhemoglobin (HbO) (Strangman et al., 2002a).

Fig. 2 illustrates the preprocessing pipeline. After applying the modified Beer-Lambert law to resolve the 36 channels seen in Fig. 2b from the source-detector pairings in Fig. 2a, band-pass frequency filters were applied to extract traditional hemodynamic bands (< 0.2 Hz) (Strangman et al., 2002a; Pinti et al., 2019) or the proposed increased hemodynamic band (0.2–0.6 Hz). A wavelet filtering approach using a Daubechies 5 wavelet, nine level decomposition was used to further reduce instrumentation noise such as movement in the oxy- and deoxyhemoglobin signals (Molavi and Dumont, 2012). Detrending and pre-whitening with an AR(1) model was then performed to remove temporally structured noise in the signal (Huppert, 2016). The signals were separated into epochs of 30 s (comprising the 15 s task phase, and the 15 s rest phase), with each channel at each activity block being referenced to the mean of the previous 5 s. We exclude the 2 second instruction phase before the task-phase from the analysis to reduce

potential confabulation of task-evoked responses from processing the instructions. Total hemoglobin was computed as the addition of both Hb and HbO. The total hemoglobin is important considering the integration of the concentrations may present unique temporal dynamics revealing distinct nonlinear temporal effects. For all three signals, entropies and mean estimates were extracted then passed into the 1st and 2nd level analysis. Before second level analysis, the results were spatially smoothed to increase the sensitivity of random effects from detector locations (Tak et al., 2016). A Gaussian smoothing kernel with a full-width-half-max of 1.5 cm (half the source-detector separation distance) was used as seen in Fig. 1b.

## 2.4. Entropy analysis

To analyze nonlinearity and regularity of the fNIRS signals, we exploited sample entropy (SampEn) (Richman and Moorman, 2000). SampEn is a method to calculate the entropy of a dynamical system in its phase space. In other words, it assesses how much information it takes to characterize the dynamics of a system. It shares similarities with methods like approximate entropy (ApEn) (Pincus, 1991), except the SampEn algorithm alleviates biases that are introduced in ApEn by not considering self-matches in its calculation of the correlation integral used in the definition of entropy in dynamical systems (Delgado-Bonal and Marshak, 2019b; Sinai, 1959). There are a whole slew of other entropy estimation methods that could have been used instead, however we have previously performed a study comparing a whole battery of entropy estimation methods for assessing fNIRS signals, with the conclusion that SampEn provides similar results as the other entropy estimation methods (Ghouse et al., 2020).

A delay-time $\tau$ and embedding dimension $m$ are needed to reconstruct manifolds using delay-coordinates. $\tau$ was selected as the first zero of the autocorrelation, while $m$ was found using the false nearest neighbors approach (Abarbanel et al., 1993).

For calculating SampEn, radius $R = 0.2 \times \sigma_x$ was used as the threshold to determine whether states were neighbors, where $\sigma_x$ is the standard deviation of the fNIRS time series (Delgado-Bonal and Marshak, 2019a). The particular equation describing how to calculate SampEn is:

<figure><img src="image_6.png" /><figcaption>Fig. 1. Illustration of spatial analysis of fNIRS signals. (a) shows the position of the optodes used to derive the 36 fNIRS channels seen in (b). (b) shows an example of the full-width-half-max of the smoothing Gaussian kernel prior to second level group analysis (Tak et al., 2016) in the color green.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Fig. 2. Overview of the analysis pipeline used for each fNIRS signal in the dataset, which includes computation of hemoglobin concentrations, band pass filtering, motion artifact correction, removal of linear trends and serial correlations, total hemoglobin calculation, state space reconstructions using delay-coordinate embeddings, estimation of sample entropy or mean value for each trial, first level analysis using median value over trials, and spatial smoothing to improve sensitivity in group analysis results (Tak et al., 2016).</figcaption></figure>

$$ \text{SampEn} = - \log \left( \frac{\sum_{i=1, i \neq j}^{N-m} \frac{1}{N} \text{ number of } |X_i^{m+1} - X_j^{m+1}| \langle R \rangle}{\sum_{i=1, i \neq j}^{N-m} \frac{1}{N} \text{ number of } |X_i^m - X_j^m| \langle R \rangle} \right) \quad (6) $$

X in this equation is a state space reconstructed using a time series. The superscript denotes the embedding dimension while the subscript denotes the state index, for which there exists N states.

## 2.5. Representational similarity analysis

Representational similarity analysis is a method used to compare multivariate data such as different brain data types (Kriegeskorte et al., 2008). This paper calculates similarities between tasks using the distance correlation in order to construct the representational similarity matrices for fast or slow wave (Székely et al., 2007; Geerligs et al., 2016).

In brief, given random vector X and random vector Y with dimensions $\mathbb{R}^p$ and $\mathbb{R}^q$ respectively, distance correlation evaluates independence by integrating the distance between the random vectors in conjunction with a weighting function $w(t, s) = (c_p c_q |t|^1 |s|^q)^{-1}$, where $c_d = \frac{\pi^{\frac{1}{2}d}}{\Gamma(\frac{1}{2}d)}$ corresponds to half the surface area of a unit sphere in the given dimensionality d. This leads to the following statistic:

$$ \tilde{Y}^{-2}(X, Y) = \int_{\mathbb{R}^{p+q}} |\phi_{X,Y}(t, s) - \phi_x(t)\phi_y(s)|^2 w(t, s) dt ds \quad (7) $$

$$ \mathcal{R}^2(X, Y) = \frac{\tilde{Y}^{-2}(X, Y)}{\tilde{Y}^{-2}(X, X) \tilde{Y}^{-2}(Y, Y)} $$

As the dimensions of the random vectors, p and q, approach infinity, the statistic converges to a Student's t-distribution and can be approximated as such for hypothesis testing (Székely and Rizzo, 2013). Given $a_{ij} = |x_i - x_j|$ and $b_{ij} = |y_i - y_j|$, where i and j represent the ith and jth observations of x or y, the sample covariance can be estimated as:

$$ \begin{aligned} A_{ij} &= a_{ij} - \frac{1}{N} \sum_{n=1}^{N} a_{in} - \frac{1}{N} \sum_{k=1}^{N} a_{kj} + \frac{1}{N^2} \sum_{n=1}^{N} \sum_{k=1}^{N} a_{kn} \\ B_{ij} &= b_{ij} - \frac{1}{N} \sum_{n=1}^{N} b_{in} - \frac{1}{N} \sum_{k=1}^{N} b_{kj} + \frac{1}{N^2} \sum_{n=1}^{N} \sum_{k=1}^{N} b_{kn} \end{aligned} \quad (8) $$

$$ \tilde{Y}^{-2}(x, y) = \frac{1}{N^2} \sum_{i=1}^{N} \sum_{j=1}^{N} A_{ij} B_{ij} $$

In our application, random vectors X and Y have equivalent dimensionality $\mathbb{R}^3$ as we are looking at multivariate correlations between oxy, deoxy and total hemoglobin. Furthermore, due to there being 4 experimental conditions (baseline, mental arithmetic, left- and right-hand imagery), the resulting representational dissimilarity matrix (RDM) is a 4 × 4 matrix. Representational dissimilarity outputs were obtained for each detector location for each measure (entropies or mean value). Extracting the upper triangle, we have two 29 matrices for slow and fast wave fNIRS respectively to perform statistical analysis upon.

## 2.6. Fast vs. slow wave spatial analysis

For each subject, the data obtained from each detector had a shape of N repetitions × N concentrations. By calculating the median over the N repetitions for either fast or slow wave fNIRS, we obtain a vector of size N Concentrations × N detectors for each subject. This analysis enables us to determine whether the spatial fNIRS arrays for fast and slow waves are dependent, or whether the null hypothesis that they are independent can be rejected. We expect that the slow and fast wave analysis will become more independent when there is a task inducing changes in autoregulatory activity. We also perform a Wilcoxon paired analysis at the detector level to compare which medians are significantly different in the group of subjects for each concentration.

## 2.7. Statistical analysis

In order to perform statistical analysis of the data, we exploited the Scipy statistics package in Python (Virtanen et al., 2020).
To compare RDMs between fast and slow wave fNIRS, a group analysis was conducted using a Friedman test, with a Bonferroni correction to account for multiple comparisons. A family-wise p-value of less than 0.05 was considered significant. Nonparametric Friedman statistical tests were used to avoid assumptions of normality in the dataset when performing intergroup analysis (Friedman, 1937). Then, the Bonferroni correction was applied to be as conservative as possible when correcting for family-wise error rates such that any inferences are robust to Type I errors (Tukey, 1953; Dunn, 1961).

For median analysis, rainclouds were generated from standardized Z transformations to compare the shapes of distributions between fast and slow wave analysis. A Kolmogorov-Smirnov test was applied to compare concentrations analyzed from fast and slow wave analysis, correcting for multiple comparisons. A paired sample Wilcoxon signed rank test was performed on the difference of median value for each concentration between mental arithmetic vs baseline or right-hand vs left-hand motor imagery. Detectors were retained for further analysis if any of the 3 concentrations returned significant with a false alarm rate of $α = 0.05/6 = 0.008$. For remaining detectors, distance correlation was performed between fast and slow wave fNIRS.

# 3. Results

The first set of results are from simulations of neurovascular signals generated by a hemodynamic model, evoked from events that either induce or do not induce autoregulatory changes. The second set of results are from group statistics of representational similarity analysis between fast wave and slow wave cortical distribution of activity, and median analysis of cortical activity. Standard fNIRS analysis refers to time averaged value of fNIRS activity during a task, while nonlinear fNIRS analysis refers to SampEn during the time duration of the task.

## 3.1. Simulations

Fig. 3 illustrates the results of realizations of the hypothesized autoregulatory feedback changes as according to eq. (5). The median peak frequency was at 0.02 Hz, corresponding to periodicity of 50 s, with periodicity ranges from 12 ss to 100 s. Fig. 4 illustrates the effects that these autoregulatory changes may have on the hemodynamics, where HbO is the oxyhemoglobin and HbR is the deoxyhemoglobin. Over the 100 realizations of the block design experiment, the mean entropy of the hemodynamic response in the high frequency (0.2–0.6 Hz) regime without autoregulation modulations was at 0.212

($\pm 0.015$) bits compared to 0.163 ($\pm 0.011$) bits when the modulations occur; through a paired Wilcoxon signed rank test, the SampEns were found to be significantly different ($p \ll 0.001$). On the other hand, the power of this high frequency band without autoregulation modulation was found to be 0.0368 ($\pm 0.00441$) as compared to 0.0372 ($\pm 0.00385$) during autoregulation modulation; the power bands were not significantly different according to the Wilcoxon signed rank test ($p = 0.46$).

## 3.2. Experimental fNIRS data

### 3.2.1. Representational similarity analysis

Results from the Friedman analysis comparing the upper triangle of the RDMs constructed using either fast wave or slow wave fNIRS signals according to the methods described in section 2.5 can be seen in Fig. 5, for mean and SampEn respectively. For the mean value, no detector has significantly different representations between slow and fast wave fNIRS, whereas each quadrant of the sensor space on the cortex contained significantly different detectors for SampEn. Generally, SampEn contained higher Friedman test scores than the mean analysis.

### 3.2.2. Fast vs slow wave spatial analysis

Using the methods described in 2.6, we obtained group level median analysis results for fast and slow wave fNIRS. To reiterate, this entailed calculating the median value of the mean or SampEn signal over the repetition of the task performed by the participant, and then performing a contrast analysis between the effects of mental arithmetic (MA) and baselines (BL), or left-hand motor imagery (LH0 and right-hand motor imagery (RH). A distillation of the spatial results can be seen in the raincloud Fig. 6 to visualize distributions of contrasts between the task comparisons across subjects; Kolmogorov-Smirnov tests were performed to ascertain whether the distributions were significantly different. Supplementary fig. S1 shows the spatial cortical map of contrasts between mental arithmetic (MA) and baseline (BL), or right- (RH) or left-hand (LH) motor imagery tasks.

From the raincloud plots in Fig. 6, distributions for mean results consistently appear to be overlapped when comparing slow and fast wave different in medians, whereas SampEn demonstrates a flatter distribution for slow wave analysis as compared to fast wave analysis. Motor imagery also presents a flatter distribution in SampEn when using slow wave analysis compared to fast wave analysis, however this greater variance is similarly observed in mean value analysis. Performing a Kolmogorov-Smirnov two sample test, we found that the distribution of fast and slow wave analysis was significantly different between SampEn

<figure><img src="image_8.png" /><figcaption>Fig. 3. Spectral analysis of realizations of $\tau_f$ simulated from the stochastic differential equation in eq. (5). a) Represents the average power spectral density (PSD) of $\tau_f$ as estimated by the Welch method over 100 realizations. The units of the PSD are arbitrary as $\tau_f$ is the time constant of $f_{in}$ seen in eq. (1). b) on the other hand is a histogram of the peak frequency of the PSD of $\tau_f$ over the 100 realizations.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Fig. 4.** (a) demonstrates a single realization of a band passed oxyhemoglobin power spectrum density when either autoregulatory activity is occurring with a task, or when there is no autoregulatory activity with a task. (b) displays a representative simulation of the signals in the model and the events that influence them (whether a task induces autoregulatory behavior or not).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Fig. 5.** Results of a group level Friedman test comparison between the upper triangle of the RDM for fast (0.2–0.6 Hz) and slow wave (< 0.2 Hz) fNIRS on the group level, where significance indicates at least one element in the upper triangle had a significantly different median between fast and slow wave fNIRS. Channels that are significant are marked by a thick outline on the marker.</figcaption></figure>

in oxy and deoxyhemoglobin, while in no mean analysis were they significantly different.

For mental arithmetic and baseline Wilcoxon paired analysis, mean estimates show significant detectors only in slow wave fNIRS, with a cluster in the right lateral cortical areas and detectors bilaterally in the frontal cortex. SampEn mainly determined significant detectors in the frontal right cortex and left parietal cortex in the slow wave analysis, while fast wave analysis complemented information in the right parietal cortical areas.

Left- and right-hand motor imagery Wilcoxon paired analysis demonstrated activity bilaterally in the parietal cortical areas with both SampEn and mean estimates. Both mean and SampEn demonstrated slow activity predominantly in the right parietal cortex, while fast

activity was found in the left lateral parietal cortex. Mean also contained significant detectors in the medial left parietal cortex in slow wave analysis and in the occipital regions for fast wave analysis.

As a sanity check, we also compared the neural representations between slow frequency bands (0 Hz, 0.2 Hz) and cardiac frequency bands (0.8 Hz, 3 Hz). The distance correlation of topolots of fast wave vs slow wave was not significant, indicating that the fast wave fNIRS analysis in freqs (0.2 Hz, 0.6 Hz) was not a result of systemic physiological founders (see also supplementary fig. S2).

# 4. Discussion and conclusion

In this study, we compared neural activity representations between
<figure><img src="image_11.png" /><figcaption>**Fig. 6.** Rainclouds illustrating the standardized Z transformed group level distributions of the absolute value of the difference between medians for mental arithmetic (a) and motor imagery (b) for either fast wave (0.2–0.6 Hz) or slow wave (0–0.2 Hz) fNIRS signals. The “*” represents significant difference between the fast and slow wave distribution for the given measure.</figcaption></figure>

fast wave (0.2–0.6 Hz) fNIRS and the standard slow wave (< 0.2 Hz) fNIRS signal using both linear measures of mean value and nonlinear measures derived from sample entropy. The fNIRS signals were obtained from a publicly available dataset of 29 subjects performing mental arithmetic and left-/right-hand motor imagery tasks (Shin et al., 2017). Our hypothesis was that mental arithmetic tasks, which induce cognitive stress, would provide more variation in areas of cortical activation than motor imagery when comparing slow and fast wave fNIRS using nonlinear analysis. This hypothesis was motivated by literature suggesting the nonlinear modulation of fNIRS dynamics sustained by various factors, including those induced by autonomic nervous system activity (Friston et al., 2000; Gianaros et al., 2012; Sheng and Zhu, 2018; Candia-Rivera et al., 2023). Hemodynamic models are also known to be nonlinear in nature (Friston et al., 2000; Buxton et al., 1998; Friston, 2001).

To investigate the potential effects of vasomotor property dynamics on autoregulation modulations, we first conducted a simulation study using equations (1) and (5). Our theoretical model for modulations of vasomotor tone closely matched the expected spectral profile of real-world vasomotor tone variations at frequencies below 0.2 Hz, as shown in Fig. 3. Additionally, we found that autoregulation variations in flow signals that propagate hemoglobin concentration dynamics (as depicted in Fig. 4) reduce entropy in the observed signal even with an

SNR of 0 dB. This suggests that with proper statistical power, these autoregulatory dynamics could potentially be observed in real-world data using nonlinear analysis methods. It is important to note that our empirical data analysis should not be interpreted as evidence that high frequency signal entropy is a measure of autoregulation. Rather, we hypothesize that the nonlinear activity observed in high frequency fNIRS provides complementary information that is crucial for a comprehensive characterization of cognitive states, such as those associated with stress processing during mental arithmetic tasks.

To explore potential distinctions in neural signal representations between fast and slow waves in fNIRS data, we conducted a Friedman test to compare results on the tasks of baseline, mental arithmetic, and motor imagery, as shown in Fig. 5. The mean value analysis did not reveal any significant differences in task combinations, as expected since the mean corresponds to the DC value of the signal and high pass filtering only affects attenuation of the DC value. In contrast, nonlinear SampEn analysis demonstrated significant differences between fast and slow wave fNIRS detectors in all four quadrants of the sensor space covering cortex, with a higher Friedman score across the cortex than the mean analysis. While it is unclear how filtering affects state space regularity, high pass filtering may result in lower mean entropy and higher variance than low-pass filtering (Borges et al., 2020). This effect may correspond to the significantly different SampEn representational
similarity seen across subjects.

Upon examining motor imagery results we observed changes in lat-
erality in SampEn analysis depending on the frequency band of focus,
with fast wave analysis being more sensitive to left hemisphere activity
and slow wave analysis being more sensitive to right hemisphere activity
(see Supplementary fig. S1). A similar phenomenon was seen in the
mean analysis, although an additional medial left hemispheric detector
was observed for slow waves. We expect this activity to correspond to
the somatosensory cortex, located in the central cortical regions, where
the hands of the cortical homunculus are located between the lateral and
medial portions of the cortex (Grodd et al., 2001). By considering both
frequency bands using SampEn analysis, we may have discovered its
bilateral effects. This bilateral difference between left- and right-hand
motor imagery was only found using slow wave analysis with just the
mean analysis. It is important to mention that the handedness of a
participant plays a critical role in the effect seen in motor imagery.
Right-handed individuals appear to have a stronger lateralization effect
specific to the side of the hand imagined; left-handed individuals on the
other hand have bilateral responses (Crotti et al., 2022). Considering the
set of participants used in this study were primarily right-handed, the
results presented in this study may be specific to this group. Nonetheless,
considering we do paired analysis of the same participants performing
left- and right-hand motor imagery, the resulting differences of cortical
signal effects observed when comparing fast and slow wave analyses
appear to indeed be a neurovascular result. This suggests that to fully
describe the complementary areas of activity that entropy provides, we
must consider the full frequency profile of potential hemodynamic
activity.

We also conducted a distance correlation analysis of the difference of
median maps and found that all mean maps were significantly corre-
lated, as expected. The cortical representations reflected in median
difference maps between tasks should be similar if we scale the signal by
merely a DC value, as previously discussed. However, our hypothesis
was that mental arithmetic should not be significantly correlated in
entropy analysis while motor imagery should be significant. Instead, we
found that neither was significant, although mental arithmetic vs.
baseline had a lower correlation than the motor imagery distinction.
This may imply that the nonlinear irregular effects induced by vaso-
motor mechanical dynamics are always occurring, although the degree
to which they affect the resulting neural signal varies depending on the
task. It's known that a significant portion of cerebral vascular resistance
comes from vasomotor control in arterioles (Mandeville et al., 1999).
The nonlinear transform of the constant dynamics of vasomotor control
to hemoglobin concentrations may always be present, yet their scale is
modulated by a stress task such as mental arithmetic, making their high
frequency contributions to hemoglobin concentrations more irregular in
the case of mental arithmetic. We finally remark that the proposed fast
wave analysis within (0.2 Hz, 0.6 Hz) is not affected by non-cognitive
phenomena associated with oscillations in the (0.8 Hz, 3 Hz) band
(Hoshi, 2016).

We recognize that this study may have its limitations. First, we did
not perform simultaneous fNIRS and autonomic signals analysis.
Respiration and heart rate, typically occurring within the fast-wave
band, could be confounding factors in our data analysis since a multi-
tude of regulatory processes may change during cognitive tasks. How-
ever, the heterogeneity in the response across the cortex suggests that
such oscillations may not systematically affect fNIRS series all over the
scalp. Literature suggests possible differential respiratory effects on the
scalp, indicating the need for future studies to take actual respiratory
and heart rate signals into account (Candia-Rivera et al., 2022c; Dubois
et al., 2016). Additionally, signals corresponding to blood pressure in
Mayer waves contaminate fNIRS signals within the low frequency bands
(Pinti et al., 2015; Kirilina et al., 2013), and nonlinearity evoked in
blood pressure variations could similarly reflect in high frequency
fNIRS. Hence, future studies should also monitor end-tidal CO2 as par-
tial pressure of CO2 has been shown to affect vascular dynamics (Mas

et al., 2000). Despite these limitations may introduce some uncertainty
regarding the origin of nonlinear fNIRS dynamics in the fast wave band,
our study's outcomes remain robust and reliable. This is further sub-
stantiated by the corroborative evidence provided through our
comprehensive simulation-based analysis.

In conclusion, our study shows that nonlinear analysis can detect distinct neural representations in fast (0.2–0.6 Hz) waves compared to conventional slow (<0.2 Hz) wave frequencies in fNIRS data. Traditional methods like mean estimates could not reveal such distinct representations. While the origin of neural hemodynamic activity from a behavioral stimulus may be associated with oscillations at frequencies <0.2 Hz, nonlinear neurovascular interactions may generate fNIRS oscillations at higher frequencies. Thus, when assessing effects on fNIRS, comprehensive characterizations should also consider the nonlinear properties of high-frequency band oscillations.

CRediT authorship contribution statement

**Ameer Ghouse:** Conceptualization, Methodology, Software, Writing- Original draft preparation, Formal Analysis, Visualization, Investigation, Validation. **Diego Candia Rivera:** Writing- Review and Editing, Validation. **Gaetano Valenza:** Conceptualization, Supervision, Writing- Review and Editing, Validation, Funding Acquisition.

Declaration of Competing Interest

The authors of the submitted manuscript titled "Nonlinear Neural Patterns Are Revealed In High Frequency Functional Near Infrared Spectroscopy Analysis" hereby declare they have no conflict of interests to report.

Data availability

All analysis code is available at the GitHub repository https://github.
com/lemiceterieux/fastSlowfNIRSAnalysis. and data are available at
https://doc.ml.tu-berlin.de/hBCI/.

Acknowledgments

This study was financed in part by the European Commission Hori-
zon 2020 project RHUMBO, grant agreement ID 813234, and by the
Italian Ministry of Education and Research (MIUR) in the framework of
the Forelab project (Departments of Excellence).

**Appendix A. Supporting information**

Supplementary data associated with this article can be found in the
online version at doi:10.1016/j.brainresbull.2023.110759.

References

Abarbanel, H.D.I., Brown, R., Sidorowich, J.J., Tsimring, L.S., 1993. The analysis of observed chaotic data in physical systems. Rev. Mod. Phys. 65, 1331-1392.

Barbieri, R., Scilingo, E.P., Valenza, G., 2017. Complexity and nonlinearity in cardiovascular signals. Springer International Publishing AG.

Berntson, G.G., Cacioppo, J.T., Fieldstone, A., 1996. Illusions, arithmetic, and the bidirectional modulation of vagal control of the heart. Biol. Psychol. 44 (1), 1-17.

Borges, V.S., Nepomuceno, E.G., Duque, C.A., Butusov, D.N., 2020. Some remarks about entropy of digital filtered signals. Entropy (Basel, Switzerland) 22, 365.

Buxton, R.B., Wong, E.C., Frank, L.R., 1998. Dynamics of blood flow and oxygenation changes during brain activation: the balloon model. Magn. Reson. Med. 39 (6), 855-864.

Candia-Rivera, D., Catrambone, V., Thayer, J.F., Gentili, C., Valenza, G., 2022a. Cardiac sympathetic-vagal activity initiates a functional brain-body response to emotional arousal. Proc. Natl. Acad. Sci. 119 (21), e2119599119.

Candia-Rivera, D., Catrambone, V., Barbieri, R., Valenza, G., 2022b. Functional assessment of bidirectional cortical and peripheral neural control on heartbeat dynamics: a brain-heart study on thermal stress. NeuroImage 251, 119023.
Candia-Rivera, D., Sappia, M.S., Horschig, J.M., Colier, W.N.J.M., Valenza, G., 2022c. Confounding effects of heart rate, breathing rate, and frontal fnirs on interoception. Sci. Rep. 12 (1), 20701.

Candia-Rivera, D., Norouzi, K., Ramsey, T.Z., Valenza, G., 2023. Dynamic fluctuations in ascending heart-to-brain communication under mental stress. Am. J. Physiol. -Regul., Integr. Comp. Physiol. 324 (4), R513-R525.

Carlin, J.D., Calder, A.J., Kriegeskorte, N., Nili, H., Rowe, J.B., 2011. A head view-invariant representation of gaze direction in anterior superior temporal sulcus. Curr. Biol.: CB 21, 1817-1821.

Crotti, M., Koschutnig, K., Wriesnegger, S.C., 2022. Handedness impacts the neural correlates of kinesthetic motor imagery and execution: a fmri study. J. Neurosci. Res. 100 (3), 798-826.

Cui, X., Bray, S., Reiss, A.L., 2010. Functional near infrared spectroscopy (nirs) signal improvement based on negative correlation between oxygenated and deoxygenated hemoglobin dynamics. NeuroImage 49 (4), 3039-3046.

Delgado-Bonal, A., Marshak, A., 2019a. Approximate entropy and sample entropy: a comprehensive tutorial. Entropy 21 (6).

Delgado-Bonal, A., Marshak, A., 2019b. Approximate entropy and sample entropy: a comprehensive tutorial. Entropy 21 (6), 541.

Dubois, M., Chenivesse, C., Raux, M., Morales-Robles, A., Nierat, M.-C., Garcia, G., Navarro-Sune, X., Chavez, M., Martinerie, J., Similowski, T., 2016. Neurophysiological evidence for a cortical contribution to the wakefulness-related drive to breathe explaining hypocapnia-resistant ventilation in humans. J. Neurosci. 36 (41), 10673-10682.

Dunn, O.J., 1961. Multiple comparisons among means. J. Am. Stat. Assoc. 56 (293), 52-64.

Friedman, M., 1937. The use of ranks to avoid the assumption of normality implicit in the analysis of variance. J. Am. Stat. Assoc. 32 (200), 675-701.

Friston, K., Mechelli, A., Turner, R., Price, C., 2000. Nonlinear responses in fmri: the balloon model, volterra kernels, and other hemodynamics. NeuroImage 12 (4), 466-477.

Friston, K.J., 2001. Book review: brain function, nonlinear coupling, and neuronal transients. Neuroscientist 7 (5), 406-418. PMID: 11597100.

Geerligs, L., Cam-CAN, Henson, R.N., 2016. Functional connectivity and structural covariance between regions of interest can be measured more accurately using multivariate distance correlation. NeuroImage 135, 16-31.

Ghouse, A., Nardelli, M., Valenza, G., 2020. fnirs complexity analysis for the assessment of motor imagery and mental arithmetic tasks. Entropy 22 (7).

Gianaros, P.J., Onyewuenyi, I.C., Sheu, L.K., Christie, I.C., Critchley, H.D., 2012. Brain systems for baroreflex suppression during stress in humans. Hum. Brain Mapp. 33, 1700-1716.

Glover, G.H., 1999. Deconvolution of impulse response in event-related bold fmri. NeuroImage 9 (4), 416-429.

Goldberger, A.L., Peng, C.-K., Lipsitz, L.A., 2002. What is physiologic complexity and how does it change with aging and disease? Neurobiol. Aging 23 (1), 23-26.

Grodd, W., Hülsmann, E., Lotze, M., Wildgruber, D., Erb, M., 2001. Sensorimotor mapping of the human cerebellum: fmri evidence of somatotopic organization. Hum. Brain Mapp. 13, 55-73.

Hoshi, Y., 2016. Hemodynamic signals in fnirs. Prog. Brain Res. 225, 153-179.

Huppert, T.J., 2016. Commentary on the statistical properties of noise and its implication on general linear models in functional near-infrared spectroscopy. Neurophotonics 3 (1), 1-10.

Huppert, T.J., Diamond, S.G., Franceschini, M.A., Boas, D.A., 2009. Homer: a review of time-series analysis methods for near-infrared spectroscopy of the brain. Appl. Opt. 48, D280-D298, 19340120[pmid].

Kirilina, E., Yu, N., Jelzow, A., Wabnitz, H., Jacobs, A.M., Tachtsidis, I., 2013. Identifying and quantifying main components of physiological noise in functional near infrared spectroscopy on the prefrontal cortex. Front. Hum. Neurosci. 7, 864.

Kriegeskorte, N., Mur, M., Bandettini, P., 2008. Representational similarity analysis - connecting the branches of systems neuroscience. Front. Syst. Neurosci. 2, 4.

Lemkuil, B.P., Drummond, J.C., Patel, P.M., 2013. Chapter 8 - central nervous system physiology: Cerebrovascular. In: Hemmings, H.C., Egan, T.D. (Eds.), Pharmacology and Physiology for Anesthesia. W.B. Saunders, Philadelphia, pp. 123-136.

Mandeville, J.B., Marota, J.J., Ayata, C., Zaharchuk, G., Moskowitz, M.A., Rosen, B.R., Weisskoff, R.M., 1999. Evidence of a cerebrovascular postarteriole windkessel with delayed compliance. J. Cereb. Blood Flow. Metab. 19 (6), 679-689.

Marmarelis, V.Z., 2004. Nonlinear dynamic modeling of physiological systems, vol. 10. John Wiley & Sons.

Mas, A., Sara, P., Joseph, D., Blanch, L., Baigorri, F., Artigas, A., Fernández, R., 2000. Effect of acute moderate changes in paco2 on global hemodynamics and gastric perfusion. Crit. Care Med. 28 (2), 360-365.

Molavi, B., Dumont, G.A., 2012. Wavelet-based motion artifact removal for functional near-infrared spectroscopy. Physiol. Meas. 33, 259-270.

Pincus, S.M., 1991. Approximate entropy as a measure of system complexity. Proc. Natl. Acad. Sci. 88 (6), 2297-2301.

Pinti, P., Cardone, D., Merla, A., 2015. Simultaneous fnirs and thermal infrared imaging during cognitive task reveal autonomic correlates of prefrontal cortex activity. Sci. Rep. 5 (1), 17471.

Pinti, P., Scholkmann, F., Hamilton, A., Burgess, P., Tachtsidis, I., 2019. Current status and issues regarding pre-processing of fnirs neuroimaging data: an investigation of diverse signal filtering methods within a general linear model framework. Front. Hum. Neurosci. 12, 505, 30687038[pmid].

Platen, E., Bruti-Liberati, N., 2010. Numerical solution of stochastic differential equations with jumps in finance, vol. 64. Springer Science & Business Media.

Richman, J.S., Moorman, J.R., 2000. Physiological time-series analysis using approximate entropy and sample entropy. Am. J. Physiol. -Heart Circ. Physiol. 278 (6), H2039-H2049. PMID: 10843903.

Santosa, H., Zhai, X., Fishburn, F., Huppert, T., 2018. The fnirs brain analyzir toolbox. Algorithms 11 (5), 73.

Sauer, T., Yorke, J.A., Casdagli, M., 1991. Embedology. J. Stat. Phys. 65, 579-616.

Sheng, Y., Zhu, L., 2018. The crosstalk between autonomic nervous system and blood vessels. Int. J. Physiol., Pathophysiol. Pharmacol. 10, 17-28.

Shin, J., vonLühmann, A., Blankertz, B., Kim, D., Jeong, J., Hwang, H., Müller, K., 2017. Open access dataset for eeg.nirs single-trial classification. IEEE Trans. Neural Syst. Rehabil. Eng. 25, 1735-1745.

Y.G. Sinai, On the notion of entropy of a dynamical system, 1959.

Strangman, G., Boas, D.A., Sutton, J.P., 2002a. Non-invasive neuroimaging using near-infrared light. Biol. Psychiatry 52 (7), 679-693.

Strangman, G., Culver, J.P., Thompson, J.H., Boas, D.A., 2002b. A quantitative comparison of simultaneous bold fmri and fnirs recordings during functional brain activation. NeuroImage 17 (2), 719-731.

Sunagawa, K., Kawada, T., Nakahara, T., 1998. Dynamic nonlinear vago-sympathetic interaction in regulating heart rate. Heart Vessels 13 (4), 157-174.

Székely, G.J., Rizzo, M.L., 2013. The distance correlation t-test of independence in high dimension. J. Multivar. Anal. 117, 193-213.

Székely, G.J., Rizzo, M.L., Bakirov, N.K., 2007. Measuring and testing dependence by correlation of distances. Ann. Stat. 35 (6), 2769-2794.

Tak, S., Uga, M., Flandin, G., Dan, I., Penny, W., 2016. Sensor space group analysis for fnirs data. J. Neurosci. Methods 264, 103-112.

Tukey, J.W., 1953. The problem of multiple comparisons. Mult. Comp.

Vikner, T., Eklund, A., Karalija, N., Malm, J., Riklund, K., Lindenberger, U., Bäckman, L., Nyberg, L., Wählin, A., 2021. Cerebral arterial pulsatility is linked to hippocampal microvascular function and episodic memory in healthy older adults. J. Cereb. Blood Flow. Metab. 41 (7), 1778-1790.

Virtanen, P., Gommers, R., Oliphant, T.E., Haberland, M., Reddy, T., Cournapeau, D., Burovski, E., Peterson, P., Weckesser, W., Bright, J., van der Walt, S.J., Brett, M., Wilson, J., Millman, K.J., Mayorov, N., Nelson, A.R.J., Jones, E., Kern, R., Larson, E., Carey, C.J., Polat, I., Feng, Y., Moore, E.W., VanderPlas, J., Laxalde, D., Perktold, J., Cimrman, R., Henriksen, I., Quintero, E.A., Harris, C.R., Archibald, A.M., Ribeiro, A. H., Pedregosa, F., van Mulbregt, P., SciPy 1.0 Contributors, 2020. SciPy 1.0: fundamental algorithms for scientific computing in python. Nat. Methods 17, 261-272.

Welch, P., 1967. The use of fast fourier transform for the estimation of power spectra: a method based on time averaging over short, modified periodograms. IEEE Trans. Audio Electro 15 (2), 70-73.

Yücel, M.A., Lühmann, A.v., Scholkmann, F., Gervain, J., Dan, I., Ayaz, H., Boas, D., Cooper, R.J., Culver, J., Elwell, C.E., et al., 2021. Best practices for fnirs publications. Neurophotonics 8 (1), 012101.