A ROBUST AND PRACTICAL METHOD TO SEPARATE PERIODIC

SIGNALS FROM MEG DATA USING SECOND ORDER STATISTICS

Hidekazu Fukai

Department of Engineering, Gifu University, Yanagido 1-1, Gifu-shi, Gifu-ken, Japan

Keywords:

Blind source separation, Magnetoencephalography, Periodic signals, Second-order statistics.

Abstract:

The analyses of recordings of magnetoencephalography (MEG) and other imaging techniques may require

the separation of periodic signals from the observed signals. Blind source separation (BSS) is widely used

for the separation of speciﬁc signals these days. Though several algorithms based on the BSS scheme for

the separation of periodic signals have been proposed, they usually assume the system to be well-posed,

satisfactory results often cannot be obtained for practical recordings. In this study, we show that a method

based on the joint approximate diagonalization of correlation matrices with several time delays (JADCM)

is robust and good results can be obtained by choosing the time delays carefully, especially in practical ill-

posed situations such as signal separation from MEG recordings. The performance of the proposed method is

compared with that of Periodic BSS and JADCM using the conventional parameter set.

1 INTRODUCTION

The blind source separation (BSS) scheme is a pow-

erful tool for the analysis of recordings of magnetoen-

cepharography (MEG) and other imaging techniques

to separate artifacts, noise, and brain signals. The po-

tential of BSS for use in MEG data analysis has been

conﬁrmed in many studies (Hyv¨arinen et al., 2001).

Generally, we consider a linear and instantaneous

superposition of independent source signals and mul-

tidimensional observations of the form

x(t) = A· s(t)

for the BSS problem, where x(t) = (x

(t),· ·· ,x

(t))

denote the observations, s(t) = (s

(t), ··· , s

(t))

de-

note the unknown source signals, and the mixing ma-

trix A is the transfer function between sources and

sensors. BSS involves identifying A and retrieving

the source signals s(t) without a priori information

about the mixing matrix A.

Many types of BSS methods have been proposed

thus far (Hyv¨arinen et al., 2001) (Theis and Inouye,

2006). Each of them require additional conditions for

separation, so that has each suitable data types. We

have to choose the methods and parameters carefully

depending on the entities to be separated or extracted.

MEG data to which the BSS scheme is applied

have several features such that (c1) most source sig-

nals have time structures; (c2) the number of source

signals is larger than that of observations, i.e., the

problem is ill-posed; and (c3) some sources are pe-

riodic.

The second-orderstatistics is useful for signal sep-

aration in case each source signal has speciﬁc auto-

correlation and no cross-correlations. In this case,

correlation matrices of source signals for arbitrary

time delays are diagonal. Furthermore, if a correla-

tion matrix for a time delay has distinct eigenvalues,

we can separate the source signals by the diagonal-

ization of the correlation matrix (Belouchrani et al.,

1997) (Ziehe and M¨uller, 1998). Several methods for

the separation of periodic signals have been proposed

(Barros and Cichocki, 2001) (Jafari et al., 2006).

Joint diagonalizations of several correlation ma-

trices with different time delays are effective for the

case of correlation matrices having the same or al-

most the same eigenvalues (Belouchrani et al., 1997)

(Ziehe and M¨uller, 1998). Several methods based on

the joint approximated diagonalizationhave been pro-

posed (Theis and Inouye, 2006) and applied to MEG

data analysis (Ziehe et al., 2000).

In the joint diagonalization process, we have to

choose the values and the total number of time delays

for the correlation matrices. If the problem is well-

posed, as assumed in most proposed methods, the se-

lection of the time delays does not critically affect the

results of the separation. On the other hand, in the

case of practical application to MEG data, the number

372

Fukai H..

A ROBUST AND PRACTICAL METHOD TO SEPARATE PERIODIC SIGNALS FROM MEG DATA USING SECOND ORDER STATISTICS.

DOI: 10.5220/0003296103720377

In Proceedings of the International Conference on Bio-inspired Systems and Signal Processing (BIOSIGNALS-2011), pages 372-377

ISBN: 978-989-8425-35-5

 2011 SCITEPRESS (Science and Technology Publications, Lda.)

of source signals is larger than that of observations,

so that the problem is ill-posed (c2). In this case,

the information has already been lost during obser-

vation, and generally, we can only get approximated

source signals from any decomposing process. There

is no matrix that diagonalizes all the correlation ma-

trices simultaneously. The approximated results vary

considerably depending on the selection of the values

of time delays, but very few reports have discussed

the selection of the values of time delays (Tang et al.,

2005).

In order to observe and separate speciﬁc signals,

we can use the features of the signals for the selec-

tion of parameters of the method. In this paper, we

propose practical and effective values and numbers of

time delays for the extraction of periodic signals (c3)

from MEG data.

In this paper, we treat 60 Hz power supply noise

as an example of separate and eliminate periodic sig-

nals from MEG data. The BSS-type method we

use is based on the approximately joint diagonal-

ization of several time-delayed correlation matrices

without noise term (Ziehe and M¨uller, 1998). We

use the Jacobi-like joint diagonalization method (Car-

doso and Souloumiac, 1996). The performance of

the proposed method is compared with that of con-

ventional usage of JADCM and Periodic BSS (Jafari

et al., 2006), which is a sequential algorithm based on

second-order statistical information.

2 METHODS

2.1 Experimental Procedure

MEG data were obtained by using a whole-

head helmet-shaped 64-channel SQUID sensor array

(Model 100, CTF Systems Inc.) in a magnetically

shielded room. Each MEG sensor has an axial gra-

diometer that approximately measures the magnetic

ﬁeld in the normal direction. One healthy volunteer

(32-year-old male) participated this study. After ex-

plaining the nature of the study, informed consent

was obtained from him. The experimental procedures

were in accordance with the Declaration of Helsinki.

The subject had no history of neurological and psy-

chiatric disorders. During recording, the subject was

seated relaxed under a helmet-shaped device. The

sampling ratio was 625 Hz.

2.2 Blind Separation using

Second-order Statistics

The problem of blind source separation (BSS)

consists of recovering unknown source signals

s(t) = (s

(t),. .. ,s

(t))

from observations x(t) =

(t),. .. ,x

(t))

, generated by the mixing model

x(t) = As(t)

where the mixing operator A ∈ R

m×n

is a ﬁxed but

unknown matrix and t is discrete time index (Ziehe

and M¨uller, 1998).

Generally, we suppose well-posed problem, i.e.,

the matrix A has full column rank. In the case where

the number of sources is equal to the number of mea-

surements m = n, the model becomes exactly deter-

mined problem. The separating system then estimates

the original source signals according to

y(t) = Bx(t) = BAs(t) = ΛPs(t) ∝ s(t),

where y(t) = (y

(t),. .. ,y

(t))

represents the recov-

ered sources. Matrices B, Λ and P ∈ R

m×m

are sep-

arating, scaling and permutation matrix respectively.

Because we cannot know the variances of sources, we

assume the correlation matrix of s(t) as unit matrix

and each s(t) has mean zero without losing generality.

In the joint diagonalization scheme, we consider cor-

relation matrices C

(τ) of time-lagged mixed signals

x(t),

(τ)

def

= E{x(t)x

(t + τ)}

where E{ ·} denotes statistical expectation over t and

τ is a time-shift parameter. Practically the expecta-

tion is computed from the available data with a ﬁnite

sample size N as a sample average

(τ

) =

N−τ

∑

t=1

(x(t)x

(t + τ

)).

We see that the correlation of x(t) is related to the

correlation of s(t) according to

(τ) = E{(As(t))(As(t + τ))

}

= AE{ s(t)s(t + τ)

= AC

(τ)A

due to the linearity of the expectation operator and the

mixing model. Here, all cross-correlation terms of s

which are the off-diagonal elements of C

(τ) are zero

for independent signals and thus C

(τ) is a diagonal

matrix. If A is invertible, this can also be written as

(τ)B

= C

(τ)

where the matrix B = A

−1

diagonalizes C

(τ).

A ROBUST AND PRACTICAL METHOD TO SEPARATE PERIODIC SIGNALS FROM MEG DATA USING

SECOND ORDER STATISTICS

373

We take two steps to identify B, i.e., whitening

step and joint-diagonalization step.

The whitening step is achieved by applying to x(t)

a whitening matrix W satisfying:

E{Wx(t)x(t)

} = WC

(0)W

= WAA

= I.

In other words, whitening step correspond to principal

component analysis (PCA) procedure which removes

correlations between the observations x(t). This fol-

lows that if W is a whitening matrix, then there ex-

ists a orthogonal matrix U such that WA = U for

any whitening matrix W. With the orthogonal ma-

trix U, the mixture matrix A can be determined as

A = W

−1

U. The whitening step reduces the determi-

nation of the matrix A to that of a orthogonal matrix

The next joint-diagonalization step is achieved

by ﬁnding U as a rotetion matrix which diagonalize

the time-delayed correlation matrices E{Wx(t)x(t +

)

} simultaneously. Here, any rotating proce-

dure does not change the correlation matrix of zero

time delay.

Though only one additional time delayed corre-

lation matrix is required theoretically (Tong et al.,

1991)(Molgedey and Schuster, 1994), we may fail to

ﬁnd the diagonalization matrix B if C

(τ) does not

have distinct eigenvalues. This can be avoided by

simultaneous diagonalization of multiple C

(τ

),1 ≤

k ≤ K (Belouchrani et al., 1997)(Ziehe and M¨uller,

1998). We use Jacobi like algorithm proposed by

Cardoso and Souloumiac (Cardoso and Souloumiac,

1996) based on several Givens rotations to ﬁnd U.

2.3 Selection of Time Delays for

Periodic Sources on Ill-posed

Problem

In most theoretical discussions on the BSS problem,

we treat exactly determined problem in which we

can get exact solutions. In the joint-diagonalization

scheme, not much focus on the selection of time

delays is required if the problems are well-posed.

Theoretically, in case each source signal s(t) has its

auto-correlations and the correlation matrix for the

time delay has distinct eigenvalues, only one time de-

lay is required. Several time delays and their joint-

diagonalizations are effective in case the correlation

matrices have the same or almost the same eigen-

values (Belouchrani et al., 1997)(Ziehe and M¨uller,

1998).

However, the problems are ill-posed in most prac-

tical cases, especially in MEG recordings, and only

approximated solutions can be obtained. There is

no orthogonal matrix that makes all the correlation

matrices with the time delay diagonal on the prac-

tical ill-posed situations. In this case, the joint-

diagonalization algorithm minimizes the sum of the

off-diagonal values to obtain the target orthogonal

matrix.

Generally, in ill-posed problems, the separation

performance depends on the trade-off between the

source signals. In the case of joint-diagonalization,

the results vary depending on the selection of the set

of the time delay τ. If speciﬁc sources have to be ex-

tracted, we have to set the number and value of time

delays depending on the condition and achive good

separation performance.

In case the target signals have speciﬁc features in

the auto-correlation function, we can select the values

of time delays depending on that. If the signals are

periodic with cycle T, the auto-correlation function

has peaks on the time delay kT, k = 1, 2, · ·· and we

can get larger evaluation values to minimize the peaks

by setting τ = kT, k = 1, 2, · ··.

Furthermore, if brain signals are not the target to

be separated, it is important to choose the time delays

to avoid the effect of brain signals. Because almost all

brain signals have prominent auto-correlation values

within 1 s of the time delay, it is effective to choose

time delays over 1 s to avoid the effect of brain sig-

nals.

In conclusion, for the separation and elimination

of power supply noise of 60 Hz, we set the time delay

as τ = kT, k = 61, 62, ··· , and T = 1/60 (s). The total

number of τ may be chosen such that it is sufﬁcient for

convergence.

3 RESULTS

As an example of the separation of periodic signals

from MEG data, we tried to separate 60 Hz power

supply noise and eliminate them from the data.

Figure 1 shows auto-correlation functions of sep-

arated signals (SS). We selected ﬁve signals from sep-

arated 64 signals, labeled SS1 to SS4, which have

remarkable peaks showing time delays associated

with 60 Hz cycle and its higher harmonics on auto-

correlation functions. The peak values become con-

stant over 0.05 s, which approximately indicates the

noise originated power supply has little temporal ﬂuc-

tuations. On the other hand, the signals originate from

the brain activities have temporal ﬂuctuations, so that

the values of auto-correlations decrease after a short

period (SS9). This difference is a useful factor for the

separation. The time delay parameter used here was

τ = kT, k = 61,62,·· · ,80, T = 1/60 (s) for 60 Hz

power supply noise. This parameter set was chosen

BIOSIGNALS 2011 - International Conference on Bio-inspired Systems and Signal Processing

374

㪄㪈

㪄㪇㪅㪌

㪇

㪇㪅㪌

㪈

㪇㪅㪌

㪈

㪇

㪇㪅㪌

㪈

㪇㪅㪇㪇㪅㪈㪇㪅㪉㪇㪅㪊㪇㪅㪋㪇㪅㪌

㪇

㪇㪅㪌

㪈

㪪㪪㪈

㪪㪪㪊

㪪㪪㪉

㪪㪪㪋

㪇

㪘㫌㫋㫆㪄㪺㫆㫉㫉㪼㫃㪸㫋㫀㫆㫅

㪛㪼㫃㪸㫐㩷㫋㫀㫄㪼㩷㩿㫊㪀

㪇㪅㪇㪇㪅㪈㪇㪅㪉㪇㪅㪊㪇㪅㪋㪇㪅㪌

㪪㪪㪐

㪇㪅㪇㪇㪅㪈㪇㪅㪉㪇㪅㪊㪇㪅㪋㪇㪅㪌

㪄㪈

㪇

㪈

Figure 1: Auto-correlation functions of separated signals

related to power supply noise (SS1–SS4) and brain signal

(SS9).

to utilize peaks of auto-correlation functions of target

signals and circumvent the effect of auto-correlation

values of signals originating from brain activities.

Figure 2 shows the performance indexes of sep-

aration to estimate a sufﬁcient number of correlation

matrices. The abscissas are total number of correla-

tion matrices, used for joint-diagonalization, the num-

ber of τ. The performance index is deﬁned as aver-

age of auto-correlation values of each separated sig-

nal y(t) on the peaks correspond to 60Hz. The index

values are calculated as

− K

∑

k=K

N−kT

∑

t=1

(t)y

(t +kT)), j = 1,··· ,64,

where (K

) was set as (61,80) after auto-

correlation functions originating from other biologi-

cal signals decrease. All index values calculated for

64 components are plotted.

Using conventional τ set, τ = 1, 2, ··· , we need

over 10 correlation matrices until we get converged

results (Fig.2 (1)). On the other hand, the index val-

ues converge at the total number of τ is about 3 or 4

with proposed τ set, τ = kT,k = 61, 62, ·· · , 80. This

smaller total number of correlation matrices for joint-

diagonalization has an advantage in calculation time.

㪇㪉㪋㪍㪏㪈㪇㪈㪉㪈㪋㪈㪍㪈㪏㪉㪇

㪫㫆㫋㪸㫃㩷㫅㫌㫄㪹㪼㫉㩷㫆㪽㩷㱠

㪧㪼㫉㪽㫆㫉㫄㪸㫅㪺㪼㩷㪠㫅㪻㪼㫏

㪇

㪇㪅㪈

㪇㪅㪉

㪇㪅㪊

㪇㪅㪋

㪇㪅㪌

㪇㪅㪍

㪇㪅㪎

㪇㪅㪏

㪇㪅㪐

㩿㪉㪀㩷τ㩷㪔㩷㫂㪫㪃㩷㫂㩷㪔㩷㪍㪈㪃㩷㪍㪉㪃㩷㪅㪅㪅㪃㩷㪏㪇㩷㩿㫇㫉㫆㫇㫆㫊㪼㪻㩷㫇㪸㫉㪸㫄㪼㫋㪼㫉㩷㫊㪼㫋㪀

㩿㪈㪀㩷τ㩷㪔㩷㪈㪃㪉㪃㩷㪅㪅㪅㩷㩷㩿㪺㫆㫅㫍㪼㫅㫋㫀㫆㫅㪸㫃㩷㫇㪸㫉㪸㫄㪼㫋㪼㫉㩷㫊㪼㫋㪀

㪇㪉㪋㪍㪏㪈㪇㪈㪉㪈㪋㪈㪍㪈㪏㪉㪇

㪇

㪇㪅㪈

㪇㪅㪉

㪇㪅㪊

㪇㪅㪋

㪇㪅㪌

㪇㪅㪍

㪇㪅㪎

㪇㪅㪏

㪇㪅㪐

Figure 2: Performance Indexes of separation improved as

total number of correlation matrix to diagonalize increase.

Furthermore, the score of the performance index af-

ter convergence is also better on proposed τ set (Fig.2

(2)).

The separation performanceis also compared with

Periodic BSS (Jafari et al., 2006). Figure 3 shows

the performanceindexesof separated signals obtained

by Periodic BSS (dashed), by JADCM with a con-

ventional parameter set, τ = 1, 2, ··· (dotted), and

by JADCM with proposed parameter set τ = kT, k =

61, 62, ··· , 80 (solid). The separated signals are

sorted according to the performance index, and on

10 from 64 signals are shown. In the case of an ill-

posed system, the proposed method shows good per-

formance.

To evaluate the noise reduction performance of

the proposed method, we reconstructed the original

MEG sensor signals from output signals after reject-

ing four components that were related to power sup-

ply noise corresponding to SS1 to SS4 (Fig.1). Figure

3 (1) shows the power spectral density function of the

original row data from MEG channel L42. Promi-

nent peaks can be conﬁrmed at frequencies close to

those of the power supply and also at higher harmon-

ics. Even after noise reduction using Periodic BSS

(Fig.3 (2)) nor JADCM with the ordinary parameter

set, τ = 1, ··· , 20, (Fig.3 (3)) the peaks are still ob-

served. After noise reduction using the proposed pa-

rameter set, τ = kT, k = 61, 62, ·· · , 80, peaks related

to power supply noise are hardly found.

A ROBUST AND PRACTICAL METHOD TO SEPARATE PERIODIC SIGNALS FROM MEG DATA USING

SECOND ORDER STATISTICS

375

㪈㪉㪊㪋㪌㪍㪎㪏㪐㪈㪇

㪇㪅㪇

㪇㪅㪈

㪇㪅㪉

㪇㪅㪊

㪇㪅㪋

㪇㪅㪌

㪇㪅㪍

㪇㪅㪎

㪇㪅㪏

㪧㪼㫉㪽㫆㫉㫄㪸㫅㪺㪼㩷㪠㫅㪻㪼㫏

㪪㪪㩷㫅㫌㫄㪹㪼㫉

Figure 3: Performance indexes of SSs separated by Peri-

odic BSS (dashed), JADCM with conventional parameter

set τ = 1, 2, ··· , 20 (dotted), and by JADCM with proposed

parameter set τ = kT,k = 61, 62, · · · , 80 (solid). The SSs

are sorted according to the performance indexes.

㪧㫆㫎㪼㫉㩷㪪㫇㪼㪺㫋㫉㪸㫃㩷㪛㪼㫅㫊㫀㫋㫐㩷㩿㪻㪙㪆㪟㫑㪀

㪇㪌㪇㪈㪇㪇㪈㪌㪇㪉㪇㪇㪉㪌㪇㪊㪇㪇

㪄㪍㪇

㪄㪋㪇

㪄㪉㪇

㪇

㪝㫉㪼㫈㫌㪼㫅㪺㫐㩷㩿㪟㫑㪀

㩿㪈㪀㩷㩷㪦㫉㫀㪾㫀㫅㪸㫃㩷㫉㫆㫎㩷㪻㪸㫋㪸

㩿㪊㪀㩷㩷㪡㫆㫀㫅㫋㩷㪛㫀㪸㪾㩷㩿㱠㪔㩷㪈㪃㩷㪉㪃㩷㵺㪃㩷㪉㪇㪀

㩿㪋㪀㩷㩷㪡㫆㫀㫅㫋㩷㪛㫀㪸㪾㩷㩿㱠㪔㩷㫂㪫㪃㩷㩷㫂㪔㪍㪈㪃㩷㪍㪉㪃㩷㵺㪃㩷㪏㪇㪀

㪇㪌㪇㪈㪇㪇㪈㪌㪇㪉㪇㪇㪉㪌㪇㪊㪇㪇

㩿㪉㪀㩷㩷㪧㪼㫉㫀㫆㪻㫀㪺㩷㪙㪪㪪

㪄㪍㪇

㪄㪋㪇

㪄㪉㪇

㪇

㪄㪍㪇

㪄㪋㪇

㪄㪉㪇

㪇

㪄㪍㪇

㪄㪋㪇

㪄㪉㪇

㪇

Figure 4: Power spectral density of MEG channel L42.

(1) original row data, (2) after noise reduction using Pe-

riodic BSS, (3) after noise reduction using JADCM with

conventional parameter set τ = 1, 2, · ·· , 20, and (4) after

noise reduction using JADCM with proposed parameter set

τ = kT, k = 61, 62, · · · , 80.

4 CONCLUSIONS

Several methods to separate periodic signals from

multichannel recordings have been proposed thus far.

However, they usually assume the problem to be ex-

actly determinable as a result of which their perfor-

mance in practical signal processings is not good.

In this report, we proposed a robust and prac-

tical to those method to extract periodic signals of

an ill-posed system. This method is based on the

jointdiagonalization of several time-delayed correla-

tion matrices. We showed the importance of select-

ing time delay parameters for the extraction of spe-

ciﬁc signals. We also veriﬁed the efﬁciency of the

proposed method in reducing power supply noise in

MEG recordings. We compared results with those of

Periodic BSS method and JADCM with the conven-

tional parameter set.

Because the power supply noise is artiﬁcial, it has

little ﬂuctuations in its period and auto-correlation

functions are constant even for long time delays. The

same feature is also found in periodic brain responses,

indicating that the method is effective.

REFERENCES

Barros, A. and Cichocki, A. (2001). Extraction of speciﬁc

signals with temporal structure. Neural Computation,

13:1995–2003.

Belouchrani, A., Abed-Meraim, K., Cardoso, J., and

Moulines, E. (1997). A blind source separation tech-

nique using second-order statistics. IEEE Trans. Sig-

nal Processing, 45:434–443.

Cardoso, J. and Souloumiac, A. (1996). Jacobi angles for

simultaneous diagonalization. SIAM J. Math. Anal.

Appl., 17:161–164.

Hyv¨arinen, A., Karhunen, J., and Oja, E. E. (2001). In-

dependent Component Analysis. John Wiley & Sons

Inc., New York.

Jafari, M., Wang, W., Chambers, J., Hoya, T., and Cichocki,

A. (2006). Sequential blind source separation based

exclusively on second-order statistics developed for a

class of periodic signals. IEEE Trans. Signal Process.,

54(3):1028–1040.

Molgedey, L. and Schuster, H. (1994). Separation of a mix-

ture of independent signals using time delayed corre-

lations. Phys. Rev. Let., 72:3634–3637.

Tang, A., Liu, J., and Sutherland, M. (2005). Recovery of

correlated neuronal sources from eeg: The good and

bad ways of using sobi. NeuroImage, 28:507–519.

Theis, F. and Inouye, Y. (2006). On the use of joint diago-

nalization in blind signal processing. In Proc. ISCAS,

pages 3586–3589.

BIOSIGNALS 2011 - International Conference on Bio-inspired Systems and Signal Processing

376

Tong, L., Liu, R., Soon, C., and Huang, Y. (1991). In-

determinacy and identiﬁability of blind identiﬁcation.

IEEE Trans. Circuit and Sys., 38:499–509.

Ziehe, A. and M¨uller, K. (1998). Tdsep-and efﬁcient

algorithm for blind separation using time structure.

In Proc. Int. Conf. on Artiﬁcial Neural Networks

(ICANN’98), pages 675–680.

Ziehe, A., M¨uller, K., Nolte, G., Mackert, B., and Curio,

G. (2000). Artifact reduction in magnetoneurogra-

phy based on time-delayed second-order correlations.

IEEE Trans. Biomed. Eng., 47:75–87.

A ROBUST AND PRACTICAL METHOD TO SEPARATE PERIODIC SIGNALS FROM MEG DATA USING

SECOND ORDER STATISTICS

377