Partial Least Squares-discriminant Analysis for the Prediction of Hemodynamic Changes Using Near Infrared Spectroscopy

  • cc icon

    Using continuous wave near-infrared spectroscopy, we measured time-resolved concentration changes of oxy-hemoglobin and deoxy-hemoglobin from the primary motor cortex following finger tapping tasks. These data were processed using partial least squares-discriminant analysis (PLS-DA) to develop a prediction model for a brain-computer interface. The tasks were composed of a series of finger tapping for 15 sec and relaxation for 45 sec. The location of the motor cortex was confirmed by the anti-phasic behavior of the oxy- and deoxy-hemoglobin changes. The results were compared with those obtained using the hidden Markov model (HMM) which has been known to produce the best prediction model. Our data imply that PLS-DA makes better judgments in determining the onset of the events than HMM.


    Near-infrared spectroscopy , Brain-computer interface , Partial least squares-discriminant analysis , Hidden Markov model


    Near-infrared spectroscopy (NIRS) measures temporal and spatial fluctuations of oxy- and deoxy-hemoglobin concen-trations along with total blood volume from the cortical area [1, 2]. It provides a unique tool in neuroscience, since brain activities can be detected noninvasively and potentially in mobile conditions [3, 4]. Applications of NIRS extend from a simple measurement of tissue oxygenation to large-scale imaging correlated with cognitive and motor cortex responses [5, 6].

    Brain-computer interfaces (BCIs) stem from communi-cation and control technology with which a user can facilitate interactions with the outside world using specific signals originating from the brain [7, 8]. Several modalities have been tested for potential BCI applications, including electro-encephalography (EEG) [9], electrocorticography (ECoG) [10], functional magnetic resonance imaging (fMRI) [11, 12] and NIRS [13-15]. These modalities can be categorized as invasive (ECoG) and non-invasive (EEG, fMRI and NIRS) methods. Each of these modalities has advantages, such as high temporal resolution (EEG, ECoG, and NIRS) and the ability to perform imaging of the deep brain areas (fMRI). However, these cannot be easily applied to BCI due to shortcomings, such as low spatial resolution (EEG), invasi-veness (ECoG), and large size and low temporal resolution (fMRI). NIRS has a number of merits to be used for BCI since it provides high temporal resolution, real-time detection, and mobility. Also, since optical signals do not interfere with electrical signals, it can be easily combined with other systems such as EEG [16].

    Regardless of which modality is used, a prediction model is needed so that the brain signal can be interpreted. So far, a number of algorithms have been developed to build a prediction model for NIRS-based BCI such as linear discriminant analysis (LDA) [6], the support vector machine (SVM) [13], and the hidden Markov model (HMM) [13, 17].

    The HMM was modified from the original Markov chain theory by Baum in the 1960s [18], and it has been used for a number of applications, including speech recognition studies [19]. It has also been successfully used as a prediction model for NIRS-based BCI. Researchers compared the SVM and the HMM to evaluate an NIRS prediction model during finger tapping and motor imagery studies and concluded that HMM performs better than SVM [13]. HMM has also been used to classify mental arithmetic and music imagery signals in a separate study [20]. PLS-DA is a combined version of the partial least square (PLS) regression method and the LDA [6, 21]. Recently, applications of the PLS-DA have been expanded to behavioral analysis for positron emission tomography (PET) [22], the spatiotemporal study of auditory and visual perception [23], and structure seeking[17] using fMRI. The LDA was used to classify motor imagery signals using hemodynamic responses obtained by NIRS [24].

    In this study, we measured the hemodynamic responses for finger tapping in the left hemisphere motor cortex using NIRS and developed an analytic tool for classification between the finger tapping task and relaxation. Using the PLS-DA method, we developed the prediction model to discriminate TASK and REST states and compared our results with those obtained using the HMM.


       2.1. Subjects

    Seven volunteers, 3 men and 4 women, with a mean age of 24 (± 5.5 years) partici-pated in this study. All participants were right-handed and had black hair. They had no history of any neuropathologic conditions, nor were they taking any specific medicines that might interfere with the study.

       2.2. Experimental Procedure

    NIRS activities were recorded from the left hemisphere primary motor cortex located in the C3 area [25]. The C3 position in the brain is shown in Fig. 1. Fig. 1(a) shows the top view of the brain, and Fig. 1(b) indicates the left view of the face. A custom-made headset accommodated eight light source pairs (785 and 830 nm) and two detectors. The source and detector distances were maintained at about 3.5-4 cm, which corresponds to a penetration depth of 1 to 2 cm in the adult head [26]. The light sources are delivered using optical fibers with diameters of 400 ㎛, and the reflected lights are coupled with two detectors using 3 mm-diameter fiber bundles. The optical fibers are arranged

    as shown in Fig. 2(a). With eight sources and two detectors, the total number of channels is 16, and the channels with large signal fluctuations were chosen for analysis.

    Before the experiment, every subject was asked to sit on a chair with hands resting upon each knee. The experi-mental protocol is shown in Fig. 2(b), which was composed of a series of TASK and REST sessions. The TASK session was composed of 15 second-long finger tapping periods using all five fingers with frequency of 3-4 Hz, after which the subject was totally relaxed for 45 sec (REST session). These series of sessions lasted for 20 minutes. During the REST session, the subjects are asked to remain motionless.

       2.3. Calculation of Hemoglobin Concentration

    To remove excessive high frequency noise and low frequency drift in hemodynamic responses, we adopted a moving average filter with every 3000 points divided from raw ac intensity data. The signal was then low- and high-pass filtered with cutoff frequencies of 0.5 and 0.01 Hz.

    A modified Beer-Lambert law was used to calculate the oxy-hemoglobin (ΔHbO2) and the deoxy-hemoglobin (ΔHbr) concentrations as in Eq. 1.


    Where OD is the optical density, IFinal is the measured light density, IInitial is the initial light intensity, ε represents the extinction coefficients of the oxy and deoxy-hemoglobin, B is the differential path-length factor, and L is the source-detector separation length at the given wavelength (785 and 830 nm). The change in concentration of oxy-hemoglobin (ΔHbO2) and deoxy-hemoglobin (Hbr) were calculated using the following equations:




       2.4. Hidden Markov Model

    The hidden Markov model (HMM) [18, 19] is a well-known statistical tool for modeling time series data. It contains discrete states, S, and emits an observation vector, ot, at every time t. It uses two sets of probability vectors, including the state transition probability distribution A = {aij} and the observation probability distribution, B = {bj(ot)} at each state j = 1, …, S. In general, each observation probability distribution is modeled by Gaussian mixtures [13]. To determine these parameters, we adopted the classification algorithm defined in a previous study [17]. So-called “train and test” was performed using the hemodynamic responses data from each participant; the whole data set was divided into two sub-sets: one training set and one test set. During the training process, the state transition probability distribution (A), obser-vation probability distribution (B), and initial state probability distribution (π) were optimized using the expectation-maximi-zation algorithm [13]. The log-likelihood (LL) was determined both from the TASK and REST session data using the forward-backward algorithm [13]. In HMM, the higher LL value indicates the test data lean towards the specific model (TASK or REST). We performed cross-validation to determine the classification accuracy between the training and test data set. The actual calculation was done using a public toolbox written in Matlab.

       2.5. Partial Least Squares-Discriminant Analysis (PLS-DA)

    PLS-DA is one of the classification methods and is based on the PLS regression algorithm [27]. PLS regression is a method for modeling a linear relationship between a set of output variables (chemical concentration), y, and a set of input variables (regressors), X, which can be expressed as y = Xb + b0. Here X is an n × p matrix, containing the values of the predictor values at n data points, b = (b1, …, bp)T is a p-dimensional column vector which contains regression coefficients; the T denotes the transpose of a vector or matrix.The regression coefficients, b, was defined as




    In this study, we noted y as TASK (1) and REST (-1), determined the regression model using principal components, and computed classification accuracy using cross-validation. We used the PLS_toolbox (Eigenvector Research, Inc., WA, USA) for analysis.


    In one set of experiments, the TASK session was continued for a whole minute to observe hemodynamic changes (Fig. 3). The maximum peak was observed at the middle of duration(around 30-40 sec). The low p-value (< 0.0001) between oxy-HB and deoxy-HB in the TASK session indicates strong anti-correlation.

    Out of 16 channels with 8 sources and two detectors, only the 8 channels connecting the closest source-detector combinations were used for our experiment. Only statistically significant signals were counted for our BCI study. A typical signal fluctuation is shown in Fig. 4. The signal indicates that a rapid reaction from the finger tapping action is reported after 2~5 seconds latency. Also, oxy-Hb started to decrease immediately when the finger tapping ceased.

    In Fig. 5, an example of the 8-channel signals is shown. The onset of the TASK session was initiated at time zero in each graph. Only the signals from channels 4 and 5 were consistent with expectations. Signals from other channels

    do not seem to reflect the expected hemodynamic changes. There are two possible reasons for this low success rate; the diffusely reflected light from the cortical area does not pass through the motor cortex of the participant, or a sufficient amount of light was not received on the detector to create significant signals. For many participants in general, the light transmitted through the scalp, skull, and cortex area was

    too weak to generate a signal on the detector. However, approximately 3~5 cm separation between a source and a detector is known to be sufficient to cover the motor cortex region [28, 29]. A prediction model based on HMM and PLS-DA was developed to classify TASK and REST signals from raw oxy-Hb data. The developed model was used to predict the intention of the participants when the oxy-Hb signal with different time segments (15, 10, and 5 sec) was used as input data. The results are shown in Table. 1 and 2 and visualized in Fig 6.

    As anticipated, in most cases, the prediction result scored the highest when longer time segments were used. More interestingly, the PLS-DA and HMM-based models produced similar predictions, while PLS-DA created a slightly better result when longer time segments were used. When the 15 sec data sets from the prediction model were compared, the success rate was 90.21% for PLS-DA and 85.78% for HMM. Additionally, Table. 1 and 2 imply that longer time segments do not always guarantee a higher success rate. For the case of HMM results, the success rate for the 10 sec segment (79.76%) was slightly worse than the 5 sec segment success rate (80.29%).

    At this stage, it is not fair to conclude that PLS-DA is superior to HMM in constructing the prediction model using

    the hemodynamics data. Even if the HMM is frequently used in obtaining time series data pattern recognition, the time segmented data used to construct the model in this study probably was not large enough to develop a robust model. On the other hand, since PLS-DA can obtain the lowest predicted residual sum of squares (PRESS) using iterative methods, such as nonlinear iterative partial least squares (NIPALS) or the simple partial least squares (SIMPLS) algorithm [30], the effect of short time segmentation may not have been as pronounced.

    It is well known that the strength of the finger tapping

    significantly affects the degree of the hemodynamic changes [31]. Therefore, each model used for prediction of the participant’s intention using NIRS has to be individually calibrated. In our study, the data from a given participant were used to construct the prediction model for that individual.

    Recent BCI studies using NIRS have shown great potential in a number of neuroscience studies, such as a mirror study in the motor cortex [5], a beverage preference test in the prefrontal lobe [20], and a word speller application [22]. Since the hemodynamic signal fluctuation is caused not only by certain specific stimuli but by other natural causes, a robust model based on PLS-DA may give a higher success rate for BCI applications using NIRS.

  • 1. Boas D. A, Gaudette T, Strangman G, Cheng X, Marota J. J, Mandeville J. B (2001) “The accuracy of near infrared spectroscopy and imaging during focal changes in cerebral hemodynamics” [Neuroimage] Vol.13 P.76-90 google
  • 2. Jobsis F. F (1977) “Noninvasive infrared monitoring of cerebral and myocardial oxygen sufficiency and circulatory parameters” [Science] Vol.198 P.1264-1267 google doi
  • 3. Muehlemann T, Haensse D, Wolf M (2008) “Wireless miniaturized in-vivo near infrared imaging” [Opt. Express] Vol.16 P.10323-10330 google doi
  • 4. Muehlemann T, Reefmann N, Wechsler B, Wolf M, Gygax L (2011) “In vivo functional near-infrared spectroscopy measures mood-modulated cerebral responses to a positive emotional stimulus in sheep” [Neuroimage] Vol.54 P.1625-1633 google doi
  • 5. Shimada S (2009) “Modulation of motor area activity by the outcome for a player during observation of a baseball game” [PLoS One] Vol.4 P.e8034 google doi
  • 6. Tai K, Chau T (2009) “Single-trial classification of NIRS signals during emotional induction tasks: towards a corporeal machine interface” [J. Neuroeng. Rehabil.] Vol.6 P.39 google doi
  • 7. Wolpaw J. R, Birbaumer N, Heetderks W. J, McFarland D. J, Peckham P. H, Schalk G, Donchin E, Quatrano L. A, Robinson C. J, Vaughan T. M (2000) “Brain-computer interface technology: a review of the first international meeting” [IEEE Trans. Rehabil. Eng.] Vol.8 P.164-173 google doi
  • 8. Daly J. J, Wolpaw J. R (2008) “Brain-computer interfaces in neurological rehabilitation” [Lancet Neurol.] Vol.7 P.1032-1043 google doi
  • 9. Lotte F, Congedo M, Lecuyer A, Lamarche F, Arnaldi B (2007) “A review of classification algorithms for EEG-based brain-computer interfaces” [J. Neural Eng.] Vol.4 P.R1-R13 google doi
  • 10. Brunner P, Ritaccio A. L, Emrich J. F, Bischof H, Schalk G (2009) “A brain-computer interface using event-related potentials (Erps) and electrocorticographic signals (Ecog) in humans” [Epilepsia] Vol.50 P.389-389 google
  • 11. Weiskopf N, Mathiak K, Bock S. W, Scharnowski F, Veit R, Grodd W, Goebel R, Birbaumer N (2004) “Principles of a brain-computer interface (BCI) based on real-time functional magnetic resonance imaging (fMRI)” [IEEE Trans. Biomed. Eng.] Vol.51 P.966-970 google doi
  • 12. Sitaram R, Caria A, Veit R, Gaber T, Rota G, Kuebler A, Birbaumer N (2007) “FMRI brain-computer interface: a tool for neuroscientific research and treatment” [Intell. Neurosci.] Vol.2007 P.25487 google
  • 13. Sitaram R, Zhang H, Guan C, Thulasidas M, Hoshi Y, Ishikawa A, Shimizu K, Birbaumer N (2007) “Temporal classifi-cation of multichannel near-infrared spectroscopy signals of motor imagery for developing a brain-computer interface” [Neuroimage] Vol.34 P.1416-1427 google doi
  • 14. Coyle S, Ward T, Markham C, McDarby G (2004) “On the suitability of near-infrared (NIR) systems for next-generation brain-computer interfaces” [Physiol. Meas.] Vol.25 P.815-822 google doi
  • 15. Coyle S. M, Ward T. E, Markham C. M (2007) “Brain-computer interface using a simplified functional near-infrared spectroscopy system” [J. Neural Eng.] Vol.4 P.219-226 google doi
  • 16. Roche-Labarbe N, Zaaimi B, Berquin P, Nehlig A, Grebe R, Wallois F (2008) “NIRS-measured oxy- and deoxy-hemoglobin changes associated with EEG spike-and-wave discharges in children” [Epilepsia] Vol.49 P.1871-1880 google doi
  • 17. Power S. D, Falk T. H, Chau T (2010) “Classification of prefrontal activity due to mental arithmetic and music imagery using hidden Markov models and frequency domain near-infrared spectroscopy” [J. Neural Eng.] Vol.7 P.026002 google doi
  • 18. Ghahramani Z (2001) “An introduction to hidden Markov models and Bayesian networks” [Int. J. Pattern Recogn.] Vol.15 P.9-42 google doi
  • 19. Rabiner L. R (1989) “A tutorial on hidden Markov-models and selected applications in speech recognition” [Proc. IEEE] Vol.77 P.257-286 google doi
  • 20. Luu S, Chau T (2009) “Decoding subjective preference from single-trial near-infrared spectroscopy signals” [J. NeuralEng.] Vol.6 P.016003 google doi
  • 21. Stahle L, Wold S (1987) “Partial least squares analysis with cross-validation for the two-class problem: a Monte Carlo study” [J. Chemometr.] Vol.1 P.185-196 google doi
  • 22. McIntosh A. R, Bookstein F. L, Haxby J. V, Grady C. L (1996) “Spatial pattern analysis of functional brain imagesusing partial least squares” [Neuroimage] Vol.3 P.143-157 google doi
  • 23. McIntosh A. R, Chau W, Protzner A. B (2004) “Spatiotemporalanalysis of event-related fMRI data using partial least squares” [Neuroimage] Vol.23 P.764-775 google doi
  • 24. Holper L, Biallas M, Wolf M (2009) “Task complexity relates to activation of cortical motor areas during uni- and bimanual performance: a functional NIRS study” [Neuroimage] Vol.46 P.1105-1113 google doi
  • 25. Homan R. W, Herman J, Purdy P (1987) “Cerebral location of international 10-20 system electrode placement” [Electro-encephalogr.Clin. Neurophysiol.] Vol.66 P.376-382 google doi
  • 26. Strangman G, Boas D. A, Sutton J. P (2002) “Non-invasive neuroimaging using near-infrared light” [Biol. Psychiat.] Vol.52 P.679-693 google doi
  • 27. Krishnan A, Williams L. J, McIntosh A. R, Abdi H (2011) “Partial least squares (PLS) methods for neuroimaging: a tutorial and review” [Neuroimage] Vol.56 P.455-475 google doi
  • 28. Franceschini M. A, Joseph D. K, Huppert T. J, Diamond S. G, Boas D. A (2006) “Diffuse optical imaging of the whole head” [J. Biomed. Opt.] Vol.11 P.054007 google doi
  • 29. Franceschini M. A, Boas D. A (2004) “Noninvasive measure-ment of neuronal activity with near-infrared optical imaging” [Neuroimage] Vol.21 P.372-386 google doi
  • 30. Wold S, Sjostrom M, Eriksson L (2001) “PLS-regression: a basic tool of chemometrics” [Chemometr. Intell. Lab.] Vol.58 P.109-130 google doi
  • 31. Nambu I, Osu R, Sato M. A, Ando S, Kawato M, Naito E (2009) “Single-trial reconstruction of finger-pinch forces from human motor-cortical activation measured by near-infrared spectroscopy (NIRS)” [Neuroimage] Vol.47 P.628-637 google doi
  • [FIG. 1.] Location of the motor cortex in the human brain.
    Location of the motor cortex in the human brain.
  • [FIG. 2.] Arrangements of the light sources (blue circles) and the detectors (red rectangles) (A) and the experimental protocol for Task/Rest sessions lasting for 20 minutes (B).
    Arrangements of the light sources (blue circles) and the detectors (red rectangles) (A) and the experimental protocol for Task/Rest sessions lasting for 20 minutes (B).
  • [FIG. 3.] The anti-phasic nature of oxy- and deoxy-Hb fluctuation during the 1-minute Task session. The red line is the oxy-Hb and the blue one is the deoxy-Hb. Vertical green lines indicate the start and stop time of the finger tapping.
    The anti-phasic nature of oxy- and deoxy-Hb fluctuation during the 1-minute Task session. The red line is the oxy-Hb and the blue one is the deoxy-Hb. Vertical green lines indicate the start and stop time of the finger tapping.
  • [FIG. 4.] A typical example of the oxy-Hb and deoxy-Hb fluctuation during a series of TASK and REST sessions.
    A typical example of the oxy-Hb and deoxy-Hb fluctuation during a series of TASK and REST sessions.
  • [FIG. 5.] Hemodynamic signals from all eight channels. Oxy-Hb (red) deoxy-Hb (blue) and total-Hb (green). The units of the y-axis is arbitrary.
    Hemodynamic signals from all eight channels. Oxy-Hb (red) deoxy-Hb (blue) and total-Hb (green). The units of the y-axis is arbitrary.
  • [TABLE 1.] Test results with HMM in 15 10 and 5 sec data with 7 participants
    Test results with HMM in 15 10 and 5 sec data with 7 participants
  • [TABLE 2.] Prediction results with PLS-DA in 15 10 and 5 sec data with 7 participants
    Prediction results with PLS-DA in 15 10 and 5 sec data with 7 participants
  • [FIG. 6.] Results of the PLS-DA for 15 sec (A 92.5%) 10 sec (B 83%) and 5 sec (C 79%) data. Predicted value 1 indicates the Task and -1 indicates the Rest. The actual values produced using the PLS-DA are classified based on their polarity.
    Results of the PLS-DA for 15 sec (A 92.5%) 10 sec (B 83%) and 5 sec (C 79%) data. Predicted value 1 indicates the Task and -1 indicates the Rest. The actual values produced using the PLS-DA are classified based on their polarity.