对一组信号使用 pwelch：一些问题 (Matlab)

Question

我想在一组信号上使用 pwelch，但我有一些问题。

首先，假设我们有 32 个持续时间为 30 秒的 (EEG) 信号。采样频率为 fs=256 samples/sec，因此每个信号的长度为 7680。我想使用 pwelch 来估计这些信号的功率谱密度 (PSD)。

问题 1: 基于 pwelch 的 documentation、

pxx = pwelch(x) returns the power spectral density (PSD) estimate, pxx, of the input signal, x, found using Welch's overlapped segment averaging estimator. When x is a vector, it is treated as a single channel. When x is a matrix, the PSD is computed independently for each column and stored in the corresponding column of pxx.

但是，如果调用pwelch如下

% ch_signals: 7680x32; one channel signal per each column
[pxx,f] = pwelch(ch_signals);

结果 pxx 的大小为 1025x1，而不是我预期的 1025x32，因为文档指出如果 x 是矩阵，则 PSD 是针对每一列独立计算的并存储在pxx对应的列中。

问题二： 假设我克服了这个问题，并且独立计算了每个信号的 PSD（通过将 pwelch 应用于 ch_signals 的每一列），我想知道这样做的最佳方法是什么。假设信号是一个 30 秒的信号，采样频率为 fs=256，我应该如何调用 pwelch（使用什么参数？）以使 PSD 有意义？

问题 3： 如果我需要将我的 32 个信号中的每一个拆分为 windows 并将 pwech 应用到每个 windows，会怎样是最好的方法？假设我想将每个 30 秒信号拆分为 3 秒的 windows，重叠 2 秒。我应该如何为每个 windows 调用 pwelch？

Answer 1

这是一个例子，就像你的情况一样，

结果表明，该算法指示的信号频率恰到好处。

矩阵的每一列，y 都是正弦曲线以检查其工作原理。

windows是3秒重叠2秒，

Fs = 256;                       
T = 1/Fs;                    
t = (0:30*Fs-1)*T;           
y = sin(2 * pi * repmat(linspace(1,100,32)',1,length(t)).*repmat(t,32,1))';
for i = 1 : 32
[pxx(:,i), freq] = pwelch(y(:,i),3*Fs,2*Fs,[],Fs); %#ok
end
plot(freq,pxx);
xlabel('Frequency (Hz)'); 
ylabel('Spectral Density (Hz^{-1})');

对一组信号使用 pwelch：一些问题 (Matlab)

Using pwelch to a set of signals: some questions (Matlab)

matlab

spectral-density