Detection and classification of three-class initial dips from prefrontal cortex

In this paper, the use of initial dips using functional near-infrared spectroscopy (fNIRS) for brain-computer interface (BCI) is investigated. Features and window sizes for detecting initial dips are also discussed. Three mental tasks including mental arithmetic, mental counting, and puzzle solving are performed in obtaining fNIRS signals from the prefrontal cortex. Vector-based phase analysis method combined with a threshold circle, as a decision criterion, are used to detect the initial dips. Eight healthy subjects participate in experiment. Linear discriminant analysis is used as a classifier. To classify initial dips, five features (signal mean, peak value, signal slope, skewness, and kurtosis) of oxy-hemoglobin (HbO) and four different window sizes (0~1, 0~1.5, 0~2, and 0~2.5 sec) are examined. It is shown that a combination of signal mean and peak value and a time period of 0~2.5 sec provide the best average classification accuracy of 57.5% for three classes. To further validate the result, three-class classification using the conventional hemodynamic response (HR) is also performed, in which two features (signal mean and signal slope) and 2~7 sec window size have yielded the average classification accuracy of 65.9%. This reveals that fNIRS-based BCI using initial dip detection can reduce the command generation time from 7 sec to 2.5 sec while the classification accuracy is a bit sacrificed from 65.9% to 57.5% for three mental tasks. Further improvement can be made by using deoxy hemoglobin signals in coping with the slow HR problem. © 2016 Optical Society of America OCIS codes: (170.2655) Functional monitoring and imaging; (300.0300) Spectroscopy; (070.5010) Pattern recognition; (200.3050) Information processing. References and links 1. R. D. Frostig, E. E. Lieke, D. Y. Ts’o, and A. Grinvald, “Cortical functional architecture and local coupling between neuronal activity and the microcirculation revealed by in vivo high-resolution optical imaging of intrinsic signals,” Proc. Natl. Acad. Sci. U.S.A. 87(16), 6082–6086 (1990). 2. A. Grinvald, R. D. Frostig, R. M. Siegel, and E. Bartfeld, “High-resolution optical imaging of functional brain architecture in the awake monkey,” Proc. Natl. Acad. Sci. U.S.A. 88(24), 11559–11563 (1991). 3. D. Malonek and A. Grinvald, “Interactions between electrical activity and cortical microcirculation revealed by imaging spectroscopy: implications for functional brain mapping,” Science 272(5261), 551–554 (1996). 4. R. S. Menon, S. Ogawa, X. Hu, J. P. Strupp, P. Anderson, and K. Uğurbil, “BOLD based functional MRI at 4 Tesla includes a capillary bed contribution: echo-planar imaging correlates with previous optical imaging using intrinsic signals,” Magn. Reson. Med. 33(3), 453–459 (1995). 5. E. Yacoub and X. Hu, “Detection of the early negative response in fMRI at 1.5 Tesla,” Magn. Reson. Med. 41(6), 1088–1092 (1999). 6. T. Q. Duong, D. S. Kim, K. Uğurbil, and S. G. Kim, “Spatiotemporal dynamics of the BOLD fMRI signals: toward mapping submillimeter cortical columns using the early negative response,” Magn. Reson. Med. 44(2), 231–242 (2000). 7. X. Hu and E. Yacoub, “The story of the initial dip in fMRI,” Neuroimage 62(2), 1103–1108 (2012). 8. R. B. Buxton, E. C. Wong, and L. R. Frank, “Dynamics of blood flow and oxygenation changes during brain activation: the balloon model,” Magn. Reson. Med. 39(6), 855–864 (1998). 9. R. B. Buxton, K. Uludağ, D. J. Dubowitz, and T. T. Liu, “Modeling the hemodynamic response to brain activation,” Neuroimage 23(Suppl 1), S220–S233 (2004). Vol. 8, No. 1 | 1 Jan 2017 | BIOMEDICAL OPTICS EXPRESS 367 #277775 Journal © 2017 http://dx.doi.org/10.1364/BOE.8.000367 Received 29 Sep 2016; revised 20 Nov 2016; accepted 12 Dec 2016; published 19 Dec 2016 10. J. B. Mandeville, J. J. Marota, C. Ayata, M. A. Moskowitz, R. M. Weisskoff, and B. R. Rosen, “MRI measurement of the temporal evolution of relative CMRO2 during rat forepaw stimulation,” Magn. Reson. Med. 42(5), 944–951 (1999). 11. G. Jasdzewski, G. Strangman, J. Wagner, K. K. Kwong, R. A. Poldrack, and D. A. Boas, “Differences in the hemodynamic response to event-related motor and visual paradigms as measured by near-infrared spectroscopy,” Neuroimage 20(1), 479–488 (2003). 12. T. Kato, “Principle and technique of NIRS imaging for human brain FORCE: fast-oxygen response in capillary event,” Proc. ISBET 1270, 85–90 (2004). 13. N. K. Logothetis, H. Guggenberger, S. Peled, and J. Pauls, “Functional imaging of the monkey brain,” Nat. Neurosci. 2(6), 555–562 (1999). 14. M. Jones, J. Berwick, D. Johnston, and J. Mayhew, “Concurrent optical imaging spectroscopy and laser-Doppler flowmetry: the relationship between blood flow, oxygenation, and volume in rodent barrel cortex,” Neuroimage 13(6), 1002–1015 (2001). 15. X. Hu, T. H. Le, and K. Uğurbil, “Evaluation of the early response in fMRI in individual subjects using short stimulus duration,” Magn. Reson. Med. 37(6), 877–884 (1997). 16. T. Ernst and J. Hennig, “Observation of a fast response in functional MR,” Magn. Reson. Med. 32(1), 146–149 (1994). 17. D. S. Kim, T. Q. Duong, and S. G. Kim, “High-resolution mapping of iso-orientation columns by fMRI,” Nat. Neurosci. 3(2), 164–169 (2000). 18. M. Watanabe, A. Bartels, J. H. Macke, Y. Murayama, and N. K. Logothetis, “Temporal jitter of the BOLD signal reveals a reliable initial dip and improved spatial resolution,” Curr. Biol. 23(21), 2146–2150 (2013). 19. M. Ferrari and V. Quaresima, “A brief review on the history of human functional near-infrared spectroscopy (fNIRS) development and fields of application,” Neuroimage 63(2), 921–935 (2012). 20. N. Naseer and K.-S. Hong, “fNIRS-based brain-computer interfaces: a review,” Front. Hum. Neurosci. 9, 3 (2015). 21. D. A. Boas, C. E. Elwell, M. Ferrari, and G. Taga, “Twenty years of functional near-infrared spectroscopy: Introduction for the special issue,” Neuroimage 85(Pt 1), 1–5 (2014). 22. M. R. Bhutta, K.-S. Hong, B.-M. Kim, M. J. Hong, Y.-H. Kim, and S.-H. Lee, “Note: three wavelengths nearinfrared spectroscopy system for compensating the light absorbance by water,” Rev. Sci. Instrum. 85(2), 026111 (2014). 23. H. Santosa, M. J. Hong, S.-P. Kim, and K.-S. Hong, “Noise reduction in functional near-infrared spectroscopy signals by independent component analysis,” Rev. Sci. Instrum. 84(7), 073106 (2013). 24. N. Liu, X. Cui, D. M. Bryant, G. H. Glover, and A. L. Reiss, “Inferring deep-brain activity from cortical activity using functional near-infrared spectroscopy,” Biomed. Opt. Express 6(3), 1074–1089 (2015). 25. T. Akiyama, T. Ohira, T. Kawase, and T. Kato, “TMS orientation for NIRS-functional motor mapping,” Brain Topogr. 19(1-2), 1–9 (2006). 26. G. R. Wylie, H. L. Graber, G. T. Voelbel, A. D. Kohl, J. DeLuca, Y. Pei, Y. Xu, and R. L. Barbour, “Using covariations in the Hb signal to detect visual activation: a near infrared spectroscopic imaging study,” Neuroimage 47(2), 473–481 (2009). 27. K. Yoshino and T. Kato, “Vector-based phase classification of initial dips during word listening using nearinfrared spectroscopy,” Neuroreport 23(16), 947–951 (2012). 28. N. Oka, K. Yoshino, K. Yamamoto, H. Takahashi, S. Li, T. Sugimachi, K. Nakano, Y. Suda, and T. Kato, “Greater activity in the frontal cortex on left curves: a vector-based fNIRS study of left and right curve driving,” PLoS One 10(5), e0127594 (2015). 29. K. Yoshino, N. Oka, K. Yamamoto, H. Takahashi, and T. Kato, “Correlation of prefrontal cortical activation with changing vehicle speeds in actual driving: a vector-based functional near-infrared spectroscopy study,” Front. Hum. Neurosci. 7, 895 (2013). 30. K.-S. Hong and N. Naseer, “Reduction of delay in detecting initial dips from functional near-infrared spectroscopy signals using vector-based phase analysis,” Int. J. Neural Syst. 26(3), 1650012 (2016). 31. F. Matthews, B. A. Pearlmutter, T. E. Ward, C. Soraghan, and C. Markham, “Hemodynamics for brain computer interfaces,” IEEE Signal Process. Mag. 25(1), 87–94 (2008). 32. J. R. Wolpaw, N. Birbaumer, D. J. McFarland, G. Pfurtscheller, and T. M. Vaughan, “Brain-computer interfaces for communication and control,” Clin. Neurophysiol. 113(6), 767–791 (2002). 33. D. J. McFarland and J. R. Wolpaw, “Brain-computer interfaces for the operation of robotic and prosthetic devices,” Adv. Comput. 79, 169–187 (2010). 34. L. F. Nicolas-Alonso and J. Gomez-Gil, “Brain computer interfaces, a review,” Sensors (Basel) 12(12), 1211– 1279 (2012). 35. A. Ortiz-Rosario and H. Adeli, “Brain-computer interface technologies: from signal to action,” Rev. Neurosci. 24(5), 537–552 (2013). 36. K.-S. Hong and H.-D. Nguyen, “State-space models of impulse hemodynamic responses over motor, somatosensory, and visual cortices,” Biomed. Opt. Express 5(6), 1778–1798 (2014). 37. S. D. Power, A. Kushki, and T. Chau, “Towards a system-paced near-infrared spectroscopy brain-computer interface: differentiating prefrontal activity due to mental arithmetic and mental singing from the no-control state,” J. Neural Eng. 8(6), 066004 (2011). Vol. 8, No. 1 | 1 Jan 2017 | BIOMEDICAL OPTICS EXPRESS 368


Introduction
The initial dip is an early small decrease of the concentration change of oxygenated hemoglobin (ΔHbO) on the locus of neuronal activity after the presentation of stimuli [1][2][3][4][5][6][7][8][9][10][11][12]. Since its first reveal [1], a number of invasive/non-invasive studies showed its occurrence prior to the increase of blood oxygenation [13][14][15][16]. The initial dip is closely related to the neuronal metabolism and is a direct indication of the generation of a neuronal command. The initial dip is found to be faster than the conventional hemodynamic response (HR) and spatially more specific [17,18]. The aim of this study is to measure the initial dips in association with three kinds of mental tasks from the prefrontal cortex, and classifies the tasks using the initial dip related features to be introduced. The three mental tasks are mental arithmetic (MA), mental counting (MC), and puzzle solving (PS).
Function near-infrared spectroscopy (fNIRS) is a non-invasive method to measure the concentration changes of oxygenated hemoglobin (ΔHbO) and deoxygenated hemoglobin (ΔHbR) by emitting near-infrared light (650~1000 nm) into the brain and measuring the photons that passed through the brain tissues by detectors placed on the scalp with appropriate distances. The principle of fNIRS has been described in detail in [19][20][21][22][23][24]. In fNIRS, to the best of our knowledge, the first study that observed the delayed response by approximately 2 sec (in term of initial dip) was done by Jasdzewski et al. [11]: The initial dip was found in the visual and motor cortex. In 2004, another study by Kato indicated that a fast oxygen response is followed by a neuronal activity [12]. Kato measured an active oxygen exchange in capillaries, which has been named to fast-oxygen response in capillary event (FORCE) related to the neuronal responses. Akiyama and coauthors incorporated transcranial magnetic stimulation (TMS) to the detection of an initial dip in the motor cortex, using NIRS imaging, to enhance the spatial specificity [25]: They revealed distinctive biphasic responses in the cortical oxygenation in the center of the primary motor cortex during the motor task, and observed an early response occurring within 1 to 3 sec after the task initiation. Wylie et al. also found the rising of HR at approximately 2 sec in the visual cortex during a contrastreversing checkerboard presentation [26]. Efforts related to a systematic detection of initial dips were result in the vector-based phase analysis method [27][28][29]. Yoshino and Kato provided an in-depth study on the classification of initial dips [27]. In their work, five types of initial dips were classified using the vector-based phase analysis: The dips were divided into canonical dips (phases I, II), hypoxic-hyperemic dips (phase III), and hypoxic-ischemic dips (phases IV, V) on the basis of phases in which the initial dips occur. However, in their study the time required for the detection of initial dip was not clear. Hong and Naseer introduced a threshold criterion on the phase diagram in detecting initial dips, and proposed the use of an auto-regressive moving average with exogenous signals (ARMAX) model to predict q-stepahead occurrence of such dips in the measured fNIRS signals [30]: The authors were able to reduce the time lag in detecting the initial dips to 0.9 sec. In fNIRS imaging for brain-computer interface (BCI), variation detection of hemodynamic responses (i.e. increase/decrease in ΔHbO/ΔHbR) is the main focus for neuronal-command interpretation [31]. BCI provides a communication link between the hardware and the brain signals by translating the brain signals into meaningful commands to actuate the hardware [32-35]. The main issue associated with using fNIRS signals for BCI applications is the inherent 2 sec time delay from the neuronal activation [11,30,36]. For this reason, researchers in fNIRS community used various features in 0~5, 2~7, 0~10, 0~15, 0~17, and 0~20 sec time windows to classify the HRs from the same or different brain regions using multi-class classification algorithms [37][38][39][40][41][42][43][44][45][46][47]. However, the use of HRs will result in timedelayed command generations due to the delays in HRs. If the detection of initial dips from the same or different brain regions is plausible, a quicker generation of commands is now possible in the BCI framework as shown in Fig. 1.
In this paper, we investigate the feasibility of classifying three different mental activities using the features associated with the initial dips measured from the prefrontal cortex. To detect the initial dips, the vector-based phase analysis like Kato et al. is adopted [27][28][29]. However, a new threshold circle with a radius defined by the peak value of either ΔHbO or ΔHbR signals during the resting state is investigated, which is in contrast with the work of Hong and Naseer in that max (ΔHbO 2 + ΔHbR 2 ) 1/2 during the resting state was used [30]. For classifying the initial dips related to the mental tasks, five features of HbO are tested: the signal mean, peak magnitude (delta), signal slope, skewness, and kurtosis of the HbO signal during 0~1, 0~1.5, 0~2, and 0~2.5 sec windows, respectively, are investigated. Linear discriminant analysis (LDA) is used as a classifier for training and testing. Revealing the obtained result briefly, the signal mean and peak value (delta) during 0~2.5 sec window provided the best three-class classification accuracies of 57.5%. Also, we have investigated the initial dip characteristics in different brain regions: prefrontal, motor, somatosensory, and visual cortices. To the best of our knowledge, this is the first study investigating the classification of initial dips in various brain regions including the prefrontal cortex. Finally, to prepare the cases when initial dips are not detected, a secondary classification loop by using the conventional HRs has been set, see Fig. 1. In this case, as features for classification, the signal mean and the signal slope during the 2~7 sec window are to be used.  when ΔHbO comes back to the baseline value after decreasing, t o denotes the time instance at which the signal reaches its negative peak value, and δ (delta) is the (negative) peak value of ΔHbO at t o . After the initial stage, t 1 denotes the time for the HR to reach its peak magnitude, a 1 is the peak value, and a 2 is the magnitude of undershoot, see the definitions in Hong and Nguyen [36]. In this study, the parameters of initial dip (t init , t o , δ) are utilized for the comparison of initial dips from the prefrontal, motor, somatosensory, and visual cortices. Also, δ is used as a feature for classification.

Subjects
A total of 8 male subjects (mean age: 27.5 ± 4.5) participated in the experiment. The subjects were MS and Ph.D. students, Pusan National University. All subjects were healthy and normal or corrected-to-normal vision. They had no history of any neurological or visual disorder. A detailed description on the experimental procedure was given to all the subjects prior to the experiment, and a consent was taken from all the subjects before conducting the experiment. The experiment was conducted in accordance with the latest declaration of Helsinki [48].

Experimental paradigms
In this study, we have conducted three experiments: prefrontal, motor/somatosensory, and visual (see Figs. 3 and 4). In Experiment 1, the mental arithmetic (MA), mental counting (MC), and a puzzle solving (PS) tasks related to the prefrontal cortex were investigated. In Experiment 2, the right-hand index finger tapping (RHFT) and right-hand index finger poking (RHFP) tasks associated with the motor and somatosensory cortices were performed. Finally, Experiment 3 contains checkerboard presentation activities for the visual cortex. The subjects were seated on a comfortable chair. They were instructed to relax and avoid body movement as much as possible. The stimuli with associated beep sounds were presented on a monitor screen placed in front of the subject to avoid the onset delay. The experimental paradigm for the three experiments is shown in Fig. 3. Experiment 1 consists of three sessions of the prefrontal mental tasks with a resting state period of 120 sec before each session. The resting period of 120 sec was used to avoid HR fluctuation as much as possible. Each session contains five trials of 30 sec, in which 10 sec is used for task and 20 sec is used for rest.
Mental arithmetic (MA) task: During the MA task, the subjects were instructed to perform backward subtraction of a two-digit number (between 15 to 25) from a three-digit number, successively, by subtracting a new number from the previous answer [37, 46,49,50]. A series of mental arithmetic tasks were shown on the monitor screen, for instance, 550 -15 = ?, 535 -19 = ?, 516 -21 = ?, etc.
Mental counting (MC) task: The subjects were instructed to start backward counting from a three-digit number shown on the monitor screen during the task period [42,[50][51][52].
Puzzle solving task (PS): In the PS task, four pieces of Pusan National University logo were presented to the participants in a random order. The subjects were instructed to reconstruct the Pusan National university logo by flipping and rotating the pieces by their thinking [50,53]. Experiment 2 consists of two sessions and each session has 5 trials. In the first session, the subjects were instructed to start right-hand index finger tapping during the task period. In the second session, poking on the subject's right-hand index finger was done. In Experiment 3, the subjects were advised to focus on the 8 Hz checkerboard flashing on the monitor during the task period. It is noted that the data of Experiments 2 and 3 were used only for the comparison of the initial dips in the prefrontal, motor/somatosensory, and visual cortices.

Signal acquisition and preprocessing
The frequency domain fNIRS system (ISS Imagent, ISS Inc., USA) was used to acquire the brain signals from the prefrontal, motor/somatosensory, and visual cortices. The system utilizes near-infrared light of two wavelengths (690 nm, 830 nm) to determine the concentration changes of both HbO/HbR. Figure 4 shows the electrode placement for signal acquisition in the prefrontal, motor/somatosensory, and visual cortices, in which individual brain areas were thoroughly covered by using 8 emitters and 2 detectors. For instance, the detectors were placed on both sides of FPz for Experiment 1, and C3 for Experiment 2 in accordance to the 10-20 International System for the prefrontal and motor/somatosensory cortices respectively. For the visual cortex (Experiment 3), the detectors were placed near to O1 and O2 locations. A total of 16 channels were made using emitter-detector combination. The distances between emitter-detector pairs and reference points (FPz, C3, O1, O2) are shown in Fig. 4. In acquiring fNIRS signals, a sampling rate of 15.625 Hz was used.
The raw intensity data of ΔHbO and ΔHbR for all the channels were obtained with the ISS Imagent data acquisition and analysis software (ISS-Boxy). Then, they were converted into ΔHbO and ΔHbR using their official software, ISS-Boxy, with extinction coefficients ε HbO = 0.95 µM 1 cm 1 , ε HbR = 4.93 µM 1 cm 1 for 690 nm wavelength and ε HbO = 2.135 µM 1 cm 1 , ε HbR = 1.791 µM 1 cm 1 for 830 nm wavelength in the modified Beer-Lambert law (MBLL) [54]. The converted ΔHbO and ΔHbR are normally contaminated with physiological noises. To remove physiological noises related to respiration and cardiac signals, a 4th order Butterworth low-pass filter with a cutoff frequency of 0.15 Hz was used [23, 39,45]. The lowfrequency drift from the data were then minimized by using a 4th order Butterworth high-pass filter with a cutoff frequency of 0.033 Hz. The cutoff frequency of the high-pass filter was selected as such because the longest possible time period between the onset of two consecutive conditions is 30 sec. In this study, a trial consists of a 10 sec task and 20 sec resting period, and therefore the stimulation frequency was approximately 1/30 sec = 0.033 Hz as in [55,56].

Vector-based phase analysis with a threshold circle
The detection of initial dips was done by employing the vector-based phase analysis with a threshold circle as a decision criterion [30], in which two independent vectors defined by ΔHbO and ΔHbR signals are orthogonally utilized [27][28][29]. Besides ΔHbO and ΔHbR, two other components are defined as follows.
The use of a threshold circle as a decision criteria incorporated in the vector-based phase analysis helps to minimize possible misidentifications of initial dips. The radius of a threshold circle for each channel is determined by detecting the peak value in ΔHbO or ΔHbR signal (i.e., max(HbO) or max(HbR)) in the final 60 sec of the initial resting period (i.e., 120 sec). MATLAB TM function max was used to determine the peak value in ΔHbO and ΔHbR signal. The maximum value among them was then selected. The initial dip instance is concluded when the magnitude of the signal breaks out the threshold circle while its phase lies within the dip phases (i.e., phases 1 to 5).

Features extraction and classification
The channels in which initial dips were detected by the vector-based phase analysis were averaged for the given task. Four different time windows (i.e., 0~1, 0~1.5, 0~2, and 0~2.5 sec) were investigated for comparison of the classification accuracies of initial dips. The signal mean, delta, signal slope, skeweness, and kurtosis values were computed for each window. For this, mean, min, polyfit, skew, and kurt functions available in MATLAB TM were used, in which the min function was used to find the delta value. Then, the extracted features of each subject were normalized between 0 and 1 by the following equation [40,49]: where n yR  represent the original values of the features, ' y is the rescaled values between 0 and 1, max(y) is the maximum value, and min(y) is the minimum value. The rescaled features were then classified using LDA. LDA is a linear classifier that intends to find a separating hyperplane by maximizing the distance between the classes' means and minimizing the within-class variance to classify diverse classes of data [20,49]. LDA is used in various BCI studies due to its ease of use and fast execution time. To determine the classification accuracies, 5-runs of 5-fold cross-validation were used. The 5-fold cross-validation randomly breaks the data into 5 equal sets and used 4 sets for training and 1 set for testing. The process was repeated for five times and the mean accuracy was taken. Furthermore, to measures the performance of the classifier, the true positive rate (TPR), true negative rate (TNR), and false positive rate (FPR) were calculated for three-class classification as follows [57]: , where TP, TN, FP, and FN denote true positive, true negative, false positive, and false negative respectively. The values of TP, TN, FP, and FN were calculated from the confusion matrix [57]. A Matlab TM function, classperf, was used to create the confusion matrix. Finally, the receiver operating characteristic (ROC) curve was drawn using the TPR and FPR for the MA, MC, and PS tasks. Figure 6 shows the phase portrait analysis of the mental arithmetic task (Trial 1) of Subject 1. It can be seen that the initial dips were detected in Channels 2, 4, 8, 9, 10, and 11. Table 1 summarizes the channels in which initial dips were detected for three tasks and eight subjects.

Results
To further confirm whether the channels showing initial dips were really active or not, the tvalues were computed by using the MATLAB function, robustfit [39,45,58]. This function compares the averaged HbO of each channel with the designed HR and returns the t-and pvalues. The criteria for concluding activeness were i) t-value > critical t-value (t crt ) and ii) pvalue < 0.05. In this study, t crt was set to 1.648, as computed from the degree of freedom of the data (i.e., trial period = 30 sec, the number of data points N = 30 × 15.625 = 468, N-1 = 467) and the statistical significance level (i.e., 0.05 for one-tailed test). The t-values of the initial dip detected channels in the vector phase analysis were higher than t crt . Then, the initial dip detected channels were averaged for each given task. Figure 7 shows the averaged HbOs over all 8 subjects and 5 trials, and their standard deviation for individual MA, MC and PS task. The shaded areas along the mean values represent their standard errors. First, for the conventional HRs, the three-class classification of three tasks was performed by using the signal mean and signal slope, as features, for 2~7 sec window, which is tabulated in Table 2. The average classification accuracy of 64.9% was obtained. It is noted that the overall classification accuracies are significantly higher than the chance level (i.e. 33.3%). Particularly, the classification accuracies of Subs. 1 and 8 are higher than the recommended classification accuracy needed for BCI (i.e., 70%). Second, Table 3 shows the averaged classification accuracies based upon the signals during the initial dip period of 0~2.5 sec window, in which two feature combinations out of five features (mean, slope, delta, skewness, and kurtosis) are compared. It is noted that the combination of signal mean and delta gives a significantly higher classification accuracy, 57.5%, than others. It is seen that the classification accuracies of other combinations still show a higher value than the chance level of the three-class classification. Table 4 shows the individual subjects' classification accuracies in 0~1, 0~1.5, 0~2, and 0~2.5 sec windows using the signal mean and the delta as features.
From Table 4, it is seen that the 0~2.5 sec window gives the best classification accuracy for each subject. There is a slight difference in the accuracies of 0~2 and 0~2.5 sec windows, but Sub. 6 in 0~2 sec window case showed almost the chance level. However, in the case of 0~2.5 sec window, all the classification accuracies were far higher than the chance level.
Particularly, Subs. 1 and 8 showed higher classification accuracy than required for BCI (i.e., 70%), whereas those of Subs. 3 and 7 were higher than 60%. Table 5 provides a comparison of the classification accuracies based upon initial dip signals of 0~2.5 sec window and those based upon HRs of 2~7 sec window.    Figures 8 and 9 show the ROC curves for the MA, MC, and PS tasks for HRs and initial dips in 2~7 sec and 0~2.5 sec windows, respectively. The TPR and FPR were calculated for each subject using the confusion matrix. The ROC curves were then drawn for multi-class classification using the averaged TPR and FPR over all subjects for each task. The area under the curve (AUC) was calculated by trapezoidal method. The MATLAB TM function trapz were used. For HR-based classification, the AUC were 0.80, 0.72, and 0.70 for MA, MC, and PS tasks, respectively. The AUC for the initial dip based classification for MA, MC, and PS tasks were 0.75, 0.64, and 0.60, respectively. In both cases, the LDA classifier shows the performance above the random guess line for all tasks.

Comparison of initial dips in various brain cortices
The initial dips in the motor/somatosensory and visual cortices were also detected by using the vector phase analysis method with proper threshold circles. For each task, the signals that have depicted initial dips were averaged over trials and subjects. Figure 10 shows the averaged HbOs of MA, MC, PS, RHFT, RHFP, and checkerboard tasks and their standard deviations.
The characteristic parameters (see Fig. 2, Section 2.1) were calculated from the averaged HbOs in Fig. 10. Those parameters for individual tasks are collected in Table 6. The time to the peak magnitude in the initial dip stage, t o , is almost the same for all tasks. However, there is a slight difference in the duration of initial dip (t init ) among the tasks. In the somatosensory cortex, it is noteworthy to mention that kind of initial dips were detected in a few trials, but it can be said that there is no initial dip phenomenon in the somatosensory cortex (the averaging gives almost zero). Further investigation is needed for this issue. Finally, the initial dip magnitude of the RHFT task (i.e., 0.89529) was higher than any other tasks.

Discussion
In this paper, we have utilized the vector-based phase analysis method and the threshold circle criterion to detect the occurrence of initial dips from the prefrontal cortex during the MA, MC, and PS tasks. The advantage of using a threshold circle in vector-phase analyses is that it helps to minimize a false identification of initial dips from the fNIRS signals [30]. In this paper, the threshold value for each channel was set to the maximum peak value of either ΔHbO or ΔHbR during the resting period just ahead of the task. The use of a peak value of {ΔHbO, ΔHbR} instead of max(ΔHbO 2 + ΔHbR 2 ) 1/2 during the resting state can reduce the radius of a threshold circle, and consequently enables an earlier detection of initial dips. The novelties of this paper exist in the demonstration of the use of initial dips for BCI, suggestion of initial dip features (signal mean, delta), and determination of window size (2.5 sec) for classification of three mental tasks in the prefrontal cortex, which allows a quicker command generation than any conventional HR-based BCI. The early fNIRS studies based upon HRs had to use a large window size: For instance, 0~5, 2~7, 0~10, 0~15, 0~17, and 0~20 sec [37][38][39][40][41][42][43][44][45][46][47]. In this case, due to the inherent delay of 2 sec of initial dip, the 2~7 sec window was shown to provide the best classification performance [46]. Also, the combination of signal mean and signal slope in two class classification based upon HRs for various prefrontal tasks was shown to generate the best accuracy in classification [20,46,51,59]. In Fig. 10, the fNIRS signals of MA tasks were stronger than MC and PS tasks. Naturally, the MC tasks show weaker signals due to that they do not require heavy cognitive load. This finding is consistent with the previous studies [42,51].
For the purpose of comparing the proposed method (i.e., initial dip based classification) and the conventional HR-based method, we performed a similar analysis using the HR signals as well. Adopting the tips from the previous studies in this case, the window size was chosen to be 2~7 sec and the features used were signal mean and signal slope during the 2~7 sec period. Upon these, the average classification accuracy for MA, MC, and PS tasks was 64.9%, see Table 5. Individual subjects' classification accuracies were all higher than the chance level (i.e., 33.3%). Typically, the classification accuracies of Subs. 1 and 8 were higher than the required classification accuracy for BCI (i.e., 70%).
Regarding the HR based classification, the work of Power et al. for 3-class classification (mental arithmetic, mental singing, and rest) achieved an average classification accuracy of 56.2% [38]. Similarly, Herf et al. reported an average ternary classification accuracy of 50.3% in differentiating different levels of workload from the prefrontal measurements [60]. Hong et al. obtained average three-class classification accuracy of 75.6% using the prefrontal and  motor cortex signals [46]. Schuldo and Chau achieved an average classification accuracy of 71.7% for verbal fluency, stroop, and rest tasks from the prefrontal and parietal cortices [61]. It is noted that two brain regions and distinctive tasks were used in [46] and [61], whereas only the prefrontal cortex with similar mental tasks were examined in our work, yielding 64.9%. For the initial dip classification, we investigated five different features (mean, slope, delta, skewness, kurtosis) and four different window sizes (0~1, 0~1.5, 0~2, and 0~2.5 sec). If using only two features, the combination of signal mean and delta provided the average classification accuracy of 57.5% using 0~2.5 sec window, which is the best. Recalling that several previous studies indicate that the initial dip peaks occur at around 2 sec and complete around 4 sec [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18], the window size of 0~2.5 sec is the smallest size that have ever been used in fNIRS studies. Also, since only 64.9% with 7 sec delay was obtained from the HR case, 57.5% with 2.5 sec delay using initial dips is a significant improvement. We also observed that the initial dips of MA, MC and PS tasks occurred in the form of hypoxicischemic type (i.e., decreased ΔHbO), which is consistent with the previous studies [27,28]. The average initial dip magnitude of MA tasks was the largest among the three tasks. Discussing individual subjects' performance, Subs. 1 and 8 achieved a sufficiently higher classification accuracy required for BCI applications. Also, Subs. 3 and 7 showed a classification accuracy higher than 60%. The overall classification accuracy of most subjects was greater than the chance level.
It is also noted that the initial dip phenomenon did not occur in some trials. Therefore, the accuracy was evaluated only for those HRs in which the initial dip appeared. Several factors can be discussed for this issue, for example, caffeine reduces the possibility of detection of initial dip [62,63]. Also, variation of classification accuracy among subjects can be due to individual differences [64].
The parameters of the initial dips detected from the prefrontal cortex were compared with those of motor/somatosensory, and visual cortices, see Table 6. This shows that the vector phase analysis with a threshold circle can successfully detect initial dips in various brain regions. However, in the case of somatosensory, the initial dip was not clearly identified. When performing averaging, the parametric values were indistinctive, which is consistent with [36]. The initial dips in most brain regions were found peaking during approximately 2-2.5 sec, which is in line with the previous fNIRS and fMRI studies [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18]65]. The completion time of initial dip period varies from 3.5 to 4.5 sec. The magnitudes were also varied in different brain regions. Therefore, a further investigation is needed to fully investigate optimal features to distinguish the initial dips for different brain regions.
The following limitations of this study are discussed as follows. In this work, only HbO signal was investigated for examining the signal mean, delta, signal slope, skewness, and kurtosis. Therefore, other features including HbR, COE, and HbT should be investigated for further improvement of initial dip classification accuracy. The second limitation of our study is the insufficient number of subjects. The results would be more acceptable if using a large number than eight. The third is that only a linear classifier (LDA) was used in this work. A more advanced classifier can improve the classification accuracy.

Conclusions
In this paper, a systematic method for using initial dips appearing in fNIRS signals upon brain activities for BCI applications was demonstrated. Three mental tasks including mental arithmetic, mental counting, and puzzle solving tasks in association with the prefrontal cortex were performed. The vector-based phase analysis method with a peak HbO value obtained from the resting period were used in detecting the initial dips occurrence for the three mental tasks. While the average classification accuracy of 64.9% after 7 sec was obtained by using the conventional hemodynamic response method, the proposed method has achieved 57.5% after 2.5 sec by using the initial dip parameters for three classes. The results are encouraging and show a great potential in reduction of window size of fNIRS-based BCI within 0~2.5 sec, resulting in a significant temporal resolution improvement.