Enhancing Classification Performance of Functional Near-Infrared Spectroscopy- Brain–Computer Interface Using Adaptive Estimation of General Linear Model Coefficients

In this paper, a novel methodology for enhanced classification of functional near-infrared spectroscopy (fNIRS) signals utilizable in a two-class [motor imagery (MI) and rest; mental rotation (MR) and rest] brain–computer interface (BCI) is presented. First, fNIRS signals corresponding to MI and MR are acquired from the motor and prefrontal cortex, respectively, afterward, filtered to remove physiological noises. Then, the signals are modeled using the general linear model, the coefficients of which are adaptively estimated using the least squares technique. Subsequently, multiple feature combinations of estimated coefficients were used for classification. The best classification accuracies achieved for five subjects, for MI versus rest are 79.5, 83.7, 82.6, 81.4, and 84.1% whereas those for MR versus rest are 85.5, 85.2, 87.8, 83.7, and 84.8%, respectively, using support vector machine. These results are compared with the best classification accuracies obtained using the conventional hemodynamic response. By means of the proposed methodology, the average classification accuracy obtained was significantly higher (p < 0.05). These results serve to demonstrate the feasibility of developing a high-classification-performance fNIRS-BCI.


INTRODUCTION
A brain-computer interface (BCI) system bypasses the peripheral nervous system and provides means of communication for patients suffering from motor disabilities or in a persistent vegetative state using devices, such as robotic arms or other prostheses (Wolpaw et al., 2002). The brain signals are acquired either invasively or non-invasively. Although the quality of brain signals acquired using invasive methods is better than those using non-invasive methods, their acquisition entails extensive surgical risk (Wester et al., 2009). With non-invasive methods, on the other hand, there is no such risk. Non-invasive techniques include electroencephalography (EEG) (Wolpaw et al., 2002;Pfurtscheller et al., 2003;Salvaris and Sepulveda, 2010;Cong et al., 2011Cong et al., , 2015Jin et al., 2011Jin et al., , 2014Jin et al., , 2015Choi, 2013;Chen et al., 2015), functional magnetic resonance imaging (fMRI) (Enzinger et al., 2008;Sorger et al., 2009), and functional near-infrared spectroscopy (fNIRS) (Ferrari et al., 1985;Kato et al., 1993;Coyle et al., 2004Coyle et al., , 2007Naito et al., 2007;Naseer and Hong, 2013;Naseer et al., 2014;Noori et al., 2017). Over the course of the past decade, fNIRS-based BCI systems have been the focus of considerable research interest and discussion due to their portability, affordable cost and better temporal resolution relative to fMRI. Moreover, compared with the EEG system, they offer better spatial resolution and a superior signal-to-noise ratio (Hu et al., 2012;. In general, fNIRS has evolved into a neuroimaging technique that has contributed to ground-breaking advances in the understanding of human brain functionality (Irani et al., 2007;Aqil et al., 2012b;Ferrari and Quaresima, 2012;Hong and Nguyen, 2014;Hong and Naseer, 2016;Hong and Santosa, 2016). fNIRS utilizes near-infrared (NI) light within the 650-1000 nm wavelength range to measure changes in the concentrations of oxygenated and deoxygenated hemoglobin [Δc HbO (t) and Δc HbR (t)] according to the modified Beer-Lamberts Law (Delpy et al., 1988;Villringer et al., 1993;Hoshi et al., 1994;Hoshi and Tamura, 1997). Since the introduction of the principle of NI spectroscopy by Jobsis (1977), fNIRS has been used effectively for functional and structural brain imaging as well as for BCI purposes (Naseer and Hong, 2015a,b;Nguyen et al., 2016;Zafar and Hong, 2017). The first step in fNIRS-BCI is to acquire signals from a suitable mental task. Over the past decade, the mental tasks used by fNIRS-BCI researchers have been motor imagery (MI), mental rotation (MR), mental arithmetic, music imagery, and letter padding (Zhang et al., 2011a;Ayaz et al., 2013Ayaz et al., , 2014Khan et al., 2014;Khan and Hong, 2015). In this study, we used right-hand MI and MR as the brain activity. Generation of control commands for fNIRS-based BCI systems proceeds according to the following conventional steps: first, acquisition of the desired signals; second, removal of motion artifacts and physiological noises; third, extraction of significant information (features), usually from the hemodynamic signals' physical properties; fourth and finally, classification of the extracted features preparatory to generation of the desired control commands. Researchers have devoted considerable efforts to the improvement and enhancement of classification accuracies for fNIRS-BCI, specifically by use of different features and classifiers (Ayaz et al., 2013(Ayaz et al., , 2014Naseer and Hong, 2013;Naseer et al., 2014;Noori et al., 2016;Qureshi et al., 2016;Khan and Hong, 2017). In this paper, we propose that features be extracted from the estimated coefficients of the general linear model (GLM).
The GLM methodology was first employed by Abdelnour and Huppert (2009) in a fNIRS-based BCI study. Since that time, multiple GLM-based fNIRS studies have been performed for noise removal and brain mapping (Hu et al., 2010;Zhang et al., 2011bZhang et al., , 2012Aqil et al., 2012a;Kamran and Hong, 2013). Abdelnour and Huppert (2009) have proposed the use of filter coefficients obtained by Kalman filtering as the features for classification. They assumed that different brain activities will produce different filter coefficients, using which different signals can be classified. Similarly, recursive least square estimation (Aqil et al., 2012a), and wavelet transform (Khoa and Nakagawa, 2008;Abibullaev et al., 2011;Abibullaev and An, 2012) have also been used for brain mapping using GLM. In this study, GLM is used with least square to estimate filter coefficient values. Afterward these values are used to extract features. To the best of our knowledge, this is the first work that uses filter coefficient values to extract statistical features that can be used for classification. In the proposed methodology, signals are acquired from the left motor cortex of the brain for right-hand MI (clenching of the right hand) and rest tasks, whereas MR (rotation of rectangular box) and rest signals are acquired from the prefrontal cortex; these signals are filtered to remove physiological noises and the GLM coefficients are extracted using the least squares estimation (LSE) technique; the feature values of these coefficients are then fed to support vector machine (SVM) for classification. The motivation of using GLM-based features for fNIRS data came from Abdelnour and Huppert (2009). They showed promising results using beta (β) values extracted from GLM as features. In this study, authors have used GLM with least square to estimate β values. Afterward β values are used to extract features in order to calculate classification accuracies.

Experimental Procedure
Subjects A total of 10 subjects participated in the experiments. Five subjects performed MI (right-hand clenching) versus rest, whereas the other five performed MR (rotation of rectangular box) versus rest. The reason for introducing two different experiments was to establish generalization of the proposed methodology. The subjects were each seated on a comfortable chair in front of a display screen and asked to restrict their body movements as much as possible during the experiment. Verbal consent was obtained from all of the subjects after explaining the experimental paradigm in detail. The subjects had little or no previous experience of fNIRS recording. This work was approved by the Institutional Review Board of Pusan National University. All experiments were conducted in accordance with the ethical standards encoded in the latest Declaration of Helsinki. The complete details of the baseline system (conventional methodology) can be found in Naseer et al. (2016a,b).

Motor Imagery
The first 20 s was the rest period, required in order to set up the baseline condition; it was followed by 20 s of a right-hand MI task (clenching of the right hand), followed by another 20 s rest period that allowed the signals to return to their baseline values before the start of the next trial. This pattern was repeated 11 times; the total duration of experiment for each subject, therefore, was 440 s. During the MI task, the subjects were asked to imagine clenching of their right hand with a self-paced frequency of around 1 Hz; during the rest period, they were asked to relax.

Mental Rotation
Similar to MI task in MR task the first 20 s was the rest period, required in order to set up the baseline condition; it was followed by 10 s of object rotation (rotation of rectangular box), followed by another 20 s rest period that allowed the signals to return to their baseline values before the start of the next trial. This pattern was repeated 10 times; the total duration of experiment for each subject, therefore, was 330 s. During the MR task, the subjects were asked to imagine rotation of a rectangular box; during the rest period, they were asked to relax.

Signal Acquisition
In order to acquire fNIRS signals from the right motor cortex of the brain, a multi-channel continuous-wave imaging system (dynamic near-infrared optical tomography; two wavelengths: 760 and 830 nm; NIRx Medical Technologies, NY, USA) with a sampling rate of 1.81 Hz was employed. The acquired signals' light intensities were first converted to Δc HbO (t) and Δc HbR (t) using the modified Beer-Lamberts Law (1) where ΔA(t; λ j ) (j = 1, 2) is the unit-less absorbance (optical density) variation of a light emitter of wavelength λ j , a HbX (λ j ) is the extinction coefficient of HbX (HbO and HbR) in μM −1 mm −1 , d is the unit-less differential path length factor, and l is the distance (in millimeters) between the emitter and detector. As shown in Figure 1A, four emitters and five detectors were positioned over the left motor cortex of the brain for right-hand MI task. Figure 1B shows eight emitters and three detectors placed on prefrontal cortex of the brain region in order to acquire signals for MR, the distance between each emitter-detector pair was of 3 cm. This emitter-detector distance is in accordance with the literature (McCormick et al., 1992;Gratton et al., 2006). In order to remove physiological noises (heartbeat, respiration) from the obtained signals, the Butterworth filter of order four was used with a cut-off frequency of 0.6 Hz; for removal of low-oscillation Mayer waves, a high-pass filter with a cut-off frequency of 0.01 Hz was used (Naseer and Hong, 2015a,b).

General Linear Model
The GLM has been very widely utilized by researchers of fNIRS-BCI systems in order to identify brain-activation patterns for multiple cognitive tasks (Abdelnour and Huppert, 2009;Hu et al., 2010;Zhang et al., 2011bZhang et al., , 2012Aqil et al., 2012a;Kamran and Hong, 2013). The GLM-based methods were developed initially for fMRI-based functional brain mapping. To analyze fMRI data, GLM methodology has been developed to explain the timeline blood oxygenation level dependent signal. Currently, they are frequently used in fNIRS studies. The GLM defines measured data in the form of a linear combination of several variables and an error term. The observation of hemodynamic changes can be expressed as where the y vector represents the measured data (in fNIRS, the vector is the observed time-series of the hemodynamic response), G is the design matrix obtained by convolving the canonical hemodynamic response with the experimental box-car function (Ye et al., 2009), β is the set of coefficients for the functional response that we want to estimate, and e is the error term. The vital part of the model function of a GLM is the box-car function, which reflects the temporal structures of the experimental paradigm and is convolved with the canonical hemodynamic response function ( Ye et al., 2009). As physiological noises had already been removed using the Butterworth filter, only one explanatory variable (the design matrix) was used to extract the β values.

Least Squares Estimation
Least squares estimation is used to estimate the β values from the GLM. The time-course values predicted by the model are obtained by linear combination of the predictorŝ In order to achieve a good fit, the β values should be close to the predicted values that are as close as possible to the measured values y. Thus, the system of equations should be rearranged as Although the GLM methodology does not estimate β values, it can be applied to minimize the sum of squared error values by using where e ′ e shows the vector notation for the sum of squares. Utilizing LSE, the β weights minimizing the square error values are obtained by The resulting matrix (G −1 G) −1 plays an important role in the calculation of the β values. The remaining term on the right side, G ′ y, evaluates a vector containing as many elements as predictors. Figure 2 plots the Δc HbO (t) signals for the MI task and rest period with their corresponding adaptively estimated β values for subject 5.

Feature Extraction and Classification
In this study, the statistical properties of the β values were used as the features. Signal peak (SP), signal skewness (SSk), signal mean (SM), signal variance (SV), signal kurtosis (SK), and signal slope (SS) were extracted from the β values obtained by LSE. The SSk values were determined by measuring the asymmetry of the signal values around the mean relative to a normal distribution: where N is the number of observations and Y i represents the β values. The variance is calculated as follow: where μ is the mean value of Y. The kurtosis is computed as follows: The SS is calculated using the polyfit function in MATLAB®. The SP values, which measure the peaks of signals, were determined using MATLAB®max function. These features were calculated across all 12 channels for the MI and MR. All of the feature values were scaled between 0 and 1 using the equation where x ∈ R n represents the original feature values, x ′ denotes the rescaled feature values between 0 and 1, max(x) is the largest value, and min(x) is the smallest value. After extracting the features from the β values, SVM was used to classify the MI and MR tasks (Naseer et al., 2016b). SVM maximizes the margins between classes by creating hyperplanes that minimize the cost function where w T , x i ∈ R 2 and b ∈ R 1 , ||w|| 2 = w T w, C is the trade-off parameter between the error and the margin, ξ i is the measure of the training data, and z i is the class label for the i-th sample. The most significant advantage of SVM is that it can be used as a linear as well as a non-linear classifier; in fact, in this study, a thirddegree polynomial kernel function was used with C = 0.5. Tenfold cross-validation was utilized to extract the classification accuracies for the MI and MR tasks versus rest periods. Moreover, in order to measure classification performance, recall and precision were calculated for both paradigms as follows: where TP, FP, and FN denote true positive, false positive, and false negative, respectively. These values were calculated from confusion matrix (Fawcett, 2006).

RESULTS
Multiple feature combinations were used in order to extract significant classification accuracies for proposed and conventional methodologies. The classification accuracies obtained for the five subjects using the proposed method for MI versus rest were 79.5, 83.7, 82.6, 81.4, and 84.1% using SM and SSk, whereas those for MR versus rest were and 85.5, 85.2, 87.8, 83.7, and 84.8% using SP and SSk. To establish the superiority of the proposed method over the previous methods, the classification accuracies using the conventional hemodynamic response feature also were calculated. Figures 3A,B provides a schematic of the conventional and proposed methodology for fNIRS-based BCI study. Furthermore, the classification accuracies obtained for the five subjects using the conventional method for MI versus rest were 60.4, 78.9, 70.4, 68.9, and 54.4% using SM and SP, whereas those for MR versus rest were and 66.7, 73.0, 72.2, 68.5, and 63.3% using SM and SP. Tables 1 and 2 list the classification accuracies, precisions, and recalls of all subjects using the proposed methodology and the conventional method, for all possible two-feature combinations for MI versus rest task, respectively. Tables 3 and 4 list the classification accuracies, precisions, and recalls of all subjects using the proposed methodology and the conventional method, for all possible twofeature combinations for MR versus rest task, respectively. The results show that in MI task the optimal feature combinations that yielded best classification accuracies were "SM and SSk" and "SM     and SP" for beta values and conventional hemodynamic response, respectively. In MR task, the optimal feature-combination that yielded best classification accuracies were "SP and SSk" and "SM and SP" for beta values and conventional hemodynamic response, respectively. In order to ensure that the data are normally distributed Kolmogorov-Smirnov method was applied, the significant value was found to be greater than 0.05 which shows normal distribution of the data. These high classification accuracies of the proposed method relative to the conventional method were statistically verified by a statistical significance test (the Student's t-test): the p-values obtained by performing t-test on the subjectwise accuracy scores was less than 0.05, which confirmed the statistical significance of the proposed methodology's superior performance for both tasks.

DISCUSSION
In previous studies, researchers have focused their efforts on enhancing the classification performance of multiple mental tasks in order to generate commands effective for control of external devices or for communication with patients suffering from amyotrophic lateral sclerosis, locked in syndrome, or other physical disabilities. However, distinct BCI signals for a specific mental task were unsuitable for classification, even when using current advanced methods. Previously, Tai and Chau (2009); Khan and Hong (2015); Naseer et al. (2016a);Naseer et al. (2016b) have used features extracted directly from hemodynamic response in order to acquire classification accuracies. In this study, a novel methodology that proceeds by adaptive estimation of GLM coefficients and extraction of the classification performances of MI versus rest and MR versus rest task were developed and evaluated. The results indicated enhanced classification performance as compared with a conventional hemodynamic-response-based fNIRS-BCI. Moreover, the proposed methodology can enhance classification performance if a user is not able to generate distinct brain signals for a specific mental task. The GLM methodology has been frequently employed to analyze time-series fMRI data: Abdelnour and Huppert (2009) first used the GLM in an fNIRS study in order to minimize physiological noises; soon thereafter, Hu et al. (2010) developed a novel online data analysis scheme using the GLM and Kalman estimator to reduce physiological noises for finger-tapping experiments; Aqil et al. (2012a) presented an online brain-imaging framework for finger-tapping tasks using GLM and a recursive least squares estimation method; Zhang et al. (2011bZhang et al. ( , 2012 tested multiple recursive algorithms for removal of physiological noises and, thereby, extraction of better neuronrelated concentration changes in observed fNIRS data. All of these studies used the GLM for the removal of physiological noises and demonstrated brain-activation mapping for multiple cognitive tasks. However, the GLM coefficients, as estimated using LSE, have not been used as features for classification. The difference in classification accuracies is possible since in the conventional method we use statistical features obtained directly from HbO signals; whereas in the proposed method, we use statistical features obtained from β values. Multiple feature combinations have been used in order to determine optimal feature combination, which yields best classification accuracies, using proposed and conventional methodologies for both mental tasks. It was found that in MI task the optimal feature-combinations that yielded best classification accuracies were "SM and SSk" and "SM and SP" for beta values and conventional hemodynamic response, respectively. In MR task, the optimal feature-combination that yielded best classification accuracy were "SP and SSk" and "SM and SP" for beta values and conventional hemodynamic response, respectively. The proposed method has shown improved overall classification accuracies as compared to conventional methodology. This study showed that there is a significant difference between the classification accuracies of the proposed and conventional methodologies: the result is improved by an average of 18.4% for MI versus rest and 16.7% for MR versus rest using the proposed method. Moreover, it was found that features extracted from proposed methodology are statistically significant from conventional methodology for both paradigms. It should be noted that in Naseer et al. (2016a), the best twoand three-feature combinations yielded accuracies of more than 90% for a seemingly very similar classification task. In this work, the accuracies obtained are in the range of 70%. These differences in the accuracies might be attributed to different recording conditions and different mental tasks. It is observed that the signal quality in mental arithmetic tasks is better as compared to MI and MR tasks. This could be attributed to user training as well. The subjects used in Naseer et al. (2016a) were regular fNIRS-based BCI users. All subjects in this paper had little or no experience of fNIRS recording/BCI training. The effect of using current method on the data from Naseer et al. (2016a) can be evaluated in future works.
This study has some limitations. The first is that only six features were used for classification. The combination of several other statistical features acquired from β values also should be utilized as features, as, thereby, classification performance could be further enhanced. The second limitation is that only SVM was used as the classifier. The positive effects of several other classifiers, however, have been seen. As shown in Naseer et al. (2016b), classification accuracies acquired using artificial neural networks (ANNs) are better than those acquired using SVM and, therefore, ANN's can be considered for classification in future studies. The third limitation is that the proposed methodology is complex as compared to conventional method since an extra step of calculating general linear model coefficients is involved. This will increase the computational cost as well. The fourth limitation of this study is the fact that only two mental tasks for each paradigm were considered, which restricts this study to a two-class BCI problem. Certainly in any case, it can be upgraded to a multi-class BCI problem in a further study.
In conclusion, we present a novel methodology for enhanced classification accuracy of two-class fNIRS-based BCI. The hemodynamic signals of five subjects were modeled using the GLM and the beta values estimated by LSE were used to extract the features for classification. The classification accuracies obtained using the proposed methodology were significantly higher than those obtained using conventional hemodynamic-response-based features. These results indeed show enhanced classification performance relative to the conventional methodology and represent a step forward in the important task of making fNIRS-based BCIs more accurate and reliable.

ETHICS STATEMENT
This work was approved by the Institutional Review Board of Pusan National University. All experiments were conducted in accordance with the ethical standards encoded in the latest Declaration of Helsinki.

AUTHOR CONTRIBUTIONS
NQ conceived this study and was involved in the experiments, data processing, and writing of the manuscript. FN was involved in the experiments and data analysis. HN and RK were involved in data analysis, rechecking of results, and revision. SS helped in revision of the manuscript. NN was involved in the writing of the manuscript and supervised the research.