【问题标题】:Extracting EEG Components from Signal within MATLAB在 MATLAB 中从信号中提取 EEG 分量
【发布时间】:2011-08-19 14:30:36
【问题描述】:

我在 MATLAB 中有一个简单的 EEG 信号,如下图所示。而我想要的是根据下表提取脑电图的成分。

  • Delta - 高达 4 Hz;
  • Theta - 4 -> 8 Hz
  • 阿尔法 - 8 -> 13 赫兹
  • 测试版 - 13 -> 30 赫兹
  • 伽玛 - 30 -> 100 赫兹

在第一次尝试解决这个问题时,我尝试使用 MATLAB 中的“fdatool”构建带通滤波器,以提取分量“theta”信号,但没有成功。

附上使用“fdatool”获得的过滤器的代码。

function Hd = filt_teta
%FILTROPARA TETA Returns a discrete-time filter object.

%
% M-File generated by MATLAB(R) 7.9 and the Signal Processing Toolbox 6.12.
%
% Generated on: 05-May-2011 16:41:40
%

% Butterworth Bandpass filter designed using FDESIGN.BANDPASS.

% All frequency values are in Hz.
Fs = 48000;  % Sampling Frequency

Fstop1 = 3;           % First Stopband Frequency
Fpass1 = 4;           % First Passband Frequency
Fpass2 = 7;           % Second Passband Frequency
Fstop2 = 8;           % Second Stopband Frequency
Astop1 = 80;          % First Stopband Attenuation (dB)
Apass  = 1;           % Passband Ripple (dB)
Astop2 = 80;          % Second Stopband Attenuation (dB)
match  = 'stopband';  % Band to match exactly

% Construct an FDESIGN object and call its BUTTER method.
h  = fdesign.bandpass(Fstop1, Fpass1, Fpass2, Fstop2, Astop1, Apass, ...
                      Astop2, Fs);
Hd = design(h, 'butter', 'MatchExactly', match);

有什么建议可以解决这个问题吗?

谢谢大家

【问题讨论】:

    标签: matlab signal-processing


    【解决方案1】:

    一种更简单的方法可能是简单地采用 FFT 并将除您可能感兴趣的特定范围之外的频率分量归零,然后逆 FFT 以返回到时域。

    Keep in mind that you'll have to zero out the positive frequency and negative frequency to maintain that the signal in the frequency domain is conjugate symmetric about the 0 frequency. 如果不这样做,在计算逆 FFT 时会得到一个复信号。

    编辑: 例如下面的代码,在时域中产生两个正弦曲线,一个对应的 DFT(用 FFT 计算),然后去除一个峰值。

    t = 0:0.01:0.999;
    x = sin(t*2*pi*4) + cos(t*2*pi*8);
    subplot(2,2,1);
    plot(x)
    title('time domain')
    subplot(2,2,2);
    xf = fft(x);
    plot(abs(xf))
    title('frequency domain');
    subplot(2,2,3);
    xf(9) = 0; xf(93) = 0;  % manual removal of the higher frequency
    plot(abs(xf));
    title('freq. domain (higher frequency removed)');
    subplot(2,2,4);
    plot(ifft(xf));
    title('Time domain (with one frequency removed)')
    

    有几点需要注意。 DFT 中的频域有几个不同的范围: DC 偏移(常数),即 0 频率;一个正频率范围,它是(对于长度为 N 的原始信号)从 1 到 N/2 的条目;负频率范围是从 N/2 到 N-1 的条目;请注意,这不是错字,最高频率(N/2 处的那个)是重复的,并且对于正频率和负频率来说是相同的值。 (有些人使用fftshift 来表示这是人类可能会绘制的,但这实际上只是为了看起来/理解。)

    至于要删除哪些频率,您必须自己弄清楚,但我可以给您一个提示。最高频率(在频率位置 N/2 处)将对应于您的系统可表示的最高频率,即 fs/2,其中 fs 是您的采样率。您可以相应地缩放以确定要否定的那些。

    如果你没有正确否定相应的负频率,你会在反 fft 时得到一个虚信号。

    最后一条评论。这种方法只有在您可以提前获得所有信号的情况下才有效,因为您需要对整个信号使用 DFT。如果您想实时执行此操作,则需要像以前一样创建某种过滤器。

    【讨论】:

    • 嗨,Chris A.,感谢您的回答!您能否在 sn-p MATLAB 代码中给我一些指导,我如何“将特定范围以外的频率分量归零”?我对 FFT 和信号处理不是很熟悉:S
    【解决方案2】:

    如果过滤器的长度没有任何限制,请为过滤器选择更锐利的边缘。如果我在你的鞋子里,我会构建不同的滤波器(低通和高通)并在傅里叶变换中处理结果,以查看任何高频或低频与频率范围混合。 1)构建低通,提取增量 2) 为 theta alpha beta 构建带通 3)构建高通滤波器,提取伽马。

    【讨论】:

    • 嗨赫菲斯托斯,感谢您的回答!锐化滤镜的边缘是什么意思?也许降低阻带频率的值?
    • 你使用 Fstop1 = 3; ,但您的信号范围为 0 到 4 kHz。因此,您会在边缘丢失一些信号。例如,我会选择 3.98。它使过滤器更大。但是,精度会更好,并且差异可能很大。
    猜你喜欢
    • 2020-02-28
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2020-08-29
    • 1970-01-01
    • 2018-08-01
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多