Matlab / Python中的低通滤波器去除运动噪声

Question

我记录一个信号(皮肤电导)随着时间的推移，也就是说，我有一个时间序列。不幸的是，信号受到运动的影响。在用户指南中，我现在阅读了以下内容：

对数据采用低通滤波器去除高频噪声，高频噪声可归因于运动伪影和其他噪声成分。由于GSR响应的性质变化缓慢，可以在不影响感兴趣的数据的情况下使用1-5 Hz的截止频率。

如何在Matlab或Python中将这种截止频率的低通滤波器应用于我的时间序列？

Answer 1

在Matlab中，最简单的方法可能是使用butter命令(查找数字滤波器的系数)，然后是filter命令：https://www.mathworks.com/help/signal/ref/butter.html

butter(n, Wn)将接受您想要设计的过滤器的顺序(阶数越高，滤波器在行为上越理想)和归一化频率w。通常，为了正确设计滤波器，您需要知道原始数据采样的频率(Hz)。

由于您只需要消除基线低频噪声，请使用高通滤波器。

示例：

rawdata=rand(1000,1); %randomly generated data

fs=1000; %frequency at which your data was sampled, assuming 1000 Hz here

fc=5; %cut-off frequency (you mentioned 1-5 Hz in your question so I used 5)

Wn= fc/fs; % Normalized frequency

filt_order=8; % arbitrarily chose an 8th order filter

[b,a]=butter(filt_order,Wn, 'high'); 
% butter creates the coefficients of your digital filter
% the optional parameter 'high' specifies a high-pass filter

filtered_data=filter(b,a,rawdata);

Answer 2

设计一个过滤器可能是一个相当复杂的过程，没有人能给你一个蝙蝠的答案。过滤器设计是一种从经验中获得的技能，但是在获得这种经验之前，重要的是要对过滤器的基本知识有一个扎实的理解。您至少应该熟悉FIR和IIR滤波器(利弊)、什么是低/高/带通滤波器以及采样的基本知识。

在您的例子中，我假设您已经将时间序列存储在内存中的某个地方，因此您不需要在测量时过滤数据。在这种情况下，我会选择一个FIR滤波器，然后用Matlab的函数filtfilt对数据进行正向和后向滤波，这样你的滤波信号就不会有相位失真。这种类型的过滤只有在数据已经收集后处理时才有可能，因为这是非因果过滤，在这里您知道数据点在某个数据点之前和之后的值。然而，如果在记录时进行过滤，则不会知道未来的值。

通常通过迭代过程选择滤波器阶数和截止频率。要记住的一件事是，通过增加过滤器的顺序，您将提高过滤器的性能(更陡峭的滚动斜率)，这会带来更高的计算要求。

下面是一个简短的例子，信号由两个频率分别为3赫兹和25赫兹的正弦波组成。滤波的目的是设计一个低通滤波器，以消除25赫兹频率分量。我想展示两种不同过滤器的性能，一种是45级，另一种是70级。

% sampling frequency [Hz]
Fs = 1000;
% sampling period [s]
Ts = 1/Fs;

% time vector [s]
t = 0:Ts:3;

% Signal A frequency [Hz] and amplitude
amplitude_A = 5;
f_A = 3;
sig_A = amplitude_A.*sin(2*pi*f_A.*t);

% Signal A frequency [Hz] and amplitude
amplitude_B = 2;
f_B = 25;
sig_B = amplitude_B.*sin(2*pi*f_B.*t);

% sum of the two signals
X = sig_A + sig_B;
plot(t,X,'LineWidth',2);
hold on;
grid on;

% Low pass filter at ~15Hz
% normalized frequency wn
wn = 15/Fs;
b_1 = fir1(45,wn);
b_2 = fir1(70,wn);

% filtered signal
a = 1; % fir filter does not have poles (transfer function denominator = 1)
Y_1 = filtfilt(b_1,a,X);
plot(t,Y_1,'-r','LineWidth',2);
Y_2 = filtfilt(b_2,a,X);
plot(t,Y_2,'-k','LineWidth',2);
xlabel('Time (s)');
ylabel('Amplitude (V)');
legend('original signal','filtered order 45',...
    'filtered order 70');

这是输出：

​

​

在图中可以看到，低阶滤波器没有对高频分量进行完全滤波，但这两个滤波器都衰减了信号的幅度。通过调整截止频率和滤波器阶数，可以提高输出。这也是一个简单的过滤例子，如果这两个组件在频率上更接近，那么要获得一个好的结果就会更具挑战性。

对于你的情况，我假设皮肤电导(GSR)比运动引起的伪影有更低的频率(低于1Hz)，所以运动是一个更高的频率成分，你想摆脱它。我认为您应该仔细观察原始信号图，并尝试确定运动引起的变化，然后对信号进行滤波并检查该运动分量是否减少，继续更改滤波参数，直到您对结果满意为止。另外(如果需要的话)，你应该了解信号的频谱(执行FFT)，这可以帮助你选择截止频率。如果建议的截止频率在1-5赫兹的某个地方，我会从它开始，看看在这个范围内的频率是否能做到这一点。

在一篇文章中解释过滤几乎是不可能的，但是仅仅为了让你开始，了解一些基本知识是很好的。请记住，在过滤器的设计中，通常没有完美的解决方案，因为每个过滤器都以某种方式扭曲了原始数据。