非均匀采样信号的频谱分析

非均匀采样信号

非均匀采样信号经常出现在汽车工业、通信以及医学和天文学等领域。非均匀采样可能是由于传感器不完善、时钟不匹配或事件触发现象造成的。

频谱内容的计算和研究是信号分析的重要部分。传统的频谱分析方法，如周期图和韦尔奇方法，要求输入信号均匀采样。当采样非均匀时，可以对信号进行重采样或插值到均匀采样网格上。然而，这会给频谱增加不需要的伪影，并可能导致分析误差。

更好的替代方法是使用隆布-斯卡格尔方法，该方法直接处理非均匀采样，因此不需要重采样或插值。该算法已在 plomb 函数中实现。

缺失数据的信号的频谱分析

假设有一个温度监控系统，其中微控制器记录房间的温度，并每 15 分钟将此读数传输到存储度数的基于云的服务器。众所周知，Internet 连接中的小故障会阻止基于云的系统接收微控制器发送的一些读数。此外，在测量期间，微控制器的电池会至少耗尽一次，导致存在大采样间隔。

加载温度读数和对应的时间戳。

load('nonuniformdata.mat','roomtemp','t1')

figure
plot(t1/(60*60*24*7),roomtemp,'LineWidth',1.2)

grid
xlabel('Time (weeks)')
ylabel('Temperature (\circF)')

确定信号是否均匀采样的简便方法是获取连续采样时间间隔的直方图。

绘制以分钟为单位的采样间隔（时间差）直方图。仅包括存在采样的点。

tAtPoints = t1(~isnan(roomtemp))/60;
TimeIntervalDiff = diff(tAtPoints);

figure
hist(TimeIntervalDiff,0:100)
grid
xlabel('Sampling intervals (minutes)')
ylabel('Occurrences')
xlim([10 100])

与预期相符，大多数测量间隔为大约 15 分钟。然而，相当多的事件具有大约 30 到 45 分钟的采样间隔，这对应于缺失一个或两个连续采样。这会导致信号采样不均匀。此外，直方图中表示高频次的条形周围都呈现一些抖动。这可能与 TCP/IP 延迟有关。

使用隆布-斯卡格尔方法计算和可视化信号的频谱内容。为了帮助更好地可视化频谱，假设有最高 0.02 mHz 的频率，等效于每周约 13 个周期。

[Plomb,flomb] = plomb(roomtemp,t1,2e-5,'power');

figure
plot(flomb*60*60*24*7,Plomb)
grid
xlabel('Frequency (cycles/week)')
ylabel('Power (dBW)')

频谱显示主要周期性是每周 7 个周期和每周 1 个周期。这是可以理解的，因为数据来自以七天为日历周期的温控建筑物。显示每周 1 个周期的峰值的频谱线表明建筑物中的温度遵循以一周为单位的周期，周末温度较低，一周的工作日温度较高。每周包含 7 个周期的频谱线表明，还存在夜间温度较低、白天温度较高的以日为单位的周期。

非等间距采样信号的频谱分析

心率变异性 (HRV) 信号按时间表示心跳之间的生理变化，通常采用非均匀采样，因为人体心率不是恒定的。HRV 信号来自心电图 (ECG) 读数。

HRV 信号的采样点位于 ECG 的 R 峰值时间。每个点的振幅计算为连续 R 峰值之间时间差的倒数，并放在第二个 R 峰值的瞬时位置。

% Load the signal, the timestamps, and the sample rate
load('nonuniformdata.mat','ecgsig','t2','Fs')

% Find the ECG peaks
[pks,locs] = findpeaks(ecgsig,Fs, ...
    'MinPeakProminence',0.3,'MinPeakHeight',0.2);

% Determine the RR intervals
RLocsInterval = diff(locs);

% Derive the HRV signal
tHRV = locs(2:end);
HRV = 1./RLocsInterval;

% Plot the signals
figure
a1 = subplot(2,1,1); 
plot(t2,ecgsig,'b',locs,pks,'*r')
grid
a2 = subplot(2,1,2);
plot(tHRV,HRV)
grid
xlabel(a2,'Time(s)')
ylabel(a1,'ECG (mV)')
ylabel(a2,'HRV (Hz)')

R 峰值之间的变化间隔导致 HRV 数据中的采样时间非均匀性。考虑信号的峰值位置，并以秒为单位绘制其间隔直方图。

figure
hist(RLocsInterval)

grid
xlabel('Sampling interval (s)')
ylabel('RR distribution')

HRV 频谱中感兴趣的典型频带是：

这些频带大致限定了对 HRV 有贡献的独特生物调节机制的频率范围。这些频带中的任何波动都有生物学意义。

使用 plomb 计算 HRV 信号的频谱。

figure
plomb(HRV,tHRV,'Pd',[0.95, 0.5])

虚线表示 95% 和 50% 的检测概率。这些阈值测量峰值的统计意义。频谱显示上述所有三个感兴趣频带的峰值。然而，只有位于 VLF 内 23.2 mHz 处的峰值显示出 95% 的检测概率，其他峰值的检测概率小于 50%。低于 40 mHz 处的峰值被认为源于长期调节机制，如体温调节系统和激素因素。