Fractal and Multifractal Properties of Electrographic Recordings of Human Brain Activity: Toward Its Use as a Signal Feature for Machine Learning in Clinical Applications

The quantification of brain dynamics is essential to its understanding. However, the brain is a system operating on multiple time scales, and characterization of dynamics across time scales remains a challenge. One framework to study such dynamics is that of fractal geometry; and currently there exist several methods for the study of brain dynamics using fractal geometry. We aim to highlight some of the practical challenges of applying fractal geometry to brain dynamics—and as a putative feature for machine learning applications, and propose solutions to enable its wider use in neuroscience. Using intracranially recorded electroencephalogram (EEG) and simulated data, we compared monofractal and multifractal methods with regards to their sensitivity to signal variance. We found that both monofractal and multifractal properties correlate closely with signal variance, thus not being a useful feature of the signal. However, after applying an epoch-wise standardization procedure to the signal, we found that multifractal measures could offer non-redundant information compared to signal variance, power (in different frequency bands) and other established EEG signal measures. We also compared different multifractal estimation methods to each other in terms of reliability, and we found that the Chhabra-Jensen algorithm performed best. Finally, we investigated the impact of sampling frequency and epoch length on the estimation of multifractal properties. Using epileptic seizures as an example event in the EEG, we show that there may be an optimal time scale (i.e., combination of sampling frequency and epoch length) for detecting temporal changes in multifractal properties around seizures. The practical issues we highlighted and our suggested solutions should help in developing robust methods for the application of fractal geometry in EEG signals. Our analyses and observations also aid the theoretical understanding of the multifractal properties of the brain and might provide grounds for new discoveries in the study of brain signals. These could be crucial for the understanding of neurological function and for the developments of new treatments.

1. Introduction

Brain dynamics are non-linear and are often considered as one of the most complex natural phenomena, involving several different and interacting temporal scales. For example, fast electric activity, slower chemical reactions, and even slower diffusive processes have been observed in the brain. Interestingly, brain dynamics have also been characterized as “scale-free” (Stam and de Bruin, 2004; Fraiman and Chialvo, 2012), meaning that certain signal properties stay preserved across different time scales. To describe and quantify such time scale invariant dynamics, the framework of fractal geometry is often applied (Werner, 2010).

Fractals have two specific properties: they consist of parts that are similar to the whole—termed self-similarity, and they have a fractional Hausdorff-Besicovitch dimension, also called fractal dimension (FD) (Mandelbrot, 1982; Feder, 1988; Falconer, 2003). Fractal geometry has been applied to the study of temporal dynamics, such as human brain dynamics in health (Lutzenberger et al., 1995; Pereda et al., 1998; Eke et al., 2000, 2002; Linkenkaer-Hansen et al., 2001; Bullmore et al., 2003, 2009; Gong et al., 2003; Acharya et al., 2005; Bassett et al., 2006, 2010; Hsu et al., 2007; Van De Ville et al., 2010; Papo et al., 2017) and disease (Esteller et al., 1999; Gómez et al., 2009; Zappasodi et al., 2014), providing intriguing results. For example, FD has been shown to vary prior to and during epileptic seizures (Esteller et al., 1999).

Objects adequately characterized by a single fractal dimension are referred to as monofractals. However, the fractal formalism has to be extended to capture certain phenomena that cannot be described by a single fractal dimension; these are called multifractals (Stanley et al., 1999). Multifractal objects can be conceived as decomposable into different subsets or parts, each characterized by its own distinct fractal dimension. The subsets are more precisely described as different statistical moments, and a multifractal is an object where the fractal dimension depends on the statistical moment being examined (Mukli et al., 2015). Thus multifractal objects are often described by a spectrum, showing the subsets/statistical moments and their corresponding fractal dimensions. Some natural phenomena exhibit multifractal patterns in space, for example, turbulence (Meneveau and Sreenivasan, 1987; Chhabra and Jensen, 1989), soil composition (Miranda et al., 2006; Zeleke and Si, 2006; Vázquez et al., 2008; Paz-Ferreiro et al., 2010a,b); and in time, for example heart beat patterns (Ivanov et al., 1999), and human physical activity (França et al., 2019).

There is also considerable evidence that brain dynamics are multifractal (Suckling et al., 2008; Ihlen and Vereijken, 2010; Ciuciu, 2012; Zorick and Mandelkern, 2013; Zhang et al., 2015; Papo et al., 2017; Xue and Bogdan, 2017; Racz et al., 2018). At the very least, additional statistical moments appear to be required, to characterize such dynamics (Fraiman and Chialvo, 2012). Furthermore, it is known that interacting processes with different time scales, similar to those observed in the brain, can generate multifractal patterns (Argoul et al., 1989; Suckling et al., 2008).

To measure the multifractal spectrum in brain dynamics, Multifractal Detrended Fluctuation Analysis (MF-DFA) (Kantelhardt et al., 2002) is currently the most used approach (Ihlen, 2012; Zhang et al., 2015). However, more advanced and potentially more stable estimation techniques have been proposed, such as the Multifractal Detrended Moving Average (Xu et al., 2017), and Chhabra-Jensen approaches (Chhabra and Jensen, 1989). These techniques, to our knowledge, however, have not yet been evaluated with brain signals. In addition, there are several parameter choices to be made for the purpose of the analysis. For example, to capture time-varying changes in multifractal properties, the epoch length and sampling frequency have to be chosen. These parameters may impact the multifractal estimation (Eke et al., 2002), but, to date, have not been studied systematically in the context of brain dynamics.

The biggest gap in the literature so far, however, is how multifractal properties relate to existing time series signal measures of brain dynamics (e.g., variance of the signal, band power, etc.). A major concern is that complex methods of analysis may not offer a significant advance over simpler, already established methods—this is crucial to a putative feature in machine learning applications, e.g., seizure prediction or detection (Mormann et al., 2007; Freestone et al., 2015; Brinkmann et al., 2016; Baldassano et al., 2017; Karoly et al., 2017; Kuhlmann et al., 2018a,b; Varatharajah et al., 2018). For example, in the analysis of the electroencephalogram of epileptic seizures, complex methods were found to actually reproduce patterns detected by simpler measures such as variance of the signal (Martinerie et al., 2003; McSharry et al., 2003). It is therefore, essential to understand how the (mono- and multi-) fractal measures relate to more traditional measures, and if new features can be obtained from the signal by applying a mono- or multi-fractal formalism.

To summarize, there is a knowledge gap in three critical areas: (1) which (multi)fractal characterization methodology is best suited for brain signals? (2) what are the optimal estimation parameters (e.g., in terms of recording epoch length) of potentially time varying multifractal properties? (3) what is the relationship between (multi)fractal properties and more traditional and established time series signal measures? To address these questions, we chose to analyse intracranially recorded human electroencephalography (icEEG) data, due to its high temporal resolution and high signal to noise ratio.

2. Materials and Methods

To address the questions above, we will first outline four experiments. We will then provide details on monofractal and and multifractal estimation methods, and also show how time series data with known mono- and multifractal properties can be generated to test the performance of the estimation methods. To test the multifractal measures on real-life brain signals, we then applied our analysis on human intracranial EEG. Thus, finally, we will summarize the EEG datasets used in this work.

The original scripts used in this work are available in https://github.com/yujiangwang/MultiFractalEEG. In addition, the following software packages were used: MATLAB; R (R Core Team, 2017); and ggplot2, R.matlab, reshape2, PerformanceAnalytics, and RColorBrewer (Wickham, 2007, 2009; Neuwirth, 2014; Peterson and Carl, 2014; Bengtsson, 2016).

2.1. Experiments

2.1.1. Experiment 1: Monofractal Estimation With Respect to Changing Signal Variance

Estimation of monofractal properties has been applied to EEG signals in the past with varying and often contrasting results (Esteller et al., 1999; Li et al., 2005). A particular concern is that complex measures may simply reflect simple properties of the signal (Martinerie et al., 2003; McSharry et al., 2003). Hence, in our first analysis, we focus on the relationship between monofractal measures and signal variance. For this, we used a simulated monofractal time series (termed fractional Brownian motion, or short fBm) with its standard deviation modulated by a modified ramp function.

The fBm was simulated with a Hurst exponent H = 0.7 and a modulating function M (described in more detail later and in Appendix A in Supplementary Material) and split into 1,800 1,024-sample epochs. We estimated the monofractal dimension of this simulated signal using the Higuchi and DFA methods. To assess the impact of signal variance, we have also tested the effect of epoch-based standardization. To ensure that our effects were not simply an artifact generated by the fBm, we also repeated the analysis on one exemplary icEEG data segment.

2.1.2. Experiment 2: Multifractal Estimation Stability

In order to evaluate the stability of multifractal properties in time, we generated a time series exhibiting stable multifractal properties over time using the p-Model. The time series was then evaluated using an epoch-based approach with the three estimators: MF-DFA, MF-DMA, and Chhabra-Jensen. The stability of the estimator can then simply be assessed as the temporal variability of its output.

2.1.3. Experiment 3: Multifractal Estimation of Human EEG and Its Potential Added Value

To assess whether the chosen multifractal metrics contribute any non-redundant features about the signal in addition to more established signal metrics, we analyzed human EEG signals recorded intracranially. Again, we used an epoch-based approach, and we compared the multifractal metrics to a number of conventional signal metrics (mean, standard deviation, line length, bandpower) on each epoch. The similarity between signal features was evaluated using Pearson correlation and Mutual Information (Guyon and Elisseeff, 2003) (the code is available at https://github.com/robince/gcmi) (Ince et al., 2017). Furthermore, monofractal metrics were also included in this comparison, to further demonstrate the advantages in applying a multifractal over monofractal approaches.

2.1.4. Experiment 4: Impact of Sampling Frequency and Epoch Length on Multifractal Estimation of Human EEG

Finally, we also evaluated the impact of the multifractal estimation parameters in the characterization of a seizure. We used intracranial EEG signals recorded from four patients undergoing pre-surgical planning, the signals were originally sampled at 5,000 Hz. For this analysis, down-sampled versions were evaluated with epochs of different sizes. To assess the effect of sampling frequency, we down-sampled the signal to 4,000, 3,000, 2,500, 2,000, 1,000, 800, 750, 600, 500, 400, 300, and 250 Hz. For each sampling frequency, we evaluated different epoch sizes (1,024, 2,048, 4,096, 8,192, and 16,384 points).

We defined a difference in multifractal spectrum width (Δα^†) during the seizure compared to the background as the effect size (Cohen's D) between the ictal and interictal periods:

$\begin{array}{cc}D=\frac{<\Delta {\alpha }_{ictal}^{†}>-<\Delta {\alpha }_{interictal}^{†}>}{s(\Delta {\alpha }_{interictal}^{†})}& (1)\end{array}$

where < Δα^† > represents the mean and s denotes standard deviation.

2.2. Fractal Dimension Estimation

To estimate the monofractal properties from a time series, we used two established estimation approaches: Higuchi method (Higuchi, 1988) and Detrended Fluctuation Analysis (Peng et al., 1994). These methods are widely applied in the literature and aim to capture the features of a time series in a single scaling exponent.

Mandelbrot (1982) defined fractals as self-similar structures with fractal dimensions (FD) that are between their topological and embedding dimensions T and E, and an established relationship of FD+H = E, where E = T+1, T = 1 in the case of a time series, and H is the Hurst exponent. When assuming this self-similarity, we can measure both FD and H in our EEG time series as alternative ways of estimating the fractal dimension.

However, we also note that more generally, the fractal dimension FD and the Hurst exponent H do not necessarily reflect the same property of the time series (see Gneiting and Schlather, 2004 for more details). Indeed, we empirically tested the relationship between our estimated FD and H for an example EEG time series and found that FD and H correlate with ρ = −0.8, and their empirical relationship is FD = −0.86H+2.74. For our application in EEG time series, we conclude that FD and H measure two related, but slightly different signal properties [FD: a measure of roughness, H: a measure of long memory dependency (Gneiting and Schlather, 2004)]. Note that this of course also depends on how the FD and H are estimated exactly. Nevertheless, for our paper, we will apply these two established methods to estimate FD and H, respectively. We will assess the properties of both methods in the context of EEG, to demonstrate that our conclusions generalize to both types of measures.

2.2.1. Higuchi Method

The Higuchi method measures the fractal dimension FD of a time series. It consists of constructing series with elements of an original time series and measuring their lengths (Higuchi, 1988). Given a time series with N time points X(1), X(2), …, X(N), the Equation (2) shows a rule for reconstructing smaller time series with elements of the original recording. The lengths of the time series can be assessed according to Equation (3). The brackets ⌊⌋ represent Gauss' notation, i.e., the rounded integer of the division (Higuchi, 1988). The variable d represents a down-sampling factor of the original time series.

$\begin{array}{c}\begin{array}{l}X(m),X(m+d),X(m+2d),\dots ,X\left(m+\left[\frac{N-m}{d}\right]d\right)\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}where}\end{array}\\ m=1,2,\dots ,d& (2)\end{array}$

$L_m(d)=\frac{\left\{{\displaystyle {\sum }_{i=1}^{[(N-m)/d]}\left|X(m+id)-X(m+(i-1)d)\right|\frac{N-1}{\lfloor (N-m)/d\rfloor d}}\right\}}{d}(3)$

If the average curve length < L_m(d) >_m over d sets follows a power law, according to Equation (4), the time series has scaling properties, with a fractal dimension FD_Hig.

$<L(d)>\propto d^{-F{D}_{Hig}}(4)$

The routine used in the estimation of Higuchi fractal dimension FD is available at https://uk.mathworks.com/matlabcentral/fileexchange/50290-higuchi-and-katz-fractal-dimension-measures.

2.2.2. Detrendred Fluctuation Analysis

The Detrended Fluctuation Analysis (DFA) method is an alternative method (Peng et al., 1994, 1995), which estimates the Hurst exponent H in time series data instead of the fractal dimension.

The method consists of the following steps: Initially the time series with N time points X(1), X(2), …, X(N) is integrated as follows:

$y(k)={\displaystyle \sum _{i=1}^k}(X(i)-\langle X\rangle )(5)$

Where X(i) represents the i−th element of the time series and < X > denotes the mean over the whole recording. The second step consists of dividing the time series into N_l windows of length l, then the mean square root of the integrated series is subtracted from the local trend, in every window (Peng et al., 1995), as shown in Equation (6).

$\begin{array}{cc}\begin{array}{c}\begin{array}{c}\begin{array}{c}F(l)=\sqrt{\frac{1}{N_l}{\displaystyle \sum _{k=1}^{N_l}}{[y(k)-y_l(k)]}^2}\end{array}\end{array}\end{array}& (6)\end{array}$

The local trend (y_l(k)) is obtained from a linear regression over the time series in the window, and number N_l represents the total number of windows. In the following step, Equation (6) is obtained for several window lengths (l). The relation between F(l) and l is described by a power law, according to Equation (7), where H is the Hurst exponent.

$F(l)\propto l^H(7)$

The code used here is available in the Physionet repository (https://www.physionet.org/physiotools/dfa/) (Peng et al., 1995; Goldberger et al., 2000).

2.3. Multifractal Spectrum Estimation

In this section, we describe three multifractal spectrum estimators: Multifractal Detrended Moving Average (Gu and Zhou, 2010), Multifractal Detrended Fluctuation Analysis (Kantelhardt et al., 2002; Ihlen and Vereijken, 2010; Ihlen, 2012), and Chhabra-Jensen (Chhabra and Jensen, 1989), as these are the most established methods used in the literature.

Multifractal properties are represented as spectra (Figure 1), where essentially the fractal scaling properties, or more precisely Hausdorff dimensions [often noted as f(α)], are measured over a range of different singularities (α). Formally, the singularity spectrum is a function that describes the Hausdorff dimension of subsets of the time series X(t) with a specific Hölder exponent, according to:

$\begin{array}{cc}f(\alpha )=D_F\{X(t_s),H(X(t_s))=\alpha \}& (8)\end{array}$

Essentially, f(α) is the Hausdorff dimension (D_F) of the subset (t_s) of the time series X(t_s) that has a the Hölder exponent α (van den Berg, 1999; Murcio et al., 2015). A definition of the Hausdorff dimension is available in Appendix D in Supplementary Material.

FIGURE 1

Figure 1. Multifractal singularity spectrum with a characteristic parabolic shape. The spectrum width (Δα) and height (Δf) measures are indicated by the arrows.

To characterize the function, or singularity spectrum f(α), usually, the width (Δα) and height (Δf)—differences of maximum and minimum values of α and f(α), respectively—of the spectrum are used. Δα indicates the range of singularities present in a time series, this is also the most commonly used measure of how multifractal a time series is. The spectrum height Δf indicates the range of Hausdorff dimensions present in the time series. See Figure 1 for an exemplary singularity spectrum plot.

2.3.1. MF-DMA

Multifractal Detrended Moving Average (MF-DMA) is one of the most commonly used methods for the estimation of multifractal measures. The method of calculation consists of the following steps (Gu and Zhou, 2010): Given time series X(t) with time points X(1), X(2), …, X(N), the cumulative sum time series is obtained:

$\begin{array}{cc}y(t)={\displaystyle \sum _{t=1}^N}X(t)& (9)\end{array}$

We then calculate the moving average over time windows of length l:

$\begin{array}{cc}\tilde{y}(t)=\frac{1}{l}{\displaystyle \sum _{z=0}^{l-1}}y(t-z)& (10)\end{array}$

A detrended version of the signal is obtained by the subtraction:

$\begin{array}{cc}ϵ(i)=y(i)-\tilde{y}(i)& (11)\end{array}$

The resulting series is then divided in N_l disjoint sets of points of size l and a root-mean-square function is obtained for each set ν via:

$F_{\nu }(l)={\left\{\frac{1}{l}{\displaystyle \sum _{i=1}^l}{ϵ}_{\nu }^2(i)\right\}}^{\frac{1}{2}}(12)$

A generalized qth-order overall fluctuation function can be obtained from:

$\begin{array}{cc}{F}_q(l)={\left\{\frac{1}{N_l}{\displaystyle \sum _{\nu =1}^{N_l}}{F}_{\nu }{(l)}^q\right\}}^{\frac{1}{q}}q\ne 0& (13)\end{array}$

and

$\begin{array}{cc}\ln F_0(l)=\frac{1}{N_l}{\displaystyle \sum _{\nu =1}^{N_l}}\ln F_{\nu }(l)\text{ for }q=0& (14)\end{array}$

It is possible to find a power-law relationship between F^q(l) and the window length, or scale l by:

$F^q(l)\propto l^{\alpha (q)}(15)$

The multifractal “mass exponent” (Biswas and Cresswell, 2012) can be defined as:

$\begin{array}{cc}\tau (q)=q\alpha (q)-D_f& (16)\end{array}$

where D_f is the fractal dimension of the support measure. For a single-channel time series, D_f = 1. The spectrum, f(α), can be obtained with a Legendre transform (Gu and Zhou, 2010):

$\begin{array}{cc}\alpha (q)=\frac{d\tau (q)}{dq}& (17)\end{array}$

$\begin{array}{cc}f(q)=q\alpha -\tau (q)& (18)\end{array}$

It is important to note that the Legendre transform is known to cause problems in multifractal spectra derivations if some heterogeneities are present in the signal, as has been reported elsewhere (Chhabra and Jensen, 1989; Mukli et al., 2015).

2.3.2. MF-DFA

The Multifractal Detrended Fluctuation Analysis (MF-DFA) method is essentially a generalization of the DFA approach (Kantelhardt et al., 2002; Ihlen, 2012). The time series is first rebuilt according to Equation (5).

It is then divided into $N_l=\frac{N}{l}$ non-overlapping epochs ν of length l. The variance of the detrended series is calculated as follows:

$F_{\nu }^2(l)=\frac{1}{l}{\displaystyle \sum _{k=1}^n}{(y((\nu -1)l+1)-y_v(k))}^2(19)$

where y_ν represents the fitting in the epoch ν obtained via linear regression. The overall q-th order fluctuation functions can be obtained as:

$F_q(l)={\left\{\frac{1}{N_l}{\displaystyle \sum _{\nu =1}^{N_l}}{(F_{\nu }^2(l))}^{\frac{q}{2}}\right\}}^{\frac{1}{q}}(20)$

A log-log plot of F_q(l) vs. l for different values of q should present a linear curve defined by the power law in Equation (15). Similarly to the MF-DMA method, the multifractal scaling exponent can be defined as in Equation (16) and the spectrum f(α) can be determined in the same way as in the MF-DMA approach.

2.3.3. Chhabra-Jensen

Multifractal spectra can be obtained in a more direct way, without the need for the Legendre transform using the Chhabra-Jensen (CJ) method (Chhabra and Jensen, 1989; Miranda et al., 2006; Zeleke and Si, 2006; Vázquez et al., 2008; Paz-Ferreiro et al., 2010a,b; Murcio et al., 2015; Xu et al., 2017; França et al., 2019). Considering a time series as a distribution over time, the approach consists of calculating a family of generalized measures by covering the time series with windows. These are probabilistic measures with an emphasis factor q that accentuates different singularities depending on its value. More singular regions are emphasized by q>1 whereas less singular regions will have a higher weight with q < 1 (Chhabra and Jensen, 1989).

First, we define:

$\begin{array}{cc}{\mu }_i(q,l)=\frac{P_i{(l)}^q}{{\displaystyle {\sum }_j{P}_j{(l)}^q}}& (21)\end{array}$

where P_i(l) represents the cumulative probability of a window i. l corresponds to the size of the window in which the generalized measures are obtained. The window epochs are indexed by the variables i and j. Then the multifractal spectra can be obtained directly from:

$\begin{array}{cc}\alpha (q)=\underset{l\to 0}{\lim }\frac{{\displaystyle {\sum }_i{\mu }_i(q,l)\log P_i(l)}}{\log l}& (22)\end{array}$

and

$\begin{array}{cc}f(q)=\underset{l\to 0}{\lim }\frac{{\displaystyle {\sum }_i{\mu }_i(q,l)\log {\mu }_i(q,l)}}{\log l}& (23)\end{array}$

A numerical approximation to the equations above is provided by the measures Mα and Mf functions in Equations (24) and (25).

$\begin{array}{cc}M\alpha ={\displaystyle \sum _i}{\mu }_i(q,l)\log P_i(l)& (24)\end{array}$

$\begin{array}{cc}Mf={\displaystyle \sum _i}{\mu }_i(q,l)\log {\mu }_i(q,l)& (25)\end{array}$

α and f(q) can then be obtained as the slopes by regressing these two measures against the scales l: Mα~l and Mf~l.

The algorithmic summary of the Chhabra-Jensen method consists of the following steps:

• The algorithm has as input the time series, a range of q values to which the spectrum will be evaluated, and window sizes l that vary in a dyadic scale.

• The time series is divided into non-overlapping epochs of length l and the generalized measures are estimated according to Equation (21).

• The measures Mα and Mf are obtained from the generalized measures.

• α and f(q) in Equations (22) and (23), respectively, are obtained with a linear regression procedure: log(Mα) is regressed against −log(l) and log(Mf) is regressed against −log(l), they give α and f respectively as the slopes.

• A rejection criterion is also used, where all q exponent values with R² < 0.9 in the linear regression are not considered.

The code used in this study to calculate the multifractal spectrum is available at: https://github.com/lucasfr/chhabra-jensen. A flow-chart diagram of the algorithm is included in the repository above and in Appendix Figure E1 in Supplementary Material.

2.4. Data

2.4.1. Simulating Fractal Time Series: Modulated Fractional Brownian Motion

To fully test methods of estimating the monofractal dimension from time series, we computationally produced time series that are known to be fractal (used for Experiment 1). We generated fractional Brownian motion (fBm) (Mandelbrot and Van Ness, 1968) profiles/time series using a novel modified version of the Wood-Chan or circulant embedding approach (Kroese and Botev, 2015; Shevchenko, 2015) that allow us to change the variance of the signal over time, in order to evaluate its influence on the fractal estimation. Our modulated fBM approach uses a modulating function, M(t), which produces a signal that has an amplitude varying over time. The details of fBm and our Modulated fBm (ModfBm) are described in Appendix A in Supplementary Material. The fBm time series was simulated with Hurst exponent H = 0.7; the value was chosen due to its persistent features, i.e., it generates a time series with memory. The modulating function M(t) used to modify the variance of the signal over time (see also described in Appendix A in Supplementary Material) is shown in Figure 2C. Using this method, we generated time series to evaluate the impact of variance change on monofractal estimators.

FIGURE 2

Figure 2. Impact of the signal standard deviation on monofractal scaling exponent estimation. (A) Modulation of the standard deviation of the time series over time; (B) Time series simulated using fractional Brownian motion based on modulation in (A); (C) Standard deviation of the simulated signal in (B). (D) Monofractal dimension obtained with the Higuchi method from signal without epoch-based standardization. (E) Monofractal dimension obtained with the Higuchi method from epoch-based standardized signal. (F) Hurst exponent obtained with the DFA method from signal without epoch-based standardization. (G) Hurst exponent obtained with the DFA method from epoch-based standardized signal.

Note that there are alternative methods to generate monofractal time series (Davies and Harte, 1987; Eke et al., 2002; Mukli et al., 2015; Nagy et al., 2017). However, as our aim was not to compare generative models of monofractal time series, but rather simply demonstrate that the effects we observe in EEG signals could be more general. We chose the above mentioned approaches as example demonstrations.

2.4.2. Simulating Multifractal Time Series: p-Model

Similarly to the fBm, we also used a computational procedure to generate time series that are known to be multifractal (for Experiment 2) based on the p-model, which was developed to reproduce features observed in turbulence experiments known to have multifractal properties (Meneveau and Sreenivasan, 1987). This is a simple model, having a single fraction p₁ as its only input and is often mentioned in literature (Meneveau and Sreenivasan, 1987, 1991; Lipa and Buschbeck, 1989; She and Leveque, 1994; Consolini et al., 1996; Davis et al., 1997; Sreenivasan and Antonia, 1997; Kestener and Arneodo, 2003; Zhou, 2008; Pechlivanidis and Arheimer, 2015). Briefly the algorithm works as follows: From an interval of length L and height ϵ_L = c (is a constant), we create two segments of length L/2. Based on the input parameter p₁, it is possible to establish a second fraction in which a second parameter will be given by p₂ = 1−p₁. The heights of each interval will thus be given by y = 2p₁ϵ_L, and y = 2p₂ϵ_L, respectively. This procedure is repeated for each remaining segment, selecting left or right for p₁ randomly (Meneveau and Sreenivasan, 1987).

We employed the p-model in the simulation of a time series profile with multifractal properties to be evaluated by different estimation methods. It was generated with a code available at http://www2.meteo.uni-bonn.de/staff/venema/themes/surrogates/pmodel/ (Davis et al., 1997; Venema et al., 2006). Using this algorithm, we generated time series to evaluate the performance of different multifractal estimators with p = 0.4. The value was rounded (for simplicity) from the figure used elsewhere (p = 0.375) (Davis et al., 1997).

2.4.3. Human EEG Data

Intracranial EEG data segments extracted from recordings in patients undergoing evaluation for epilepsy surgery were used for Experiments 3 and 4. In order to evaluate the effect of EEG signal variance change on multifractal properties (Experiment 3), we specifically looked for one recording, where the signal variance changes dramatically over time. One such recording was found in one patient (male, 28 years old, temporal lobe epilepsy, recorded at the National Hospital for Neurology and Neurosurgery (NHNN) (UCLH NHS Foundation Trust, Queen Square, London, UK), patient ID: “NHNN1”) near one seizure event. We used a 60-min recording segment around the epileptic seizure for our analysis. The seizure onset and offset were marked by expert clinicians, independent of this research project. Note that we used this segment specifically due to the dramatic change in signal variance, which actually occurred before the seizure and evolves over about 15 min. We do not make conclusions about the seizure event itself at this stage, but rather use this recording as an example to illustrate a technical point about multifractal property estimation from EEG.

To analyse the possible changes in multifractal properties during seizures (Experiment 4), we used a different dataset: Intracranial EEG from four subjects were retrieved from the ieeg.org repository (http://www.ieeg.org/) (Wagenaar et al., 2013): “I001_P005_D01,” “I001_P034_D01,” “I001_P010_D01,” and “Study 040.” These subjects were chosen due to the high sampling rate of their recordings (5 kHz), as we evaluated the impact of sampling frequency on multifractal properties. We extracted a 15-min segment around every seizure in each patient for further analysis. In Experiment 4, we performed the multifractal analysis on channels that were marked as seizure onset channels. We show the results for one patient in the main figure and the results for the remaining three patients are shown in Appendix B in Supplementary Material. Further information on the recordings is available in Appendix F in Supplementary Material.

The anonymized data analyzed in this study were recorded in patients undergoing evaluation for epilepsy surgery. iEEG.org portal provided EEG data and ethical approval for analyzing the data was provided by Mayo Clinic IRB (Brinkmann et al., 2009, 2016).

For NHNN data, the subject gave informed written consent, and the study was approved by the Joint Research Ethics Committee of the NHNN (UCLH NHS Foundation Trust) and UCL Queen Square Institute of Neurology, Queen Square, London, UK.

2.5. Pre-processing and Analysis of Time Series

Unless stated otherwise, we have applied the same pre-processing and analysis parameters to the computationally generated time series and the human EEG recordings and performed the fractal and multifractal estimations on 1,024-sample epochs. In Experiments 1 and 2, we were specifically interested in the effect of signal variance on the (multi) fractal estimation, and for comparison we also subjected the EEG signal to a standardization procedure, as follows:

$\begin{array}{cc}{x}^{\prime }=\frac{X-<X>}{s}& (26)\end{array}$

where < X> is the epoch mean and s the epoch standard deviation of the time series X, resulting in a time series with zero mean and unit standard deviation in each epoch.

The Chhabra-Jensen method requires as input a distribution function over the domain of positive real numbers, which is incompatible with EEG data which contain positive and negative values. Hence, we propose the use of a sigmoid-transformation here (Equation 27) to map the time series onto positive values, in order to apply the Chhabra-Jensen method. Example sigmoid functions and correspondingly transformed EEG signal are shown in Appendix Figure E2 in Supplementary Material.

$\begin{array}{cc}\sigma (X)=\frac{1}{1+e^{vX}}& (27)\end{array}$

The parameter v was chosen based on its effect on the estimated multifractal width for three types of time series: icEEG (NHNN1-channel 1), surrogate EEG (temporally shuffled values of the original time series from NHNN1-channel 1) and a simulated random series (with the same mean and variance), across the range v = [0.1, 2.0] in steps of 0.1. To find the optimal value for the parameter v, we needed to balance the trade-off between the three series in terms of presenting the most distinct Δα values (Appendix Figure E3A in Supplementary Material), while showing minimum distortion on the recording, or maximum correlation with the original time series (Appendix Figure E3B in Supplementary Material). We chose v = 1 as an acceptable trade-off point. Finally, to compare multifractal properties to classical EEG frequency band power, we used the following definitions for the classical EEG frequency bands: δ (0.5–4 Hz), θ (4–8 Hz), α (8–15 Hz), β (15–30 Hz), and γ (30–60 Hz).

3. Results

3.1. Experiment 1: Monofractal Estimation With Respect to Changing Signal Variance

We evaluated the relationship between monofractal measures and signal variance using a simulated time series based on fractional Brownian motion (fBm), where its signal variance is modulated by a modified ramp function. The modulation function is shown in Figure 2A and resulting time series in Figure 2B. The standard deviation of the generated time series indeed tracks the shape of the modulating function (Figure 2C).

We estimated the monofractal dimension of this simulated signal using two standard methods: Higuchi and DFA. We observe that both methods appear to be affected by the changing signal variance (Figures 2D,F). Furthermore, the effect persists even after epoch-based standardization (Figures 2E,G): the monofractal properties and standard deviation correlate with ρ = 1.00 and ρ = 0.99 for the Higuchi and DFA methods, respectively. A similar effect was observed for a real icEEG recording that contained changes in signal variance over time (Appendix Figure E4 in Supplementary Material).

In conclusion, monofractal properties derived for each epoch from DFA and Higuchi methods (with, or without signal standardization) correlate highly with the signal standard deviation of the epoch. Therefore, in epoch-based approaches (e.g., for application such as detecting or predicting epileptic seizures), the monofractal properties cannot be regarded as a new useful EEG feature of an epoch that is not redundant to standard deviation of the epoch. Thus we turn our attention to multifractal properties of the signal next.

3.2. Experiment 2: Multifractal Estimation Stability

In the following, we will denote the epoch-wise estimates of multifractal width Δα and height Δf (and Δα^† and Δf^† for the measure of the epoch-based standardized time series).

This experiment was designed to assess the reliability of the different multifractal estimation methods over time. In other words, if the multifractal properties of the time series remain constant over different epochs, then we expect the multifractal estimation method to show the same output over these different epochs. Note that the accuracy of these methods (i.e., the method outputting the expected multifractal measures of a predefined multifractal object with known multifractal properties) has been demonstrated elsewhere (Chhabra and Jensen, 1989; Kantelhardt et al., 2002; Gu and Zhou, 2010).

Figure 3 shows the simulated signal by the p-model and the outputs of the three multifractal spectral estimation methods. In all cases, the magnitude of Δα^† and Δf^† were clearly different from zero. The (Δα^†, Δf^†) output variances over time for the MF-DFA, MF-DMA, and Chhabra-Jensen estimation methods were: (0.018, 0.18), (4.17e-4, 0.0028), and (2.3e-30, 6.5e-30), respectively. In addition, the MF-DFA output violated the theoretical topological limit of Δf^† = 1, again indicating problems in the MF-DFA method, potentially due to the inversion of multifractal spectrum (Mukli et al., 2015). As the Chhabra-Jensen method shows the lowest variance over time (i.e., most reliable/stable), it will be our multifractal analysis method of choice for the remainder of this work.

FIGURE 3

Figure 3. Comparison of three multifractal spectrum estimation methods (MF-DFA, MF-DMA, and Chhabra-Jensen) for p-Model simulated time series. (A) Time series simulated for p = 0.4. (B) Estimated multifractal spectra width Δα^† and (C) height Δf^†.

3.3. Experiment 3: Multifractal Estimation of Human EEG and Its Potential Added Value

Next, we evaluated the relationship between multifractal signal properties and other widely used conventional EEG measures (such as signal variance). Figure 4 shows the results of the multifractal spectrum and conventional measures in comparison. The pattern of multifractal spectrum width without epoch-based standardization (Δα) reflects the signal variance closely, in contrast to the estimate for the epoch-based standardized signal (Δα^†). Finally, signal line length also shows a very different temporal profile from Δα^†. A similar figure showing the variation of Δf and Δf^† metrics is available in Appendix Figure E5 in Supplementary Material.

FIGURE 4

Figure 4. Temporal dynamics of multifractal spectrum width compared with conventional measures for human intracranial EEG. (A) Intracranial EEG segment containing a seizure (onset and offset marked with red vertical lines). Note that this recording was chosen because it showed a dramatic change in signal variance during non-seizure periods, not because of any seizure related properties. (B) Variation of multifractal spectrum width without epoch-wise standardization (Δα). (C) Multifractal spectrum width based on epoch-wise standardized time series (Δα^†). (D) Standard deviation in each epoch. (E) Line length in each epoch. Black line: moving average of each measure.

Figure 5 shows the quantification of similarities of the signals in Figure 4 through a correlation analysis. In summary, a high degree of correlation is present between the signal standard deviation, multifractal spectrum width (Δα), and detrended fluctuation analysis (monofractal approach) both with and without epoch-based standardization. We found that standardization reduces the correlation between Δα and the standard variation from ρ = 0.86 (for Δα ) to ρ = −0.14 (for Δα^†). We also note that Δα is highly correlated with DFA and DFA^† estimates (ρ = 0.74 and ρ = 0.71, respectively) while it is markedly reduced for Δα^† (|ρ| < 0.3). The analysis based on the mutual information (Ince et al., 2017) rather than correlation showed a similar pattern (Appendix Figure E6 in Supplementary Material).

FIGURE 5

Figure 5. Correlation between multifractal spectrum and conventional EEG measures for human icEEG data (from Figure 4). The diagonal of the matrix shows the distribution for each measure across epochs. The lower triangle contains the scatter plots for each pair of measures across epochs. The upper triangle shows the Pearson correlation value for each pair of measure, where the size of the font additionally corresponds to the correlation coefficient to provide an additional visual cue.

The relationships of the multifractal properties and specific EEG frequency band power are shown in Figure 6. In summary, the correlation values between the multifractal measures Δα^†, Δf^†, and signal power in the classical EEG bands are low (|ρ| < 0.3). A supplementary analysis of EEG time series data containing different sleep stages (which are known to be dominated by specific frequencies) shows similar results (see Appendix C in Supplementary Material). Based on these results, we focused on Δα^† (using epoch-wise standardization of the time series) in the subsequent analysis.

FIGURE 6

Figure 6. Comparison of multifractal measures with classical spectral band power. Scatter plot matrix comparing both standardized multifractal spectrum width and height (Δα^† and Δf^†) with the δ, θ, α, β, and γ average band power in each epoch. Each scatter point is derived from a single epoch of the time series. The diagonal of the matrix features the histograms for each measure. The lower triangle contains the scatter plots for each pair of measures. The upper triangle shows the Pearson correlation for each pair of measure, where the size of the font additionally corresponds to the correlation coefficient to provide an additional visual cue. The icEEG data underlying this figure is shown in Figure 4A.

3.4. Experiment 4: Impact of Sampling Frequency and Epoch Length on Multifractal Estimation of Human EEG

The variation of the multifractal spectrum width Δα^† for different combinations of epoch sizes and sampling frequencies is shown in Figure 7A. On visual inspection, it is clear that there are some combinations of epoch size and sampling frequency that show a clear increase of Δα^† during the ictal period (marked by the red lines). To quantify this effect, Figure 7B shows the Cohen's effect size D of the ictal vs. interictal Δα^† distributions plotted against epoch duration (in seconds). In this plot, we included 15 different sampling frequencies, and also data from three different EEG channels (all in the seizure onset zone). A peak in D can be seen at about 1 s (across all sampling frequencies), indicating that the change in Δα^† during a seizure can be best captured when using 1 s epochs (regardless of sampling frequency). This effect was not found for the sampling frequency or epoch length separately. Similar results for additional patients are shown on Supplementary materials (Appendix B in Supplementary Material).

FIGURE 7

Figure 7. Influence of EEG sampling frequency and epoch length on multifractal spectrum width around and during an epileptic seizure. (A) Multifractal spectrum width (Δα^†) in a 15-min intracranial EEG segment containing one seizure (onset and offset marked by the red lines). The signal was initially sampled at 5,000 Hz. Each column shows Δα^† for 5,000, 2,500, 500, and 100 Hz sampling rates. Different epoch sizes were used ranging from 1,024 to 16,384 samples (in each row). (B) Relationship of effect size D (between the interictal and ictal distribution of Δα^†) and epoch duration in seconds (obtained by dividing the number of sampling points by the sampling rate of the signal). Channel 1 is the data shown in (A). The solid line represents a LOESS curve fitting of the data points, with formula “y ~ x.” The data used for this figure is obtained from for subject “I001_P005_D01” around seizure 1. Channel 1: ADMacro_01. Channel 2: ADMacro_02. Channel 3: ADMacro_03.

4. Discussion

In this study, we have explored the monofractal and multifractal properties of human EEG recordings and used simulated data to test the performance of fractal property estimation methods. Although mono- and multi-fractal approaches have been widely employed in the study of physiological signals in humans (Ivanov et al., 1999; Stanley et al., 1999; Hu et al., 2004, 2009; Costa et al., 2017; França et al., 2019), we have demonstrated that the monofractal dimension may be capturing a similar signal feature as the signal variance. When using standardization to remove the effect of signal variance, we demonstrated that multifractal measures (estimated by the Chhabra-Jensen method) capture information not contained in widely used conventional signal measures, making it a viable feature for machine learning in clinical EEG applications. Finally, using epileptic seizure as an example, we showed that the epoch length can significantly impact the detection of time-varying effects in multifractal properties, suggesting the need for data- and application-specific optimization.

4.1. Methodological Considerations

One of our key observations is that monofractal estimators are tightly correlated with signal variance—even following epoch-wise standardization, whereas multifractal properties following epoch-wise standardization are no longer tightly correlated with signal variance. This may appear to be a curious and non-intuitive observation that, to our knowledge, has not been reported before.

To interpret this observation, it is worth noting the relationship between monofractal and multifractal analyses. Essentially, in multifractal analysis, at the point for which q = 2, the corresponding f(α) is the so-called correlation dimension, which is an alternative way of estimating the monofractal dimension (Murcio et al., 2015). The relationship between monofractal dimension and signal variance has been established and explained before (Cannon et al., 1997). By the same token, signal variance also affects higher statistical moments (q>2 or q < −2). However, when analysing the exact effect of variance on the multifractal spectrum (Appendix Figure E9 in Supplementary Material), we observe that the variance particularly impacts the multifractal spectrum width and height, but maintains an almost constant value of f(α) for q = 2. This explains why epoch-wise standardization does not impact monofractal dimension but does impact multifractal spectrum width and height. The mono- and multifractal properties we are investigating here are essentially describing different properties of the multifractal spectrum. Note that through our standardization procedure, we do not abolish “multifractality,” but only its dependence on signal variance. Future work has to show mathematically the exact reason for this observation, although intuitively it is understandable that the standardization procedure (a linear transformation of the signal) changes the q = 2 moment least and affects higher moment more.

We further observed that the Chhabra-Jensen method is the most reliable out of the three multifractal estimation methods. As was pointed out in the original publication (Chhabra and Jensen, 1989), this is most likely due to the fact that the Chhabra-Jensen method avoids a Legendre transform that the other methods require. The Legendre transformation requires smoothing of the D_q curve and can lead to errors. For further advantages of the Chhabra-Jensen method, the reader is referred to the original publication (Chhabra and Jensen, 1989). A recent development, FMF method (Mukli et al., 2015; Nagy et al., 2017), may be an alternative to the approach proposed in this study.

Finally, our analysis highlighted the importance of choosing an adequate epoch size given a sampling frequency, in order to study events such as epileptic seizures. However, our study was based on the analysis of ictal vs. interictal epochs, i.e., a hard separation that may not represent continuous phenomena accurately. Future work should take into account that multifractal properties may be continuously changing over time (a striking example is shown in Appendix Figure E7 in Supplementary Material), and an explicitly time based approach may be needed. Along similar lines, our finding of a optimal time scale may be due to the non-stationary nature of the multifractal properties. Further theoretical work may have to develop a temporally resolved multifractal estimator, in order to fully understand this aspect.

4.2. Implications for the Understanding of Brain Activity and Brain Generators

Previous studies reported that the brain is characterized by critical dynamics (Eguíluz et al., 2005; Chialvo, 2010, 2012; Racz et al., 2018). This characteristic, found from microscopic spatial scales (such as neuronal networks) (Beggs and Plenz, 2003, 2004) to the whole-brain level (Eguíluz et al., 2005), is thought to facilitate the storage and processing of information. It has been further suggested that more than one scaling exponent would be necessary to properly characterize the brain's critical dynamics (Suckling et al., 2008; Ihlen and Vereijken, 2010; Ciuciu, 2012; Fraiman and Chialvo, 2012; Zorick and Mandelkern, 2013; Papo, 2014; Zhang et al., 2015; Papo et al., 2017; Racz et al., 2018), as departures from the power-law pattern have been frequently observed in brain signals. Hence, it has been proposed that using additional, higher-order statistical moments can better characterize such data (Fraiman and Chialvo, 2012). In this work, we contribute a complementary observation: while monofractal measures of EEG appeared to essentially follow the slow changes of signal variance, multifractal characterization is capable of revealing new information.

In terms of generative processes that can produce monofractal properties, it has been suggested that a property called Self-Organized Criticality (SOC) (Bak et al., 1987) may play an essential role. SOC describes the capacity of a system to evolve naturally into a critical state (a state in which a minimum perturbation could lead to events of all sizes). Such phenomena display power-law distributions and fractal properties as signatures (Bak and Paczuski, 1995). An example process that displays SOC is the so-called single avalanche or Bak–Tang–Wiesenfeld model (also known as Abelian sandpile model) (Bak et al., 1987). SOC behavior has been linked to physiological control mechanisms, such as in human heart rate variability (Goldberger et al., 2002). Similar to SOC, a related regime—termed non-classical SOC—is thought to give rise to multifractal properties (Lovejoy and Schertzer, 2007). The analysis and understanding of the non-classical SOC is, however, still under development.

In this context, our multifractal spectral analyses of human EEG data suggest that cerebral phenomena should not be modeled by a single avalanche model (classical SOC), in agreement with findings in a previous study (Fraiman and Chialvo, 2012). Moreover, it is hypothesized that brain dynamics are non-ergodic (Bianco et al., 2007), i.e., display preferential states and depends on previous states (Papo, 2014), which are all properties of multifractal processes (Lovejoy and Schertzer, 2007). Thus, multifractal analyses could provide a new paradigm for studying brain function and structure, as previously suggested in other studies of normal (Suckling et al., 2008; Ihlen and Vereijken, 2010; Ciuciu, 2012; Zorick and Mandelkern, 2013; Papo, 2014; Papo et al., 2017; Racz et al., 2018) and pathological brain activity (Zhang et al., 2015). Furthermore, generative processes displaying multifractal properties could help understanding the observed multifractal changes on a mechanistic level.

4.3. On the Detection of Brain State Transitions in Health and Disease

We want to emphasize that the conclusions from our work are drawn on the basis that slow changes in signal fractal features can be captured by using an epoch-wise feature extraction procedure. It is also from a feature redundancy perspective that we argue for the need of multifractal approaches over monofractal measures. We do not dispute the usefulness of monofractal measures in other general applications. In our work, we essentially performed a feature selection procedure using correlation and mutual information (Guyon and Elisseeff, 2003). We evaluated how different signal feature compare on an epoch-wise basis. Feature selection is crucial to obtain faster and cost effective models, and avoids overfitting of the available data. It might also help achieving a deeper insight into the nature of the studied phenomena (Blum and Langley, 1997; Liu et al., 1998; Guyon and Elisseeff, 2003; Liu and Yu, 2005; Saeys et al., 2007).

A fundamental observation in our work is that an optimal time scales may exist for specific physiological processes (such as epileptic seizures) in terms of their multifractal dynamics (Figure 7 and Appendix B in Supplementary Material). This result suggests that, at least in an epoch-based study, for any given epileptic seizure in a given patient, the variety of scaling exponents (Δα) will depend on the length of the epoch analyzed. This is further supported by similar findings in monofractal analysis (Eke et al., 2002). The implications of this observation are that certain scaling exponents will only exist in specific time scales and the diversity of scaling exponents will depend on the duration of the epoch. These results suggest the potential need for “tuning,” i.e., potentially having to find the characteristic time for every studied phenomenon. If this is indeed the case, a temporally resolved (not epoch-based) multifractal method should be developed in future to adequately characterize brain dynamics.

Furthermore, the slow temporal changes in multifractal dynamics need to be characterized in a systematic way. Using epileptic seizures as an example, Appendix Figure E7 in Supplementary Material shows that dramatic changes in multifractal properties can sometimes be seen before an epileptic seizure. This observation requires further investigation to address questions such as: are all epileptic seizures characterized by pre-ictal changes in multifractal properties? Do other physiological processes, such as sleep, influence this finding? To answer these questions, we will most likely also need well-characterized experimental conditions, where seizures can be triggered in a controlled manner.

Finally, it is well-recognized that epileptic seizures are spatio-temporal processes (see e.g., Wang et al., 2014, 2017), and our current approach of only focusing in the temporal aspect in one location will need to be expanded. Data-driven unsupervised approaches, such as dimensionality reduction, may help summarize spatial aspects. Additionally, the challenge will be to develop a spatio-temporal multifractal analysis approach that can also deal with the challenges of low spatial sampling resolution in EEG recordings.

4.4. Outlook

Our work has highlighted several challenges that need to be considered when analysing multifractal properties of EEG signals; namely choice of the appropriate estimation method, estimation parameters, and the influence of the time series variance on signal features. We have suggested some solutions to these problems, such as the used of the Chhabra-Jensen approach combined with an epoch-wise standardization approach, which has shown potential capabilities as a signal feature for machine learning applications. We have also highlighted possible process-specific challenges. In terms of epileptic seizures, future work is required to analyse a larger number of patients in order to draw firmer conclusions on the potential clinical relevance of multifractal analyses. Furthermore, the study of mechanistic generative models of EEG may shed light on why those multifractal changes occur. For example, a generative process of potential interest could feature a modified version of Bak–Tang–Wiesenfeld model (Bak et al., 1987).

4.5. Summary

In this paper, we have analyzed the monofractal and multifractal properties of human EEG recordings. We have shown that monofractal estimates are influenced by the standard deviation of the time series, thus not capturing features beyond signal variance. For multifractal estimation, we have shown that the Chhabra-Jensen approach is the most stable, and we have developed a method of signal pre-processing to remove the influence caused by the variance of the signal. Using the suggested approach, the multifractal estimates do not correlate with traditional EEG measures, thus yielding additional information about the signal and being a relevant signal feature. Finally, our results also indicate a preferential time scale to identify differences in multifractal properties between ictal and interictal state recordings in patients with epilepsy.

Data Availability Statement

The datasets analyzed for this study can be found in https://github.com/yujiangwang/MultiFractalEEG.

Author Contributions

LF, LL, MW, and YW: conceptualization, project administration, resources, and writing of the original draft. LF: data curation, funding acquisition, investigation, and visualization. LF and YW: formal analysis and validation. JM, LF, LL, MW, and YW: methodology. JM, LF, and YW: software. LL, MW, and YW: supervision. JM, LF, LL, MW, ML, NS, and YW: writing of the review, and editing.

Funding

LF was supported by a grant from Brazilian National Council for Scientific and Technological Development (CNPq) (206907/2014-1). YW gratefully acknowledges the Wellcome Trust (208940/Z/17/Z).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors acknowledge the use of the UCL Legion High Performance Computing Facility (Legion@UCL), and associated support services, in the completion of this work.

The authors would like to thank Benjamin H. Brinkmann, Joost Wagenaar, and Hoameng Ung from iEEG.org, and Dr. Peter Taylor for helpful comments throughout the project.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphys.2018.01767/full#supplementary-material

References

Acharya, R., Faust, O., Kannathal, N., Chua, T., and Laxminarayan, S. (2005). Non-linear analysis of EEG signals at various sleep stages. Comput. Methods Prog. Biomed. 80, 37–45. doi: 10.1016/j.cmpb.2005.06.011

CrossRef Full Text | Google Scholar

Argoul, F., Arnéodo, A., Grasseau, G., Gagne, Y., Hopfinger, E. J., and Frisch, U. (1989). Wavelet analysis of turbulence reveals the multifractal nature of the Richardson cascade. Nature 338, 51–53. doi: 10.1038/338051a0

CrossRef Full Text | Google Scholar

Bak, P., and Paczuski, M. (1995). Complexity, contingency, and criticality. Proc. Natl. Acad. Sci. U.S.A. 92, 6689–6896. doi: 10.1073/pnas.92.15.6689

PubMed Abstract | CrossRef Full Text | Google Scholar

Bak, P., Tang, C., and Wiesenfeld, K. (1987). Self-organized criticality: an explanation of the 1/ f noise. Phys. Rev. Lett. 59, 381–384. doi: 10.1103/PhysRevLett.59.381

PubMed Abstract | CrossRef Full Text | Google Scholar

Baldassano, S. N., Brinkmann, B. H., Ung, H., Blevins, T., Conrad, E. C., Leyde, K., et al. (2017). Crowdsourcing seizure detection: algorithm development and validation on human implanted device recordings. Brain 140, 1680–1691. doi: 10.1093/brain/awx098

PubMed Abstract | CrossRef Full Text | Google Scholar

Bassett, D. S., Greenfield, D. L., Meyer-Lindenberg, A., Weinberger, D. R., Moore, S. W., and Bullmore, E. T. (2010). Efficient physical embedding of topologically complex information processing networks in brains and computer circuits. PLoS Comput. Biol. 6:e1000748. doi: 10.1371/journal.pcbi.1000748

PubMed Abstract | CrossRef Full Text | Google Scholar

Bassett, D. S., Meyer-Lindenberg, A., Achard, S., Duke, T., and Bullmore, E. (2006). Adaptive reconfiguration of fractal small-world human brain functional networks. Proc. Natl. Acad. Sci. U.S.A. 103, 19518–19523. doi: 10.1073/pnas.0606005103

PubMed Abstract | CrossRef Full Text | Google Scholar

Beggs, J. M., and Plenz, D. (2003). Neuronal avalanches in neocortical circuits. J. Neurosci. 23, 11167–11177. doi: 10.1523/JNEUROSCI.23-35-11167.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

Beggs, J. M., and Plenz, D. (2004). Neuronal avalanches are diverse and precise activity patterns that are stable for many hours in cortical slice cultures. J. Neurosci. 24, 5216–5229. doi: 10.1523/JNEUROSCI.0540-04.2004

PubMed Abstract | CrossRef Full Text | Google Scholar

Bengtsson, H. (2016). R.matlab: Read and Write Mat Files and Call Matlab From Within R. R package version 3.6.1.

Bianco, S., Ignaccolo, M., Rider, M. S., Ross, M. J., Winsor, P., and Grigolini, P. (2007). Brain, music, and non-Poisson renewal processes. Phys. Rev. E 75:061911. doi: 10.1103/PhysRevE.75.061911

PubMed Abstract | CrossRef Full Text | Google Scholar

Biswas, A. P. H., and Si, C. (2012). “Application of multifractal and joint multifractal analysis in examining soil spatial variation: a review,” in Fractal Analysis and Chaos in Geosciences, ed S.-A. Ouadfeul (InTech), 109–138. doi: 10.5772/51437

CrossRef Full Text | Google Scholar

Blum, A. L., and Langley, P. (1997). Selection of relevant features and examples in machine learning. Artif. Intell. 97, 245–271.

Google Scholar

Brinkmann, B. H., Bower, M. R., Stengel, K. A., Worrell, G. A., and Stead, M. (2009). Large-scale electrophysiology: acquisition, compression, encryption, and storage of big data. J. Neurosci. Methods 180, 185–192. doi: 10.1016/j.jneumeth.2009.03.022

PubMed Abstract | CrossRef Full Text | Google Scholar

Brinkmann, B. H., Wagenaar, J., Abbot, D., Adkins, P., Bosshard, S. C., Chen, M., et al. (2016). Crowdsourcing reproducible seizure forecasting in human and canine epilepsy. Brain 139, 1713–1722. doi: 10.1093/brain/aww045

PubMed Abstract | CrossRef Full Text | Google Scholar

Bullmore, E., Barnes, A., Bassett, D. S., Fornito, A., Kitzbichler, M., Meunier, D., et al. (2009). Generic aspects of complexity in brain imaging data and other biological systems. NeuroImage 47, 1125–1134. doi: 10.1016/j.neuroimage.2009.05.032

PubMed Abstract | CrossRef Full Text | Google Scholar

Bullmore, E., Fadili, J., Breakspear, M., Salvador, R., Suckling, J., and Brammer, M. (2003). Wavelets and statistical analysis of functional magnetic resonance images of the human brain. Stat. Methods Med. Res. 12, 375–399. doi: 10.1191/0962280203sm339ra

PubMed Abstract | CrossRef Full Text | Google Scholar

Cannon, M. J., Percival, D. B., Caccia, D. C., Raymond, G. M., and Bassingthwaighte, J. B. (1997). Evaluating scaled windowed variance methods for estimating the Hurst coefficient of time series. Physica A. 241, 606–626. doi: 10.1016/S0378-4371(97)00252-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Chhabra, A., and Jensen, R. V. (1989). Direct determination of the f(α) singularity spectrum. Phys. Rev. Lett. 62, 1327–1330. doi: 10.1103/PhysRevLett.62.1327

PubMed Abstract | CrossRef Full Text | Google Scholar

Chialvo, D. R. (2010). Emergent complex neural dynamics. Nat. Phys. 6, 744–750. doi: 10.1038/nphys1803

CrossRef Full Text | Google Scholar

Chialvo, D. R. (2012). Critical brain dynamics at large scale. Crit. Neural Syst. 2013, 1–22. doi: 10.1002/9783527651009.ch3

CrossRef Full Text | Google Scholar

Ciuciu, P. (2012). Scale-free and multifractal time dynamics of fMRI signals during rest and task. Front. Physiol. 3:186. doi: 10.3389/fphys.2012.00186

PubMed Abstract | CrossRef Full Text | Google Scholar

Consolini, G., Marcucci, M. F., and Candidi, M. (1996). Multifractal structure of auroral electrojet index data. Phys. Rev. Lett. 76, 4082–4085. doi: 10.1103/PhysRevLett.76.4082

PubMed Abstract | CrossRef Full Text | Google Scholar

Costa, I. d. S., Gamundí, A., Miranda, J. G. V., França, L. G. S., De Santana, C. N., and Montoya, P. (2017). Altered functional performance in patients with fibromyalgia. Front. Hum. Neurosci. 11:14. doi: 10.3389/fnhum.2017.00014

PubMed Abstract | CrossRef Full Text | Google Scholar

Davies, R. B., and Harte, D. S. (1987). Tests for Hurst effect. Biometrika 74, 95–101. doi: 10.1093/biomet/74.1.95

CrossRef Full Text | Google Scholar

Davis, A., Marshak, A., Cahalan, R., and Wiscombe, W. (1997). The landsat scale break in stratocumulus as a three-dimensional radiative transfer effect: implications for cloud remote sensing. J. Atmos. Sci. 54, 241–260. doi: 10.1175/1520-0469(1997)054<0241:TLSBIS>2.0.CO;2

CrossRef Full Text | Google Scholar

Eguíluz, V. M., Chialvo, D. R., Cecchi, G. A., Baliki, M., and Apkarian, A. V. (2005). Scale-free brain functional networks. Phys. Rev. Lett. 94:018102. doi: 10.1103/PhysRevLett.94.018102

PubMed Abstract | CrossRef Full Text | Google Scholar

Eke, A., Hermán, P., Bassingthwaighte, J. B., Raymond, G., Percival, D. B., Cannon, M., et al. (2000). Physiological time series: distinguishing fractal noises from motions. Pflügers Archiv Eur. J. Physiol. 439, 403–415. doi: 10.1007/s004249900135

PubMed Abstract | CrossRef Full Text | Google Scholar

Eke, A., Herman, P., Kocsis, L., and Kozak, L. R. (2002). Fractal characterization of complexity in temporal physiological signals. Physiol. Meas. 23, R1–R38. doi: 10.1088/0967-3334/23/1/201

PubMed Abstract | CrossRef Full Text | Google Scholar

Esteller, R., Vachtsevanos, G., Echauz, J., Henry, T., Pennell, P., Epstein, C., et al. (1999). “Fractal dimension characterizes seizure onset in epileptic patients,” in 1999 IEEE International Conference on Acoustics, Speech, and Signal Processing. Proceedings. ICASSP99 (Cat. No.99CH36258), Vol. 4, 2343–2346.

Google Scholar

Falconer, K. (2003). Fractal Geometry. Chichester: John Wiley & Sons, Ltd.

Google Scholar

Feder, J. (1988). Fractals. Boston, MA: Springer US.

Google Scholar

Fraiman, D., and Chialvo, D. R. (2012). What kind of noise is brain noise: anomalous scaling behavior of the resting brain activity fluctuations. Front. Physiol. 3:307. doi: 10.3389/fphys.2012.00307

PubMed Abstract | CrossRef Full Text | Google Scholar

França, L. G. S., Montoya, P., and Miranda, J. G. V. (2019). On multifractals: a non-linear study of actigraphy data. Phys. A 514, 612–619. doi: 10.1016/j.physa.2018.09.122

CrossRef Full Text | Google Scholar

Freestone, D. R., Karoly, P. J., Peterson, A. D. H., Kuhlmann, L., Lai, A., Goodarzy, F., et al. (2015). Seizure prediction: science fiction or soon to become reality? Curr. Neurol. Neurosci. Rep. 15:73. doi: 10.1007/s11910-015-0596-3

CrossRef Full Text | Google Scholar

Gneiting, T., and Schlather, M. (2004). Stochastic models that separate fractal dimension and the hurst effect. SIAM Rev. 46, 269–282. doi: 10.1137/S0036144501394387

CrossRef Full Text | Google Scholar

Goldberger, A. L., Amaral, L. A. N., Glass, L., Hausdorff, J. M., Ivanov, P. C., Mark, R. G., et al. (2000). PhysioBank, physioToolkit, and physioNet : components of a new research resource for complex physiologic signals. Circulation 101, e215–e220. doi: 10.1161/01.CIR.101.23.e215

PubMed Abstract | CrossRef Full Text | Google Scholar

Goldberger, A. L., Amaral, L. A. N., Hausdorff, J. M., Ivanov, P. C., Peng, C., and Stanley, H. E. (2002). Fractal dynamics in physiology: alterations with disease and aging. Proc. Natl. Acad. Sci. U.S.A. 99(Suppl. 1), 2466–2472. doi: 10.1073/pnas.012579499

PubMed Abstract | CrossRef Full Text | Google Scholar

Gómez, C., Mediavilla, Á., Hornero, R., Abásolo, D., and Fernández, A. (2009). Use of the Higuchi's fractal dimension for the analysis of MEG recordings from Alzheimer's disease patients. Med. Eng. Phys. 31, 306–313. doi: 10.1016/j.medengphy.2008.06.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Gong, P., Nikolaev, A. R., and van Leeuwen, C. (2003). Scale-invariant fluctuations of the dynamical synchronization in human brain electrical activity. Neurosci. Lett. 336, 33–36. doi: 10.1016/S0304-3940(02)01247-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu, G.-F., and Zhou, W.-X. (2010). Detrending moving average algorithm for multifractals. Phys. Rev. E 82:011136. doi: 10.1103/PhysRevE.82.011136

PubMed Abstract | CrossRef Full Text | Google Scholar

Guyon, I., and Elisseeff, A. (2003). An introduction to variable and feature selection. J. Mach. Learn. Res. 3, 1157–1182. Available online at: http://www.jmlr.org/papers/v3/guyon03a.html

Google Scholar

Higuchi, T. (1988). Approach to an irregular time series on the basis of the fractal theory. Phys. D 31, 277–283. Guyon and Elisseeff, 2003

Google Scholar

Hsu, W.-Y., Lin, C.-C., Ju, M.-S., and Sun, Y.-N. (2007). Wavelet-based fractal features with active segment selection: application to single-trial EEG data. J. Neurosci. Methods 163, 145–160. doi: 10.1016/j.jneumeth.2007.02.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Hu, K., Ivanov, P. C., Chen, Z., Hilton, M. F., Stanley, H., and Shea, S. A. (2004). Non-random fluctuations and multi-scale dynamics regulation of human activity. Phys. A Stat. Mech. Appl. 337, 307–318. doi: 10.1016/j.physa.2004.01.042

PubMed Abstract | CrossRef Full Text | Google Scholar

Hu, K., Van Someren, E. J. W., Shea, S. A., and Scheer, F. A. J. L. (2009). Reduction of scale invariance of activity fluctuations with aging and Alzheimer's disease: involvement of the circadian pacemaker. Proc. Natl. Acad. Sci. U.S.A. 106, 2490–2494. doi: 10.1073/pnas.0806087106

PubMed Abstract | CrossRef Full Text | Google Scholar

Ihlen, E. A. F. (2012). Introduction to multifractal detrended fluctuation analysis in Matlab. Front. Physiol. 3, 1–18. doi: 10.3389/fphys.2012.00141

PubMed Abstract | CrossRef Full Text | Google Scholar

Ihlen, E. A. F., and Vereijken, B. (2010). Interaction-dominant dynamics in human cognition: beyond 1/fα fluctuation. J. Exp. Psychol. 139, 436–463. doi: 10.1037/a0019098

CrossRef Full Text | Google Scholar

Ince, R. A., Giordano, B. L., Kayser, C., Rousselet, G. A., Gross, J., and Schyns, P. G. (2017). A statistical framework for neuroimaging data analysis based on mutual information estimated via a gaussian copula. Hum. Brain Mapp. 38, 1541–1573. doi: 10.1002/hbm.23471

PubMed Abstract | CrossRef Full Text | Google Scholar

Ivanov, P. C., Amaral, L. A. N., Goldberger, A. L., Havlin, S., Rosenblum, M. G., Struzik, Z. R., et al. (1999). Multifractality in human heartbeat dynamics. Nature 399, 461–465. doi: 10.1038/20924

PubMed Abstract | CrossRef Full Text | Google Scholar

Kantelhardt, J. W., Zschiegner, S. A., Koscielny-Bunde, E., Havlin, S., Bunde, A., and Stanley, H. (2002). Multifractal detrended fluctuation analysis of nonstationary time series. Phys. A 316, 87–114. doi: 10.1016/S0378-4371(02)01383-3

CrossRef Full Text | Google Scholar

Karoly, P. J., Ung, H., Grayden, D. B., Kuhlmann, L., Leyde, K., Cook, M. J., et al. (2017). The circadian profile of epilepsy improves seizure forecasting. Brain 140, 2169–2182. doi: 10.1093/brain/awx173

PubMed Abstract | CrossRef Full Text | Google Scholar

Kestener, P., and Arneodo, A. (2003). Three-dimensional wavelet-based multifractal method: the need for revisiting the multifractal description of turbulence dissipation data. Phys. Rev. Lett. 91:194501. doi: 10.1103/PhysRevLett.91.194501

PubMed Abstract | CrossRef Full Text | Google Scholar

Kroese, D. P., and Botev, Z. I. (2015). “Spatial process simulation,” Lecture Notes in Mathematics, ed V. Schmidt (Cham: Springer International Publishing), 369–404.

Google Scholar

Kuhlmann, L., Karoly, P., Freestone, D. R., Brinkmann, B. H., Temko, A., Barachant, A., et al. (2018a). Epilepsyecosystem.org: crowd-sourcing reproducible seizure prediction with long-term human intracranial EEG. Brain 141, 2619–2630. doi: 10.1093/brain/awy210

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuhlmann, L., Lehnertz, K., Richardson, M. P., Schelter, B., and Zaveri, H. P. (2018b). Seizure prediction-ready for a new era. Nat. Rev. Neurol. 14, 618–630. doi: 10.1038/s41582-018-0055-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X., Polygiannakis, J., Kapiris, P., Peratzakis, A., Eftaxias, K., and Yao, X. (2005). Fractal spectral analysis of pre-epileptic seizures in terms of criticality. J. Neural Eng. 2, 11–16. doi: 10.1088/1741-2560/2/2/002

PubMed Abstract | CrossRef Full Text | Google Scholar

Linkenkaer-Hansen, K., Nikouline, V. V., Palva, J. M., and Ilmoniemi, R. J. (2001). Long-range temporal correlations and scaling behavior in human brain oscillations. J. Neurosci. 21, 1370–1377. doi: 10.1523/JNEUROSCI.21-04-01370.2001

PubMed Abstract | CrossRef Full Text | Google Scholar

Lipa, P., and Buschbeck, B. (1989). From strong to weak intermittency. Phys. Lett. B 223, 465–469.

Google Scholar

Liu, H., Motoda, H., and Dash, M. (1998). “A monotonic measure for optimal feature selection,” in European Conference on Machine Learning, 101–106. doi: 10.1007/BFb0026678

CrossRef Full Text | Google Scholar

Liu, H., and Yu, L. (2005). Toward integrating feature selection algorithms for classification and clustering. IEEE Trans. Knowl. Data Eng. 17, 491–502. doi: 10.1109/TKDE.2005.66

CrossRef Full Text | Google Scholar

Lovejoy, S., and Schertzer, D. (2007). Scaling and multifractal fields in the solid earth and topography. Nonlin. Process. Geophys. 14, 465–502. doi: 10.5194/npg-14-465-2007

CrossRef Full Text | Google Scholar

Lutzenberger, W., Preissl, H., and Pulvermüller, F. (1995). Fractal dimension of electroencephalographic time series and underlying brain processes. Biol. Cybern. 73, 477–482. doi: 10.1007/BF00201482

PubMed Abstract | CrossRef Full Text | Google Scholar

Mandelbrot, B. B. (1982). The Fractal Geometry of Nature, Vol. 1. New York, NY: WH freeman.

Google Scholar

Mandelbrot, B. B., and Van Ness, J. W. (1968). Fractional brownian motions, fractional noises and applications. SIAM Rev. 10, 422–437. doi: 10.1137/1010093

CrossRef Full Text | Google Scholar

Martinerie, J., Le Van Quyen, M., Baulac, M., and Renault, B. (2003). Reply to "Prediction of epileptic seizures: are nonlinear methods relevant?” Nat. Med. 9, 242–242. doi: 10.1038/nm0303-242

CrossRef Full Text | Google Scholar

McSharry, P. E., Smith, L. A., and Tarassenko, L. (2003). Reply to "Prediction of epileptic seizures: are nonlinear methods relevant?” Nat. Med. 9, 241–242. doi: 10.1038/nm0303-2

CrossRef Full Text | Google Scholar

Meneveau, C., and Sreenivasan, K. R. (1987). Simple multifractal cascade model for fully developed turbulence. Phys. Rev. Lett. 59, 1424–1427. doi: 10.1103/PhysRevLett.59.1424

PubMed Abstract | CrossRef Full Text | Google Scholar

Meneveau, C., and Sreenivasan, K. R. (1991). The multifractal nature of turbulent energy dissipation. J. Fluid Mech. 224:429. doi: 10.1017/S0022112091001830

CrossRef Full Text | Google Scholar

Miranda, J., Montero, E., Alves, M., Paz González, A., and Vidal Vázquez, E. (2006). Multifractal characterization of saprolite particle-size distributions after topsoil removal. Geoderma 134, 373–385. doi: 10.1016/j.geoderma.2006.03.014

CrossRef Full Text | Google Scholar

Mormann, F., Andrzejak, R. G., Elger, C. E., and Lehnertz, K. (2007). Seizure prediction: the long and winding road. Brain 130, 314–333. doi: 10.1093/brain/awl241

PubMed Abstract | CrossRef Full Text | Google Scholar

Mukli, P., Nagy, Z., and Eke, A. (2015). Multifractal formalism by enforcing the universal behavior of scaling functions. Phys. A 417, 150–167. doi: 10.1016/j.physa.2014.09.002

CrossRef Full Text | Google Scholar

Murcio, R., Masucci, A. P., Arcaute, E., and Batty, M. (2015). Multifractal to monofractal evolution of the London street network. Phys. Rev. E 92:062130. doi: 10.1103/PhysRevE.92.062130

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagy, Z., Mukli, P., Herman, P., and Eke, A. (2017). Decomposing multifractal crossovers. Front. Physiol. 8:533. doi: 10.3389/fphys.2017.00533

PubMed Abstract | CrossRef Full Text | Google Scholar

Neuwirth, E. (2014). Rcolorbrewer: Colorbrewer Palettes. R package version 1.1-2.

Papo, D. (2014). Functional significance of complex fluctuations in brain activity: from resting state to cognitive neuroscience. Front. Syst. Neurosci. 8:112. doi: 10.3389/fnsys.2014.00112

PubMed Abstract | CrossRef Full Text | Google Scholar

Papo, D., Goñi, J., and Buldú, J. M. (2017). Editorial: on the relation of dynamics and structure in brain networks. Chaos Interdisc. J. Nonlin. Sci. 27:047201. doi: 10.1063/1.4981391

PubMed Abstract | CrossRef Full Text | Google Scholar

Paz-Ferreiro, J., Miranda, J., and Vidal Vázquez, E. (2010b). Multifractal analysis of soil porosity based on mercury injection and nitrogen adsorption. Vadose Zone J. 9:325. doi: 10.2136/vzj2009.0090

CrossRef Full Text | Google Scholar

Paz-Ferreiro, J., Vázquez, E. V., and Miranda, J. (2010a). Assessing soil particle-size distribution on experimental plots with similar texture under different management systems using multifractal parameters. Geoderma 160, 47–56. doi: 10.1016/j.geoderma.2010.02.002

CrossRef Full Text | Google Scholar

Pechlivanidis, I. G., and Arheimer, B. (2015). Large-scale hydrological modelling by using modified PUB recommendations: the India-HYPE case. Hydrol. Earth Syst. Sci. 19, 4559–4579. doi: 10.5194/hess-19-4559-2015

CrossRef Full Text | Google Scholar

Peng, C., Buldyrev, S. V., Havlin, S., Simons, M., Stanley, H. E., and Goldberger, A. L. (1994). Mosaic organization of DNA nucleotides. Phys. Rev. E 49, 1685–1689. doi: 10.1103/PhysRevE.49.1685

PubMed Abstract | CrossRef Full Text | Google Scholar

Peng, C., Havlin, S., Stanley, H. E., and Goldberger, A. L. (1995). Quantification of scaling exponents and crossover phenomena in nonstationary heartbeat time series. Chaos Interdisc. J. Nonlin. Sci. 5, 82–87. doi: 10.1063/1.166141

PubMed Abstract | CrossRef Full Text | Google Scholar

Pereda, E., Gamundi, A., Rial, R., and González, J. (1998). Non-linear behaviour of human EEG: fractal exponent versus correlation dimension in awake and sleep stages. Neurosci. Lett. 250, 91–94. doi: 10.1016/S0304-3940(98)00435-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Peterson, B. G., and Carl, P. (2014). Performanceanalytics: Econometric Tools for Performance and Risk Analysis. R package version 1.4.3541.

R Core Team (2017). R: A Language and Environment for Statistical Computing (Vienna).

Racz, F. S., Mukli, P., Nagy, Z., and Eke, A. (2018). Multifractal dynamics of resting-state functional connectivity in the prefrontal cortex. Physiol. Meas. 39:024003. doi: 10.1088/1361-6579/aaa916

PubMed Abstract | CrossRef Full Text | Google Scholar

Saeys, Y., Inza, I., and Larranaga, P. (2007). A review of feature selection techniques in bioinformatics. Bioinformatics 23, 2507–2517. doi: 10.1093/bioinformatics/btm344

PubMed Abstract | CrossRef Full Text | Google Scholar

She, Z.-S., and Leveque, E. (1994). Universal scaling laws in fully developed turbulence. Phys. Rev. Lett. 72, 336–339. doi: 10.1103/PhysRevLett.72.336

PubMed Abstract | CrossRef Full Text | Google Scholar

Shevchenko, G. (2015). “Fractional Brownian motion in a nutshell,” in International Journal of Modern Physics: Conference Series, Vol. 36. doi: 10.1142/S2010194515600022

CrossRef Full Text | Google Scholar

Sreenivasan, K. R., and Antonia, R. A. (1997). The phenomenology of small-scale turbulence. Annu. Rev. Fluid Mech. 29, 435–472. doi: 10.1146/annurev.fluid.29.1.435

CrossRef Full Text | Google Scholar

Stam, C. J., and de Bruin, E. A. (2004). Scale-free dynamics of global functional connectivity in the human brain. Hum. Brain Mapp. 22, 97–109. doi: 10.1002/hbm.20016

PubMed Abstract | CrossRef Full Text | Google Scholar

Stanley, H., Amaral, L., Goldberger, A., Havlin, S., Ivanov, P., and Peng, C. (1999). Statistical physics and physiology: monofractal and multifractal approaches. Phys. A Stat. Mech. Appl. 270, 309–324. doi: 10.1016/S0378-4371(99)00230-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Suckling, J., Wink, A. M., Bernard, F. A., Barnes, A., and Bullmore, E. (2008). Endogenous multifractal brain dynamics are modulated by age, cholinergic blockade and cognitive performance. J. Neurosci. Methods 174, 292–300. doi: 10.1016/j.jneumeth.2008.06.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Van De Ville, D., Britz, J., and Michel, C. M. (2010). EEG microstate sequences in healthy humans at rest reveal scale-free dynamics. Proc. Natl. Acad. Sci. U.S.A. 107, 18179–18184. doi: 10.1073/pnas.1007841107

PubMed Abstract | CrossRef Full Text | Google Scholar

van den Berg, J. C. (ed.). (1999). Wavelets in Physics. Cambridge: Cambridge University Press.

Google Scholar

Varatharajah, Y., Berry, B., Cimbalnik, J., Kremen, V., Van Gompel, J., Stead, M., et al. (2018). Integrating artificial intelligence with real-time intracranial EEG monitoring to automate interictal identification of seizure onset zones in focal epilepsy. J. Neural Eng. 15:046035. doi: 10.1088/1741-2552/aac960

PubMed Abstract | CrossRef Full Text | Google Scholar

Vázquez, E. V., Ferreiro, J. P., Miranda, J. G. V., and González, A. P. (2008). Multifractal analysis of pore size distributions as affected by simulated rainfall. Vadose Zone J. 7, 500–511. doi: 10.2136/vzj2007.0011

CrossRef Full Text | Google Scholar

Venema, V., Bachner, S., Rust, H. W., and Simmer, C. (2006). Statistical characteristics of surrogate data based on geophysical measurements. Nonlin. Process. Geophys. 13, 449–466. doi: 10.5194/npg-13-449-2006

CrossRef Full Text | Google Scholar

Wagenaar, J. B., Brinkmann, B. H., Ives, Z., Worrell, G. A., and Litt, B. (2013). “A multimodal platform for cloud-based collaborative research,” in 2013 6th International IEEE/EMBS Conference on Neural Engineering (NER) (IEEE), 1386–1389. doi: 10.1109/NER.2013.6696201

CrossRef Full Text | Google Scholar

Wang, Y., Goodfellow, M., Taylor, P. N., and Baier, G. (2014). Dynamic mechanisms of neocortical focal seizure onset. PLoS Comput. Biol. 10:e1003787. doi: 10.1371/journal.pcbi.1003787

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Y., Trevelyan, A. J., Valentin, A., Alarcon, G., Taylor, P. N., and Kaiser, M. (2017). Mechanisms underlying different onset patterns of focal seizures. PLoS Comput. Biol. 13:e1005475. doi: 10.1371/journal.pcbi.1005475

PubMed Abstract | CrossRef Full Text | Google Scholar

Werner, G. (2010). Fractals in the nervous system: conceptual implications for theoretical neuroscience. Front. Physiol. 1:15. doi: 10.3389/fphys.2010.00015

PubMed Abstract | CrossRef Full Text | Google Scholar

Wickham, H. (2007). Reshaping data with the reshape package. J. Stat. Softw. 21, 1–20. doi: 10.18637/jss.v021.i12

CrossRef Full Text | Google Scholar

Wickham, H. (2009). ggplot2: Elegant Graphics for Data Analysis. New York, NY: Springer-Verlag.

Google Scholar

Xu, H. C., Gu, G. F., and Zhou, W. X. (2017). Direct determination approach for the multifractal detrending moving average analysis. Phys. Rev. E 96:052201. doi: 10.1103/PhysRevE.96.052201

PubMed Abstract | CrossRef Full Text | Google Scholar

Xue, Y., and Bogdan, P. (2017). Reliable multi-fractal characterization of weighted complex networks: algorithms and implications. Sci. Rep. 7:7487. doi: 10.1038/s41598-017-07209-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Zappasodi, F., Olejarczyk, E., Marzetti, L., Assenza, G., Pizzella, V., and Tecchio, F. (2014). Fractal dimension of EEG activity senses neuronal impairment in acute stroke. PLoS ONE 9:e100199. doi: 10.1371/journal.pone.0100199

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeleke, T. B., and Si, B. C. (2006). Characterizing scale-dependent spatial relationships between soil properties using multifractal techniques. Geoderma 134, 440–452. doi: 10.1016/j.geoderma.2006.03.013

CrossRef Full Text | Google Scholar

Zhang, Y., Zhou, W., and Yuan, S. (2015). Multifractal analysis and relevance vector machine-based automatic seizure detection in intracranial EEG. Int. J. Neural Syst. 25:1550020. doi: 10.1142/S0129065715500203

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, W.-X. (2008). Multifractal detrended cross-correlation analysis for two nonstationary signals. Phys. Rev. E 77:066211. doi: 10.1103/PhysRevE.77.066211

PubMed Abstract | CrossRef Full Text | Google Scholar

Zorick, T., and Mandelkern, M. A. (2013). Multifractal detrended fluctuation analysis of human EEG: preliminary investigation and comparison with the wavelet transform modulus maxima technique. PLoS ONE 8:e68360. doi: 10.1371/journal.pone.0068360

PubMed Abstract | CrossRef Full Text | Google Scholar