Nonlinear and machine learning analyses on high-density EEG data of math experts and novices

Current trend in neurosciences is to use naturalistic stimuli, such as cinema, class-room biology or video gaming, aiming to understand the brain functions during ecologically valid conditions. Naturalistic stimuli recruit complex and overlapping cognitive, emotional and sensory brain processes. Brain oscillations form underlying mechanisms for such processes, and further, these processes can be modified by expertise. Human cortical functions are often analyzed with linear methods despite brain as a biological system is highly nonlinear. This study applies a relatively robust nonlinear method, Higuchi fractal dimension (HFD), to classify cortical functions of math experts and novices when they solve long and complex math demonstrations in an EEG laboratory. Brain imaging data, which is collected over a long time span during naturalistic stimuli, enables the application of data-driven analyses. Therefore, we also explore the neural signature of math expertise with machine learning algorithms. There is a need for novel methodologies in analyzing naturalistic data because formulation of theories of the brain functions in the real world based on reductionist and simplified study designs is both challenging and questionable. Data-driven intelligent approaches may be helpful in developing and testing new theories on complex brain functions. Our results clarify the different neural signature, analyzed by HFD, of math experts and novices during complex math and suggest machine learning as a promising data-driven approach to understand the brain processes in expertise and mathematical cognition.

Continuous brain imaging data, which is collected over a long time span during naturalistic stimuli, enables the application of data-driven analyses (Cantlon, 2020 [2]; Zhang et al., 2021 [4]).Machine learning (ML) analyses may assist in generating new hypotheses about the underlying task-relevant brain processes, especially in the naturalistic context.In such contexts, several low and high-level overlapping brain processes occur simultaneously (Nastase et al., 2020 [3]).Due to the overlapping nature of several brain processes, extension of the neuroscientific theories formulated based on reductionist and simplified study designs is both challenging and questionable (Cantlon, 2020 [2]).Novel methodologies in analyzing naturalistic data are required and data-driven intelligent approaches form a good candidate for developing and testing new theories on the brain functions in the real world (Nastase et al., 2020 [3]).
Recent developments in ML are already applied in healthcare and extend to several fields: spike detection in epilepsy, dementia prediction, and mental health and sleep stage classification (Singh et al., 2022 [10]).These data-driven methods aim to transform healthcare delivery and to change the trajectory of brain health by addressing brain care earlier in the lifespan (Singh et al., 2022 [10]).For example, recent advances utilizing ML, specifically techniques with Brain-Computer Interfaces (BCI), help stroke patients either restore neurologic pathways or communicate with an electronic prosthetic (Cervera et al., 2018 [11]; Baniqued et al., 2021 [12]).On the other hand, ML may help in diagnostics of conditions like stroke or Alzheimer through the detection of disease-specific EEG biomarkers (Karthik et al., 2020 [13]; Meghdadi et al., 2021 [14]).
In addition to the applications for prediction and diagnostics in healthcare, ML for brain imaging has application possibilities in the contexts of learning and education (Bavelier and Green, 2019 [7]; Cantlon, 2020 [2]).For decades, scientists have studied the brain processes during cognitive tasks, like mathematics or language.These studies have brought valuable knowledge on the domaingeneral brain functions of working memory, attention, and solving strategies (e.g.De Smedt et al., 2009 [15]; Kulasingham et al., 2021 [16]; Wang et al., 2020 [17]) and domain-specific brain functions on numeric and verbal processing (e.g.Amalric and Dehaene 2016 [18], 2019 [19]).Some studies have focused to understand healthy development and expertise (Jeon et al., 2019 [20]; Zhang et al., 2015 [21]), whereas others bring insights on disrupted development and learning deficits (Klados et al. 2017 [22]; Rubinsten, 2015 [23]).Neuroscientific studies made in learning sciences have not yet utilized ML in the data analysis.However, ML has potential to be used in data-driven hypothesis formation of the brain functions underlying expertise development or learning deficits, and for real-time adaptive feedback in learning and focused attention (Kefalis et al., 2020 [24]; Hunkin et al., 2021 [25]).
Previous studies show differences in the brain functions of math experts and novices during short and simple math tasks (e.g.Grabner et al., 2007 [26]).Such differences are associated with brain functions modified through expertise, such as rote learning and strategy selection for solving the tasks at hand (Grabner and De Smedt, 2012 [27]; Hinault and Lemaire, 2016 [28]).However, a few second simple math tasks, which are used traditionally as stimuli in studies on math expertise, seldomly create enough of continuous brain imaging data for which to successfully apply the ML methods.
Despite some ML algorithms are designed to evaluate raw EEG data (da Silva Lourenco et al., 2021 [29]), several studies which focus on the comparison of brain states have preprocessed the data before ML classifications.The brain, as many biological systems, behaves in a nonlinear manner.Nonlinear behavior of biological systems is characterized by a high degree of variability in the time domain (nonstationarity) and randomness that could be attributed to the interaction of internal and external factors influencing the organism (Glass, 2001 [30]; Eke et al., 2002 [31]).Engagement with complex math recruits several cognitive brain processes which overlap with sensory and emotional processes (Suarez-Pellicioni et al., 2016 [32]; Wang et al., 2015 [33]).Therefore, the EEG data collected during such cognitively challenging task is likely highly complex, and therefore, a potentially optimal way to process such data includes an analysis which is suitable for nonlinear systems.
Cognitively challenging tasks create a brain state which is clearly different from those of relaxed states (Finn, 2021 [34]).Fractal dimension is a highly sensitive measure in the detection of hidden information contained in physiological time series (Klonowski [35], 2002; Raghavendra and Dutt, 2010 [36]) and is shown to vary depending on the brain state.An often-used nonlinear measure for signal analysis is Higuchi's fractal dimension (HFD) which is a measure of signal complexity in the time domain (Higuchi, 1988 [37]; Spasic et al., 2008 [38]).Previous studies utilizing such methods classified successfully different sleep stages and detected the difference in the brain state during drowsiness and wakefulness (Inoye et al., 1994 [39]; Šušmáková and Krakovská, 2008 [40]).HFD showed the most robust results and seems to be superior to other FD methods for EEG signals (Solhjoo and Nasrabadi, 2005 [41]; Ahmadlou et al., 2012 [42]).
This study investigated the neural signature of math expertise with a relatively robust nonlinear analysis, HFD, and explored a new paradigm by applying ML to EEG data collected from math experts and novices when they engaged with long and complex math demonstrations.Such math demonstrations with a duration up to one minute form a part of the current trend in investigating the brain with naturalistic stimuli.Our aim was to describe the EEG data during advanced mathematical cognition with a nonlinear method and evaluate whether the neural signature of math experts and novices differ in a way which is detectable with artificial intelligence.We hypothesized that the experts' and novices' brain functions during long math tasks differ in signal complexity detectable with HFD, which further, can be classified by a ML model.

Participants
Thirty-four math experts (bachelor and master students in math or mathrelated disciplines, like physics or engineering) and thirty-five math novices (no university-level math studies) participated in the experiment.However, eleven participants from the group of math experts and twelve participants from the novice group were discarded from the data analysis because their EEG data was too noisy, or some of the relevant data was missing due to malfunctioning EEG amplifier.Therefore, in the group of math experts, there were 22 participants (5 female and 17 male), and in the novice group, 22 participants (7 female and 15 male).The background of the participants was screened by a math questionnaire.
The age of the participants ranged from 19 to 24 years (mean 21.0 years) among math experts and from 19 to 35 years (mean 23.8 years) among novices.All participants in both groups were right-handed.No participants reported hearing loss nor history of neurological illnesses.The experiment protocol was conducted in accordance with the Declaration of Helsinki and approved by the Executive Board of ETH Zurich after a review by the ETH Zurich Ethics Commission.All participants provided written informed consent.

Task design
Participants watched 16 math demonstrations.After each demonstration they were asked three self-evaluation reflections to which they answered by pressing a button in a 4-button response box.Each set of trials consisted of four excerpts of the same presentation style (symbolic or geometric), and these sets were presented in a pseudo random order via a monitor.The pseudo randomization defined the presentation order (symbolic first or geometric first).However, each participant saw the same four math demonstrations presented in both symbolic and geometric form before seeing them in the other form.
Each math demonstration consisted of several slides, from 4 up to 12 slides (6.9 slides on average) depending on the complexity of each demonstration.The total duration of math demonstrations varied from 13 seconds to 68 seconds (33.1 seconds on average).The timing of each slide was the same for all the participants.The duration of each slide was defined according to an online screening in which 25 math experts and 25 math novices watched the math demonstrations slides and auto-regulated the following slide with a button press.The participants who attended the online screening did not attend the actual EEG experiment.The duration of each slide in the EEG experiment was the average time the participants spent on each slide during the online screening.In the online screening, there was no statistically significant difference between experts and novices in the duration of time spent on each slide.

Data acquisition
The stimuli were presented to the participants with the MATLAB via Psych-Toolbox.The experimenter launched the playback of the presentation program after which participant could navigate to the math demonstrations by a button press once they had read the instruction slides on the screen.The total length of the experiment material was approximately 15 minutes.
The data were recorded using Ant Neuro eego mylab electrode caps with active 128 EEG channels 1 .
Four external electrodes placed below, above and on the left side of the left eye and on the right side of the right eye.The offsets of the active electrodes were kept below 30 mv at the beginning of the measurement, and the data were collected with a sampling rate of 2048 Hz.A timestamp (trigger) was marked into to EEG data at the beginning of each slide of the math presentations.The triggers were sent wirelessly via Lab Streaming Layer2 .

Data pre-processing
The EEG data of all the participants were first preprocessed with EEGLAB (version 2019.1;Delorme & Makeig, 2004 [43]).The reference was set as the average of all the EEG electrodes.The data were high-pass filtered at 0.5 Hz and low-pass filtered at 40 Hz.Finite impulse response (FIR) filtering, based on the firls (least square fitting of FIR coefficients) MATLAB function, was used as a filter for all the data.Then, the data were treated with independent component analysis (ICA) decomposition with the runica algorithm of EEGLAB (Delorme and Makeig, 2004 [43]) to detect and remove artefacts related to eye movements and blinks.ICA decomposition gives as many spatial signal source components as there are channels in the EEG data.Typically, one to four ICA components related to the eye artefacts were removed.Noisy EEG data channels for some participants were interpolated.

Feature extraction 2.5.1. Higuchi Fractial Dimension (HFD)
The EEG time-series has a duration between 10-20 minutes, resulting in a large data size per sample.Hence, feature extraction is necessary to capture relevant information.The extracted feature are then used to draw conclusions regarding the relevance of each brain area for mathematical calculations.For this purpose the fractal dimension (FD) [44] for each sample is calculated and is used to measure the complexity of the signal.A simple pattern that is repeating continuously can become a very complex series which is the basis for the fractal constructs.A fractal is a shape that retains its structural detail despite scaling and is the reason why complex objects can be described with the help of fractal dimension.One variant of FD, the Higuchi's fractal dimension, [37] has its roots in chaos theory and has been successfully applied as a complexity in various domains of signal processing.It has been shown to be a good numerical solution to nonlinear signals [45].The speed, accuracy, and cost of applying the HFD method for research and medical diagnosis make it stand out from the widely used linear methods [46].Among the different FD algorithms Higuchi's method [45] demonstrates to be a more accurate option for EEG signals, since it is accurate for stationary and non-stationary signals.
Say X is an EEG signal of length T is the length of a time window on which we calculate a HFD value.Following [37], we calculate HFD as: A new signal x k m is constructed from X, with window size N where m = (1, 2, ..., k) denotes the starting point and k = (1, 2, ..., k max ) the interval size: L m (k) describes the length of the curve of x k m for every k given m: where is the normalization factor.Length L(k) is defined by the average of the k lengths: HFD is the slope of the best fitted curve between all the data points of timeseries X for a given time window N for for k = (1, 2, ..., k max ) between log(1/k) and log L(k): It is possible to calculate HFD for the whole signal (T = N ).However, this is not recommended if the signal is nonstationary.In such cases the HFD value does not represent the true measure, and division into windows (or segments) is advised.In [47], Accardo and colleagues have shown on synthetic fractal signals that Higuchi's algorithm is more efficient, faster, more accurate and able to estimate fractal dimension for short segments, compared to Maragos and Sun's algorithm proposed in [48].

Hyperparameter tuning
An important hyperparameter that requires finetuning is k max .There is no agreed methodology to optimize this parameter [49].As per equation 3, HFD is summed up to k max , therefore increasing k max will lead to an increase in HFD.A poor choice of k max will result in uninformative HFD, thus, it has to be carefully tuned.
We propose the following methodology to identify the best value for k max : 1. We compute the HFD values as per equation 4 for a wide range of k max values, i.e., k max ∈ 2, 5, 20, 100, 150, 200, 400 over all subjects and presentations.2. We identify the k max at which the difference (equation 5) between HFD values of significant and non-significant channels is maximized.Significance/nonsignificance is assessed by taking the maximum/minimum HFD value across all electrodes for a subject.Here, the minimum value is understood as the baseline fractal dimension and is therefore subtracted from the maximum value, which is the complexity of the relevant channels.We base this requirement on the assumption that certain EEG regions are more relevant than others for the mathematical tasks.Hence, there will be a difference in HFD values and we want to select the k max that maximizes this difference.3. The k max value that satisfies requirement 2) and 3) is chosen to compute the HFD values for further analyses and for the machine learning classification.

HFD features analyses
Estimating HFD values for each channel of each participant allows to investigate which brain areas are most active while performing mathematical tasks.Since HFD values have no physical interpretation, a relative comparison between two different groups is performed.
First, a comparison between experts and novices is investigated, by taking the average of all HFD values of the expert group and the novice group and subtracting them from each other: where j ∈{1,...11} is the index of experts and novices, respectively, k ∈{1,...16} is the index of presentations and i ∈{1,...129} is the index of EEG channels.
A one-sided t-test is calculated, testing whether there is a significant difference between the two groups.A visual heatmap of the difference between experts and novices based on equation 5 is mapped onto the head for better qualitative interpretation.
Subsequently a more fine grained analysis is performed by comparing the difference between expert and novice for algebraic and geometric separately: where k A and k G ∈{1,...8} is the index of the algebraic and geometric presentations respectively.

Machine learning classification
We posit the question if a prediction can be made whether a new subject is a novice or an expert based on EEG recordings while performing mathematical tasks.We frame this problem as a two-class classification task.To understand and interpret the outcome of the machine learning classifiers, care needs to be taken while generating the classification dataset and splitting it into training and testing sets.
We first define the classification-dataset as a collection of subject-presentation pairs (e.g.Expert1-Presentation1A etc.).Together with the 16 presentations, the full dataset include 704 samples, i.e., subject-presentation pairs.Subsequently, we calculate either a unique HFD value per EEG channel, meaning that each sample consists of 124 HFD features, or divide the EEG signals of total length T into non-overlapping windows of length N and calculate a HFD value for each window leading to (T/N)*124 HFD features.To be noted that the channels "VEOGL", "HEOGL", "HEOGR", "VEOGU", "HEART" are discarded, since they do not record brain signals but eye movement and cardiac activity.
Since this work is the first in the literature to attempt an automatic classification of mathematical cognitive behavior, we propose three different cases of dataset splitting, illustrated in Figure 1: 1. Subject-presentation pairs: We randomly split all 704 samples without considering whether a sample is coming from different subjects.This means that the samples from the same subject can either be entirely in the training set or in the validation set, or partially in the training and in the validation set.With case 1, we verify if the machine learning (ML) classifier is able to discern between the 22 experts and novices present in the dataset based on a single mathematical presentation.With case 2, we validate the ML classifier on new subjects of which data it has never seen before.The former is a relatively easier classification task, but necessary as a first proof-of-concept, whereas the latter tackles the most challenging problem of inter-subject variability common to all biomedical data.With case 3, we analyze whether a prediction can be made based on samples coming from a single presentation.By training a separate classifier for each presentation, we can compare the classification accuracy among the presentations and draw insights about which mathematical presentation is more suitable for discerning between math novices and experts.
For cases 1 and 2 we calculate a single HFD value per EEG channel throughout the whole duration of the presentations.This choice is motivated by the fact that all presentations, of different recording lengths, belong to the same dataset on which a machine learning classifier is trained on and, in general, the classifiers require a fixed numbers of features.This is no longer an issue for case 3, because each sub-dataset consists of data from a single presentation of fixed length.Hence, we can increase the granularity and use a non-overlapping moving window of length N to calculate the HFD value in equation 4 for each   window.More precisely, a HFD value is calculated every N seconds of the duration of the presentation HF D 1:N , ..., HF D t:t+N with t being time steps.This allows to analyze the temporal evolution of the presentation and draw conclusions regarding the classification differences.We test several values of N, i.e., 5, 8, 11 seconds.
Once the datasets are prepared, we proceed with classifiers training using the scikit-learn Python package.We investigate several ML algorithms including Nearest Neigbours, Linear SVM, Decision Tree and Adaboost.We first optimize the classifiers by tuning the hyperparameters under case 1, i.e., subjectpresentation level.Once the optimal parameters are found, we keep them for case 2 and 3.The various ML algorithm tested are summarized in their corresponding parameters ranges.Once the best performing ML algorithm has been identified, we further optimize it with a grid-search algorithm.Given the small sample size, 10 fold cross-validation (90 percent training/ 10 percent validation set) has been applied with a fixed seed.

Results
As described in the introduction, extracting the neural signature of math experts and novices requires careful features extraction via the HFD method.calculate the HFD correctly, hyperparameter k max requires finetuning.Therefore, section 3.1 presents the optimization results of hyperparameter k max .
Based on the extracted HFD features, experts and novices are compared in section 3.2 giving insights which brain region is relevant for performing mathematical tasks.Finally, based on the features, classification results between experts and novices are shown in section 3.3.

Optimal k max
Figure 2 shows the value of HFD for all subjects averaging over all channels for different values of k max .HFD is steadily increasing but starts to plateau at a value of 100. Figure 3 shows the difference between the maximum and minimum HFD values for different k max with accordance to equation 5.It can be observed that the difference in HFD value corresponding to k max reaches a peak at 20 and 100 and progressively declines with increasing k max .Based on the fact that HFD is plateauing at k max equal to 100 and the largest difference between the maximum and minimum HFD values is also found at the same value, k max = 100 is used for further analysis.

HFD feature analyses
Figure 4 shows the difference between the average HFD values between experts and novices, for the top 10 channels that present the highest difference between expert and novices.All top 10 channels are statistically significant under p = 0.05 constraint.All channels are depicted in form of a heatmap in Figure 5.The dark blue shaded areas indicate the highest positive difference between expert and novices.
The subsequent more finegrained analysis comparing the difference between expert and novice for algebraic and geometric is shown in Figure 6 given equation 6.Although there are differences between algebraic and geometric presentations, none of them is statistically different under p-value 0.05 hypothesis.

Average over all Algebraic presentations
Average over all Geometric presentations Figure 6: HFDmax-HFD min calculated as average for all Algebraic vs Geometric presentations for all channels as defined in equation 6

Expert/Novice classification
Table 2 summarizes the classification results between expert and novices.On a subject-presentation split the accuracy reaches 97% demonstrating that it is possible to automatically classify between math experts and math novices based on their electroencephalogram (EEG) signals while watching math demonstrations because the ML model can successfully learn each subject's brainwaves signatures.However, when we split the training and test sets on a subject level, meaning that we increase the difficulty of the task by introducing inter-subject variability that is well-known to be challenging in biosignals classification, i.e., the trained model is validated on new subjects whose data it has never seen before, the accuracy falls to 66%.

Subject
So far the results are shown by considering all presentations for each subject, i.e., the calculated HFD features for all presentations are concatenated for the final classification stage.We suspect that the poor classification accuracy could be partially caused by some of the presentations that might perform poorly.Hence, we perform presentation-specific classification on subject level and the classification accuracy improves up to 79% (presentation 7A).
Figure 7 and Figure 8 show the HFD values when window size of 8 seconds is applied for the presentation with the highest (presentation 7A) and the lowest (presentation 4G) classification accuracy.The difference in classification accuracy may be explained through a better separation between Experts and Novices.

Discussion
Advantages of ML for brain research include the data driven approach which enables generation of hypotheses about underlying brain processes in rest or in active engagement with a cognitive or emotional task.Such underlying processes are sometimes impossible to detect by experts' observations.ML also enables explorations of new paradigms with respect to their neurophysiological signatures (Lemm et al., 2011).One of such new paradigms is naturalistic study design which aims to understand the brain during real-life tasks, like when solving complex math.
Our novel approach on applying ML to EEG data recorded in math experts and novices during complex math encourages to expand the usage of data driven brain imaging methods from healthcare to education.Our approach utilizing nonlinear HFD, which measures signal complexity, was reliable in describing the data by systematically detecting the difference in the neural signature of math experts and novices with a 98% cross-validation accuracy.However, the results gained with ML discriminative algorithm were mixed and showed 50-80 percent classification accuracy when tested with unseen subjects.
Nonlinear fractal dimension methods seem ideal for tracing fluctuations in biological systems, including the brain, which are nonlinear by nature.HFD is a measure of signal complexity in the time domain (Higuchi, 1988 [37]; Spasic et al., 2008 [38]) and has been successfully applied for brain state analysis of EEG in sleep, drowsiness, and wakefulness (Inouye et al., 1994 [39]; Klonowski et al., 2005 [35]; Peiris et al., 2006).Our results gained with HFD show a difference in the neural signature between math experts and novices during long and complex math tasks with a high classification accuracy.These results encourage to use the HFD method in detecting subtle differences in the brain states, like those of math experts and novices, which go beyond the more drastic differences in the brain states during the levels of arousal, like sleep stages, or drowsiness and wakefulness.
Despite the successful classification to experts and novices based on HFD was relatively stable for the entire dataset, the ML model adapted poorly to unseen subjects, and we could not overcome the overfitting and high generalization error caused by inter-subject variability.The most important reason for such a poor generalization is that our dataset is incorrigibly small to be divided to the training and test sets on a subject level.In healthcare, big data platforms are being formed increasingly (Eickhoff et al., 2016; Zbontar et al., 2019), and it is important to take similar steps to create large and clearly labeled open data pools for educational neurosciences.
Our small dataset may function reasonably well for method development of data-driven approaches, since the differences between math demonstrations are statistically significant especially over several frontal electrodes showing higher frontal signal complexity in math novices in comparison to experts.Cognitively, these results may indicate novices' stronger recruitment of domain-general processes in comparison to experts, which is in line with previous literature (Amalric and Dehaene, 2016 [18]; Wang et al., 2020 [17]).
Some studies have investigated for the connection between nonlinear FD methods and linear oscillation analyses over delta, theta and alpha bands.These studies show a dependence between the nonlinear and linear methods and suggest that the most reliable results are gained when combining nonlinear and linear methods to classify different brain states (Acharya et al., 2005 [18]; Šušmáková and Krakovská, 2008 [40]).Since combination of nonlinear and linear methods seem to bring the most robust classification results, we could combine the HFD and oscillation analyses and feed the combined information to a machine learning model.Our novel analysis with machine learning utilized only fractal dimension; however, we report on other papers the brain oscillations for the same dataset (Formaz et al., unpublished data; Poikonen et al., 2022 [9]).
Another interesting way to deepen the analysis of our dataset was to break the temporal data stream to segments.With a larger dataset and statistical power, time points during which the neural signatures of math experts and novices differ the most could potentially be found.This data-driven approach may have practical implications after detecting whether the cortical functions of experts and novices differ the most at the beginning, at the end, or at some other time point during the math demonstrations.With our dataset, ML algorithm showed 50-80 percent classification accuracy for unseen subjects when breaking the data to a temporal stream.Such a high variation may be explained by a small dataset, or by a combination several features related to the length, content, and difficulty level of the math demonstrations.
Understanding which parts of the math demonstrations to emphasize when teaching complex math may be helpful in supporting students' development towards math expertise.Such time-dependent information may be hard to collect with questionnaires or other behavioral measures, and therefore, brainoriginated data-driven methods may be the only way to access such information in the context of learning.Further, these ML models could be used to create learning contexts in which adaptive feedback is given to adjust to the individual needs of a learner or those of a specific group during collaborative learning, building on the previous examples like BCI applications for post-stroke motor rehabilitation, or relatively simple neurofeedback applications for focused attention or working memory (Cervera et al., 2018 [11]; Kefalis et al., 2020;Hunkin et al., 2021 [25]).Simple options for BCI interventions for the math demonstrations used in our study might be to adjust the velocity of presenting new information, or by scaffolding the learning process via instructions or remarks depending on the EEG signal of the learner.

Limitations
Our novel paradigm combining mathematical cognition, cortical activity and ML is exploratory in nature and we recognize the following limitations.First, the most drastic limitation is the small dataset in use.The straightforward way around it would be to increase significantly the amount of data, e.g., by at least doubling the number of participants.The more data the better we can estimate the real data distribution of the general population.The second limitation is related to the classes chosen for the ML classification.We chose to compare two groups of participants during the same cognitive task.Other strategy for a small dataset would be to explore individual differences, for example, by aiming to classify the data excerpts of resting state and cognitively active state for each participant.Earlier studies show that differentiation of brain states for an individual participant during simple sensory tasks is rather robust whereas the generalizations of the cortical activation patterns across a group of participants, and during complex cognitive tasks, is challenging.However, such individual brain state classification would not give us hardly any insights for the expert-novice differences during mathematical cognition.As the third limitation to consider, when preprocessing, we chose to band-pass filter the data with a bandwidth of 0.5-40 Hz due to the contamination of the data with the 50 Hz line noise.HFD is associated with changes in delta, theta and alpha oscillations which all were included in our analysis.However, also gamma oscillation is known to be important during cognitive tasks (Herrmann et al., 2004), and it has been connected to HFD.Due to bandpass filtering chosen, gamma activity is not included in our analysis.Based on previous literature, HFD seems the most stable fractal dimension methods (Kesic and Spasic, 2016 [45]).However, as the fourth limitation of our study, is the general criticism for the HFD that it has a short margin of scale which may give the same complexity number to signals with only subtle differences.For detecting the possibly small differences in the cortical activity of math experts and novices, some other method with more detailed scale may be more suitable.Fifth, for the cross-validation, different models could be compared to find a model with ideal complexity which balances between overfitting of an unnecessarily complex model and simple model's inability to adapt to the details of the complex cognitive data.Ideally for ML algorithms, each sample (e.g.EEG data collected during each math demonstration) would have the same number of data points (e.g. the same duration).However, it is difficult to realize in practice due to different duration it takes to solve different naturalistic math tasks.In the future, research of brain processes during abstract cognition might be conducted, for example, within a video game context, in which the duration is easier to match to be the same over all the rounds played.
Conclusions The present study used a unique paradigm to compare neural correlates of math experts and novices while solving naturalistic math demonstrations.Overcoming limitations of previous studies with reductionist stimuli and linear EEG analysis methods, the brain functions during abstract cognition were measured with a high-density EEG during long and complex math demonstrations and analyzed with a relatively rigor nonlinear method, HFD.Our results indicated that math novices have a higher signal complexity measure with HFD than experts over several frontal electrodes suggesting a stronger engagement of domain-general brain functions.Further, we explored ML algorithms for classifying math experts and novices based on their neural signature.These results were promising but we also acknowledge the inevitably small dataset we had in use for consistent results.We encourage taking example from brain imaging databases created in healthcare for a creation of a similar database for educational neuroscience.In the future, application possibilities for such a database and deep learning lay in data-driven theory formation for normal and disrupted learning and development, and adaptive feedback systems for learning contexts.
Data and code availability statement The dataset analyzed in the present study as well as scripting and plotting code are available from the corresponding authors via email on request. CRediT

2 .
Subject-specific: We split the dataset on the level of subjects, meaning that all subject-presentation pairs of the same subject are either in the training or validation set.3. Presentation-specific: We deal with each presentation as a separate machine learning task.In other words, we divide the full dataset into subdatasets, each of which consists in a single presentation, and perform the training and testing procedure on each of the sub-datasets.

Figure 2 :
Figure 2: HFD value, averaged across all channels all subjects and presentations for different values of kmax.

Figure 3 :
Figure 3: HFD value, difference between HFD values between the maximum and minimum of all channels averaged across all subject and presentations for different values of kmax.

Figure 4 :
Figure 4: Top 10 channels with the highest difference between their HFD values.Asterisk * indicates that the average value of HFD expert is statistically different than HFD novice under p=0.05 threshold for that specific channel.

Figure 5 :
Figure 5: Heatmap of HFD difference between Expert and Novices.Darker blue color indicates the brain areas where the positive differences between experts and novices are the highest.

Table 1 :
Machine Learning algorithms used for classification between experts and novices

Table 2 :
Classification results between experts and novices based on different classification algorithms for Subject-presentation pairs, Subject-specific and Presentation-specific split.All results are based 10 fold cross validation and averaged over 3 random seeds.