EEG 信号的带宽答案

【问题标题】：Bandwidth of an EEG signalEEG 信号的带宽
【发布时间】：2014-01-10 18:47:36
【问题描述】：

我正在尝试在 Python 中对 EEG 信号执行 FFT，然后根据带宽确定它是 alpha 还是 beta 信号。它看起来不错，但结果图与它们应该的完全不同，频率和幅度值不是我所期望的。任何帮助表示赞赏，这是代码：

from scipy.io import loadmat
import scipy
import numpy as np
from pylab import *
import matplotlib.pyplot as plt

eeg = loadmat("eeg_2013.mat");
eeg1=eeg['eeg1'][0]
eeg2=eeg['eeg2'][0]
fs = eeg['fs'][0][0]
fft1 = scipy.fft(eeg1)
f = np.linspace (fs,len(eeg1), len(eeg1), endpoint=False)
plt.figure(1)
plt.subplot(211)
plt.plot (f, abs (fft1))
plt.title ('Magnitude spectrum of the signal')
plt.xlabel ('Frequency (Hz)')
show()
plt.subplot(212)
fft2 = scipy.fft(eeg2)
f = np.linspace (fs,len(eeg2), len(eeg2), endpoint=False)
plt.plot (f, abs (fft2))
plt.title ('Magnitude spectrum of the signal')
plt.xlabel ('Frequency (Hz)')
show()

还有情节：

【问题讨论】：

有一个潜在的问题，取决于输入数据，或者它们是如何归一化的，因为 EEG 信号通常有 10 到 50 Hz 的范围，你的 1 到 9 kHz，这把输入信号？你把它带到哪里？
如果您提供数据链接可能会有所帮助。
关于数据图，单位会很有帮助。
您没有指定与您的采样频率匹配的单位，但是，猜测单位是赫兹（即 200 Hz），那么根据 Nyquest 定理，您无法解析任何高于100 赫兹。
似乎您的 FFT 发生了偏移（0 在您的图中约为 4500Hz）。使用fftfreq 获取真实频率。

标签： python fft canopy

【解决方案1】：

为了得到 fft 频率的数组，你应该使用fftfreq；它为您提供了一组频率用作横坐标：

from scipy.fftpack import fftfreq

eeg = loadmat("eeg_2013.mat");
eeg1=eeg['eeg1'][0]
eeg2=eeg['eeg2'][0]
fs = eeg['fs'][0][0]
fft1 = scipy.fft(eeg1)
f=fftfreq(eeg1.size,1/fs)

抱歉，我无法在真实条件下测试此代码，因为您没有发布数据样本，但我希望这应该可以工作。

关于如何确定带宽，据我了解，您想获得基频。有不同的方法，无论您的信号是否嘈杂，或多或少复杂，......在您的情况下，您只想知道基频 f0 是否在 8-13Hz（alpha）或 13-30Hz（beta );一种非常简单的方法是计算 fft 在 8-13Hz 范围内的最大值：fft1[(f>8) & (f<13)].max()，如果大于 1000，则为 alpha 波，否则为 beta。如果您的信号不太相似，请发布一些不同类型样本的示例以及您将获得的结果，以便我们尝试更复杂的算法。

【讨论】：

【解决方案2】：

如果您的采样频率为fs 并且您有N=len(eeg1) 样本，那么fft 过程当然会返回一个N 值数组。其中第一个N/2对应频率范围0..fs/2，后半部分频率对应镜像频率范围-fs/2..0。对于真实输入信号，镜像一半只是正半部分的复共轭，因此在进一步分析中可以忽略（但在反 fft 中不能）。

所以本质上，你应该格式化

f=linspace(0,N-1,N)*fs/N

编辑：甚至更简单，只需对初始代码进行最少的更改

f = np.linspace (0,fs,len(eeg1), endpoint=False)

所以f 的范围从0 到fs 之前，忽略输出中fft 结果的后半部分：

plt.plot( f(0:N/2), abs( fft1(0:N/2) ) )

补充：可以用fftshift交换两半，那么正确的频率范围是

f = np.linspace (-fs/2,fs/2,len(eeg1), endpoint=False)

【讨论】：