A Robust Causal Brain Network Measure and Its Application on Ictal Electrocorticogram Analysis of Drug-Resistant Epilepsy

Measuring causal brain network is a significant topic for exploring complex brain functions. While various data-driven algorithms have been proposed, they still have some drawbacks such as ignoring time non-separability, cumbersome parameter settings, and poor robustness. To solve these deficiencies, we developed a novel framework: “time-shift permutation cross-mapping, TPCM,” integrating steps of: <xref rid="deqn1" ref-type="disp-formula">(1)</xref> delayed improved phase-space reconstruction (DIPSR), <xref rid="deqn2" ref-type="disp-formula">(2)</xref> rank transformation of embedding vectors’ distances, <xref rid="deqn3" ref-type="disp-formula">(3)</xref> cross-mapping with a fitting estimation, and <xref rid="deqn4" ref-type="disp-formula">(4)</xref> causality quantification using multi-delays. Based on the synthetic models and comparison with baseline methods, numerical validation results demonstrate that TPCM significantly improves the robustness for data length with or without noise interference, and achieves the best quantification accuracy in detecting time delay and coupling strength, with the highest determination coefficient (<inline-formula> <tex-math notation="LaTeX">${R}^{{2}}={0}.{96}\text {)}$ </tex-math></inline-formula> of fitting verse coupling parameters. The developed TPCM was finally applied to ictal electrocorticogram (ECoG) analysis of patients with drug-resistant epilepsy (DRE). A total of 17 patients with DRE were included into the retrospective study. For 8 patients undergoing successful surgeries, the causal coupling strength (<inline-formula> <tex-math notation="LaTeX">$0.58~\pm ~0.20$ </tex-math></inline-formula>) within epileptogenic zone network is significantly higher than those suffering failed surgeries (<inline-formula> <tex-math notation="LaTeX">$0.38~\pm ~0.16$ </tex-math></inline-formula>) with <inline-formula> <tex-math notation="LaTeX">${P}\lt {0}. {001}$ </tex-math></inline-formula> through Mann-Whitney-U-test. Therefore, the epileptic brain network measured by TPCM is a credible biomarker for predicting surgical outcomes. These findings additionally confirm TPCM’s superior performance and promising potential to advance precision medicine for neurological disorders.


A Robust Causal Brain Network Measure and Its Application on Ictal Electrocorticogram Analysis of Drug-Resistant Epilepsy
Yalin Wang , Wentao Lin , Hong Peng , Ligang Zhou , Wei Chen , Senior Member, IEEE, and Bin Hu, Fellow, IEEE Abstract-Measuring causal brain network is a significant topic for exploring complex brain functions.While various data-driven algorithms have been proposed, they still have some drawbacks such as ignoring time non-separability, cumbersome parameter settings, and poor robustness.To solve these deficiencies, we developed a novel framework: "time-shift permutation crossmapping, TPCM," integrating steps of: (1) delayed improved phase-space reconstruction (DIPSR), (2) rank transformation of embedding vectors' distances, (3) cross-mapping with a fitting estimation, and (4) causality quantification using multi-delays.Based on the synthetic models and comparison with baseline methods, numerical validation results demonstrate that TPCM significantly improves the robustness for data length with or without noise interference, and achieves the best quantification accuracy in detecting time delay and coupling strength, with the highest determination coefficient (R 2  = 0.96) of fitting verse coupling parameters.The developed TPCM was finally applied to ictal electrocorticogram (ECoG) analysis of patients with drug-resistant epilepsy (DRE).A total of 17 patients with DRE were included into the retrospective study.For 8 patients undergoing successful surgeries, the causal coupling strength (0.58 ± 0.20) within epileptogenic zone network is significantly higher than those suffering failed surgeries (0.38 ± 0.16) with P < 0.001 through Mann-Whitney-U-test.Therefore, the epileptic brain

I. INTRODUCTION
A. Background C AUSAL interaction is a fundamental problem in natural sciences, engineering, and socioeconomic systems [1], [2], [3], [4], [5].Herein, the detailed dynamics of complex systems are usually unknown, while only their generated time series are observable.Correctly identifying causal interactions among these variables is indispensable, as it provides quantitative insights to study the internal dynamics of a target system [6], [7].In neuroscience, specifically, assessing directed connectivity and constructing topological causality mapping of the brain network are crucial for exploring brain function [8], [9], [10].Causality coupling analysis has been successfully applied in brain disease classification, pathological mechanism investigation, precise diagnosis, and treatment, etc. [11], [12], [13], [14], [15], [16].The Granger causality (GC) is currently the most successful technology, providing a data-driven framework to compute causality based on their predictability relationship [17].The basic GC is a bivariate measure through an autoregressive model, with some later variants extending to multivariate systems using multivariate autoregressive models (MVAR) and proposing conditional Granger causality (cGC) algorithm [18], [19], [20].Transfer entropy (TE) is another development in nonlinear causality estimation by calculating conditional mutual information differences [21].Together, TE and its variants, partial transfer entropy (pTE), symbolic transfer entropy (STE), frequency domain transfer spectral entropy, etc., facilitate nonlinear causal discovery among variables [20], [22], [23], [24].However, all these methods, GC, TE, and their variants, are based on the assumption of independent separation of system variables and lead to inevitable information loss in the computation process.In realistic coupling systems, time-nonseparability (TNS) is prevalent among system variables [7], [25], [26].
By introducing nonlinear attractor reconstruction, convergent cross-mapping (CCM) has been recently developed and implemented for causality estimation in non-separable dynamic systems [7].Nowadays, CCM-inspired techniques have gained widespread attention in many fields [25], [27], [28], [29], [30], [31].For example, it has been used to analyze functional connectivity during the resting state for discovering default brain networks [32].Brain-heart interactions are also detected by CCM measure, exploring how the human's cardiovascular system is coupled with the central nervous system [28].CCM is further applied to explore the underlying pathological mechanisms of certain neurological disorders, facilitating their precision medicine approaches [29], [30], [31].

B. Problems Statement
TNS is invited by Takens's embedding theory, which reveals that in a deterministic nonlinear dynamic system, the information of the entire dynamic system could be generally embedded into only one single variable and thus could be reconstructed by the observation data of that variable [7], [26].When predicting the dynamic system, removing some variables' information from others is usually impossible [26].For example, the human brain is an extremely complex networked dynamics system, where the activities of neuronal populations are closely interconnected, including physical synaptic connections and functional coupling (in fact, the intrinsic laws are still poorly understood), which together constitute a complex, sophisticated brain function characterized by typical non-separability.TNS indicates that in such a complex system, the effect variable is unlikely to be fully separated from the cause variables; otherwise, it will lead to a failed estimation [30], [33], [34], [35].All traditional measures, including GC-, TE-based algorithms ignore TNS premise, leading to failed quantification.They follow the premise of separation of variables: variable X is regarded as cause effect on Y if the predictability of Y declines when X is removed from the universe of all possible causative variables [7].The rigorous concept of TNS can also be found in [26].
Herein, we give a numerical example to illustrate it.Consider a typical 2-species logistic coupling model [7]: where C is the coupling coefficient, and an explicit causality Y → X is defined.The TE and GC results are represented in Fig. 1.Moreover, significance tests are performed to determine the existence of a causal relationship between variables [8], [23], [36].Herein, we computed the 100 surrogate causalities Y ′ → X , where Y ′ is surrogate series of Y after random disordering [37], [38].When data length is 1000, for GC and TE calculation, there are no significantly difference between raw causality Y → X and 100 surrogate causalities Y ′ → X , thus they fail to identify the true causality.In TE and GC calculation, X and Y are considered as two mutually separated variables, and Y is removed from the set of hypothetical causal variables without reducing X s predictability; therefore, they generate an incorrect estimation: Y does not cause X .
CCM starts with phase space reconstruction (PSR) from Taken's embedding theorem.For a multivariate complex system, each variable's embedding dimension and delay need to be considered simultaneously, and comprehensively to characterize system dynamics more accurately with low redundancy.Considering a time series X = x (t) , t ∈ [1, L], the corresponding reconstruction manifold M X is constructed as: where X m,δ i = {x (i) , x (i + δ) , . . ., x (i + (m − 1) δ)} is with an embedding dimension of m and a delay coefficient of δ.Various methods are used to determine parameters m and δ, however they often yield different results, leading to instable causality results.And the process of determining the parameters is also tedious and time-consuming.These shortcomings regarding the parameter settings inspired our study to introduce a novel parameter-free phase space reconstruction method.
Referring to CCM, based on the reconstruction manifold M X and M Y , a cross-mapped estimate of Y (t), denoted by Ŷ (t) |M X is generated to compute causality coupling from X to Y .It introduces exponentially weighted distances from nearby points on a manifold M X by false nearest neighbors (FNNs) estimation.However, CCM leads to poor robustness because short series cannot produce a refined prediction Ŷ (t) |M X by FNNs information.Consider the system modeled by Eq. ( 1) with data length ranging from L = 200 to L = 4000 samples.As the L increases, the CCM results grow up from around 0.22 to 0.8.It implies that for the same dynamics system, the observed series's length can lead to large differences in causal quantification and even false estimates.Meanwhile, the distributions of observed raw series Y and predicted series Ŷ |M X by cross-mapping between manifolds when data length L = 200 and L = 4000 are represented in Fig. 2. Obviously, the CCM results of short and long series vary greatly, and short series lead to unrefined cross-prediction due to the lack of state space information via FNNs.
In conclusion, data-driven causality measures still face the following limits: (1) traditional GC, TE, and their variants lead to poor effectiveness ignoring the TNS properties of complex systems; (2) recent CCM-inspired measures still suffer from the limits of cumbersome parameter setting, poor robustness, and poor effectiveness in detecting coupling strength and time delay.

C. Contributions and Structure of This Study
Finally, we apply the TPCM method to a real-world neuroelectrophysiological dataset.As a retrospective study, 17 DRE patients with totally 55 seizures were included.8 patients underwent successful surgeries, while other 9 patients underwent failed surgeries.Based on the TPCM calculation and statistical analysis, we examined multi-frequency ECoG, and characterized the causal coupling dynamics in the epileptic brain network involving surgical outcomes.Meanwhile, a comparative study with other baseline methods is also included.The specific contributions of this study can be summarized as follows.
(1) We introduced a novel causality algorithm, called"TPCM".It is a parameter-free measure that eliminates the time-consuming and inconsistent parameter settings, meanwhile maintaining excellent quantification accuracy.
(2) TPCM also achieves the best quantification robustness.It proposes incorporating the normalized rank ordering of pairwise distances in two manifolds instead of the traditional FNNs prediction [39], [40].
(3) TPCM method was also applied to analyze DRE's ECoG brain network.It subsequently can be found that the causal coupling within EZs were closely related to the surgical outcomes, which were expected to be a more favorable biomarker for predicting surgical outcomes of DRE.
The remainder of this paper is structured as follows.Section II explains the proposed TPCM.Section III introduces synthetic datasets used for validation experiments.Section IV presents the numerical experiment results.Section V gives TPCM's pilot application on the ECoG analysis.Finally, Section VI discusses and concludes this study.

II. METHOD
Inspired by a recent study about an embedding method, we introduce a novel parameter-free algorithm to generate a delayed attractor manifold in this study.Considering a dynamics system with time series To estimate the causality X → Y , set the time delay as τ .The delayed sub-series are generated as: Based on the delayed sub-series X τ , the attractor manifold M τ X of variate X is generated as follows.The corresponding embedding vectors are calculated with constant embedding parameters δ = 1 and m = 2: x (3) . . .
where X m,δ i = {x (i) , x (i + δ) , . . ., x (i + (m − 1) δ)} Then the covariance matrix C can be obtained as: C has been performed eigendecomposition as: where is a square matrix with each column being eigenvectors, is a diagonal matrix with the diagonal elements being the eigenvalues.Finally, project M τ X to its principal directions and obtain the final attractor manifold: where The above steps, Eq. ( 4) ∼ (7), consider the time delay effect in causal coupling.They do not require the choice of embedding parameters, avoiding the tedious process of parameter selection and quantification bias resulting from different embedding parameter settings.The above process is termed "delayed improved phase-space reconstruction (DIPSR)".The series Y = y(t) can also be processed by DIPSR to obtain the attractor manifold: Then we calculate pairwise Euclidean distances between all the points in the delayed attractor manifold Mτ X and Mτ Y : Sparsify the distance matrices E X τ and E Y τ , remove each distance matrix's trivial elements surrounding the central diagonal, and then the pairwise Euclidean distances are retained: where U pper T riangular is the upper triangular matrix, k is the first diagonal meets that std Sort the remained distances E X τ and E Y τ in ascending order and obtain the index ranks, I X τ and I Y τ .Then, normalize these retained distances to range (0,1] by the corresponding ranks: where R = 1 : length (I ) length (I ).Get the cross-mapping of normalized ranks by other variables, Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.which has been proven a more reliable indicator of manifold structure [40]: Obtain the square error (SE) between the true ranks and cross-sorted ranks: where R ≤ T , T is an adaptive threshold that meets: . Where, Y τ ′ is the first derivative.
Normalize the square error S E (R X τ →Y τ ) by the expected error of randomly distributed ranks: , where uncorrelated ranks are defined as To create a smoother estimate of asymmetric divergence between the manifolds, compute the cumulative average of N S E (R X τ →Y τ ) as Finally, by fitting the polynomial function [N S E (R X τ →Y τ )], the TPCM causality quantification with time delay τ is obtained by y-intercept of fitting function: This y-axis intercept quantifies the causal coupling strength because it captures the behavior of the error function [N S E (R X τ →Y τ )] when distance reaches minimal, revealing the topology for testing the smooth estimation of asymmetric divergence between two manifolds.A series of time delays are set as {τ 1 , τ 2 , . . .τ i , . . .τ max } and the corresponding delayed causality results are computed as: T PC M τ i (X → Y ).Based on the findings of existing studies, in data-driven causality estimation, when the delay τ i is equal to the true delay of a coupled system, the corresponding calculated causality result will reach maximum.Therefore, in this study, when the TPCM curve reaches an extremum value point, the delay index at this point is the system's true delay.The final causality X → Y can also be obtained: Flowchart of the proposed TPCM is shown in Fig. 3.

III. NUMERICAL EXPERIMENT
To verify the TPCM's performances, we carried out a series of numerical experiments.All operation codes run on a PC with Intel Core i5-9400F @ 2.9 GHz, 24.0 GB RAM, Windows 10 operating system, and Matlab R2020a (The MathWorks Inc.).
The first synthetic dataset (1) in numerical experiments is a discrete-time multivariate logistic coupling model [34], [41], [42]: where r x = 3.78, r y = 3.77, r z = 3.75, β x,y , β y,x and β z,y are coupling parameters controlling causal connection strength.ε X , ε Y and ε Z are zero-mean white Gaussian noise, and their standard deviation σ will be referred to as the system's noise level.The system noise in either variable in Eq. ( 16) propagates to later values of the variable and to another variable through coupling terms.τ 1 , τ 2 and τ 3 are intrinsic time delays.
The second synthetic dataset (2) is a nonlinear strongly-coupled dynamics system, modeled as [8] and [37]: Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

A. Numerical Experiment Results of Dataset (1)
Based on dataset (1), we first test the TPCM's robustness for data length.We considered two coupling scenarios, unidirectional and bidirectional causality connections, and then performed 100 implementations by randomizing the model's initial values at each data length condition.The results calculated by the baseline methods and developed TPCM are shown in Fig. 5.In addition, we also computed the surrogate causality X ′ → Y and Y ′ → X , where X ′ is surrogate data after randomly disordering the raw series X .The TPCM results present the best robustness and highest stability from short to long series.In contrast, traditional methods show poor robustness and high relative variation, meaning that these methods require sufficiently long time series when quantifying causality, otherwise leading to misleading results.TE, pTE, GC and cGC measures cannot effectively detect bidirectional causality, demonstrating these methods' failure in nonlinear non-separability systems.Coefficients of variation (CV) are also calculated for all results, defined as the mean value of 100 realizations divided by the standard deviation.CV is an effective index that can accurately quantify causality estimates' consistency.The CV Noise interference is another important consideration for causal estimation algorithms.However, previous studies have not reported the algorithms' robustness disturbed by system noise.In this study, the TPCM's robustness is further explored by comprehensively considering system noise interference, using the dataset (1) defined by Eq. ( 16), where σ X = σ Y = σ Z ∈ {0.001, 0.002, 0.003, 0.004, 0.005}.Other numerical experimental settings are detailed Fig. 4.
The causal coupling of each dynamics model is calculated using different algorithms, and 100 quantification results by randomizing initial value of each model are obtained separately, with the mean ± std as shown in Fig. 7.The results indicate that as noise's effect increases, the TPCM results gradually decrease.For coupling X → Y, when the noise intensity increases from σ = 0.001 to σ = 0.005, the average result of the causal quantification TPCM(X → Y) gradually decreases from 0.85 to 0.68.Similarly, the average result of TPCM(Y → X) decreases from 0.6 to 0.12.These findings align with expectation that the noise injected into the dynamical system disrupts the causally driven structure, and the noise is more destructive to the weaker coupled system.Moreover, the TPCM results are consistent for each noise condition over the entire data length range.For example, for σ X = 0.001, TPCM(X → Y) results always remain near 0.9, and TPCM(Y → X) results remain near 0.6.In summary, the developed TPCM method maintains the optimal robustness under noise interference.The last numerical experiment using dataset (1) is time delay detection in the coupling system.Traditional CCM algorithm cannot detect time delays, even lead to false estimates.Recently, a time delay version inspired by CCM has been developed, termed "CCM delay " herein for a comparative study.
Based on dataset (1), three typical delayed causality scenarios are generated, i.e., unidirectional coupling system, bidirectional causality with different delays, and bidirectional causality with the same delay.Other specific settings are set, as shown in Fig. 4. For comparison purposes, the TPCM and CCM delay methods are included in the calculations, and the results are shown in Fig. 8.
The unidirectional coupling model is set with three delays τ 2 = 2, 4, 6, and the TPCM(X → Y) result curves in Fig. 8 peaked at time points 2, 4, and 6, respectively.Therefore, TPCM method accurately detected the true coupling delay.In contrast, the peak point of CCM delay (X → Y) curves are lagging the actual delay.Namely, CCM delay cannot accurately detect the causality delay.The reasons for this result may include the existence of autoregressive lag order in raw series, the estimation error of the simplex projection, etc.For bidirectional model, the TPCM method can also accurately detect the accurate delays τ 1 = 0,τ 2 = 2, 4, 6, while the detection accuracy of CCM delay is limited, with the results larger than actual delays.Finally, for the third delayed causality scenario, TPCM can also accurately detect all delays and guarantee the quantification accuracy of both causal strengths.In the present numerical experiments, the causal strength determined by the coupling coefficient relationship (β yx = 0.08 >β x y = 0.02) is correctly quantified: T PC M(X → Y ) > T PC M(Y → X ).In summary, TPCM achieves the best performance in delayed causality detection.

B. Numerical Experiment Results of Dataset (2)
Dataset ( 2) is a widely used synthetic model to test nonlinear strongly coupled systems verified in TE and GC-based Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

TABLE I PATIENT DEMOGRAPHICS AND CLINICAL INFORMATION 9. Robustness results for data length in estimating dataset (2).
In each validation condition, all methods carry out 100 implementations by randomizing the synthetic model's initial values.
measures.Herein, we first test the TPCM's robustness using dataset (2).The results in Fig. 9 show that the TPCM achieves the best robustness.TPCM results always remain stable when series's length increases from 200 to 4000 samples; meanwhile, it successfully quantifies two different causal strengths, X → Y,C 1 = 0.9 and Y → X,C 2 = 0.3.
The CCM, TE, and pTE methods are also able to detect the correct causality, but these results are less stable, such as the CCM(X → Y) result curve gradually increases from about 0.4 to 0.8, and the relative change rate reaches 0.8−0.→ Y), pTE(X ′ → Y), with quantification accuracy of 65% and 31% respectively.Similarly, the pTE method suffers from incorrect estimation when data lengths below 1000.As for GC and cGC, they both fail to measure nonlinear causality, showing no significant difference between surrogate causalities and estimated causalities.
By varying the coupling parameters of synthetic model, we further test the TPCM method's quantification accuracy.When parameter C 1 increases from 0.1 to 0.9 with a step of 0.1, the causality results are calculated by TPCM and other baseline methods.Theoretically, the causality results will  rise with increasing parameter C 1 linearly.The mean value of 100 implementations is sketched versus the parameter C 1 along with its standard deviation as an error bar.The causality results against C 1 along with the fitted curve are plotted in Fig. 10.The coefficients of determination (R 2 ) are calculated to check the fitting, which indicates the precision of coupling estimation.
It can be found that the TPCM method can accurately identify causal relationships with an R 2 of 0.96, which is higher than other methods.The coefficients of determination for the CCM and TE methods were R 2 = 0.94 and R 2 = 0.95, respectively, which also confirmed the good quantification performance of these two methods to some extent.However, when the coupling parameters are low, i.e., for weaker causal relationships, all methods except TPCM suffer from false estimation.As marked by the dashed box in Fig. 10, the statistical results of CCM(X → Y) and the shuffled causal estimates CCM(X ′ → Y) partially overlap, which means that some of the results of the CCM will be judged as no causality, i.e., the CCM cannot infer the real existing causal coupling X → Y, and when C = 0.1, the CCM's causal quantification accuracy are only 75% and 78% when C = 0.2.Similarly, when the parameter C = 0.1, the TE method also suffers from low quantification accuracy.And the and cGC methods are invalid in these nonlinear causal systems.This result is consistent with previous studies [7].Therefore, in summary, the developed TPCM achieves the highest quantification accuracy.
Finally, time delays of dataset (2) are set τ 1 = 3, 4, 5 to validate the TPCM's quantitative accuracy in delayed models.Fit the coupling parameter C and detected delayed causality (i.e., the maximum point of the causal coupling curve for each coupling parameter condition) and calculate the corresponding R 2 .The results are shown in Fig. 11.In each C condition, TPCM can accurately detect both time delay and coupling strength, while the delays detected by the CCM delay lag behind the true delays.In addition, TPCM's R 2 results are 0.93, 0.93 and 0.92, respectively, confirming the TPCM's optimal quantification performance in delayed systems.

Epilepsy is one of the most common neurological disorders, affecting approximately 70 million individuals worldwide.
There is approximately one-third of epilepsy patients who are resistant to antiepileptic medications [45], called "drug resistant epilepsy (DRE)".For patients with DRE, neurosurgery is the most effective treatment.Successful surgery requires complete removal or disconnection of the epileptogenic zone (EZ), while current surgical success rates are extremely low.Only approximately 50% of patients with DRE who undergo surgery achieve seizure freedom.Surgical failure can lead to potential consequences such as permanent nerve damage, persistent seizures, etc.Therefore, predicting surgical outcomes has attracted considerable attention among clinicians.
As shown in Fig. 12, clinically annotated epileptogenic zone (EZ) is surgically resected, (or thermally ablated).A successful Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
surgery entails the complete removal of the actual epileptogenic zone, whereas a failed surgical treatment may leave residual epileptogenic zones.Additionally, the ictal ECoG signals sampled from areas in contact with the EZ do not exhibit significant differences.Quantitative analysis of individual ECoG channels in isolation poses challenges in effectively advancing precision medicine for DRE [45].
Epilepsy is growingly considered as a network-related disorder [45].Network neuroscience has constituted a more powerful methodology in precision medicine for DRE [46], [47], [48], [49].Causal brain network, defined as the directed influence (or causal driving effect) that one brain region or neuron exerts over another [2], has become a useful quantitative method.Seizure generation is believed to be a process driven by the epileptic network; thus, understanding the causal connectivity within epileptogenic zone (EZ, a clinically annotated region which are surgically resected), is crucial in revealing the triggering mechanism of epileptic seizures originating.Characterizing causal connectivity in epileptic network would considerably facilitate accurate prediction of DRE surgical outcomes.
In this paper, we apply the TPCM method to a retrospective study on ECoG epileptic network.A public dataset [45], totally 17 DRE patients from two clinical centers with 55 seizures, were included in this retrospective cohort study.During ECoG detection, each DRE patient will have a number of seizures, usually lasting tens of seconds to several minutes.For the 8 successful surgeries, there are 29 seizures were recorded by ECoG.For the 9 failed surgeries, there are totally 26 seizures.
At all centers, data were recorded using Nihon Kohden or Natus acquisition system with a typical sampling rate of 1000 Hz.For each patient, the clinically annotated EZ is defined as the anatomical area to be treated (i.e., resected).Surgical outcomes were classified by epileptologists using the Engel Surgical Outcome Scale classification system.Successful surgical outcomes are defined as free of disabling seizures (Engel class E I) and failed outcomes as not free of disabling seizures (Engel classes, E II ∼ E IV) at 12+ months postoperation.Furthermore, the number of seizures per patient varies from 1 to 6.The detailed demographics and clinical information are listed in Tab.I.

B. Causal Brain Network Results
Using the TPCM method, we calculated the causal coupling between each two preprocessed ECoG in the EZ regions.To more intuitively present the specific analysis of this study, we firstly give the full-band ECoG causal network calculated by our TPCM method for two typical DRE patients, namely JH108 with failed surgery, and Pt2 with successful surgery.JH108 have 8 EZ contacts, including rdi1∼rdi4, rsi1∼rsi4, and Pt2 includes EZ contacts of PST1∼PST4, AST1∼AST2, MST1∼MST2.It can be found from Fig. 13 that the causal strength of successful surgery is significantly higher than that of failed surgery.It provides a preliminary finding of the difference in EZ brain network between the successful and failed groups.
Similarly, we analyzed all patients' causal EZ networks using the proposed TPCM method, including all frequency bands ECoG: Full, δ, θ, α, β, γ , HF.All the causal networks of 17 DRE patients are divided into two groups: the successful surgery group and the failed surgery group.Corresponding statistical violin plots are shown in Fig. 14   Other baseline methods, CCM, TE, pTE, GC and cGC were also applied to calculate the same ECoG recordings.Their specific analysis is exactly the same as the TPCM results.The P values obtained by Mann-Whitney-U-test when successful surgery group Vs.failed surgery group are listed in Tab.II.By comparison, the TPCM results demonstrate the most significant statistical difference, confirming its best performance in predicting DRE's surgery outcomes.In addition, the CCM method achieves a sub-optimum with P = 5.02e − 04.

VI. DISCUSSION AND CONCLUSION A. Implications for Causal Brain Network Measure
Quantifying causality coupling from observed time series is a significant topic in studying complex physiology system.The commonly used methods still suffer from the following defects, cumbersome parameter setting, poor robustness, and accuracy.These limits seriously affect further development in study complex physiology system.This study develops a novel TPCM algorithm to solve those problems and improve quantification performances.The four steps, DIPSR, rank transformation of embedding vectors' distances, crosswith a fitting estimation, and multi-delays causality quantification, together construct TPCM computation.Two datasets in this paper are widely used in many previous studies, and are provenly reliable.The baseline methods for comparable studies are the most mainstream data-driven algorithms at present, including CCM, TE, pTE, GC, cGC, etc. Synthetic datasets confirm the effectiveness and superiority of developed TPCM method.In each numerical experiment condition, we perform 100 calculations by randomizing the synthetic models' initial values, and retain the statistical value (mean±std) of 100 values as the result of one certain experimental condition.Such scheme eliminates accidental errors and ensuring the reliability of the numeric verification.
In addition, the computation cost is also of interest.In fact, we aim to achieve the best causal quantification without significantly increasing the running time.Here we list the running time of all methods when computing both dataset (1) and (2) modeled by Eq. ( 16) and (17).The two datasets are unified as unidirectional coupling X → Y , with the data length L of 1000 and 4000 samples, respectively.Other parameters are set as follows, Eq. ( 16): β x,y = 0, β y,x = 0.08, τ 1 = τ 2 = σ = 0; Eq. ( 17): τ 1 = τ 2 = σ = 0, C1 = 0.9, C2 = 0.3.To eliminate chance errors, the calculations are repeated 100 times for each model, and running times are summed up in the Tab.III.Compared with the methods based on phase space reconstruction (CCM, TE and pTE), the running time of TPCM is significantly reduced.Although TPCM takes longer to run compared to the Granger-based methods, this trade-off is acceptable given the superior performance of TPCM.
Some potential efforts may be made to promote further TPCM-related studies.The current study only experimented on limited datasets (1) ∼ (2).Applying TPCM method to various synthetic datasets, e.g., coupling Lorenz system [50], neural mass model [23], [50], etc., is necessary.Meanwhile, more numerical verification indexes, e.g., sensitivity and specificity are also suggested.Frequency-domain expansion of traditional CCM method, termed "FDCCM", has been shown better general applicability and resilient to noise in measuring causal network [30].In future studies, the frequency-domain TPCM method will be developed tentatively.

B. Implications for ECoG Epileptic Network
Based on the network neuroscience framework, causal brain network analysis has demonstrated its advantages in DRE Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
research [51], where neurophysiological signals, fMRI, scalp EEG and iEEG are most frequently analyzed.As for the DRE study about their surgical outcomes, previous publications of iEEG analysis did not adequately account for causal interactions between multiple brain areas, especially within the clinically annotated EZs closely related to DRE's treatment [52], [53], [54], [55], [56].Systematic analysis of iEEG causal brain networks in all frequency bands within EZ regions is warranted, potentially facilitating DRE precision medicine.In this retrospective study, we established the DRE's causal epileptic network involving multi-frequency ECoG and EZ regions.Quantitative results using the developed TPCM together with statistical analysis have been performed.We found that the causal brain network strength in the EZ regions can clearly distinguish surgical success from failure.For DRE patients, there represents an anomalous network dynamic in epileptic focus.For patients with failed surgery, the resected EZ is only part of or completely not the true epileptic focus.
Comparing the causal network results of two groups (successful and failed surgery), typical frequency-band network, such as HF, have stable statistical significance, and the EZ network strength of success group is significantly higher than that of failure group.To some extent, it also supports the unique effect of high-frequency oscillations (HFOs) on identifying epileptic focus.Research on DRE's HFOs analysis is currently very manifold [53], [57], [58].It is worth noting that the ECoG calculation of the full frequency-band shows the largest statistical difference, meaning that neuroelectrophysiological activity of each frequency band may be involved in the triggering mechanism of epileptic seizures.Therefore, the full frequency-band analysis is more effective than single feature band (δ, θ, α, β, γ , etc.,).
Finally, future work will recruit more DRE subjects and may extend to multicenter clinical studies, including more types of DRE (CC1 ∼ CC4), and other neuroelectrophysiological recordings, such as stereoelectroencephalogram (SEEG).Moreover, the brain network dynamics will be combined with artificial intelligence techniques to more precisely predicting surgical outcomes for DRE patients.

C. Conclusion
In conclusion, to improve the data-driven causality algorithm's universality, robustness and accuracy, we propose a novel almost-parameter-free framework, TPCM, under the premise of time non-separability of complex systems.Meanwhile, the developed TPCM eliminates the cumbersome parameter settings and simplifies the computational process.Numerical experiments using two synthetic datasets have tested the TPCM's superior performances.The TPCM is finally applied to ECoG analysis of patients with DRE.The EZ network's causal strength of successful surgeries is significantly higher than that of failed surgeries with P = 0.00036 through Mann-Whitney-U-test.These findings further confirm TPCM's promising potential to advance precision medicine for neurological disorders.

Fig. 1 .
Fig. 1.Numerical examples of failed estimation in time non-separability system using TE and GC method.The error bars indicate the standard deviation of every 100 calculations by randomizing the initial value.

Fig. 2 .
Fig. 2. Numerical examples of the CCM method indicate the poor robustness; The error bars in the graph indicate the standard deviation of 100 calculations.

Fig. 4 .
Fig. 4. (a) Example time series of dataset (1) ∼ (2), and (b) validation schemes based on the baseline methods, CCM, CCMdelay, TE, pTE, GC, cGC.Note that, in each simulation, 100 implementations for each specific condition are performed by randomizing the model's initial values to avoid chance errors.

Fig. 5 .Fig. 6 .
Fig. 5.Robustness results for data length in causal estimation using dataset(1).In each validation condition, six methods carry out 100 implementations by randomizing the logistic model's initial values and then recording the mean value and standard deviation of the 100 implementation results.The surrogate causalities are also calculated and labeled as X ′ → Y and Y ′ →X.

Fig. 7 .
Fig. 7. Robustness results of different methods interfered by systematic noise: (a) the causality results of X → Y and (b) the causality results of Y → X.The error bars in the graph indicate the standard deviation of 100 calculations.

Fig. 8 .
Fig. 8. Time delay detection results of TPCM and CCM methods in unidirectional and bidirectional logistic models.The shading area indicates a standard deviation of every 100 realizations.

4 0. 4
×100% = 100%, TE(X → Y) gradually rises from 0.06 to 0.28, with a relative change rate of 0.28−0.060.06 ×100% = 367%, pTE(X → Y) gradually increases from 0.02 to 0.135, and the relative change rate of 0.135−0.020.02 ×100% = 575%.In addition, for shorter series, such as for a length of 200, the error bars of the TE(X → Y) and pTE(X → Y) curves are overlapped with the results of the random shuffle substitution data TE(X ′

Fig. 10 .
Fig. 10.Estimated causality results using six different methods.The solid curves are mean values with error bars versus different parameters C 1 ∈ [0.1,0.9], and dashed lines indicate the fitted curves.The coefficients of determination (R 2 ) are also recorded in each figure.

Fig. 12 .
Fig. 12. Schematic of successful and failed surgeries, and sample ictal ECoG signals in clinically annotated EZ.

Fig. 13 .
Fig. 13.The sample results of full-band ECoG (0.5 ∼ 500 Hz) causal network calculated by the TPCM method.(a) a failed surgery of patient JH108; (b) a successful surgery of patient Pt2.

Fig. 14 .
Fig. 14.Statistical violin plots of causal EZ networks from all subjects with 55 seizures.For comparison, in each frequency-band plot, there are successful surgeries (the left violin plot filled with saturated color) and failed surgeries (the right violin plot filled with unsaturated color).The marker " * " indicates the P value through Mann-Whitney-U-test. * : P < 0.05; * * : P < 0.01; * * * : P < 0.001.

TABLE III RUNNING
TIMES OF ALL METHODS (SECOND)