Spatiotemporal Dynamical Analysis of Brain Activity During Mental Fatigue Process

Mental fatigue is a common phenomenon with implicit and multidimensional properties. It brings dynamic changes in functional brain networks. However, the challenging problem of false positives appears when the connectivity is estimated by electroencephalography (EEG). In this article, we propose a novel framework based on spatial clustering to explore the sources of mental fatigue and functional activity changes caused by them. To suppress false positive observations, spatial clustering is implemented in brain networks. The nodes extracted by spatial clustering are registered back to the functional magnetic resonance imaging (fMRI) source space to determine the sources of mental fatigue. The wavelet entropy of EEG in a sliding window is calculated to find the temporal features of mental fatigue. Our experimental results show that the extracted nodes correspond to the fMRI sources across different subjects and different tasks. The entropy values on the extracted nodes demonstrate clearer staged decreasing changes (deactivation). Additionally, the synchronization among the extracted nodes is stronger than that among all the nodes in the deactivation stage. The initial time of the strong synchronized deactivation is consistent with the subjective fatigue time reported by the subjects themselves. It means the synchronization and deactivation correspond to the subjective feelings of fatigue. Therefore, this functional activity pattern may be caused by sources of mental fatigue. The proposed framework is useful for a wide range of prolonged functional imaging and fatigue detection studies.


I. INTRODUCTION
M ENTAL fatigue represents a kind of neurophysiologic state caused by prolonged cognitive load and easily leads to decrease in mental and physical performance capacity [3], [4]. It is multidimensional in nature [5]. From the perspective of information flow, Desmond and Hancock [6] made a distinction between active and passive fatigue states that reflect different styles of perceptual-motor response requirements. From the perspective of psychophysiology, the state change has subjective and objective manifestations, which include an increased feeling of tiredness, an increased propensity toward less analytic information processing, and alterations in brain activity [4], [7]. Therefore, fatigue has implications for many aspects, not only in the researches but also in daily life [4]. In the studies of cognition and psychology, the commonly used stimulus sequences certainly will induce fatigue more or less [1], [3], [8]. In the workplace, mental fatigue has been found to have a strong association with errors and inadequate performance adjustments after errors [4], [9], [10]. In the transportation system, mental fatigue gives rise to road crashes [11]. It causes 10%-20% of all traffic accidents every year all over the world. Despite these multifaceted effects, still little is known about the brain mechanisms underlying mental fatigue and previous studies did not embody the implicit information of fatigue from multidimensions [5], [9], [12]. For example, where are the cortical generators of mental fatigue (in spatial dimension) and how do the source activities in our brains change (in time dimension), when we feel tired?
Changes in functional connectivity across mental states can provide richer information about the human cognition than simpler univariate approaches [13]. Recently, several studies have reported that in adults with chronic fatigue syndrome, subjective reports of fatigue are correlated with decreased functional connectivity between the salience network and the default mode network during resting state [14]- [16]. Sun et al. [13] also revealed the degraded performance of small-world features in a mentally demanding test of sustained attention, providing further support for the presence of a reshaped global topology in cortical connectivity networks under fatigue state. Langner et al. [17] attributed the decrease of cerebral blood flow in the network of areas to a depletion of neural resources, or an inability to retrieve these resources related to mental fatigue.
In spite of the existing contributions, significant challenges in the functional connectivity analysis of time series 2379-8920 c 2020 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission.
remain [18]. Traditional approaches of functional connectivity analysis, which can be used to explore the brain mechanisms underlying mental fatigue, focus on the long time behavior of at most a few dynamical variables, characterizing either different single connections or the system's average behavior [19]. In addition, the connectivity analysis from sensor space data is generally not reliable and robust. This issue is even serious for the measure of electroencephalography (EEG) which records brain electrical activities noninvasively at a millisecond time scale. Because of signal mixing, which translates to volume conduction in the case of EEG recordings and to signal leakage in source reconstructed EEG data, a spatially widespread group of sensors detects the activity of any single neuronal source [20]. This means that the same source may be detected by multiple sensors and the same sensor may always pick up multiple sources, resulting in false positives such as spurious interactions (referred to as ghost interactions). Even though source reconstruction can be used to reduce signal mixing and elucidate the likely neuroanatomical sources of the EEG signals, so far no source reconstruction approach can yield an unambiguous estimate of the source topography due to ill-posed nature of the EEG inverse problem [21]- [23]. The residual signal mixing in source space, signal leakage, is quantitatively dependent on the source-reconstruction method of choice but qualitatively characteristic to all such methods [23].
Since EEG source reconstruction has the inevitable ill-posed problem, we explore a solution through dual-modal imaging. With high spatial resolution, the functional magnetic resonance imaging (fMRI) is a powerful tool to probe and monitor human full-brain activity for cognition. On the basis of fMRI localization precision in brain imaging space, Han et al. [24] developed a dense individualized and common connectivitybased cortical landmarks (DICCCOL) system to map and annotate the large-scale functional networks. They proved that it could provide more comprehensive, accurate, and reliable Regions of Interests (ROIs). Nevertheless, it is impractical to acquire large-scale task-based fMRI data for the same group of subjects because of the cost and time constraints [24]. An interesting observation from the results in [25] is that structural connection patterns and the same functional ROI obtained by fMRI are reasonably consistent across different subjects. Furthermore, Zhang et al. [26] established a framework to predict ROIs from the diffusion tensor imaging (DTI) data of new subjects according to the correspondence relationship of the dual-modal imaging results. It avoided the large-scale taskbased fMRI acquisition and the average prediction error was around 3.94 mm.
This motivates us to develop an EEG connectivity analysis approach based on spatial clustering to explore the sources of mental fatigue and functional activity changes caused by them. The EEG acquisition is adopted because it is suitable for dynamic analysis of mental fatigue in the prolonged task with millisecond time resolution. The sources of mental fatigue are determined by the correspondence relationship of EEG (scalp electrical activity) and fMRI (source localization) results. To suppress the false positive observations of EEG coupling, we consider that all the links connected to a sensor (node) is the basic unit of connectivity analysis and not a single link or the node itself. The same source may be detected by multiple sensors. However, different basic analysis units contain different information of the source, because every time we consider the relationship between the signals from a sensor and the other sensors. The basic analysis units can be clustered spatially over time with the progress of fatigue to acquire more precise determination of local areas (represented by nodes) modulated by fatigue. Thereafter, a basic analysis unit will be referred to a cluster. To determine the sources of mental fatigue, we attempt to register the nodes extracted by spatial clustering back to the fMRI neuroimages. The fMRI neuroimages (alternative sources) come from the previous fMRI studies. The Mahalanobis distance is utilized to measure the matching degree of two modality results from different subjects and even different tasks. We will show the consistent performance between the results from the two modalities under the same dominant factor (i.e., fatigue).
The signal flowchart of the proposed methodology is shown in Fig. 1. The kernel of the methodology is spatial clustering. Through spatial clustering, the spatiotemporal analysis of mental fatigue is implemented. After EEG preprocessing, brain networks were constructed based on EEG connectivity estimation. The nodes of the discriminative clusters were extracted through spatial clustering to reduce the influence of false positives. In the spatial dimension, the sources localized by the previous fMRI studies related to mental fatigue are introduced as alternative sources. The nodes extracted by spatial clustering are used as the fiducial points to realize the coregistration with the alternative sources. The fMRI sources with the shortest Mahalanobis distance to the fiducial points are determined as the sources of mental fatigue. In the temporal dimension, the EEG features based on wavelet entropy are extracted from the nodes after spatial clustering. Finally, the temporal nonlinear dynamic features on the extracted nodes were linked to the neuroimages so that mental fatigue can be examined over time on task. The data-driven results are verified regarding the functional significance of the spatiotemporal dynamics and the subjectivity of mental fatigue.
This remainder of this article is organized as follows. In Section II, we describe the experimental details and EEG data (see Sections II-A and II-B), elaborate functional connectivity analysis and spatial clustering applied to suppress signal mixing and leakage to find the fiducial points (see Section II-C), illustrate temporal feature extraction according to the fiducial points to examine mental fatigue over time on task (see Section II-D), and discuss the spatiotemporal imaging used for fMRI coregistration based on the fiducial points to evaluate mental fatigue in both time and space domains (see Section II-E). The results of the study are presented in Section III and discussed in Section IV. Finally, this article concludes in Section V.

A. Experiments and Data
This article was reviewed and approved by the Ethics Committee, Dalian University of Technology (protocol number: 2018-040). Written informed consents were obtained from all participants before the experiments. The electrophysiological data (EEG) were collected from twelve right-handed subjects (age range: 20-35) while performing the simulated driving task lasting for more than one hour under a monotonous highway environment produced on the mediumfidelity driving simulator of Dalian University of Technology. The subjects are all college students. The driving simulator mainly consisted of a cockpit, a host computer, a CRT monitor, an audio system, operating sensors, and a data collector. The cockpit contained the same operational components as the real vehicle, such as a steering wheel, a clutch, a foot brake, an accelerator pedal, and a hand brake. The simulator had both automatic and manual transmission (including five forward gears and one reverse gear). To reduce the movements of the subjects and make the task boring, we asked the subjects to drive using the automatic transmission. Their responses to the driving scene stimuli are the control outputs of the vehicle to avoid accidents. When the subjects felt fatigued, they selfreported this in real time by telling it to the researcher and the corresponding times were recorded by defining timestamp to the data. The temporal indication of the events was used as the criterion (critical position) from subjective fatigue [ Fig. 8(d)]. The subjects who already had the driving licenses drove in a monotonous highway environment produced on the driving simulator. There is no difference of driving environment in the simulator for subjects with different lengths of driving experience. Altogether, 10 out of 12 subjects provided the exact temporal indication of their subjective fatigue and avoided the accident in the experiments. Thence the EEG data throughout the fatigue process from these ten subjects were used. Nineteen standard electrodes (i.e., Fp1, Fp2, F7, F3, Fz, F4, F8, T3, C3, Cz, C4, T4, T5, P3, Pz, P4, T6, O1, and O2) were attached to the scalps following the International 10-20 system to collect the subjects' EEG data. The EEG's sampling frequency is 200 Hz.

B. Data Preprocessing
Since some studies have demonstrated the high correlation between alpha activity and fatigue [27], [28], we used wavelet decomposition with four levels to obtain the subband EEG waves (i.e., alpha). Wavelet transform is an effective tool to analyze the various components of a nonstationary signal due to its attractive properties, such as good local representation in both time and frequency domain and multirate filtering (differentiating the signals having various frequencies) [29]. The decomposition based on the wavelet transform at each level output approximate and detailed components corresponding to scale and wavelet coefficients, respectively. The detailed component at the fourth level was roughly within the alpha range (8)(9)(10)(11)(12)(13). The discrete Meyer (dmey) was selected as a mother wavelet because of its resemblance to the EEG signal in our case.
As the decomposed wavelet coefficients reflect detailed information, the signal reconstruction in alpha band has a certain suppression effect on low-frequency artifacts, such as large body movements and respiration [see Fig. 2 Besides, the wavelet coefficients also represent the correlation between the signal and the basic functions in terms of the high-frequency content. The signal will generate high amplitude coefficients at places where the artifacts [e.g., the blinks in the red boxes of Fig. 2(b)] are present. We utilized the following threshold technique to remove them. It has proven to be effective in the research about driver fatigue analysis [30], [31] where C j denotes the wavelet coefficient at the jth level of wavelet decomposition and T j is its threshold. If any of the coefficients is greater than this threshold, it is halved to reduce the impact of artifacts. Then, the wavelet-corrected EEG signals are obtained by reconstruction from the processed coefficients. To reduce the effect of the artifacts further, high amplitude coefficients actually need to be processed recursively. In this analysis, the wavelet coefficients with artifacts were halved no more than twice taking into account the efficiency and effectiveness of the algorithm. The EEG signals from all channels were preprocessed in the same way as mentioned above. The comparison of the EEG signals before and after preprocessing is shown in Fig. 2. Fig. 2(a) gives a typical 10-s original EEG signal disturbed by a variety of artifacts. Fig. 2(b) shows alpha waves extracted by wavelet decomposition. The relatively large fluctuations which may be caused by body movements were apparently eliminated, whereas a few of artifacts are still left (e.g., the artifacts probably caused by blinks in the red boxes). Fig. 2(c) shows the result of the wavelet-based artifact removal. Compared to the original signal, the corrected signal is not only free from the dominated noise contamination but also retains the main trend. It is reflected by the alpha rhythm nonlinear oscillations.

C. Source Analysis Model
Until now, all the EEG source localization methods are developed based on the generative model of EEG data [32], [33]. According to Maxwell's equations, the EEG signals measured at K EEG electrodes are linear functions of M dipole amplitudes on the cortex (i.e., sources), represented by primary current density. The generative model of EEG data is given by where (t) ∈ R K is the EEG signal measured at K EEG electrodes at time t. J(t) ∈ R 3M is the activity of N sources with 3-D orientation. L ∈ R K×3M is a lead matrix summarizing the propagation of M sources throughout the brain. ε(t) ∈ R K denotes the measurement noise. Given the geometry and electrical conductivity of the majority tissues in the head (between the electrodes and the sources), L can be calculated. This process is called forward modeling. Reversely, the sources J(t) can be estimated from the measured signal (t) and calculated L by inverse modeling. This reverse process is also called source reconstruction.
The EEG inverse problem of source estimation from EEG on the scalp has ill-posed nature because the solution is underdetermined for all admissible output voltages. The number of dipoles M is usually much larger than the number of electrodes K. Additionally, the solution is unstable, because it is highly sensitive to small changes in the noisy data [34]. In the present study, we utilize the fMRI coregistration to avoid the ill-posed problem. From (2), we know that there must be the related specific EEG activities on the scalp, if the sources J(t) are given, and vice versa. We do not solve this equation. The alternative sources of J(t) come from related fMRI studies. The specific EEG activity is revealed by the following functional connectivity and spatial clustering (see the next section). The sources J(t) are determined by registration back to fMRI source space after spatial clustering [see Fig. 3(a)].

D. Functional Connectivity and Spatial Clustering
As shown in Fig. 1, the key method which can be used to integrate spatial and temporal information of fatigue is divided into two parts: 1) brain network constructing and 2) spatial clustering on the constructed brain networks.
Representing the complex system brain as a network requires identification of nodes and edges [35]. Here, we assign EEG electrodes to the nodes of the brain networks. The adjacency relations among the nodes in the networks can be described by the adjacency matrix A whose element A(i, j) shows the measured edge between electrodes (nodes) i and j. Since the electrodes for EEG recording are distributed in different brain regions, the edges embody functional connectivity among these brain regions and this kind of brain network is known as a functional brain network. A simple and most commonly used measure of the functional connectivity is the cross-correlation function [36]. The correlation between EEG signals s i and s j can be calculated by the following equations: where time delay τ is set to 0. γ ij is corresponding to the element of the cross-correlation matrix R which presents in the ith row and the jth column. s i and s j represent the EEG signals from electrodes i and j, respectively. The adjacency matrix A can be obtained by the correlation matrix R. First, define the adjacency matrix A = R. To exclude self-connections of nodes, the elements on the main diagonal of A were set to zero. The other elements of A reflect functional coupling (i.e., functional connectivity) of the signals in different brain regions.
In the present study, the spatial clustering on brain networks is defined as a process of classifying a set of networks' clusters into groups in which similarity in the same group is maximized in terms of distance formed by the functional connectivity between the clusters. A cluster consists of a node and all the connections to it. We employ the number of the node to indicate the number of the cluster and consider the connections as the node's attributes.
Hence, spatial clustering is also a process of electrode selection considering multiple interaction attributes. Mathematically, a brain network with N nodes, each of which has M links (M = N − 1 for weighted network), is denoted by is a scalar represents the strength of the functional connectivity between node i and node j.
The following is based on agglomerative hierarchical clustering which establishes the foundations for inducing a hierarchical clustering from a newly represented, or newly encoded, mapping of functional connectivity matrix A.
To begin with, each node of a brain network is used to represent a cluster at the lowest level, i.e., Then, based on the distance (similarity function) calculation, two closest clusters are successively merged to reduce the number of clusters by 1 until a single cluster remains at the highest level (i.e., k = 1) where C i and C j denote two different clusters; k is the level number; Dendrogram indicates the cluster set and its element number is reduced by 1 to reach a higher level in (7); and d represents the Euclidean distance between two clusters. Equations (5) and (6) means when the Euclidean distance of C i and C j is minimal among k clusters, a new cluster C 2N +1−k is generated. At the next level, C i and C j are replaced by C 2N +1−k . C 2N +1−k consists of C i and C j representing a merging operation. The Euclidean distance can be calculated by the following expressions: where M denotes the number of attributes of the node vectors X a and X b . x m a and x m b indicate the attributes (corresponding to functional connectivity) of the two node vectors.
Finally, we cut the dendrogram to complete clustering on brain networks to obtain different node groups whose number (Z) is predetermined. As the number of the clusters contained in Dendrogram is equal to the level number k, Z clusters will be left by the Zth level cut through the dendrogram. That means we classify the nodes of brain networks into Z groups according to the elements (clusters) of Dendrogram Z . In this article, Z is set to 2 to find the node cluster in the specific brain areas associated with mental fatigue (Fig. 5). As mentioned in Section II-C, the EEG activity related to the sources J(t) has a specific topology and is different from other activities. Therefore, Z is set to 2 so that two groups of clusters are left.
The specific one is determined by fatigue information (see Figs. 3, 5, and 6). In this article, the brain network construction and spatial clustering were performed in a sliding window with a length of 30 s without overlap to obtain the evolution of the node clusters over time. The result of the spatiotemporal dynamical analysis was like a movie. Every time after spatial clustering in the sliding window, it generated a timeframe.

E. Temporal Feature Extraction Based on Wavelet Entropy
After spatial clustering, the result of functional connectivity analysis is no longer presented in a manner of the system's average behavior but optimized local information (node cluster). As the dynamical EEG activities on the extracted nodes still contain voluminous information from mentality and physicality, fatigue-related information is often implied. Normally, in order to differentiate between different levels of fatigue, a set of nonlinear EEG features need to be extracted. Among the nonlinear analysis approaches of the dynamical EEG activities, entropy-based algorithms have been useful and robust estimators for evaluating regularity or predictability [30], [31], [37]. The Shannon entropy (SE) is a measure of average uncertainty of event set. The information theory reflects the information content and complexity of the source (event set). In this article, we define the set of the EEG activities in a time window as the event set. The alpha subband after preprocessing [ Fig. 2(c)] was employed to calculate the SE values. The same sliding window as spatial clustering was used for dynamically estimating the fatigue of drivers. The equations are as follows: where V i is the amplitude of the ith sampling point. If V i > 0, V lim is defined to be the upper limit of the measurement range of the instrument (the upper limit > 0). If V i < 0, V lim is defined to be the lower limit (the lower limit < 0). If V i = 0, then E j is defined to be 0. L indicates the length of the sliding window (L = 30 × sampling length). m denotes the shift factor. N represents the number of samples. The EEG signals on the nodes gained from the spatial clustering (in the following section) were selected to extract the fatigue-related features. The alpha waves were among the earliest described functional oscillatory components in the human EEG [38], and research has supported the general notion that alpha band power is inversely related to brain activation [39], [40] and reflects deactivated or inhibited cortical processes [41]. Thus, the alpha waves were used to calculate the entropy features. The window moves forward along the experimental data. Every time when the window moves a step, we first decompose the data of the window into four levels by wavelet transform to get the alpha subband and take advantage of the wavelet-based threshold technique (1) to suppress artifacts. Then, we calculate the entropy using (11).
As the window moves to the last part of the data, we can get the entropy versus time curve [ Fig. 7(a)].

F. Spatiotemporal Imaging
A T1-weighted structural MRI (fieldtrip template) from the scanner in DICOM format was manually coregistered to the EEG coordinate system (10-20 system) through the position of the template fiducials and EEG channel positions. The MRI was then segmented to obtain a representation of the brain, including gray and white matter, and cerebrospinal fluid. A single-shell model in [42] was used to construct a volume conduction model. The brain regions (voxels) associated with mental fatigue provided in the previous fMRI studies [1], [2] were marked with blue points in the structural MRI template in Fig. 3(a).
When we focus on the specific node clusters on brain networks, spatial clustering actually plays a role in information extraction in the spatial domain for adapting to fMRI registration and dynamic feature extraction for prolonged task. To link the EEG features to the fMRI images, we consider the nodes extracted by spatial clustering [red points in Fig. 3(a)] as fiducial points to search the nearest coordinates of brain regions (voxels) associated with mental fatigue provided in the previous fMRI studies [1], [2]. These studies reported the coordinates of sources. In the same MRI coordinate system, the sources can be registered to the template MRI in this article. The coordinates from the previous fMRI results which were nearest to the coordinates of the extracted nodes were considered as the source coordinates of mental fatigue [the nodes with two colors in Fig. 3(a)]. The Mahalanobis distance was calculated to measure the coordinate closeness. Mahalanobis distance represents the covariance distance of the data. It is an effective method to calculate the similarity between two unknown samples sets [43], [44]. Fig. 3(c) shows the Mahalanobis distances between the three extracted nodes and fMRI alternative sources. The alternative sources with shortest Mahalanobis distances were determined as the sources J(t) in (2). All the alternative sources came from previous studies of mental fatigue. Thus, the coordinates of the sources J(t) were considered as the source coordinates of mental fatigue. Through the spatiotemporal combination (Fig. 3), we have not only spatial information reflected by the fMRI sources (specific clusters) but also temporal information reflected by the dynamic EEG features (wavelet entropy). The feature data are represented as a tensor (i.e., sequence of space-time matrices), whereas 2-D images (feature dynamics after spatial clustering) evolve over subject, as shown in Fig. 3(b). If let I be the tensor of the entropy-based features [WE in Fig. 3(b)], the element I(i, j, k) corresponds to the wavelet entropy value of subject k at the jth key spatial node and in the ith sliding window. According to the deactivation level [ Fig. 8(c)], the critical points of mental fatigue can be obtained by the edge detection of spatial values of the measured feature field at t − y, I(t, y), at each instant s [ Fig. 8(d)].

A. Functional Connectivity and Spatial Clustering
As mentioned above, the analysis in a sliding window generated timeframes. Each timeframe was marked as the time corresponding to the first data point in the sliding window. Fig. 4(a) shows the functional connectivity matrices in different stages (T1, T2, and T3). T1, T2, and T3 correspond to timeframes 0, 30, and 60 min, respectively. T1 indicates the initial state, which is considered as an alert state according to the experimental setup. The time of self-reported subjective fatigue from the subjects is 26 ± 9 min (mean value ± standard deviation). Thus, most subjects already began to feel tired at stage T2 (30 min). T3 (60 min) was selected to present the result of the final state of the subjects in the prolonged task. Each small square surface in Fig. 4(a) represents a connection between two nodes. The color indicates the value of the correlation coefficient. The correlation coefficients were set to zero (deep blue) if the nodes were self-connected. Since fatigue reduced the activity complexity of the whole brain, the overall functional connectivity in the final state became weaker comparing with the initial state. To measure the overall functional connectivity of each brain network in the three stages, the energy of all the elements of the functional connectivity matrix A was calculated and the results of statistical analysis are shown in Fig. 4(b). The global connectivity energy averaged across the subjects significantly decreased at the stage T3. Through statistical unpaired two-sample t-test, p < 0.05 between stages T1 and T3 revealed statistically significant differences in connectivity with a confidence interval of 95%. Nevertheless, as mentioned in Section I, the functional connectivity analysis at the sensor level involves the problem of signal mixing which may bring the fluctuations to the networks' average behavior. As shown in Fig. 4(b), the average connectivity energy at stage T2 was abnormally higher than that at the stage T1 and there was no statistically significant difference between the two stages. It assuredly brought challenges for interpreting and reporting the connectivity results, and it was ambiguous whether fatigue impacted the functional connectivity at the stage T2 according to the traditional connectivity analysis. As most subjects already began to feel tired at the stage T2 in the experiments [see Fig. 8(d)], it probably left some traces in the objective brain signals. Obviously, the fatigue information was "implied" for the traditional connectivity analysis.
To improve the information accuracy from the space domain, we extracted the significant nodes of the brain networks using spatial clustering. The significant nodes will be adopted as the fiducial points in the following spatiotemporal imaging.
Through spatial clustering based on the connectivity, the nodes of the brain networks in Fig. 5(a) were classified into two groups by using (5)-(9) and nodes Fz, C3, and P3 (the cyan nodes) were distinguished out as the sparse spatial sources whose synchronized deactivation information was retained [ Fig. 7(b) and (d)]. The dendrograms in Fig. 5(b) demonstrated the clustering processes of the brain networks. Two clusters were merged at each level. Hence, at the penultimate level corresponding to the position of the cutting line, one cluster contained only one EEG node and the other cluster contained the remaining 18 nodes. The cluster corresponding to the single node at the penultimate level was the most different one from others, as the other 18 clusters are merged. Such multiple merging or dendrogram structure change was probably caused by internal source activity change since signal mixing does not vary over time [22]. Therefore, at each stage of spatial clustering, the single node is considered as the pivotal node which is modulated by mental fatigue in the prolonged task. As shown in Fig. 5(a), the clusters with red links and cyan nodes had spatial gradient changes in the three stages not subject to the problem of spurious interactions in Fig. 4. In addition, the temporal features on the three extracted nodes also significantly improved comparing with all the network nodes and every single node (see next section).

B. Comparison With Existing Methods
We compared the performance of our proposed method with that of the commonly used method, connectivity estimation in EEG source space (see technical details in [32] and [45]). The EEG recorded on electrodes was projected back to the cortex level using the two source reconstruction methods, dynamic imaging of coherent sources (DICS) and minimum norm estimation (MNE). In the source space, connectivity was calculated by the Granger causality and phase lag index (PLI). The networks were constructed by the connectivity and the cortex nodes from the same structural MRI template used in the proposed method. The cortex nodes with the maximum number of links were determined as the source nodes. Fig. 6 shows the extracted source nodes in different stages (T1, T2, and T3). The link number in the Granger and PLI networks was normalized by dividing by maximum link number and converted to relative link number. From Fig. 6(a) to (c), the source nodes (red nodes) obtained by all the three methods under the DICS source reconstruction condition have a tendency to move to the hindbrain as driving time increases. However, when the source reconstruction condition was changed to MNE, the extracted source nodes by the connectivity methods based on Granger causality and PLI did not maintain this trend. Since the proposed method was combined with the fMRI neuroimages, it was unaffected by the EEG source reconstruction and the extracted source nodes correspond to the activity on the scalp electrode.

C. Temporal Feature Extraction Based on Wavelet Entropy
Even though the global topology in the brain networks demonstrated the degraded performance of the connectivity with the accumulation of fatigue [T3 stage in Fig. 4(a) and (b)], it was hard to know when it started just using the traditional connectivity analysis due to the lack of dynamic variables and the problem of signal mixing. Fig. 7(a) presents the dynamic entropy features of all the EEG channels. The red dashed line shows the time when the subject reports the feeling of fatigue (subjective fatigue) in the experiment. In Fig. 7(a), the features from different spatial loci of EEG recordings had a certain synchronous trend and the inflection point (turning-point) seemed to correspond to the dashed line of subjective fatigue. Nonetheless, the spatial residual fluctuations in different channels blur the inflection point position which was an important labeling information of fatigue. In Fig. 7(b), the temporal features from the significant nodes (fiducial points) obtained by spatial clustering had better synchronized changes. The synchronization among all the nodes and extracted nodes was estimated by fuzzy synchronization likelihood (FSL, see details in [46]) and presented in Fig. 7(c). Note that FSL was calculated in a moving window with a length of 5 min without overlap. Therefore, Fig. 7(c) did not include the synchronization of the features in the last 5 min of Fig. 7(b). The most points on the FSL curve of the extracted nodes were above the curve of all the nodes, especially after the subjective fatigue time, which demonstrated the effectiveness of the proposed method. As a measure of average uncertainty of event set [47], the entropy values decreasing revealed complexity decreasing of brain activities [ Fig. 7(a)]. The whole brain was less active. Consequently, the subjects were in the fatigue stage after the critical points. This was confirmed by the temporal records of subjective feelings in the experiment [ Figs. 7(a) and 8(d)].
The synchronous deactivations corresponding to less alternative brain activities means only the most important function areas will stay relatively active after the subjects are fatigued.
It is important to note that the residual fluctuations exist not only in the space domain but also in the time domain. Although the usefulness and robustness of entropy-based algorithms have been proven by the previous studies [30], [31], [37], the entropy features extracted from a single channel of EEG may still be unreliable and hard to be interpreted sometimes. In Fig. 7(b), the entropy value of C3 is larger than that of Fz from 0 to 25 min. However, in Fig. 7(d), the mean entropy value of Fz is larger than that of C3 at 0 min. This means the entropy values of Fz of most subjects except the subject corresponding to Fig. 7(b) are larger than those of C3 at 0 min. The entropy feature of C3 from 0 to 25 min has a fluctuation. Additionally, in Fig. 7(d), the mean entropy value of Fz at 60 min is larger than that at 30 min and even close to the mean entropy value at 0 min, which seems anomalous comparing with the mean entropy values of C3 and P3 with the fatigue accumulation. To avoid multiple fluctuations, the results of spatial clustering were added to the temporal features at the three stages. The mean entropy values of the extracted nodes were marked with red dots. As shown in Fig. 7(d), the mean entropy values with red dots had temporal gradient changes in the three stages not subject to the problem of temporal fluctuations. Therefore, the proposed method is also useful for the statistical analysis based on the mean entropy values to show the effect of time on the task of mental fatigue.

D. Spatiotemporal Imaging
For further source analysis in the space domain, we linked our EEG results to fMRI images. According to previous fMRIbased studies of mental fatigue [1], [2], the main regions significantly correlated with mental fatigue included: anterior cingulate [Brodmann area (BA) 24/32], middle cingulated (BA 23), left inferior frontal (BA 44), temporal (BA 21 and 22), postcentral (BA 3 and 4), and parietal (BA 7) regions. We registered between the fMRI and the EEG features based on three fiducial points (Fz, C3, and P3) to determine the space coordinates. As these sources were sparse enough and had a certain distance from each other, the unrelated ones with our experiment can be rejected easily [ Fig. 3(a)]. The coupling foci in this article were located at the anterior cingulate cortex (ACC: BA 24 and 32), postcentral gyrus (PCG: BA 3 and 4, left hemisphere), and posterior parietal cortex (PPC: BA 7, left hemisphere) regions [see Fig. 8(a)]. After the correction and projection of EEG features, the fMRI images had temporal dynamics and showed strong synchronized deactivations in a specific period of time [ Fig. 8(b)].
To further visualize the fatigue patterns, we converted the temporal features on the spatial regions into spatiotemporal images [ Fig. 8(c)]. The feature values were normalized and represented by different colors. The color gradients on the large time scale revealed the degradations of brain activity and retained the synchronization in certain periods of time. The blocks containing wide range of pixels in the similar colors are observable. With larger gradients, the more definite critical points in time were obtained, distinguishing the strong synchronized deactivations. As mentioned above, FSL was employed to quantify the intensity of synchronized deactivations. Fig. 7(c) corresponded to Fig. 8(c). As shown in Fig. 8(c), the strongest deactivation can be observed in a rough range of 30-50 min. In Fig. 7(c), the range was determined as 28.5-45.5 min. We considered the initial time [e.g., 28.5 min in Fig. 7(c)] of the range of the strongest deactivation as the physiological fatigue time of the subjects. The fatigue time of all the subjects was represented by the red asterisks presented in Fig. 8(d). Since the fatigue time was obtained by the spatiotemporal dynamical analysis of objective EEG signals, the state of the subjects was defined as the state of objective fatigue when the strongest synchronized deactivation occurred. As shown in Fig. 8(d), the objective fatigue timing points of all the subjects were corresponding to the recording points of subjective fatigue.

IV. DISCUSSION
As we know, the effects of fatigue are multifaceted [5]. For our brains, it weakens functional coupling among different areas and brings synchronously inhibitory performance of electrical activities. The results of the present study confirm the link between the proposed spatiotemporal images and mental fatigue.
As shown in Fig. 5, the clustering algorithm distinguished the special nodes with distinct connectivity attributes throughout the networks. If we project them to the time dimension, the clustered sources indicating enriched alpha activities move from the frontal region to the parietal region in the fatigue process. Mental fatigue brings variations of space distribution features. The results of the present study are consistent with the existing evidence on cognitive science. The rhythmic alpha activity increases in posterior brain regions (parietaloccipital) during attentional lapses [48], [49] and during states of drowsiness relative to states of alertness [28], [38], [50].
In the time domain, the critical points between different stages definitely have practical significance, for they can become triggers of feedback warnings. As shown in Fig. 7(c), the temporal features on the extracted nodes had better synchronized changes, especially when strong synchronized deactivation occurred in Fig. 8(c). This also proves that the compressed spatial information by clustering is effective.
In the space domain, the three critical source areas further localized by fMRI play key roles in brain circuitry that include essentially cognitive and executive functions [1], [51], [52]. The anterior source of BA 24 and 32 is located in ACC which is involved in higher cognitive functions, such as working memory, selective attention, error detection, and controlled information processing [8], [51], [53], [54]. The primary somatosensory and motor cortices (also known as M1 and S1, respectively) covering PCG (BA 3 and 4) contribute to the somatosensory input and the motor control adapting the scene changes [55], [56]. PPC (BA 7) lying between visual areas and M1 is believed to play a role in transforming the spatial information from the visual system into a motor plan [57]- [59]. All the functional activations centered on the three sources exhibited phased declines linked to the subjective feedbacks of mental fatigue. It illustrated the multiple negative effects of mental fatigue over time in the prolonged task. Centered on these essential functions, mental fatigue may have a powerful impact upon the external performances. Since the goal-directed attention was necessarily deployed with every movement in the driving task, an ACC involved in the brain's error checking and the conflictmonitoring system would trigger compensatory adjustments to avoid accidents [60], [61]. This explained why the anterior source was active again at the later stage in Fig. 8(c). In line with the previous mental fatigue analysis, the nodes detected in Fig. 5 show hemispheric asymmetry [13], [62]. The previous study [13] reported that the asymmetrical pattern of connectivity in frontoparietal regions associated with sustained attention. Using the spatial clustering algorithm, we found the sets of connections with nodes (cyan nodes in Fig. 5) which were the most different from others were distributed in the left hemisphere. In this article, all the subjects were right handed. As is well known, each hemisphere of the brain controls the opposite side of the body. Hence, the left lateralization of the extracted nodes involved in the motor, somatosensory, and posterior regions may have a relationship with the excessive use of attention and body control on the right to complete the prolonged driving task [55].
As shown in Fig. 1, the core of the proposed neuroimaging method in this article is spatial clustering and we have added dynamical variables (i.e., entropy features) to it. Note that a single link of EEG may be a spurious link and the networks' average behavior is also unstable for the fatigue investigation (see Fig. 4) because of the problem of EEG signal mixing. The proposed method considers the connections on a node together as a cluster. Even if the signals are mixed, it can be resolved in such a way that when we analyze the connectivity considering the relationship between the signals from every single sensor and the other sensors. Additionally, since signal mixing does not vary over time [22], the multiple merging and dendrogram structure changes of spatial clustering is probably caused by internal source activity change. The results (Fig. 5) of this article have proved that the clusters involve the relationship between the signals from every single sensor and the other sensors can be classified into groups by spatial clustering over time with the progress of fatigue and acquire more precise determination of local area modulated by fatigue. Thus, the proposed method makes it possible to combine the EEG results with the fMRI neuroimages to analyze the source changes related to fatigue [Figs. 3(a) and 8(b)].
As shown in Fig. 6, the proposed method is unaffected by the EEG source reconstruction. Because we did not solve (2) directly, it avoided the ill-posed problem. As mentioned in Section II-C, the solution of (2) for the EEG inverse problem is nonunique and unstable. The accuracy of the source reconstruction is affected by a diverse range of factors, including the head model errors, the source-modeling errors, and EEG noise (instrumental or biological) [32]. Therefore, the results of Granger causality and PLI are inconsistent under the two conditions of source reconstruction. In this article, the sources are determined by registration back to fMRI source space after spatial clustering. The alternative sources are from related fMRI studies. Obviously, some errors come from the similarity of the EEG and fMRI results, because the data of the two modalities were collected from different subjects and even different tasks. The different task (i.e., prolonged driving task) was selected to exclude the possibility that the consistent performance came from the same task. The dominant factor was the same (i.e., fatigue), so there was a corresponding relationship between the results of the two modalities. We used the Mahalanobis distance to measure the similarity. In Fig. 3(c), the sources with shortest distances were determined as the sources correlated with mental fatigue. It has a more consistent performance compared with other existing methods in Fig. 6. In addition, the synchronous deactivations from the extracted sources correspond to subjective fatigue in time domain [see Figs. 7(b) and 8(c)]. It is important to note that there may be coordinates from related fMRI studies that are not included in Fig. 3(a). The shortest Mahalanobis distances in Fig. 3(c) may become smaller. We do not focus on proving the accuracy of these sources from the fMRI studies. In this article, we just provide an interface method which can combine the EEG and fMRI results and determine which sources inside the brain are more relevant to mental fatigue.
The limitation of the current spatial clustering is that it calculates the distance only between the extracted nodes and the alternative sources. Only one source can be selected for the multiple sources under an extracted node. In the future, we will calculate the distance considering the connections on the extracted nodes to improve the algorithm.
Additionally, deep-learning-based studies have impressive performance in many areas [63]- [65]. Of cause, it also has the potential to investigate mental fatigue. The clustering algorithm in spatial clustering is a kind of machine learning algorithm. Another future direction is to replace the clustering algorithm by a deep learning model to improve the analysis performance of mental fatigue.
In this article, the sliding window length is selected as 30 s so that the number of samples is enough to estimate the connectivity and extract the temporal features. This parameter of sliding windows may affect the final result, but there is no standard for the length selection. Fortunately, the consistent performance is obtained by spatial clustering in Fig. 6, which can point to mental fatigue.
Because of residual temporal fluctuations, the temporal features from a single channel of EEG are also not reliable sometimes, even for the nonlinear dynamic features like entropy. In Fig. 7, the temporal fluctuations bring difficulties to interpret the results of statistical analysis, or the tradeoff between individual differences and averages is difficult. Hence, the results of spatial clustering (extracted nodes) are used as the fiducial points to suppress the problem of temporal fluctuations. As shown in Fig. 7(b) and (d), the fluctuations can be excluded in some extent, when we just focus on the periodic statistical analysis of the temporal features based on the fiducial points (the mean entropy values with red dots) and the synchronization of the temporal features of the fiducial points. The mean entropy values with red dots have temporal gradient changes in the three stages corresponding to the process of increasing fatigue.
Note that we do not reduce the residual fluctuations of the temporal features in every channel in the current study. The proposed method emphasizes the directive role of the effective information of mental fatigue. It is particularly helpful to add correct dynamical variables to the neuroimages after spatial clustering is applied. As shown in Fig. 7(c), we have quantified the synchronized deactivation from the spatiotemporal neuroimages in Fig. 8(c). The synchronized deactivation presents phase changes and the time range of the strong synchronized deactivation is obtained. The initial time of the strongest synchronized deactivation corresponds to the appearance of subjective fatigue reported by the subjects themselves [see Figs. 7(c) and 8(d)]. Therefore, the spatiotemporal dynamical analysis may help to find the source of subjective feelings of mental fatigue or in which way it is truly modulated. It may not be one or several spatial loci, but a pathway in brain.
Concurring with previous studies focusing specifically on each functional subsystem [8], [51]- [61], [66], our results indicated that ACC, PCG, and PPC served as the pivotal nodes engaged in the fatigue modulation. Through community structure analysis, the fMRI study [66] reported that frontoparietal network compensated with cognitive decline due to mental fatigue. That may explain the strong synchronized deactivation of the frontoparietal regions in this article. What is more, from ACC to PCG, other research also proved that rich anatomical connectivity with the association and motor cortices supported the operation of the higher cognitive subsystems centered on ACC [67], [68]. From the point of view of neuropharmacology, multiple neurotransmitters, e.g., dopamine and adenosine, are speculated to convey the connectivity of the pathway, which shows that regulating motor control and monitoring functions of ACC rely on dopaminergic inputs from the midbrain [9], [69]. At the same time, the accumulation of adenosine (an inhibitory neurotransmitter) has also been found in the motor area and ACC (due to a mentally fatiguing task) and explains the endurance performance reduction [4]. Similarly, the link between ACC and PCG is also indispensable to support functional implementation from visuomotor coordination to intended movement [57], [69]. Coupled with the electrophysiological connectivity in this article, the negative effects of mental fatigue on the pathway (ACC-PCG-PPC) may be attributed to the microstructure and multidimensional property.
Overall, the strong deactivations on nodes ACC, PCG, and PPC modulated by mental fatigue exhibit synchronization and gradient from temporal and spatial dimensions, respectively. By means of the pathway, mental fatigue reduces the complexity and diversity of brain activities, which reflects the multifaceted effect. In accordance with the order of functions to achieve, fatigue first may block sensory information receiving or sensory stimulus interpreting. In the halfway, it affects the motor function. At last, the overshoot of higher order cognitive responses will be increased and there may be more error responses. With the spatiotemporal attribute structure, such a pathway provides evidence of multidimensional property that should be accommodated by any mechanism model of mental fatigue.

V. CONCLUSION
In this article, we have presented a spatiotemporal mapping method based on spatial clustering for imaging mental fatigue with the advantages of EEG and fMRI. It demonstrates the consistent performance related to mental fatigue across different subjects and different tasks. During the fatigue process, the clustered sources (ACC, PCG, and PPC) show strong deactivations corresponding to subjective feelings of fatigue. The deactivations are interrelated to each other on account of synchronization and gradient from temporal and spatial dimensions, respectively. The multidimensional property of mental fatigue can be expressed as the combined effect of ACC-PCG-PPC to a certain degree.