Preserving cerebral autoregulation (CA), namely the mechanism aiming at maintaining a nearly stable cerebral blood flow (CBF) despite variations in mean arterial pressure (MAP), is fundamental in critical situations such as during the different phases of cardiac surgery (Hori et al 2015, Caldas et al 2018). The maintenance of CA is particularly important after the induction of general anesthesia, known to depress autonomic function and affect cerebrovascular control (Ebert et al 1992, Deutschman et al 1994, Strebel et al 1995, Keyl et al 2000, Dagal and Lam 2009, Porta et al 2013, Slupe and Kirsch 2018, Bari et al 2021, 2022, Manquat et al 2022, Saputo et al 2025), as well as in other critical phases of the surgical intervention, such as during cardiopulmonary bypass (Ono et al 2012, Easley et al 2013, Hori et al 2017) and at the end of the intervention, at the weaning (Joshi et al 2010). It has been suggested that an impaired CA could be linked to a higher risk of developing postoperative adverse events, such as stroke, either overt or silent (Joshi et al 2010, Ono et al 2012), acute kidney injury (Ono et al 2013a, Vaini et al 2019), or, more in general, major complications (Ono et al 2014, Lewis et al 2021, Liu et al 2022).

Given the impossibility of noninvasively measuring the CBF without knowing the diameter of the cerebral arteries, traditionally the CA is assessed by measuring cerebral blood velocity (CBv) at the level of middle cerebral artery via a transcranial Doppler (TCD) device (Aaslid et al 1982). However, TCD has some limitations related to the unfeasibility of acquiring signals in a relevant percentage of subjects, increasing with age and pathology, and to the need to have a trained operator to perform the acquisition. Near infrared spectroscopy (NIRS) allows the monitoring of regional saturation of oxygen (rSO2) from the prefrontal cortex. In situations of steady arterial saturation, metabolism and diffusivity, a measure of regional brain tissue oxygenation provides an indication of variations of mean CBv (MCBv). Remarkably, compared to measurement of MCBv using TCD probe, the acquisition of rSO2 based on NIRS is easier, more robust, less operator-dependent and does not require a time-consuming insonating phase because it requires only to stitch the probe on the patient’s forehead (Smielewski et al 1995, Sainbhi et al 2023). Consequently, NIRS has been increasingly utilized to monitor the CA in critical situations as during surgery or in intensive care unit (Steiner et al 2009, Brady et al 2010, Ono et al 2013b, Hori et al 2015, Rivera-Lara et al 2017).

Time-domain methods based on correlation analysis have been developed to assess the CA from CBv and arterial pressure (AP) signals via the mean flow index (Mx), or from rSO2 and AP, via the cerebral oximetry index (COx) (Ono et al 2013b, Liu et al 2021, Mol et al 2021, Kostoglou et al 2024). The rationale underlying methods based on the assessment of the degree of CBv-AP or CBv-rSO2 association is that, if the CA is functioning, its main aim is to limit the positive dependence of the response variable, namely the CBv or rSO2, on variations of the drive variable, namely the AP. Therefore, any significant increase in the degree of positive association well above the level setting the null correlation should be considered as a hallmark of the CA impairment. This observation was originally proposed based on squared coherence analysis (Giller 1990, Zhang et al 1997, 1998) computed in specific frequency bands typical of the functioning of the CA, namely the very low frequency (VLF, from 0.02 –0.07 Hz) and low frequency (LF, from 0.07 –0.2 Hz) (Claassen et al 2016). Despite the conjecture that an increase of the MCBv-MAP, or MCBv-rSO2, association above the significance threshold could indicate a CA impairment is mainly valid in the VLF and LF bands, this concept has been extended to time-domain indexes of correlation under the hypothesis that VLF and LF rhythms are dominant in CBv, rSO2 and AP signals after applying suitable filtering procedures (Czosnyka et al 1996, Kostoglou et al 2024).

Time-domain correlation-based methods were found useful to assess CA and its dysfunction in different clinical situations, such as during cardiac surgery, in intensive care unit or in traumatic brain injury patients (Lang et al 2003, Radolovich et al 2011, Ono et al 2012, Schramm et al 2012, Highton et al 2015, Severdija et al 2015, Rivera-Lara et al 2017, Bögli et al 2025b). However, the procedure set to check the significance of the positive MCBv-MAP and rSO2-MAP association, and consequently, to detect the situation of CA impairment, led to different reference values that have not been fully justified based on the methodological considerations. In the study that originally proposed the time-domain correlation-based method (Czosnyka et al 1996) and in a subsequent study that provided a correlation of Mx with the autoregulatory index (Czosnyka et al 2008), the Mx threshold exploited to separate favorable from unfavorable outcome could be reasonably set to 0 with positive values indicating CA impairment and negative values suggesting that the CA was preserved. This threshold has been subsequently corrected to 0.3 (Sorrentino et al 2011, Manquat et al 2022, Thudium et al 2023). However, in literature there is a disagreement among the value of the threshold that needs to be utilized given that 0.2 (Budohoski et al 2012), 0.25 (Kostoglou et al 2024), 0.4 (Ono et al 2013b, Kho et al 2022, Desebbe et al 2025) and 0.45 (Brady et al 2010, Liu et al 2021) have been utilized in Mx estimation as a likely consequence of parameter modifications in the computation of Mx. The same scatter of thresholds is present in association with the evaluation of COx and 0.3 is again one of the most utilized thresholds (Ono et al 2014, Joran et al 2020) but also 0.35 (Liu et al 2021). Even the suggestions that the threshold utilized to detect the CA dysfunction with COx should be different from that applied to Mx (Liu et al 2021, Thudium et al 2023) and the threshold applied to Mx could be exploited to decide that of COx (Ono et al 2013b) have contributed to make the rationale underlying the setting of this important parameter difficult to understand. The different setting of the threshold might have provided the basis for the contradictory results in assessing CA from the analysis of MCBv-MAP and rSO2-MAP correlation. Indeed, some studies suggested poor agreement between Mx and COx (Thudium et al 2023, Bögli et al 2025a), while others stressed their interchangeability (Brady et al 2010, Ono et al 2013b, Mol et al 2021).

In the present study we propose a short-time iterative time-domain correlation-based approach for the CA characterization. The method was applied in connection with two approaches testing on a frame-by-frame basis and individually the null hypothesis of zero or negative correlation in the planes [MAP, MCBv] and [MAP,rSO2]., namely a method based on the classical t-test checking the significant and positive departure of the Pearson correlation coefficient from 0 in linear regression (LR) analysis, and a method based on generation of a set of uncorrelated pairs of variability series preserving the distribution and power spectral density of the original variability series, but destroying their cross-correlation (Porta and Faes 2016) via iteratively-refined amplitude-adjusted Fourier transform (IAAFT) procedure (Schreiber and Schmitz 1996). Analysis was carried out in patients scheduled for cardiac surgery in basal condition (BASAL), namely at rest in supine position before general anesthesia induction, and after induction of general propofol-based anesthesia (ANESTH) before opening the chest. We tested four main hypotheses: i) the proposed approach, in association with the classical t-test applied to the Pearson correlation coefficient using Student’s t-distribution, allows a fast and robust detection of the CA impairment on an individual basis that can be utilized instead of a more time consuming surrogate approach; ii) the CA estimate obtained from the analysis of MCBv-MAP and rSO2-MAP variability interactions might exhibit similar between-condition trends at population level but their values might be unrelated at the individual level; iii) ANESTH affects cerebrovascular control by reducing the likelihood of observing positive MCBv-MAP and rSO2-MAP association; iv) the proposed MCBv-MAP approach provides similar information compared to more traditional Mx approach, while the proposed rSO2-MAP technique does not provide results in agreement with COx.

2.1. Short-time iterative time-domain correlation-based analysis to assess the CA impairment

The Pearson correlation coefficient r was computed in the plane , where n is the nth cardiac beat. The analysis was carried out over segments of 50 consecutive values of MAP and MCBv and iterated with a shift of 1 beat (i.e. with an overlap of 49 values) until the entire sequence of 250 consecutive values was covered. The choice of a moving window length of 50 consecutive values focuses the analysis on the LF band or faster times scales. For each segment, the significance of the positive MCBv-MAP correlation was tested according to two different approaches, namely LR and IAAFT techniques. The value of allowing the rejection of the null hypothesis of zero or negative correlation was retained and labeled as . The median value of and the percentage of the patterns with significant and positive MCBv-MAP correlation, labeled as , were taken as indicators of the relevance of the CA impairment. The analysis was also performed in the plane with the same parameter setting, thus computing and . To distinguish the analyses according to the strategy followed to reject the null hypothesis, we labeled the indexes as , , , , , , , and .

2.2. Detecting the CA impairment

The CA impairment was detected according to two different approaches.

The first approach was based on the classical t-test checking the zero or negative correlation applied to the Pearson correlation coefficient r in LR analysis in the plane and . Under the assumptions that the relationship between variables is linear, they are normally distributed and residuals are completely independent of each other over time, the null hypothesis was rejected with type I error probability p < 0.05, if computed over the original pairs was above the critical positive value of the Student’s t-distribution with N–2 degrees of freedom associated with the level of significance of 0.05, where N is the series length. Conversely, the alternative hypothesis, namely r exhibited a significant and positive deviation from 0, was accepted.

The second approach was based on surrogate data generation (Porta and Faes 2016). From each original [MAP, MCBv] and [MAP, rSO2] pairs we generated 100 surrogate couples. Each surrogate pair preserved the distribution of values and power spectral density of the original sequence of MAP, MCBv and rSO2 values, while they were completely uncorrelated because two independent identically-distributed random sets of Fourier phases ranging from 0 to 2π were utilized (Palus 1997). The IAAFT method was utilized to create iso-distributed iso-spectral surrogate windows of data, and the IAAFT procedure was iterated 100 times (Schreiber and Schmitz 1996). The Pearson correlation coefficient r was computed from each surrogate pair and the 95th percentile of the distribution of r was extracted. If the r computed over the original sequence pair was higher than the 95th percentile of the r calculated over surrogate couples, the null hypothesis was rejected and the alternative hypothesis of positive and significant correlation was accepted.

2.3. Assessing CA via traditional time-domain correlation-based markers

CA was estimated via traditional time-domain correlation-based analysis between CBv and AP signals and between rSO2 and AP signals by assessing respectively Mx (Czosnyka et al 1996) and COx (Brady et al 2007). We implemented the procedure reported in (Kostoglou et al 2024). Briefly, the AP, CBv and rSO2 signals were first band-pass filtered between 0.005 and 0.05 Hz, thus limiting the fluctuations of AP, CBv and rSO2 to those in the VLF band. The Pearson correlation coefficient r was computed over epochs of 30 values. Each value was obtained by averaging signals over windows of 10 s. These peculiar characteristics of Mx and COx computation limit the time scales to those in the VLF band. When r was computed over AP and CBv the index was referred to as Mx, while when r was calculated over AP and rSO2 the marker was denoted as COx. Mx and COx assessment was carried out starting from the same onset utilized for the analysis described in section 2.1. In the case the variability series of AP, CBv and rSO2 were shorter than 300 s, less than 30 averaged values were considered. Values of Mx > 0 and COx > 0, suggesting the CA impairment, were labeled as and respectively (Czosnyka et al 1996). The percentage of subjects with Mx > 0 and COx > 0 was computed and indicated as and respectively.

3.1. Experimental protocol and signal acquisition

Fifty-three patients (age: 62 ± 11 years, 39 males, 14 females) scheduled for major cardiac surgery with cardiopulmonary bypass were enrolled at the Department of Cardiac Surgery, IRCCS Policlinico San Donato, San Donato Milanese, Italy. Exclusion criteria were: non sinus rhythm at the enrollment, age under 18 years, pregnancy, or reintervention. The study was approved by the San Raffaele Hospital Ethical Review Board and authorized by IRCCS Policlinico San Donato (approval number: 06/int/2023; approval date: 15/02/2023; authorization date: 01/03/2023; ClinicalTrials registration: NCT05786274). The study was compliant with Helsinki Declaration for studies involving human subjects. Participants signed a written informed consent before participating. The sex imbalance in this study reflects the typical sex distribution in major cardiac surgery patients at IRCCS Policlinico San Donato San Donato Milanese, Italy.

Electrocardiogram (ECG) from modified lead II and invasive AP from the radial artery were acquired from the patient’s monitor, while CBv was recorded via TCD (Multi-DopX, DWL, USA) from middle cerebral artery and rSO2 through NIRS (INVOS 7100 c, Medtronic, USA) positioning the probe on the forehead of the patient. CBv and rSO2 were acquired from the same side of the patient’s head. ECG, AP, and TCD analog outputs were connected to an analog-to-digital board (National Instruments, USA) plugged into a personal computer and were synchronized with the NIRS’s serial output via a custom-made software written in LabVIEW environment (National Instruments, USA). All the signals were sampled at 400 Hz. Recording sessions lasted 10 min in BASAL and during ANESTH. Before BASAL, patients received an intramuscular injection of 0.5 mg of atropine and remifentanil at 100 µg. Then, after BASAL recording, general anesthesia was induced with propofol (intravenous bolus of 1 mg·kg−1 and continuous injection of 3 mg kg−1·h−1 for maintenance) and remifentanil (initial infusion at 0.2 µg kg−1·min−1 and maintenance of 0.05–0.5 µg kg−1·min−1). The entire anesthesia protocol was described in (Bari et al 2024). Subjects breathed spontaneously in BASAL and were mechanically ventilated at a rate of 12–16 breaths·min−1 during ANESTH with a mixture 1:1 of oxygen and air.

3.2. Beat-to-beat variability series extraction and data analysis

Heart period (HP) was computed from the ECG as the time distance between two consecutive R-wave peaks. The nth systolic peak from the AP signal was taken within the HP. The nth diastolic valley was searched after the nth systolic peak. MAP and MCBv were extracted as the integral between the (n–1)th diastolic and nth diastolic points from AP and CBv respectively divided by the interdiastolic time (Bari et al 2016). The time occurrence of the first R-wave peak defining the nth HP was utilized to sample rSO2. MAP-MCBv and MAP-rSO2 indexes were derived from 44 patients during BASAL and 28 during ANESTH. The lower number of acquired patients during ANESTH than in BASAL was due to difficulties in utilizing TCD in the surgical theater, poor quality TCD recordings, and impossibility in completing the acquisition during ANESTH due to time constraints set by surgical team. The original sex distribution in our group was maintained between the experimental conditions. Time series were corrected in the case of misdetections via a graphical interface. The impact of isolated arrhythmias was limited via linear interpolation using values unaffected by arrhythmic beats. Stability of the mean and variance was checked. Analysis was carried out over beat-to-beat series of MAP, MCBv and rSO2 recorded in BASAL and during ANESTH. Mean and variance were denoted as , , , , , , , and respectively and expressed in ms, ms2, mmHg, mmHg2, cm·s−1, cm2·s−2, %unit, and %unit2 respectively.

3.3. Statistical analysis

Normality of the distribution of the parameters was tested via the Shapiro-Wilk test. Unpaired t-test, or Mann-Whitney rank sum test when appropriate, was used to assess differences in time-domain markers of MAP, MCBv and rSO2 between BASAL and ANESTH. Two-way analysis of variance (Holm-Sidak test for multiple comparisons) was exploited to evaluate the difference between , or , derived from MAP-MCBv and MAP-rSO2 analyses within the same experimental condition (i.e. BASAL or ANESTH) and the effect of ANESTH within the same type of analysis (i.e. MAP-MCBv and MAP-rSO2). Analysis was performed separately according to the strategy followed to reject the null hypothesis of null or negative correlation (i.e. LR and IAAFT). Two-way repeated measures analysis of variance (one factor repetition, Holm-Sidak test for multiple comparisons), was exploited to test the significance of difference between , or , derived from the MAP-MCBv and MAP-rSO2 analyses within the same method to reject the null hypothesis and the impact of technique for the assessment of the threshold of significance within the same type of analysis. Two-way analysis of variance (Holm-Sidak test for multiple comparisons) was applied to compare and within the same experimental condition and to evaluate differences between experimental conditions given the same marker (i.e. or ). χ2 test was utilized to test the difference between and within the same experimental condition and the impact of ANESTH within the same marker (i.e. and ). The level of significance of the test was lowered according to the number of comparisons (i.e. 4) to account for the multiple comparison issue. After pooling together data relevant to BASAL and ANESTH sessions, the linear association was tested in the planes , , , , , , , , , , %, , , , and , . Pearson correlation coefficient r and type-I error probability p was computed. Statistical analysis was performed with a commercial statistical software (Sigmaplot v.14.5, Systat Software, San Jose, CA, USA). The level of statistical significance of all the tests was set to 0.05.

Table 1 shows time domain markers extracted from HP, MAP, MCBv, rSO2. All the parameters were reduced as expected during ANESTH due to known depression of cardiovascular and cerebrovascular controls, except for that lengthened and that remained stable.

Table 1. Time domain parameters of HP, MAP, MCBv and rSO2 variability series.

ParameterBASALANESTHμHP [ms]846 ± 151981 ± 221*σ2HP [ms2]1992 ± 2315344 ± 530*μMAP [mmHg]105 ± 1175 ± 11*σ2MAP [mmHg2]20 ± 198 ± 14*μMCBv [cm·s−1]41 ± 1527 ± 9*σ2MCBv [cm2·s−2]16 ± 195 ± 9*μrSO2 [%unit]69 ± 766 ± 13σ2rSO2 [%unit2]0.60 ± 0.290.43 ± 0.47*

BASAL = at supine resting before general anesthesia induction; ANESTH = after induction of general anesthesia and before opening the chest; HP = heart period; MAP = mean arterial pressure; MCBv = mean cerebral blood velocity; rSO2 = regional saturation of oxygen; µHP= HP mean; σ2HP = HP variance; µMAP = MAP mean; σ2MAP = MAP variance; µMCBv = MCBv mean; σ2MCBv = MCBv variance; µrSO2 =rSO2 mean; σ2 rSO2 = rSO2 variance. Results are reported as mean ± standard deviation. The symbol. *indicates p< 0.05 versus BASAL.

The vertical grouped box-and-whisker plots in figure 1 show the median of computed in the plane and in the plane in BASAL (white boxes) and during ANESTH (grey boxes). Analysis was based on LR (figure 1(a)) and on IAAFT (figure 1(b)) approaches. The decrease of during ANESTH compared to BASAL was significant only when assessed over MAP and MCBv series (figure 1(a)). was lower than only in BASAL condition (figure 1(a)). was lower than both in BASAL and during ANESTH (figure 1(b)). Both and decreased during ANESTH (figure 1(b)).

Figure 1. The vertical grouped box-and-whisker plots show the median of computed according to MCBv-MAP and rSO2-MAP analyses in BASAL (white boxes) and during ANESTH (grey boxes). The analyses are carried out according to LR (a) and IAAFT (b) techniques. The height of the box represents the distance between the first and third quartiles, with the median marked as a horizontal line, and the whiskers indicate the 5th and 95th percentiles. The symbol # indicates p< 0.05 versus MCBv-MAP association within the same experimental condition (i.e. BASAL or ANESTH), while the symbol * indicates p< 0.05 within the same type of analysis (i.e. MCBv-MAP or rSO2-MAP).

Download figure:

Standard image High-resolution image

Figure 2 has the same structure as figure 1 but it shows